ApolloChoiceModelling forum

Hi David and Stephane!

Unfortunately, I have a problem with interpreting my ICLV model results.

I would like to estimate the following structure: (1) MNL model, (2) MNL model with latent variable, (3) Mixed logit model, (4) Mixed logit model with latent variable.

I have one latent variable, and I use 17 indicator with ordered logit structure.

I used the following interaction term to estimate latent variable effect in choice model:

b_indication_value = b_indication + b_indication_shift_lambda * LV
(the indication is one of my attribute with levels yes or no)

In simple multinomial logit model, I get the following parameter estimate for indication attribute:

Estimate s.e. t.rat.(0) p(1-sided) Rob.s.e. Rob.t.rat.(0) p(1-sided)
b_indication 0.88030 0.034108 25.8090 0.0000 0.055539 15.8501 0.0000

In hybrid context, I get the following parameters:

Estimate s.e. t.rat.(0) p(1-sided) Rob.s.e. Rob.t.rat.(0) p(1-sided)
b_indication 1.02070 0.067759 15.0636 0.000000 0.079906 12.7738 0.000000
b_indication_shift_lambda 0.85525 0.041877 20.4227 0.000000 0.074643 11.4578 0.000000

In simple mixed logit model, I get the following parameters:

Estimate s.e. t.rat.(0) p(1-sided) Rob.s.e. Rob.t.rat.(0) p(1-sided)
b_indication_mu 1.40060 0.10932 12.812 0.00000 0.11695 11.9764 0.00000
b_indication_sig 2.43056 0.11569 21.010 0.00000 0.12538 19.3852 0.00000


In hybrid context, I get the following parameters:

Estimate s.e. t.rat.(0) p(1-sided) Rob.s.e. Rob.t.rat.(0) p(1-sided)
b_indication_mu 1.29933 0.12177 10.6708 0.000000 0.16258 7.9918 6.661e-16
b_indication_sig -1.92295 0.09302 -20.6731 0.000000 0.10091 -19.0567 0.000000
b_indication_shift_lambda 1.41567 0.09522 14.8667 0.000000 0.13035 10.8609 0.000000

In this case, how can I interpret the "b_indication_shift_lambda" parameter? I thought that, the higher level of the latent variable cause higher level of utility for existence of indication (relative to no indication), but how can I interpret the growth of latent variable?
What is the unit of latent variable?

Thank you very much for your help!

Best regards,
Peter

Peter

the interpretation that "the higher level of the latent variable cause higher level of utility for existence of indication" is correct.

The latent variable does not have a measurement scale that you can use for interpretation. There are two things you can do. First, you can look at the combined heterogeneity that you now have in b_indication as a result of the latent variable, and recognise that in your MRS calculation. Second, if you have a structural equation, then you can look at how much MRS change as a result of the covariates that are in the structural equation

Stephane

Hi Stephane!

Thank you for the quick and useful answer!

