Original link:tecdat.cn/?p=23050

Original source:Tuo End number according to the tribe public number

In this article, we will use R language to fit the linear mixing effects model of the data, and then visualize your results.

Linear mixed effects models are used when there are random effects, which occur when a random sample of units is measured multiple times. Measurements from the same natural group are not themselves independent random samples. Thus, these units or groups are assumed to be randomly drawn from the “population” of a group. Example cases include

  • When you divide and test each part individually (random unit).
  • When your sampling design is nested, such as within a cross-sectional quartile; A cross section within a woodland; Woodlands within districts (cross sections, woodlands and districts are random groups).
  • When you measure related individuals (families are random units).
  • When you measure the subject repeatedly (subject is random).

Linear models of mixed effects are implemented in the R command lME4 and lmerTest packages. Another option is to use the LME methods in the NMLE package. The method used in LME4 to calculate approximate degrees of freedom is more accurate than the method in the NMLE package, especially when the sample size is small.

Measure plaque length

This first data set was extracted from a paper by Griffith and Sheldon (2001, Ethology 61:987-993), in which white forehead patches of 30 male lead flycatchers on gotland Island, Sweden, were measured over two years. The plaque is important in attracting mates, but its size varies from year to year. Our goal here is to estimate plaque length in millimeters.

Read and examine data

  1. Read data from a file.
  2. Look at the first few lines of the data to see if it was read correctly.
  3. Create a pair of measurements showing each bird during the two-year study. You can try to make a lattice diagram. Is there evidence of measured variability between years?

Construct linear mixed effect model

  1. The data were modeled with linear mixing effects, and individual birds were treated as random groups. Note: Two measurements for each bird were taken in consecutive years of the study. For simplicity, the year is not included in the model. Convert it to a character or factor in R so that it is not treated as a numeric variable. Recalculate repeatability with this model as described in steps (2) and (3) below. How does the repetitive interpretation change?
  2. Extract parameter estimates (coefficients) from the saved LMER object. Check the output of random effects. What are the two sources of random variation? What do fixed effects mean?
  3. In the output, check the standard deviation of the random effect. There should be two standard deviations: “intercept” and “residual”. This is because the mixed-effect model has two sources of random variation: differences in repeated measurements within birds, and real differences in frontal patch length between birds. Which of these two sources corresponds to “(intercept)” and which corresponds to “residual”?
  4. Also check the output of the fixed effect results. The only fixed effect in the model formula is the average of all length measurements. This is called “intercept”, but is not to be confused with the intercept of random effects. The fixed effect output gives you an estimate of the mean value and the standard error of that estimate. Note how the fixed effects output provides an estimate of the mean, while the random effects output provides an estimate of the variance (or standard deviation).
  5. The variance component was extracted from the fitting model to estimate the repeatability of plaque length in each year *.
  6. Explain the repeatable measurements obtained in the previous step. If you get a repeatability of less than 1.0, what is the source of variation between measurements within individuals. Is it just measurement error?
  7. Generates a graph of residuals and fitted values. Notice anything wrong with that? There seems to be a slight positive trend. This is not an error, but the result of a “contraction” of the best linear unbiased predictor (BLUPs).

Analysis steps

Read and examine the data.


head(fly)
Copy the code

# Chart (patch ~ bird)Copy the code

# But a better way to display paired data is to display the plot(res=patch, x = year) with paired interaction plotsCopy the code

Plot (y = patch, x = factor(year), theme_classic)Copy the code

Fit a linear mixed effect model. The output of summary() shows two sources of random variation: variation between individual birds (bird intercept), and variation between repeated measurements of the same birds (residuals). Each source has an estimated variance and standard deviation. The fixed effect is just the average of all birds — another intercept.

# 2. Parameter Estimation Summary (z)Copy the code

# 5. Variance component VarCorr(z)Copy the code

1.11504^2/(1.11504^2 + 0.59833^2)Copy the code
# # 0.7764342 [1]Copy the code
# 7. Residuals and fitting values plot(z)Copy the code

Goldfish visual

