[R] Complex nonlinear model

[Murray Jorgensen]

I am running 1.3.1 on a Windows (NT 4.0) machine.

I am trying to fit a nonlinear model intended to predict crop yield from
nutrient information.

I am troubled by an error message:

> library(nls)
> simparj.fm <- nls(Y ~ Y.model(gN, MnmN, OptN, DIs,  beta, eta1, eta2,
+       Popn, Dmax, AWC, SumEp, PotYield3, Nsupply),
+       start = simparj.st, trace =T)
Error in numericDeriv(form[[3]], names(ind), env) : 
        theta should be of type character

I suspect that I am doing something stupid. My code in full is

# parameter assignments
Nmin      <- 0.8852
Nopt1     <- 16.78
gN <- 0.5511
E.nfert1 <- 0.3271
E.nfert2 <- 0.6132
Beta <- 0.8902
Dls <- 0.5378
eta1 <- 0.3791
eta2 <- 0.6332
PopStd <- 90468
beta <- Beta
DIs <- Dls
MnmN <- Nmin
OptN <- Nopt1

# simulate experimental data for predictors
nsim <- 70
Popn <- rnorm(nsim,PopStd,0.1*PopStd)
Dmax <- rnorm(nsim,140.68,47.45)
AWC <- rnorm(nsim,186.86,47.41)
SumEp <- rnorm(nsim,318.54,32.53)
PotYield3 <- rnorm(nsim,0.16180,0.01167)
Nsoil <- rnorm(nsim,94.07,34.06)
Bdfield <- rnorm(nsim,1.0590,0.1420)
Bdlab <- rnorm(nsim,0.7876,0.1169)

Nfert.broad <- runif(nsim,95.3,576.5)
Nfert.band <- runif(nsim,122,250)
broad <- rbinom(nsim,1,0.2)
Nfert.broad <- Nfert.broad*broad
Nfert.band <- Nfert.band*(1 - broad)
Nsupply<-Nsoil*Bdfield/Bdlab + Nfert.broad*E.nfert1 + Nfert.band*E.nfert2

# define model function

Y.model <- function(gN, MnmN, OptN, DIs,  beta, eta1, eta2,
      Popn, Dmax, AWC, SumEp, PotYield3, Nsupply)
      {
      Ymax<- 1-ifelse(Popn<=PopStd, eta1, eta2)*log(Popn/PopStd)
      Ymax <- Ymax*PotYield3*Popn/1000
      Ymax <- Ymax*ifelse(Dmax<=DIs*AWC, 1, 1 - beta*(Dmax -DIs*AWC)/SumEp)
      Nstar <- (Nsupply- MnmN*Ymax) / (OptN*Ymax - MnmN*Ymax)
      Nstar<-pmax ( 0,Nstar)
      Ystar<-ifelse(Nstar<1, (1 + gN*(1 - Nstar))* Nstar^(1+gN ), 1)
      Ystar<-pmax ( 0, Ystar)
      Y.model <- Ystar*Ymax
      }

# generate response variable from model

Y <- Y.model(gN, MnmN, OptN, DIs,  beta, eta1, eta2,
      Popn, Dmax, AWC, SumEp, PotYield3, Nsupply)
Y <- Y + rnorm(nsim,0,1)

# attempt to fit starting from true parameters
library(nls)
simparj.st <- c(gN, MnmN, OptN, DIs,  beta, eta1, eta2)
simparj.fm <- nls(Y ~ Y.model(gN, MnmN, OptN, DIs,  beta, eta1, eta2,
      Popn, Dmax, AWC, SumEp, PotYield3, Nsupply),
      start = simparj.st, trace =T)




Dr Murray Jorgensen      http://www.stats.waikato.ac.nz/Staff/maj.html 
Department of Statistics, University of Waikato, Hamilton, New Zealand 
Email: maj at waikato.ac.nz                            Fax +64-7 838 4155
Phone +64-7 838 4773 home phone +64-7 856 6705  Mobile +64-21 139 5862

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._