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@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.
______________________________________________ R-help@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.
