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.

Logistic distribution of a random parameter

Ask questions about model specifications. Ideally include a mathematical explanation of your proposed model.
Post Reply
bokapatsila
Posts: 21
Joined: 28 Jul 2021, 02:41

Logistic distribution of a random parameter

Post by bokapatsila »

Good morning Apollo forum moderators. Can you please help me to reflect in code the assumption that a random error term follows a logistic distribution and not a normal one? For example, how that can be coded in example #24 for "eta" (meaning that eta follows a logistic distribution). Thank you in advance!
stephanehess
Site Admin
Posts: 974
Joined: 24 Apr 2020, 16:29

Re: Logistic distribution of a random parameter

Post by stephanehess »

Hi

you would need to start with uniform draws in apollo_draws instead of Normal draws, and then use the appropriate formula to turn these into logistic draws inside apollo_randCoeff.

Something like this:

Code: Select all

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="halton", 
  interNDraws=100,          
  interUnifDraws=c("eta")
)

### Create random parameters
apollo_randCoeff=function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["LV"]] = gamma_reg_user*regular_user + gamma_university*university_educated + gamma_age_50*over_50 + 2*b*atanh(2*eta-1)+m
  
  return(randcoeff)
}
where b and m are the parameters of the logistic distribution, and you need to think about normalisation

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
bokapatsila
Posts: 21
Joined: 28 Jul 2021, 02:41

Re: Logistic distribution of a random parameter

Post by bokapatsila »

Thank you!
bokapatsila
Posts: 21
Joined: 28 Jul 2021, 02:41

Re: Logistic distribution of a random parameter

Post by bokapatsila »

stephanehess wrote: 12 Nov 2021, 09:40 Hi

you would need to start with uniform draws in apollo_draws instead of Normal draws, and then use the appropriate formula to turn these into logistic draws inside apollo_randCoeff.

Something like this:

Code: Select all

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="halton", 
  interNDraws=100,          
  interUnifDraws=c("eta")
)

### Create random parameters
apollo_randCoeff=function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["LV"]] = gamma_reg_user*regular_user + gamma_university*university_educated + gamma_age_50*over_50 + 2*b*atanh(2*eta-1)+m
  
  return(randcoeff)
}
where b and m are the parameters of the logistic distribution, and you need to think about normalisation

Stephane
Dear Stephane,

I have another question regarding the example that you kindly provided. In the code that you've shared, is b a scale parameter, and m a location parameter? When you speak about normalisation, do you mean that b and m should be kept fixed in apollo_fixed?

In the model below I'm trying to classify my respondents into latent classes using a threshold for a latent variable. While the model seems to be running without errors, I'm not fully certain that my code is correct. Can you please take a look at it?

Here's the code:

Code: Select all

# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS                       ####
# ################################################################# #

### Clear memory
rm(list = ls())

### Load libraries
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="LatentClass_LV1_N1.4_imp",
  modelDescr ="ILVLC with 1 LV",
  indivID    ="UNID",
  panelData = FALSE,
  mixing     = TRUE,
#workInLogs = TRUE,
  nCores     = 3)

# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS                     ####
# ################################################################# #

setwd("D:/Data")

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

# ################################################################# #
#### DEFINE MODEL PARAMETERS                                     ####
# ################################################################# #

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(zeta_aconC = 1, 
              zeta_sconC = 1, 
              zeta_bothC = 1,
              zeta_seatC = 1,
              zeta_offpeakC = 1, 
              zeta_altC = 1, 
              tau_acon_1C =-2, 
              tau_acon_2C =-1, 
              tau_acon_3C = 1, 
              tau_acon_4C = 2,
              tau_scon_1C =-2, 
              tau_scon_2C =-1, 
              tau_scon_3C = 1, 
              tau_scon_4C = 2,
              tau_both_1C =-2, 
              tau_both_2C =-1, 
              tau_both_3C = 1, 
              tau_both_4C = 2,
              tau_seat_1C =-2, 
              tau_seat_2C =-1, 
              tau_seat_3C = 1, 
              tau_seat_4C = 2,
              tau_offpeak_1C =-2, 
              tau_offpeak_2C =-1,
              tau_offpeak_3C = 1, 
              tau_offpeak_4C = 2,
              tau_alt_1C =-2, 
              tau_alt_2C =-1,
              tau_alt_3C = 1, 
              tau_alt_4C = 2,

              zeta_sconX_A = 1,
              zeta_sconX_B = -1,

              tau_scon_1X =-2, 
              tau_scon_2X =-1,
              tau_scon_3X = 1, 
              tau_scon_4X = 2,

              gamma_LV1_female  = 1, 
              gamma_LV1_fulltime  = 1,
              gamma_LV1_incomelow  = 1,
              gamma_LV1_age2534  = 1,
              gamma_LV1_edbach = 1,
              gamma_LV1_kids = 1,
              gamma_LV1_car = 1,
              gamma_LV1_wave2  = 1,
              
              b1 = 1,
              m1 = 0,
              thre1 = 1)
              