Cronly – Dillon and Muntz (1965; J. Exp. Biol 42: 481-493) to measure color vision in goldfish by visual motor response. Here, we will fit the data, including the full wavelength of the test. Each of the five fish was tested at all wavelengths in a random order. A high sensitivity value indicates that the fish can detect low light intensity. An important feature of the visuomotor response is that fish are not used to it, and a measurement of visual sensitivity at one wavelength is unlikely to have an impact on a later measurement at another wavelength.

Read and examine data

  1. Read the data in the file and look at the first few lines to make sure it was read correctly.
  2. Interaction diagrams were used to compare the reactions of individual fish at different wavelengths of light.
  3. What type of experimental design is used? * This will determine the linear mix model to use when fitting the data.

Construct linear mixed effect model

  1. Fit a linear mixed effects model to the data. This can be done with lmer(). Find “malformed fit”, “Boundary (singular) fit: see? IsSingular”
  2. Plot the fitted (predicted) values **. The difference between the predicted and observed values for each fish represents the residual.
  3. What assumptions did you make in (1)? Create a graph of residuals with fit values to check one of these assumptions.
  4. Extract parameter estimates from the saved LMER object. Check the results of fixed effects. The coefficients given are the same as the interpretation of categorical variables using LM analysis.
  5. Check the output of random effects. Again, two sources of random error appear in our mixed effects model. What are they? Which of these corresponds to “(intercept)” in the output, and which corresponds to “residual”? Note that the estimated standard deviation of one of the sources of change in this data set is very small. This is the reason behind the malformation fit information. It is unlikely that the variance between fish is really zero, but this data set is so small that low variance estimates can occur due to sampling errors.
  6. Generate model-based estimates of average sensitivity for each wavelength.
  7. Are there significant differences between wavelengths? Generate anOVA tables for LMER objects. What effect is being tested here, random or fixed? Interpret anOVA results.

* This is a “by subject” repeated measurement design, as each fish is measured once under each experiment. It is essentially the same as a random complete block design (treating each fish as a “block”).

* Visualization is preferred because the data and fit values are plotted. Notice how similar the predicted values are from fish to fish. This suggests that the estimated differences between individual fish in this study were very small.

*** In general, only fixed effects are tested in anOVA tables. It is possible to use invalid assumptions with no variance in test random effects.

Analysis steps

Read and examine the data.

x <- read.csv("fish.csv", 
        stringsAsFactors = FALSE)
head(x)
Copy the code

Fit a linear mixed effect model.

The model assumes that the residuals of all fitted values are normally distributed and the variances are equal. The method also assumes that the random intercepts between individual fish are normally distributed. The method also assumes random sampling of groups (fish), with no effect on the measurements of the same fish.

# # 1. Fit the mixed effect model.Copy the code
## boundary (singular) fit: see ? isSingularCopy the code

# 2. This plots the fit values for each fish individually. vis(z)Copy the code

 

# 3. Test hypothesis plot(z)Copy the code

# 4. Extract parameter estimates Summary (z)Copy the code

# 6. Model-based mean sensitivity Estimation Means (Z)Copy the code

# 7. ANOVACopy the code

Concentration of phenolic substances in Yarrow

The project experimentally investigated the effects of fertilization and grazing in the boreal forest ecosystem of a national park (Krebs, C.J., Boutin, s. & Boonstra, R., Eds (2001a) Ecosystem Dynamics of the Boreal Forest.Kluane Project. Oxford University Press, New York), and present data from a study of the effects of plant resources and herbivores on the defensive chemistry of underlying plant species.

Each of 16 5×5 meter plots was randomly assigned to one of four experiments. 1) Fencing to exclude herbivores; 2) Fertilize with n-P-K fertilizer; 3) Fence and fertilize; 4) Untested controls. Then, each of the 16 plots was divided in two. One side of each plot (randomly selected) was continuously tested during the 20-year study. The other half of each plot was tested for the first ten years and then returned to its untested state. The data analysed here recorded the concentration of phenols (a rough measure of plant defense compounds) in Achillea millefolium, a common herb in the plot. It is measured in milligrams of tannic acid equivalent per gram of dry weight.

