To make your example reproducible you have to provide the data somehow; I am pretty sure nprint doesn't effect the solution, but if it does this would be a bug and I would appreciate a reproducible report.
The example in nls.lm is a little complicated in order to show how to use an analytical expression for the gradient (possible to provide in the argument jac); since you don't need this, note that your residual function can be simplified into something along the lines of p <- list(a=-0.003, b=0.13, c=0.50, E=400) optim.f <- function(p, ST, SM, R, q) { res <- R -(p$a*SM^I(2)+p$b*SM+p$c)*exp(p$E*((1/(283.15-227.13))-(1/(ST+273.15-227.13)) abserr <- abs(res) qnum <- quantile(abserr, probs=q, na.rm=T) res[abserr > qnum] <- 0 res } Rmodel <- nls.lm(par = p, fn = optim.f, ST = ST, SM = SM, R = R, q = q) On Thu, 8 May 2008, Fernando Moyano wrote: > I solved the problem arising from using certain quantile values simply > by printing the iterations with the argument nprint. This gave me > correct estimates. I have no idea why. > > ----- Original Message ---- > From: elnano <[EMAIL PROTECTED]> > To: r-help@r-project.org > Sent: Thursday, May 8, 2008 5:43:31 PM > Subject: Re: [R] function in nls argument > > > I've basically solved the problem using the nls.lm function from the > minpack.lm (thanks Katharine) with some modifications for ignoring residuals > above a given percentile. This is to avoid the strong influence of points > which push my modeled vs. measured values away from the 1:1 line. > I based it on the example given for nls.lm. Here it is: > > R # soil respiration data > ST <- ST [!is.na(R)] # soil temeprature data. Had to remove na to make > nls.lm work > SM <- SM [!is.na(R)] # soil moisture data > R <- R [!is.na(R)] > q <- 0.95 # quantile > > p <- c("a"=-0.003, "b"=0.13, "c"=0.50, E=400) # model parameters with > estimated values > > Rf <- function(ST, SM, a, b, c, E) > { > expr <- > expression((a*SM^I(2)+b*SM+c)*exp(E*((1/(283.15-227.13))-(1/(ST+273.15-227.13))))) > eval(expr) > } > > optim.f <- function(p, ST, SM, R, Rfcall) > { > res <- R - do.call("Rfcall", c(list(ST = ST, SM = SM), as.list(p))) # > the nls.lm example divides this by sqrt(R), I don't know why. I removed > that. > abserr <- abs(res) > qnum <- quantile(abserr, probs=q, na.rm=T) # calculate the "q" > quantile of the absolute errors > res[abserr > qnum]=0 # convert > residuals above qnum to 0 > res > } > > Rmodel<- nls.lm(par = p, fn = optim.f, Rfcall = Rf, ST = ST, SM = SM, R > = R) > > summary(Rmodel) > > The only error I still get is when using lower q values. A q of around 0.95 > or less (depending on the dataset) gives me completely wrong parameter > estimates resulting in negative predicted values. Maybe someone has a > suggestion here. > > Fernando > > > > Katharine Mullen wrote: > > > > The error message means that the gradient (first derivative of residual > > vector with respect to the parameter vector) is not possible to work with; > > calling the function qr on the gradient multiplied by the square root of > > the weight vector .swts (in your case all 1's) fails. > > > > If you want concrete advice it would be helpful to provide the commented, > > minimal, self-contained, reproducible code that the posting guide > > requests. what are the values of ST04, SM08b, ch2no, and tower? > > > > General comments: If your goal is to minimize sum( (observed - > > predicted)^2), the function you give nls to minimize (optim.fun in your > > case) should return the vector (observed - predicted), not the scalar sum( > > (observed - predicted)^2). You may want to see the nls.lm function in > > package minpack.lm for another way to deal with the problem. > > > > On Wed, 7 May 2008, Fernando Moyano wrote: > > > >> Greetings R users, maybe there is someone who can help > >> me with this problem: > >> > >> I define a function "optim.fun" and want as output the > >> sum of squared errors between predicted and measured > >> values, as follows: > >> > >> optim.fun <- function (ST04, SM08b, ch2no, a, b, d, E) > >> { > >> predR <- > >> (a*SM08b^I(2)+b*SM08b+d)*exp(E*((1/(283.15-227.13))-(1/(ST04+273.15-227.13)))) > >> abserr <- abs(ch2no-predR) > >> qnum <- quantile(abserr, probs=0.95, na.rm=T) > >> > >> is.na(abserr) <- (abserr > qnum) > >> errsq <- sum(abserr^2, na.rm=T) > >> errsq > >> } > >> > >> Then I want to optimize parameters a,b,d and E as to > >> minimize the function output with the following: > >> > >> optim.model<-nls(nulldat ~ optim.fun(ST04, SM08b, > >> ch2no, a, b, d, E), data=tower, > >> start=list(a=-0.003,b=0.13,d=0.50, E=400), na.action = > >> na.exclude ) > >> > >> were nulldat=0 > >> At this point I get the following error message: > >> > >> Error in qr.default(.swts * attr(rhs, "gradient")) : > >> NA/NaN/Inf in foreign function call (arg 1) > >> Warning messages: > >> 1: In if (na.rm) x <- x[!is.na(x)] else if > >> (any(is.na(x))) stop("missing values and NaN's not > >> allowed if 'na.rm' is FALSE") ... : > >> the condition has length > 1 and only the first > >> element will be used > >> (this warning is repeated 12 times) > >> > >> Question: what does the error mean? What am I doing > >> wrong? > >> Thanks a bunch. > >> Nano > >> Jen, Germany > >> Max Planck for Biogeochemistry > >> > >> > >> > >> > >> > >> ____________________________________________________________________________________ > >> > >> [[elided Yahoo spam]] > >> > >> ______________________________________________ > >> 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. > > > > > > -- > View this message in context: > http://www.nabble.com/function-in-nls-argument-tp17108100p17127514.html > Sent from the R help mailing list archive at Nabble.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. > > > > > ____________________________________________________________________________________ > > [[elided Yahoo spam]] > > ______________________________________________ > 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.