My data consist of two columns - time and cumulative number like below:
time <- c(1:14)
cum.num <- c(20, 45, 99, 195, 301, 407, 501, 582, 679, 753, 790, 861, 1011, 1441)
My non linear function is:
B/(B*C*exp(-A*B*time) + 1)
My objective is to model my data using non linear method. I tried the following:
m1 < -nls(cum.num ~ B/((B*C)*exp(-A*B*time) + 1)
I have tried several initial values but got the following error:
Error in nls(cum.vul ~ B/((B * C) * exp(-A * B * time) + 1),
start = list(B = 60, : singular gradient
When this happens you generally have to do a bit of detective work to understand the parameters of your function and get a rough estimate of the values by looking at plots.
time <- c(1:14)
cum.num <- c(20, 45, 99, 195, 301, 407, 501, 582, 679,
753, 790, 861, 1011, 1441)
1/(1+x)
). It will be easier to visualize the data and estimate parameters if we invert the y
value, so that we have 1/cum.num ~ C*exp(-A*B*time) + 1/B
.plot(time,1/cum.num,log="y",xlim=c(0,14),ylim=c(0.0005,0.5))
(plotting on a log-y scale, and extending the y-axis limits, based on what I discovered below ...)
1/B ~ 0.001
or B ~ 1000
.C + 1/B = C + 0.001
. Looking at the graph, we have C ~ 0.5
1/(A*B)
is the characteristic scale of decrease (i.e. the time for an e-fold decrease). It looks like the e-folding time ~ 1 (from t=1 to t=2) so 1/(A*B) ~ 1
so A ~ 0.001
Use these starting values to fit:
m1 <- nls(cum.num ~ B/((B*C)*exp(-A*B*time) + 1),
start=list(A=0.001,B=1000,C=0.5))
Seems to work. Plot the predicted values:
tpred <- seq(0,14,length=101)
cpred <- predict(m1,newdata=data.frame(time=tpred))
par(las=1,bty="l")
plot(time,cum.num)
lines(tpred,cpred)