微信公众号搜"智元新知"关注
微信扫一扫可直接关注哦!

使用 R 中的逻辑函数拟合 sigmoid 曲线

如何解决使用 R 中的逻辑函数拟合 sigmoid 曲线

我有一个遵循 sigmoid 曲线的数据,我想拟合一个逻辑函数来为每个参与者提取三个(或两个)参数。我在网上找到了一些方法,但我不确定哪个是正确的选项。

This tutorial 说明您应该像这样使用 nls() 函数

fitmodel <- nls(y~a/(1 + exp(-b * (x-c))),start=list(a=1,b=.5,c=25))
 ## get the coefficients using the coef function 
params=coef(fitmodel)

...您显然需要起始值来找到最合适的值 (?)。

然后 this post 解释说要获得起始值,您可以使用“自启动模型可以为您估计好的起始值,因此您不必指定它们”:

fit <- nls(y ~ SSlogis(x,Asym,xmid,scal),data = data.frame(x,y))

但是在其他地方我也读到你应该使用 SSlogis 函数来拟合逻辑函数。请有人确认这两个步骤是否是最好的方法?还是应该使用从以前的类似数据中提取的值作为起始值?

另外,如果我根本不想让渐近线定义逻辑函数,我该怎么办?

谢谢!

解决方法

没有最好的方法,但 SSlogis 确实无需设置起始值,而如果您指定公式,您可以更好地控制参数化。

如果问题真的是如何将 a 固定在预定级别,这里是值 1,而不重写公式,然后在运行 a 之前设置 nls 并从开头省略它价值。

a <- 1
fo <- y ~ a / (1 + exp(-b * (x-c)))
nls(fo,start = list(b = 0.5,c = 25))

或者,这将 a=1 代入公式 fo 给出 fo2,而无需自己重写公式。

fo2 <- do.call("substitute",list(fo,list(a = 1)))
nls(fo2,c = 25))
,

作为@G。 Grothendieck 写道,没有通用的“最佳方式”,它总是取决于您的特定目标。使用 SSLogis 是个好主意,因为您不需要指定起始值,但定义自己的函数更灵活。请参见以下示例,其中我们使用启发式方法自行导出起始值,而不是手动指定它们。然后我们拟合了一个逻辑模型,作为一个小小的奖励,我们拟合了具有显式滞后阶段的 Baranyi 增长模型。

# time (t)
x <- c(0,2,4,6,8,10,12,14,16,18,20)

# Algae cell counts (Mio cells per ml)
y <- c(0.88,1.02,1.43,2.79,4.61,7.12,6.47,8.16,7.28,5.67,6.91)

## we now plot the data linearly and logarithmically
## the layout function is another way to subdivide the plotting area
nf <- layout(matrix(c(1,3,3),byrow = TRUE),respect = TRUE)
layout.show(nf) # this shows how the plotting area is subdivided

plot(x,y)
plot(x,log(y))

## we see that the first points show the steepest increase,## so we can estimate a start value of the growth rate
r <- (log(y[5]) - log(y[1])) / (x[5] - x[1])
abline(a=log(y[1]),b=r)

## this way,we have a heuristics for all start parameters:
## r:  steepest increase of y in log scale
## K:  maximum value
## N0: first value

## we can check this by plotting the function with the start values
f <- function(x,r,K,N0) {K /(1 + (K/N0 - 1) * exp(-r *x))}
plot(x,y,pch=16,xlab="time (days)",ylab="algae (Mio cells)")
lines(x,f(x,r=r,K=max(y),N0=y[1]),col="blue")

pstart <- c(r=r,N0=y[1])
aFit   <- nls(y ~ f(x,N0),start = pstart,trace=TRUE)

x1 <- seq(0,25,length = 100)
lines(x1,predict(aFit,data.frame(x = x1)),col = "red")
legend("topleft",legend = c("data","start parameters","fitted parameters"),col = c("black","blue","red"),lty = c(0,1,1),pch = c(16,NA,NA))

summary(aFit)
(Rsquared <- 1 - var(residuals(aFit))/var(y))

## =============================================================================
## Approach with Baranyi-Roberts model
## =============================================================================

## sometimes,a logistic is not good enough. In this case,use another growth
## model
baranyi <- function(x,N0,h0) {
  A <- x + 1/r * log(exp(-r * x) + exp(-h0) - exp(-r * x - h0))
  y <- exp(log(N0) + r * A - log(1 + (exp(r * A) - 1)/exp(log(K) - log(N0))))
  y
}

pstart <- c(r=0.5,K=7,N0=1,h0=2)
fit2   <- nls(y ~ baranyi(x,h0),trace=TRUE)

lines(x1,predict(fit2,col = "forestgreen",lwd=2)

legend("topleft","logistic model","Baranyi-Roberts model"),"red","forestgreen"),NA))

Sigmoidal curve

版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。