AGSCL
commented
11 months ago I decided to calculate them manually applying the formula to transform odds ratios to probabilities in multinomial logistic regressions (in this example I had 5 classes and two terms: "(Intercept)" and "outcome1"):
`# Use map_df to loop through each call class and apply the summarise logic
Long, S. and Freese, J. (2014) Regression Models for Categorical Dependent Variables Using Stata. 3rd Edition, Stata Press, College Station.
https://www.ncbi.nlm.nih.gov/pmc/articles/PMC9041638/ #eq 3
$\frac{exp(\alpha+\sum{k=1}^{K}\beta{mk}X{ik})}{1+exp(\alpha+\sum{k=1}^{K}\beta{mk}X{ik})}$
https://www3.nd.edu/~rwilliam/stats3/Mlogit1.pdf
require(glca) df_best_model_glca_w_distal_outcome<- data.frame(coef(best_model_glca_w_distal_outcome)) %>% rownames_to_column("term") %>% janitor::clean_names() %>% rownames_to_column("rowname") %>% gather(key = "key", value = "value", -rowname) %>% spread(key = "rowname", value = "value") %>% setnames(as.character(unlist(tail(., 1)))) %>% slice(-n()) %>% dplyr::mutate(term=strsplit(sub('^(.*?.?_.?)_(.*)$', '\1,\2', term), ',')) %>% separate(col=term,into = c("prefix", "suffix"), sep = ", ", extra = "merge") %>% dplyr::mutate(across(c("prefix", "suffix"), ~gsub('\(|"|\)', "", .)))
df_best_model_glca_w_distal_outcome2<-
df_best_model_glca_w_distal_outcome%>%
dplyr::filter(suffix == "coefficient" | suffix == "std_error") %>%
pivot_wider(names_from=suffix, values_from = c("(Intercept)", "outcome1")) %>%
dplyr::mutate_at(2:5, ~as.numeric(.)) %>%
dplyr::mutate(
lower_log_or_int = (Intercept)_coefficient - 1.96 (Intercept)_std_error,
upper_log_or_int = (Intercept)_coefficient + 1.96 outcome1_std_error,
lower_log_or_comp = outcome1_coefficient - 1.96 outcome1_std_error,
upper_log_or_comp = outcome1_coefficient + 1.96 outcome1_std_error) %>%
dplyr::rename("int_coef"="(Intercept)_coefficient", "int_std_error"="(Intercept)_std_error", "comp_coef"="outcome1_coefficient","comp_std_error"="outcome1_std_error") %>%
dplyr::select(prefix,#t_coef int_std_error comp_coef comp_std_error
int_coef, int_std_error, comp_coef, comp_std_error,
lower_log_or_int, upper_log_or_int, lower_log_or_comp, upper_log_or_comp)
List of call classes
call_classes <- unique(df_best_model_glca_w_distal_outcome2$prefix)
result2 <- map_df(call_classes, ~ { df_best_model_glca_w_distal_outcome2 %>% dplyr::mutate(int_comp_coef=int_coef+comp_coef) %>% dplyr::summarise(call_class = .x, nom = sum(ifelse(prefix == .x, exp(int_comp_coef), 0)), den = sum(exp(int_comp_coef)), prob = nom/(1+den)) })
result2_lo <- map_df(call_classes, ~ { df_best_model_glca_w_distal_outcome2 %>% dplyr::mutate(int_comp_coef=lower_log_or_int+lower_log_or_comp) %>% dplyr::summarise(call_class = .x, nom = sum(ifelse(prefix == .x, exp(int_comp_coef), 0)), den = sum(exp(int_comp_coef)), prob = nom/(1+den)) })
result2_hi <- map_df(call_classes, ~ { df_best_model_glca_w_distal_outcome2 %>% dplyr::mutate(int_comp_coef=upper_log_or_int+upper_log_or_comp) %>% dplyr::summarise(call_class = .x, nom = sum(ifelse(prefix == .x, exp(int_comp_coef), 0)), den = sum(exp(int_comp_coef)), prob = nom/(1+den)) })
result2b <- map_df(call_classes, ~ { df_best_model_glca_w_distal_outcome2 %>% dplyr::mutate(int_comp_coef=int_coef) %>% dplyr::summarise(call_class = .x, nom = sum(ifelse(prefix == .x, exp(int_comp_coef), 0)), den = sum(exp(int_comp_coef)), prob = nom/(1+den)) })
result2b_lo <- map_df(call_classes, ~ { df_best_model_glca_w_distal_outcome2 %>% dplyr::mutate(int_comp_coef=lower_log_or_int) %>% dplyr::summarise(call_class = .x, nom = sum(ifelse(prefix == .x, exp(int_comp_coef), 0)), den = sum(exp(int_comp_coef)), prob = nom/(1+den)) })
result2b_hi <- map_df(call_classes, ~ { df_best_model_glca_w_distal_outcome2 %>% dplyr::mutate(int_comp_coef=upper_log_or_int) %>% dplyr::summarise(call_class = .x, nom = sum(ifelse(prefix == .x, exp(int_comp_coef), 0)), den = sum(exp(int_comp_coef)), prob = nom/(1+den)) })
df_lca40522_probs<- cbind.data.frame( est=c(result2$prob, 1/(1+unique(result2$den))), lo=c(result2_lo$prob, 1/(1+unique(result2_lo$den))), hi=c(result2_hi$prob, 1/(1+unique(result2_hi$den))), est_int=c(result2b$prob, 1/(1+unique(result2b$den))), lo_int=c(result2b_lo$prob, 1/(1+unique(result2b_lo$den))), hi_int=c(result2b_hi$prob, 1/(1+unique(result2b_hi$den))))*100
df_lca40522_probs %>% knitr::kable("markdown", caption="Tabla de probabilidades de LCA con resultado distal (int: no interrumpe)", digits=1) `
I'm using the glca package to estimate a latent class model with covariates. In the table of covariate coefficients, the reference value for the regiao variable is set to CENTRO-OESTE, but I would like to set it to SUDESTE.
However, I can't seem to find any argument in the
glcafunction to set the reference levels for each covariate. I checked the package documentation but there's no mention of how to do this.Here is an example of the coefficient table: