Multilevel ModelsCorn ExampleComputationlmerMultilevel StructuresBret LargetDepartments of Botany and of StatisticsUniversity of Wisconsin—MadisonApril 17, 20081 / 15Data DescriptionWe consider a subset of a larger data set on corn grown on the islandAntigua.The response variable we consider is the harvest weight (harvwt) perplot (units unknown).There are eight sites with eight separate plots within each site wherethe corn is grown under the same treatment conditions.We can ask if the site has an effect on the harvest weight.In a standard regression framework, we could analyze the data as aone-way ANOVA with eight fixed parameters for the expected values(an intercept which is the mean of a reference group and sevendifferences in means between the other groups and the reference) anda single plot-level source of error.In a multilevel model, we can have covariates and error associatedwith the plot level and separate covariates and error associated withthe site level.Multilevel Models Corn Example 2 / 15ModelsStandard ANOVA model:yi= β1+ β2· 1(site 2) + · · · + β8· 1(site 8) + eiwhere i = 1, . . . , 64 indexes the observation.Iei∼ iid N(0, σ2);Iβj, j = 1, . . . , 8 and σ2are fixed parameters.In a multilevel model, we may haveyi= αj[i]+ eiwhere i = 1, . . . , 64 indexes the observation and j[i] = 1, . . . , 8indicates which of the eight sites contains the ith observation.Iαj, j = 1, . . . , 8, ∼ N(µα, σ2α) are random effects for the sites;Iei∼ iid N(0, σ2);Iµα, σα, and σ2are fixed and unknown.Notice here that we have a regression model for the response, and alsoa regression model for the coeficients of the first regression model.Multilevel models include sources of variation at more than one level.Multilevel Models Corn Example 3 / 15Data> corn = read.table("corn.txt", header = T)> summary(corn)site block ears harvwtDBAN : 8 I :16 Min. :13.00 Min. :1.280LFAN : 8 II :16 1st Qu.:37.75 1st Qu.:2.935NSAN : 8 III:16 Median :43.00 Median :4.300ORAN : 8 IV :16 Mean :41.22 Mean :4.292OVAN : 8 3rd Qu.:46.00 3rd Qu.:5.442TEAN : 8 Max. :58.00 Max. :7.530(Other):16Multilevel Models Corn Example 4 / 15Plot of Datasiteharvwt246DBAN LFAN NSAN ORAN OVAN TEAN WEAN WLAN●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●Multilevel Models Corn Example 5 / 15Standard Regression Model> corn.lm = lm(harvwt ~ site, data = corn)> display(corn.lm)lm(formula = harvwt ~ site, data = corn)coef.est coef.se(Intercept) 4.89 0.31siteLFAN -0.68 0.44siteNSAN -2.79 0.44siteORAN 2.03 0.44siteOVAN -0.05 0.44siteTEAN -1.85 0.44siteWEAN 0.64 0.44siteWLAN -2.04 0.44---n = 64, k = 8residual sd = 0.87, R-Squared = 0.77Multilevel Models Corn Example 6 / 15Multilevel Model> corn.lmer = lmer(harvwt ~ (1 | site), data = corn)> display(corn.lmer)lmer(formula = harvwt ~ (1 | site), data = corn)coef.est coef.se(Intercept) 4.29 0.56Error terms:Groups Name Std.Dev.site 1.55Residual 0.87---number of obs: 64, groups: site, 8AIC = 192.9, DIC = 190.3deviance = 189.6Multilevel Models Corn Example 7 / 15Comparing ModelsDiscuss the different parameter estimates on the board.Multilevel Models Corn Example 8 / 15Linear Mixed Effects Models using lmerThe most recently developed R package for fitting linear models withrandom effects is in the library lme4.The function to use instead of lm is named lmer.A model formula with a random effect in lmer differs from lm byincluding a term of the form (a | b) where a is a model matrix(often the intercept 1) for the scope of the random effect and b is thegroup to which the random effect applies.Computation lmer 9 / 15Residual Plotfitted(corn.lmer)residuals(corn.lmer)−3−2−10122 3 4 5 6 7●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●Computation lmer 10 / 15Compare Fitted Values> means = with(corn, sapply(split(harvwt, site), mean))> fitted.lm = with(corn, sapply(split(fitted(corn.lm), site),+ mean))> fitted.lmer = with(corn, sapply(split(fitted(corn.lmer), site),+ mean))> signif(rbind(means, fitted.lm, fitted.lmer), 3)DBAN LFAN NSAN ORAN OVAN TEAN WEAN WLANmeans 4.88 4.21 2.09 6.92 4.83 3.04 5.53 2.84fitted.lm 4.88 4.21 2.09 6.92 4.83 3.04 5.53 2.84fitted.lmer 4.86 4.21 2.17 6.82 4.81 3.08 5.48 2.90The overall mean is 4.29.The multilevel model shrinks the estimates toward the overall mean.Computation lmer 11 / 15Other ExamplesMultilevel models are often used in these situations.IRepeated measures — when a single individual is measured multipletimes, it is often appropriate to model two levels of variation, one forindividuals and one for measurements.ISplit-plot designs — in agricultural or ecological studies, it is often thecase that sites are broken into plots and possibly subplots. Variablescan be measured at the site, plot, subplot, or individual measurementlevel.IMultilevel models are also appropriate for non-nested variables. Forexample, measurements could be clustered by year and by site if asingle site is measured over multiple years.Computation lmer 12 / 15Summary of Classical RegressionPrediction for continuous and discrete outcomes;Fitting nonlinear relationships using transformations;Inclusion of categorical predictors with indicator random variables;Modeling interactions;Causal inference.Computation lmer 13 / 15Motivations for Multilevel ModelsAccounting for both individual and group level variation in estimatinggroup-level effects.Modeling individual level regression coefficients.Estimation of effects for subgroups.Computation lmer 14 / 15When is it worth fitting multilevel models?If the group size is small, there may not be much data to estimaterandom effects and there is little to gain.The complexity of multilevel models is greater than classicalregression. The added complexity is often worthwhile, but perhapsnot when there are only a small number (say less than five)individuals in a group.Computation lmer 15 /