R version 2.13.0 alpha (2011-03-18 r54865)
Copyright (C) 2011 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
Platform: i686-pc-linux-gnu (32-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

  Natural language support but running in an English locale

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

[Previously saved workspace restored]

> library(bbmle)
Loading required package: stats4
Loading required package: numDeriv
Loading required package: lattice
Loading required package: MASS
> 
> ## source("~/lib/R/pkgs/bbmle/pkg/R/mle.R
> 
> set.seed(1002)
> X <- rexp(1000, rate = 0.0001)
> f <- function(X, rate) {
+   if (rate<0) cat("rate<0: ",rate,"\n")
+   -sum(dexp(X, rate = rate, log = TRUE))
+ }
> if (FALSE) {
+   ## L-BFGS-B violates bounds, and gets stuck at lower bound
+   m <- mle2(minuslogl = f,
+             data = list(X = X),
+             start = list(rate = 0.01),
+             method = "L-BFGS-B",
+             control = list(trace = 1, parscale = 1e-4),
+             lower = list(rate = 1e-9))
+ 
+   profile(m, std.err=0.0001) ## finds new optimum
+ 
+   fsc <- function(X, rate) {
+     -sum(dexp(X, rate = rate*1e-4, log = TRUE))
+   }
+   msc <- mle2(minuslogl = fsc,
+             data = list(X = X),
+             start = list(rate = 100),
+             method = "L-BFGS-B",
+             control = list(trace = 1),
+             lower = list(rate = 1e-5))
+ 
+   ## does it work if we scale by hand?
+   ##   no, identical problem
+ }
> 
> ## works fine with a better starting point
> m <- mle2(minuslogl = f,
+           data = list(X = X),
+           start = list(rate = 0.001),
+           method = "L-BFGS-B",
+           control = list(trace = 1,parscale=1e-4),
+               lower = list(rate = 1e-9))
iter    0 value 12490.509482
final  value 10188.014396 
converged
> vcov(m)
             rate
rate 1.045669e-11
> confint(m)
       2.5 %       97.5 % 
9.605011e-05 1.087274e-04 
> 
> 
> ## works OK despite warnings about 1-dimensional opt. with N-M
> (m0 <- mle2(minuslogl = f,
+           data = list(X = X),
+           start = list(rate = 0.01),
+           method = "Nelder-Mead",
+           control = list(trace = 1, parscale = 1e-4)))
  Nelder-Mead direct search function minimizer
function value for initial parameters = 102397.310635
  Scaled convergence tolerance is 0.00152584
Stepsize computed as 10.000000
BUILD              2 112081.214500 102397.310635
EXTENSION          4 102397.310635 83062.026096
EXTENSION          6 83062.026096 44638.317097
HI-REDUCTION       8 63791.280079 44638.317097
REFLECTION        10 44638.317097 25773.036188
HI-REDUCTION      12 35146.785125 25773.036188
REFLECTION        14 25773.036188 16686.969324
HI-REDUCTION      16 21171.111238 16686.969324
REFLECTION        18 16686.969324 12490.509482
HI-REDUCTION      20 14529.847885 12490.509482
REFLECTION        22 12490.509482 10738.853151
HI-REDUCTION      24 11555.789799 10738.853151
REFLECTION        26 10738.853151 10209.598576
HI-REDUCTION      28 10415.334346 10209.598576
LO-REDUCTION      30 10209.598576 10191.680210
HI-REDUCTION      32 10191.680210 10190.329749
HI-REDUCTION      34 10190.329749 10188.037612
HI-REDUCTION      36 10188.497339 10188.037612
HI-REDUCTION      38 10188.089444 10188.037612
HI-REDUCTION      40 10188.037612 10188.018175
HI-REDUCTION      42 10188.018175 10188.016446
HI-REDUCTION      44 10188.016446 10188.014462
Exiting from Nelder Mead minimizer
    46 function evaluations used

Call:
mle2(minuslogl = f, start = list(rate = 0.01), method = "Nelder-Mead", 
    data = list(X = X), control = list(trace = 1, parscale = 1e-04))

Coefficients:
        rate 
0.0001022949 

Log-likelihood: -10188.01 
Warning message:
In optim(par = 0.01, fn = function (p)  :
  one-diml optimization by Nelder-Mead is unreliable: use optimize
> vcov(m0)
             rate
rate 1.046414e-11
> 
> confint(m0)
       2.5 %       97.5 % 
9.604965e-05 1.087271e-04 
> confint(m0,method="quad")
       2.5 %       97.5 % 
9.595477e-05 1.086351e-04 
> ## very similar (good quadratic surface, not surprising)
> 
> m1 <- mle2(minuslogl = f,
+           data = list(X = X),
+           start = list(rate = 0.01),
+           method = "BFGS",
+           control = list(trace = 1, parscale = 1e-4))
initial  value 102397.310635 
rate<0:  -0.08679214 
rate<0:  -0.009358428 
rate<0:  -0.5831727 
rate<0:  -0.1117319 
rate<0:  -0.01744372 
rate<0:  -0.07719334 
rate<0:  -0.01430754 
rate<0:  -0.001730383 
rate<0:  -0.08426903 
rate<0:  -0.01622577 
rate<0:  -0.002617114 
final  value 10188.014408 
converged
There were 11 warnings (use warnings() to see them)
> 
> 
> ## gets stuck? will have to investigate ...
> m2 <- mle2(minuslogl = f,
+            data = list(X = X),
+            start = list(rate = 0.01),
+            optimizer = "optimize",
+            lower=1e-9,upper=0.1)
> 
> vcov(m2)
             rate
rate 1.407176e-11
> 
> proc.time()
   user  system elapsed 
  0.856   0.196   1.065