Important: Read this before posting to this forum

  1. This forum is for questions related to the use of Apollo. We will answer some general choice modelling questions too, where appropriate, and time permitting. We cannot answer questions about how to estimate choice models with other software packages.
  2. There is a very detailed manual for Apollo available at http://www.ApolloChoiceModelling.com/manual.html. This contains detailed descriptions of the various Apollo functions, and numerous examples are available at http://www.ApolloChoiceModelling.com/examples.html. In addition, help files are available for all functions, using e.g. ?apollo_mnl
  3. Before asking a question on the forum, users are kindly requested to follow these steps:
    1. Check that the same issue has not already been addressed in the forum - there is a search tool.
    2. Ensure that the correct syntax has been used. For any function, detailed instructions are available directly in Apollo, e.g. by using ?apollo_mnl for apollo_mnl
    3. Check the frequently asked questions section on the Apollo website, which discusses some common issues/failures. Please see http://www.apollochoicemodelling.com/faq.html
    4. Make sure that R is using the latest official release of Apollo.
  4. If the above steps do not resolve the issue, then users should follow these steps when posting a question:
    1. provide full details on the issue, including the entire code and output, including any error messages
    2. posts will not immediately appear on the forum, but will be checked by a moderator first. This may take a day or two at busy times. There is no need to submit the post multiple times.

CALCULATION ISSUE - Log-likelihood calculation fails at values close to the starting values!

Ask questions about errors you encouunter. Please make sure to include full details about your model specifications, and ideally your model file.
Post Reply
thoscha
Posts: 1
Joined: 12 Mar 2023, 11:13

CALCULATION ISSUE - Log-likelihood calculation fails at values close to the starting values!

Post by thoscha »

Hi there,

I am estimating a mode choice model with data from an SP experiment on mobility coins. First step was to estimate a simple MNL in WTP space with two components.
Model run by ge75mum using Apollo 0.2.9 on R 4.3.1 for Windows.
www.ApolloChoiceModelling.com

Model name : 4.1_mnl_wtp_space.R
Model description : mnl in wtp, No income individuals removed. Preliminary model. Includes shifts (multipliers) for purposes
Model run at : 2023-07-12 00:53:30.562108
Estimation method : bgw
Model diagnosis : Relative function convergence
Optimisation diagnosis : Maximum found
hessian properties : Negative definitive
maximum eigenvalue : -1.61998
Number of individuals : 1053
Number of rows in database : 11172
Number of modelled outcomes : 11172
mocono : 5448
mocosi : 5724

Number of cores used : 1
Model without mixing

LL(start) : -9359.71
LL (whole model) at equal shares, LL(0) : -12822.54
LL (whole model) at observed shares, LL(C) : -12146.8
LL(final, whole model) : -9349.34
Rho-squared vs equal shares : 0.2709
Adj.Rho-squared vs equal shares : 0.2664
Rho-squared vs observed shares : 0.2303
Adj.Rho-squared vs observed shares : 0.2261
AIC : 18812.68
BIC : 19189.05

LL(0,mocono) : -6283.17
LL(final,mocono) : -4602.63
LL(0,mocosi) : -6539.37
LL(final,mocosi) : -4746.71

Estimated parameters : 57
Time taken (hh:mm:ss) : 00:03:2.82
pre-estimation : 00:00:43.12
estimation : 00:00:8.31
initial estimation : 00:00:6.57
estimation after rescaling : 00:00:1.73
post-estimation : 00:02:11.39
Iterations : 11
initial estimation : 10
estimation after rescaling : 1

Unconstrained optimisation.

