Как заставить lm () выводить t-тест сварного шва в R

#r

#r

Вопрос:

Мне говорили, и я видел примеры, когда линейная модель и t-тест в основном являются одним и тем же тестом, только t-тест представляет собой специализированную линейную модель с фиктивно закодированными предикторами. Есть ли способ заставить вывод lm выводить те же значения t, значения p, доверительные интервалы и стандартную ошибку, что и обычная t.test функция в r, где значение по умолчанию для var.equal аргумента равно FALSE ?

Например, прямо сейчас выходные данные lm и t.test прямо сейчас отличаются

 data("mtcars")
#these outputs below give me different values 
summary(lm(mpg ~ am, mtcars))
t.test(mpg ~ am, mtcars)
  

Я хочу, чтобы у lm были те же значения, что и у функции t.test., которая является t-тестом Уэлча. Как бы я это сделал?

Комментарии:

1. stats.stackexchange.com/a/144480/20833 показано, как выполнить t-тест Уэлча с помощью lme4::lmer ; с stats::lm это невозможно.

Ответ №1:

Во-первых, существует отличный пост о перекрестной проверке Как регрессия, t-тест и ANOVA являются версиями общей линейной модели? это дает много справочной информации о взаимосвязи между t-тестом, линейной регрессией и ANOVA.

По сути, p-значение из t-теста соответствует p-значению параметра наклона в линейной модели.

В вашем случае вам нужно сравнить

 t.test(mpg ~ am, mtcars, alternative = "two.sided", var.equal = T)
#
#   Two Sample t-test
#
#data:  mpg by am
#t = -4.1061, df = 30, p-value = 0.000285
#alternative hypothesis: true difference in means is not equal to 0
#95 percent confidence interval:
# -10.84837  -3.64151
#sample estimates:
#mean in group 0 mean in group 1
#       17.14737        24.39231

fit <- lm(mpg ~ as.factor(am), mtcars)
summary(fit)
#
#Call:
#lm(formula = mpg ~ as.factor(am), data = mtcars)
#
#Residuals:
#    Min      1Q  Median      3Q     Max
#-9.3923 -3.0923 -0.2974  3.2439  9.5077
#
#Coefficients:
#               Estimate Std. Error t value Pr(>|t|)
#(Intercept)      17.147      1.125  15.247 1.13e-15 ***
#as.factor(am)1    7.245      1.764   4.106 0.000285 ***
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
#Residual standard error: 4.902 on 30 degrees of freedom
#Multiple R-squared:  0.3598,   Adjusted R-squared:  0.3385
#F-statistic: 16.86 on 1 and 30 DF,  p-value: 0.000285
  

Обратите внимание, что p-значения совпадают.

Два комментария:

  1. as.factor(am) превращается am в категориальную переменную
  2. Чтобы соответствовать предположениям линейной модели (где термин ошибки epsilon ~ N(0, sigma^2) ), нам нужно использовать t.test with var.equal = T , который предполагает, что дисперсия одинакова для измерений из обеих групп.
  3. Разница в знаке значения t обусловлена другим определением опорного уровня «категоризированный» am .

Чтобы получить те же групповые средние значения в линейной модели, мы можем удалить перехват

 lm(mpg ~ as.factor(am) - 1, mtcars)
#
#Call:
#lm(formula = mpg ~ as.factor(am) - 1, data = mtcars)
#
#Coefficients:
#as.factor(am)0  as.factor(am)1
#         17.15           24.39
  

Комментарии:

1. Вместо того, чтобы показывать, как выполнить t-тест Уэлча с помощью lm этого ответа, показано, как заставить t.test вести себя подобно тесту в summary(lm()) (а именно, t-тесту Стьюдента).

Ответ №2:

Предположение линейной регрессии заключается в том, что остатки обычно распределяются со средним значением 0 и постоянной дисперсией. Следовательно, ваш t.test и сводка регрессии будут иметь согласованные результаты, только если вы предполагаете, что отклонения равны.