### 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
apollo_fixed = c("b1", "m1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="halton", 
  interNDraws=50,          
  interUnifDraws=c("eta1"),      
  interNormDraws=c(), 
  
  intraDrawsType='',
  intraNDraws=0,          
  intraUnifDraws=c(),     
  intraNormDraws=c()      
)

### Create random parameters
apollo_randCoeff=function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["LVC1"]] = gamma_LV1_female * FEMALE + gamma_LV1_wave2 * WAVE2 +
    gamma_LV1_fulltime * EMPLOY_FULL + gamma_LV1_incomelow * INCOME_LOW + gamma_LV1_age2534 * AGE_2534 + 
    gamma_LV1_edbach * EDUCATION_BACH + gamma_LV1_kids * HOUSEHOLD_CHILD_CLEAN_N + gamma_LV1_car * CAR_B +
    2*b1*atanh(2*eta1-1)+m1 #Error term follows a logistic distribution (to match the logit in allocation function)
  
  randcoeff[["LVX1"]] = gamma_LV1_female * FEMALE + gamma_LV1_wave2 * WAVE2 + 
    gamma_LV1_fulltime * EMPLOY_FULL + gamma_LV1_incomelow * INCOME_LOW + gamma_LV1_age2534 * AGE_2534 + 
    gamma_LV1_edbach * EDUCATION_BACH + gamma_LV1_kids * HOUSEHOLD_CHILD_CLEAN_N + gamma_LV1_car * CAR_B

  return(randcoeff)
}

# ################################################################# #
#### DEFINE LATENT CLASS COMPONENTS                              ####
# ################################################################# #

apollo_lcPars=function(apollo_beta, apollo_inputs){
  lcpars = list()
  
  ### Define lists of parameters for each class
  ###                                    classA   classB  ...
  lcpars[["zeta_sconX"]]  = list(zeta_sconX_A, zeta_sconX_B)
  
  ### Class allocation probabilities
  ### These are the probabilities of a binary logit model
  classAlloc_settings = list(
    V = list(A = LVC1-thre1,
             B = 0)
  )
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  
  return(lcpars)
}

