匠人府分享 http://blog.sciencenet.cn/u/meiweipingg

博文

非线性回归(迭代法)及其两种拟合曲线:y=a+b*exp(c*x)

已有 4545 次阅读 2016-11-11 14:10 |个人分类:R语言|系统分类:科研笔记|关键词:在R中运行非线性回归(迭代法)并添加拟合曲线:y=a+b*exp(c*x)

参考文献:

1.  薛毅,陈立萍.统计建模与R软件,第六章回归分析. 清华大学出版社

2.  决定系数R^2与调整决定系数adjR^2的相互转换计算方法

###############################################################

# ara  #导入数据ara,包含amountaraSL.mm.两个指标

ara <- read.csv("~/ara.csv")

View(ara)

attach(ara)

plot(amountara~SL.mm.)

#---------------建立所需的拟合模型------------------

# y=a+b*exp(c*x)  


mod1<-lm(amountara~SL.mm.)  # 线性 y~a+b*x

mod.ara<-nls(amountara~a+b*exp(c*SL.mm.),start= list(a=0.1,b=0.1,c=0.1))

# 非线性1y~a+b*expc*x


summary(mod1)   #导出both residual standard error and R^2

#Residual standard error: 0.1518 on 110 degrees of freedom# Multiple R-squared:  0.8093, Adjusted R-squared:  0.8076

summary(mod.ara)   #导出 residual standard errorno R^2

# Residual standard error: 0.1234 on 109 degrees of freedom # Number of iterations to convergence: 12


AIC(mod1)

# [1] -100.4275

AIC(mod.ara)       # AIC值越小,拟合的模型就越好

# [1] -145.7978# 上述AIC结果显示模型 mod.ara更好


#------------------next作图,并添加拟合曲线-----------------


x.ara<-seq(min(SL.mm.),max(SL.mm.),len=length(SL.mm.))  #  therange of x-axis

# N.B. 此处SL.mm.不可以使用 ara$SL.mm.,否则后面拟合曲线会出错,下同

# 或者

x.ara2<-seq(3,9,len=112) # 直接定义用于拟合的x轴的取值范围,如3 and 9 for the range of x-axis


y1<-predict(mod1,data.frame(SL.mm.=x.ara))   # lm

y.ara<-predict(mod.ara,data.frame(SL.mm.=x.ara,data=ara))    # nls

# N.B. 上述SL.mm.不可以使用ara$SL.mm.,否则后面拟合曲线会出错,下同

plot(amountara~SL.mm.,main="ara")


lines(x.ara,y1,col="black")   #线性

lines(x.ara,y.ara,col="red")     #非线性


#---------------如果想看所有迭代次数拟合曲线,如下-----------

#y.ara<-predict(mod.ara,data.frame(SL.mm.=x.ara,data=ara))    # nls,只画出1条最优拟合曲线

y.ara.iterations<-predict(mod.ara,data.frame(SL.mm.,data=ara))# Number of iterations to convergence: 12

# 不同之处在于:未将SL.mm.定义给 x.ara

plot(amountara~SL.mm.,main="ara")

lines(SL.mm.,y.ara.iterations,col="red")     #非线性,绘制所有迭代拟合曲线

# iterations means 迭代

#-------------------所有迭代次数的拟合曲线,结束---------------

#mod.ara<-nls(amountara~a+b*exp(c*SL.mm.),start =list(a=0.1,b=0.1,c=0.1))

# summary(mod.ara) 拟合结果显示 # -0.3695  0.1232 0.3186 分别对应上行的参数abc


text(5,1.2,"y=-0.3695+0.1232*exp(0.3186*SL.mm.)")


#######--------next 手动计算nls非线性拟合的R^2---------#######

#-----------how to get R^2 byresidual standard error????

#-----------the residual standard error σˆ =sqrt(SSE/(n-k))

# k 为参数个数,即y=a+b*X中的ab

# SSE= (residual standard error^2)*(n-k)

# R^2= 1- SSE/[var(amountara)*(n-1)]

# R^2= 1- [(residual standarderror^2)*(n-k)]/[var(amountara)*(n-1)]


# var(amountara)   #计算y的方差,不是标准差

sum.ara<-summary(mod.ara)

# sum.ara$sigma  # sigmais the value of "residual standard error"

# sum.ara$df[2]  # 即为 n-k

# sum.ara$df[1]  # 即为 k,参数个数,mod2k=3,即abc

# length(ara$SL.mm.) # 即为 n,样本数即112


R.square= 1-((sum.ara$sigma^2)*(sum.ara$df[2]))/(var(amountara)*(length(ara$SL.mm.)-1))

R.square

# [1] 0.8750826   # 非线性拟合结果的R^2=0.875082比线性拟合结果R^2=0.8093


sum.ara

# over






http://blog.sciencenet.cn/blog-651374-1014133.html

上一篇:在R中运行metaMDS
下一篇:x~y1+y2+...的情况下,ggplot2作柱状图及误差棒

2 张海权 杨正瓴

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备14006957 )

GMT+8, 2018-9-19 03:18

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部