ApolloChoiceModelling forum

Forum for users of the Apollo choice modelling software
http://www.apollochoicemodelling.com/forum/

Calculating logsum from the estimated model

http://www.apollochoicemodelling.com/forum/viewtopic.php?t=338

Page 1 of 1

Calculating logsum from the estimated model

Posted: 06 Oct 2021, 18:41
by hrrichard
Hi Stephane and David,

I am working on a project related to accessibility. The model has been estimated, and I would like to calculate/obtain the logsum for each individual in the dataset using the estimated parameters of the model.

Could you please advise how to do this?

Thank you for your help!

Richard

Re: Calculating logsum from the estimated model

Posted: 07 Oct 2021, 17:19
by dpalma
Hi Richard,

There is no pre-coded function to extract the logsum from a logit model in Apollo. However, you can easily calculate it by calculating the deterministic utilities (V) using the estimated parameters, and applying the logsum formula to them.

Cheers
David

Re: Calculating logsum from the estimated model

Posted: 07 Oct 2021, 19:11
by hrrichard
Hi David,

Thank you for your advice! I worked out a code below, could you please advise if it makes sense to you:

The deterministic utilities from the model are as follows:

Code: Select all

V = list()
  V[['MMB']]  = asc_MMB  + tt_1*mm.tt1 + tt_2*mm.tt2 + tc*mm.tc + cont_1*(mm.cont==0) + cont_2*(mm.cont==1) + cont_3*(mm.cont==2) +
    cont_4*(mm.cont==3) + alig_1*(mm.aligm==0) + alig_2*(mm.aligm==1) + alig_3*(mm.aligm==2)
  V[['MBUS']]  = asc_MB  + tt_1*mb.tt1 + tt_2*mb.tt2 + tc*mb.tc + cont_1*(mb.cont==0) + cont_2*(mb.cont==1) + cont_3*(mb.cont==2) +
    cont_4*(mb.cont==3) + alig_1*(mb.aligb==0) + alig_2*(mb.aligb==1) + alig_3*(mb.aligb==2) 
  V[['TAXI']]  = asc_TX + tt_t*tx.ttt + tc_t*tx.ttc
  V[['NO']] = att*mm.attr + ratt_1*(mm.rattr==0) + ratt_2*(mm.rattr==1) + ratt_3*(mm.rattr==2) + ratt_4*(mm.rattr==3)
The logsum calculation using the model estimates would be then:

Code: Select all

V_MMB  = model$estimate['asc_MMB']  + model$estimate['tt_1']*database$mm.tt1 + model$estimate['tt_2']*database$mm.tt2 + model$estimate['tc']*database$mm.tc + model$estimate['cont_1']*(database$mm.cont==0) + model$estimate['cont_2']*(database$mm.cont==1) + model$estimate['cont_3']*(database$mm.cont==2) +
    model$estimate['cont_4']*(database$mm.cont==3) + model$estimate['alig_1']*(database$mm.aligm==0) + model$estimate['alig_2']*(database$mm.aligm==1) + model$estimate['alig_3']*(database$mm.aligm==2)
V_MB  = model$estimate['asc_MB']  + model$estimate['tt_1']*database$mm.tt1 + model$estimate['tt_2']*database$mm.tt2 + model$estimate['tc']*database$mm.tc + model$estimate['cont_1']*(database$mm.cont==0) + model$estimate['cont_2']*(database$mm.cont==1) + model$estimate['cont_3']*(database$mm.cont==2) +
    model$estimate['cont_4']*(database$mm.cont==3) + model$estimate['alig_1']*(database$mm.aligm==0) + model$estimate['alig_2']*(database$mm.aligm==1) + model$estimate['alig_3']*(database$mm.aligm==2)
V_Taxi  = model$estimate['asc_TX'] + model$estimate['tt_t']*database$tx.ttt + model$estimate['tc_t']*database$tx.ttc
V_NO = model$estimate['att']*database$mm.attr + model$estimate['ratt_1']*(database$mm.rattr==0) + model$estimate['ratt_2']*(database$mm.rattr==1) + model$estimate['ratt_3']*(database$mm.rattr==2) + model$estimate['ratt_4']*(database$mm.rattr==3)

Logsum= exp(V_MMB)+exp(V_MB)+exp(V_Taxi)+exp(V_NO)
Thank you for your advice!

Cheers,
Richard

Re: Calculating logsum from the estimated model

