[R] Fix for nls bug???

Keith Jewell k.jewell at campden.co.uk
Tue Aug 5 18:19:05 CEST 2008


Hi All,

I've hit a problem using nls. I think it may be a restriction in the 
applicability of nls and I may have found a fix, but I've been wrong before.

This example is simplified to the essentials. My real application is much 
more complicated.

Take a function of matrix 'x' with additional arguments:
matrix 'aMat' whose values are _not_ to be determined by nls
vector 'Coeffs' whose vales _are_ to be determined.
For simplicity, this isn't a selfStart function with an 'initial' attribute, 
but that doesn't change things.

Myfunc<-function(x, aMat, Coeffs)
{
#
# result = quadratic response in x with
# terms selected by aMat
#
aMat[aMat!=0] <- Coeffs
rowSums((x%*%aMat)%*%t(x))
}

If aMat is passed in by name (e.g. aMat = bMat) nls fails.
e.g.
#
# data frame with some noise
DF <- data.frame(x1 = runif(20, 1, 20), x2=runif(20, 1, 20))
DF$y <- 1 +DF$x1 +DF$x2 +DF$x1*DF$x2 +DF$x1^2 + DF$x2^2 + rnorm(20)
#
# matrix to pass in as aMat
bMat <- matrix(c(1,1,0,0), 2, 2)
#
# and nls fails
nls(y ~ Myfunc(cbind(x1, x2), bMat, aVec), DF, start=list(aVec=c(1,2)))
#
# pass in the same matrix other than by name and it works
nls(y ~ Myfunc(cbind(x1, x2), matrix(c(1,1,0,0), 2, 2), aVec), DF, 
start=list(aVec=c(1,2)))

I think the problem lies in this line in nls

  for (var in varNames[!varIndex]) mf[[var]] <- eval(as.name(var), data, 
env)

This adds values for some named arguments (bMat) as columns of the data 
frame. The problem is that generally they don't have the same number of 
rows. I've made it work for my example by replacing that line with this 
line, which adds values for those arguments to the data frame as parameters 
rather than as a column

  attributes(mf)[["parameters"]] <- 
c(attributes(mf)[["parameters"]],lapply(varNames[!varIndex], function(var) 
eval(as.name(var), data, env)) )

Problem is, I really don't know nls internals enough to be sure I haven't 
broken something.
And anyway, if this is really an improvement I ought to share it, but don't 
know how.

Or I could have totally the wrong end of the stick...

Comments, corrections and advice are welcome.

Thanks in advance,

Keith Jewell
-----------------------
> sessionInfo()
R version 2.7.0 (2008-04-22)
i386-pc-mingw32

locale:
LC_COLLATE=English_United Kingdom.1252;LC_CTYPE=English_United 
Kingdom.1252;LC_MONETARY=English_United 
Kingdom.1252;LC_NUMERIC=C;LC_TIME=English_United Kingdom.1252

attached base packages:
[1] stats     graphics  grDevices datasets  tcltk     utils     methods 
base

other attached packages:
[1] xlsReadWrite_1.3.2 svSocket_0.9-5     svIO_0.9-5         R2HTML_1.58
[5] svMisc_0.9-5       svIDE_0.9-5

loaded via a namespace (and not attached):
[1] tools_2.7.0 VGAM_0.7-7



More information about the R-help mailing list