[R] NLS question:Quadratic plus plateau fit

[Jose A. Hernandez]

Dear R colleagues:

Am trying to fit a simple NL model to determine Economical Optimum Nitrogen 
Rates.
The segmented (quadratic + plateau) model only works with some y's, in some 
cases I get a "singular gradient" error.
I'll appreciate any ideas in how to solve the singular gradient error.

Thanks,

Jose

# The following code works using yield2 in the nls model but not using yield.

# Economical constants, nitrogen price in dollars per lb and corn price in 
dollars per bushel
nprice <- 0.17
cprice <- 2.25
ratio <- nprice/cprice

# Example data,

nrate <- c(0,60,90,120,150,180)
yield2 <- c(161.7,187.1,188.5,196.6,196.0,196.5)
yield <- c(163.4,178.1,179.6,178.2,184.4,184.5)

data.1 <- data.frame(nrate = nrate, yield = yield)
plot(data.1)

qp.nls.fit<- nls(yield ~ (b0 + b1*nrate + b2*I(nrate^2))*(nrate <= x0)
                 + (b0 + b1*x0 + b2*I(x0^2))*(nrate > x0),
                 data=data.1,
                 start=list(b0=125, b1=0.5, b2=-0.001, x0=135),
                 trace=T)
qp.nls.fit
summary(qp.nls.fit)
c.qp.fit <- coefficients(qp.nls.fit)
attach(as.list(c.qp.fit))

yld.x0 <- b0 + b1*x0 + b2*x0^2
maxn <- -b1/(2*b2)
eonr <- (ratio - b1)/(2*b2)
eonr2 <- round(eonr, digits = 0)
yldmaxn <- b0 + b1*maxn + b2*maxn^2
yldeonr <- b0 + b1*eonr + b2*eonr^2
yldeonr2 <- round(yldeonr, digits = 0)

plot(yield ~ nrate,
     pch = 16,
     main = "Quadratic plus Plateau Model",
     xlab = expression(paste("Nitrogen rate [lbs ac" ^-1,"]")),
     ylab = expression(paste("Corn yield [bu ac"^-1,"]")))
curve(predict(qp.nls.fit, data.frame(nrate = x)), add = T)
text(100,175, paste("EONR=",eonr2))
text(100,173, paste("Yield at EONR=",yldeonr2))


-- 
Jose A. Hernandez
Ph.D. Candidate
Precision Agriculture Center

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