R не пропускает базовый уровень взаимодействия факторов с числовым значением, если основной эффект числового значения исходит от функции

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

Сгенерируйте (коррелированные) данные игрушек:

set.seed(123)
test<-data.frame(x=rnorm(100,1,.5),z=factor(sample(c('a','b','c'),100,T)))
test$y<-.3*test$x+0*(test$z=='a')-.07*(test$z=='b')-.15*(test$z=='c')+rnorm(100,0,.1)

Запустите линейную модель:

> lm(y ~ x + z, test)
Call:
lm(formula = y ~ x + z, data = test)

Coefficients:
(Intercept)            x           zb           zc  
    0.02453      0.27484     -0.08279     -0.12868

Выглядит неплохо. Первый факторный уровень «а» опущен, как и должно быть. Теперь включите взаимодействие между числовым значением x и фактором z:

> lm(y ~ x + z + z:x, test)
Call:
lm(formula = y ~ x + z + z:x, data = test)

Coefficients:
(Intercept)            x           zb           zc         x:zb         x:zc  
   0.037008     0.262650    -0.134938    -0.118896     0.049068    -0.009225 
        lm(y ~ poly(x,2) + z:x, test)

Все по-прежнему в порядке. Теперь используйте функцию «poly», чтобы добавить квадратичное преобразование x:

> lm(y ~ poly(x, 2) + z + z:x, test)

Call:
lm(formula = y ~ poly(x, 2) + z + z:x, data = test)

Coefficients:
(Intercept)  poly(x, 2)1  poly(x, 2)2           zb           zc         za:x         zb:x         zc:x  
    0.33928      1.23017     -0.18029     -0.15478     -0.15574     -0.02749      0.04165           NA

И вот оно. Вместо того, чтобы исключить первый уровень z 'a' из термина взаимодействия, он включается вместе с двумя другими уровнями. Теперь za:x имеет «псевдоним», потому что модель, конечно, была бы единственной с включением всех трех уровней факторов. Это плохо, потому что не работают такие функции, как 'vif' из пакета 'car':

> vif(lm(y ~ poly(x,2) + z + z:x, test))
Error in vif.lm(lm(y ~ poly(x, 2) + z + z:x, test)) : 
  there are aliased coefficients in the model

Я пробовал такие вещи, как y ~ poly(x,2) + z + z:poly(x,1) или y ~ poly(x,2) + z + relevel(z, ref='a'):x, но ничего не казалось работать. Это ошибка или кто-то может объяснить этот результат? Есть ли способ избежать этой проблемы и по-прежнему использовать функциональные возможности формулы так, как я предполагал? Спасибо.


r interaction formula regression
person Nima    schedule 10.05.2015    source источник

Ответы (1)


arrow_upward
1
arrow_downward

Поскольку формулы позволяют использовать любую функцию, R не может узнать, какие функции вернут значения, равные другим значениям, уже включенным в уравнение. Специального кодирования для poly() не существует.

Если вы хотите просто включить термин x и x^2, вы можете сделать

lm(formula = y ~ x + I(x^2) + z + z:x, data = test)

избегать использования poly() вместе. Вы просто должны быть более осторожными в построении формулы.

person MrFlick    schedule 11.05.2015
comment
Разве y ~ x + I(x^2) + z + z:x и y ~ x + I(x^2) + z + z:x не разные? - person Jaehyeon Kim; 11.05.2015
comment
Я почти уверен, что вы набрали одно и то же дважды, но если вы имеете в виду, что y~poly(x,2) и y~x+I(x^2) разные, то да, потому что poly() по умолчанию создает ортогональные значения, но это то же самое, что и y~poly(x,2,raw=TRUE) - person MrFlick; 11.05.2015
comment
@MrFlick: должен ли R знать, какие значения вернет функция? Почему включение или исключение x влияет на формат z:x? Почему R не ведет себя последовательно, т.е. е. используя первый уровень z в качестве базы? I(x^2) работает в этом случае, но poly() был просто простым примером других вещей, которые я хотел бы сделать, таких как bs()... - person Nima; 11.05.2015
comment
Это связано с идентифицируемостью модели и чрезмерной спецификацией при использовании кодирования опорного уровня. Когда в модели присутствуют и x, и x:z, R распознает избыточность и удаляет параметр, не поддающийся оценке. Когда вы скрываете x в вызове функции, R не может сделать эту очистку за вас. Если вы не разбираетесь в статическом моделировании, вы можете вместо этого задать этот вопрос на странице Cross Validated. - person MrFlick; 11.05.2015
comment
@MrFlick: я понимаю статистический фон и знаком с математическими ограничениями линейных моделей. Я просто надеялся, что есть способ сказать R, как выполнить эту «очистку», когда задействованы функции, а встроенный механизм не может удалить избыточные векторы. Спасибо, что помогли мне понять, что это не так. - person Nima; 11.05.2015