[R-pkg-devel] Strange behaviour of function called from within a function

John Fox j|ox @end|ng |rom mcm@@ter@c@
Sat Aug 13 00:54:40 CEST 2022


Dear John,

This is a scoping issue due to how model formulas are evaluated by nls() 
and other modeling functions. Try this:

tw2 <- function(formula, data, start, control, trace, weights) {
   firstcoef <- c(b1=199, b2=50, b3=0.3)
   cat("firstcoef:\n")
   print(firstcoef)
   cat("weights:"); print(weights)
   data$weights <- weights
   secondw<-nls(formula, data, firstcoef, control, algorithm=NULL, TRUE, 
weights=weights)
   secondw
}

 > tw2(hobbsu, weeddf, st2, control=list(), trace=TRUE, weights=wts)
firstcoef:
    b1    b2    b3
199.0  50.0   0.3
weights: [1] 0.5000000000 0.2500000000 0.1250000000 0.0625000000 
0.0312500000 0.0156250000
  [7] 0.0078125000 0.0039062500 0.0019531250 0.0009765625 0.0004882812 
0.0002441406
0.3342481   (3.98e+00): par = (199 50 0.3)
0.01987378  (1.14e-01): par = (199.2897 50.14155 0.3139337)
0.01961683  (2.88e-04): par = (199.8481 50.2582 0.3133964)
0.01961683  (1.97e-06): par = (199.8595 50.26061 0.3133932)
Nonlinear regression model
   model: weed ~ b1/(1 + b2 * exp(-b3 * tt))
    data: data
       b1       b2       b3
199.8595  50.2606   0.3134
  weighted residual sum-of-squares: 0.01962

Number of iterations to convergence: 3
Achieved convergence tolerance: 1.975e-06

(Your second attempt worked because you stuck the variable weights in 
the global environment.)

I hope that this does what you want.

John

On 2022-08-12 6:40 p.m., J C Nash wrote:
> I have been trying to sort out an odd issue for the last few days in 
> upgrading package
> nlsr (for nonlinear least squares). It generally does much better in 
> finding solutions from
> poor starts than nls(), but it does not generate the full nls returned 
> object. I've written a
> wrapper wrapnlsr() that calls my nlsr functions then uses the resulting 
> parameters in nls()
> to get the nls object. But up to now, neither weights nor subset are 
> allowed. I've not got
> to subset, but have been trying weights, with very strange results. The 
> full example is too
> much code for the list, so I've faked the nlsr output as "firstcoef" in 
> the example below.
> What is troubling is that a hand-run example "works", but using a 
> function to do it fails
> when there are weights.
> 
> The nls.R code referenced in the example around line 570 may suggest 
> that its handling the
> weights may be the culprit here. I don't understand why
>      eval(substitute(weights), data, environment(formula))
> is needed for a fixed numeric vector of weights (as required by the 
> manual).
> 
> Any help welcome.
> 
> John Nash
> 
> The example:
> 
> # wterr.R -- reproduce nls error when called from function
> rm(list=ls()) # probably not needed, but clear the workspace
> 
> tnow <- function(formula, data, start, control, trace) {
>    firstcoef <- c(b1=199, b2=50, b3=0.3) # as if from another solver
>    cat("firstcoef:\n")
>    print(firstcoef)
>    cat("No weights\n")
>    second<-nls(formula, data, firstcoef, control, algorithm=NULL, TRUE)
>    second
> }
> 
> tw <- function(formula, data, start, control, trace, weights) {
>    firstcoef <- c(b1=199, b2=50, b3=0.3)
>    cat("firstcoef:\n")
>    print(firstcoef)
>    cat("weights:"); print(weights)
>    secondw<-nls(formula, data, firstcoef, control, algorithm=NULL, TRUE, 
> weights=weights)
>    secondw
> }
> 
> weed <- c(5.308, 7.24, 9.638, 12.866, 17.069, 23.192, 31.443,
>            38.558, 50.156, 62.948, 75.995, 91.972)
> tt <- 1:12
> weeddf <- data.frame(tt, weed)
> hobbsu <- weed ~ b1/(1+b2*exp(-b3*tt))
> st2 <- c(b1=200, b2=50, b3=0.3)
> wts <- 0.5^tt # a straight scaling comes via wts <- rep(0.01, 12)
> # From this level, we can do the calculation
> firstcoef <- c(b1=199, b2=50, b3=0.3)
> 
> # Does NOT work from the function tw with weights, but unweighted OK.
> # Is the issue line 570/571 of ./R-4.2.1/src/library/stats/R/nls.R ??
> tn2<-tnow(hobbsu, weeddf, st2, control=list(), trace=TRUE)
> tn2
> tt2<-tw(hobbsu, weeddf, st2, control=list(), trace=TRUE, weights=wts)
> tt2
> 
> # But we can get it running directly
> formula<-hobbsu
> data<-weeddf
> control<-list()
> weights<-wts
> second<-nls(formula, data, firstcoef, control, algorithm=NULL, TRUE)
> second
> secondw<-nls(formula, data, firstcoef, control, algorithm=NULL, TRUE, 
> weights=weights)
> secondw
> # And if we rerun the functions AFTER succeeding here, tw now works
> tn2<-tnow(hobbsu, weeddf, st2, control=list(), trace=TRUE)
> tn2
> tt2<-tw(hobbsu, weeddf, st2, control=list(), trace=TRUE, weights=wts)
> tt2
> 
> ______________________________________________
> R-package-devel using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-package-devel
-- 
John Fox, Professor Emeritus
McMaster University
Hamilton, Ontario, Canada
web: https://socialsciences.mcmaster.ca/jfox/



More information about the R-package-devel mailing list