library(readxl)
library(dplyr)
BASE <- read_excel("BD_articuloaceptabilidad_autotomaVPH.xlsx") #cargue base datos
BASE1 <- BASE %>% select(edad, tipo_afiliacion, vacuna_vph, meses_lastscreen, n_partos, resultado_vph) #seleccionar
BASE1 <- as.data.frame(BASE1) #manejo como dataframe
#str(BASE1)
BASE1$tipo_afiliacion <- factor(BASE1$tipo_afiliacion, levels = c(1,2,3), labels = c("Contributivo", "Subsidiado", "Otro"))
BASE1$vacuna_vph <- factor(BASE1$vacuna_vph, levels = c(1,2), labels = c("Si", "No"))
BASE1$resultado_vph <- replace(BASE1$resultado_vph, BASE1$resultado_vph == 3, NA)
BASE1$resultado_vph <- factor(BASE1$resultado_vph, levels = c(1,2), labels = c("Detectado", "No detectado"))
library(gtsummary)
tabla1 <- BASE1 %>%
select(edad, tipo_afiliacion, vacuna_vph, meses_lastscreen, n_partos, resultado_vph) %>%
tbl_summary(
by = resultado_vph,
type = list(
all_continuous() ~ "continuous2",
all_categorical() ~ "categorical"
),
label = list(
edad ~ "Edad",
tipo_afiliacion ~ "Afiliacion",
vacuna_vph ~"Vacunacion vph",
meses_lastscreen ~ "Tiempo en meses último tamizaje",
n_partos ~ "Paridad"
),
statistic = list(
all_continuous2() ~ c(
"{mean} ({sd})",
"{median} ({p25}, {p75})",
"{min} - {max}"
),
all_categorical() ~ "{n} ({p}%)"
),
digits = list(
all_continuous() ~ 2,
all_categorical() ~ 2
),
missing = "ifany",
missing_text = "Sin dato"
) %>%
add_p(
test = list(
all_continuous() ~ "kruskal.test",
all_categorical() ~ "chisq.test"
),
pvalue_fun = ~ style_pvalue(.x, digits = 3)
) %>%
add_overall() %>%
modify_header(label ~ "**Variable**") %>%
bold_labels() %>%
bold_p(t = 0.05) %>%
modify_caption("**Tabla 1**")
tabla1
| Variable | Overall N = 8691 |
Detectado N = 2261 |
No detectado N = 6431 |
p-value2 |
|---|---|---|---|---|
| Edad | 0.263 | |||
| Mean (SD) | 43.91 (9.87) | 43.42 (10.24) | 44.09 (9.73) | |
| Median (Q1, Q3) | 42.00 (36.00, 51.00) | 41.00 (35.00, 50.00) | 43.00 (36.00, 51.00) | |
| Min - Max | 24.00 - 77.00 | 25.00 - 77.00 | 24.00 - 69.00 | |
| Sin dato | 3 | 0 | 3 | |
| Afiliacion | 0.030 | |||
| Contributivo | 129.00 (15.58%) | 28.00 (13.02%) | 101.00 (16.48%) | |
| Subsidiado | 697.00 (84.18%) | 185.00 (86.05%) | 512.00 (83.52%) | |
| Otro | 2.00 (0.24%) | 2.00 (0.93%) | 0.00 (0.00%) | |
| Sin dato | 41 | 11 | 30 | |
| Vacunacion vph | 0.731 | |||
| Si | 98.00 (13.80%) | 27.00 (14.84%) | 71.00 (13.45%) | |
| No | 612.00 (86.20%) | 155.00 (85.16%) | 457.00 (86.55%) | |
| Sin dato | 159 | 44 | 115 | |
| Tiempo en meses último tamizaje | 0.168 | |||
| Mean (SD) | 49.45 (435.50) | 32.26 (28.32) | 55.37 (504.69) | |
| Median (Q1, Q3) | 25.82 (16.89, 39.39) | 25.07 (13.96, 37.06) | 25.89 (17.74, 40.90) | |
| Min - Max | 0.36 - 12,015.80 | 0.56 - 234.09 | 0.36 - 12,015.80 | |
| Sin dato | 108 | 31 | 77 | |
| Paridad | 0.214 | |||
| Mean (SD) | 2.39 (1.47) | 2.31 (1.52) | 2.41 (1.45) | |
| Median (Q1, Q3) | 2.00 (2.00, 3.00) | 2.00 (1.00, 3.00) | 2.00 (2.00, 3.00) | |
| Min - Max | 0.00 - 10.00 | 0.00 - 8.00 | 0.00 - 10.00 | |
| Sin dato | 4 | 1 | 3 | |
| 1 n (%) | ||||
| 2 Kruskal-Wallis rank sum test; Pearson’s Chi-squared test | ||||
BASE2 <- BASE %>% select(encuesta_post_toma, como_se_sintio, razones_no_satisfaccion___1, razones_no_satisfaccion___2, razones_no_satisfaccion___3, razones_no_satisfaccion___4, razones_no_satisfaccion___5, otra_no_satisfaccion, confianza_prueba, toma_prueba)
#lapply(BASE2, table)
BASE2$emocion <- ifelse(
BASE2$razones_no_satisfaccion___1 == 1 |
BASE2$razones_no_satisfaccion___2 == 1 |
grepl("inseguridad", BASE2$otra_no_satisfaccion, ignore.case = TRUE),
2, 1)
BASE2$tecnico <- ifelse(
BASE2$razones_no_satisfaccion___3 == 1 |
BASE2$razones_no_satisfaccion___4 == 1 |
grepl("instrumento", BASE2$otra_no_satisfaccion, ignore.case = TRUE),
2, 1
)
BASE2$confianza <- ifelse(BASE2$confianza_prueba == 1, 1, 2)
BASE2$acept_autotoma <- ifelse(BASE2$toma_prueba == 1, 1,2)
library(dplyr)
library(ggplot2)
BASE3 <- as.data.frame(BASE[, c("resultado_vph", "como_se_sintio",
"razones_no_satisfaccion___1", "razones_no_satisfaccion___2",
"razones_no_satisfaccion___3", "razones_no_satisfaccion___4",
"razones_no_satisfaccion___5", "otra_no_satisfaccion",
"confianza_prueba", "toma_prueba")])
BASE3$como_se_sintio <- factor(BASE3$como_se_sintio, levels = c(1,2), labels = c("Bien", "No Bien"))
BASE3$resultado_vph <- replace(BASE3$resultado_vph, BASE3$resultado_vph == 3, NA)
BASE3$resultado_vph <- factor(BASE3$resultado_vph, levels = c(1,2), labels = c("Detectado", "No detectado"))
BASE3$barrera_emo <- factor(ifelse(
BASE3$razones_no_satisfaccion___1 == 1 |
BASE3$razones_no_satisfaccion___2 == 1 |
grepl("inseguridad", BASE3$otra_no_satisfaccion, ignore.case = TRUE),
"Si", "No"))
BASE3$barrera_tec <- factor(ifelse(
BASE3$razones_no_satisfaccion___3 == 1 |
BASE3$razones_no_satisfaccion___4 == 1 |
grepl("instrumento", BASE3$otra_no_satisfaccion, ignore.case = TRUE),
"Si", "No"
))
BASE3$confianza_prueba <- ifelse(BASE2$confianza_prueba == 1, 1, 2)
BASE3$confianza_prueba <- factor(BASE3$confianza_prueba, levels = c(1,2), labels = c("Si", "No"))
BASE3$toma_prueba <- factor(BASE3$toma_prueba, levels = c(1,2,3), labels = c("Misma mujer", "Profesional salud", "Sin preferencia"))
BASE3 <- BASE3 %>%
filter(!is.na(toma_prueba)) #870 --> 856
BASE3 <- BASE3 %>%
filter(!is.na(resultado_vph)) #856 --> 855
#BASE2$consumo_diario_std_3cat[is.na(BASE1$consumo_diario_std_3cat)] <- "Consumo Ligero"
datos_sankey <- as.data.frame(BASE3[, c("resultado_vph", "como_se_sintio", "barrera_emo", "barrera_tec", "confianza_prueba", "toma_prueba")])
#dput(names(BASE3))
library(forcats)
datos_sankey <- datos_sankey %>%
rename(VPH = resultado_vph, Sintio_toma = como_se_sintio, Barrera_emocional = barrera_emo, Barrera_tecnica = barrera_tec, Confianza = confianza_prueba, Preferencia_prueba = toma_prueba)
BASE_transformado <- datos_sankey %>% group_by(VPH, Sintio_toma, Barrera_emocional, Barrera_tecnica, Confianza, Preferencia_prueba) %>% summarise(Freq = n()) %>% ungroup()
library(easyalluvial)
alluvial_wide( data = datos_sankey,
max_variables = 8, fill_by = 'last_variable',
col_vector_flow = c("skyblue", "salmon", "#52b788"))
library(gt)
BASE_transformado %>%
arrange(desc(Freq)) %>%
gt() %>%
tab_header(title = "Conteos") %>%
cols_label_with(fn = toupper)
| Conteos | ||||||
| VPH | SINTIO_TOMA | BARRERA_EMOCIONAL | BARRERA_TECNICA | CONFIANZA | PREFERENCIA_PRUEBA | FREQ |
|---|---|---|---|---|---|---|
| No detectado | Bien | No | No | Si | Misma mujer | 297 |
| Detectado | Bien | No | No | Si | Misma mujer | 87 |
| No detectado | Bien | No | No | Si | Profesional salud | 87 |
| No detectado | Bien | No | No | Si | Sin preferencia | 81 |
| Detectado | Bien | No | No | Si | Profesional salud | 45 |
| No detectado | No Bien | Si | No | Si | Profesional salud | 34 |
| No detectado | Bien | No | No | No | Profesional salud | 32 |
| Detectado | Bien | No | No | Si | Sin preferencia | 25 |
| No detectado | Bien | Si | No | Si | Profesional salud | 17 |
| No detectado | No Bien | Si | No | No | Profesional salud | 17 |
| Detectado | Bien | No | No | No | Profesional salud | 11 |
| No detectado | Bien | No | No | No | Misma mujer | 11 |
| No detectado | Bien | Si | No | Si | Misma mujer | 11 |
| No detectado | No Bien | Si | No | Si | Misma mujer | 9 |
| Detectado | No Bien | Si | No | No | Profesional salud | 8 |
| Detectado | No Bien | Si | No | Si | Profesional salud | 7 |
| No detectado | No Bien | No | Si | Si | Profesional salud | 6 |
| No detectado | No Bien | Si | No | Si | Sin preferencia | 5 |
| No detectado | No Bien | Si | Si | No | Profesional salud | 5 |
| Detectado | Bien | Si | No | Si | Misma mujer | 4 |
| Detectado | Bien | Si | No | No | Profesional salud | 4 |
| Detectado | No Bien | Si | No | Si | Misma mujer | 4 |
| No detectado | Bien | Si | No | Si | Sin preferencia | 4 |
| No detectado | No Bien | Si | No | No | Misma mujer | 4 |
| No detectado | Bien | No | No | No | Sin preferencia | 3 |
| Detectado | Bien | No | No | No | Misma mujer | 2 |
| Detectado | Bien | Si | No | Si | Profesional salud | 2 |
| Detectado | No Bien | Si | No | Si | Sin preferencia | 2 |
| Detectado | No Bien | Si | No | No | Misma mujer | 2 |
| Detectado | No Bien | Si | Si | No | Profesional salud | 2 |
| No detectado | Bien | No | Si | Si | Misma mujer | 2 |
| No detectado | Bien | No | Si | No | Profesional salud | 2 |
| No detectado | Bien | Si | No | No | Misma mujer | 2 |
| No detectado | Bien | Si | No | No | Profesional salud | 2 |
| No detectado | No Bien | No | No | No | Profesional salud | 2 |
| No detectado | No Bien | No | Si | No | Profesional salud | 2 |
| No detectado | No Bien | Si | Si | Si | Misma mujer | 2 |
| Detectado | Bien | No | No | No | Sin preferencia | 1 |
| Detectado | Bien | Si | No | Si | Sin preferencia | 1 |
| Detectado | No Bien | No | No | Si | Profesional salud | 1 |
| Detectado | No Bien | No | No | No | Profesional salud | 1 |
| Detectado | No Bien | No | Si | No | Profesional salud | 1 |
| Detectado | No Bien | Si | Si | Si | Misma mujer | 1 |
| Detectado | No Bien | Si | Si | Si | Profesional salud | 1 |
| No detectado | Bien | Si | No | No | Sin preferencia | 1 |
| No detectado | Bien | Si | Si | Si | Misma mujer | 1 |
| No detectado | No Bien | No | No | Si | Sin preferencia | 1 |
| No detectado | No Bien | No | Si | Si | Misma mujer | 1 |
| No detectado | No Bien | Si | No | No | Sin preferencia | 1 |
| No detectado | No Bien | Si | Si | Si | Profesional salud | 1 |
tabla_acept <- BASE3 %>%
select(como_se_sintio, barrera_emo, barrera_tec,
confianza_prueba, toma_prueba, resultado_vph) %>%
tbl_summary(
by = resultado_vph,
label = list(
como_se_sintio ~ "¿Cómo se sintió?",
barrera_emo ~ "Barrera emocional",
barrera_tec ~ "Barrera técnica",
confianza_prueba ~ "Confianza en la prueba",
toma_prueba ~ "Preferencia de toma"
),
statistic = all_categorical() ~ "{n} ({p}%)",
missing = "ifany",
missing_text = "Sin dato"
) %>%
add_p(pvalue_fun = ~ style_pvalue(.x, digits = 3)) %>%
add_overall() %>%
bold_labels() %>%
bold_p(t = 0.05) %>%
modify_caption("**Aceptabilidad de autotoma según resultado VPH**")
tabla_acept
| Characteristic | Overall N = 8551 |
Detectado N = 2121 |
No detectado N = 6431 |
p-value2 |
|---|---|---|---|---|
| ¿Cómo se sintió? | 0.955 | |||
| Bien | 735 (86%) | 182 (86%) | 553 (86%) | |
| No Bien | 120 (14%) | 30 (14%) | 90 (14%) | |
| Barrera emocional | 0.970 | |||
| No | 701 (82%) | 174 (82%) | 527 (82%) | |
| Si | 154 (18%) | 38 (18%) | 116 (18%) | |
| Barrera técnica | 0.443 | |||
| No | 828 (97%) | 207 (98%) | 621 (97%) | |
| Si | 27 (3.2%) | 5 (2.4%) | 22 (3.4%) | |
| Confianza en la prueba | 0.454 | |||
| Si | 739 (86%) | 180 (85%) | 559 (87%) | |
| No | 116 (14%) | 32 (15%) | 84 (13%) | |
| Preferencia de toma | 0.178 | |||
| Misma mujer | 440 (51%) | 100 (47%) | 340 (53%) | |
| Profesional salud | 290 (34%) | 83 (39%) | 207 (32%) | |
| Sin preferencia | 125 (15%) | 29 (14%) | 96 (15%) | |
| 1 n (%) | ||||
| 2 Pearson’s Chi-squared test | ||||
library(poLCA)
set.seed(2026)
Class <- cbind(como_se_sintio, emocion, confianza, acept_autotoma) ~ 1
lca_2 <- poLCA(Class, data = BASE2, nclass = 2, nrep = 10, maxiter = 1000, graphs = FALSE)
## Model 1: llik = -1440.129 ... best llik = -1440.129
## Model 2: llik = -1440.129 ... best llik = -1440.129
## Model 3: llik = -1440.129 ... best llik = -1440.129
## Model 4: llik = -1440.129 ... best llik = -1440.129
## Model 5: llik = -1440.129 ... best llik = -1440.129
## Model 6: llik = -1440.129 ... best llik = -1440.129
## Model 7: llik = -1440.129 ... best llik = -1440.129
## Model 8: llik = -1440.129 ... best llik = -1440.129
## Model 9: llik = -1440.129 ... best llik = -1440.129
## Model 10: llik = -1440.129 ... best llik = -1440.129
## Conditional item response (column) probabilities,
## by outcome variable, for each class (row)
##
## $como_se_sintio
## Pr(1) Pr(2)
## class 1: 0.191 0.809
## class 2: 1.000 0.000
##
## $emocion
## Pr(1) Pr(2)
## class 1: 0.1309 0.8691
## class 2: 0.9646 0.0354
##
## $confianza
## Pr(1) Pr(2)
## class 1: 0.6320 0.3680
## class 2: 0.9132 0.0868
##
## $acept_autotoma
## Pr(1) Pr(2)
## class 1: 0.1856 0.8144
## class 2: 0.5843 0.4157
##
## Estimated class population shares
## 0.1733 0.8267
##
## Predicted class memberships (by modal posterior prob.)
## 0.1787 0.8213
##
## =========================================================
## Fit for 2 latent classes:
## =========================================================
## number of observations: 856
## number of estimated parameters: 9
## residual degrees of freedom: 6
## maximum log-likelihood: -1440.129
##
## AIC(2): 2898.258
## BIC(2): 2941.028
## G^2(2): 43.33248 (Likelihood ratio/deviance statistic)
## X^2(2): 42.29904 (Chi-square goodness of fit)
##
lca_3 <- poLCA(Class, data = BASE2, nclass = 3, nrep = 10, maxiter = 1000, graphs = FALSE)
## Model 1: llik = -1420.693 ... best llik = -1420.693
## Model 2: llik = -1420.692 ... best llik = -1420.692
## Model 3: llik = -1420.694 ... best llik = -1420.692
## Model 4: llik = -1419.552 ... best llik = -1419.552
## Model 5: llik = -1420.694 ... best llik = -1419.552
## Model 6: llik = -1420.694 ... best llik = -1419.552
## Model 7: llik = -1420.694 ... best llik = -1419.552
## Model 8: llik = -1420.693 ... best llik = -1419.552
## Model 9: llik = -1419.551 ... best llik = -1419.551
## Model 10: llik = -1419.552 ... best llik = -1419.551
## Conditional item response (column) probabilities,
## by outcome variable, for each class (row)
##
## $como_se_sintio
## Pr(1) Pr(2)
## class 1: 0.2111 0.7889
## class 2: 0.8729 0.1271
## class 3: 0.9975 0.0025
##
## $emocion
## Pr(1) Pr(2)
## class 1: 0.0000 1.0000
## class 2: 0.9943 0.0057
## class 3: 0.9668 0.0332
##
## $confianza
## Pr(1) Pr(2)
## class 1: 0.6448 0.3552
## class 2: 0.5419 0.4581
## class 3: 0.9677 0.0323
##
## $acept_autotoma
## Pr(1) Pr(2)
## class 1: 0.2037 0.7963
## class 2: 0.0000 1.0000
## class 3: 0.6716 0.3284
##
## Estimated class population shares
## 0.1553 0.1247 0.72
##
## Predicted class memberships (by modal posterior prob.)
## 0.1612 0.0736 0.7652
##
## =========================================================
## Fit for 3 latent classes:
## =========================================================
## number of observations: 856
## number of estimated parameters: 14
## residual degrees of freedom: 1
## maximum log-likelihood: -1419.551
##
## AIC(3): 2867.102
## BIC(3): 2933.634
## G^2(3): 2.17683 (Likelihood ratio/deviance statistic)
## X^2(3): 2.134992 (Chi-square goodness of fit)
##
## ALERT: iterations finished, MAXIMUM LIKELIHOOD NOT FOUND
##
library(tidyverse)
fit_table <- tibble(
Clases = c(2, 3),
LogLik = c(lca_2$llik, lca_3$llik),
AIC = c(lca_2$aic, lca_3$aic),
BIC = c(lca_2$bic, lca_3$bic),
Chi2 = c(lca_2$Chisq, lca_3$Chisq),
df = c(lca_2$resid.df, lca_3$resid.df)
)
library(gt)
fit_table %>%
gt() %>%
tab_header(title = "Comparación de modelos LCA") %>%
fmt_number(columns = c(LogLik, AIC, BIC, Chi2), decimals = 2) %>%
tab_footnote("Menor AIC y BIC indican mejor ajuste. Entropía > 0.80 indica buena clasificación.") %>%
tab_style(
style = cell_fill(color = "#d4edda"),
locations = cells_body(rows = BIC == min(BIC))
)
| Comparación de modelos LCA | |||||
| Clases | LogLik | AIC | BIC | Chi2 | df |
|---|---|---|---|---|---|
| 2 | −1,440.13 | 2,898.26 | 2,941.03 | 42.30 | 6 |
| 3 | −1,419.55 | 2,867.10 | 2,933.63 | 2.13 | 1 |
| Menor AIC y BIC indican mejor ajuste. Entropía > 0.80 indica buena clasificación. | |||||