# ################################################################# #
#### 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()
  PClass4 = list()
  
  ### Likelihood of indicators
  ol_settings1C = list(outcomeOrdered=AGREE_BOTH_CROWD_BC_R, 
                      V=zeta_bothC*LVC1, 
                      tau=list(tau_both_1C, tau_both_2C, tau_both_3C, tau_both_4C))
  ol_settings2C = list(outcomeOrdered=AGREE_CONCERN_BC_R, 
                      V=zeta_aconC*LVC1, 
                      tau=list(tau_acon_1C, tau_acon_2C, tau_acon_3C, tau_acon_4C))
  ol_settings3C = list(outcomeOrdered=AGREE_SEAT_BC_R, 
                      V=zeta_seatC*LVC1, 
                      tau=list(tau_seat_1C, tau_seat_2C, tau_seat_3C, tau_seat_4C))
  ol_settings4C = list(outcomeOrdered=STATE_CONCERNED_R, 
                      V=zeta_sconC*LVC1, 
                      tau=list(tau_scon_1C, tau_scon_2C, tau_scon_3C, tau_scon_4C))
  ol_settings5C = list(outcomeOrdered=AGREE_OFFPEAK_BC_R, 
                      V=zeta_offpeakC*LVC1, 
                      tau=list(tau_offpeak_1C, tau_offpeak_2C, tau_offpeak_3C, tau_offpeak_4C))
  ol_settings6C = list(outcomeOrdered=AGREE_ALT_BC_R, 
                      V=zeta_altC*LVC1, 
                      tau=list(tau_alt_1C, tau_alt_2C, tau_alt_3C, tau_alt_4C))
  
  P[["indic_bothC"]]     = apollo_ol(ol_settings1C, functionality)
  P[["indic_aconC"]]     = apollo_ol(ol_settings2C, functionality)
  P[["indic_seatC"]]      = apollo_ol(ol_settings3C, functionality)
  P[["indic_sconC"]]      = apollo_ol(ol_settings4C, functionality)
  P[["indic_offpeakC"]]     = apollo_ol(ol_settings5C, functionality)
  P[["indic_altC"]]      = apollo_ol(ol_settings6C, functionality)

  ### Categorical LV_1 - Indicator 4  

  ### Define settings for OL model component that are generic across classes
  ol_settings4X = list(
    outcomeOrdered = STATE_CONCERNED_R, 
    tau            = list(tau_scon_1X, tau_scon_2X, tau_scon_3X, tau_scon_4X),
    coding         = c(1,2,3,4,5))
  ### Loop over classes
  S = 2 # number of classes
  for(s in 1:S){
    ### Class-specific utilities
    ol_settings4X$V = zeta_sconX[[s]]*LVX1
    
    ### Within-class choice probabilities using OL model
    label = paste0("sconClass",s)
    PClass4[[label]] = apollo_ol(ol_settings4X, functionality)
    
    ### Take product across observation for same individual
    #PClass[[label]] = apollo_panelProd(PClass[[label]], apollo_inputs, functionality)
  }
  ### Mix the probabilities from each class
  lc_settings4   = list(inClassProb=PClass4, classProb=pi_values)
  P[["indic_sconX"]] = apollo_lc(lc_settings4, apollo_inputs, functionality) 

  
  ### Likelihood of the whole model
  P = apollo_combineModels(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)
}

# ################################################################# #
#### ESTIMATE SETTINGS                                           ####
# ################################################################# #

estimate_settings = list(maxIterations = 250, scaling = TRUE)

# ################################################################# #
#### MODEL ESTIMATION                                            ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs, estimate_settings)
And here's the output

Code: Select all

Model name                       :LatentClass_LV1_N1.4_imp
Model description                : ILVLC with 1 LV
Model run at                     : 2021-12-01 13:58:16
Estimation method                : bfgs
Model diagnosis                  : successful convergence 
Number of individuals            : 1201
Number of rows in database       : 1201
Number of modelled outcomes      : 8407
                     indic_bothC : 1201
                     indic_aconC : 1201
                     indic_seatC : 1201
                     indic_sconC : 1201
                  indic_offpeakC : 1201
                      indic_altC : 1201
                     indic_sconX : 1201

Number of cores used             :  3 
Number of inter-individual draws : 50 (halton)
WARNING: Inter-individual draws were used
         without a panel data structure.

LL(start)                        : -16489.55
LL(0, whole model)               : Not applicable
LL(C, whole model)               : Not applicable
LL(final, whole model)           : -12185.15
Rho-square (0)                   : Not applicable
Adj.Rho-square (0)               : Not applicable
AIC                              :  24460.3 
BIC                              :  24776.96 

LL(0,indic_bothC)                : Not applicable
LL(final,indic_bothC)            : -1737.703
LL(0,indic_aconC)                : Not applicable
LL(final,indic_aconC)            : -1928.372
LL(0,indic_seatC)                : Not applicable
LL(final,indic_seatC)            : -1912.865
LL(0,indic_sconC)                : Not applicable
LL(final,indic_sconC)            : -1779.392
LL(0,indic_offpeakC)             : Not applicable
LL(final,indic_offpeakC)         : -1912.602
LL(0,indic_altC)                 : Not applicable
LL(final,indic_altC)             : -1827.499
LL(0,indic_sconX)                : Not applicable
LL(final,indic_sconX)            : -1777.838