Posted: 07 Oct 2021, 20:21
by dpalma
Hi Richard,

If your model is an MNL, then you would only be missing the log:

Code: Select all

Logsum= log( exp(V_MMB)+exp(V_MB)+exp(V_Taxi)+exp(V_NO) )
Cheers
David

Re: Calculating logsum from the estimated model

Posted: 08 Oct 2021, 08:23
by hrrichard
Hi David,

Thanks for your swift reply! Our estimated model has a more complex preference structure than the MNL, accounting for random, deterministic and latent-variable induced heterogeneity. I include the code of the model below. Could you please help with suggesting how to calculate the logsum when we account for different sources of heterogeneity?

Code: Select all

# ################################################################# #
#### DEFINE RANDOM COMPONENTS                                    ####
# ################################################################# #

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType = "pmc",
  interNDraws    = 1000,
  interUnifDraws = c(),
  interNormDraws = c("draws_LV_Tr_perc",
                     "draws_tt_1","draws_tt_2","draws_tc","draws_tc_t","draws_tt_t",
                     "draws_att","draws_ratt_2","draws_ratt_3","draws_ratt_4",
                     "draws_alig_2","draws_alig_3",
                     "draws_cont_2","draws_cont_3","draws_cont_4"),
  intraDrawsType = "halton",
  intraNDraws    = 0,
  intraUnifDraws = c(),
  intraNormDraws = c()
)

### Create random parameters
apollo_randCoeff = function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  
  randcoeff[["LV_Tr_perc"]] =  draws_LV_Tr_perc + gamma_Repeat_vis_LV_Transport_perception*Repeat_vis +
    gamma_Long_stay_LV_Transport_perception*four_nights_long
 
  randcoeff[["tt_1"]] = tt_1_mu +  tt_1_sig * draws_tt_1+gamma_tt_1_long_stay*four_nights_long +gamma_tt_1_Repeat_vis*Repeat_vis
  randcoeff[["tt_2"]] = tt_2_mu +  tt_2_sig * draws_tt_2+gamma_tt_2_long_stay*four_nights_long +gamma_tt_2_Repeat_vis*Repeat_vis
  randcoeff[["tc"]] = tc_mu +  tc_sig * draws_tc+gamma_tc_long_stay*four_nights_long +  gamma_tc_Repeat_vis*Repeat_vis
  randcoeff[["tt_t"]] = tt_t_mu +  tt_t_sig * draws_tt_t+gamma_tt_t_long_stay*four_nights_long +gamma_tt_t_Repeat_vis*Repeat_vis
  randcoeff[["tc_t"]] = tc_t_mu +  tc_t_sig * draws_tc_t+gamma_tc_t_long_stay*four_nights_long +gamma_tc_t_Repeat_vis*Repeat_vis
  randcoeff[["att"]] = att_mu +  att_sig * draws_att+gamma_att_long_stay*four_nights_long + gamma_att_Repeat_vis*Repeat_vis
  randcoeff[["ratt_2"]] = ratt_2_mu +  ratt_2_sig * draws_ratt_2+gamma_ratt_2_long_stay*four_nights_long +  gamma_ratt_2_Repeat_vis*Repeat_vis
  randcoeff[["ratt_3"]] = ratt_3_mu +  ratt_3_sig * draws_ratt_3+gamma_ratt_3_long_stay*four_nights_long +  gamma_ratt_3_Repeat_vis*Repeat_vis
  randcoeff[["ratt_4"]] = ratt_4_mu +  ratt_4_sig * draws_ratt_4+gamma_ratt_4_long_stay*four_nights_long +  gamma_ratt_4_Repeat_vis*Repeat_vis
  randcoeff[["cont_2"]] = cont_2_mu +  cont_2_sig * draws_cont_2+gamma_cont_2_long_stay*four_nights_long +  gamma_cont_2_Repeat_vis*Repeat_vis
  randcoeff[["cont_3"]] = cont_3_mu +  cont_3_sig * draws_cont_3+gamma_cont_3_long_stay*four_nights_long +  gamma_cont_3_Repeat_vis*Repeat_vis
  randcoeff[["cont_4"]] = cont_4_mu +  cont_4_sig * draws_cont_4+gamma_cont_4_long_stay*four_nights_long +  gamma_cont_4_Repeat_vis*Repeat_vis
  randcoeff[["alig_2"]] = alig_2_mu +  alig_2_sig * draws_alig_2+gamma_alig_2_long_stay*four_nights_long +  gamma_alig_2_Repeat_vis*Repeat_vis
  randcoeff[["alig_3"]] = alig_3_mu +  alig_3_sig * draws_alig_3+gamma_alig_3_long_stay*four_nights_long +  gamma_alig_3_Repeat_vis*Repeat_vis
  
  return(randcoeff)
}