Estimates:
Estimate s.e. t.rat.(0) Rob.s.e. Rob.t.rat.(0)
asc_pt_base 0.000000 NA NA NA NA
asc_c_base -1.254459 0.374202 -3.3524 0.678908 -1.84776
asc_c_work 0.000000 NA NA NA NA
asc_c_errand 0.532274 0.194014 2.7435 0.325964 1.63293
asc_c_leisure -0.453211 0.197837 -2.2908 0.347810 -1.30304
asc_c_female 0.000000 NA NA NA NA
asc_c_male 0.299004 0.061825 4.8363 0.108138 2.76502
asc_c_divers 1.618516 0.632749 2.5579 0.811620 1.99418
asc_c_age 0.040300 0.012786 3.1519 0.023521 1.71337
asc_c_age_2 -4.4839e-04 1.3464e-04 -3.3302 2.4721e-04 -1.81377
asc_b_base -3.850031 0.319654 -12.0444 0.535656 -7.18751
asc_b_work 0.000000 NA NA NA NA
asc_b_errand 0.203590 0.149451 1.3623 0.261974 0.77714
asc_b_leisure 0.106287 0.148497 0.7157 0.255373 0.41620
asc_b_female 0.000000 NA NA NA NA
asc_b_male 0.175134 0.055085 3.1793 0.090286 1.93976
asc_b_divers 0.146116 0.439456 0.3325 1.130482 0.12925
asc_b_age 0.058708 0.010796 5.4381 0.017788 3.30043
asc_b_age_2 -5.9628e-04 1.1541e-04 -5.1666 1.9361e-04 -3.07982
asc_w_base -2.587003 0.654354 -3.9535 1.193279 -2.16798
asc_w_work 0.000000 NA NA NA NA
asc_w_errand 0.299264 0.412286 0.7259 0.854242 0.35033
asc_w_leisure 0.072829 0.421915 0.1726 0.857031 0.08498
asc_w_female 0.000000 NA NA NA NA
asc_w_male -0.177903 0.103643 -1.7165 0.184262 -0.96549
asc_w_divers 2.212046 0.449991 4.9158 0.711920 3.10715
asc_w_age 0.052144 0.020034 2.6027 0.032396 1.60957
asc_w_age_2 -3.2780e-04 2.0570e-04 -1.5936 3.4243e-04 -0.95727
s_mocono 1.000000 NA NA NA NA
s_mocosi 0.972736 0.029385 33.1034 0.032661 29.78236
b_cost_base -0.150005 0.011513 -13.0297 0.018597 -8.06591
elast_cost_dist -0.717950 0.053069 -13.5287 0.107998 -6.64782
b_mcost_rembudget75days15 -0.036301 0.009344 -3.8850 0.010686 -3.39720
b_mcost_rembudget75days20 -0.013583 0.008767 -1.5493 0.009071 -1.49732
b_mcost_rembudget75days25 -0.008614 0.009316 -0.9246 0.011301 -0.76224
b_mcost_rembudget50days15 -0.080285 0.014135 -5.6800 0.017842 -4.49975
b_mcost_rembudget50days20 -0.011427 0.012140 -0.9413 0.012592 -0.90748
b_mcost_rembudget50days25 -0.013329 0.012619 -1.0563 0.014097 -0.94551
b_mcost_rembudget25days15 -0.080614 0.011172 -7.2158 0.013679 -5.89342
b_mcost_rembudget25days20 -0.054174 0.010026 -5.4033 0.012696 -4.26700
b_mcost_rembudget25days25 -0.058214 0.011242 -5.1783 0.014191 -4.10215
elast_mcost_dist -0.782566 0.086440 -9.0533 0.101831 -7.68494
elast_mcost_inc -0.255097 0.094963 -2.6863 0.116134 -2.19657
v_tt_pt_base 0.400504 0.050325 7.9584 0.078902 5.07597
v_tt_c_base 0.289455 0.066391 4.3599 0.102945 2.81174
v_tt_b_base 0.407471 0.047868 8.5124 0.076805 5.30525
v_tt_w_base 0.535308 0.079314 6.7492 0.156265 3.42565
elast_vtt_dist 0.546130 0.069166 7.8960 0.134232 4.06855
mult_vtt_pt_leisure 0.910207 0.126314 7.2059 0.183174 4.96910
mult_vtt_pt_shopping 0.627084 0.107013 5.8599 0.159394 3.93419
mult_vtt_c_leisure 0.429602 0.188720 2.2764 0.310528 1.38346
mult_vtt_c_shopping 0.458857 0.191549 2.3955 0.292903 1.56659
mult_vtt_w_leisure 0.740412 0.111332 6.6505 0.211563 3.49972
mult_vtt_w_shopping 0.757512 0.114156 6.6358 0.235887 3.21134
mult_vtt_b_leisure 0.865357 0.087197 9.9241 0.129897 6.66187
mult_vtt_b_shopping 0.693743 0.079515 8.7247 0.125469 5.52918
v_aet_pt 0.116845 0.031697 3.6863 0.064021 1.82512
v_freq_pt 0.061717 0.030121 2.0490 0.046273 1.33377
v_trans_pt 0.705802 0.156862 4.4995 0.205693 3.43134
b_mgain_b_base 0.208319 0.086838 2.3989 0.099356 2.09670
b_qual_b_bad 0.000000 NA NA NA NA
b_qual_b_med 0.762972 0.061527 12.4007 0.069484 10.98061
b_qual_b_vgood 0.942095 0.062026 15.1886 0.067958 13.86289
b_weather_c_rainy 0.000000 NA NA NA NA
b_weather_c_sunny -0.373262 0.068852 -5.4212 0.083909 -4.44843
b_weather_b_rainy 0.000000 NA NA NA NA
b_weather_b_sunny 2.830420 0.073280 38.6249 0.096785 29.24438
b_weather_w_rainy 0.000000 NA NA NA NA
b_weather_w_sunny 2.104389 0.111653 18.8475 0.165460 12.71845

