On Thu, 24 Jan 2008, peter salzman wrote: > Dear list, > > i'm trying to test if a linear combination of coefficients of glm is equal > to 0. For example : > class 'cl' has 3 levels (1,2,3) and 'y' is a response variable. We want to > test H0: mu1 + mu2 - mu3 =0 where mu1,mu2, and mu3 are the means for each > level.
Look at linear.hypothesis() in package "car". ## data (reduce treatment to three levels) data("OrchardSprays") os <- subset(OrchardSprays, treatment %in% c("A", "B", "C")) os$treatment <- factor(os$treatment) ## model and test fm <- lm(decrease ~ 0 + treatment, data = os) coef(fm) linear.hypothesis(fm, "treatmentA + treatmentB - treatmentC = 0") See the accompanying documentation for an overview which extractor functions (such as coef, vcov, etc.) are used to compute the test. hth, Z > for me, the question is how to get the covariance matrix of the estimated > parameters from glm. but perhaps there is a direct solution in one of the > packages. > > i know how to solve this particular problem (i wrote it below) but i'm > curious about the covariance matrix of coefficient as it seems to be > important. > > the R code example : > ### > nObs <- 10 > cl <- as.factor( sample(c(1,2,3),nObs,replace=TRUE) ) > y <- rnorm(nObs) > > model <- glm(y ~ cl) > b <- model$coefficients > H <- c(1,1,-1) # we want to test H0: Hb = 0 > > ### the following code will NOT run unless we can compute covModelCoeffs > > #the mean of Hb is > mu = H %*% model$coefficients > #the variance is HB is > var = H %*% covModelCoeffs %*% t(H) > > p.val <- 2 * pnorm( -abs(mu), mean=0, sd=sqrt(var),lower.tail = TRUE) > > > how do i get the covariance matrix of the estimated parameters ? > > thanks, > peter > > P.S. the simple solution for this particular problem: > > ## get the mean for each level > muV <- by(y,cl,mean) > ## get the variance for each level > varV <- by(y,cl,var) > > ## the mean of Hb is > muHb <- H %*% muV > ## because of independence, the variance of Hb is > varHb <- sum(varV) > > ## the probability of error, so-called p-value: > p.val <- 2 * pnorm( -abs(muHb), mean=0, sd=sqrt(varHb),lower.tail = TRUE) > > thanks again, > peter > > > -- > Peter Salzman, PhD > Department of Biostatistics and Computational Biology > University of Rochester > > [[alternative HTML version deleted]] > > ______________________________________________ > R-help@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > > ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.