tips:rma_vs_lm_lme_lmer
Differences
This shows you the differences between two versions of the page.
Both sides previous revisionPrevious revisionNext revision | Previous revision | ||
tips:rma_vs_lm_lme_lmer [2022/04/22 12:00] – Wolfgang Viechtbauer | tips:rma_vs_lm_lme_lmer [2023/11/14 08:00] (current) – Wolfgang Viechtbauer | ||
---|---|---|---|
Line 1: | Line 1: | ||
===== A Comparison of the rma() and the lm(), lme(), and lmer() Functions ===== | ===== A Comparison of the rma() and the lm(), lme(), and lmer() Functions ===== | ||
- | In essence, the commonly | + | Commonly |
==== Equal-Effects Model ==== | ==== Equal-Effects Model ==== | ||
Let's start with an equal-effects model. As an example, consider the data from the meta-analysis by Molloy et al. (2014) on the relationship between conscientiousness and medication adherence. For each study, we can compute the r-to-z transformed correlation coefficient and corresponding sampling variance with: | Let's start with an equal-effects model. As an example, consider the data from the meta-analysis by Molloy et al. (2014) on the relationship between conscientiousness and medication adherence. For each study, we can compute the r-to-z transformed correlation coefficient and corresponding sampling variance with: | ||
+ | |||
<code rsplus> | <code rsplus> | ||
library(metafor) | library(metafor) | ||
Line 30: | Line 31: | ||
16 Wiebe & Christensen 1997 65 0.040 0.0400 0.0161 | 16 Wiebe & Christensen 1997 65 0.040 0.0400 0.0161 | ||
</ | </ | ||
+ | |||
The '' | The '' | ||
We can now fit an equal-effects model to these data with: | We can now fit an equal-effects model to these data with: | ||
+ | |||
<code rsplus> | <code rsplus> | ||
res.ee <- rma(yi, vi, data=dat, method=" | res.ee <- rma(yi, vi, data=dat, method=" | ||
Line 42: | Line 45: | ||
I^2 (total heterogeneity / total variability): | I^2 (total heterogeneity / total variability): | ||
H^2 (total variability / sampling variability): | H^2 (total variability / sampling variability): | ||
- | + | ||
- | Test for Heterogeneity: | + | Test for Heterogeneity: |
Q(df = 15) = 38.1595, p-val = 0.0009 | Q(df = 15) = 38.1595, p-val = 0.0009 | ||
Model Results: | Model Results: | ||
- | estimate | + | estimate |
- | 0.1252 | + | 0.1252 |
--- | --- | ||
Line 56: | Line 59: | ||
Now let's try to fit the same model with the '' | Now let's try to fit the same model with the '' | ||
+ | |||
<code rsplus> | <code rsplus> | ||
res.lm <- lm(yi ~ 1, weights = 1/vi, data=dat) | res.lm <- lm(yi ~ 1, weights = 1/vi, data=dat) | ||
Line 65: | Line 69: | ||
Weighted Residuals: | Weighted Residuals: | ||
- | Min 1Q Median | + | Min 1Q Median |
- | -3.1919 -1.0719 | + | -3.1919 -1.0719 |
Coefficients: | Coefficients: | ||
- | Estimate Std. Error t value Pr(> | + | Estimate Std. Error t value Pr(>|t|) |
(Intercept) | (Intercept) | ||
--- | --- | ||
Line 76: | Line 80: | ||
Residual standard error: 1.595 on 15 degrees of freedom | Residual standard error: 1.595 on 15 degrees of freedom | ||
</ | </ | ||
+ | |||
Two things are of note here: | Two things are of note here: | ||
- The estimated intercept ('' | - The estimated intercept ('' | ||
- | - The standard error (of the estimated intercept) is different (and hence, the test statistic and p-value also differs). | + | - The standard error (of the estimated intercept) is different (and hence, the test statistic and p-value also differ). |
The reason for this discrepancy is that the model fitted by the '' | The reason for this discrepancy is that the model fitted by the '' | ||
This can be demonstrated by extracting the estimated error variance from the '' | This can be demonstrated by extracting the estimated error variance from the '' | ||
+ | |||
<code rsplus> | <code rsplus> | ||
rma(yi, vi*sigma(res.lm)^2, | rma(yi, vi*sigma(res.lm)^2, | ||
Line 103: | Line 109: | ||
Signif. codes: | Signif. codes: | ||
</ | </ | ||
+ | |||
Now these results are exactly the same as those obtained by the '' | Now these results are exactly the same as those obtained by the '' | ||
Therefore, if we want to obtain the same standard error from '' | Therefore, if we want to obtain the same standard error from '' | ||
+ | |||
<code rsplus> | <code rsplus> | ||
coef(summary(res.ee)) | coef(summary(res.ee)) | ||
Line 119: | Line 127: | ||
[1] 0.01699805 | [1] 0.01699805 | ||
</ | </ | ||
+ | |||
+ | This is now the exact same standard error of the pooled estimate from an equal-effects model. | ||
==== Random-Effects Model ==== | ==== Random-Effects Model ==== | ||
A random-effects model can be fitted to the same data using the '' | A random-effects model can be fitted to the same data using the '' | ||
+ | |||
<code rsplus> | <code rsplus> | ||
res.re <- rma(yi, vi, data=dat) | res.re <- rma(yi, vi, data=dat) | ||
Line 135: | Line 146: | ||
H^2 (total variability / sampling variability): | H^2 (total variability / sampling variability): | ||
- | Test for Heterogeneity: | + | Test for Heterogeneity: |
Q(df = 15) = 38.1595, p-val = 0.0009 | Q(df = 15) = 38.1595, p-val = 0.0009 | ||
Model Results: | Model Results: | ||
- | estimate | + | estimate |
- | 0.1499 | + | 0.1499 |
--- | --- | ||
Line 148: | Line 159: | ||
Now let's try to fit the same model with the '' | Now let's try to fit the same model with the '' | ||
+ | |||
<code rsplus> | <code rsplus> | ||
library(nlme) | library(nlme) | ||
Line 156: | Line 168: | ||
<code output> | <code output> | ||
Linear mixed-effects model fit by REML | Linear mixed-effects model fit by REML | ||
- | Data: dat | + | Data: dat |
AIC | AIC | ||
-8.781569 -6.657418 7.390784 | -8.781569 -6.657418 7.390784 | ||
Line 167: | Line 179: | ||
Variance function: | Variance function: | ||
| | ||
- | | + | |
- | Fixed effects: yi ~ 1 | + | Fixed effects: yi ~ 1 |
Value Std.Error DF t-value p-value | Value Std.Error DF t-value p-value | ||
(Intercept) 0.1415557 0.03071906 16 4.608073 | (Intercept) 0.1415557 0.03071906 16 4.608073 | ||
Standardized Within-Group Residuals: | Standardized Within-Group Residuals: | ||
- | | + | |
- | -1.1596768 -0.6903414 | + | -1.1596768 -0.6903414 |
Number of Observations: | Number of Observations: | ||
Line 181: | Line 193: | ||
Alternatively, | Alternatively, | ||
+ | |||
<code rsplus> | <code rsplus> | ||
library(lme4) | library(lme4) | ||
- | res.lmer <- lmer(yi ~ 1 + (1 | study), weights = 1/vi, data=dat, | + | res.lmer <- lmer(yi ~ 1 + (1 | study), weights = 1/vi, data=dat, |
| | ||
summary(res.lmer) | summary(res.lmer) | ||
Line 192: | Line 205: | ||
Data: dat | Data: dat | ||
Weights: 1/vi | Weights: 1/vi | ||
- | Control: | + | Control: |
lmerControl(check.nobs.vs.nlev = " | lmerControl(check.nobs.vs.nlev = " | ||
REML criterion at convergence: | REML criterion at convergence: | ||
- | Scaled residuals: | + | Scaled residuals: |
- | Min 1Q Median | + | Min 1Q Median |
- | -1.1597 -0.6903 | + | -1.1597 -0.6903 |
Random effects: | Random effects: | ||
| | ||
- | | + | |
- | | + | |
Number of obs: 16, groups: | Number of obs: 16, groups: | ||
Line 211: | Line 224: | ||
(Intercept) | (Intercept) | ||
</ | </ | ||
+ | |||
Note that some of the default checks need to be switched off before '' | Note that some of the default checks need to be switched off before '' | ||
A couple things are of note here: | A couple things are of note here: | ||
- The estimated intercept differs from the model coefficient obtained earlier with the '' | - The estimated intercept differs from the model coefficient obtained earlier with the '' | ||
- | - The standard error is also different (and hence, the test statistic and p-value also differs). | + | - The standard error is also different (and hence, the test statistic and p-value also differ). |
- The estimated standard deviation of the random effects also differs ('' | - The estimated standard deviation of the random effects also differs ('' | ||
- | These discrepancies | + | These discrepancies |
To illustrate this, we can again factor in that constant into the sampling variances and refit the model with '' | To illustrate this, we can again factor in that constant into the sampling variances and refit the model with '' | ||
+ | |||
<code rsplus> | <code rsplus> | ||
rma(yi, vi*sigma(res.lme)^2, | rma(yi, vi*sigma(res.lme)^2, | ||
Line 232: | Line 247: | ||
H^2 (total variability / sampling variability): | H^2 (total variability / sampling variability): | ||
- | Test for Heterogeneity: | + | Test for Heterogeneity: |
Q(df = 15) = 24.0532, p-val = 0.0642 | Q(df = 15) = 24.0532, p-val = 0.0642 | ||
Model Results: | Model Results: | ||
- | estimate | + | estimate |
- | 0.1416 | + | 0.1416 |
--- | --- | ||
Signif. codes: | Signif. codes: | ||
</ | </ | ||
+ | |||
Now these results match exactly what was obtained with the '' | Now these results match exactly what was obtained with the '' | ||
Line 249: | Line 265: | ||
==== Meta-Regression ==== | ==== Meta-Regression ==== | ||
- | The examples above show models without any covariates/ | + | The examples above show models without any covariates/ |
The possibility of using both an additive and multiplicative component simultaneously was not discussed by Thompson and Sharp (1999). This is in fact what the model fitted with the '' | The possibility of using both an additive and multiplicative component simultaneously was not discussed by Thompson and Sharp (1999). This is in fact what the model fitted with the '' | ||
Line 255: | Line 271: | ||
==== More Complex Models ==== | ==== More Complex Models ==== | ||
- | Note that the same issues apply when fitting more complex meta-analytic models (e.g., multivariate/ | + | Note that the same issue arises |
==== S-Plus Version of lme() ==== | ==== S-Plus Version of lme() ==== | ||
- | The '' | + | The '' |
<code rsplus> | <code rsplus> | ||
res.lme <- lme(yi ~ 1, random = ~ 1 | study, weights = varFixed(~ vi), control=lmeControl(sigma = 1), data=dat) | res.lme <- lme(yi ~ 1, random = ~ 1 | study, weights = varFixed(~ vi), control=lmeControl(sigma = 1), data=dat) | ||
Line 266: | Line 283: | ||
<code output> | <code output> | ||
Linear mixed-effects model fit by REML | Linear mixed-effects model fit by REML | ||
- | Data: dat | + | Data: dat |
- | AIC BIC | + | AIC BIC |
-10.44653 -9.03043 7.223265 | -10.44653 -9.03043 7.223265 | ||
Random effects: | Random effects: | ||
| | ||
- | (Intercept) Residual | + | (Intercept) Residual |
StdDev: | StdDev: | ||
Variance function: | Variance function: | ||
| | ||
- | | + | |
- | Fixed effects: yi ~ 1 | + | Fixed effects: yi ~ 1 |
- | Value Std.Error DF t-value p-value | + | Value Std.Error DF t-value p-value |
(Intercept) 0.1499171 0.03155994 16 4.750234 | (Intercept) 0.1499171 0.03155994 16 4.750234 | ||
Standardized Within-Group Residuals: | Standardized Within-Group Residuals: | ||
- | | + | |
| | ||
Line 289: | Line 306: | ||
Number of Groups: 16 | Number of Groups: 16 | ||
</ | </ | ||
+ | |||
These are the exact same results as obtained earlier with the '' | These are the exact same results as obtained earlier with the '' | ||
- | **Update:** The R version of the '' | + | **Update:** The R version of the '' |
==== Summary ==== | ==== Summary ==== |
tips/rma_vs_lm_lme_lmer.txt · Last modified: 2023/11/14 08:00 by Wolfgang Viechtbauer