[Rd] Bug in profile.nls with algorithm = "plinear"

[Benjamin Christoffersen]

Dear sirs

It seems like there is a bug in `profile.nls` with `algorithm =
"plinear"` when a matrix is supplied on the right hand side. Here is
the bug and a potential fix

#####
# example where profile.nls does not work with `plinear` but does with
# `default`
require(graphics)
set.seed(1)
DNase1 <- subset(DNase, Run == 1)
x <- rnorm(nrow(DNase1))

f1 <- nls(density ~ b1/(1 + exp((xmid - log(conc))/scal)) + b2 / x,
          data = DNase1, start = list(xmid = 0, scal = 1, b1 = 1, b2 = 1))
coef   (f1)
#R>      xmid      scal        b1        b2
#R> 1.5055308 1.0455951 2.3707962 0.0006887
confint(f1)
#R> Waiting for profiling to be done...
#R>          2.5%    97.5%
#R> xmid  1.323034 1.727130
#R> scal  0.974952 1.123036
#R> b1    2.193703 2.600172
#R> b2   -0.001597 0.002978

f2 <- nls(density ~ cbind(1/(1 + exp((xmid - log(conc))/scal)), x),
          data = DNase1, start = list(xmid = 0, scal = 1),
          algorithm = "plinear")
coef   (f2)
#R>    xmid     scal    .lin1   .lin.x
#R> 1.461636 1.028726 2.323707 0.008807
confint(f2) # this fails
#R> Waiting for profiling to be done...
#R> Error in attr(ans, "gradient")[c(TRUE, TRUE, TRUE, TRUE, TRUE,
TRUE, TRUE,  :
#R>  (subscript) logical subscript too long
traceback()
# [output output abbreviated]

#####
# It seems to fail due to the dimension in gradCall
environment(f1$m$setPars)$gradCall
#R> attr(ans, "gradient")[c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE,
#R> TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE),
#R>    c(TRUE, TRUE, TRUE, TRUE), drop = FALSE]
(ca <- environment(f2$m$setPars)$gradCall)
#R> attr(ans, "gradient")[c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE,
#R> TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE),
#R>    c(TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE,
#R>    TRUE, TRUE, TRUE, TRUE, TRUE, TRUE), c(TRUE, TRUE)]

f2$m$setVarying
# [output output abbreviated]
environment(f2$m$setVarying)$getRHS.varying
#R> function ()
#R> {
#R>    ans <- getRHS.noVarying()
#R>    attr(ans, "gradient") <- eval(gradCall)
#R>    ans
#R> }
#R> <bytecode: 0x000000000ce31768>
#R> <environment: 0x000000000cbdbad0>
dim(attr(environment(f2$m$setVarying)$getRHS.noVarying(), "gradient"))
#R> [1] 16  2  2
sapply(ca[3:5], length)
#R> [1] 16 16  2

#####
# It seems like the `nlsModel.plinear` in `stats/R/nls.R` should be
#
#   # ...
#   for(i in 1:marg + 1L)
#      gradSetArgs[[i]] <- rep_len(TRUE, dimGrad[i-1])
#   # ...
#   gradCall <-
#     switch(length(gradSetArgs) - 1L,
#            call("[", gradSetArgs[[1L]], drop = FALSE, gradSetArgs[[2L]]),
#            call("[", gradSetArgs[[1L]], drop = FALSE,
gradSetArgs[[2L]], gradSetArgs[[3L]]),
#            call("[", gradSetArgs[[1L]], drop = FALSE,
gradSetArgs[[2L]], gradSetArgs[[3L]],
#                 gradSetArgs[[4L]]),
#            call("[", gradSetArgs[[1L]], drop = FALSE,
gradSetArgs[[2L]], gradSetArgs[[3L]],
#                 gradSetArgs[[4L]], gradSetArgs[[5L]]))

# this seems to fix the problem.
source("R/nls_mod.R") # `nls_mod` is assigned in `R/nls_mod.R`
f3 <- nls_mod(density ~ cbind(1/(1 + exp((xmid - log(conc))/scal)), x),
              data = DNase1, start = list(xmid = 0, scal = 1),
              algorithm = "plinear")
confint(f3)
#R> Waiting for profiling to be done...
#R>       2.5% 97.5%
#R> xmid 1.3169 1.635
#R> scal 0.9667 1.096

#####
# The example on the help page still gives the same
ca <- quote(nls(density ~ 1/(1 + exp((xmid - log(conc))/scal)),
                data = DNase1,
                start = list(xmid = 0, scal = 1),
                algorithm = "plinear"))
confint(eval(ca))
#R> Waiting for profiling to be done...
#R>       2.5% 97.5%
#R> xmid 1.3225 1.677
#R> scal 0.9748 1.114

ca[[1]] <- quote(nls_mod)
confint(eval(ca))
#R> Waiting for profiling to be done...
#R>       2.5% 97.5%
#R> xmid 1.3225 1.677
#R> scal 0.9748 1.114

#####
# session info
sessionInfo()
#R> R version 3.5.0 (2018-04-23)
#R> Platform: x86_64-w64-mingw32/x64 (64-bit)
#R> Running under: Windows >= 8 x64 (build 9200)
#R>
#R> Matrix products: default
#R>
#R> locale:
#R> [1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United
States.1252    LC_MONETARY=English_United States.1252
#R> [4] LC_NUMERIC=C                           LC_TIME=English_United
States.1252
#R>
#R> attached base packages:
#R> [1] stats     graphics  grDevices utils     datasets  methods   base
#R>
#R> loaded via a namespace (and not attached):
#R> [1] MASS_7.3-49     compiler_3.5.0  tools_3.5.0     packrat_0.4.8-1