非线性拟合
先画图，再观测以提供a和b的初始值
#根据均匀分布,正态分布生成随机数
x <- seq(from=0,to = 50,by = 1)
y <- runif(n=1,min = 5,max = 15)*exp(-runif(n = 1,min = 0.01,max = 0.05)*x)+rnorm(n = 51,mean = 0,sd = 0.5)
plot(x,y)
a_start <- 14
b_start <- -log(2/a_start)/40
拟合
m <- nls(y ~ a*exp(-b*x),start = list(a=a_start,b=b_start))
coef(m)#查看拟合的参数
cor(y,predict(m))#查看残差平方和
cor(y,predict(m))#查看相关系数
重新画图和添加标注及公式
plot(x,y)
lines(x,y = predict(m),col='blue')
legend('topright',
legend = c('Data','Non-linear fit'),
pch = c(1,NA),
lty = c(NA,1),
col = c('black','blue'),
bty = "n")
eq <- substitute(italic(y) == a*italic(e)^{-b*italic(x)}*","~italic(RSE)==cor,
list(a = format(unname(coef(m)[1]), digits = 3),
b = format(unname(coef(m)[2]), digits = 2),
cor = format(summary(m)$sigma, digits = 3)))
text(x = 20,y = 10,labels = eq,adj = 0,col="blue")
ggplot绘图
data1=data.frame(x,y)
data2=data.frame(x,predict(m))
ggplot()+
geom_point(aes(x,y),data=data1,size=2,col=brewer.pal(8,"Set3")[4])+
geom_line(aes(x,predict.m.),data=data2,col=brewer.pal(8,"Set3")[1],lty=1,size=1.2)+
geom_text(aes(x = 25, y = 10),label=as.character(as.expression(eq)),parse = T,color="black")
