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