For solution one:
So, for example if I use WTP-space transformation (V_indication=b_indication+b_indication_shift_lambda*LV -> U=b_cost*(V_indication*Indication_alti +....), is it enough to interpret "b_indication" and ignore the "b_indication_shift_lambda" parameter as well as highlight the combined heterogeneity what "b_indication" contains?

For solution two:
Yes, I have the following structural equation: gamma_city*City + gamma_highedu*HigherEducation + gamma_above60*Above60. I got the following parameter estimates: gamma_city: -0.30, gamma_highedu: -0.20, gamma_above60: 0.24
So in this case, should I use for example the covariate of gamma_city (-0.30) to calculate and interpret WTP by the following: b_indication_shift_lambda*-0.30/b_cost, where the result means the WTP of city resident relative to non-city resident for indication (of course in addition to all other factors are considered fixed in the structural equation).

Would this interpretation be correct?

Thank you,
Peter

Peter

if you can show your entire code, then I can help you.

Thanks

Stephane

Hi Stephane,

I copy the code used to estimate the ICLV model.

rm(list = ls())

library(apollo)

apollo_initialise()

apollo_control = list(
modelName ="MLM_based_ICLV",
modelDescr ="MLM_based_ICLV",
indivID ="ID",
mixing = TRUE,
nCores = 4
)

database = read.csv("database.csv",header=TRUE)

apollo_beta = c(asc_alt1 = 0,
asc_alt2 = 0,
asc_alt3 = 0,
asc_nc = 0,
b_brand_mu = 0,
b_brand_sig = 0,
b_indication_mu = 0,
b_indication_sig = 0,
b_method_mu = 0,
b_method_sig = 0,
b_price_mu = -3,
b_price_sig = 0,
b_indication_shift_lambda = 1,
gamma_city = 0,
gamma_highedu = 0,
gamma_above60 = 0,
zeta_et_1 = 1,
zeta_et_2 = 1,
zeta_et_3 = 1,
zeta_et_4 = 1,
zeta_et_5 = 1,
zeta_et_6 = 1,
zeta_et_7 = 1,
zeta_et_8 = 1,
zeta_et_9 = 1,
zeta_et_10 = 1,
zeta_et_11 = 1,
zeta_et_12 = 1,
zeta_et_13 = 1,
zeta_et_14 = 1,
zeta_et_15 = 1,
zeta_et_16 = 1,
zeta_et_17 = 1,
tau_et_1_1 =-2,
tau_et_1_2 =-1,
tau_et_1_3 = 1,
tau_et_1_4 = 2,
tau_et_2_1 =-2,
tau_et_2_2 =-1,
tau_et_2_3 = 1,
tau_et_2_4 = 2,
tau_et_3_1 =-2,
tau_et_3_2 =-1,
tau_et_3_3 = 1,
tau_et_3_4 = 2,
tau_et_4_1 =-2,
tau_et_4_2 =-1,
tau_et_4_3 = 1,
tau_et_4_4 = 2,
tau_et_5_1 =-2,
tau_et_5_2 =-1,
tau_et_5_3 = 1,
tau_et_5_4 = 2,
tau_et_6_1 =-2,
tau_et_6_2 =-1,
tau_et_6_3 = 1,
tau_et_6_4 = 2,
tau_et_7_1 =-2,
tau_et_7_2 =-1,
tau_et_7_3 = 1,
tau_et_7_4 = 2,
tau_et_8_1 =-2,
tau_et_8_2 =-1,
tau_et_8_3 = 1,
tau_et_8_4 = 2,
tau_et_9_1 =-2,
tau_et_9_2 =-1,
tau_et_9_3 = 1,
tau_et_9_4 = 2,
tau_et_10_1 =-2,
tau_et_10_2 =-1,
tau_et_10_3 = 1,
tau_et_10_4 = 2,
tau_et_11_1 =-2,
tau_et_11_2 =-1,
tau_et_11_3 = 1,
tau_et_11_4 = 2,
tau_et_12_1 =-2,
tau_et_12_2 =-1,
tau_et_12_3 = 1,
tau_et_12_4 = 2,
tau_et_13_1 =-2,
tau_et_13_2 =-1,
tau_et_13_3 = 1,
tau_et_13_4 = 2,
tau_et_14_1 =-2,
tau_et_14_2 =-1,
tau_et_14_3 = 1,
tau_et_14_4 = 2,
tau_et_15_1 =-2,
tau_et_15_2 =-1,
tau_et_15_3 = 1,
tau_et_15_4 = 2,
tau_et_16_1 =-2,
tau_et_16_2 =-1,
tau_et_16_3 = 1,
tau_et_16_4 = 2,
tau_et_17_1 =-2,
tau_et_17_2 =-1,
tau_et_17_3 = 1,
tau_et_17_4 = 2)

apollo_fixed = c("asc_alt1")

apollo_draws = list(
interDrawsType="mlhs",
interNDraws=500,
interUnifDraws=c(),
interNormDraws=c("eta","draws_brand","draws_indication","draws_method","draws_price")
)

apollo_randCoeff=function(apollo_beta, apollo_inputs){
randcoeff = list()

randcoeff[["LV"]] = gamma_city*city + gamma_highedu*highedu + gamma_above60*above60 + eta
randcoeff[["b_brand"]] = ( b_brand_mu + b_brand_sig * draws_brand )
randcoeff[["b_indication"]] = ( b_indication_mu + b_indication_sig * draws_indication )
randcoeff[["b_method"]] = ( b_method_mu + b_method_sig * draws_method )
randcoeff[["b_price"]] = -exp( b_price_mu + b_price_sig * draws_price )

return(randcoeff)
}

apollo_inputs = apollo_validateInputs()

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

apollo_attach(apollo_beta, apollo_inputs)
on.exit(apollo_detach(apollo_beta, apollo_inputs))

P = list()

b_indication_value = b_indication + b_indication_shift_lambda * LV


op_settingse1 = list(outcomeOrdered = et_1,
V = zeta_et_1*LV_e,
tau = c(tau_et_1_1, tau_et_1_2, tau_et_1_3, tau_et_1_4),
rows = (situation==1),
componentName = "indic_et_1")
op_settingse2 = list(outcomeOrdered = et_2,
V = zeta_et_2*LV_e,
tau = c(tau_et_2_1, tau_et_2_2, tau_et_2_3, tau_et_2_4),
rows = (situation==1),
componentName = "indic_et_2")
op_settingse3 = list(outcomeOrdered = et_3,
V = zeta_et_3*LV_e,
tau = c(tau_et_3_1, tau_et_3_2, tau_et_3_3, tau_et_3_4),
rows = (situation==1),
componentName = "indic_et_3")
op_settingse4 = list(outcomeOrdered = et_4,
V = zeta_et_4*LV_e,
tau = c(tau_et_4_1, tau_et_4_2, tau_et_4_3, tau_et_4_4),
rows = (situation==1),
componentName = "indic_et_4")
op_settingse5 = list(outcomeOrdered = et_5,
V = zeta_et_5*LV_e,
tau = c(tau_et_5_1, tau_et_5_2, tau_et_5_3, tau_et_5_4),
rows = (situation==1),
componentName = "indic_et_5")
op_settingse6 = list(outcomeOrdered = et_6,
V = zeta_et_6*LV_e,
tau = c(tau_et_6_1, tau_et_6_2, tau_et_6_3, tau_et_6_4),
rows = (situation==1),
componentName = "indic_et_6")
op_settingse7 = list(outcomeOrdered = et_7,
V = zeta_et_7*LV_e,
tau = c(tau_et_7_1, tau_et_7_2, tau_et_7_3, tau_et_7_4),
rows = (situation==1),
componentName = "indic_et_7")
op_settingse8 = list(outcomeOrdered = et_8,
V = zeta_et_8*LV_e,
tau = c(tau_et_8_1, tau_et_8_2, tau_et_8_3, tau_et_8_4),
rows = (situation==1),
componentName = "indic_et_8")
op_settingse9 = list(outcomeOrdered = et_9,
V = zeta_et_9*LV_e,
tau = c(tau_et_9_1, tau_et_9_2, tau_et_9_3, tau_et_9_4),
rows = (situation==1),
componentName = "indic_et_9")
op_settingse10 = list(outcomeOrdered = et_10,
V = zeta_et_10*LV_e,
tau = c(tau_et_10_1, tau_et_10_2, tau_et_10_3, tau_et_10_4),
rows = (situation==1),
componentName = "indic_et_10")
op_settingse11 = list(outcomeOrdered = et_11,
V = zeta_et_11*LV_e,
tau = c(tau_et_11_1, tau_et_11_2, tau_et_11_3, tau_et_11_4),
rows = (situation==1),
componentName = "indic_et_11")
op_settingse12 = list(outcomeOrdered = et_12,
V = zeta_et_12*LV_e,
tau = c(tau_et_12_1, tau_et_12_2, tau_et_12_3, tau_et_12_4),
rows = (situation==1),
componentName = "indic_et_12")
op_settingse13 = list(outcomeOrdered = et_13,
V = zeta_et_13*LV_e,
tau = c(tau_et_13_1, tau_et_13_2, tau_et_13_3, tau_et_13_4),
rows = (situation==1),
componentName = "indic_et_13")
op_settingse14 = list(outcomeOrdered = et_14,
V = zeta_et_14*LV_e,
tau = c(tau_et_14_1, tau_et_14_2, tau_et_14_3, tau_et_14_4),
rows = (situation==1),
componentName = "indic_et_14")
op_settingse15 = list(outcomeOrdered = et_15,
V = zeta_et_15*LV_e,
tau = c(tau_et_15_1, tau_et_15_2, tau_et_15_3, tau_et_15_4),
rows = (situation==1),
componentName = "indic_et_15")
op_settingse16 = list(outcomeOrdered = et_16,
V = zeta_et_16*LV_e,
tau = c(tau_et_16_1, tau_et_16_2, tau_et_16_3, tau_et_16_4),
rows = (situation==1),
componentName = "indic_et_16")
op_settingse17 = list(outcomeOrdered = et_17,
V = zeta_et_17*LV_e,
tau = c(tau_et_17_1, tau_et_17_2, tau_et_17_3, tau_et_17_4),
rows = (situation==1),
componentName = "indic_et_17")

P[["indic_et_1"]] = apollo_op(op_settingse1, functionality)
P[["indic_et_2"]] = apollo_op(op_settingse2, functionality)
P[["indic_et_3"]] = apollo_op(op_settingse3, functionality)
P[["indic_et_4"]] = apollo_op(op_settingse4, functionality)
P[["indic_et_5"]] = apollo_op(op_settingse5, functionality)
P[["indic_et_6"]] = apollo_op(op_settingse6, functionality)
P[["indic_et_7"]] = apollo_op(op_settingse7, functionality)
P[["indic_et_8"]] = apollo_op(op_settingse8, functionality)
P[["indic_et_9"]] = apollo_op(op_settingse9, functionality)
P[["indic_et_10"]] = apollo_op(op_settingse10, functionality)
P[["indic_et_11"]] = apollo_op(op_settingse11, functionality)
P[["indic_et_12"]] = apollo_op(op_settingse12, functionality)
P[["indic_et_13"]] = apollo_op(op_settingse13, functionality)
P[["indic_et_14"]] = apollo_op(op_settingse14, functionality)
P[["indic_et_15"]] = apollo_op(op_settingse15, functionality)
P[["indic_et_16"]] = apollo_op(op_settingse16, functionality)
P[["indic_et_17"]] = apollo_op(op_settingse17, functionality)



V = list()
V[['alt1']] = asc_alt1 + b_brand * (brand_alt1 == 1) + b_indication_value * (indication_alt1 == 1) + b_method * (method_alt1 == 1) + b_price * price_alt1
V[['alt2']] = asc_alt2 + b_brand * (brand_alt2 == 1) + b_indication_value * (indication_alt2 == 1) + b_method * (method_alt2 == 1) + b_price * price_alt2
V[['alt3']] = asc_alt3 + b_brand * (brand_alt3 == 1) + b_indication_value * (indication_alt3 == 1) + b_method * (method_alt3 == 1) + b_price * price_alt3
V[['alt4']] = asc_nc


mnl_settings = list(
alternatives = c(alt1=1, alt2=2, alt3=3, alt4=4),
avail = list(alt1=av_alt1, alt2=av_alt2, alt3=av_alt3, alt4=av_alt4),
choiceVar = choice,
V = V,
componentName= "choice"
)

P[["choice"]] = apollo_mnl(mnl_settings, functionality)

P = apollo_combineModels(P, apollo_inputs, functionality)

P = apollo_panelProd(P, apollo_inputs, functionality)

P = apollo_avgInterDraws(P, apollo_inputs, functionality)

P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}


