Полиномиальные контрасты для регрессии

17

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

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

весряTезнак равно52,7870+14,2587Икс-0,9680Икс2-0,1554Икс3,

где Икс должен принимать значения , , или 4 в соответствии с различным уровнем переменной интервала.2 31234

Это верно? И если да, то какова была цель полиномиальных контрастов?

Пиппо
источник
7
Нет, эти коэффициенты предназначены для ортогональных полиномиальных терминов: вы написали модель для необработанных полиномиальных терминов. Замените Икс , Икс2 и Икс3 значениями L , Q и С соответственно (из справочной таблицы).
Scortchi - Восстановить Монику
1
Уважаемый @ Scortchi, спасибо за ваш ответ. Я думаю, чтобы понять, что вы имеете в виду, но тогда я не совсем понял, как работают эти ортогональные полиномиальные термины. : P
Пиппо
1
Что примечательно, то, что у вас есть, не совсем подходящая модель. Вы либо нуждаетесь в гигантской «шляпе» над записью (или E [write]), что означает прогнозируемое значение записи или ожидаемое значение записи; или вам нужно «+ e» в конце, чтобы указать остатки.
gung - Восстановить Монику
@ Scortchi Что такое или как найти «справочную таблицу»?
Антони Пареллада
2
@AntoniParellada: это таблица на странице, на которую ссылается ОП: ats.ucla.edu/stat/r/library/contrast_coding.htm#ORTHOGONAL . и получил contr.polyв R.
Scortchi - Восстановить Монику

Ответы:

29

Напомним (и в случае, если гиперссылки OP потерпят неудачу в будущем), мы смотрим на набор данных hsb2как таковой:

   id     female race ses schtyp prog read write math science socst
1  70        0    4   1      1    1   57    52   41      47    57
2 121        1    4   2      1    3   68    59   53      63    61
...
199 118      1    4   2      1    1   55    62   58      58    61
200 137      1    4   3      1    2   63    65   65      53    61

которые можно импортировать сюда .

Мы превращаем переменную readв и упорядоченную / порядковую переменную:

hsb2$readcat<-cut(hsb2$read, 4, ordered = TRUE)
(means = tapply(hsb2$write, hsb2$readcat, mean))
 (28,40]  (40,52]  (52,64]  (64,76] 
42.77273 49.97849 56.56364 61.83333 

Теперь мы все настроены на запуск обычного ANOVA - да, это R, и у нас в основном есть непрерывная зависимая переменная writeи пояснительная переменная с несколькими уровнями readcat. В R мы можем использоватьlm(write ~ readcat, hsb2)


1. Генерация контрастной матрицы:

Упорядоченная переменная имеет четыре разных уровня readcat, поэтому у нас будет контрастов.N-1знак равно3

table(hsb2$readcat)

(28,40] (40,52] (52,64] (64,76] 
     22      93      55      30 

Для начала давайте взглянем на деньги и посмотрим на встроенную функцию R:

contr.poly(4)
             .L   .Q         .C
[1,] -0.6708204  0.5 -0.2236068
[2,] -0.2236068 -0.5  0.6708204
[3,]  0.2236068 -0.5 -0.6708204
[4,]  0.6708204  0.5  0.2236068

Теперь рассмотрим, что происходило под капотом:

scores = 1:4  # 1 2 3 4 These are the four levels of the explanatory variable.
y = scores - mean(scores) # scores - 2.5

Yзнак равно[-1,5,-0,5,0,5,1,5]

seq_len (n) - 1знак равно[0,1,2,3]

n = 4; X <- outer(y, seq_len(n) - 1, "^") # n = 4 in this case

[1-1,52,25-3,3751-0,50,25-0,12510,50,250,12511,52,253,375]

Что там произошло? outer(a, b, "^")поднимает элементы aк элементам b, таким образом , что первые результаты столбцов из операций, , ( - 0,5 ) 0 , 0,5 0 и 1,5 0 ; второй столбец из ( - 1,5 ) : 1 ,(-1,5)0(-0,5)00,501,50(-1,5)1 , 0,5 1 и 1,5 1 ; третье из ( - 1,5 ) 2 = 2,25(-0,5)10,511,51(-1,5)2знак равно2,25, , 0,5 2 = 0,25 и 1,5 2 = 2,25 ; и четвертый, ( - 1,5 ) 3 = - 3,375 , ( - 0,5 ) 3 = - 0,125 , 0,5 3 = 0,125 и 1,5 3 = 3,375 .(-0,5)2знак равно0,250,52знак равно0,251,52знак равно2,25(-1,5)3знак равно-3,375(-0,5)3знак равно-0,1250,53знак равно0,1251.53=3.375

