Comment effectuer un test d'hypothèse linéaire sur les coefficients de régression avec une matrice de covariance en grappes?

Question

Je suis intéressé par le calcul des estimations et des erreurs standard pour des combinaisons linéaires de coefficients après une régression linéaire dans R. Par exemple, supposons que j'ai la régression et le test:

data(mtcars) 
library(multcomp) 
lm1 <- lm(mpg ~ cyl + hp, data = mtcars) 
summary(glht(lm1, linfct = 'cyl + hp = 0')) 

Cette estimera la valeur de la somme de les coefficients sur cyl et hp, et fournissent l'erreur standard basée sur la matrice de covariance produite par lm.

Mais, je suppose que veux regrouper mes erreurs types, sur une troisième variable:

data(mtcars) 
library(multcomp) 
library(lmtest) 
library(multiwayvcov) 
lm1 <- lm(mpg ~ cyl + hp, data = mtcars) 
vcv <- cluster.vcov(lm1, cluster = mtcars$am) 
ct1 <- coeftest(lm1,vcov. = vcv) 

ct1 contient le regroupement pour mon SEs par am. Cependant, si je tente d'utiliser l'objet ct1 dans glht, vous obtenez une erreur indiquant

Erreur dans modelparm.default (modèle, ...): pas de méthode « coef » pour « modèle » trouvé! Un conseil sur la façon de faire l'hypothèse linéaire avec la matrice de covariance de la variance clusterisée?

Merci!

Answer 1

glht(ct1, linfct = 'cyl + hp = 0') ne fonctionnera pas, car ct1 n'est pas un objet glht et ne peut pas être forcé par as.glht. Je ne sais pas s'il existe un paquet ou une fonction existante pour le faire, mais ce n'est pas un travail difficile à faire nous-mêmes. La petite fonction qui suit fait:

LinearCombTest <- function (lmObject, vars, .vcov = NULL) { 
    ## if `.vcov` missing, use the one returned by `lm` 
    if (is.null(.vcov)) .vcov <- vcov(lmObject) 
    ## estimated coefficients 
    beta <- coef(lmObject) 
    ## sum of `vars` 
    sumvars <- sum(beta[vars]) 
    ## get standard errors for sum of `vars` 
    se <- sum(.vcov[vars, vars])^0.5 
    ## perform t-test on `sumvars` 
    tscore <- sumvars/se 
    pvalue <- 2 * pt(abs(tscore), lmObject$df.residual, lower.tail = FALSE) 
    ## return a matrix 
    matrix(c(sumvars, se, tscore, pvalue), nrow = 1L, 
     dimnames = list(paste0(paste0(vars, collapse = " + "), " = 0"), 
         c("Estimate", "Std. Error", "t value", "Pr(>|t|)"))) 
    } 

Ayons un test:

data(mtcars) 
lm1 <- lm(mpg ~ cyl + hp, data = mtcars) 
library(multiwayvcov) 
vcv <- cluster.vcov(lm1, cluster = mtcars$am) 

Si nous laissons .vcov non spécifiée dans LinearCombTest, il est le même que multcomp::glht:

LinearCombTest(lm1, c("cyl","hp")) 

#    Estimate Std. Error t value  Pr(>|t|) 
#cyl + hp = 0 -2.283815 0.5634632 -4.053175 0.0003462092 

library(multcomp) 
summary(glht(lm1, linfct = 'cyl + hp = 0')) 

#Linear Hypotheses: 
#    Estimate Std. Error t value Pr(>|t|)  
#cyl + hp == 0 -2.2838  0.5635 -4.053 0.000346 *** 

Si nous fournir une covariance, il fait ce que vous voulez:

LinearCombTest(lm1, c("cyl","hp"), vcv) 

#    Estimate Std. Error t value Pr(>|t|) 
#cyl + hp = 0 -2.283815 0.7594086 -3.00736 0.005399071 

---

Remarque

LinearCombTest est mis à jour à Get p-value for group mean difference without refitting linear model with a new reference level, où l'on peut tester toute combinaison avec des coefficients combinaison alpha:

alpha[1] * vars[1] + alpha[2] * vars[2] + ... + alpha[k] * vars[k] 

plutôt que la somme

vars[1] + vars[2] + ... + vars[k]