title: BFGS Fails to Converge
body: <p>The model I'm working on is a multinomial logit choice model. It's a very specific dataset so other existing MNLogit libraries don't fit with my data. </p>
  
  <p>So basically, it's a very complex function which takes 11 parameters and returns a loglikelihood value. Then I need to find the optimal parameter values that can minimize the loglikelihood using scipy.optimize.minimize. </p>
  
  <p>Here are the problems that I encounter with different methods: </p>
  
  <p>'Nelder-Mead’: it works well, and always give me the correct answer. However, it's EXTREMELY slow. For another function with a more complicated setup, it takes 15 hours to get to the optimal point. At the same time, the same function takes only 1 hour on Matlab using fminunc (which uses BFGS by default)</p>
  
  <p>‘BFGS’: This is the method used by Matlab. It works well for any simply functions. However, for the function that I have, it always fails to converge and returns 'Desired error not necessarily achieved due to precision loss.’. I've spent lots of time playing around with the options but still failed to work. </p>
  
  <p>'Powell': It quickly converges successfully but returns a wrong answer. The code is printed below (x0 is the correct answer, Nelder-Mead works for whatever initial value), and you can get the data here: <a href="https://www.dropbox.com/s/aap2dhor5jyxy94/data.csv" rel="nofollow">https://www.dropbox.com/s/aap2dhor5jyxy94/data.csv</a></p>
  
  <p>Thanks!</p>
  
  <pre><code>import pandas as pd
  import numpy as np
  from scipy.optimize import minimize
  
  # https://www.dropbox.com/s/aap2dhor5jyxy94/data.csv
  df = pd.read_csv('data.csv', index_col=0)
  dfhh = df.hh
  B = df.ix[:,'b0':'b4'].values # NT*5
  P = df.ix[:,'p1':'p4'].values # NT*4
  F = df.ix[:,'f1':'f4'].values # NT*4
  SDV = df.ix[:,'lagb1':'lagb4'].values
  
  def Li(x):
      b1 = x[0] # coeff on prices
      b2 = x[1] # coeff on features
      a = x[2:7] # take first 4 values as alpha
      E = np.exp(a + b1*P + b2*F) # (1*4) + (NT*4) + (NT*4) build matrix (NT*J) for each exp()
      E = np.insert(E, 0, 1, axis=1) # (NT*5)
      denom = E.sum(1)
      return -np.log((B * E).sum(1) / denom).sum()
  
  
  x0 = np.array([-32.31028223, 0.23965953, 0.84739154, 0.25418215,-3.38757007,-0.38036966])
  np.random.seed(0)
  x0 = x0 + np.random.rand(6)
  
   minL = minimize(Li, x0, method='Nelder-Mead',options={'xtol': 1e-8, 'disp': True})
  # minL = minimize(Li, x0, method='BFGS')
  # minL = minimize(Li, x0, method='Powell', options={'xtol': 1e-12, 'ftol': 1e-12})
  print minL
  </code></pre>
  
  <p>Update: 03/07/14 Simpler Version of the Code
  Now Powell works well with very small tolerance, however the speed of Powell is slower than Nelder-Mead in this case. BFGS still fails to work. </p>