Latent variable inside delta function
Posted: 18 Feb 2023, 11:41
Hello, Professor.
I would like to run a mediation analysis with a latent variable as a mediator. Consider, for example, the following model:
Independent variable (X): treatment variable in which subjects receive or not a benefit;
Mediator (M): a latent variable (ex., car lover) using ordered measurement model for indicators estimated with Apollo's HCM framework;
Dependent variable (Y): choice between mode A or B.
Other variables: price, time etc.
In SEM packages such as Lavaan in R, we have to estimate extra parameters to demonstrate that M mediates the effect of X on Y. For example, I could specify the following model in Lavaan to test a mediation:
model <- '
# direct effect
Y ~ c*X
# mediator
M ~ a*X
Y ~ b*M
# indirect effect (a*b)
ab := a*b
# total effect
total := c + (a*b)
'
1- Do you see any problem in estimating indirect effects (a*b) with apollo_deltaMethod? Or does it sound good?
2- Is it possible to run the apollo_deltaMethod with a latent variable? I mean, my mediator is a latent variable, so I would need to use it inside the apollo_deltaMethod. However, I received the following error when tried:
Please find below the code I am working on.
Thank you very much for your time and help.
Best, Rafael.
I would like to run a mediation analysis with a latent variable as a mediator. Consider, for example, the following model:
Independent variable (X): treatment variable in which subjects receive or not a benefit;
Mediator (M): a latent variable (ex., car lover) using ordered measurement model for indicators estimated with Apollo's HCM framework;
Dependent variable (Y): choice between mode A or B.
Other variables: price, time etc.
In SEM packages such as Lavaan in R, we have to estimate extra parameters to demonstrate that M mediates the effect of X on Y. For example, I could specify the following model in Lavaan to test a mediation:
model <- '
# direct effect
Y ~ c*X
# mediator
M ~ a*X
Y ~ b*M
# indirect effect (a*b)
ab := a*b
# total effect
total := c + (a*b)
'
Since Apollo has a function to calculate standard errors for derived measures using Delta method either, I think I could do this mediation analysis there. Thus, I would have two questions:. The := operator ‘defines’ new parameters which take on values that are an arbitrary function of the original model parameters. [...] By default, the standard errors for these defined parameters are computed by using the so-called Delta method. As with other models, bootstrap standard errors can be requested...
1- Do you see any problem in estimating indirect effects (a*b) with apollo_deltaMethod? Or does it sound good?
2- Is it possible to run the apollo_deltaMethod with a latent variable? I mean, my mediator is a latent variable, so I would need to use it inside the apollo_deltaMethod. However, I received the following error when tried:
Code: Select all
deltaMethod_settings = list(
expression = c(
effect = "LV_sec*gamma_sec_treat"
)
)
Error in apollo_deltaMethod(model, deltaMethod_settings) :
The expression effect includes variables that are not parameters from the model: LV_sec.Thank you very much for your time and help.
Best, Rafael.
Code: Select all
# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS ####
# ################################################################# #
### Initialise code
apollo_initialise()
### Set core controls
apollo_control = list(
modelName = "Hybrid_with_OL",
modelDescr = "Hybrid choice model with SP data, using ordered measurement model for indicators",
indivID = "responseid",
mixing = TRUE,
nCores = 30,
outputDirectory = "output"
)
# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #
### Loading data from package
### if data is to be loaded from a file (e.g. called data.csv),
### the code would be: database = read.csv("data.csv",header=TRUE)
database = base
### for data dictionary, use ?apollo_drugChoiceData
table(base$choice, useNA = "ifany")
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(# Choice parameters
asc_sec = -1,
b_price_angel_ref = -3,
b_price_angel_lv = -3,
b_price_sec_ref = -3,
b_price_sec_lv = -3,
b_time = 0,
gamma_sec_treat = -1,
gamma_ang_treat = 1,
gamma_sec_inc17 = 1,
gamma_sec_inc29 = 1,
gamma_sec_inc41 = 1,
gamma_ang_inc17 = 1,
gamma_ang_inc29 = 1,
gamma_ang_inc41 = 1,
gamma_sec_depen = 1,
gamma_ang_depen = 1,
# gamma_sec_life = 1,
# gamma_ang_life = 1,
# Measurement equations parameters
# sec
zeta_secTrust1 = 1, zeta_secTrust2 = 1, zeta_secTrust3 = 1,
tau_secTrust1_1 =-3, tau_secTrust1_2 =-2, tau_secTrust1_3 =-1, tau_secTrust1_4 =1, tau_secTrust1_5 =2, tau_secTrust1_6 =3,
tau_secTrust2_1 =-3, tau_secTrust2_2 =-2, tau_secTrust2_3 =-1, tau_secTrust2_4 =1, tau_secTrust2_5 =2, tau_secTrust2_6 =3,
tau_secTrust3_1 =-3, tau_secTrust3_2 =-2, tau_secTrust3_3 =-1, tau_secTrust3_4 =1, tau_secTrust3_5 =2, tau_secTrust3_6 =3,
# angel
zeta_angTrust1 = 1, zeta_angTrust2 = 1, zeta_angTrust3 = 1,
tau_angTrust1_1 =-3, tau_angTrust1_2 =-2, tau_angTrust1_3 =-1, tau_angTrust1_4 =1, tau_angTrust1_5 =2, tau_angTrust1_6 =3,
tau_angTrust2_1 =-3, tau_angTrust2_2 =-2, tau_angTrust2_3 =-1, tau_angTrust2_4 =1, tau_angTrust2_5 =2, tau_angTrust2_6 =3,
tau_angTrust3_1 =-3, tau_angTrust3_2 =-2, tau_angTrust3_3 =-1, tau_angTrust3_4 =1, tau_angTrust3_5 =2, tau_angTrust3_6 =3)
### Vector with names (in quotes) of parameters to be kept fixed at their starting value in apollo_beta, use apollo_beta_fixed = c() if none
# fix one zeta parameter to 1 per latent variable, but then estimate the standard deviation of the latent variable
# http://www.apollochoicemodelling.com/forum/viewtopic.php?f=11&t=23&p=433&hilit=parameter+latent+variable#p433
apollo_fixed = c("zeta_secTrust1", "zeta_angTrust1")
### Read in starting values for at least some parameters from existing model output file
# apollo_beta = apollo_readBeta(apollo_beta, apollo_fixed, "Hybrid_with_OL", overwriteFixed=FALSE)
# ################################################################# #
#### DEFINE RANDOM COMPONENTS ####
# ################################################################# #
### Set parameters for generating draws
apollo_draws = list(
interDrawsType = "MLHS", #Halton Sobol
interNDraws = 500,
interUnifDraws=c(),
interNormDraws = c("eta_sec", "eta_ang"),
intraDrawsType="",
intraNDraws=0,
intraUnifDraws=c(),
intraNormDraws=c()
)
### Create random parameters
apollo_randCoeff=function(apollo_beta, apollo_inputs){
randcoeff = list()
randcoeff[["LV_sec"]] = eta_sec + gamma_sec_treat*treat + gamma_sec_inc17*(baselineprice==17) + gamma_sec_inc29*(baselineprice==29) + gamma_sec_inc41*(baselineprice==41) + gamma_sec_depen*(depenscaled>0) #+ gamma_sec_life*(life_ins==1)
randcoeff[["LV_ang"]] = eta_ang + gamma_ang_treat*treat + gamma_ang_inc17*(baselineprice==17) + gamma_ang_inc29*(baselineprice==29) + gamma_ang_inc41*(baselineprice==41) + gamma_ang_depen*(depenscaled>0) #+ gamma_ang_life*(life_ins==1)
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()
### Likelihood of indicators
ol_settings1 = list(outcomeOrdered = sec_tru1,
V = zeta_secTrust1*LV_sec,
tau = list(tau_secTrust1_1, tau_secTrust1_2, tau_secTrust1_3, tau_secTrust1_4, tau_secTrust1_5, tau_secTrust1_6),
componentName ='Indicator_1')
ol_settings2 = list(outcomeOrdered = sec_tru2,
V = zeta_secTrust2*LV_sec,
tau = list(tau_secTrust2_1, tau_secTrust2_2, tau_secTrust2_3, tau_secTrust2_4, tau_secTrust2_5, tau_secTrust2_6),
componentName ='Indicator_2')
ol_settings3 = list(outcomeOrdered = sec_tru3,
V = zeta_secTrust3*LV_sec,
tau = list(tau_secTrust3_1, tau_secTrust3_2, tau_secTrust3_3, tau_secTrust3_4, tau_secTrust3_5, tau_secTrust3_6),
componentName ='Indicator_3')
ol_settings4 = list(outcomeOrdered = ang_tru1,
V = zeta_angTrust1*LV_ang,
tau = list(tau_angTrust1_1, tau_angTrust1_2, tau_angTrust1_3, tau_angTrust1_4, tau_angTrust1_5, tau_angTrust1_6),
componentName ='Indicator_4')
ol_settings5 = list(outcomeOrdered = ang_tru2,
V = zeta_angTrust2*LV_ang,
tau = list(tau_angTrust2_1, tau_angTrust2_2, tau_angTrust2_3, tau_angTrust2_4, tau_angTrust2_5, tau_angTrust2_6),
componentName ='Indicator_5')
ol_settings6 = list(outcomeOrdered = ang_tru3,
V = zeta_angTrust3*LV_ang,
tau = list(tau_angTrust3_1, tau_angTrust3_2, tau_angTrust3_3, tau_angTrust3_4, tau_angTrust3_5, tau_angTrust3_6),
componentName ='Indicator_6')
P[["indic_secTrust1"]] = apollo_ol(ol_settings1, functionality)
P[["indic_secTrust2"]] = apollo_ol(ol_settings2, functionality)
P[["indic_secTrust3"]] = apollo_ol(ol_settings3, functionality)
P[["indic_angTrust1"]] = apollo_ol(ol_settings4, functionality)
P[["indic_angTrust2"]] = apollo_ol(ol_settings5, functionality)
P[["indic_angTrust3"]] = apollo_ol(ol_settings6, functionality)
### Likelihood of choices
### List of utilities: these must use the same names as in mnl_settings, order is irrelevant
V = list()
asc_sec_value = asc_sec #+ asc_sec_shift*treat
b_price_angel_value = b_price_angel_ref*exp(b_price_angel_lv*LV_ang)
b_price_sec_value = b_price_sec_ref *exp(b_price_sec_lv *LV_sec)
V[["angel"]] = b_price_angel_value*(ang_price/10) + b_time*ang_time
V[["sec"]] = asc_sec_value + b_price_sec_value *(sec_price/10) + b_time*sec_time
### Define settings for MNL model component
mnl_settings = list(
alternatives = c(angel=1, sec=0),
avail = list(angel=1, sec=1),
choiceVar = choice,
utilities = V,
componentName = "choice"
)
### Compute probabilities for MNL model component
P[["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)
}
# ################################################################# #
#### MODEL ESTIMATION ####
# ################################################################# #
### Optional: calculate LL before model estimation
# apollo_llCalc(apollo_beta, apollo_probabilities, apollo_inputs)
### Estimate model
model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)
apollo_modelOutput(model)