Search code examples
rnls

Fit an exponential curve using nls with a custom data frame in R

I have been through many questions on SO but I am not able to find a fix for my problem in related answers. I have a table like this stored as PAO1.data:

Name         P          AO
Prog1        0.654      59.702
Prog2        0.149      49.595
Prog3        0.505      50.538
Prog4        0.777      59.954
Prog5        0.237      49.611
Prog6        0.756      50.630
Prog7        0.560      118.014
Prog8        0.015      53.779
Prog9        0.789      68.096
Prog10       0.825      79.558

I tried to use nls to fit an exponential curve to the data.

df = data.frame(PAO1.data)
p = df$P
ao = df$AO
RMSE <- function(fit, act){
    sqrt(mean((fit - act)^2))
}

expmodel = nls(formula = p ~ exp(ao), data = df, start = list(ao = 0.01))
fit1 = fitted.values(expmodel)
err1 = RMSE(fit1, p)
plot(ao, p)
lines(ao, predict(expmodel))
print(err1)

When I try to run it, I get these warning messages while creating expmodel:

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

At the same time, I get this error at lines():

Error in xy.coords(x, y) : 'x' and 'y' lengths differ
Calls: lines -> lines.default -> plot.xy -> xy.coords

Another question I read on SO with the 'lengths differing' error actually had differing lengths in x and y. However, my x and y (here ao and p) have exactly the same number of values.

Please note, it is not likely that an exponential curve will actually be a good fit, but I am trying out several different models and I want to know how to use nls properly so I can do the same with other models.

Some related curve-fitting questions:

This question says that the starting data was the key. The smallest value for AO in my table is 0.015 and I chose 0.01, which is close enough in my opinion. This question asks about nls and the answer is given using a polynomial using lm. I specifically need to know how to use nls for a lot of complex models in the future and this will not work for me. This question looked promising but I am not able to find the problem in my code by looking at that question and answer - I have similar statements in my code too.

How do I do it?

Edit:

Here are screenshots of the solutions posted in the comments by Roland: (The actual dataset is bigger)

After changing the call to nls to expmodel = nls(formula = p ~ exp(beta * ao), data = df, start = list(beta = 0.01))

After changing the call to nls

After sorting the AO values with lines(sort(ao), predict(expmodel))

After sorting the AO values

Solution 
  • df <- read.table(text = "Name         P          AO
Prog1        0.654      59.702
                 Prog2        0.149      49.595
                 Prog3        0.505      50.538
                 Prog4        0.777      59.954
                 Prog5        0.237      49.611
                 Prog6        0.756      50.630
                 Prog7        0.560      118.014
                 Prog8        0.015      53.779
                 Prog9        0.789      68.096
                 Prog10       0.825      79.558", header = TRUE)

#use correct syntax:
expmodel = nls(formula = P ~ exp(beta * AO), data = df, start = list(beta = 0.01))

plot(P ~ AO, data = df)
#you could use lines after sorting, but this is more convenient:
curve(predict(expmodel, newdata = data.frame(AO = x)), from = 49, to = 120, add = TRUE)

    resulting plot

    Obviously this is not a good model for your data. As you know the exponential function goes through (0,1). You should consider adding an intercept.