The optim() function has trouble calculating derivatives at/near the boundary, because it is using a simplistic finite-difference formula centered on the parameter. optimx::optimr() may work better.
-pd > On 26 Apr 2020, at 09:02 , Abby Spurdle <spurdl...@gmail.com> wrote: > > I ran your example. > It's possible that it's another bug in the optim function. > > Here's the optim call (from within fitdistr): > > stats::optim(x = c(1, 4, 1, 2, 3, 1, 1, 1, 2, 2, 2, 2, 1, 1, > 1, 1, 1, 4, 4, 3, 1, 2, 2, 1, 1, 3, 1, 1, 1, 4, 1, 1, 1, 1, 1, #more lines... > 1, 4, 1, 1, 1, 5, 5, 5, 4, 5, 2, 5, 5, 5, 5, 3, 3, 5, 4, 5, 2, #removed... > 4, 5, 5), topn = 5, lower = lwr, par = list(alpha = 1.010652, > beta = 1.929018), fn = function(parm, ...) -sum(log(dens(parm, ...))), > hessian = TRUE, method = "L-BFGS-B") > > And here's dens: > > function (parm, x, ...) > densfun(x, parm[1], parm[2], ...) > > I can't see any reason why it should call dens with parm[1] < lower[1]. > > On Sun, Apr 26, 2020 at 5:50 PM Abby Spurdle <spurdl...@gmail.com> wrote: >> >> I haven't run your example. >> I may try tomorrow-ish if no one else answers. >> >> But one question: Are you sure the "x" and "i" are correct in your function? >> It looks like a typo... >> >> >> On Sun, Apr 26, 2020 at 2:14 PM Rolf Turner <r.tur...@auckland.ac.nz> wrote: >>> >>> >>> For some reason fitdistr() does not seem to be passing on the "..." >>> argument "lower" to optim() in the proper manner, and as result >>> falls over. >>> >>> Here is my example; note that data are attached in the file "x.txt". >>> >>> dhse <- function(i,alpha,beta,topn) { >>> x <- seq(0,1,length=topn+2)[-c(1,topn+2)] >>> p <- dbeta(x,alpha,beta) >>> if(any(!is.finite(p))) browser() >>> (p/sum(p))[i] >>> } >>> >>> lwr <- rep(sqrt(.Machine$double.eps),2) >>> par0 <- c(alpha=1.010652,beta=1.929018) >>> x <- dget("x.txt") >>> fit <- MASS::fitdistr(x,densfun=dhse,topn=5,start=as.list(par0), >>> lower=lwr) >>> >>> The browser() in dhse() allows you to see that alpha has gone negative, >>> taking a value: >>> >>>> alpha >>>> -0.001999985 >>> >>> Continuing causes fitdistr() to fall over with the error message: >>> >>>> Error in stats::optim(x = c(1, 4, 1, 2, 3, 1, 1, 1, 2, 2, 2, 2, 1, 1, : >>>> non-finite finite-difference value [1] >>> >>> If I eschew using fitdistr() and "roll-my-own" as follows: >>> >>> foo <- function(par,x,topn){-sum(log(dhse(i=x,alpha=par[1], >>> beta=par[2], >>> topn=topn)))} >>> >>> fit <- optim(par0,fn=foo,method="L-BFGS-B",lower=lwr,topn=5,x=x) >>> >>> then optim() returns a result without complaint. >>> >>> Am I somehow messing up the syntax for fitdistr()? >>> >>> cheers, >>> >>> Rolf Turner >>> >>> P. S. I've tried supplying the "method" argument, method="L-BFGS-B" >>> explicitly to fitdistr(); doesn't seem to help. >>> >>> R.T. >>> >>> -- >>> Honorary Research Fellow >>> Department of Statistics >>> University of Auckland >>> Phone: +64-9-373-7599 ext. 88276 >>> ______________________________________________ >>> R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see >>> 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@r-project.org mailing list -- To UNSUBSCRIBE and more, see > 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. -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Office: A 4.23 Email: pd....@cbs.dk Priv: pda...@gmail.com ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see 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.