2013-02-23 39 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表示,总是绘制数据。