model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)

L = apollo_probabilities(apollo_beta, apollo_inputs, functionality="estimate")

L[1:30]

IDs=unique(database$ID)
failures=IDs[is.na(log(L))]
database$choice[database$ID%in%failures]
failures=IDs[is.na(log(L))]
database$choice[database_ID%in%failures]
database$et_1[database$ID%in%failures&database$situation==1]
database$et_2[database$ID%in%failures&database$situation==1]
database$et_3[database$ID%in%failures&database$situation==1]
database$et_4[database$ID%in%failures&database$situation==1]
database$et_5[database$ID%in%failures&database$situation==1]
database$et_6[database$ID%in%failures&database$situation==1]
database$et_7[database$ID%in%failures&database$situation==1]
database$et_8[database$ID%in%failures&database$situation==1]
database$et_9[database$ID%in%failures&database$situation==1]
database$et_10[database$ID%in%failures&database$situation==1]
database$et_11[database$ID%in%failures&database$situation==1]
database$et_12[database$ID%in%failures&database$situation==1]
database$et_13[database$ID%in%failures&database$situation==1]
database$et_14[database$ID%in%failures&database$situation==1]
database$et_15[database$ID%in%failures&database$situation==1]
database$et_16[database$ID%in%failures&database$situation==1]
database$et_17[database$ID%in%failures&database$situation==1]

