Message-ID: <9a39963d-9d48-806c-81f6-5f6174da5c2f@mcmaster.ca>
Date: 2022-08-12T22:54:40Z
From: John Fox
Subject: [R-pkg-devel] Strange behaviour of function called from within a function
In-Reply-To: <4312_1660344060_27CMexOO029722_81078bd7-d89f-ab7e-b8ec-793744ef2971@gmail.com>
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 at 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/