Линейная модель в R не подходит должным образом

Question

Я знаю, что заголовок не указывает в точности, что я имею в виду, поэтому позвольте мне объяснить это здесь. Я работаю над набором данных, который состоит из урожайности пшеницы для определенного типа пшеницы (A, B, C, D). Теперь моя проблема при подборе линейной модели заключается в том, что я пытаюсь подобрать:

lm1 = lm(yield ~ type), при этом R фиксирует первый тип пшеницы (A) и отмечает его как глобальный перехват, а затем оценивает влияние все остальные виды по доходности. Я знаю, что могу подобрать такую ​​линейную модель: lm2 = lm(yield ~ 0 + type), которая даст мне оценки влияния каждого типа на доходность, однако на самом деле я хочу видеть своего рода комбинацию этих двух из них.

Есть ли возможность подобрать линейную модель в R st lm3 = lm(yield ~ GlobalIntercept + type), где GlobalIntercept будет представлять общий отрезок моей линейной модели, и тогда я смогу увидеть влияние каждого типа пшеницы на этот общий отрезок. Так же, как в первой модели, но на этот раз мы оценили влияние всех типов пшеницы (A, B, C, D) на общую урожайность.

Answer 1

Вопросы к SO должны включать минимальный воспроизводимый пример данных - см. Инструкции в верхней части страницы тегов r . Поскольку вопрос не включал это, мы предоставим его на этот раз с использованием встроенного набора данных InsectSprays, который поставляется с R.

Вот несколько подходов:

1) lm / contr.sum / dummy.coef Попробуйте использовать contr.sum контрасты суммы к нулю для фактора spray и посмотрите на фиктивные коэффициенты. Это расширит коэффициенты, чтобы включить все 6 уровней фактора spray в этом примере:

fm <- lm(count ~ spray, InsectSprays, contrasts = list(spray = contr.sum))
dummy.coef(fm)
## Full coefficients are 
##                                                                           
## (Intercept):          9.5                                                  
## spray:                  A         B         C         D         E         F
##                  5.000000  5.833333 -7.416667 -4.583333 -6.000000  7.166667

sum(dummy.coef(fm)$spray)  # check that coefs sum to zero
## [1] 0

2) tapply Если каждый уровень имеет одинаковое количество строк в данных например, как в случае с InsectSprays, где каждый уровень имеет 12 строк, тогда мы можем взять среднее значение для каждого уровня и затем вычесть пересечение (которое является общим средним). Это не работает, если набор данных несбалансирован, т. Е. Если разные уровни имеют разное количество строк. Обратите внимание, как приведенные ниже вычисления дают тот же результат, что и (1).

mean(InsectSprays$count)  # intercept
## [1] 9.5

with(InsectSprays, tapply(count, spray, mean) - mean(count))
##         A         B         C         D         E         F 
##  5.000000  5.833333 -7.416667 -4.583333 -6.000000  7.166667 

3) aov / model.tables Мы также можем использовать aov с model.tables следующим образом:

fm2 <- aov(count ~ spray, InsectSprays)
model.tables(fm2)
## Tables of effects
##
##  spray 
## spray
##      A      B      C      D      E      F 
##  5.000  5.833 -7.417 -4.583 -6.000  7.167 

model.tables(fm2, type = "means")
## Tables of means
## Grand mean
##    
## 9.5 
##
##  spray 
## spray
##      A      B      C      D      E      F 
## 14.500 15.333  2.083  4.917  3.500 16.667 

4) emmeans Мы можем использовать lm, за которым следует emmeans, например:

library(emmeans)

fm <- lm(count ~ spray, InsectSprays)
emmeans(fm, "spray")
##  spray emmean   SE df lower.CL upper.CL
##  A      14.50 1.13 66   12.240    16.76
##  B      15.33 1.13 66   13.073    17.59
##  C       2.08 1.13 66   -0.177     4.34
##  D       4.92 1.13 66    2.656     7.18
##  E       3.50 1.13 66    1.240     5.76
##  F      16.67 1.13 66   14.406    18.93
##
## Confidence level used: 0.95 

Answer 2

Согласно предоставленной вами информации, я могу сделать вывод, что вы моделируете доходность как линейную функцию типа, имеющего четыре категории. Вы ожидаете иметь перехват помимо коэффициентов каждого из типов. В этом нет смысла.

Вы прогнозируете доходность на основе номинальной переменной. Если вы хотите иметь регрессию с перехватом, вам необходимо иметь переменную-предиктор с источником. Свойство номинальной переменной состоит в том, что она не имеет происхождения. Происхождение означает, что нулевое значение для предиктора. У номинальной переменной не может быть происхождения. Другими словами, перехват (с непрерывной переменной-предиктором) означает значение зависимой переменной y, когда значение предиктора равно нулю (в вашем случае категория типа равна нулю, что практически невозможно). Вот почему ваша модель берет одну из категорий как ссылочную и вычисляет для нее точку пересечения. Изменения переменной y, когда категория отличается от эталонной категории, задаются коэффициентами.