modelOutput_settings = list()
modelOutput_settings$printPVal=TRUE
apollo_modelOutput(model,modelOutput_settings)

saveOutput_settings = list()
saveOutput_settings$printPVal=TRUE
apollo_saveOutput(model,modelOutput_settings)


Peter

Hi

you don't need WTP space for this, but you also can't ignore the b_indication_shift_lambda in computing WTP.

In particular, in your model, you have that that marginal utility of (indication_alt1 == 1) is equal to b_indication_value

This in turn is equal to:
which is:
And then you have that: WTP= b_indication_value/b_price

what you could do after estimation is to use this code:
This is the WTP that takes into account the role of the LV in your model

Stephane

Hi Stephane,

thank you very much for your help!

I just tried to use the WTP-space transformation (it worked perfectly for the other models) in order to reach a consistent approach for all models analyzed (to avoid future criticism from reviewer).

Based on your suggestion, in this case I should use the following code:

##Interaction term:
V_indication_value = V_indication + b_indication_shift_lambda * LV

##Utility function:
V[['alt1']] = asc_alt1 + b_price*(V_brand * (brand_alt1 == 1) + V_indication_value * (indication_alt1 == 1) + V_method * (method_alt1 == 1) + price_alt1)

##Unconditionals
unconditionals=apollo_unconditionals(model,apollo_probabilities,apollo_inputs)
WTP_indication=unconditionals[["V_indication"]]+model$estimate["b_indication_shift_lambda"]*unconditionals[["LV"]]


