# NOT RUN {
# load data (CLSI EP05-A2 Within-Lab Precision Experiment)
data(dataEP05A2_1)
# perform (V)variance (C)component (A)nalysis (also comute A-matrices)
res <- anovaVCA(y~day/run, dataEP05A2_1)
# get confidence intervals for total and error (VC, SD, CV)
VCAinference(res)
# additionally request CIs for all other VCs; default is to constrain
# CI-limits to be >= 0
# first solve MME
res <- solveMME(res)
VCAinference(res, VarVC=TRUE)
# now using Satterthwaite methodology for CIs
VCAinference(res, VarVC=TRUE, ci.method="satt")
# request unconstrained CIs
VCAinference(res, VarVC=TRUE, constrainCI=FALSE)
# additionally request Chi-Squared Tests of total and error, default
# is that claim values are specified as variances (claim.type="VC")
VCAinference(res, total.claim=4.5, error.claim=3.5)
# perform Chi-Squared Tests, where claim-values are given as SD,
# compare p-values to former example
VCAinference(res, total.claim=sqrt(4.5), error.claim=sqrt(3.5), claim.type="SD")
# now using Satterthwaite methodology for CIs
VCAinference(res, total.claim=sqrt(4.5), error.claim=sqrt(3.5),
claim.type="SD", ci.method="satt")
# now add random error to example data forcing the ANOVA-estimate of the
# day-variance to be negative
set.seed(121)
tmpData <- dataEP05A2_1
tmpData$y <- tmpData$y + rnorm(80,,3)
res2 <- anovaVCA(y~day/run, tmpData)
# call 'VCAinference' with default settings
VCAinference(res2)
# extract components of the returned 'VCAinference' object
inf <- VCAinference(res2, total.claim=12)
inf$ConfInt$VC$OneSided # one-sided CIs for variance components
inf$ConfInt$VC$TwoSided # two-sided CI for variance components
inf$ChiSqTest
# request CIs for all VCs, default is to exclude CIs of negative VCs (excludeNeg=TRUE)
# solve MMEs first (or set MME=TRUE when calling anovaVCA)
res2 <- solveMME(res2)
VCAinference(res2, VarVC=TRUE)
# request CIs for all VCs, including those for negative VCs, note that all CI-limits
# are constrained to be >= 0
VCAinference(res2, VarVC=TRUE, excludeNeg=FALSE)
# request unconstrained CIs for all VCs, including those for negative VCS
# one has to re-fit the model allowing the VCs to be negative
res3 <- anovaVCA(y~day/run, tmpData, NegVC=TRUE, MME=TRUE)
VCAinference(res3, VarVC=TRUE, excludeNeg=FALSE, constrainCI=FALSE)
### use the numerical example from the CLSI EP05-A2 guideline (p.25)
data(Glucose)
res.ex <- anovaVCA(conc~day/run, Glucose)
### also perform Chi-Squared tests
### Note: in guideline claimed SD-values are used, here, claimed variances are used
VCAinference(res.ex, total.claim=3.4^2, error.claim=2.5^2)
# load another example dataset and extract the "sample_1" subset
data(VCAdata1)
sample1 <- VCAdata1[which(VCAdata1$sample==1),]
# generate an additional factor variable and random errors according to its levels
sample1$device <- gl(3,28,252)
set.seed(505)
sample1$y <- sample1$y + rep(rep(rnorm(3,,.25), c(28,28,28)),3)
# fit a crossed-nested design with main factors 'lot' and 'device'
# and nested factors 'day' and 'run' nested below, also request A-matrices
res1 <- anovaVCA(y~(lot+device)/day/run, sample1)
# get confidence intervals, covariance-matrix of VCs, ...,
# explicitly request the covariance-matrix of variance components
# solve MMEs first
res1 <- solveMME(res1)
inf1 <- VCAinference(res1, VarVC=TRUE, constrainCI=FALSE)
inf1
# print numerical values with more digits
print(inf1, digit=12)
# print only parts of the 'VCAinference' object (see \code{\link{print.VCAinference}})
print(inf1, digit=12, what=c("VCA", "VC"))
# extract complete covariance matrix of variance components
# (main diagonal is part of standard output -> "Var(VC"))
VarCovVC <- vcovVC(inf1$VCAobj)
round(VarCovVC, 12)
# use by-processing and specific argument-values for each level of the by-variable
data(VCAdata1)
fit.all <- anovaVCA(y~(device+lot)/day/run, VCAdata1, by="sample", NegVC=TRUE)
inf.all <- VCAinference(fit.all, total.claim=c(.1,.75,.8,1,.5,.5,2.5,20,.1,1))
print.VCAinference(inf.all, what="VC")
# }
Run the code above in your browser using DataLab