I'm trying to run a simple Latent Class model in Apollo 0.2.5 but am running into an error when computing posterior class membership probabilities with apollo_lcConditionals. The error I get reads
"Error in data.frame(ID = unique(database[, apollo_control$indivID]), conditionals) :
arguments imply differing number of rows: 1274, 5499".
1274 is the number of individual IDs in the dataset, and 5499 is the number of rows. So the issue seems to be that "conditionals" is computed for every row in the dataset, when apollo_lcConditionals needs one row per individual. I'm not sure why that's happening though.
My code is below in case helpful, though I can't share the data.
Thanks so much!
Code: Select all
# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS ####
# ################################################################# #
### Clear memory
rm(list = ls())
### Load Apollo library
library(apollo)
library(data.table)
### Initialise code
apollo_initialise()
apollo_control = list(
modelName ="Apollo_example_20",
modelDescr ="LC model of student participation",
indivID ="household_anon",
nCores = 4,
analyticGrad = FALSE
)
# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #
dt = fread(...)
cols <- grep("hh_allswab|hrf",names(dt),value = T)
database <- dt[term == "E21" & consent == T & wifi == 1 & between(n_hh,4,7) & week %in% c(26:31),
lapply(.SD,function(x) na.omit(unique(x))),.SDcols = cols,.(household_anon,week)]
database <- database[order(household_anon,week)]
database <- na.omit(database[,.(household_anon,week,hh_allswab,hh_allswab_l1,sex_hrf)])
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
### Vector of parameters, including any that are kept fixed in estimation
apollo_beta = c(delta_a = 0,
delta_b = 1,
delta_c = 1,
delta_d = 1,
alpha_a = 0.5,
alpha_b = 0.5,
alpha_c = 0.5,
alpha_d = 0.5,
rho_a = 0.4,
rho_b = 0,
rho_c = 0,
rho_d = 0)
### 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("delta_a")
# ################################################################# #
#### DEFINE LATENT CLASS COMPONENTS ####
# ################################################################# #
# apollo_lcPars=function(apollo_beta, apollo_inputs){
# lcpars = list()
# lcpars[["alpha"]] = list(alpha_a, alpha_b, alpha_c, alpha_d)
# lcpars[["rho"]] = list(rho_a, rho_b, rho_c, rho_d)
#
# V=list()
# V[["class_a"]] = 0
# V[["class_b"]] = delta_b
# V[["class_c"]] = delta_c
# V[["class_d"]] = delta_d
#
#
# mnl_settings = list(
# alternatives = c(class_a=1, class_b=2, class_c=3, class_d=4),
# avail = 1,
# choiceVar = NA,
# V = V
# )
# lcpars[["pi_values"]] = apollo_mnl(mnl_settings, functionality="raw")
#
# lcpars[["pi_values"]] = apollo_firstRow(lcpars[["pi_values"]], apollo_inputs)
#
# return(lcpars)
# }
apollo_lcPars=function(apollo_beta, apollo_inputs){
lcpars = list()
lcpars[["alpha"]] = list(alpha_a, alpha_b, alpha_c, alpha_d)
lcpars[["rho"]] = list(rho_a, rho_b, rho_c, rho_d)
lcpars[["pi_values"]] = list(class_a = exp(delta_a)/(exp(delta_a) + exp(delta_b) + exp(delta_c) + exp(delta_d)),
class_b = exp(delta_b)/(exp(delta_a) + exp(delta_b) + exp(delta_c) + exp(delta_d)),
class_c = exp(delta_c)/(exp(delta_a) + exp(delta_b) + exp(delta_c) + exp(delta_d)),
class_d = exp(delta_d)/(exp(delta_a) + exp(delta_b) + exp(delta_c) + exp(delta_d)))
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()
### Define settings for MNL model component that are generic across classes
mnl_settings = list(
alternatives = c(alt1=0, alt2=1),
avail = list(alt1=1, alt2=1),
choiceVar = hh_allswab
)
### Loop over classes
s=1
while(s<=4){
### Compute class-specific utilities
V=list()
V[['alt1']] = 0
V[['alt2']] = alpha[[s]] + rho[[s]]*hh_allswab_l1
mnl_settings$V = V
mnl_settings$componentName = paste0("Class_",s)
### Compute within-class choice probabilities using MNL model
P[[paste0("Class_",s)]] = apollo_ol(mnl_settings, functionality)
### Take product across observation for same individual
P[[paste0("Class_",s)]] = apollo_panelProd(P[[paste0("Class_",s)]], apollo_inputs ,functionality)
s=s+1}
### Compute latent class model probabilities
lc_settings = list(inClassProb = P, classProb=pi_values)
P[["model"]] = apollo_lc(lc_settings, apollo_inputs, functionality)
### Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
}
# ################################################################# #
#### MODEL ESTIMATION ####
# ################################################################# #
#apollo_beta=apollo_searchStart(apollo_beta, apollo_fixed,apollo_probabilities, apollo_inputs)
#apollo_outOfSample(apollo_beta, apollo_fixed,apollo_probabilities, apollo_inputs)
### Estimate model
model = apollo_estimate(apollo_beta, apollo_fixed,
apollo_probabilities, apollo_inputs,
estimate_settings=list(writeIter=FALSE))
### Show output in screen
apollo_modelOutput(model)
### Save output to file(s)
apollo_saveOutput(model)
# ################################################################# #
##### POST-PROCESSING ####
# ################################################################# #
### Print outputs of additional diagnostics to new output file (remember to close file writing when complete)
sink(paste(model$apollo_control$modelName,"_additional_output.txt",sep=""),split=TRUE)
# ----------------------------------------------------------------- #
#---- OUT OF SAMPLE TESTING ----
# ----------------------------------------------------------------- #
apollo_outOfSample(apollo_beta, apollo_fixed,
apollo_probabilities, apollo_inputs)
# ----------------------------------------------------------------- #
#---- BOOTSTRAP ESTIMATION ----
# ----------------------------------------------------------------- #
apollo_bootstrap(apollo_beta, apollo_fixed,apollo_probabilities, apollo_inputs)
# ----------------------------------------------------------------- #
#---- POSTERIOR ANALYSIS ----
# ----------------------------------------------------------------- #
### Compute unconditionals
unconditionals = apollo_lcUnconditionals(model,apollo_probabilities,apollo_inputs)
vtt_class_a = unconditionals[["beta_tt"]][[1]]/unconditionals[["beta_tc"]][[1]]
vtt_class_b = unconditionals[["beta_tt"]][[2]]/unconditionals[["beta_tc"]][[2]]
vtt_unconditional = unconditionals[["pi_values"]][[1]]*vtt_class_a +
unconditionals[["pi_values"]][[2]]*vtt_class_b
### Compute conditionals
conditionals = apollo_lcConditionals(model,apollo_probabilities, apollo_inputs)
summary(conditionals)
summary(as.data.frame(unconditionals[["pi_values"]]))
vtt_conditional=conditionals[,1]*vtt_class_a+conditionals[,2]*vtt_class_b
summary(vtt_unconditional)
summary(vtt_conditional)
### Take first value of covariates for each person
commute_n = apollo_firstRow(database$commute, apollo_inputs)
car_availability_n = apollo_firstRow(database$car_availability, apollo_inputs)
### Compute posterior values for covariates
post_commute=colSums(commute_n*conditionals)/colSums(conditionals)
post_car_availability=colSums(car_availability_n*conditionals)/colSums(conditionals)
post_commute
post_car_availability
# ----------------------------------------------------------------- #
#---- switch off writing to file ----
# ----------------------------------------------------------------- #
if(sink.number()>0) sink()