Visual data

  1. Read data from a file.
  2. Check the first few lines of data. The experiment was presented as a single variable with four levels (rather than as two variables, fence and fertilizer, in a model designed with 2×2 factors). Duration indicates whether half of the field was tested for the full 20 years, or whether the experiment stopped after 10 years. The variable “CH” is the concentration of phenols in Yarrow.
  3. Draw a graph to illustrate the concentration of phenols in Yarrow for different experiments and duration categories. There are not many data points in the combination of each experiment and duration level, so it may be better to plot a bar chart by group rather than a box chart by group.
  4. Add line segments to connect pairs of points.

Fit a linear mixed effect model

  1. What type of experimental design was used? * This will determine the fitting of the linear mixed model of the data.
  2. In the absence of interaction between experiment and duration, linear mixed model fitting was performed on the data. Logarithms of phenols were used as the dependent variable because logarithmic conversion improves the fitting of the data to the assumptions of the linear model.
  3. Visual model fitting of data. Separate panels by duration (if xvar is experiment) or experiment (if xvar is duration). Visreg () does not preserve pairings, but allows you to check for residuals.
  4. The model fit is now repeated, but this time involving the interaction between experiment and duration. Visualize the fit between model and data. What is the most obvious difference between the two model fits, one with interaction and the other without? Describe what is “allowed” in a model that includes interaction items, but not in a model that does not. Decide, which model best fits the data?
  5. Use diagnostic diagrams to examine a key assumption of the linear hybrid model of the model including interaction terms.
  6. Use the fit model object to estimate the parameters (including interactions) of the linear model. Notice that there are now many coefficients in the fixed effects table.
  7. In the output from the previous step, you will see two quantities of “std.dev” under the “random effects” TAB. Explain what these quantities refer to.
  8. To estimate the model fit mean for all fixed effect combinations.
  9. Generate anOVA tables for fixed effects. Which items are statistically significant?
  10. By default, lmerTest tests model items using the sum of squares of Type 3, rather than in order (Type 1). Repeat the ANOVA table with type 1. Did it make a difference? 支那

* Experiments used a chunking design, in which the entire block was randomly assigned to different experiments, and then different levels of the second experiment (duration) were assigned to half of the block.

* There should be no difference as the design is completely balanced.

Analysis steps

Read and review the data.

A good strategy is to sort the categories of experiments and put the control group first. This will make the output of the linear model more useful.

# 2. Check head(x)Copy the code

First of all, Factor (treat,levels=c("cont","exc","fer","bo") plot(data = x, y = log(phe), x = treat, fill = dura, color = dura)Copy the code

# 4. Plot pairs of data on multiple panels (data = x,y = log(ach, x = dur, fill = dur, col = dur)Copy the code

Fit a linear mixed effect model. Fixed effects are “experiments” and “duration”, while “chunks” are random effects. When fitting interactions, the size of the difference between experimental levels will vary between duration levels.

Since random effects are also present (blocks), the coefficient table will show variance estimates for two sources of random variation. One is the variance of the residual of the fitting model. The second is the variance between (random) block intercepts.

# 2. Fit mixed effects model - no interactionCopy the code
# 3. Visualization VIS (Z)Copy the code

# 4. Including interaction and visual again z.i nt < - lmer (log (phen. Ach) ~ * duration treatment + 1 | (plot), data = x) vis (z.i nt, overlay = TRUE)Copy the code

# 5. Graph to test homogeneity of variance (and normality) plot(z)Copy the code

# 6. Coefficient Summary (z)Copy the code

# 8. Mean value of model fitting means(z, data = x)Copy the code

The default of anOVA (Z) # lmerTest is the sum of squares of 3 classes.Copy the code

Anova (z, type = 1)Copy the code

Most welcome insight

1. Lmer mixed linear regression model based on R language

2. Exploring LME4 Generalized Linear Mixing Model (GLMM) and Linear Mixing Model (LMM) with Rshiny

3. Practical case of linear mixed effect model of R language

4. Practical case of linear mixed effect model of R language 2

5. Practical case of linear mixed effect model of R language

6. Partially folded Gibbs sampling of Linear mixed-effects Models

7. A mixed-effect model of R language LME4 to study teacher’s popularity

8. Har-rv model based on Mixed data sampling (MIDAS) regression in R language predicts GDP growth

9. Hierarchical linear model HLM using SAS, Stata, HLM, R, SPSS and Mplus