Estimated parameters             :  45
Time taken (hh:mm:ss)            :  00:31:14.68 
     pre-estimation              :  00:00:27.47 
     estimation                  :  00:24:53.57 
     post-estimation             :  00:05:53.64 
Iterations                       :  335  
Min abs eigenvalue of Hessian    :  1e-05 

Unconstrained optimisation.

These outputs have had the scaling used in estimation applied to them.
Estimates:
                       Estimate        s.e.   t.rat.(0)    Rob.s.e. Rob.t.rat.(0)
zeta_sconC            -0.516301    0.044179    -11.6866    0.057098       -9.0424
zeta_aconC            -1.367395    0.105542    -12.9560    0.123505      -11.0716
zeta_bothC            -0.936651    0.066140    -14.1616    0.078279      -11.9655
zeta_seatC            -0.693456    0.049631    -13.9723    0.053084      -13.0634
zeta_offpeakC         -0.576208    0.045995    -12.5276    0.058799       -9.7997
zeta_altC             -0.608091    0.049704    -12.2341    0.063638       -9.5554
tau_acon_1C           -2.735325    0.169862    -16.1032    0.185532      -14.7431
tau_acon_2C           -0.905383    0.111129     -8.1471    0.112286       -8.0632
tau_acon_3C            0.702853    0.107679      6.5273    0.109086        6.4431
tau_acon_4C            2.280086    0.151013     15.0986    0.161528       14.1157
tau_scon_1C           -3.183680    0.140997    -22.5798    0.145596      -21.8666
tau_scon_2C           -1.759633    0.089108    -19.7472    0.092471      -19.0291
tau_scon_3C           -0.177871    0.067357     -2.6407    0.068057       -2.6136
tau_scon_4C            1.269287    0.078478     16.1739    0.081167       15.6379
tau_both_1C           -3.990899    0.193115    -20.6659    0.199527      -20.0018
tau_both_2C           -2.344303    0.123126    -19.0399    0.127686      -18.3599
tau_both_3C           -0.660522    0.086547     -7.6319    0.086505       -7.6356
tau_both_4C            0.829748    0.087774      9.4533    0.088464        9.3795
tau_seat_1C           -2.153858    0.104168    -20.6768    0.104823      -20.5476
tau_seat_2C           -0.979777    0.079791    -12.2793    0.079203      -12.3705
tau_seat_3C            0.322617    0.073798      4.3716    0.073773        4.3731
tau_seat_4C            1.604270    0.090116     17.8024    0.091638       17.5065
tau_offpeak_1C        -1.989205    0.096152    -20.6882    0.102495      -19.4078
tau_offpeak_2C        -0.935259    0.074985    -12.4727    0.075837      -12.3324
tau_offpeak_3C         0.248404    0.069138      3.5929    0.069207        3.5893
tau_offpeak_4C         1.596185    0.086389     18.4768    0.088816       17.9718
tau_alt_1C            -0.675723    0.073267     -9.2227    0.075722       -8.9238
tau_alt_2C             0.253764    0.070448      3.6022    0.069945        3.6281
tau_alt_3C             1.232414    0.080909     15.2320    0.080954       15.2236
tau_alt_4C             2.300133    0.107494     21.3978    0.112074       20.5233
zeta_sconX_A         -56.710067   79.067208     -0.7172   50.526818       -1.1224
zeta_sconX_B         235.966040  310.127577      0.7609  112.277078        2.1016
tau_scon_1X           -2.878808    0.281470    -10.2278    0.460261       -6.2547
tau_scon_2X           -1.520953    0.257625     -5.9037    0.449470       -3.3839
tau_scon_3X            0.143124    0.241754      0.5920    0.418130        0.3423
tau_scon_4X            1.950833    0.278134      7.0140    0.470956        4.1423
gamma_LV1_female       0.002473    0.003350      0.7384    0.001499        1.6498
gamma_LV1_fulltime     0.003476    0.004622      0.7520    0.001162        2.9928
gamma_LV1_incomelow    0.003467    0.004670      0.7424    0.001600        2.1669
gamma_LV1_age2534      0.002807    0.003796      0.7394    0.001270        2.2095
gamma_LV1_edbach     8.6371e-04    0.001634      0.5287    0.001444        0.5980
gamma_LV1_kids      -4.2224e-04    0.001333     -0.3168    0.001353       -0.3120
gamma_LV1_car          0.002165    0.003054      0.7089    0.001515        1.4291
gamma_LV1_wave2        0.003599    0.004836      0.7441    0.002138        1.6831
b1                     1.000000          NA          NA          NA            NA
m1                     0.000000          NA          NA          NA            NA
thre1                 -0.733647    0.373854     -1.9624    0.582837       -1.2588