apollo_probabilities
----------------------
function(apollo_beta, apollo_inputs, functionality="estimate"){

### Attach inputs and detach after function exit
apollo_attach(apollo_beta, apollo_inputs)
on.exit(apollo_detach(apollo_beta, apollo_inputs))

### Create list of probabilities P
P = list()

### List of utilities: these must use the same names as in cnl_settings, order is irrelevant
V = list()

# ascs with shifts
asc_pt = asc_pt_base
asc_c = asc_c_base +
asc_c_work * (purpose == 1) + asc_c_errand * (purpose == 2) + asc_c_leisure * (purpose == 3) +
asc_c_female * (p_gender == "weiblich") + asc_c_male * (p_gender == "männlich") + asc_c_divers * (p_gender == "divers") +
asc_c_age * p_age + asc_c_age_2 * p_age**2
asc_b = asc_b_base +
asc_b_work * (purpose == 1) + asc_b_errand * (purpose == 2) + asc_b_leisure * (purpose == 3) +
asc_b_female * (p_gender == "weiblich") + asc_b_male * (p_gender == "männlich") + asc_b_divers * (p_gender == "divers") +
asc_b_age * p_age + asc_b_age_2 * p_age**2
asc_w = asc_w_base +
asc_w_work * (purpose == 1) + asc_w_errand * (purpose == 2) + asc_w_leisure * (purpose == 3) +
asc_w_female * (p_gender == "weiblich") + asc_w_male * (p_gender == "männlich") + asc_w_divers * (p_gender == "divers") +
asc_w_age * p_age + asc_w_age_2 * p_age**2

# b_cost_base = cost sensitivity at mean distance, elast_cost_dist < 0, decreasing cost elasticity in relation to distance
b_cost = (b_cost_base) *
((trip_dist/mean(trip_dist)) ^ elast_cost_dist)

#b_mcost = b_mcost_base

b_mcost = (b_mcost_rembudget75days15 * (cbudget == 75) * (days_into_month == 15) +
b_mcost_rembudget75days20 * (cbudget == 75) * (days_into_month == 20) +
b_mcost_rembudget75days25 * (cbudget == 75) * (days_into_month == 25) +
b_mcost_rembudget50days15 * (cbudget == 50) * (days_into_month == 15) +
b_mcost_rembudget50days20 * (cbudget == 50) * (days_into_month == 20) +
b_mcost_rembudget50days25 * (cbudget == 50) * (days_into_month == 25) +
b_mcost_rembudget25days15 * (cbudget == 25) * (days_into_month == 15) +
b_mcost_rembudget25days20 * (cbudget == 25) * (days_into_month == 20) +
b_mcost_rembudget25days25 * (cbudget == 25) * (days_into_month == 25)) *
((trip_dist/mean(trip_dist)) ^ elast_mcost_dist) *
((p_inc_lin/mean(p_inc_lin)) ^ elast_mcost_inc)

b_mgain_b = (b_mgain_b_base)

v_tt_pt = v_tt_pt_base *
((trip_dist/mean(trip_dist)) ^ elast_vtt_dist) *
(1*(purpose == 1) + mult_vtt_pt_shopping*(purpose == 2) + mult_vtt_pt_leisure * (purpose == 3))

v_tt_c = v_tt_c_base *
((trip_dist/mean(trip_dist)) ^ elast_vtt_dist) *
(1*(purpose == 1) + mult_vtt_c_shopping*(purpose == 2) + mult_vtt_c_leisure * (purpose == 3))

v_tt_b = v_tt_b_base *
((trip_dist/mean(trip_dist)) ^ elast_vtt_dist) *
(1*(purpose == 1) + mult_vtt_b_shopping*(purpose == 2) + mult_vtt_b_leisure * (purpose == 3))

v_tt_w = v_tt_w_base *
((trip_dist/mean(trip_dist)) ^ elast_vtt_dist) *
(1*(purpose == 1) + mult_vtt_w_shopping*(purpose == 2) + mult_vtt_w_leisure * (purpose == 3))

# mocono
# reference alternative
V[["pt_mocono"]] = asc_pt +
b_cost * (cost_pt + v_tt_pt * tt_pt + v_aet_pt * aett_pt + v_freq_pt * freq_pt + v_trans_pt * trans_pt)
# other alternatives
V[["c_mocono"]] = asc_c +
b_cost * (cost_c + v_tt_c * tt_c) +
b_weather_c_rainy * (weather == "rainy") + b_weather_c_sunny * (weather == "sunny")
V[["b_mocono"]] = asc_b +
b_cost * (v_tt_b * tt_b) +
b_qual_b_bad * (quality_b == "bad") + b_qual_b_med * (quality_b == "medium") + b_qual_b_vgood * (quality_b == "very good") +
b_weather_b_rainy * (weather == "rainy") + b_weather_b_sunny * (weather == "sunny")
V[["w_mocono"]] = asc_w +
b_cost * (v_tt_w * tt_w) +
b_weather_w_rainy * (weather == "rainy") + b_weather_w_sunny * (weather == "sunny")

# mocosi
# reference alternative
V[["pt_mocosi"]] = asc_pt +
b_cost * (cost_pt + v_tt_pt * tt_pt + v_aet_pt * aett_pt + v_freq_pt * freq_pt + v_trans_pt * trans_pt) +
b_mcost * ccost_pt
# other alternatives
V[["c_mocosi"]] = asc_c +
b_cost * (cost_c + v_tt_c * tt_c) +
b_mcost * ccost_c +
b_weather_c_rainy * (weather == "rainy") + b_weather_c_sunny * (weather == "sunny")
V[["b_mocosi"]] = asc_b +
b_cost * (v_tt_b * tt_b) +
b_mgain_b * cgain_b +
b_qual_b_bad * (quality_b == "bad") + b_qual_b_med * (quality_b == "medium") + b_qual_b_vgood * (quality_b == "very good") +
b_weather_b_rainy * (weather == "rainy") + b_weather_b_sunny * (weather == "sunny")
V[["w_mocosi"]] = asc_w +
b_cost * (v_tt_w * tt_w) +
b_weather_w_rainy * (weather == "rainy") + b_weather_w_sunny * (weather == "sunny")

# Compute probabilities for the mocono mode choice part of the data using MNL model
mnl_settings_mocono = list(
alternatives = c(pt_mocono = 1, c_mocono = 2, b_mocono = 3, w_mocono = 4),
avail = list(pt_mocono = pt_availability, c_mocono = car_availability, b_mocono = bicycle_availability, w_mocono = walk_availability),
choiceVar = CHOICE,
utilities = list(pt_mocono = s_mocono * V[["pt_mocono"]],
c_mocono = s_mocono * V[["c_mocono"]],
b_mocono = s_mocono * V[["b_mocono"]],
w_mocono = s_mocono * V[["w_mocono"]]),
rows = (exp_nr == 0)
)

P[["mocono"]] = apollo_mnl(mnl_settings_mocono, functionality)

# Compute probabilities for the RP mode choice part of the data using MNL model
mnl_settings_mocosi = list(
alternatives = c(pt_mocosi = 5, c_mocosi = 6, b_mocosi = 7, w_mocosi = 8),
avail = list(pt_mocosi = pt_availability, c_mocosi = car_availability, b_mocosi = bicycle_availability, w_mocosi = walk_availability),
choiceVar = CHOICE,
utilities = list(pt_mocosi = s_mocosi * V[["pt_mocosi"]],
c_mocosi = s_mocosi * V[["c_mocosi"]],
b_mocosi = s_mocosi * V[["b_mocosi"]],
w_mocosi = s_mocosi * V[["w_mocosi"]]),
rows = (exp_nr == 1)
)

P[["mocosi"]] = apollo_mnl(mnl_settings_mocosi, functionality)

# Combined model
P = apollo_combineModels(P, apollo_inputs, functionality, components = c("mocono","mocosi"))

# Take product across observation for same individual
P = apollo_panelProd(P, apollo_inputs, functionality)

# Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
All good so far. Now I would like to estimate an MMNL with random parameters for 4 different travel times variables (walk, bike, car and PT) as well as one for the cost parameter. see here:

Code: Select all

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "sobol",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws_tt_pt", "draws_tt_c", "draws_tt_b","draws_tt_w","draws_tc")
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  # Log normal
  randcoeff[["v_tt_c"]]  = exp( mu_log_v_tt_c_base  + sigma_log_v_tt_c_base  * draws_tt_c  ) *
    ((trip_dist/mean(trip_dist)) ^ elast_vtt_dist)  *
    (1*(purpose == 1) + mult_vtt_c_shopping*(purpose == 2) + mult_vtt_c_leisure * (purpose == 3))
  
  randcoeff[["v_tt_w"]]  = exp( mu_log_v_tt_w_base  + sigma_log_v_tt_w_base  * draws_tt_w  ) *
    ((trip_dist/mean(trip_dist)) ^ elast_vtt_dist)  *
    (1*(purpose == 1) + mult_vtt_w_shopping*(purpose == 2) + mult_vtt_w_leisure * (purpose == 3))
  
  randcoeff[["v_tt_b"]]  = exp( mu_log_v_tt_b_base  + sigma_log_v_tt_b_base  * draws_tt_b  ) *
    ((trip_dist/mean(trip_dist)) ^ elast_vtt_dist)  *
    (1*(purpose == 1) + mult_vtt_b_shopping*(purpose == 2) + mult_vtt_b_leisure * (purpose == 3))
  
  randcoeff[["v_tt_pt"]] = exp( mu_log_v_tt_pt_base + sigma_log_v_tt_pt_base * draws_tt_pt ) *
    ((trip_dist/mean(trip_dist)) ^ elast_vtt_dist)  *
    (1*(purpose == 1) + mult_vtt_pt_shopping*(purpose == 2) + mult_vtt_pt_leisure * (purpose == 3))
  
  # Negative log normal
  randcoeff[["b_cost"]]  = -exp( mu_b_cost_base +  sigma_b_cost_base  * draws_tc ) *
    ((trip_dist/mean(trip_dist)) ^ elast_cost_dist)
  
  return(randcoeff)
}
For the vtt's we assume a lognormal whereas for the cost coefficient we assume a negative lognormal. These are the starting values:

