Qual método de comparação múltipla usar para um modelo mais antigo: lsmeans ou glht?

15

Estou analisando um conjunto de dados usando um modelo de efeitos mistos com um efeito fixo (condição) e dois efeitos aleatórios (participante devido ao design do sujeito e ao par). O modelo foi gerado com o lme4pacote: exp.model<-lmer(outcome~condition+(1|participant)+(1|pair),data=exp).

Em seguida, realizei um teste de razão de verossimilhança desse modelo em relação ao modelo sem o efeito fixo (condição) e tenho uma diferença significativa. Existem três condições no meu conjunto de dados, por isso quero fazer uma comparação múltipla, mas não tenho certeza de qual método usar . Encontrei várias perguntas semelhantes no CrossValidated e em outros fóruns, mas ainda estou bastante confuso.

Pelo que vi, as pessoas sugeriram usar

1. O lsmeanspacote - lsmeans(exp.model,pairwise~condition)que me fornece a seguinte saída:

condition     lsmean         SE    df  lower.CL  upper.CL
 Condition1 0.6538060 0.03272705 47.98 0.5880030 0.7196089
 Condition2 0.7027413 0.03272705 47.98 0.6369384 0.7685443
 Condition3 0.7580522 0.03272705 47.98 0.6922493 0.8238552

Confidence level used: 0.95 

$contrasts
 contrast                   estimate         SE    df t.ratio p.value
 Condition1 - Condition2 -0.04893538 0.03813262 62.07  -1.283  0.4099
 Condition1 - Condition3 -0.10424628 0.03813262 62.07  -2.734  0.0219
 Condition2 - Condition3 -0.05531090 0.03813262 62.07  -1.450  0.3217

P value adjustment: tukey method for comparing a family of 3 estimates 

2. O multcomppacote de duas maneiras diferentes - usando mcp glht(exp.model,mcp(condition="Tukey"))resultando em

     Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts


Fit: lmer(formula = outcome ~ condition + (1 | participant) + (1 | pair), 
    data = exp, REML = FALSE)

Linear Hypotheses:
                             Estimate Std. Error z value Pr(>|z|)  
Condition2 - Condition1 == 0  0.04894    0.03749   1.305    0.392  
Condition3 - Condition1 == 0  0.10425    0.03749   2.781    0.015 *
Condition3 - Condition2 == 0  0.05531    0.03749   1.475    0.303  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1
(Adjusted p values reported -- single-step method)

e usando lsm glht(exp.model,lsm(pairwise~condition))resultando em

Note: df set to 62

     Simultaneous Tests for General Linear Hypotheses

Fit: lmer(formula = outcome ~ condition + (1 | participant) + (1 | pair), 
    data = exp, REML = FALSE)

Linear Hypotheses:
                             Estimate Std. Error t value Pr(>|t|)  
Condition1 - Condition2 == 0 -0.04894    0.03749  -1.305   0.3977  
Condition1 - Condition3 == 0 -0.10425    0.03749  -2.781   0.0195 *
Condition2 - Condition3 == 0 -0.05531    0.03749  -1.475   0.3098  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1
(Adjusted p values reported -- single-step method)

Como você pode ver, os métodos fornecem resultados diferentes. Esta é minha primeira vez trabalhando com R e estatísticas, então algo pode estar errado, mas eu não saberia onde. Minhas perguntas são:

Quais são as diferenças entre os métodos apresentados? Li numa resposta a perguntas relacionadas que se trata dos graus de liberdade ( lsmeansvs. glht). Existem algumas regras ou recomendações quando usar qual, ou seja, o método 1 é bom para esse tipo de conjunto de dados / modelo etc.? Qual resultado devo relatar? Sem conhecer melhor, eu provavelmente reportaria o valor p mais alto que consegui jogar pelo seguro, mas seria bom ter um motivo melhor. obrigado

schvaba986
fonte

Respostas:

17

Não é uma resposta completa ...

A diferença entre glht(myfit, mcp(myfactor="Tukey"))e os dois outros métodos é que dessa maneira usa uma estatística "z" (distribuição normal), enquanto os outros usam uma estatística "t" (distribuição de Student). A estatística "z" é igual à estatística "t" com um grau infinito de liberdade. Este método é assintótico e fornece valores de p menores e intervalos de confiança mais curtos do que os outros. Os valores p podem ser muito pequenos e os intervalos de confiança podem ser muito curtos se o conjunto de dados for pequeno.

Quando executo, lsmeans(myfit, pairwise~myfactor)a seguinte mensagem aparece:

Loading required namespace: pbkrtest

Isso significa que lsmeans(para um lmermodelo) usa o pbkrtestpacote que implementa o método Kenward & Rogers para os graus de liberdade da estatística "t". Este método pretende fornecer melhores valores de p e intervalos de confiança do que o assintótico (não há diferença quando o grau de liberdade é grande).

Agora, sobre a diferença entre lsmeans(myfit, pairwise~myfactor)$contrastse glht(myfit, lsm(pairwise~factor), acabei de fazer alguns testes e minhas observações são as seguintes:

  • lsmé uma interface entre o lsmeanspacote e o multcomppacote (consulte ?lsm)

  • para um design equilibrado, não há diferença entre os resultados

  • para um projeto desequilibrado, observei pequenas diferenças entre os resultados (os erros padrão e a relação t)

Infelizmente, não sei qual é a causa dessas diferenças. Parece apenas lsmchamadas lsmeanspara obter a matriz de hipóteses lineares e os graus de liberdade, mas lsmeansusa uma maneira diferente de calcular os erros padrão.

Stéphane Laurent
fonte
Obrigado pela resposta detalhada! Perdi completamente a diferença na estatística do teste ... Você mencionou que os valores podem ser muito pequenos e os ICs muito estreitos para o método assintótico. Meu conjunto de dados consiste em ~ 30 participantes, então acho que vou me ater à estatística t. Quando você diz que o método Kenward & Rogers leva a melhores valores de p, você quer dizer mais preciso ou menor? Portanto, as diferenças são devidas a diferenças nos métodos de cálculo df e SE e não devido ao uso incorreto de um deles com o meu modelo, se eu entendi corretamente. Existe uma maneira de escolher o "melhor" método aqui?
schvaba986
11
(Eu sou o desenvolvedor do pacote lsmeans ) lsmeansusa o pacote pbkrtest, que fornece (1) cálculos de Kenward-Rogers df e (2) uma matriz de covariância ajustada com viés reduzido nas estimativas. Se você definir primeiro lsm.options(disable.pbkrtest=TRUE), a lsmeanschamada com adjust="mvt"produzirá os mesmos resultados que glht, exceto por pequenas diferenças devido ao algoritmo aleatório usado pelos dois pacotes para a distribuição t multivariada.
Russ Lenth
3
No entanto, sugiro o ajuste "mvt" sem desativar o teste de pbkr, devido ao ajuste de viés e ao fato de que sem df, os valores assintóticos (z) assumem essencialmente df infinito, produzindo assim valores de P irrealisticamente baixos.
Russ Lenth
3
A propósito, o summarymétodo glhtpermite vários métodos de teste de redução, além do ajuste de multiplicidade padrão de uma etapa (ICs simultâneos). Em um ponto completamente diferente, se você tiver mais de um fator, lsmpoderá criar os tipos usuais de comparação com bastante facilidade, embora mcpnão possa fazê-lo.
Russ Lenth