Could this be correct?

Thank you and sorry for the many questions,
Peter

Yes, that seems okay

Thank you so much for your help Stephane!

Dear Stephane/David,

I have a question related to this: When you look at how a latent parameter affects WTP in WTP-space using the unconditional estimates in Apollo, is there a simple way to estimate confidence interval of WTP for a given shift in the latent variable?

Let me try to explain my question more thoroughly:

For example:

Say you have the following interaction term between WTP for an attribute and a latent variable (LV):

##Interaction term:
V_att_value = V_att + b_indication_shift_LV * LV

Then, as written previosly in this thread, I can use the unconditional to look at how WTP changes over different levels of LV by using the following command:

##Unconditionals
unconditionals=apollo_unconditionals(model,apollo_probabilities,apollo_inputs)
WTP_LV=unconditionals[["V_indication"]]+model$estimate["b_indication_shift_LV"]*unconditionals[["LV"]]

Say that LV is normalised to have mean of 0 and standard deviation of 1. Then, mean WTP if LV decreases by 1 standard deviation becomes (if I am not mistaken):

mean(WTP_LV[unconditionals[["LV"]]<=-1])

So my question is whether there is an integrated way in Apollo to get a confidence interval for this value, or perhaps another way to get the confidence interval? (Sorry if this is a silly question).

Thanks
Per