Summary of class allocation for LC model component indic_sconX:
            Mean prob.
sconClass1      0.6215
sconClass2      0.3785
Thank you so much for your input and feedback!
stephanehess
Site Admin
Posts: 974
Joined: 24 Apr 2020, 16:29

Re: Logistic distribution of a random parameter

Post by stephanehess »

Hi

with normalisation, I was indeed referring to normalising the variance of the latent variable.

However, having had a look at your model, I'm not quite sure why you're doing what you're doing here.

1. Why do you want the LV to have a logitstic error?
2. What is LVX1? This doesn't seem to have a random part in it, so should not be in randcoeff

Stephane
--------------------------------
Stephane Hess
www.stephanehess.me.uk
bokapatsila
Posts: 21
Joined: 28 Jul 2021, 02:41

Re: Logistic distribution of a random parameter

Post by bokapatsila »

Thank you, Stephane, for your prompt response.

In my model, I'm trying to apply direct categorization of individuals on the basis of latent variables approach as proposed by Bahamonde-Birke and Ortúzar here - https://bit.ly/3IeYJPQ. To answer your questions:

1. I want LVC1 to have a logistic error to match it with the binary logit used for class allocation
2. LVC1 is a latent variable accounting for awareness, while LVX1 is a discretized latent variable accounting for individuals with high awareness. You're correct, it doesn't have to be in the randcoeff, so I took it out.

I'm not interested in the choice the various classes are making, but whether the classes change over time (I have a binary variable accounting for the survey wave), that's why I'm using one of the indicators (STATE_CONCERNED_R) as a choice component too.

Given that, in your opinion, am getting to where I want with my code? Am I correct in my understanding that b is a scale parameter, and m is a location parameter, and the way I normalize them?

Attached is the code where I took LVX1 out of randcoeff.

Code: Select all

# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS                       ####
# ################################################################# #

### Clear memory
rm(list = ls())

### Load libraries
library(apollo)

### Initialise code
apollo_initialise()

### Set core controls
apollo_control = list(
  modelName  ="LatentClass_LV1_N1.4_imp",
  modelDescr ="ILVLC with 1 LV",
  indivID    ="UNID",
  panelData = FALSE,
  mixing     = TRUE,
#workInLogs = TRUE,
  nCores     = 3)

# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS                     ####
# ################################################################# #

setwd("D:/Data")

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

# ################################################################# #
#### DEFINE MODEL PARAMETERS                                     ####
# ################################################################# #

### Vector of parameters, including any that are kept fixed in estimation
apollo_beta=c(zeta_aconC = 1, 
              zeta_sconC = 1, 
              zeta_bothC = 1,
              zeta_seatC = 1,
              zeta_offpeakC = 1, 
              zeta_altC = 1, 
              tau_acon_1C =-2, 
              tau_acon_2C =-1, 
              tau_acon_3C = 1, 
              tau_acon_4C = 2,
              tau_scon_1C =-2, 
              tau_scon_2C =-1, 
              tau_scon_3C = 1, 
              tau_scon_4C = 2,
              tau_both_1C =-2, 
              tau_both_2C =-1, 
              tau_both_3C = 1, 
              tau_both_4C = 2,
              tau_seat_1C =-2, 
              tau_seat_2C =-1, 
              tau_seat_3C = 1, 
              tau_seat_4C = 2,
              tau_offpeak_1C =-2, 
              tau_offpeak_2C =-1,
              tau_offpeak_3C = 1, 
              tau_offpeak_4C = 2,
              tau_alt_1C =-2, 
              tau_alt_2C =-1,
              tau_alt_3C = 1, 
              tau_alt_4C = 2,

              zeta_sconX_A = 1,
              zeta_sconX_B = -1,

              tau_scon_1X =-2, 
              tau_scon_2X =-1,
              tau_scon_3X = 1, 
              tau_scon_4X = 2,

              gamma_LV1_female  = 1, 
              gamma_LV1_fulltime  = 1,
              gamma_LV1_incomelow  = 1,
              gamma_LV1_age2534  = 1,
              gamma_LV1_edbach = 1,
              gamma_LV1_kids = 1,
              gamma_LV1_car = 1,
              gamma_LV1_wave2  = 1,
              
              b1 = 1,
              m1 = 0,
              thre1 = 1)
              

