Hi All.
I'd like to send some control arguments to the nls function when
performing a nlsList analysis.
I'm fitting a power model to some grouped data and would like to
impose lower bounds on the estimates using the "port" algorithm.
Obtaining the lower bound constraint works fine with a direct call to
nls for a single level of the grouping variable. ?However, the bounds
aren't imposed when performing the nlsList analysis across the levels
of the grouping variable. ?Any idea why this isn't working?
# Generate example data ##########
trial ?<- 1:100
result <- list()
for (i in 1:3) {
?min <- rnorm(max(trial),250,5)
?dif <- rnorm(max(trial),500,5)
?p ? <- rnorm(max(trial),-.5,.1)
?e ? <- rnorm(max(trial),mean=0,sd=5)
?y ? <- min + dif*(trial)^p + e
?result[[i]] <- data.frame(y,trial,id=i)
}
newdf ? <-do.call('rbind',result)
df.gr ? ? <- groupedData( y ~ trial | id, data=newdf)
### Single unit analysis
........................................................................
### The boundary condition on the dplt parameter is enforced! ..........
df.one <- subset(df.gr,id==1)
nls(y~SSpowrDplt(trial,min,dplt,dif,p),data=df.one,algorithm="port",lower=c(0.0,0.0,0.0,-10))
...... example output.......>Nonlinear regression model
> ?model: ?y ~ SSpowrDplt(trial, min, dplt, dif, p)
> ?data: ?df.one
> ? ?min ? ?dplt ? ? dif ? ? ? p
>247.052 ? 0.000 491.965 ?-0.462
> residual sum-of-squares: 234322
> Algorithm "port", convergence message: relative convergence (4)
##################################
### nlsList analysis
........................................................................................
### Boundary condition on dplt is not enforced
.............................................
Lfit.nls ? ?<-
nlsList(y~SSpowrDplt(trial,min,dplt,dif,p),data=df.gr,control=list(algorithm='port',lower=c(0.0,0.0,0.0,-10),maxiter=100))
...... example output.......>Call:
> ?Model: y ~ SSpowrDplt(trial, min, dplt, dif, p) | id
> ? Data: df.gr
>Coefficients:
> ? min
> ?Estimate Std. Error ? t value ? ? Pr(>|t|)
>1 276.2354 ? 16.16609 17.087337 1.086442e-44
>2 257.0127 ? 20.43564 12.576694 3.390686e-30
>3 206.4017 ? 29.01315 ?7.114075 7.354863e-12
> ? dplt
> ? ? Estimate Std. Error ? ?t value ?Pr(>|t|)
>1 -0.06579982 0.03848086 -1.7099365 0.0951222
>2 -0.01694362 0.04161933 -0.4071093 0.6786473
>3 ?0.08981518 0.04636532 ?1.9371199 0.0528957
> ? dif
> ?Estimate Std. Error ?t value ? ? Pr(>|t|)
>1 477.5049 ? 21.89002 21.81382 6.679439e-62
>2 488.7171 ? 22.11908 22.09482 4.466288e-66
>3 552.7105 ? 25.04206 22.07129 9.215344e-65
> ? p
> ? ?Estimate Std. Error ? t value ? ? Pr(>|t|)
>1 -0.5455936 0.06262040 -8.712713 7.615265e-16
>2 -0.4839114 0.06074282 -7.966560 1.307734e-14
>3 -0.4059903 0.05455864 -7.441355 9.297527e-13
>Residual standard error: 27.43384 on 888 degrees of freedom
#####################################################
If you look at the structure of Lfit.nls, it looks like the control
arguments are passed correctly.
str(Lfit.nls)
List of 3
?$ 1:List of 6
?..$ control ? ?:List of 7
?.. ..$ maxiter ?: num 100
?.. ..$ tol ? ? ?: num 1e-05
?.. ..$ minFactor: num 0.000977
?.. ..$ printEval: logi FALSE
?.. ..$ warnOnly : logi FALSE
?.. ..$ algorithm: chr "port"
?.. ..$ lower ? ?: num [1:4] 0 0 0 -10
If it helps, here's the selfStart function that I'm using....
powrDpltInit <-
function(mCall, LHS, data) {
? xy ? ? <- sortedXyData(mCall[["x"]],LHS,data)
? min.s ?<- min(y)
? dif.s ?<- max(y)-min(y)
? dplt.s <- 0.5
? p.s ? ?<- -.20
? value ?<- c(min.s, dplt.s, dif.s, p.s)
? names(value) <-
mCall[c("min","dplt","dif","p")]
? value
}
SSpowrDplt<-selfStart(~min + dplt*x + dif*x^p,initial=powrDpltInit,
parameters=c("min","dplt","dif","p"))
Thanks for your help!
Rick DeShon
