jpw
jpw

Reputation: 41

Curve-fitting with nls() in R

I have some data that I fit a curve to a certain formula.

To do this, I use the nls-function, like so:

fitmodel <- nls(y ~ a+b/(1+exp(-((x-c)/d))),
    data = combined,
    start=list(a=200,b=2000, c=80, d=10.99),      
    trace=TRUE)

This works, but gives the warnings

"1: In min(x) : no non-missing arguments to min; returning Inf
2: In max(x) : no non-missing arguments to max; returning -Inf"

When I draw a function with the determined parameters in my plot, it does not fit at all, like here:

jpw <- coef(fitmodel)

logfit <- function(x,jpw) {jpw[1] + jpw[2]/(1+exp(-((x-jpw[3])/jpw[4])))
}

points(logfit(x, jpw),type='l')

bad fit

My friend tried to fit this data in another program. It finds the same parameters, and the function the other program draws fits the data beautifully. Also, it is really easy to find parameters manually that make the curve follow the data well.

Where did I mess up? I am a beginner, so it might be something stupid.

Thank you in advance!

Edited: Data file

Upvotes: 4

Views: 4460

Answers (2)

jeremycg
jeremycg

Reputation: 24945

Your problem is your plotting, you are only giving one value to points, so it is using that as the y, and defaulting to one value is one unit on the x axis (If you look at your original plot, you can see it ends at 439, which is the number of points you have).

You can fix this by plotting with x too:

plot(combined$V1~combined$V3)
points(x,logfit(x,jpw), type = 'l')

enter image description here

Upvotes: 2

Marco Sandri
Marco Sandri

Reputation: 24262

It would be interesting to see your dataset.
Anyway, here is a working example. Hope this can help you.

logfit <- function(x,jpw) {
   jpw[1] + jpw[2]/(1+exp(-((x-jpw[3])/jpw[4])))
}

jpw <- c(-2,1,0,.5)
x <- runif(100, -3, 3)
y <- logfit(x, jpw)+rnorm(100, sd=0.01)
df <- data.frame(x,y)

curve(logfit(x,jpw),from=-3,to=3, ,type='l')
points(x,y)

enter image description here

fitmodel <- nls(y ~ a + b/(1+exp(-((x-c)/d))),
    data = df,
    start=list(a=1, b=2, c=1, d=1),      
    trace=TRUE)

fitmodel

The output is:

Nonlinear regression model
  model: y ~ a + b/(1 + exp(-((x - c)/d)))
   data: df
        a         b         c         d 
-1.999901  1.002425  0.006527  0.498689 
 residual sum-of-squares: 0.009408

Number of iterations to convergence: 6 
Achieved convergence tolerance: 1.732e-06

Here I use the data given by @jpw.

df <- dget(file="data.txt")
names(df) <- c("y","v2","x")    
fitmodel <- nls(y ~ a + b/(1+exp(-((x-c)/d))),
    data = df,
    start=list(a=200,b=2000, c=80, d=10.99),      
    trace=TRUE)
summary(fitmodel)

The estimated parameters are:

Formula: y ~ a + b/(1 + exp(-((x - c)/d)))

Parameters:
   Estimate Std. Error t value Pr(>|t|)    
a  231.6587     2.8498   81.29   <2e-16 ***
b 1893.0646     6.3528  297.99   <2e-16 ***
c  151.5405     0.2016  751.71   <2e-16 ***
d   17.2068     0.1779   96.72   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 37.19 on 473 degrees of freedom

Number of iterations to convergence: 10 
Achieved convergence tolerance: 3.9e-06

And now I plot the results.

plot(df$x, df$y)
jpw.est <- coef(fitmodel)
curve(logfit(x,jpw.est), from=0, to=300, col="red", lwd=2, add=T)

enter image description here

Upvotes: 1

Related Questions