2013-02-23 38 views
3

我有兩個變量散點圖,比如這個:如何估計R中散點圖的最佳擬合函數?

x<-c(0.108,0.111,0.113,0.116,0.118,0.121,0.123,0.126,0.128,0.131,0.133,0.136) 

y<-c(-6.908,-6.620,-5.681,-5.165,-4.690,-4.646,-3.979,-3.755,-3.564,-3.558,-3.272,-3.073) 

,我想找到更符合這兩個變量之間的關係的功能。

準確地說我想比較三種型號的配件:linear,exponentiallogarithmic

我在考慮爲每個函數擬合我的值,計算每種情況下的可能性並比較AIC值。

但我真的不知道如何或從哪裏開始。任何可能的幫助將非常感謝。

非常感謝您提前。

Tina。

+0

您是否嘗試過使用'rgp'包進行符號迴歸?如果您包含一些示例數據,我們可以嘗試一下。更多細節在這裏:http://www.rsymbolic.org/projects/rgp/wiki/Symbolic_Regression – Ben 2013-02-23 16:07:58

+2

我們必須這麼基本嗎?你讀過數據了嗎?你有沒有做過任何探索性的地塊?你至少知道如何使用lm'包裝來配合線性模型?我們有點卡在水平沒有多一點... – Spacedman 2013-02-23 16:26:07

+0

非常感謝你,我已經添加了一個例子,我非常瞭解R中的基礎知識,但在擬合比迴歸更復雜的模型時我是新的。 – user18441 2013-02-23 17:00:06

回答

4

以下是比較五個模型的示例。由於前兩款車型的形式,我們可以使用lm來獲得良好的初始值。 (請注意,不應使用y的不同變換的模型進行比較,因此我們不應使用lm1lm2作爲比較模型,但僅用於起始值。)現在爲前兩個運行nls。在這兩個模型之後,我們在x中嘗試不同程度的多項式。幸運的是,lmnls使用一致的AIC定義(儘管它不一定是真的,其他R模型擬合函數具有一致的AIC定義),所以我們可以使用lm多項式。最後我們繪製前兩個模型的數據和擬合。

AIC越低越好,所以nls1最好跟着lm3.2後跟nls2

lm1 <- lm(1/y ~ x) 
nls1 <- nls(y ~ 1/(a + b*x), start = setNames(coef(lm1), c("a", "b"))) 
AIC(nls1) # -2.390924 

lm2 <- lm(1/y ~ log(x)) 
nls2 <- nls(y ~ 1/(a + b*log(x)), start = setNames(coef(lm2), c("a", "b"))) 
AIC(nls2) # -1.29101 

lm3.1 <- lm(y ~ x) 
AIC(lm3.1) # 13.43161 

lm3.2 <- lm(y ~ poly(x, 2)) 
AIC(lm3.2) # -1.525982 

lm3.3 <- lm(y ~ poly(x, 3)) 
AIC(lm3.3) # 0.1498972 

plot(y ~ x) 

lines(fitted(nls1) ~ x, lty = 1) # solid line 
lines(fitted(nls2) ~ x, lty = 2) # dashed line 

enter image description here

增加了幾個模型,隨後固定起來,並改變符號。此外,爲了跟進Ben Bolker的評論,我們可以用AICcmodavg軟件包中的AICc代替上述各處的AIC

+1

它可能是值得考慮的AICc這個小的數據集... – 2013-02-23 21:37:17

+0

非常感謝你! – user18441 2013-02-24 19:44:21

7

我會通過explantory地塊開始,像這樣:

x<-c(0.108,0.111,0.113,0.116,0.118,0.121,0.123,0.126,0.128,0.131,0.133,0.136) 
y<-c(-6.908,-6.620,-5.681,-5.165,-4.690,-4.646,-3.979,-3.755,-3.564,-3.558,-3.272,-3.073) 
dat <- data.frame(y=y,x=x) 
library(latticeExtra) 
library(grid) 
xyplot(y ~ x,data=dat,par.settings = ggplot2like(), 
     panel = function(x,y,...){ 
     panel.xyplot(x,y,...) 
     })+ 
    layer(panel.smoother(y ~ x, method = "lm"), style =1)+ ## linear 
    layer(panel.smoother(y ~ poly(x, 3), method = "lm"), style = 2)+ ## cubic 
    layer(panel.smoother(y ~ x, span = 0.9),style=3) + ### loeess 
    layer(panel.smoother(y ~ log(x), method = "lm"), style = 4) ## log 

enter image description here

看起來像你需要一個立方體的模型。

summary(lm(y~poly(x,3),data=dat)) 

Residual standard error: 0.1966 on 8 degrees of freedom 
Multiple R-squared: 0.9831, Adjusted R-squared: 0.9767 
F-statistic: 154.8 on 3 and 8 DF, p-value: 2.013e-07 
+0

+1這非常好,AIC值如何?在ggplot中探索平滑器的方法如下:http://www.ats.ucla.edu/stat/r/faq/smooths.htm – Ben 2013-02-23 18:33:29

+0

非常感謝你,我有安裝網格包的問題,​​我猜猜這是你的意思:http://www.stat.auckland.ac.nz/~paul/grid/grid.html(我有一個mac)。 – user18441 2013-02-23 18:44:43

+0

是的。保羅·穆雷爾的網格(保佑他)。無需安裝它,只需加載它,它就像在你給的鏈接中提到的那樣與R一起分發。 – agstudy 2013-02-23 18:48:03

0

您可以先閱讀Box和Cox關於轉換的經典論文。他們討論如何比較轉換和如何在一組或一系列潛在轉換中找到有意義的轉換。對數轉換和線性模型是Box-Cox系列的特例。

而@agstudy表示,總是繪製數據。