I need some help with running a simple latent class model with 2 classes. The model gives the following error: "WARNING: Estimation succeeded but the BHHH matrix is singular. Apollo will not attempt to compute the covariance matrix." I am not certain what the issue is - its either typographic or down to the priors used. Currently, I am using the priors from the MNL output for class 1 and 0 priors for class 2.
I have been able to run the MNL, EC, MXL and MXL with correlated errors successfully. The MXL models are why I have chosen to explore latent class - there is significant heterogeneity in the model. Particularly in the asc and mam_rev# categorical variables.
DCE info: 24 choice tasks (in 3 blocks of 8), 3 alternatives (2 appointments & 1 opt out (not attend)), 5 attributes (1 dummy with 4 levels(mam_rev#), 4 continuous(can_det, inap_rec, over, wait)). The asc is coded as 1 for either appointment and 0 for the opt out.
Any advice would be much appreciated - I have run the analysis in STATA(lclogit package) and tried using the results as priors in apollo but I am getting the same issue.
Code: Select all
# ################################################################# #
#### LOAD LIBRARY AND DEFINE CORE SETTINGS ####
# ################################################################# #
## CLEAR MEMORY
rm(list = ls())
## LOAD PACKAGES
library(haven)
library(apollo)
library(dplyr)
library(parallel)
# CHECK AVAILABLE CORES
num_cores <- parallel::detectCores()
print(num_cores)
## INITIALISE CODE
apollo_initialise()
## SET CORE CONTROLS
apollo_control = list(
modelDescr = "PreferAI: LC_MNL preference space model",
modelName = "PreferAI_lc_model_simp",
indivID = "id",
outputDirectory = "output",
nCores = 4
)
# ################################################################# #
#### LOAD DATA AND APPLY ANY TRANSFORMATIONS ####
# ################################################################# #
# Load in data from cleaned STATA dataset (panel data with 24 rows per respondent (3 alternatives * 8 choice tasks))
LongAnalysisSet <- read_dta("LongAnalysisSet.dta")
## Ensure data is ordered by id and choice task number
LongAnalysisSet <- LongAnalysisSet[order(LongAnalysisSet$id, LongAnalysisSet$csno, LongAnalysisSet$alt), ]
## Adjust column ordering to aid data exploration
ChoiceData <-LongAnalysisSet %>% select(id, csno, ncs, alt, choice, mam_rev1, mam_rev2, mam_rev3, mam_rev4, can_det, inap_rec, over, wait)
head(ChoiceData)
## Alter dataset to Apollo wide format (1 row per choice task)
longToWide_settings <- list(
altColumn = "alt",
altSpecAtts = c("mam_rev2", "mam_rev3", "mam_rev4", "can_det", "inap_rec", "over", "wait"),
choiceColumn = "choice",
idColumn = "id",
obsColumn = "ncs"
)
## Reshaped data
database = apollo_longToWide(ChoiceData,longToWide_settings)
## Rename new choice variable back to choice for clarity
head(database$choice_new)
database <- database %>%
dplyr::rename(choice = choice_new)
# ################################################################# #
#### DEFINE MODEL PARAMETERS ####
# ################################################################# #
apollo_beta <- c(
# Means
# Alternative specific constant
mu_v_asc_a = -1.4,
mu_v_asc_b = 4.29,
# Categorical mammography review variables (excluding mam_rev1 base)
mu_v_mam_rev2_a = 0.25,
mu_v_mam_rev3_a = 0.14,
mu_v_mam_rev4_a = -0.36,
mu_v_mam_rev2_b = 0,
mu_v_mam_rev3_b = 0,
mu_v_mam_rev4_b = 0,
# Continuous variables
mu_v_cd_a = 0.05,
mu_v_ir_a = -0.01,
mu_v_ov_a = -0.16,
mu_v_wt_a = -0.12,
mu_v_cd_b = 0,
mu_v_ir_b = 0,
mu_v_ov_b = 0,
mu_v_wt_b = 0,
delta_a=0,
delta_b=0
)
apollo_fixed <- c()
# ################################################################# #
#### DEFINE LATENT CLASS COMPONENTS ####
# ################################################################# #
apollo_lcPars = function(apollo_beta, apollo_inputs){
lcpars = list()
lcpars[["asc"]] = list(mu_v_asc_a, mu_v_asc_b)
lcpars[["mr2"]] = list(mu_v_mam_rev2_a, mu_v_mam_rev2_b)
lcpars[["mr3"]] = list(mu_v_mam_rev3_a, mu_v_mam_rev3_b)
lcpars[["mr4"]] = list(mu_v_mam_rev4_a, mu_v_mam_rev4_b)
lcpars[["cd"]] = list(mu_v_cd_a,mu_v_cd_b)
lcpars[["ir"]] = list(mu_v_ir_a,mu_v_ir_b)
lcpars[["ov"]] = list(mu_v_ov_a,mu_v_ov_b)
lcpars[["wt"]] = list(mu_v_wt_a,mu_v_wt_b)
V=list()
V[["class_a"]] = delta_a
V[["class_b"]] = delta_b
classAlloc_settings = list(
classes = c(class_a=1, class_b=2),
utilities = V
)
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()
# Define settings for MNL model component that are generic across classes
mnl_settings = list(
alternatives = c(A = 1, B = 2, OO = 3),
choiceVar = choice[,1]
)
### Loop over classes
for(s in 1:2){
### Compute class-specific utilities
V=list()
V = list()
V[["A"]] = asc[[s]] +
mr2[[s]]* mam_rev2_1 +
mr3[[s]]* mam_rev3_1 +
mr4[[s]]* mam_rev4_1 +
cd[[s]] * can_det_1 +
ir[[s]] * inap_rec_1 +
ov[[s]] * over_1 +
wt[[s]] * wait_1
V[["B"]] = asc[[s]]+
mr2[[s]]* mam_rev2_2 +
mr3[[s]]* mam_rev3_2 +
mr4[[s]]* mam_rev4_2 +
cd[[s]] * can_det_2 +
ir[[s]] * inap_rec_2 +
ov[[s]] * over_2 +
wt[[s]] * wait_2
V[["OO"]] = 0 # Opt-out
mnl_settings$utilities = V
### Compute within-class choice probabilities using MNL model
P[[paste0("Class_",s)]] = apollo_mnl(mnl_settings, functionality)
### Take product across observation for same individual
P[[paste0("Class_",s)]] = apollo_panelProd(P[[paste0("Class_",s)]], apollo_inputs ,functionality)
}
### 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 ####
# ################################################################# #
lc_model_simp = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)
# ################################################################# #
#### MODEL OUTPUTS ####
# ################################################################# #
# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO SCREEN) ----
# ----------------------------------------------------------------- #
apollo_modelOutput(lc_model_simp)
# ----------------------------------------------------------------- #
#---- FORMATTED OUTPUT (TO FILE, using model name) ----
# ----------------------------------------------------------------- #
apollo_saveOutput(lc_model_simp)