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.

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.

Here are the problems that I encounter with different methods:

'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)

‘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.

'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: https://www.dropbox.com/s/aap2dhor5jyxy94/data.csv


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])
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

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.

  • $\begingroup$ This takes three seconds on my computer. $\endgroup$ – k20 Mar 7 '14 at 19:14
  • $\begingroup$ Yes. I know. I was saying " For another function with a more complicated setup, it takes 15 hours to get to the optimal point". This sample code is very simple. I'm trying to understand why it can't be solve using BFGS, which is much faster than Nelder-Mead. $\endgroup$ – Titanic Mar 7 '14 at 20:53
  • $\begingroup$ would it be helpful to post a code that takes some amount of time intermediate between 3 seconds and 15 hours? $\endgroup$ – k20 Mar 7 '14 at 21:00
  • $\begingroup$ maybe try L-BFGS-B counterintuitively its implementation may be better in some ways than that of BFGS even without limited memory and even without bounds $\endgroup$ – k20 Mar 7 '14 at 21:06
  • $\begingroup$ Sure. Here you are: dropbox.com/s/petdctq0he4n9z2/Logit2.py $\endgroup$ – Titanic Mar 8 '14 at 3:30

Your Answer

By clicking "Post Your Answer", you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.