[R] nls() and lines()

Mon Jul 18 15:38:49 CEST 2011

Provided, of course, that I alter the lines for different data sets
and data frames, the code to plot a line derived from nls() onto a
plot works with no problems.

Here's an example:

Year NOP
2002   6
2003   8
2004  11
2005  19
2006  26
2007  25

mod1 <- nls(NOP~alpha*exp(beta*Year), data=aic,
start=list(alpha=1e-278, beta=0.3205), trace=T,
nls.control(maxiter=30000, minFactor=0.000005))
plot(NOP~Year, data=aic, pch=19, ylab="Number of papers")
mod1a=seq(2002, 2007, by=.0001)
lines(mod1a, predict(mod1, list(Year = mod1a)))

I've been using this code for several years not to get models from
nls() onto plot and I've never had an issue with it until the dataset
I referenced in my initial email.

Thanks for your assistance.

SR

Steven H. Ranney

http://www.steven-ranney.com
http://stevenranney.blogspot.com

On Mon, Jul 18, 2011 at 2:55 AM, Peter Ehlers <ehlers at ucalgary.ca> wrote:
> On 2011-07-17 17:37, Steven Ranney wrote:
>>
>> All -
>>
>> I'm having an issue with trying to plot a model derived from nls()
>> onto a simple plot.  I have included a sample data set and the code
>> that I've been using.
>>
>>    year month day       date location mileage  cost gallon      cpg
>>   mpg          x
>> 2009     1   4   1/4/2009      BZN  124585 19.39  14.37 1.349339
>> 10.71677 2009-01-04
>> 2009     1  15  1/15/2009      BZN  124888  23.2  16.12 1.439206
>> 18.79653 2009-01-15
>> 2009     1  27  1/27/2009      BZN  125133 21.51  14.35 1.498955
>> 17.07317 2009-01-27
>> 2009     2  16  2/16/2009      BZN  125429 27.96  15.54 1.799228
>> 19.04762 2009-02-16
>> 2009     2  27  2/27/2009      BZN  125702 26.82  14.27 1.879467
>> 19.13104 2009-02-27
>> 2009     3  19  3/19/2009      BZN  125941 24.38  12.91 1.888459
>> 18.51278 2009-03-19
>> 2009     4   9   4/9/2009      BZN  126260 32.59  16.30 1.999387
>> 19.57055 2009-04-09
>> 2009     4  28  4/28/2009      BZN  126587 34.58  16.79 2.059559
>> 19.47588 2009-04-28
>> 2009     5  17  5/17/2009      BZN  126925 35.78  16.57 2.159324
>> 20.39831 2009-05-17
>> 2009     5  27  5/27/2009      BZN  127240 35.57  15.01 2.369753
>> 20.98601 2009-05-27
>> 2009     6   7   6/7/2009      BZN  127590 40.99  16.60 2.469277
>> 21.08434 2009-06-07
>> 2009     6  21  6/21/2009      BZN  127910 41.52  15.64 2.654731
>> 20.46036 2009-06-21
>> 2009     7  21  7/21/2009      BZN  128264 43.37  16.67 2.601680
>> 21.23575 2009-07-21
>> 2009     8  11  8/11/2009      BZN  128618 42.68  16.42 2.599269
>> 21.55907 2009-08-11
>> 2009     8  27  8/27/2009      BZN  128947 43.12  16.60 2.597590
>> 19.81928 2009-08-27
>> 2009     9  21  9/21/2009      BZN  129255 40.44  15.56 2.598972
>> 19.79434 2009-09-21
>> 2009    10   1  10/1/2009      BZN  129541 38.55  14.83 2.599461
>> 19.28523 2009-10-01
>> 2009    10  11 10/11/2009      BZN  129806 36.65  14.10 2.599291
>> 18.79433 2009-10-11
>> 2009    10  22 10/22/2009      BZN  130027 30.18  11.61 2.599483
>> 19.03531 2009-10-22
>> 2009    11   9  11/9/2009      BZN  130358 43.19  16.62 2.598676
>> 19.91576 2009-11-09
>> 2009    11  22 11/22/2009      BZN  130631 39.23  15.09 2.599735
>> 18.09145 2009-11-22
>> 2009    12   5  12/5/2009      BZN  130950 44.43  17.09 2.599766
>> 18.66589 2009-12-05
>> 2009    12  30 12/30/2009      BZN  131239 42.14  16.70 2.523353
>> 17.30539 2009-12-30
>>
>> After converting my dates into R-usable dates:
>>
>> #convert my dates to R-usable dates
>> x<- strptime(date, format="%m/%d/%Y")
>> x
>> mileage<- cbind(mileage, x)
>>
>> I plot the data and model mpg as a function of date.  In the nls()
>> statement, I convert x back to a numeric value so that I can conduct
>> the regression:
>>
>> plot(mpg~x, data=mileage[year==2009,], ylab="Miles per gallon",
>> xlab="2009", yaxs="i", ylim=c(10,30))
>> nls.2009<-
>> nls(mpg~(alpha*(as.numeric(x)^2))+(bravo*as.numeric(x))+(charlie),
>> data=mileage[year==2009,], start=list(alpha=-2e-14, bravo=5e-5,
>> charlie=-31407),
>>   trace=T, na.action=na.omit,
>> nls.control(minFactor=0.000000000000000000001))
>> plot(mpg~x, data=mileage[year==2009,])
>>   modb=seq(min(as.numeric(x)), max(as.numeric(x)), by=10000)
>>   lines(modb, predict(nls.2009, lines(as.numeric(x)=modb)))
>>
>> Unfortunately, when I run the final line of this code, I get the
>> following:
>>
>> Error: unexpected '=' in "  lines(modb, predict(nls.2009,
>> lines(as.numeric(x)="
>>
>> In other similar analyses, I've been able to plot an nls() model using
>> this exact code--altered of course according to information--but here
>> I'm at a loss.  I'm certain it has something to do with the
>> lines(...as.numeric(x)) value I'm trying to plot, but I can't figure
>> out what I'm doing wrong.
>
> That last line of code doesn't look right to me. The arguments
> that you need to supply to predict() are 'object' and 'newdata',
> where 'newdata' must have the appropriate form. Unless you have
> your own function lines(), I don't think that lines(as.numeric(x)=modb)
> would qualify as newdata.
>
> It's usually a bad idea to shove too much stuff into a single command
> and a good idea to use str() often.
>
> This 'exact' code worked in the past?
>
> Peter Ehlers
>
>>
>> The model is fine, but it's the plotting of the model that escapes me.
>>
>> I'm running R version 2.12.1 on a Windows 7 machine.
>>
>> Thanks for your help -
>>
>> Steven H. Ranney
>>
>> http://stevenranney.blogspost.com
>> http://www.steven-ranney.com
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>
>