Page 1 of 1

apollo_lcConditionals error (Latent Class model)

Posted: 20 Sep 2021, 12:36
by vivrc
Hi,

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()

Re: apollo_lcConditionals error (Latent Class model)

Posted: 20 Sep 2021, 13:27
by stephanehess
Hi

this is not an issue we've seen before, and there's nothing a priori wrong with your code. However, it's difficult to debug it without access to the data. Is there a way for you to share even part of the data with us by e-mail, we will then delete it after fixing your issue

Stephane

Re: apollo_lcConditionals error (Latent Class model)

Posted: 20 Sep 2021, 15:28
by vivrc
stephanehess wrote: 20 Sep 2021, 13:27 Hi

this is not an issue we've seen before, and there's nothing a priori wrong with your code. However, it's difficult to debug it without access to the data. Is there a way for you to share even part of the data with us by e-mail, we will then delete it after fixing your issue

Stephane
Hi Stephane, just to conclude the conversation we had by email. The solution, weirdly, seems to be to add two lines to the code as below. If I comment out the lines where I write.csv and then fread "database", I get the error.

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",
  panelData = T,
  nCores     = 4,
  analyticGrad = F
)

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

dt = fread("...",header=TRUE)
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)])

write.csv(database,"...")
 
database <- fread("...")

# ################################################################# #
#### 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)
  
  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_mnl(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()