Далее мы делаем ортонормированное разложение этой матрицы и берем компактное представление Q ( ). Некоторые из внутренних функций функций, используемых в факторизации QR в R, используемых в этом посте, объясняются здесь далее .QRc_Q = qr(X)$qr

[202.500.52.23604.5840.50.447200.50.8940.92961.342]

... из которых мы сохраняем только диагональ ( z = c_Q * (row(c_Q) == col(c_Q))). Что лежит в диагонали: просто «нижние» записи части разложения Q R. Только? ну нет ... Оказывается, диагональ верхней треугольной матрицы содержит собственные значения матрицы!RQR

Далее мы вызываем следующую функцию: raw = qr.qy(qr(X), z)результат, который может быть воспроизведен «вручную» двумя операциями: 1. Превращение компактной формы , т. Е. В Q , преобразование, которое может быть достигнуто с помощью , и 2. Выполнение умножение матриц Q z , как в .Qqr(X)$qrQQ = qr.Q(qr(X))QzQ %*% z

Важно отметить, что умножение на собственные значения R не изменяет ортогональность векторов составляющих столбцов, но, учитывая, что абсолютное значение собственных значений появляется в порядке убывания сверху вниз слева направо, умножение Q z будет иметь тенденцию к уменьшению значения в полиномиальных столбцах высшего порядка:QRQz

Matrix of Eigenvalues of R
     [,1]      [,2] [,3]      [,4]
[1,]   -2  0.000000    0  0.000000
[2,]    0 -2.236068    0  0.000000
[3,]    0  0.000000    2  0.000000
[4,]    0  0.000000    0 -1.341641

Сравните значения в более поздних векторах столбцов (квадратичных и кубических) до и после операций факторизации и с незатронутыми первыми двумя столбцами.QR

Before QR factorization operations (orthogonal col. vec.)
     [,1] [,2] [,3]   [,4]
[1,]    1 -1.5 2.25 -3.375
[2,]    1 -0.5 0.25 -0.125
[3,]    1  0.5 0.25  0.125
[4,]    1  1.5 2.25  3.375


After QR operations (equally orthogonal col. vec.)
     [,1] [,2] [,3]   [,4]
[1,]    1 -1.5    1 -0.295
[2,]    1 -0.5   -1  0.885
[3,]    1  0.5   -1 -0.885
[4,]    1  1.5    1  0.295

Наконец, мы называем (Z <- sweep(raw, 2L, apply(raw, 2L, function(x) sqrt(sum(x^2))), "/", check.margin = FALSE))превращение матрицы rawв ортонормированные векторы:

Orthonormal vectors (orthonormal basis of R^4)
     [,1]       [,2] [,3]       [,4]
[1,]  0.5 -0.6708204  0.5 -0.2236068
[2,]  0.5 -0.2236068 -0.5  0.6708204
[3,]  0.5  0.2236068 -0.5 -0.6708204
[4,]  0.5  0.6708204  0.5  0.2236068

"/"col.xi2(i) apply(raw, 2, function(x)sqrt(sum(x^2)))2 2.236 2 1.341(ii)(i)

R4contr.poly(4)

[0.67082040.50.22360680.22360680.50.67082040.22360680.50.67082040.67082040.50.2236068]

(sum(Z[,3]^2))^(1/4) = 1z[,3]%*%z[,4] = 0scores - mean123


2. Какие контрасты (столбцы) вносят существенный вклад в объяснение различий между уровнями в объясняющей переменной?

Мы можем просто запустить ANOVA и посмотреть на резюме ...

summary(lm(write ~ readcat, hsb2))

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  52.7870     0.6339  83.268   <2e-16 ***
readcat.L    14.2587     1.4841   9.607   <2e-16 ***
readcat.Q    -0.9680     1.2679  -0.764    0.446    
readcat.C    -0.1554     1.0062  -0.154    0.877 

