weibull mle with optim in R I am fitting a weibull distribution to some data. Let us say we have the data in warpbreaks, that is, warpbreaks$breaks in R.
My code is as follows:
#Weibull distribution
#f(x) = (a/b) (x/b)^(a-1) exp(- (x/b)^a)
#shape= a and scale= b

ll.weibull<-function(dat,par){
  a=exp(par[1])
  b=exp(par[2])

  ll=-sum(dweibull(dat,a,scale=b,log=T))

}

a=0.4
b=0.4
par=c(a,b)
dat=warpbreaks$breaks


weibull.optim<-optim(par=par,fn=ll.weibull,dat=dat)
weibull.optim #E(X) = b Γ(1 + 1/a)
#check mean:
exp(weibull.optim$par[2])*gamma(1+1/(exp(weibull.optim$par[1])))
rm(a,b,par,dat)

Now, this works fine. Even when we vary the initial values in a and b, we seem to get the same values for the parameters after optimalization.
But if we use another data set, namely the second column of the following data set, with all values of 0 removed,
http://www.uio.no/studier/emner/matnat/math/STK4540/h14/assignment/stormcntpevd.txt
, the initial values plays too much of an importance. I do not seem to be able to optimize.
 A: Suggestions:


*

*try establishing starting values based on the method of moments, or using the survreg package to get optimal starting values. library(survival); survreg(Surv(dat)~1) for your data set gives me "scale=2.2, shape=12", which I think corresponds to shape=1/2.2, scale=exp(12), but you should check).  (You may not want to use survreg even to get initial parameters, but I think the method of moments should give you similar results.) 

*if the initial parameters differ greatly in magnitude try setting the parscale option in optim.  


The second doesn't appear to be necessary: 
weibull.optim <- optim(par=c(log(1/2.2),12),fn=ll.weibull,dat=dat)

seems to work fine for me with your suggested data.  
My guess is that you just didn't get close enough with your initial values (you didn't say what you tried).
A: I don't have enough knowledge about your optimization problem. However in the past when a local optimizer has failed to provide a optimal solution, one strategy I have used is to run a global optimizer such as genetic algorithm, differential evolution or simulated annealing and then use the output from the global optimizer as a starting point for a local optimizer. The only caveat is you need to specify bounds for the global optimizer which I'm assuming you can easily specify with a good knowledge of your optimization problem. See below, I have used differential evoution + optim function R and gave me good results. Typically DEoptim would give near global and optim would then search for global solution. In this case DEoptim gave an excellent solution. The advantage of running a local + global is that in case if you have multiple solution, you would be able to find them all in few iterations.  Hope this helps. I modified the function so I can use it in the DEoptim package, let me know if I interpreted correctly.
Reference:
http://www.sciencedirect.com/science/article/pii/S0377221702004010
#Weibull distribution
#f(x) = (a/b) (x/b)^(a-1) exp(- (x/b)^a)
#shape= a and scale= b

library(DEoptim)

ll.weibull<-function(par){
  a<-exp(par[1])
  b<-exp(par[2])

  ll=-sum(dweibull(warpbreaks$breaks,a,scale=b,log=T))

}

## Optimize using a global optimzer suh as GA or Differential Evolution to 
## near global optimal 

lower <- c(0.001,0.001)
upper <- c(5,5)

outDEoptim <- DEoptim(ll.weibull, lower, upper,DEoptim.control(NP = 50,itermax = 200,trace=FALSE))
outDEoptim

par = c(outDEoptim$optim$bestmem)

## Use near optimal from Global optimzer as a starting point to get global
## solution

weibull.optim<-optim(par=par,fn=ll.weibull)
weibull.optim #E(X) = b Γ(1 + 1/a)


## Global optimizer Soution
global.par = c(outDEoptim$optim$bestmem)
global.obj = outDEoptim$optim$bestval
global.par
global.obj

## Global + Local Optimizer
loc.glob.par = c(weibull.optim$par)
    loc.global.obj =weibull.optim$value
loc.glob.par
loc.global.obj

Outputs from both optimizers are identical
global.par
     par1      par2 
0.8252952 3.4622785 
global.obj
[1] 211.6973

loc.glob.par
     par1      par2 
0.8252952 3.4622785 
loc.global.obj
[1] 211.6973

