Search code examples
rnls

R nls singular gradient to find a sigmoid function

I am trying to find an optimal sigmoid functional that adjust to my data (at the end of this post). But there is an error: Error in nls(cumulativo ~ f(eixox, phi1, phi2, phi3), start = st, data = data, : singular gradient

Any suggestions ?

library("ggplot2")

data<-structure(list(cumulativo = c(2, 3, 17, 191, 819, 1699, 2679, 
                          3907, 5535, 7254, 9226, 11543, 13809, 15542, 16852, 17709, 18246, 
                          18661, 18976, 19256, 19412, 19539, 19639), eixox = 1994:2016), 
      class = "data.frame", row.names = c(NA, -23L))
plot(cumulativo~eixox, data=data)

st <- list(phi1=20000,phi2=-5,phi3=.0005)
f <- function(x,phi1,phi2,phi3) {phi1/(1 + exp(-phi3 * x - phi2))}
curvaS<-nls(cumulativo~f(eixox,phi1,phi2,phi3),start=st,data=data,trace=TRUE)

dvd.csv

Solution

  • You have an error in your function definition. It should be

    f <- function(x, phi1, phi2, phi3) {phi1/(1 + exp(-phi3 * (x - phi2)))}

    Where phi1 is the upper bound, phi2 is the midpoint of the sigmoid curve, and phi3 is the rate. Note the extra parentheses so that phi2 is subtracted from x and then multiplied by -phi3. Now pick reasonable starting values and run nls:

    st <- list(phi1=20000, phi2=2005, phi3=.5)
curvaS <- nls(cumulativo~f(eixox, phi1, phi2, phi3), start=st, data=data, trace=TRUE)
# 20466691 :  20000.0  2005.0     0.5
# 1334673 :  19669.7851882  2004.2327533     0.4406048
# 902806.6 :  19566.0810794  2004.1449741     0.4639131
# 901808.4 :  19578.7102128  2004.1498061     0.4637724
# 901808.4 :  19578.7652076  2004.1498401     0.4637683
curvaS
# Nonlinear regression model
#   model: cumulativo ~ f(eixox, phi1, phi2, phi3)
#    data: data
#       phi1       phi2       phi3 
# 19578.7652  2004.1498     0.4638 
#  residual sum-of-squares: 901808
# 
# Number of iterations to convergence: 4 
# Achieved convergence tolerance: 0.000003139

    Finally, plot:

    plot(cumulativo~eixox, data=data)
pred <- predict(curvaS)
lines(data$eixox, pred)

    Plot