aov.sufficient
Analysis of variance from sufficient statistics for groups.
Analysis of variance from sufficient statistics for groups.
For each group, we need the factor level, the response mean, the
within-group standard deviation, and the sample size.
The correct ANOVA table is produced. The residuals are fake.
The generic vcov
and summary.lm
don't work for the
variance of the regression coefficients in this case.
Use vcov.sufficient
.
- Keywords
- htest
Usage
aov.sufficient(formula, data = NULL,
projections = FALSE, qr = TRUE, contrasts = NULL,
weights = data$n, sd = data$s,
...)vcov.sufficient(object, ...)
Arguments
Value
- For
aov.sufficient
, an object of class c("aov", "lm"). Forvcov.sufficient
, a function that returns the covariance matrix of the regression coefficients.
Note
The residuals are fake. They are all identical and equal to the MLE
standard error (sqrt(SumSq.res/df.tot)
). They give the right
ANOVA table. They may cause confusion or warnings in other programs.
The standard errors and t-tests of the coefficients are not calculated
by summary.lm
.
Using the aov
object from aov.sufficient
in glht
requires the vcov.
and df
arguments.
See Also
Examples
## This example is from Hsu and Peruggia
## This is the R version
## See ?mmc.mean for S-Plus
if.R(s={},
r={
pulmonary <- read.table(hh("datasets/pulmonary.dat"), header=TRUE,
row.names=NULL)
names(pulmonary)[3] <- "FVC"
names(pulmonary)[1] <- "smoker"
pulmonary$smoker <- factor(pulmonary$smoker, levels=pulmonary$smoker)
row.names(pulmonary) <- pulmonary$smoker
pulmonary
pulmonary.aov <- aov.sufficient(FVC ~ smoker,
data=pulmonary)
summary(pulmonary.aov)
pulmonary.mca <- glht(pulmonary.aov,
linfct=mcp(smoker="Tukey"),
df=pulmonary.aov$df.residual,
vcov.=vcov.sufficient)
old.omd <- par(omd=c(.03,1,0,1))
plot(pulmonary.mca)
par(old.omd)
pulm.lmat <- cbind("npnl-mh"=c( 1, 1, 1, 1,-2,-2), ## not.much vs lots
"n-pnl" =c( 3,-1,-1,-1, 0, 0), ## none vs light
"p-nl" =c( 0, 2,-1,-1, 0, 0), ## {} arbitrary 2 df
"n-l" =c( 0, 0, 1,-1, 0, 0), ## {} for 3 types of light
"m-h" =c( 0, 0, 0, 0, 1,-1)) ## moderate vs heavy
dimnames(pulm.lmat)[[1]] <- row.names(pulmonary)
if.R(r={pulm.lmat.glht <- rbind(Int=0, pulm.lmat[-1,])
pulm.lmat.glht},
s={})
pulm.lmat
pulmonary.mmc <- glht.mmc(pulmonary.aov,
linfct=mcp(smoker="Tukey"),
df=pulmonary.aov$df.residual,
vcov.=vcov.sufficient,
lmat=pulm.lmat.glht,
focus.lmat=pulm.lmat,
calpha=attr(confint(pulmonary.mca)$confint,"calpha"))
old.omd <- par(omd=c(.03,.95,0,1))
plot(pulmonary.mmc, print.mca=TRUE, print.lmat=FALSE)
## tiebreaker plot, with contrasts ordered to match MMC plot,
## with all contrasts forced positive and with names also reversed,
## and with matched x-scale.
plot.matchMMC(pulmonary.mmc$mca)
## orthogonal contrasts
plot(pulmonary.mmc, print.lmat=TRUE, col.lmat.signif='blue', col.iso='gray')
plot.matchMMC(pulmonary.mmc$lmat)
## pairwise and orthogonal contrasts on the same plot
plot(pulmonary.mmc, print.mca=TRUE, print.lmat=TRUE,
col.mca.signif='red', col.lmat.signif='blue', col.iso='gray',
lty.lmat.not.signif=2)
par(old.omd)
})