### 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
apollo_fixed = c("b1", "m1")

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

### Set parameters for generating draws
apollo_draws = list(
  interDrawsType="halton", 
  interNDraws=50,          
  interUnifDraws=c("eta1"),      
  interNormDraws=c(), 
  
  intraDrawsType='',
  intraNDraws=0,          
  intraUnifDraws=c(),     
  intraNormDraws=c()      
)

### Create random parameters
apollo_randCoeff=function(apollo_beta, apollo_inputs){
  randcoeff = list()
  
  randcoeff[["LVC1"]] = gamma_LV1_female * FEMALE + gamma_LV1_wave2 * WAVE2 +
    gamma_LV1_fulltime * EMPLOY_FULL + gamma_LV1_incomelow * INCOME_LOW + gamma_LV1_age2534 * AGE_2534 + 
    gamma_LV1_edbach * EDUCATION_BACH + gamma_LV1_kids * HOUSEHOLD_CHILD_CLEAN_N + gamma_LV1_car * CAR_B +
    2*b1*atanh(2*eta1-1)+m1 #Error term follows a logistic distribution (to match the logit in allocation function)
  
  return(randcoeff)
}

# ################################################################# #
#### DEFINE LATENT CLASS COMPONENTS                              ####
# ################################################################# #

apollo_lcPars=function(apollo_beta, apollo_inputs){
  lcpars = list()
  
  ### Define lists of parameters for each class
  ###                                    classA   classB  ...
  lcpars[["zeta_sconX"]]  = list(zeta_sconX_A, zeta_sconX_B)
  
  ### Class allocation probabilities
  ### These are the probabilities of a binary logit model
  classAlloc_settings = list(
    V = list(A = LVC1-thre1,
             B = 0)
  )
  lcpars[["pi_values"]] = apollo_classAlloc(classAlloc_settings)
  
  
  return(lcpars)
}

