[R] Quadratic curve on a loop

[Jose A. Hernandez]

Dear R community, once again I request your generous help on an R issue 
I cannot seem to figure out.

I am trying to do some analysis for multiple sites, in my example I am 
using only two, but the same will done on many sites. It's just ONE 
line, so I hope somebody can give a me a hand by email.

Look at the loop to make the plot x ~ y. If the curve function is 
removed from the loop everything works fine, but when I tried to add the 
curve to the plot, I get the following error:

Error in xy.coords(x, y) : x and y lengths differ
In addition: Warning message:
'newdata' had 101 rows but variable(s) found have 6 rows

I know this error should be a hint about how to proceed but I am stuck 
figuring ways to make the curve to print in the plots.

Any hints on how to solve this issue would be greatly appreciated.

Sincerely,

Jose

PS.
 > version
          _
platform i386-pc-mingw32
arch     i386
os       mingw32
system   i386, mingw32
status
major    2
minor    0.1
year     2004
month    11
day      15
language R



###################################################################
# some data
site <- c(1,1,1,1,1,1,2,2,2,2,2,2)
nrate <- c(0,60,90,120,150,180,0,60,90,120,150,180)
yield <- 
c(161.7,187.1,188.5,196.6,196.0,196.5,160.7,186.1,189.5,194.6,195.0,198.5)
site <- as.factor(site)

# file name
fname <- levels(site)
i <- nchar(fname) == 1
fname[i] <- paste("0", fname[i], sep = "")
fname <- paste("site", fname, ".eps", sep = "")

###################################################################
###################################################################

for (i in 1:length(levels(site))) {
     j <- site == levels(site)[i]
     x <- nrate[j]
     y <- yield[j]
     out <- try(lm(y ~ x + I(x^2)), silent = TRUE)
     plot(y ~ x,
          pch = 16,
          xlab = expression(paste("Nitrogen rate [lbs ac" ^-1,"]")),
          ylab = expression(paste("Corn yield [bu ac"^-1,"]")))
          curve(predict(out, data.frame(nrate = x)), add = T)
          title(paste("Site", i))
          # dev.copy2eps(file = fname[i], horizontal = TRUE)
}
###################################################################
###################################################################



-- 
Jose A. Hernandez
Department of Soil, Water, and Climate
University of Minnesota
1991 Upper Buford Circle
St. Paul, MN 55108

Ph. (612) 625-0445, Fax. (612) 625-2208