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.
>
> 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.
See ?vcov .
BTW, help.search("covariance matrix") finds it.
>
> 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
>
>
>
--
Brian D. Ripley, [EMAIL PROTECTED]
Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/
University of Oxford, Tel: +44 1865 272861 (self)
1 South Parks Road, +44 1865 272866 (PA)
Oxford OX1 3TG, UK Fax: +44 1865 272595
______________________________________________
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.
