Full text download link: http://tecdat.cn/?p=23681
What is the difference between the linear mixed effects model and the linear model we already know(Click “Read the original text” at the end of the article to get the completeCode data)?
Related videos
A linear mixed model (sometimes called a “multilevel model” or a “hierarchical model”, depending on the context) is a regression model that takes into account both (1) the independent variable of interest (such as lm( )) the variation explained – fixed effects, and (2) the variation not explained by the independent variable of interest – random effects. Since the model includes a mixture of fixed and random effects, it is called a mixed model. These random effects essentially give the error term structure.
The definitions of fixed and random effects may vary, so be aware of your interpretation in the literature; however, for most purposes, if data are collected from all levels of interest, you can treat a variable as Fixed effects factors (e.g. gender: male/female, condition: easy/medium/hard, dose: low/high) if the variable has a bunch of possible levels but you only have a random set (e.g. subjects, stimuli Although these samples will have some specificity, you generally don’t care about them, and the goal is to make generalizations about the broader population (such as all people, all settings, all classrooms).
Example
Let’s say you’re interested in language, or more specifically, in the relationship between pitch and politeness. You ask your subjects to respond to hypothetical scenarios (IV, within-subjects) that are either formal situations that require a polite attitude (e.g., giving a professor an excuse for being late) or something more informal occasions (e.g., explaining to a friend why you are late) and measure their pitch (DV). Each subject was given a list of all scenarios, so each subject gave multiple polite or informal responses. You also note the gender of each subject (IV, between subjects), as this is another important effect on vocal tone.
Of the linear models we have seen so far, we will build such a model.
Tone = Politeness + Gender + ?
The last term here is our error term. This error term represents the deviation from our predictions due to “random” factors that we cannot control in the experiment.
For this kind of data, since each subject gave multiple responses (“repeated measures” design), we can see that this would violate the important independence assumption in linear modeling: same subject Multiple reactions cannot be considered independent of each other. In our scenario, each person’s tone of voice is slightly different, which will be a specific factor affecting all responses of the same subject, making these different responses interdependent (correlated) rather than independent.
Random effects
The way we’re going to handle this is to add a random effect to the subject. This allows us to account for this non-independence by assuming a different “baseline” pitch value for each subject. Thus, Subject 1 might have had an average pitch of 233 Hz across different utterances, while Subject 2 might have had an average pitch of 210 Hz. In our model, we account for individual differences in these tones via random effects across subjects.
Let’s analyze some data as an example.
table(subject)
Visualize the data.
qplot(condition, pitch, facets = . ~ subject)
Subject “F#” is a female subject. Object “M#” is a male object. You’ll notice right away that men’s voices are lower than women’s (this is to be expected). But beyond that, within groups of men and women, you see a lot of individual differences, with some people having relatively high gender values and others having relatively low gender values.
Correlation of samples from the same subject
Another way to put it is that within subjects, there are correlations between pitches under different conditions. Let’s visualize it.
Modeling individual means with random intercepts
We can model these individual differences by assuming a different random intercept for each participant; each participant is assigned a different intercept value (i.e., a different mean tone), and the mixed model is essentially You estimate these intercepts.
Looking back at our model, our previous formula was.
Tone = intercept + politeness + gender + ?
Our updated formula looks like this.
Tone = intercept + politeness + gender + (1|individual) + ?
“(1|subject) ” is the R syntax for random intercepts. What this sentence means is “assuming the intercept is different for each subject”…and “1” represents the intercept here. You can think of this formula as telling your model that it should expect multiple responses per subject, and that these responses will depend on each subject’s baseline level. This effectively solves the problem of non-independence caused by multiple responses from the same subject.
Note that this formula still contains a general error term ?. This is necessary because even if we account for variations in each subject, there will still be “random” differences between different pitches of the same subject.
Click on the title to view previous issues
R language LME4 mixed effects model to study the popularity of teachers
Swipe left or right to see more
01
02
03
04
Have an idea of the average across different participants in different conditions.
aggregate(pitch ~ subject, FUN = "mean")
Now using lmer(), we can estimate the mean for each participant. To do this, we will include a random intercept for each subject and then look at the estimated intercept.
coef(lmer(pitch ~ (1 | subject))
#Fixed effect + main body of random effect \['(intercept)'\] + subject
Note that the estimates are quite close to the actual mean pitch. We can see that the actual mean pitch across subjects is the estimate (Intercept), while the standard deviation of the mean pitch across subjects is the criterion for the random effects. Poor(Std.Dev).
# Use original data mean
## \[1\] 193
sd
## \[1\] 63.47
# Use estimated intercepts for each subitem mean(subject\[1\]\[,'(Intercept)'\])
## \[1\] 193
sd
## \[1\] 62.4
# This is also the summary in the model output summary(res1)
Include fixed effects
Since we predict that the condition of the hypothetical state (“informal” vs. “polite attitude”) will affect pitch (perhaps the pitch will be higher in the informal state), in addition to the gender of the subject (female pitch may be higher), let us incorporate these conditions into the model while also accounting for a random intercept for each subject (let the intercept vary across subjects).
lmer(tone~polite + gender + (1|individual))
Pitch j=intercept + intercept j + state + gender
Pitch A=intercept + change in A-intercept [_Intercept_ _Shift_] + state + gender
ggplot(x=condition, y=mean)
Here we will also compare status and gender so we can see the effect of condition on the “average” between women and men, and gender on “informal” and “impolite” average effect.
It can be seen that our average pitch is 192.88, the pitch of the informal state is higher than that of the polite state, b=9.71, t=3.03, and the pitch of women is higher than that of men, b=54.10, t=5.14. We can use a rule of thumb that if t is greater than 2, it is probably significant.
More model information
We can use many different types of information to evaluate models. This information can be useful when comparing models. A useful measure is the AIC, which is Bias + 2?(p + 1), where p is the number of parameters in the model (here, we break the parameters down, so 1 is the estimate The residuals, p are all other parameters, e.g. fixed effects coefficients + estimated variances of random effects, etc.). A lower AIC is better because higher bias means the model doesn’t fit the data well. Since AIC increases as p increases, AIC is penalized for more parameters.
deviance = -2*logLikelihood; deviance
## \[1\] 789.5
# Calculate AIC by hand p = 4 # Parameter = 3 (fixed effect) + 1 (random effect deviance + 2*(p + 1) # Total parameters = 4 + 1 used to estimate residuals
## \[1\] 799.5
AIC(res2)
## \[1\] 799.5
Extract all coefficients
ggplot(aes(x=factor(subject), y=mean_pitch))
coef(res2)
Here we can see that this model produces a separate intercept for each subject, in addition to a parameter estimate/slope for condition and sex that is constant across subjects. From here, we can try to estimate the average pitch of a given sample based on these coefficients. For example, let’s try to estimate the mean of subjects’ F1 (xˉ=232.0357) using their estimated intercept and the effect of being female.
179.3003 + 0*(9.7) + 1*(54.10244)
## \[1\] 233.4
pretty close
Now it is M3 (xˉ=168.9786):
220.3196 + 0*(9.7) + -1*(54.10244)
## \[1\] 166.2
not bad
Random slope
In the above model, we assume that the impact of polite attitude is the same for all subjects, therefore, the coefficient of polite attitude is one. However, the impact of polite attitudes may be different for different subjects; that is, there may be an interaction of polite attitudes among subjects. For example, we might expect that some people will be more polite and others less polite in situations where etiquette is required. What we need, therefore, is a random slope model in which not only subjects are allowed to have different intercepts, but also their effects on politeness are allowed to have different slopes (i.e., different effects of state on pitch).
Let’s start visualizing the data.
Based on this plot, does it look like the slope is unstable across subjects?
Now add random slopes.
coef(res3)
anova(res2, res3, refit=FALSE)
Increasing the random slope for each subject took up another 2 degrees of freedom and did not significantly improve the model fit, χ2(2)=0.024, P=0.99. Note that df=2 because we include both the slope variance and the correlation between the intercept and slope. Looking at the AIC values, more complex models have higher AIC values, so we want to use less complex (more concise) models. If we look at the estimated slopes, we can see that they are very similar across different samples. Therefore, we do not need to include conditional random slopes in the model.
However, others would argue that we should keep our models as comprehensive as possible and should include random slopes even if they do not improve the model!
Testing significance
While there is some debate over whether you should get p-values for lmer() models (e.g., this one; most of the debate revolves around how to calculate dfs), you can get an approximation of df (and therefore p-values) using the {lmerTest} package.
Get P value
summary(res3b)
Compare model output to SS/Kenward-Roger appox
anova
anova(res2b)
Model comparison
On the other hand, some argue that model comparison using likelihood ratio tests is a better way to test whether a parameter is significant. That is, if adding the parameter to your model significantly improves the fit of the model, then the parameter should be included in the model.
The likelihood ratio test essentially tells us how much more likely the data is under a more complex model than under a simpler model (these models need to be nested!):
The distribution of D is approximately χ2, with degrees of freedom df2-df1. We can either “manually” do this calculation, or just use the anova() function directly!
anova(res4, res4b)
# Manual calculation dev0 <- -2*logLik(res4) # Simpler deviation model devdiff <- (dev0-dev1) # Deviation difference dfdiff # Difference of parameters (using dfs).
Here we compare two nested models, one without conditions and the other with conditions. Through model comparison, we concluded that adding the condition to our model is necessary because it significantly improves the fit of the model, χ2(1)=8.79, P<0.01.
REML and ML
Let us start with a statistical model specifying (i) fixed effects and (ii) a variety of random effects for the normally distributed variation and covariance.
-
In ML (maximum likelihood) estimation, we calculate the logarithm (likelihood) (LL) of the data for arbitrarily chosen parameter values in groups (i) and (ii) above. Then, we look for parameter values that maximize L (or minimize -L). These optimal parameter values are called ML parameter estimates. The maximum value of L (-2 times) is called the bias of the model. For some purposes, such as describing data, we focus on ML parameter estimates; for other purposes, such as model comparison, we focus on bias. When comparing models with different fixed effects, you should use ML, and you must include lmer(, REML=FALSE). Additionally, if you are comparing an lm() and lmer() model (i.e. testing whether it is necessary to use any random effects), you should also use ML estimation.
-
In REML (REstricted ML) estimation, our main interest is in estimating random effects, not fixed effects. Imagine that we restrict our parameter space by setting the fixed effects parameters in group (i) above to some reasonable values. In this restricted space, we look for random effect parameter values in set (ii) that maximize the LL value of the data; while paying attention to the maximum LL value. Then repeat this process multiple times. Then the fixed effect parameter values, the estimated values of the random effect parameter and the maximum value of LL are averaged. This averaging method can obtain REML parameter estimates and REML deviation values. Because this procedure pays little attention to the fixed effects parameters, it should not be used to compare models with different fixed effects structures. You should use this method when comparing models with different random effects. To do this, you should get into the habit of including lmer(, REML=TRUE) and using anova(, refit=FALSE) when running model comparisons.
Example
Test whether random effects are necessary
dev1 <- -2*logLik(res5) dev0 <- -2*logLik(res5b) devdiff <- as.numeric(dev0-dev1); devdiff
# Compare AICs AIC(res5b, res5)
AICtab(res5b, res5) # AIC
Use anova to directly compare lm and nlme models
It appears that, yes, including a random intercept within subjects makes sense, χ2(1) = 19.51, P < 0.001.
Here, the χ2 distribution is not always a good approximation of the null distribution (too conservative when testing some random effects, and not conservative enough when testing some fixed effects).
Scene effect
Different scenes may elicit different “pitch” values; therefore, the tone of a particular scene may be relevant between different subjects, or even within a subject, in polite and informal situations. of. We can think of this as a random effects model.
ggplot(d, aes(x=scenario, y=pitch)
anova
coef(res4)
Similar to the random intercepts across subjects, now we have the average pitch level for each scenario.
So, what about the slope change for each scene?
ggplot(byscenario, aes(x=condition, y=pitch))
anova
There is no improvement in the model. This suggests that our schemes may be similar in terms of differences between informal and polite conditions.
Hybrid model description
There are several important things to say here. You may ask “Which random slopes should I specify?” or even “Are random slopes necessary at all?”
Conceptually, it makes perfect sense to include random slopes and random intercepts together. After all, you could argue that people respond differently to experimental manipulations! Likewise, you could argue that people respond differently to experimental manipulations. Likewise, you can always assume that the effect of an experimental manipulation will be different for all items in the experiment.
In the above model, the whole point of our study is to say something about politeness. We are not interested in sex differences, but they are well worth controlling for. This is why we have random slopes (by subject and item) for the effect of politeness, but not gender. In other words, we only modeled by-subject and by-item changes in the effect of polite attitudes on pitch.
Everything discussed in the context of linear models applies directly to mixed models. Therefore, you also have to worry about collinearity and outliers. You also have to worry about homoskedasticity (equal variances) and potential lack of normality.
Independence, as the most important assumption, requires a special word. One of the main reasons we moved to mixed models instead of just linear models was to account for non-dependencies in our data. However, mixture models can still violate independence. If you are missing important fixed or random effects. So, for example, if we analyze our data with a model that does not include the random effect “subject”, then our model does not “know” that each subject has multiple responses. This amounts to a violation of the independence assumption. Therefore, choose your fixed and random effects carefully to address non-independence issues.
Some other notes.
If your dependent variable is…
-
Continuous: Linear regression model using mixed effects
-
Binary: Logistic regression model using mixed effects
The function lmer is used to fit linear mixture models, and the function glmer is used to fit generalized (non-Gaussian) linear mixture models.
Click “Read original text” at the end of the article
Get full text complete information.
This article is selected from “R Language Uses Linear Mixed Effects (Multi-level/Level/Nested) Model to Analyze the Relationship between Tone Level and Polite Attitude”.
The data analyzed in this article is shared to the membership group. Scan the QR code below to join the group!
Click on the title to view previous issues
R language to build and visualize mixed effect model mixed effect model
R language linear mixed effects model practical case
R language uses latent class mixed effects model (Latent Class Mixed Model, LCMM) to analyze Alzheimer’s age data
R language Bayesian generalized linear mixed (multi-level/level/nested) model GLMM, logistic regression analysis of data on influencing factors of education grade repetition
R language for estimating multivariate labeled latent process mixed effects models (lcmm) to analyze cognitive processes in psychological testing
R language factorial experimental design nlme fits nonlinear mixed model to analyze nitrogen application level in organic agriculture
R language nonlinear mixed effects NLME model (fixed effect & random effect) study on the kinetics of the anti-asthmatic drug theophylline
R language uses a linear mixed effects (multi-level/level/nested) model to analyze the relationship between pitch and politeness.
R language LME4 mixed effects model to study the popularity of teachers
R language nlme, nlmer, lme4 uses (non-) linear mixed model non-linear mixed model to analyze algae data examples
R language mixed linear model, multi-level model, regression model analysis of student average grade GPA and visualization
R language linear mixed effects model (fixed effect & random effect) and interactive visualization 3 cases
R language uses lme4 multi-level (mixed effects) generalized linear model (GLM) and logistic regression to analyze education grade retention survey data
R language linear mixed effects model practical case
R language mixed effects logistic regression model to analyze lung cancer data
How to analyze depressive symptoms using latent class mixed effects model (LCMM) in R language
Research on diagnostic accuracy of Bayesian hierarchical mixture model based on copula in R language
R language to build and visualize mixed effect model mixed effect model
R language LME4 mixed effects model to study the popularity of teachers
R language linear mixed effects model practical case
R language uses Rshiny to explore lme4 generalized linear mixed model (GLMM) and linear mixed model (LMM)
Research on diagnostic accuracy of Bayesian hierarchical mixture model based on copula in R language
How R language solves the problem of Singular fit in linear mixed models
lmer mixed linear regression model based on R language
R language uses WinBUGS software to build a hierarchical (hierarchical) Bayesian model for academic ability tests
R language hierarchical linear model case
R language uses WinBUGS software to build a hierarchical model for the Academic Aptitude Test (SAT)
Hierarchical Linear Model HLM using SAS, Stata, HLM, R, SPSS and Mplus
R language uses WinBUGS software to build a hierarchical (hierarchical) Bayesian model for academic ability tests
Multilevel (hierarchical) linear models in SPSS Multilevel linear models to study plastic surgery data
Estimating HLM multi-level (hierarchical) linear model models with SPSS