Сгенерируйте K-матрицу для двухфакторного дисперсионного анализа с взаимодействием анализа Тьюки с multcomp :: glht () - PullRequest
Новая
       104

Сгенерируйте K-матрицу для двухфакторного дисперсионного анализа с взаимодействием анализа Тьюки с multcomp :: glht ()

0 голосов
/ 28 мая 2020

Я пытаюсь понять, как реализовать одновременный вывод для обобщенных линейных моделей с помощью функции multcomp :: glht (). В частности, я хотел бы провести полный анализ Тьюки для примера двустороннего ANOVA из здесь . Они проводят postho c анализ Тьюки для модели с взаимодействием.

mod <- lm(breaks ~ wool * tension, data = warpbreaks)

  • шерсть - это двухуровневый фактор (A, B)
  • натяжение - трехуровневый фактор (L, M, H)

В примере виньетки они проверяют разницу в средних значениях уровней натяжения на каждом уровне шерсти. Однако мне интересно научиться искать различия между всеми возможными комбинациями уровней. Следуя примеру, с таким кодом: 

tmp <- expand.grid(tension = unique(warpbreaks$tension),
                    wool = unique(warpbreaks$wool))
X <- model.matrix(~ wool * tension, data = tmp)
glht(mod, linfct = X)
Tukey <- contrMat(table(warpbreaks$tension), "Tukey")
K1 <- cbind(Tukey, matrix(0, nrow = nrow(Tukey), ncol = ncol(Tukey)))
rownames(K1) <- paste(levels(warpbreaks$wool)[1], rownames(K1), sep = ":")
K2 <- cbind(matrix(0, nrow = nrow(Tukey), ncol = ncol(Tukey)), Tukey)
rownames(K2) <- paste(levels(warpbreaks$wool)[2], rownames(K2), sep = ":")
K <- rbind(K1, K2)
colnames(K) <- c(colnames(Tukey), colnames(Tukey))
summary(glht(mod, linfct = K %*% X))

они получат что-то вроде следующего: 

Simultaneous Tests for General Linear Hypotheses

Fit: lm(formula = breaks ~ wool * tension, data = warpbreaks)

Linear Hypotheses:
             Estimate Std. Error t value Pr(>|t|)   
A:M - L == 0 -20.5556     5.1573  -3.986  0.00131 **
A:H - L == 0 -20.0000     5.1573  -3.878  0.00187 **
A:H - M == 0   0.5556     5.1573   0.108  0.99996   
B:M - L == 0   0.5556     5.1573   0.108  0.99996   
B:H - L == 0  -9.4444     5.1573  -1.831  0.30795   
B:H - M == 0 -10.0000     5.1573  -1.939  0.25535   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)

Итак, как мне получить правильную матрицу контрастности, чтобы я делал сравнения вроде :

  • A: M - B: M == 0
  • A: M - B: H == 0

Я знаю, как контраст матрица K должна быть, чтобы я мог разобраться вручную. Однако это всего лишь пример для ознакомления с пакетом. Мой настоящий ANOVA имеет 5-уровневый фактор и еще 10-уровневый фактор, поэтому выполнение его вручную было бы такой болью. Спасибо

1 Ответ

0 голосов
/ 29 мая 2020

Я нашел решение с помощью пакета emmeans. Функция lsm, похоже, делает то, что я искал. В частности: 

library(emmeans)        # for lsm
model_glht <- glht(mod1, lsm(pairwise ~ wool : tension), by = NULL)
summary(model_glht)

Наконец-то я получаю нужные сравнения: 

     Simultaneous Tests for General Linear Hypotheses

Fit: lm(formula = breaks ~ wool * tension, data = warpbreaks)

Linear Hypotheses:
               Estimate Std. Error t value Pr(>|t|)    
A,L - B,L == 0  16.3333     6.4767   2.522   0.1285    
A,L - A,M == 0  20.5556     6.3052   3.260   0.0216 *  
A,L - B,M == 0  15.7778     6.4135   2.460   0.1465    
A,L - A,H == 0  20.0000     6.5399   3.058   0.0369 *  
A,L - B,H == 0  25.7778     5.8918   4.375   <0.001 ***
B,L - A,M == 0   4.2222     4.1239   1.024   0.9002    
B,L - B,M == 0  -0.5556     4.2877  -0.130   1.0000    
B,L - A,H == 0   3.6667     4.4746   0.819   0.9591    
B,L - B,H == 0   9.4444     3.4589   2.730   0.0809 .  
A,M - B,M == 0  -4.7778     4.0239  -1.187   0.8294    
A,M - A,H == 0  -0.5556     4.2225  -0.132   1.0000    
A,M - B,H == 0   5.2222     3.1261   1.671   0.5384    
B,M - A,H == 0   4.2222     4.3826   0.963   0.9210    
B,M - B,H == 0  10.0000     3.3391   2.995   0.0431 *  
A,H - B,H == 0   5.7778     3.5759   1.616   0.5738    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
(Adjusted p values reported -- single-step method)

Тем не менее, после стольких работ над виньетками мне все еще любопытно найти матрицу K с помощью базы или инструменты glht. Кто-нибудь знает как?

Спасибо

...