2012-12-15 63 views
4

假設我想擬合具有二階(正交)多項式的線性迴歸模型,然後預測響應。以下是第一個模型(M1)多項式迴歸無意義預測

x=1:100 
y=-2+3*x-5*x^2+rnorm(100) 
m1=lm(y~poly(x,2)) 
prd.1=predict(m1,newdata=data.frame(x=105:110)) 

現在的代碼,讓我們嘗試相同的模式,但不是使用$聚(X,2)$,我將利用其列,如:

m2=lm(y~poly(x,2)[,1]+poly(x,2)[,2]) 
prd.2=predict(m2,newdata=data.frame(x=105:110)) 

我們來看看m1和m2的總結。

> summary(m1) 

Call: 
lm(formula = y ~ poly(x, 2)) 

Residuals: 
    Min  1Q Median  3Q  Max 
-2.50347 -0.48752 -0.07085 0.53624 2.96516 

Coefficients: 
       Estimate Std. Error t value Pr(>|t|)  
(Intercept) -1.677e+04 9.912e-02 -169168 <2e-16 *** 
poly(x, 2)1 -1.449e+05 9.912e-01 -146195 <2e-16 *** 
poly(x, 2)2 -3.726e+04 9.912e-01 -37588 <2e-16 *** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.9912 on 97 degrees of freedom 
Multiple R-squared:  1,  Adjusted R-squared:  1 
F-statistic: 1.139e+10 on 2 and 97 DF, p-value: < 2.2e-16 

> summary(m2) 

Call: 
lm(formula = y ~ poly(x, 2)[, 1] + poly(x, 2)[, 2]) 

Residuals: 
    Min  1Q Median  3Q  Max 
-2.50347 -0.48752 -0.07085 0.53624 2.96516 

Coefficients: 
        Estimate Std. Error t value Pr(>|t|)  
(Intercept)  -1.677e+04 9.912e-02 -169168 <2e-16 *** 
poly(x, 2)[, 1] -1.449e+05 9.912e-01 -146195 <2e-16 *** 
poly(x, 2)[, 2] -3.726e+04 9.912e-01 -37588 <2e-16 *** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.9912 on 97 degrees of freedom 
Multiple R-squared:  1,  Adjusted R-squared:  1 
F-statistic: 1.139e+10 on 2 and 97 DF, p-value: < 2.2e-16 

所以m1和m2基本相同。現在讓我們看看預測prd.1和prd.2

> prd.1 
     1   2   3   4   5   6 
-54811.60 -55863.58 -56925.56 -57997.54 -59079.52 -60171.50 

> prd.2 
     1   2   3   4   5   6 
    49505.92 39256.72 16812.28 -17827.42 -64662.35 -123692.53 

Q1:爲什麼prd.2與prd.1有顯着不同?

Q2:如何使用模型m2獲得prd.1?

+1

不是一個答案,但足夠高的值總是嚇到我了... –

+1

這根本不是問題。我們可以用$ y = -2 + 3 * x-5 * x^2 + x^5 + rnorm(100,15)$和R平方減少到95%來改變$ y $,但問題依然存在預測。 – 2012-12-15 20:07:19

+0

第一個模型的結果看起來像是一個病態的矩陣。預測只是從第一個模型估計的無意義係數開始。 –

回答

8

m1是正確的做法。 m2正在進入一個痛苦的整個世界......

m2做預測,模型需要知道它被安裝在一組正交基函數,以便它使用相同的基礎功能外推新數據值。比較:poly(1:10,2)[,2]poly(1:12,2)[,2] - 前十個值不一樣。如果你明確地使用poly(x,2)來匹配模型,那麼predict理解所有這些,並做正確的事情。

您需要做的是確保您的預測位置使用與首先用於創建模型相同的基函數集進行轉換。您可以使用predict.poly這個(注意我打電話給我的解釋變量x1x2使其容易搭配的名字了)的R平方(0.99東西)

px = poly(x,2) 
x1 = px[,1] 
x2 = px[,2] 

m3 = lm(y~x1+x2) 

newx = 90:110 
pnew = predict(px,newx) # px is the previous poly object, so this calls predict.poly 

prd.3 = predict(m3, newdata=data.frame(x1=pnew[,1],x2=pnew[,2])) 
+0

非常感謝您的回覆。這完全回答了我的問題。我問第二個問題的原因是:假設我們擬合一個5階多項式,如m4 = lm(y〜poly(x,5))。然後在擬合之後,我們想通過刪除兩個項來擬合一個新模型(m5):即2次和4次多項式。最後用這個最終模型做一些預測(m5)。我不能用m4來做到這一點。但是,這可以通過使用您提到的參數和模型m3來完成。 – Stat