How do I do a Variance Components Analysis? Example Problem Given: A gage study is run with three randomly selected OPERATORS, two randomly selected FIXTURES and five randomly selected PARTS. All possible combinations are studied (full factorial or fully crossed experiment). The response is some measurement take for each of the 30 combinations. How do I use R to find an estimate for the contribution each factor makes to the overall measurement variation? In other words how do I do a Variance Components analysis in R? The model is Expected(Var(Measurements)) = Var(due to OPERATORS) + Var(due to FIXTURES) + Var(due to PARTS)+Var(due to NOISE) -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
VarComp in R
3 messages · Richard Scott, Douglas Bates, Robert Scott Hoy
Richard Scott <rscott at rochester.rr.com> writes:
How do I do a Variance Components Analysis? Example Problem Given: A gage study is run with three randomly selected OPERATORS, two randomly selected FIXTURES and five randomly selected PARTS. All possible combinations are studied (full factorial or fully crossed experiment). The response is some measurement take for each of the 30 combinations. How do I use R to find an estimate for the contribution each factor makes to the overall measurement variation? In other words how do I do a Variance Components analysis in R? The model is Expected(Var(Measurements)) = Var(due to OPERATORS) + Var(due to FIXTURES) + Var(due to PARTS)+Var(due to NOISE)
You can do this with lme but, as I mentioned in our private correspondence, lme is better designed for nested random effects than for crossed random effects like this. You need to use the pdIdent form of the variance-covariance matrix on the matrix of indicator columns, created with formulas like ~ PARTS - 1, and combine those with pdBlocked. Do you really have a fully crossed design with random effects? It seems strange to pick the PART at random then reuse that particular part 6 times. Anyway, here is how you would do it.
gage <- data.frame(OPERATORS = gl(3, 1, 30), FIXTURES = gl(2, 3, 30),
+ PARTS = gl(5, 6, 30))
summary(gage)
OPERATORS FIXTURES PARTS
1:10 1:15 1:6
2:10 2:15 2:6
3:10 3:6
4:6
5:6
gage[1:10,]
OPERATORS FIXTURES PARTS 1 1 1 1 2 2 1 1 3 3 1 1 4 1 2 1 5 2 2 1 6 3 2 1 7 1 1 2 8 2 1 2 9 3 1 2 10 1 2 2
gage$Grp <- gl(1,1,30) # trivial factor representing a single group # create the simulated random effects Op <- rnorm(3, 0, 1); Fix <- rnorm(2, 0, 2); Par <- rnorm(5, 0, 0.5) attach(gage) gage$Meas <- Op[OPERATORS] + Fix[FIXTURES] + Par[PARTS] + rnorm(30, 12, 0.3) detach() gr <- groupedData(Meas ~ 1 | Grp, gage) fm <- lme(Meas ~ 1, data = gr,
+ random = pdBlocked(list(pdIdent(~ OPERATORS - 1), + pdIdent(~ FIXTURES - 1), + pdIdent(~ PARTS - 1))))
fm
Linear mixed-effects model fit by REML
Data: gr
Log-restricted-likelihood: -13.95199
Fixed: Meas ~ 1
(Intercept)
12.73601
Random effects:
Composite Structure: Blocked
Block 1: OPERATORS1, OPERATORS2, OPERATORS3
Formula: ~OPERATORS - 1 | Grp
Structure: Multiple of an Identity
OPERATORS1 OPERATORS2 OPERATORS3
StdDev: 0.5424059 0.5424059 0.5424059
Block 2: FIXTURES1, FIXTURES2
Formula: ~FIXTURES - 1 | Grp
Structure: Multiple of an Identity
FIXTURES1 FIXTURES2
StdDev: 1.723611 1.723611
Block 3: PARTS1, PARTS2, PARTS3, PARTS4, PARTS5
Formula: ~PARTS - 1 | Grp
Structure: Multiple of an Identity
PARTS1 PARTS2 PARTS3 PARTS4 PARTS5 Residual
StdDev: 0.3006655 0.3006655 0.3006655 0.3006655 0.3006655 0.245659
Number of Observations: 30
Number of Groups: 1
VarCorr(fm)
Grp = pdIdent(OPERATORS - 1), pdIdent(FIXTURES - 1), pdIdent(PARTS - 1)
Variance StdDev
OPERATORS1 0.29420413 0.5424059
OPERATORS2 0.29420413 0.5424059
OPERATORS3 0.29420413 0.5424059
FIXTURES1 2.97083530 1.7236111
FIXTURES2 2.97083530 1.7236111
PARTS1 0.09039977 0.3006655
PARTS2 0.09039977 0.3006655
PARTS3 0.09039977 0.3006655
PARTS4 0.09039977 0.3006655
PARTS5 0.09039977 0.3006655
Residual 0.06034835 0.2456590
Douglas Bates bates at stat.wisc.edu Statistics Department 608/262-2598 University of Wisconsin - Madison http://www.stat.wisc.edu/~bates/ -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
Hi. I've just installed R on Mac OS X, and when I try to run the shell script I get the message dyld: /Sci_Math_Apps/R/locl/bin/R/bin/R.bin can't open library:libz.A.dylib.1.1.3 (No such file or directory, errno = 2) I searched for the file "libz.A.dylib.1.1.3" but couldn't find it. I found a lot of dylib files in various places but couldn't find it. Any ideas? Thanks, Rob -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._