... чтобы увидеть линейное влияние readcatна write, так что исходные значения (в третьем фрагменте кода в начале сообщения) можно воспроизвести как:

coeff = coefficients(lm(write ~ readcat, hsb2))
C = contr.poly(4)
(recovered = c(coeff %*% c(1, C[1,]),
               coeff %*% c(1, C[2,]),
               coeff %*% c(1, C[3,]),
               coeff %*% c(1, C[4,])))
[1] 42.77273 49.97849 56.56364 61.83333

... или...

enter image description here

... или намного лучше ...

enter image description here

i=1tai=0a1,,at

enter image description here

X0,X1,.Xn

Графически это гораздо проще понять. Сравните фактические средние значения по группам в больших квадратных черных блоках с предсказанными значениями и посмотрите, почему оптимально приближение прямой линии с минимальным вкладом квадратичных и кубических полиномов (с кривыми, аппроксимированными только с лессом):

enter image description here

Если бы просто для эффекта коэффициенты ANOVA были такими же большими для линейного контраста для других приближений (квадратичного и кубического), то бессмысленный график, приведенный ниже, более четко отображал бы полиномиальные графики каждого «вклада»:

enter image description here

Код здесь .

Антони Пареллада
источник
+1 Вау. Может ли этот ответ (я до сих пор его не читал) рассматривать как ответ на мой старый, забытый вопрос тоже stats.stackexchange.com/q/63639/3277 ?
ttnphns
(+1) @ttnphns: Возможно, там будет еще лучше.
Scortchi - Восстановить Монику
Просто совет: вы можете прокомментировать меня там со ссылкой на здесь; или выпустите ответ там - который я, вероятно, приму.
ttnphns
1
@ttnphns и @Scortchi Спасибо! Я потратил довольно много времени, пытаясь понять смысл этих концепций, и не ожидал особой реакции. Так что это очень позитивный сюрприз. Я думаю, что есть некоторые морщины, которые нужно сгладить в отношении объяснения qr.qy()функции, но я определенно постараюсь выяснить, смогу ли я сказать что-то минимально связное по вашему вопросу, как только у меня будет время.
Антони Пареллада
1
@ Элвис Я попытался выбрать хорошее итоговое предложение и поместить его где-то в посте. Я думаю, что это хороший момент и требует хорошего математического объяснения, но на данном этапе может быть слишком много подробностей.
Антони Пареллада
5

Я буду использовать ваш пример, чтобы объяснить, как это работает. Использование полиномиальных контрастов с четырьмя группами дает следующее.

ЕвесряTе1знак равноμ-0,67L+0,5Q-0,22СЕвесряTе2знак равноμ-0,22L-0,5Q+0,67СЕвесряTе3знак равноμ+0,22L-0,5Q-0,67СЕвесряTе4знак равноμ+0,67L+0,5Q+0,22С

Где первое уравнение работает для группы с наименьшими оценками чтения, а четвертое - для группы с лучшими показателями чтения. мы можем сравнить эти уравнения с приведенным с помощью нормальной линейной регрессии (предположим,реadя это непрерывно)

ЕвесряTеязнак равноμ+реadяL+реadя2Q+реadя3С

Обычно вместо L,Q,С вам придется β1,β2,β3и написано на первой позиции. Но это письмо напоминает полиномиальные контрасты. Так что цифры передL,Q,С на самом деле вместо реadя,реadя2,реadя3, Вы можете увидеть эти коэффициенты раньшеL have linear trend, before Q quadratic and before C cubic.

Then R estimates parameters μ,L,Q,C and gives you

μ^=52.79,L^=14.26,Q^=0.97,C^=0.16
Where μ^=14i=14Ewritei and estimated coefficients μ^,L^,Q^,C^что-то вроде оценки при нормальной линейной регрессии. Таким образом, из выходных данных вы можете увидеть, значительно ли оценочные коэффициенты отличаются от нуля, поэтому вы можете ожидать какой-либо линейный, квадратичный или кубический тренд.

В этом примере значительно ненулевой L^, Таким образом, ваш вывод может быть следующим: мы видим, что лучший результат в письменной форме линейно зависит от оценки чтения, но нет значительного квадратичного или кубического эффекта.

Фимба
источник