Code: Select all

mu_b_cost_base = -1,
sigma_b_cost_base = 0,
elast_cost_dist = 0,
mu_log_v_tt_pt_base = 2,
  sigma_log_v_tt_pt_base = 0,
  mu_log_v_tt_c_base = 2,
  sigma_log_v_tt_c_base = 0,
  mu_log_v_tt_b_base = 2,
  sigma_log_v_tt_b_base = 0,
  mu_log_v_tt_w_base = 2,
  sigma_log_v_tt_w_base = 0,
  elast_vtt_dist = 0,
  mult_vtt_pt_leisure = 1,
  mult_vtt_pt_shopping = 1,
  mult_vtt_c_leisure = 1,
  mult_vtt_c_shopping = 1,
  mult_vtt_w_leisure = 1,
  mult_vtt_w_shopping = 1,
  mult_vtt_b_leisure = 1,
  mult_vtt_b_shopping = 1,
For the other parameters, I read in the MNL's estimates as starting values. However, I run into troubles =/.

Would appreciate any help here... I assume it might be related to the starting values... but I am stuck at the moment.

Thomas
stephanehess
Site Admin
Posts: 998
Joined: 24 Apr 2020, 16:29

Re: CALCULATION ISSUE - Log-likelihood calculation fails at values close to the starting values!

Post by stephanehess »

Hi

the issue is your starting values. With lognormals, you should not start at positive values for the mean of the log of the coefficient - see http://apollochoicemodelling.com/faq.html

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
Post Reply