DOC PREVIEW
UW-Madison STAT 572 - Random Effects in R

This preview shows page 1-2 out of 7 pages.

Save
View full document
View full document
Premium Document
Do you want full access? Go Premium and unlock all 7 pages.
Access to all documents
Download any document
Ad free experience
View full document
Premium Document
Do you want full access? Go Premium and unlock all 7 pages.
Access to all documents
Download any document
Ad free experience
Premium Document
Do you want full access? Go Premium and unlock all 7 pages.
Access to all documents
Download any document
Ad free experience

Unformatted text preview:

The Big PictureRandom Effects and BalanceComputationlmerRandom Effects in RBret LargetDepartments of Botany and of StatisticsUniversity of Wisconsin—MadisonMarch 15, 2007Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 1 / 13The Big Picture Random Effects and BalanceThe Big PictureThe completely randomized design with a random effect assumes thefollowing model.yij= µ + ai+ eij, j = 1, . . . , ni, i = 1, . . . , k,where ai∼ iidN(0, σ2a) and eij∼ iidN(0, σ2e) with {ai} independentof {eij}.If the data is perfectly balanced with equal numbers of individuals ineach group, (ni= n for all i), it is possible to make modifications to afixed effects ANOVA analysis by computing a different expected meansquare (EMS) for inference.When the data is not balanced, this approach fails — sometimesapproximations work, sometimes not.An alternative regression approach is appropriate even when thedesign is imbalanced, but the scope of inference is then different.I will first show the regression approach (described in Chapter 10 ofyour book).Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 2 / 13The Big Picture Random Effects and BalanceCorrelationOne consequence of the random effects model is that observationscan be correlated.In the model we study here, observations in different groups areindependent just as in the fixed effect model, because they do notshare any random variables.For example, the first observations from groups one, y11, depends a1and e11while y21depends on a2and e21.In contrast, observations in the same group are correlated.y11and y12each depend on a1.Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 3 / 13The Big Picture Random Effects and BalanceVarianceRecall from Statistics 571 the following facts about variances ofindependent random variables X and Y and constant c.1V(X + Y ) = V(X −Y ) = V(X ) + V(Y ) (if independent).2V(cX ) = c2V(X ).With this reminder, we can recognize that σameasures the size of atypical difference between a random effect and the grand mean.A typical difference between two regression effects ispV(ai− aj) =pV(ai) + V(aj) =√2σa.Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 4 / 13Computation lmerLinear 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.Mac users will need a recent operating system, OS X 10.4.4 or laterand R 2.4.1 or later to use this library.Users of older Macs can use the older package nlme and the functionlme.lme cannot fit as rich a class of random effects models as lmer (forexample, random effects cannot be nested and you cannot usegeneralized linear models), but it will suffice for much of what we doin the course.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.Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 5 / 13Computation lmerData DescriptionWe will consider the data set ant111b from the textbook.This data set is a subset of a much larger data set on corn yields ondifferent islands in the Carribean.In the subset of data, two possibly interesting response variables areharvest weight (harvwt) and the number of ears (ears).Both harvwt and ears are measured by plot, but the units areunspecified.There are four plots within each site. The variable plot is onlymeaningful nested within a site.We will consider the site as a random effect as there are manypossible sites on each island.In the subset of the data we have, there is only a single treatment andall sites are on the same island.Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 6 / 13Computation lmerData> library(DAAG)> str(ant111b)'data.frame': 32 obs. of 9 variables:$ site : Factor w/ 8 levels "DBAN","LFAN",..: 1 2 3 4 5 6 7 8 1 2 ...$ parcel: Factor w/ 4 levels "I","II","III",..: 1 1 1 1 1 1 1 1 2 2 ...$ code : num 58 58 58 58 58 58 58 58 58 58 ...$ island: num 1 1 1 1 1 1 1 1 1 1 ...$ id : num 3 40 186 256 220 ...$ plot : num 3 4 5.5 4.5 3.5 5 7 7 15.5 15 ...$ trt : num 111 111 111 111 111 111 111 111 111 111 ...$ ears : num 43.5 40.5 20 42.5 31.5 32.5 43.5 50 46 46.5 ...$ harvwt: num 5.16 2.93 1.73 6.79 3.25 ...Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 7 / 13Computation lmerPlot of Data> attach(ant111b)> library(lattice)> fig1 = xyplot(harvwt ~ site, pch = 16, col = "blue")> print(fig1)Notice that there is considerablevariation among sites.Variation within sites (from plotto plot) appears smaller.siteharvwt234567DBAN LFAN NSAN ORAN OVAN TEAN WEAN WLAN●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 8 / 13Computation lmerFitting a Random Effects Model> library(lme4)> corn1.lmer = lmer(harvwt ~ 1 + (1 | site))> summary(corn1.lmer)Linear mixed-effects model fit by REMLFormula: harvwt ~ 1 + (1 | site)AIC BIC logLik MLdeviance REMLdeviance98.42 101.3 -47.21 95.08 94.42Random effects:Groups Name Variance Std.Dev.site (Intercept) 2.36773 1.53874Residual 0.57754 0.75996number of obs: 32, groups: site, 8Fixed effects:Estimate Std. Error t value(Intercept) 4.2917 0.5604 7.659> sig2site = as.vector(VarCorr(corn1.lmer)$site)> sigma = attributes(summary(corn1.lmer))$sigmaThe first 1 is the fixed effect.The term (1 | site) meansthat there is a random effect foreach site and this effect isnested within the intercept (thewhole model).There are two sources ofrandom variation, one for siteand one for parcel within site.The estimated variancecomponents are σ2a= 2.36773and σ2e= 0.57754.Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 9 / 13Computation lmerResidual Plot> plot(fitted(corn1.lmer), residuals(corn1.lmer))> abline(h = 0)There are no bad patterns in theresidual plot.●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●●3 4 5 6−1.5 −1.0 −0.5 0.0 0.5 1.0 1.5fitted(corn1.lmer)residuals(corn1.lmer)Statistics 572 (Spring 2007) Random Effects in R March 17, 2007 10 / 13Computation lmerFitted Values> means = sapply(split(harvwt, site), mean)> siteFit = sapply(split(fitted(corn1.lmer), site), mean)> data.frame(mean = means, fitted = siteFit)mean fittedDBAN 4.88500 4.850901LFAN


View Full Document

UW-Madison STAT 572 - Random Effects in R

Download Random Effects in R
Our administrator received your request to download this document. We will send you the file to your email shortly.
Loading Unlocking...
Login

Join to view Random Effects in R and access 3M+ class-specific study document.

or
We will never post anything without your permission.
Don't have an account?
Sign Up

Join to view Random Effects in R 2 2 and access 3M+ class-specific study document.

or

By creating an account you agree to our Privacy Policy and Terms Of Use

Already a member?