# ################################################################# #
#### 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()
  PClass4 = list()
  
  ### Likelihood of indicators
  ol_settings1C = list(outcomeOrdered=AGREE_BOTH_CROWD_BC_R, 
                      V=zeta_bothC*LVC1, 
                      tau=list(tau_both_1C, tau_both_2C, tau_both_3C, tau_both_4C))
  ol_settings2C = list(outcomeOrdered=AGREE_CONCERN_BC_R, 
                      V=zeta_aconC*LVC1, 
                      tau=list(tau_acon_1C, tau_acon_2C, tau_acon_3C, tau_acon_4C))
  ol_settings3C = list(outcomeOrdered=AGREE_SEAT_BC_R, 
                      V=zeta_seatC*LVC1, 
                      tau=list(tau_seat_1C, tau_seat_2C, tau_seat_3C, tau_seat_4C))
  ol_settings4C = list(outcomeOrdered=STATE_CONCERNED_R, 
                      V=zeta_sconC*LVC1, 
                      tau=list(tau_scon_1C, tau_scon_2C, tau_scon_3C, tau_scon_4C))
  ol_settings5C = list(outcomeOrdered=AGREE_OFFPEAK_BC_R, 
                      V=zeta_offpeakC*LVC1, 
                      tau=list(tau_offpeak_1C, tau_offpeak_2C, tau_offpeak_3C, tau_offpeak_4C))
  ol_settings6C = list(outcomeOrdered=AGREE_ALT_BC_R, 
                      V=zeta_altC*LVC1, 
                      tau=list(tau_alt_1C, tau_alt_2C, tau_alt_3C, tau_alt_4C))
  
  P[["indic_bothC"]]     = apollo_ol(ol_settings1C, functionality)
  P[["indic_aconC"]]     = apollo_ol(ol_settings2C, functionality)
  P[["indic_seatC"]]      = apollo_ol(ol_settings3C, functionality)
  P[["indic_sconC"]]      = apollo_ol(ol_settings4C, functionality)
  P[["indic_offpeakC"]]     = apollo_ol(ol_settings5C, functionality)
  P[["indic_altC"]]      = apollo_ol(ol_settings6C, functionality)

  ### Categorical LV_1 - Indicator 4  
  
    LVX1 = gamma_LV1_female * FEMALE + gamma_LV1_wave2 * WAVE2 + 
    gamma_LV1_fulltime * EMPLOY_FULL + gamma_LV1_incomelow * INCOME_LOW + gamma_LV1_age2534 * AGE_2534 + 
    gamma_LV1_edbach * EDUCATION_BACH + gamma_LV1_kids * HOUSEHOLD_CHILD_CLEAN_N + gamma_LV1_car * CAR_B

  ### Define settings for OL model component that are generic across classes
  ol_settings4X = list(
    outcomeOrdered = STATE_CONCERNED_R, 
    tau            = list(tau_scon_1X, tau_scon_2X, tau_scon_3X, tau_scon_4X),
    coding         = c(1,2,3,4,5))
  ### Loop over classes
  S = 2 # number of classes
  for(s in 1:S){
    ### Class-specific utilities
    ol_settings4X$V = zeta_sconX[[s]]*LVX1
    
    ### Within-class choice probabilities using OL model
    label = paste0("sconClass",s)
    PClass4[[label]] = apollo_ol(ol_settings4X, functionality)
    
    ### Take product across observation for same individual
    #PClass[[label]] = apollo_panelProd(PClass[[label]], apollo_inputs, functionality)
  }
  ### Mix the probabilities from each class
  lc_settings4   = list(inClassProb=PClass4, classProb=pi_values)
  P[["indic_sconX"]] = apollo_lc(lc_settings4, apollo_inputs, functionality) 

  
  ### Likelihood of the whole model
  P = apollo_combineModels(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)
}

# ################################################################# #
#### ESTIMATE SETTINGS                                           ####
# ################################################################# #

estimate_settings = list(maxIterations = 250, scaling = TRUE)

# ################################################################# #
#### MODEL ESTIMATION                                            ####
# ################################################################# #

model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs, estimate_settings)
dpalma
Posts: 190
Joined: 24 Apr 2020, 17:54

Re: Logistic distribution of a random parameter

Post by dpalma »

Hi,

Sorry for the very long delay in replying. It wasn't a simple question, so it kind of got put aside for a long while.

I haven't gone through the whole Bahamonde-Birke & Ortúzar paper, but I did notice that you don't include an error term in LVC1 because you mention this is essentially the deterministic part of the utility in the class allocation model model. This is fine, but the problem is you later use LVC1 in your measurement equations as well, without any error term.

What I guess you want to do is using the same error term of the class allocation into the measurement equations. But I a afraid you cannot do that in apollo. By using an MNL in class allocation you are integrating over the Logistic error term (Logistic and not Gumbel because it's just a binary logit), which leads to the MNL probability. But if you want to use the same error term in all those measurement equations at the same time, you should integrate over all those things at the same time, which I am guessing will not lead to a product of one MNL and multiple ordered logits, because of the interactions between them.

Of course, this is what I can gather from your last post, but I may be missing a lot because I haven't gone through the whole paper. Please let us know if this is still an issue and we can continue discussing.

Best wishes
David
bokapatsila
Posts: 21
Joined: 28 Jul 2021, 02:41

Re: Logistic distribution of a random parameter

Post by bokapatsila »

Thanks, David! Back when I asked this question, I wasn't fully aware of what I was asking about :-) In any case, we figured it all out.
Post Reply