0

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

Community
  • 1
  • 1
Kajal
  • 581
  • 11
  • 24
  • Fix your model call to `expmodel = nls(formula = P ~ exp(beta * AO), data = df, start = list(beta = 0.01))`. You are misunderstanding how formulas for `nls` are specified. You could also consider creating a selfstarting model (search Google for `SSexp`). – Roland May 23 '16 at 11:13
  • Thanks, I did what you said and it removed the 'lines differing' error, but the resulting plot is a huge jumble of straight lines - not an exponential curve. I shall look into `SSexp` later if I can, but do you know how to fix this new issue? – Kajal May 23 '16 at 11:16
  • You need to order the y and x values for the plot so that the x values are increasing. – Roland May 23 '16 at 11:17
  • Another answer also said to use `sort` to remove the jumbled lines, but it's still jumbled in a different way. I'll upload screenshots. – Kajal May 23 '16 at 11:17
  • I've uploaded screenshots now. Could you take a look? – Kajal May 23 '16 at 11:19

1 Answers1

2
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.

Roland
  • 127,288
  • 10
  • 191
  • 288
  • My actual dataset is much bigger, so I changed the first part of your code to refer to my actual table. The rest remains the same. However I do not see any line on the plot generated - it's strange, I'll try making its colour different and seeing what happens. – Kajal May 23 '16 at 11:29
  • I changed it to red by adding `col="red"` in the curves line but still no line visible.. I have no idea what is going wrong here. I have only changed the dataset – Kajal May 23 '16 at 11:37
  • Have you adjusted the `from` and `to` values? – Roland May 23 '16 at 11:37
  • Oh, thanks, I hadn't - `from` remains the same but `to`'s limit should have been 155. Still no line, though. – Kajal May 23 '16 at 11:40
  • Possibly a problem with the scale of the y-axis? Leave out `add = TRUE` to see only the curve and set `ylim` in `plot` accordingly? – Roland May 23 '16 at 11:42
  • *facepalm* Sorry, it was the `from` and `to` values. In the huge dataset I had missed out some extremely small values. It's giving an exponential line now, (although it's obviously a terrible fit), but my problem has been solved so thank you so much! – Kajal May 23 '16 at 11:44