Полиномиальная регрессия с двумя переменными с R

Я пытаюсь сделать что-то довольно простое с R, но не уверен, что у меня это хорошо получается. У меня есть набор данных, содержащий три столбца V1, V4, V5, и я хочу выполнить регрессию, чтобы получить коэффициенты Ci, j следующего многочлена двух переменных:

sum[i=0->3] sum[j=0->i] Ci,j . (V4_k)^i . (V5_k)^(3-j)

Поэтому я попробовал использовать функцию polym:

lm(V1 ~ polym(V4, V5, degree=3, raw = TRUE), data)

что дает мне следующие коэффициенты

[1]  1.048122e+04 -2.050453e+02  1.407736e+00 -3.309312e-03 -3.748650e+01  8.983050e-01 -4.308559e-03  1.834724e-01 -6.868446e-04  4.030224e-04

Теперь, если я хорошо понимаю, как мы должны построить формулу, я предположил, что следующее даст то же самое:

lm(v1 ~ V4 + V5 + I(V4 * V5) + I(V4^2 * V5) + I(V4^3 * V5) + I(V4^2 * V5^2) + I(V4^2*V5^3) + I(V4^3 * V5^2) + I(V4^3 * V5^3), data)

Но получаются разные коэффициенты:

[1]  3.130403e+03 -1.652007e+01 -1.592879e+02  3.984177e+00 -2.419069e-02  3.919910e-05  1.008657e-04  4.271893e-07 -5.305623e-07 -2.289836e-09

Не могли бы вы рассказать мне, что я делаю не так и как правильно добиться этой регрессии с помощью R?


r regression polynomials
person Tanguy A.    schedule 20.08.2014    source источник
comment
Вы проверяли порядок возврата матрицы из poly ()? Похоже, что он возвращает V4, V4 ^ 2, V4 ^ 3, V5, V5 ^ 2, V5 ^ 3, ... что отличается от вашего пользовательского уравнения. Кроме того, V4 ^ 3 * V5 - полином 4-го порядка.   -  person won782    schedule 20.08.2014

Ответы (2)


arrow_upward
1
arrow_downward

С образцами данных

dd<-data.frame(x1=rnorm(50),
   x2=rnorm(50))
dd<-transform(dd, z = 2*x1-.5*x1*x2 + 3*x2^2+x1^2 + rnorm(50))

Мы видим, что

lm(z~polym(x1,x2,degree=3, raw=T), dd)
lm(z~x1+I(x1^2)+I(x1^3)+I(x2)+I(x1*x2) + 
   I(x1^2*x2)+I(x2^2) + I(x1*x2^2) + I(x2^3), dd)

одинаковы.

Обратите внимание, что в вашем расширении есть такие термины, как

I(V4^3 * V5) + I(V4^2 * V5^2)

которые оба являются членами 4-й степени (сумма показателей равна 4), поэтому они не должны появляться в полиноме третьей степени. Так что это зависит от того, чего вы хотите. Обычно для полинома третьей степени у вас есть

sum[i=0->3] sum[j=0->3-i] Ci,j . (V4_k)^i . (V5_k)^j

так i+j<=3 всегда. Мне неясно, какой именно тип регрессии вы хотите.

person MrFlick    schedule 20.08.2014
comment
В самом деле, мое расширение было неправильным. Но как теперь объяснить тот факт, что при первом вызове (с использованием polym) я получаю 10 коэффициентов, а не только 6, как я получаю при правильном раскрытии (при удалении экспонент ›3)? - person Tanguy A.; 20.08.2014
comment
@TanguyA. Если вы посмотрите на мой пример выше, который имеет правильное расширение, есть 9 членов плюс член перехвата, так что вы получите 10 оценок общих коэффициентов. Итак, 10 - правильное число. Похоже, в расширении вам не хватает некоторых терминов. - person MrFlick; 20.08.2014
comment
Ты прав. Теперь я получаю одинаковые результаты с обоими вызовами. Таким образом, я думаю, я получаю те коэффициенты, которые ищу. Спасибо за ответы! - person Tanguy A.; 20.08.2014

arrow_upward
2
arrow_downward

Вызов polym (V4, V5) не дает вам того, что вы думаете. (Неважно, используете ли вы в этом примере poly или polym)

Давайте посмотрим на пример:

v1 <- 1:10; v2 <- 1:10
poly(v1, v2, degree=3, raw=TRUE)
      1.0 2.0  3.0 0.1 1.1  2.1 0.2  1.2  0.3
 [1,]   1   1    1   1   1    1   1    1    1
 [2,]   2   4    8   2   4    8   4    8    8
 [3,]   3   9   27   3   9   27   9   27   27
 [4,]   4  16   64   4  16   64  16   64   64
 [5,]   5  25  125   5  25  125  25  125  125
 [6,]   6  36  216   6  36  216  36  216  216
 [7,]   7  49  343   7  49  343  49  343  343
 [8,]   8  64  512   8  64  512  64  512  512
 [9,]   9  81  729   9  81  729  81  729  729
[10,]  10 100 1000  10 100 1000 100 1000 1000

Метка столбца сообщает вам степень первого и второго векторов, которые вы указали в качестве аргументов. Первые три находятся от V2 ^ 0, вторые три линейны по V2 и так далее.

Это правильно, но в вашем втором примере есть термины 4-й степени. Если вы действительно ищете термины 4-й степени, просто измените степень на 4 в вызове метода.

Если вам нужна дополнительная помощь с полиномиальной регрессией, эта статья о R- Блогеры должны быть полезны. Он показывает, как создавать модели как с I(), так и с poly, хотя я думаю, что они были просто одномерными.

person MentatOfDune    schedule 20.08.2014