# ################################################################# #
#### GROUP AND VALIDATE INPUTS                                   ####
# ################################################################# #

apollo_inputs = apollo_validateInputs()

# ################################################################# #
#### DEFINE MODEL AND LIKELIHOOD FUNCTION                        ####
# ################################################################# #

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()
  
  
  ## PERCVULNEV
  OL_settings1 = list(outcomeOrdered=Q49_1.1, 
                      V=zeta_lv_perc_1*LV_Tr_perc,
                      tau=c(tau_lv_perc_1_1,tau_lv_perc_1_2),
                      rows=(Sequence==1),
                      componentName = "indic_1")
  
  OL_settings3 = list(outcomeOrdered=Q49_3, 
                      V=zeta_lv_perc_3*LV_Tr_perc,
                      tau=c(tau_lv_perc_3_1,tau_lv_perc_3_2,tau_lv_perc_3_3,tau_lv_perc_3_4),
                      rows=(Sequence==1),
                      componentName = "indic_3")
  
  P[["indic_1"]] = apollo_ol(OL_settings1, functionality)
  P[["indic_3"]] = apollo_ol(OL_settings3, functionality)
  
  
  ### Creating interactions with latent variables
  tt_1_val=tt_1+LV_tt_1_sig*LV_Tr_perc
  tt_2_val=tt_2+LV_tt_2_sig*LV_Tr_perc
  tc_val=tc+LV_tc_sig*LV_Tr_perc
  cont_2_val = cont_2+LV_cont_2_sig*LV_Tr_perc
  cont_3_val = cont_3+LV_cont_3_sig*LV_Tr_perc
  cont_4_val = cont_4+LV_cont_4_sig*LV_Tr_perc
  alig_2_val = alig_2+LV_alig_2_sig*LV_Tr_perc
  alig_3_val = alig_3+LV_alig_3_sig*LV_Tr_perc
  tt_t_val=tt_t+LV_tt_t_sig*LV_Tr_perc
  tc_t_val=tc_t+LV_tc_t_sig*LV_Tr_perc
  
  ### List of utilities: these must use the same names as in mnl_settings, order is irrelevant
  
  V = list()
  V[['MMB']]  = asc_MMB  + tt_1_val*mm.tt1 + tt_2_val*mm.tt2 + tc_val*mm.tc + cont_1*(mm.cont==0) + cont_2_val*(mm.cont==1) + cont_3_val*(mm.cont==2) +
    cont_4_val*(mm.cont==3) + alig_1*(mm.aligm==0) + alig_2_val*(mm.aligm==1) + alig_3_val*(mm.aligm==2) 
  V[['MBUS']]  = asc_MB  + tt_1_val*mb.tt1 + tt_2_val*mb.tt2 + tc_val*mb.tc + cont_1*(mb.cont==0) + cont_2_val*(mb.cont==1) + cont_3_val*(mb.cont==2) +
    cont_4_val*(mb.cont==3) + alig_1*(mb.aligb==0) + alig_2_val*(mb.aligb==1) + alig_3_val*(mb.aligb==2) 
  V[['TAXI']]  = asc_TX + tt_t_val*tx.ttt + tc_t_val*tx.ttc 
  V[['NO']] = asc_NO + att*mm.attr + ratt_1*(mm.rattr==0) + ratt_2*(mm.rattr==1) + ratt_3*(mm.rattr==2) + ratt_4*(mm.rattr==3)
  
  
  ### Define settings for MNL model component
  mnl_settings = list(
    alternatives  = c(MMB=1, MBUS=2, TAXI=3, NO=4), 
    avail         = list(MMB=1, MBUS=1, TAXI=1, NO=1), 
    choiceVar     = Choice1,
    V             = V,
    componentName = "Mode_choice"
  )

  ### Compute probabilities for MNL model component
  P[["Mode_choice"]] = apollo_mnl(mnl_settings, functionality)
  
  ### Likelihood of the whole model
  P = apollo_combineModels(P, apollo_inputs, functionality)
  
  ### Take product across observation for same individual
  P = apollo_panelProd(P, apollo_inputs, functionality)
  
  ### Average across inter-individual draws
  P = apollo_avgInterDraws(P, apollo_inputs, functionality)
  
  ### Prepare and return outputs of function
  P = apollo_prepareProb(P, apollo_inputs, functionality)
  return(P)
}
Thank you for your advice!

