Please provide reproducible code which shows the error. > dd <- structure(list(Lum = c(0.15, 0.07, 0.1, 0.19, 0.4, 0.73, 1.2, + 1.85, 2.91, 3.74, 5.08, 6.43, 8.06, 9.84, 12, 14.2, 16.6, 0.1, + 0.1, 0.17, 0.46, 1.08, 2.22, 3.74, 5.79, 8.36, 11.6, 15.4, 19.9, + 24.6, 30.4, 36.1, 43, 49.9, 0.06, 0.06, 0.08, 0.13, 0.25, 0.43, + 0.66, 1.02, 1.46, 1.93, 2.49, 3.2, 3.96, 4.9, 5.68, 6.71, 7.93 + ), GL = as.integer(c(0, 15, 31, 47, 63, 79, 95, 111, 127, 143, + 159, 175, 191, 207, 223, 239, 255, 0, 15, 31, 47, 63, 79, 95, + 111, 127, 143, 159, 175, 191, 207, 223, 239, 255, 0, 15, 31, + 47, 63, 79, 95, 111, 127, 143, 159, 175, 191, 207, 223, 239, + 255)), Gun = structure(as.integer(c(2, 2, 2, 2, 2, 2, 2, 2, 2, + 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, + 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 + )), .Label = c("bleu", "rouge", "vert"), class = "factor")), .Names = c("Lum", + "GL", "Gun"), class = "data.frame", row.names = c("1", "2", "3", + "4", "5", "6", "7", "8", "9", "10", "11", "12", "13", "14", "15", + "16", "17", "18", "19", "20", "21", "22", "23", "24", "25", "26", + "27", "28", "29", "30", "31", "32", "33", "34", "35", "36", "37", + "38", "39", "40", "41", "42", "43", "44", "45", "46", "47", "48", + "49", "50", "51")) > > st <- list(Blev = -0.06551802, beta = c(1.509686e-05, 4.55525e-05, + 7.32272e-06), gamm = 2.51187) > > dd.nls <- nls(Lum ~ Blev + beta[Gun] * GL^gamm, + data = dd, start = st) > > confint(dd.nls) Waiting for profiling to be done... 2.5% 97.5% Blev -1.612492e-01 2.988388e-02 beta1 6.108282e-06 8.762679e-06 beta2 1.269000e-05 1.792914e-05 beta3 3.844042e-05 5.388546e-05 gamm 2.481102e+00 2.542966e+00 > R.version.string [1] "R version 2.4.0 Patched (2006-10-24 r39722)"
On 11/17/06, Ken Knoblauch <[EMAIL PROTECTED]> wrote: > Hi, > > I'm fitting a standard nonlinear model to the luminances measured > from the red, green and blue guns of a TV display, using nls. > > The call is: > > dd.nls <- nls(Lum ~ Blev + beta[Gun] * GL^gamm, > data = dd, start = st) > where st was initally estimated using optim() > > st > $Blev > [1] -0.06551802 > > $beta > [1] 1.509686e-05 4.555250e-05 7.322720e-06 > > $gamm > [1] 2.511870 > > This works fine but I received an error message when I tried to > use confint(). I thought that getting derivatives with deriv might > help but found that deriv does not automatically handle the > indexing of the beta parameter. I modified the output of deriv > from the case when the term beta is not indexed to give: > > dd.deriv2 <- function (Blev, beta, gamm, GL) > { > .expr1 <- GL^gamm > .value <- Blev + rep(beta, each = 17) * .expr1 > .grad <- array(0, c(length(.value), 5), list(NULL, c("Blev", > "beta.rouge", "beta.vert", "beta.bleu", "gamm" ))) > .grad[, "Blev"] <- 1 > .grad[1:17, "beta.rouge"] <- .expr1[1:17] > .grad[18:34, "beta.vert"] <- .expr1[1:17] > .grad[35:51, "beta.bleu"] <- .expr1[1:17] > .grad[, "gamm"] <- ifelse(GL == 0, 0, rep(beta, each = 17) * (.expr1 * > log(GL))) > attr(.value, "gradient") <- .grad > .value > } > > which is not general but > now I can get a result from confint. My question: Can deriv() > be made to handle an indexed term more automatically (elegantly)? > I think that this would become more urgent (or require more work > on my part) if I were to want the Hessian, too, for example, in > anticipation of using rms.curv as described in the on-line > complements for MASS. > > The plinear algorithm can be used for this model (it is similar in some > respects to the example given on p 218 of MASS, but the intercept > terms are indexed instead, there). > > Here are the data, if of interest: > dd > Lum GL Gun > 1 0.15 0 rouge > 2 0.07 15 rouge > 3 0.10 31 rouge > 4 0.19 47 rouge > 5 0.40 63 rouge > 6 0.73 79 rouge > 7 1.20 95 rouge > 8 1.85 111 rouge > 9 2.91 127 rouge > 10 3.74 143 rouge > 11 5.08 159 rouge > 12 6.43 175 rouge > 13 8.06 191 rouge > 14 9.84 207 rouge > 15 12.00 223 rouge > 16 14.20 239 rouge > 17 16.60 255 rouge > 18 0.10 0 vert > 19 0.10 15 vert > 20 0.17 31 vert > 21 0.46 47 vert > 22 1.08 63 vert > 23 2.22 79 vert > 24 3.74 95 vert > 25 5.79 111 vert > 26 8.36 127 vert > 27 11.60 143 vert > 28 15.40 159 vert > 29 19.90 175 vert > 30 24.60 191 vert > 31 30.40 207 vert > 32 36.10 223 vert > 33 43.00 239 vert > 34 49.90 255 vert > 35 0.06 0 bleu > 36 0.06 15 bleu > 37 0.08 31 bleu > 38 0.13 47 bleu > 39 0.25 63 bleu > 40 0.43 79 bleu > 41 0.66 95 bleu > 42 1.02 111 bleu > 43 1.46 127 bleu > 44 1.93 143 bleu > 45 2.49 159 bleu > 46 3.20 175 bleu > 47 3.96 191 bleu > 48 4.90 207 bleu > 49 5.68 223 bleu > 50 6.71 239 bleu > 51 7.93 255 bleu > > and here is the sessionInfo > R version 2.4.0 Patched (2006-11-10 r39843) > i386-apple-darwin8.8.1 > > locale: > C > > attached base packages: > [1] "stats" "graphics" "grDevices" "utils" "datasets" > [6] "methods" "base" > > other attached packages: > boot MASS lattice > "1.2-26" "7.2-29" "0.14-13" > > Thanks for any suggestions. > > best, > > Ken > > > -- > Ken Knoblauch > Inserm U371 > Institut Cellule Souche et Cerveau > Département Neurosciences Intégratives > 18 avenue du Doyen Lépine > 69500 Bron > France > tel: +33 (0)4 72 91 34 77 > fax: +33 (0)4 72 91 34 61 > portable: +33 (0)6 84 10 64 10 > http://www.lyon.inserm.fr/371/ > > ______________________________________________ > R-help@stat.math.ethz.ch mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html > and provide commented, minimal, self-contained, reproducible code. > ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.