Cheers,
Richard

Re: Calculating logsum from the estimated model

Posted: 11 Oct 2021, 15:37
by dpalma
Hi Richard,

Following de Jong et al. (2007) (section 3.2.2), to calculate the logsum with random heterorgeneity you need to average out any random heterogenity in the logsum the same way you do with the probability. The only difference is that you want the logsum at the observation level, as opposed to the individual level (as you do with the probabilities). Also, if you have income effects you need to convert the logsum to monetary units before averaging.

The code below should calculate the logsum at the observation level, assuming you don't have any income effects in the utility. I am not sure it will work, as I can't try it myself. The idea is to basically do the same things you do inside apollo_probabilities when calculating the likelihood, but instead of calculating the MNL probability, you calculate its logSum. Also, you don't want to multiply all observations from the same individual, so I removed the call to apollo_panelProd, and used functionality="prediction" to avoid errors popping out due to not calling it. I also dropped the likelihood of the indicators as they are not needed to calculate the logsum.

Please try it and let us know if there are any further issues.

Cheers
David

Code: Select all

### Attach inputs and detach after function exit
  apollo_attach(apollo_beta, apollo_inputs)
  on.exit(apollo_detach(apollo_beta, apollo_inputs))
  
  ### Creating interactions with latent variables
  tt_1_val   = tt_1+LV_tt_1_sig*LV_Tr_perc
  tt_2_val   = tt_2+LV_tt_2_sig*LV_Tr_perc
  tc_val     = tc+LV_tc_sig*LV_Tr_perc
  cont_2_val = cont_2+LV_cont_2_sig*LV_Tr_perc
  cont_3_val = cont_3+LV_cont_3_sig*LV_Tr_perc
  cont_4_val = cont_4+LV_cont_4_sig*LV_Tr_perc
  alig_2_val = alig_2+LV_alig_2_sig*LV_Tr_perc
  alig_3_val = alig_3+LV_alig_3_sig*LV_Tr_perc
  tt_t_val   = tt_t+LV_tt_t_sig*LV_Tr_perc
  tc_t_val   = tc_t+LV_tc_t_sig*LV_Tr_perc
  
  ### List of utilities: these must use the same names as in mnl_settings, order is irrelevant
  V = list()
  V[['MMB']]  = asc_MMB  + tt_1_val*mm.tt1 + tt_2_val*mm.tt2 + tc_val*mm.tc + cont_1*(mm.cont==0) + cont_2_val*(mm.cont==1) + cont_3_val*(mm.cont==2) +
    cont_4_val*(mm.cont==3) + alig_1*(mm.aligm==0) + alig_2_val*(mm.aligm==1) + alig_3_val*(mm.aligm==2) 
  V[['MBUS']]  = asc_MB  + tt_1_val*mb.tt1 + tt_2_val*mb.tt2 + tc_val*mb.tc + cont_1*(mb.cont==0) + cont_2_val*(mb.cont==1) + cont_3_val*(mb.cont==2) +
    cont_4_val*(mb.cont==3) + alig_1*(mb.aligb==0) + alig_2_val*(mb.aligb==1) + alig_3_val*(mb.aligb==2) 
  V[['TAXI']]  = asc_TX + tt_t_val*tx.ttt + tc_t_val*tx.ttc 
  V[['NO']] = asc_NO + att*mm.attr + ratt_1*(mm.rattr==0) + ratt_2*(mm.rattr==1) + ratt_3*(mm.rattr==2) + ratt_4*(mm.rattr==3)
  
  logSum <- log(Reduce('+', lapply(V, exp)))
  
  ### Average across inter-individual draws
  logSum = apollo_avgInterDraws(logSum, apollo_inputs, functionality='prediction')
  
  apollo_detach()

Re: Calculating logsum from the estimated model

Posted: 14 Oct 2021, 12:50
by hrrichard
Hi David,

Thank you so much for your help! The code worked well, and I was able to derive the logsum.

Cheers,
Richard
All times are UTC+01:00
Page 1 of 1

Powered by phpBB® Forum Software © phpBB Limited