Untitled
2025-10-04
#Buscamos proporción de datos faltantes (Missingvalues)---- NO SE DETECTAN DATOS FALTANTES
#is.na(sad[ , 25:34])
#which(is.na(sad[ , 25:34]))
#which(is.na(sad[ , 42]))#pasar variables categóricas a factor
names(sad)<- tolower(names (sad))
sad <- sad %>%
mutate(across(19:35, as.factor))
#crear cicat_factor
sad$cicat_factor <- factor(sad$cicat, levels = c(0, 1), labels = c( "No Cicatrizo", "Cicatrizo"))
table(sad$cicat_factor)##
## No Cicatrizo Cicatrizo
## 174 232#paso variables chr a factor
sad <- sad %>% mutate_if(is.character, as.factor)
#crear ampu_factor
sad$ampu_factor <- factor(sad$ampu, levels = c(0, 1), labels = c( "No Amputado", "Amputado"))
table(sad$ampu_factor)##
## No Amputado Amputado
## 348 58#crear fallecio_factor
sad$fallecio_factor <- factor(sad$fallecio, levels = c(0, 1), labels = c( "No Fallecido", "Fallecido"))
table(sad$fallecio_factor)##
## No Fallecido Fallecido
## 392 14#variable unificada
sad <- sad %>%
mutate(
estado = tolower(trimws(estado))
)
table(sad$estado)##
## ampu cerro fallecio persiste
## 58 232 14 102sad <- sad %>%
mutate(
estado_unif = case_when(
estado %in% c("ampu") ~ "amputacion",
estado %in% c("cerro") ~ "cerro_herida",
estado %in% c("murio", "fallecio") ~ "fallecio",
estado %in% c("persiste") ~ "persiste",
estado %in% c("perdido") ~ "perdido_seguimiento",
TRUE ~ NA_character_
)
)
table(sad$estado_unif)##
## amputacion cerro_herida fallecio persiste
## 58 232 14 102sad <- sad %>%
mutate(
evento_ampu = case_when(
estado_unif %in% c("cerro_herida","fallecio") ~ 2,
estado_unif %in% c("amputacion") ~ 1,
estado_unif %in% c("perdido_seguimiento", "persiste") ~ 0,
TRUE ~ NA_real_
)
)
table(sad$evento_ampu)##
## 0 1 2
## 102 58 246tiempo_evento_am <- Surv(sad$dias_out,sad$ampu)library(tableone)
#Las variables van directamente sin el nombre del dataframe
catvars = c( "sexo","tbq", "iam", "hta", "dial", "local","topo", "inf", "ede", "neu", "prof", "area", "zon", "fase", "ant_pie", "am_cont")
#continuas
vars = c("dias_ev","san_elian", "sexo","tbq", "iam", "hta", "dial", "local","topo", "inf", "ede", "neu", "prof", "area", "zon", "fase", "edad", "ant_pie", "am_cont" )
tabla_1 <- CreateTableOne(vars = vars, factorVars = catvars, data = sad)
print(tabla_1)##
## Overall
## n 406
## dias_ev (mean (SD)) 50.26 (112.85)
## san_elian (mean (SD)) 18.26 (4.09)
## sexo = M (%) 324 (79.8)
## tbq (%)
## 0 169 (41.8)
## 1 80 (19.8)
## 2 155 (38.4)
## iam = 1 (%) 83 (20.5)
## hta = 1 (%) 273 (67.4)
## dial = 1 (%) 19 ( 4.7)
## local (%)
## 1 196 (48.3)
## 2 134 (33.0)
## 3 76 (18.7)
## topo (%)
## 1 152 (37.4)
## 2 104 (25.6)
## 3 150 (36.9)
## inf (%)
## 0 86 (21.2)
## 1 32 ( 7.9)
## 2 223 (54.9)
## 3 65 (16.0)
## ede (%)
## 0 143 (35.2)
## 1 98 (24.1)
## 2 124 (30.5)
## 3 41 (10.1)
## neu (%)
## 0 4 ( 1.0)
## 1 25 ( 6.2)
## 2 341 (84.0)
## 3 36 ( 8.9)
## prof (%)
## 1 32 ( 7.9)
## 2 133 (32.8)
## 3 241 (59.4)
## area (%)
## 1 188 (46.3)
## 2 157 (38.7)
## 3 61 (15.0)
## zon (%)
## 1 242 (59.6)
## 2 90 (22.2)
## 3 74 (18.2)
## fase (%)
## 1 19 ( 4.7)
## 2 51 (12.6)
## 3 336 (82.8)
## edad (mean (SD)) 57.41 (11.09)
## ant_pie (%)
## 0 142 (35.0)
## 1 263 (64.8)
## 2 1 ( 0.2)
## am_cont = 1 (%) 23 ( 5.7)kableone(tabla_1)| Overall | |
|---|---|
| n | 406 |
| dias_ev (mean (SD)) | 50.26 (112.85) |
| san_elian (mean (SD)) | 18.26 (4.09) |
| sexo = M (%) | 324 (79.8) |
| tbq (%) | |
| 0 | 169 (41.8) |
| 1 | 80 (19.8) |
| 2 | 155 (38.4) |
| iam = 1 (%) | 83 (20.5) |
| hta = 1 (%) | 273 (67.4) |
| dial = 1 (%) | 19 ( 4.7) |
| local (%) | |
| 1 | 196 (48.3) |
| 2 | 134 (33.0) |
| 3 | 76 (18.7) |
| topo (%) | |
| 1 | 152 (37.4) |
| 2 | 104 (25.6) |
| 3 | 150 (36.9) |
| inf (%) | |
| 0 | 86 (21.2) |
| 1 | 32 ( 7.9) |
| 2 | 223 (54.9) |
| 3 | 65 (16.0) |
| ede (%) | |
| 0 | 143 (35.2) |
| 1 | 98 (24.1) |
| 2 | 124 (30.5) |
| 3 | 41 (10.1) |
| neu (%) | |
| 0 | 4 ( 1.0) |
| 1 | 25 ( 6.2) |
| 2 | 341 (84.0) |
| 3 | 36 ( 8.9) |
| prof (%) | |
| 1 | 32 ( 7.9) |
| 2 | 133 (32.8) |
| 3 | 241 (59.4) |
| area (%) | |
| 1 | 188 (46.3) |
| 2 | 157 (38.7) |
| 3 | 61 (15.0) |
| zon (%) | |
| 1 | 242 (59.6) |
| 2 | 90 (22.2) |
| 3 | 74 (18.2) |
| fase (%) | |
| 1 | 19 ( 4.7) |
| 2 | 51 (12.6) |
| 3 | 336 (82.8) |
| edad (mean (SD)) | 57.41 (11.09) |
| ant_pie (%) | |
| 0 | 142 (35.0) |
| 1 | 263 (64.8) |
| 2 | 1 ( 0.2) |
| am_cont = 1 (%) | 23 ( 5.7) |
tabla_ampu <- CreateTableOne(
vars = vars,
strata = "ampu",
factorVars = catvars,
data = sad,
addOverall = TRUE
)
print(tabla_ampu)## Stratified by ampu
## Overall 0 1 p
## n 406 348 58
## dias_ev (mean (SD)) 50.26 (112.85) 49.01 (111.88) 57.74 (119.19) 0.586
## san_elian (mean (SD)) 18.26 (4.09) 17.55 (3.78) 22.52 (3.22) <0.001
## sexo = M (%) 324 (79.8) 279 (80.2) 45 (77.6) 0.781
## tbq (%) 0.974
## 0 169 (41.8) 144 (41.6) 25 (43.1)
## 1 80 (19.8) 69 (19.9) 11 (19.0)
## 2 155 (38.4) 133 (38.4) 22 (37.9)
## iam = 1 (%) 83 (20.5) 70 (20.2) 13 (22.4) 0.837
## hta = 1 (%) 273 (67.4) 228 (65.7) 45 (77.6) 0.102
## dial = 1 (%) 19 ( 4.7) 16 ( 4.6) 3 ( 5.2) 1.000
## local (%) 0.011
## 1 196 (48.3) 171 (49.1) 25 (43.1)
## 2 134 (33.0) 120 (34.5) 14 (24.1)
## 3 76 (18.7) 57 (16.4) 19 (32.8)
## topo (%) <0.001
## 1 152 (37.4) 146 (42.0) 6 (10.3)
## 2 104 (25.6) 97 (27.9) 7 (12.1)
## 3 150 (36.9) 105 (30.2) 45 (77.6)
## inf (%) <0.001
## 0 86 (21.2) 80 (23.0) 6 (10.3)
## 1 32 ( 7.9) 28 ( 8.0) 4 ( 6.9)
## 2 223 (54.9) 197 (56.6) 26 (44.8)
## 3 65 (16.0) 43 (12.4) 22 (37.9)
## ede (%) 0.163
## 0 143 (35.2) 121 (34.8) 22 (37.9)
## 1 98 (24.1) 90 (25.9) 8 (13.8)
## 2 124 (30.5) 101 (29.0) 23 (39.7)
## 3 41 (10.1) 36 (10.3) 5 ( 8.6)
## neu (%) 0.717
## 0 4 ( 1.0) 3 ( 0.9) 1 ( 1.7)
## 1 25 ( 6.2) 20 ( 5.7) 5 ( 8.6)
## 2 341 (84.0) 293 (84.2) 48 (82.8)
## 3 36 ( 8.9) 32 ( 9.2) 4 ( 6.9)
## prof (%) 0.003
## 1 32 ( 7.9) 31 ( 8.9) 1 ( 1.7)
## 2 133 (32.8) 122 (35.1) 11 (19.0)
## 3 241 (59.4) 195 (56.0) 46 (79.3)
## area (%) <0.001
## 1 188 (46.3) 178 (51.1) 10 (17.2)
## 2 157 (38.7) 138 (39.7) 19 (32.8)
## 3 61 (15.0) 32 ( 9.2) 29 (50.0)
## zon (%) <0.001
## 1 242 (59.6) 230 (66.1) 12 (20.7)
## 2 90 (22.2) 78 (22.4) 12 (20.7)
## 3 74 (18.2) 40 (11.5) 34 (58.6)
## fase (%) 0.010
## 1 19 ( 4.7) 19 ( 5.5) 0 ( 0.0)
## 2 51 (12.6) 49 (14.1) 2 ( 3.4)
## 3 336 (82.8) 280 (80.5) 56 (96.6)
## edad (mean (SD)) 57.41 (11.09) 56.54 (11.23) 62.64 (8.57) <0.001
## ant_pie (%) 0.068
## 0 142 (35.0) 114 (32.8) 28 (48.3)
## 1 263 (64.8) 233 (67.0) 30 (51.7)
## 2 1 ( 0.2) 1 ( 0.3) 0 ( 0.0)
## am_cont = 1 (%) 23 ( 5.7) 17 ( 4.9) 6 (10.3) 0.174
## Stratified by ampu
## test
## n
## dias_ev (mean (SD))
## san_elian (mean (SD))
## sexo = M (%)
## tbq (%)
## 0
## 1
## 2
## iam = 1 (%)
## hta = 1 (%)
## dial = 1 (%)
## local (%)
## 1
## 2
## 3
## topo (%)
## 1
## 2
## 3
## inf (%)
## 0
## 1
## 2
## 3
## ede (%)
## 0
## 1
## 2
## 3
## neu (%)
## 0
## 1
## 2
## 3
## prof (%)
## 1
## 2
## 3
## area (%)
## 1
## 2
## 3
## zon (%)
## 1
## 2
## 3
## fase (%)
## 1
## 2
## 3
## edad (mean (SD))
## ant_pie (%)
## 0
## 1
## 2
## am_cont = 1 (%)kableone(tabla_ampu)| Overall | 0 | 1 | p | test | |
|---|---|---|---|---|---|
| n | 406 | 348 | 58 | ||
| dias_ev (mean (SD)) | 50.26 (112.85) | 49.01 (111.88) | 57.74 (119.19) | 0.586 | |
| san_elian (mean (SD)) | 18.26 (4.09) | 17.55 (3.78) | 22.52 (3.22) | <0.001 | |
| sexo = M (%) | 324 (79.8) | 279 (80.2) | 45 (77.6) | 0.781 | |
| tbq (%) | 0.974 | ||||
| 0 | 169 (41.8) | 144 (41.6) | 25 (43.1) | ||
| 1 | 80 (19.8) | 69 (19.9) | 11 (19.0) | ||
| 2 | 155 (38.4) | 133 (38.4) | 22 (37.9) | ||
| iam = 1 (%) | 83 (20.5) | 70 (20.2) | 13 (22.4) | 0.837 | |
| hta = 1 (%) | 273 (67.4) | 228 (65.7) | 45 (77.6) | 0.102 | |
| dial = 1 (%) | 19 ( 4.7) | 16 ( 4.6) | 3 ( 5.2) | 1.000 | |
| local (%) | 0.011 | ||||
| 1 | 196 (48.3) | 171 (49.1) | 25 (43.1) | ||
| 2 | 134 (33.0) | 120 (34.5) | 14 (24.1) | ||
| 3 | 76 (18.7) | 57 (16.4) | 19 (32.8) | ||
| topo (%) | <0.001 | ||||
| 1 | 152 (37.4) | 146 (42.0) | 6 (10.3) | ||
| 2 | 104 (25.6) | 97 (27.9) | 7 (12.1) | ||
| 3 | 150 (36.9) | 105 (30.2) | 45 (77.6) | ||
| inf (%) | <0.001 | ||||
| 0 | 86 (21.2) | 80 (23.0) | 6 (10.3) | ||
| 1 | 32 ( 7.9) | 28 ( 8.0) | 4 ( 6.9) | ||
| 2 | 223 (54.9) | 197 (56.6) | 26 (44.8) | ||
| 3 | 65 (16.0) | 43 (12.4) | 22 (37.9) | ||
| ede (%) | 0.163 | ||||
| 0 | 143 (35.2) | 121 (34.8) | 22 (37.9) | ||
| 1 | 98 (24.1) | 90 (25.9) | 8 (13.8) | ||
| 2 | 124 (30.5) | 101 (29.0) | 23 (39.7) | ||
| 3 | 41 (10.1) | 36 (10.3) | 5 ( 8.6) | ||
| neu (%) | 0.717 | ||||
| 0 | 4 ( 1.0) | 3 ( 0.9) | 1 ( 1.7) | ||
| 1 | 25 ( 6.2) | 20 ( 5.7) | 5 ( 8.6) | ||
| 2 | 341 (84.0) | 293 (84.2) | 48 (82.8) | ||
| 3 | 36 ( 8.9) | 32 ( 9.2) | 4 ( 6.9) | ||
| prof (%) | 0.003 | ||||
| 1 | 32 ( 7.9) | 31 ( 8.9) | 1 ( 1.7) | ||
| 2 | 133 (32.8) | 122 (35.1) | 11 (19.0) | ||
| 3 | 241 (59.4) | 195 (56.0) | 46 (79.3) | ||
| area (%) | <0.001 | ||||
| 1 | 188 (46.3) | 178 (51.1) | 10 (17.2) | ||
| 2 | 157 (38.7) | 138 (39.7) | 19 (32.8) | ||
| 3 | 61 (15.0) | 32 ( 9.2) | 29 (50.0) | ||
| zon (%) | <0.001 | ||||
| 1 | 242 (59.6) | 230 (66.1) | 12 (20.7) | ||
| 2 | 90 (22.2) | 78 (22.4) | 12 (20.7) | ||
| 3 | 74 (18.2) | 40 (11.5) | 34 (58.6) | ||
| fase (%) | 0.010 | ||||
| 1 | 19 ( 4.7) | 19 ( 5.5) | 0 ( 0.0) | ||
| 2 | 51 (12.6) | 49 (14.1) | 2 ( 3.4) | ||
| 3 | 336 (82.8) | 280 (80.5) | 56 (96.6) | ||
| edad (mean (SD)) | 57.41 (11.09) | 56.54 (11.23) | 62.64 (8.57) | <0.001 | |
| ant_pie (%) | 0.068 | ||||
| 0 | 142 (35.0) | 114 (32.8) | 28 (48.3) | ||
| 1 | 263 (64.8) | 233 (67.0) | 30 (51.7) | ||
| 2 | 1 ( 0.2) | 1 ( 0.3) | 0 ( 0.0) | ||
| am_cont = 1 (%) | 23 ( 5.7) | 17 ( 4.9) | 6 (10.3) | 0.174 |
print(tabla_ampu,
nonnormal = c("dias_ev", "san_elian"),
showAllLevels = TRUE,
quote = FALSE,
noSpaces = TRUE)## Stratified by ampu
## level Overall 0
## n 406 348
## dias_ev (median [IQR]) 20.00 [10.00, 45.00] 20.00 [10.00, 30.00]
## san_elian (median [IQR]) 18.00 [15.00, 21.00] 18.00 [15.00, 20.00]
## sexo (%) F 82 (20.2) 69 (19.8)
## M 324 (79.8) 279 (80.2)
## tbq (%) 0 169 (41.8) 144 (41.6)
## 1 80 (19.8) 69 (19.9)
## 2 155 (38.4) 133 (38.4)
## iam (%) 0 321 (79.5) 276 (79.8)
## 1 83 (20.5) 70 (20.2)
## hta (%) 0 132 (32.6) 119 (34.3)
## 1 273 (67.4) 228 (65.7)
## dial (%) 0 387 (95.3) 332 (95.4)
## 1 19 (4.7) 16 (4.6)
## local (%) 1 196 (48.3) 171 (49.1)
## 2 134 (33.0) 120 (34.5)
## 3 76 (18.7) 57 (16.4)
## topo (%) 1 152 (37.4) 146 (42.0)
## 2 104 (25.6) 97 (27.9)
## 3 150 (36.9) 105 (30.2)
## inf (%) 0 86 (21.2) 80 (23.0)
## 1 32 (7.9) 28 (8.0)
## 2 223 (54.9) 197 (56.6)
## 3 65 (16.0) 43 (12.4)
## ede (%) 0 143 (35.2) 121 (34.8)
## 1 98 (24.1) 90 (25.9)
## 2 124 (30.5) 101 (29.0)
## 3 41 (10.1) 36 (10.3)
## neu (%) 0 4 (1.0) 3 (0.9)
## 1 25 (6.2) 20 (5.7)
## 2 341 (84.0) 293 (84.2)
## 3 36 (8.9) 32 (9.2)
## prof (%) 1 32 (7.9) 31 (8.9)
## 2 133 (32.8) 122 (35.1)
## 3 241 (59.4) 195 (56.0)
## area (%) 1 188 (46.3) 178 (51.1)
## 2 157 (38.7) 138 (39.7)
## 3 61 (15.0) 32 (9.2)
## zon (%) 1 242 (59.6) 230 (66.1)
## 2 90 (22.2) 78 (22.4)
## 3 74 (18.2) 40 (11.5)
## fase (%) 1 19 (4.7) 19 (5.5)
## 2 51 (12.6) 49 (14.1)
## 3 336 (82.8) 280 (80.5)
## edad (mean (SD)) 57.41 (11.09) 56.54 (11.23)
## ant_pie (%) 0 142 (35.0) 114 (32.8)
## 1 263 (64.8) 233 (67.0)
## 2 1 (0.2) 1 (0.3)
## am_cont (%) 0 383 (94.3) 331 (95.1)
## 1 23 (5.7) 17 (4.9)
## Stratified by ampu
## 1 p test
## n 58
## dias_ev (median [IQR]) 30.00 [15.00, 60.00] 0.011 nonnorm
## san_elian (median [IQR]) 23.00 [20.25, 24.00] <0.001 nonnorm
## sexo (%) 13 (22.4) 0.781
## 45 (77.6)
## tbq (%) 25 (43.1) 0.974
## 11 (19.0)
## 22 (37.9)
## iam (%) 45 (77.6) 0.837
## 13 (22.4)
## hta (%) 13 (22.4) 0.102
## 45 (77.6)
## dial (%) 55 (94.8) 1.000
## 3 (5.2)
## local (%) 25 (43.1) 0.011
## 14 (24.1)
## 19 (32.8)
## topo (%) 6 (10.3) <0.001
## 7 (12.1)
## 45 (77.6)
## inf (%) 6 (10.3) <0.001
## 4 (6.9)
## 26 (44.8)
## 22 (37.9)
## ede (%) 22 (37.9) 0.163
## 8 (13.8)
## 23 (39.7)
## 5 (8.6)
## neu (%) 1 (1.7) 0.717
## 5 (8.6)
## 48 (82.8)
## 4 (6.9)
## prof (%) 1 (1.7) 0.003
## 11 (19.0)
## 46 (79.3)
## area (%) 10 (17.2) <0.001
## 19 (32.8)
## 29 (50.0)
## zon (%) 12 (20.7) <0.001
## 12 (20.7)
## 34 (58.6)
## fase (%) 0 (0.0) 0.010
## 2 (3.4)
## 56 (96.6)
## edad (mean (SD)) 62.64 (8.57) <0.001
## ant_pie (%) 28 (48.3) 0.068
## 30 (51.7)
## 0 (0.0)
## am_cont (%) 52 (89.7) 0.174
## 6 (10.3)kableone(tabla_ampu)| Overall | 0 | 1 | p | test | |
|---|---|---|---|---|---|
| n | 406 | 348 | 58 | ||
| dias_ev (mean (SD)) | 50.26 (112.85) | 49.01 (111.88) | 57.74 (119.19) | 0.586 | |
| san_elian (mean (SD)) | 18.26 (4.09) | 17.55 (3.78) | 22.52 (3.22) | <0.001 | |
| sexo = M (%) | 324 (79.8) | 279 (80.2) | 45 (77.6) | 0.781 | |
| tbq (%) | 0.974 | ||||
| 0 | 169 (41.8) | 144 (41.6) | 25 (43.1) | ||
| 1 | 80 (19.8) | 69 (19.9) | 11 (19.0) | ||
| 2 | 155 (38.4) | 133 (38.4) | 22 (37.9) | ||
| iam = 1 (%) | 83 (20.5) | 70 (20.2) | 13 (22.4) | 0.837 | |
| hta = 1 (%) | 273 (67.4) | 228 (65.7) | 45 (77.6) | 0.102 | |
| dial = 1 (%) | 19 ( 4.7) | 16 ( 4.6) | 3 ( 5.2) | 1.000 | |
| local (%) | 0.011 | ||||
| 1 | 196 (48.3) | 171 (49.1) | 25 (43.1) | ||
| 2 | 134 (33.0) | 120 (34.5) | 14 (24.1) | ||
| 3 | 76 (18.7) | 57 (16.4) | 19 (32.8) | ||
| topo (%) | <0.001 | ||||
| 1 | 152 (37.4) | 146 (42.0) | 6 (10.3) | ||
| 2 | 104 (25.6) | 97 (27.9) | 7 (12.1) | ||
| 3 | 150 (36.9) | 105 (30.2) | 45 (77.6) | ||
| inf (%) | <0.001 | ||||
| 0 | 86 (21.2) | 80 (23.0) | 6 (10.3) | ||
| 1 | 32 ( 7.9) | 28 ( 8.0) | 4 ( 6.9) | ||
| 2 | 223 (54.9) | 197 (56.6) | 26 (44.8) | ||
| 3 | 65 (16.0) | 43 (12.4) | 22 (37.9) | ||
| ede (%) | 0.163 | ||||
| 0 | 143 (35.2) | 121 (34.8) | 22 (37.9) | ||
| 1 | 98 (24.1) | 90 (25.9) | 8 (13.8) | ||
| 2 | 124 (30.5) | 101 (29.0) | 23 (39.7) | ||
| 3 | 41 (10.1) | 36 (10.3) | 5 ( 8.6) | ||
| neu (%) | 0.717 | ||||
| 0 | 4 ( 1.0) | 3 ( 0.9) | 1 ( 1.7) | ||
| 1 | 25 ( 6.2) | 20 ( 5.7) | 5 ( 8.6) | ||
| 2 | 341 (84.0) | 293 (84.2) | 48 (82.8) | ||
| 3 | 36 ( 8.9) | 32 ( 9.2) | 4 ( 6.9) | ||
| prof (%) | 0.003 | ||||
| 1 | 32 ( 7.9) | 31 ( 8.9) | 1 ( 1.7) | ||
| 2 | 133 (32.8) | 122 (35.1) | 11 (19.0) | ||
| 3 | 241 (59.4) | 195 (56.0) | 46 (79.3) | ||
| area (%) | <0.001 | ||||
| 1 | 188 (46.3) | 178 (51.1) | 10 (17.2) | ||
| 2 | 157 (38.7) | 138 (39.7) | 19 (32.8) | ||
| 3 | 61 (15.0) | 32 ( 9.2) | 29 (50.0) | ||
| zon (%) | <0.001 | ||||
| 1 | 242 (59.6) | 230 (66.1) | 12 (20.7) | ||
| 2 | 90 (22.2) | 78 (22.4) | 12 (20.7) | ||
| 3 | 74 (18.2) | 40 (11.5) | 34 (58.6) | ||
| fase (%) | 0.010 | ||||
| 1 | 19 ( 4.7) | 19 ( 5.5) | 0 ( 0.0) | ||
| 2 | 51 (12.6) | 49 (14.1) | 2 ( 3.4) | ||
| 3 | 336 (82.8) | 280 (80.5) | 56 (96.6) | ||
| edad (mean (SD)) | 57.41 (11.09) | 56.54 (11.23) | 62.64 (8.57) | <0.001 | |
| ant_pie (%) | 0.068 | ||||
| 0 | 142 (35.0) | 114 (32.8) | 28 (48.3) | ||
| 1 | 263 (64.8) | 233 (67.0) | 30 (51.7) | ||
| 2 | 1 ( 0.2) | 1 ( 0.3) | 0 ( 0.0) | ||
| am_cont = 1 (%) | 23 ( 5.7) | 17 ( 4.9) | 6 (10.3) | 0.174 |
wilcox.test(sad$dias_ev, sad$ampu)##
## Wilcoxon rank sum test with continuity correction
##
## data: sad$dias_ev and sad$ampu
## W = 163502, p-value <0.0000000000000002
## alternative hypothesis: true location shift is not equal to 0wilcox.test(sad$san_elian, sad$ampu)##
## Wilcoxon rank sum test with continuity correction
##
## data: sad$san_elian and sad$ampu
## W = 164836, p-value <0.0000000000000002
## alternative hypothesis: true location shift is not equal to 0#variables categorizadas
#Las variables van directamente sin el nombre del dataframe#evaluacion de la cantidad de eventos por tiempo de seguimiento a intervalos de a 30 dias
sad %>%
filter(ampu == 1) %>%
mutate(intervalo = cut(dias_out,
breaks = seq(0, max(dias_out), by = 30))) %>%
count(intervalo) %>%
ggplot(aes(x = intervalo, y = n)) +
geom_col(fill = "darkorange") +
labs(x = "Intervalos de tiempo (30 dias)",
y = "Eventos",
title = "Eventos por intervalo de seguimiento") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
#eventos según valoren la escala de San Elian
sad %>%
filter(ampu == 1) %>%
mutate(intervalo = cut(san_elian,
breaks = seq(0, max(san_elian), by = 1))) %>%
count(intervalo) %>%
ggplot(aes(x = intervalo, y = n)) +
geom_col(fill = "darkorange") +
labs(x = "Escala de San Elian",
y = "Eventos",
title = "Eventos por escala") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
#diferencia en el valor de la escala en pacientes cicatrizados, amputados y fallecidos
table(sad$estado_unif)##
## amputacion cerro_herida fallecio persiste
## 58 232 14 102describeBy(sad$san_elian, sad$estado_unif)##
## Descriptive statistics by group
## group: amputacion
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 58 22.5 3.22 23 22.6 2.97 15 29 14 -0.29 -0.33 0.42
## ------------------------------------------------------------
## group: cerro_herida
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 232 17.2 3.76 17 17.2 4.45 8 26 18 -0.19 -0.5 0.25
## ------------------------------------------------------------
## group: fallecio
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 14 19.3 4.75 20.5 19.3 5.19 11 27 16 -0.17 -1.32 1.27
## ------------------------------------------------------------
## group: persiste
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 102 18.2 3.58 18 18.3 2.97 9 26 17 -0.28 -0.28 0.35hist(sad$san_elian[sad$estado_unif=="amputacion"])
hist(sad$san_elian[sad$estado_unif=="cerro_herida"])
hist(sad$san_elian[sad$estado_unif=="fallecio"])
hist(sad$san_elian[sad$estado_unif=="persiste"])
shapiro.test(sad$san_elian[sad$estado_unif=="amputacion"])##
## Shapiro-Wilk normality test
##
## data: sad$san_elian[sad$estado_unif == "amputacion"]
## W = 1, p-value = 0.2shapiro.test(sad$san_elian[sad$estado_unif=="cerro_herida"])##
## Shapiro-Wilk normality test
##
## data: sad$san_elian[sad$estado_unif == "cerro_herida"]
## W = 1, p-value = 0.01shapiro.test(sad$san_elian[sad$estado_unif=="fallecio"])##
## Shapiro-Wilk normality test
##
## data: sad$san_elian[sad$estado_unif == "fallecio"]
## W = 0.9, p-value = 0.5shapiro.test(sad$san_elian[sad$estado_unif=="persiste"])##
## Shapiro-Wilk normality test
##
## data: sad$san_elian[sad$estado_unif == "persiste"]
## W = 1, p-value = 0.1qqnorm(sad$san_elian[sad$estado_unif=="amputacion"])
qqnorm(sad$san_elian[sad$estado_unif=="cerro_herida"])
qqnorm(sad$san_elian[sad$estado_unif=="fallecio"])
qqnorm(sad$san_elian[sad$estado_unif=="persiste"])
kruskal.test(sad$san_elian, sad$estado_unif)##
## Kruskal-Wallis rank sum test
##
## data: sad$san_elian and sad$estado_unif
## Kruskal-Wallis chi-squared = 76, df = 3, p-value <0.0000000000000002pairwise.wilcox.test(sad$san_elian, sad$estado_unif,"bonferroni")##
## Pairwise comparisons using Wilcoxon rank sum test with continuity correction
##
## data: sad$san_elian and sad$estado_unif
##
## amputacion cerro_herida fallecio
## cerro_herida < 0.0000000000000002 - -
## fallecio 0.1 0.5 -
## persiste 0.0000000001 0.1 1.0
##
## P value adjustment method: bonferronitapply(sad$san_elian,
sad$estado_unif,
quantile,
probs = c(0.25, 0.75),
na.rm = TRUE)## $amputacion
## 25% 75%
## 20.2 24.0
##
## $cerro_herida
## 25% 75%
## 15 20
##
## $fallecio
## 25% 75%
## 15 22
##
## $persiste
## 25% 75%
## 16 20cox_model_san <- coxph(tiempo_evento_am ~ san_elian + edad, data = sad)
summary(cox_model_san)## Call:
## coxph(formula = tiempo_evento_am ~ san_elian + edad, data = sad)
##
## n= 406, number of events= 58
##
## coef exp(coef) se(coef) z Pr(>|z|)
## san_elian 0.3424 1.4083 0.0416 8.24 <0.0000000000000002 ***
## edad 0.0361 1.0368 0.0137 2.64 0.0083 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## san_elian 1.41 0.710 1.30 1.53
## edad 1.04 0.965 1.01 1.06
##
## Concordance= 0.845 (se = 0.022 )
## Likelihood ratio test= 93.7 on 2 df, p=<0.0000000000000002
## Wald test = 82.1 on 2 df, p=<0.0000000000000002
## Score (logrank) test = 85.9 on 2 df, p=<0.0000000000000002tab_model(cox_model_san,
show.ci = TRUE, # Mostrar intervalo de confianza
show.se = TRUE, # Mostrar error estándar
transform = "exp", # Para mostrar HR en lugar de log(HR)
dv.labels = "Modelo de Cox - Datos drugs")| Modelo de Cox - Datos drugs | ||||
|---|---|---|---|---|
| Predictors | Estimates | std. Error | CI | p |
| san elian | 1.41 | 0.06 | 0.00 – Inf | <0.001 |
| edad | 1.04 | 0.01 | 0.00 – Inf | 0.008 |
| Observations | 406 | |||
| R2 Nagelkerke | 0.257 | |||
confint(cox_model_san)## 2.5 % 97.5 %
## san_elian 0.26091 0.4238
## edad 0.00928 0.0629#Test de proporcionalidad de riesgos de Schoenfeld
test_ph <- cox.zph(cox_model_san)
test_ph## chisq df p
## san_elian 3.579 1 0.059
## edad 0.315 1 0.575
## GLOBAL 4.153 2 0.125#Gráficos por variables de residuos de Schoenfeld en el tiempo
ggcoxzph(test_ph)## Warning: ggtheme is not a valid theme.
## Please use `theme()` to construct themes.
## ggtheme is not a valid theme.
## Please use `theme()` to construct themes.
CUMPLE EL SUPUESTO DE RIESGOS PROPORCIONALES
# Martingale y Deviance
ggcoxdiagnostics(cox_model_san, type = "martingale", linear.predictions = FALSE)## `geom_smooth()` using formula = 'y ~ x'
ggcoxdiagnostics(cox_model_san, type = "deviance", linear.predictions = FALSE)## `geom_smooth()` using formula = 'y ~ x'
ADECUACION EVALUADA POR DEVIANCE Y MARTINGALE
#C de harrell
summary(cox_model_san)$concordance## C se(C)
## 0.8451 0.0221#bootstrap
#Indice D de Somer (Dxy), Metricas Q D U, slope y R2 con paquete `rms`
dd <- datadist(sad); options(datadist = "dd")## Warning in datadist(sad): ...44 is constantcox_metrics <- cph(tiempo_evento_am ~ san_elian,
data = sad, x = TRUE, y = TRUE, surv = TRUE)
# Validamos métricas con bootstrap
invisible(capture.output({ #esta linea es para evitar que se imprima el output del bootstrap)
cox_rms <- validate(cox_metrics, method = "boot", B = 200)
}))
# Revisamos las metricas validadas (en la columna index.corrected)
cox_rms## index.orig training test optimism index.corrected n
## Dxy 0.668 0.6733 0.6677 0.0055 0.6622 200
## R2 0.240 0.2455 0.2396 0.0059 0.2337 200
## Slope 1.000 1.0000 0.9840 0.0160 0.9840 200
## D 0.130 0.1355 0.1303 0.0052 0.1251 200
## U -0.003 -0.0031 0.0013 -0.0045 0.0014 200
## Q 0.133 0.1386 0.1290 0.0097 0.1237 200
## g 1.673 1.7166 1.6727 0.0438 1.6289 200round(cox_rms[,5],2)## Dxy R2 Slope D U Q g
## 0.66 0.23 0.98 0.13 0.00 0.12 1.63Métrica Qué evalúa Valores esperados Interpretación C (Harrell’s C) Discriminación 0.5 – 1 0.5 = azar; 0.7 aceptable; 0.8–0.9 muy buena; >0.9 excelente Dxy Discriminación (derivado de C) 0 – 1 0 = sin discriminación; 1 = perfecta R² Ajuste global (variabilidad explicada) 0 – 1 0.1–0.2 habitual; 0.2–0.4 bueno; >0.4 alto (raro en clínica) Slope (calibración) Calibración ≈ 1 1 = ideal; <1 sobreajuste; >1 subajuste Intercepto (si se reporta) Calibración ≈ 0 0 = sin sesgo sistemático D Separación de riesgos >0 mayor = mejor separación U Error de calibración ≈ 0 0 = ideal Q Calidad global (D – U) >0 mayor = mejor modelo g (Gini / log-risk separation) Separación del riesgo >0 mayor = mejor discriminación Brier score Error de predicción 0 – 1 menor = mejor; <0.25 aceptable; <0.10 bueno AUC (ROC) Discriminación 0.5 – 1 igual que C: >0.8 bueno Optimismo (bootstrap) Sobreajuste ≈ 0 cercano a 0 = modelo estable
# Calibración. Evaluación Gráfica
invisible(capture.output({cal <- calibrate(cox_metrics, method = "boot", B = 200, u = 120)
}))
plot(cal,
subtitles = FALSE,
lwd = 2,
lty = 2,
main = "Grafico de Calibracion",
cex.main = 1.2,
cex.axis = 0.7,
cex.lab = 0.9,
xlab = "Probabilidad Predicha de Ampu",
ylab = "Amaputacion observada")
#modelo de Fine and gray
modelo_fg_san <- FGR(
formula = Hist(dias_out, evento_ampu) ~ san_elian + edad
, data = sad,
cause = 1 # Aclara que ampu es el evento de interés
)
modelo_fg_san##
## Right-censored response of a competing.risks model
##
## No.Observations: 406
##
## Pattern:
##
## Cause event right.censored
## 1 58 0
## 2 246 0
## unknown 0 102
##
##
## Fine-Gray model: analysis of cause 1
##
## Competing Risks Regression
##
## Call:
## FGR(formula = Hist(dias_out, evento_ampu) ~ san_elian + edad,
## data = sad, cause = 1)
##
## coef exp(coef) se(coef) z p-value
## san_elian 0.3442 1.41 0.0361 9.53 0.0000
## edad 0.0332 1.03 0.0108 3.08 0.0021
##
## exp(coef) exp(-coef) 2.5% 97.5%
## san_elian 1.41 0.709 1.31 1.51
## edad 1.03 0.967 1.01 1.06
##
## Num. cases = 406
## Pseudo Log-likelihood = -291
## Pseudo likelihood ratio test = 97.1 on 2 df,
##
## Convergence: TRUE# Evaluar performance predictiva con Score. AUC
score_fg_san <- Score(
object = list("FG model" = modelo_fg_san),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sad,
times = c(17:200),
metrics = "auc",
summary = "IPA",
cens.model = "km",
cause = 1
)
score_fg_san$AUC$score## model times AUC se lower upper
## <fctr> <int> <num> <num> <num> <num>
## 1: FG model 17 0.892 0.0221 0.849 0.935
## 2: FG model 18 0.896 0.0212 0.854 0.937
## 3: FG model 19 0.910 0.0184 0.874 0.946
## 4: FG model 20 0.910 0.0184 0.874 0.946
## 5: FG model 21 0.915 0.0179 0.880 0.950
## ---
## 180: FG model 196 0.848 0.0267 0.796 0.900
## 181: FG model 197 0.847 0.0268 0.794 0.899
## 182: FG model 198 0.847 0.0268 0.794 0.899
## 183: FG model 199 0.847 0.0268 0.794 0.899
## 184: FG model 200 0.846 0.0269 0.793 0.898# Calibracion global
s <- Score(list("FGR" = modelo_fg_san),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sad,
plots = "calibration",
summary = "IPA",
cause = 1)
plotCalibration(s, models = "FGR" )
tendencia a la subestimación del riesgo en los pacientes con mayor
probabilidad predicha
# Calibracion tiempo especifico
t <- Score(list("FGR" = modelo_fg_san),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sad,
plots = "calibration",
summary = "IPA",
times = c(1:180),
cause = 1)
plotCalibration(t, models = "FGR", times = 30)
plotCalibration(t, models = "FGR", times = 60 )
plotCalibration(t, models = "FGR", times = 90)
plotCalibration(t, models = "FGR", times = 130 )
plotCalibration(t, models = "FGR", times = 180 )
library(prodlim)
# Modelo Fine-Gray
mod_fg <- FGR(Hist(dias_out, evento_ampu) ~ san_elian + edad,
data = sad,
cause = 1)
# Tiempos de interés
times <- c(30, 60, 90, 120, 180)
#CALCULADORA!
# CIF según covariables
newdata0 <- data.frame(
san_elian =20,
edad = 72
)
predictRisk(mod_fg, newdata = newdata0, times = times)## [,1] [,2] [,3] [,4] [,5]
## [1,] 0.103 0.142 0.163 0.188 0.198times_seq <- 1:180
cif_base <- predictRisk(mod_fg, newdata = newdata0, times = times_seq)
plot(times_seq, cif_base[1,], type = "l",
xlab = "Días", ylab = "CIF basal",
main = "Baseline CIF")
#splines, evaluacion de posible no linealidad
library(prodlim)
library(splines)
library(splines)
spl_san <- ns(sad$san_elian, df = 3)
sad$san1 <- spl_san[,1]
sad$san2 <- spl_san[,2]
sad$san3 <- spl_san[,3]
library(riskRegression)
mod_spline <- FGR(
Hist(dias_out, evento_ampu) ~ san1 + san2 + san3 + edad,
data = sad,
cause = 1
)
mod_lin <- FGR(
Hist(dias_out, evento_ampu) ~ san_elian + edad,
data = sad,
cause = 1
)
sc <- Score(
object = list("Lineal" = mod_lin, "Spline" = mod_spline),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sad,
times = c(180),
metrics = c("auc", "brier"),
cause = 1
)
sc##
## Metric AUC:
##
## Results by model:
##
## model times AUC lower upper
## <fctr> <num> <char> <char> <char>
## 1: Lineal 180 85.0 79.9 90.2
## 2: Spline 180 85.0 79.8 90.2
##
## Results of model comparisons:
##
## times model reference delta.AUC lower upper p
## <num> <fctr> <fctr> <char> <char> <char> <num>
## 1: 180 Spline Lineal -0.0 -0.1 0.1 0.4##
## NOTE: Values are multiplied by 100 and given in %.## NOTE: The higher AUC the better.##
## Metric Brier:
##
## Results by model:
##
## model times Brier lower upper
## <fctr> <num> <char> <char> <char>
## 1: Null model 180 12.7 10.2 15.2
## 2: Lineal 180 9.4 7.3 11.5
## 3: Spline 180 9.4 7.3 11.5
##
## Results of model comparisons:
##
## times model reference delta.Brier lower upper p
## <num> <fctr> <fctr> <char> <char> <char> <num>
## 1: 180 Lineal Null model -3.3 -5.0 -1.6 0.000106
## 2: 180 Spline Null model -3.3 -5.0 -1.6 0.000108
## 3: 180 Spline Lineal 0.0 -0.0 0.0 0.932968##
## NOTE: Values are multiplied by 100 and given in %.## NOTE: The lower Brier the better.Se evaluó la posible no linealidad en la relación entre la escala de San Elián y la amputación mayor mediante modelos con splines restringidos. Sin embargo, no se observaron mejoras en la discriminación ni en el error de predicción respecto al modelo lineal, lo que sugiere que la relación puede ser adecuadamente modelada como lineal en esta cohorte. Entonces… ¿por qué tenías mala calibración arriba?
Más probable que sea por:
pocos pacientes en alto riesgo pocos eventos en ese rango variabilidad estadística no por forma funcional del predictor
#voy a hacer bootstrap para evaluar optimismo. medidas de discriminacion y calibracion corregidas
#library(riskRegression)
library(prodlim)
mod_fg <- FGR(
Hist(dias_out, evento_ampu) ~ san_elian + edad,
data = sad,
cause = 1
)
set.seed(123)
sc_boot <- Score(
object = list("Modelo" = mod_fg),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sad,
times = 180,
metrics = c("auc", "brier"),
summary = "risks",
cause = 1,
split.method = "bootcv", # ← bootstrap
B = 200, # ← número de muestras bootstrap (200–500 ideal)
conf.int = TRUE,
plots = "calibration"
)## Running crossvalidation algorithm## Fitting the models in 200 learning datasets, then predicting the risks in validation datasets## | | | 0% | | | 1% | | | 2% | |= | 3% | |= | 4% | |= | 5% | |= | 6% | |= | 7% | |== | 8% | |== | 9% | |== | 10% | |== | 11% | |== | 12% | |=== | 13% | |=== | 14% | |=== | 15% | |=== | 16% | |=== | 17% | |==== | 18% | |==== | 19% | |==== | 20% | |==== | 21% | |==== | 22% | |===== | 23% | |===== | 24% | |===== | 25% | |===== | 26% | |===== | 27% | |====== | 28% | |====== | 29% | |====== | 30% | |====== | 31% | |====== | 32% | |======= | 33% | |======= | 34% | |======= | 35% | |======= | 36% | |======= | 37% | |======== | 38% | |======== | 39% | |======== | 40% | |======== | 41% | |======== | 42% | |========= | 43% | |========= | 44% | |========= | 45% | |========= | 46% | |========= | 47% | |========== | 48% | |========== | 49% | |========== | 50% | |========== | 51% | |========== | 52% | |=========== | 53% | |=========== | 54% | |=========== | 55% | |=========== | 56% | |=========== | 57% | |============ | 58% | |============ | 59% | |============ | 60% | |============ | 61% | |============ | 62% | |============= | 63% | |============= | 64% | |============= | 65% | |============= | 66% | |============= | 67% | |============== | 68% | |============== | 69% | |============== | 70% | |============== | 71% | |============== | 72% | |=============== | 73% | |=============== | 74% | |=============== | 75% | |=============== | 76% | |=============== | 77% | |================ | 78% | |================ | 79% | |================ | 80% | |================ | 81% | |================ | 82% | |================= | 83% | |================= | 84% | |================= | 85% | |================= | 86% | |================= | 87% | |================== | 88% | |================== | 89% | |================== | 90% | |================== | 91% | |================== | 92% | |=================== | 93% | |=================== | 94% | |=================== | 95% | |=================== | 96% | |=================== | 97% | |====================| 98% | |====================| 99% | |====================| 100%
## | | | 0%## Calculating the performance metrics in 200 validation data sets## | | | 1% | | | 2% | |= | 3% | |= | 4% | |= | 5% | |= | 6% | |= | 7% | |== | 8% | |== | 9% | |== | 10% | |== | 11% | |== | 12% | |=== | 13% | |=== | 14% | |=== | 15% | |=== | 16% | |=== | 17% | |==== | 18% | |==== | 19% | |==== | 20% | |==== | 21% | |==== | 22% | |===== | 23% | |===== | 24% | |===== | 25% | |===== | 26% | |===== | 27% | |====== | 28% | |====== | 29% | |====== | 30% | |====== | 31% | |====== | 32% | |======= | 33% | |======= | 34% | |======= | 35% | |======= | 36% | |======= | 37% | |======== | 38% | |======== | 39% | |======== | 40% | |======== | 41% | |======== | 42% | |========= | 43% | |========= | 44% | |========= | 45% | |========= | 46% | |========= | 47% | |========== | 48% | |========== | 49% | |========== | 50% | |========== | 51% | |========== | 52% | |=========== | 53% | |=========== | 54% | |=========== | 55% | |=========== | 56% | |=========== | 57% | |============ | 58% | |============ | 59% | |============ | 60% | |============ | 61% | |============ | 62% | |============= | 63% | |============= | 64% | |============= | 65% | |============= | 66% | |============= | 67% | |============== | 68% | |============== | 69% | |============== | 70% | |============== | 71% | |============== | 72% | |=============== | 73% | |=============== | 74% | |=============== | 75% | |=============== | 76% | |=============== | 77% | |================ | 78% | |================ | 79% | |================ | 80% | |================ | 81% | |================ | 82% | |================= | 83% | |================= | 84% | |================= | 85% | |================= | 86% | |================= | 87% | |================== | 88% | |================== | 89% | |================== | 90% | |================== | 91% | |================== | 92% | |=================== | 93% | |=================== | 94% | |=================== | 95% | |=================== | 96% | |=================== | 97% | |====================| 98% | |====================| 99% | |====================| 100%sc_boot$AUC##
## Results by model:
##
## model times AUC lower upper
## <fctr> <num> <char> <char> <char>
## 1: Modelo 180 84.753 78.412 90.357##
## NOTE: Values are multiplied by 100 and given in %.## NOTE: The higher AUC the better.sc_boot$Brier##
## Results by model:
##
## model times Brier lower upper
## <fctr> <num> <char> <char> <char>
## 1: Null model 180 12.895 9.995 15.806
## 2: Modelo 180 9.718 7.309 12.658
##
## Results of model comparisons:
##
## model times reference delta.Brier lower upper
## <fctr> <num> <fctr> <char> <char> <char>
## 1: Modelo 180 Null model -3.177 -5.059 -1.126##
## NOTE: Values are multiplied by 100 and given in %.## NOTE: The lower Brier the better.plotCalibration(sc_boot)
La validación interna se realizó mediante bootstrap con 200
repeticiones, estimando el optimismo del modelo y obteniendo medidas de
discriminación y error de predicción corregidas.
#categorización de san elian
library(dplyr)
sad <- sad %>%
mutate(
san_elian_cat = case_when(
san_elian <= 9 ~ "leve",
san_elian >= 10 & san_elian <= 20 ~ "moderada",
san_elian >= 21 ~ "grave",
TRUE ~ NA_character_
)
)#CIF
cifivhx <- cuminc(Surv(dias_out, factor(evento_ampu)) ~ san_elian_cat, data=sad)
cifivhx## ## ── cuminc() ────────────────────────────────────────────────────────────────────## • Failure type "1"## strata time n.risk estimate std.error 95% CI
## grave 100 50 0.338 0.044 0.253, 0.424
## grave 200 18 0.358 0.045 0.271, 0.445
## grave 300 3 0.371 0.046 0.282, 0.460
## leve 100 1 0.000 0.000 NA, NA
## leve 200 0 0.000 0.000 NA, NA
## leve 300 0 0.000 0.000 NA, NA
## moderada 100 110 0.042 0.012 0.022, 0.071
## moderada 200 32 0.060 0.015 0.035, 0.095
## moderada 300 9 0.060 0.015 0.035, 0.095## • Failure type "2"## strata time n.risk estimate std.error 95% CI
## grave 100 50 0.167 0.035 0.105, 0.242
## grave 200 18 0.416 0.050 0.318, 0.510
## grave 300 3 0.578 0.052 0.470, 0.671
## leve 100 1 0.833 0.202 0.086, 0.987
## leve 200 0 0.833 0.202 0.086, 0.987
## leve 300 0 0.833 0.202 0.086, 0.987
## moderada 100 110 0.488 0.032 0.425, 0.548
## moderada 200 32 0.726 0.031 0.659, 0.781
## moderada 300 9 0.844 0.031 0.772, 0.895## • Tests## outcome statistic df p.value
## 1 63.9 2.00 <0.001
## 2 41.7 2.00 <0.001library(tidycmprsk)
library(ggsurvfit)
ggcuminc(
cifivhx,
outcome = c(1, 2),
linetype_aes = TRUE
) +
scale_color_manual(
name = "san_elian_cat",
values = c("leve" = "#56B4E9", "moderada" = "#E69F00","grave"="#a25")
) +
scale_linetype_manual(
name = "Tipo de evento",
values = c("solid", "dashed"),
labels = c(
"Evento de interes (amputacion)",
"Evento competitivo (muerte / cicatrizacion)"
)
) +
labs(
x = "Tiempo (dias)",
y = "Incidencia acumulada"
) +
theme_minimal() +
theme(
legend.position = "bottom"
)## Warning: Removed 170 rows containing missing values or values outside the scale range
## (`geom_step()`).
#Gray test
cuminc_obja <- cmprsk::cuminc(
ftime = sad$dias_out,
fstatus = sad$evento_ampu,
group = sad$san_elian_cat
)
cuminc_obja## Tests:
## stat pv df
## 1 63.9 0.000000000000013 2
## 2 41.7 0.000000000869993 2
## Estimates and Variances:
## $est
## 100 200 300
## grave 1 0.3376 0.3576 0.3710
## leve 1 0.0000 NA NA
## moderada 1 0.0419 0.0605 0.0605
## grave 2 0.1671 0.4156 0.5775
## leve 2 0.8333 NA NA
## moderada 2 0.4880 0.7255 0.8442
##
## $var
## 100 200 300
## grave 1 0.001925 0.002007 0.002108
## leve 1 0.000000 NA NA
## moderada 1 0.000153 0.000233 0.000233
## grave 2 0.001251 0.002454 0.002665
## leve 2 0.040656 NA NA
## moderada 2 0.000996 0.000966 0.000957#Gray test e incidencias acumuladas a 180 días
sad$dias_180 <- pmin(sad$dias_out, 180)
sad$evento_180 <- ifelse(sad$dias_out <= 180, sad$evento_ampu, 0)
cuminc_obja_180 <- cmprsk::cuminc(
ftime = sad$dias_180,
fstatus = sad$evento_180,
group = sad$san_elian_cat
)
cuminc_obja_180## Tests:
## stat pv df
## 1 61.9 0.0000000000000359 2
## 2 37.0 0.0000000092186401 2
## Estimates and Variances:
## $est
## 50 100 150
## grave 1 0.2826 0.3376 0.3576
## leve 1 0.0000 0.0000 NA
## moderada 1 0.0340 0.0419 0.0605
## grave 2 0.0916 0.1671 0.3447
## leve 2 0.6667 0.8333 NA
## moderada 2 0.2156 0.4880 0.6276
##
## $var
## 50 100 150
## grave 1 0.001702 0.001925 0.002007
## leve 1 0.000000 0.000000 NA
## moderada 1 0.000125 0.000153 0.000233
## grave 2 0.000700 0.001251 0.002209
## leve 2 0.051512 0.040656 NA
## moderada 2 0.000643 0.000996 0.001008mod_fg_cat <- FGR(
Hist(dias_out, evento_ampu) ~ san_elian_cat + edad,
data = sad,
cause = 1
)
mod_fg_cat##
## Right-censored response of a competing.risks model
##
## No.Observations: 406
##
## Pattern:
##
## Cause event right.censored
## 1 58 0
## 2 246 0
## unknown 0 102
##
##
## Fine-Gray model: analysis of cause 1
##
## Competing Risks Regression
##
## Call:
## FGR(formula = Hist(dias_out, evento_ampu) ~ san_elian_cat + edad,
## data = sad, cause = 1)
##
## coef exp(coef) se(coef) z p-value
## san_elian_catleve -11.5527 0.00000961 0.4618 -25.01 0.000000000000
## san_elian_catmoderada -2.0027 0.13496706 0.2968 -6.75 0.000000000015
## edad 0.0387 1.03941315 0.0114 3.38 0.000720000000
##
## exp(coef) exp(-coef) 2.5% 97.5%
## san_elian_catleve 0.00000961 104062.795 0.00000389 0.0000238
## san_elian_catmoderada 0.13496706 7.409 0.07543126 0.2414928
## edad 1.03941315 0.962 1.01638485 1.0629632
##
## Num. cases = 406
## Pseudo Log-likelihood = -305
## Pseudo likelihood ratio test = 69.1 on 3 df,
##
## Convergence: TRUE#curvas roc del modelo
pred_180 <- predictRisk(mod_fg, newdata = sad, times = 180)
pred_90 <- predictRisk(mod_fg, newdata = sad, times = 90)
pred_120 <- predictRisk(mod_fg, newdata = sad, times = 120)
library(timeROC)## Warning: package 'timeROC' was built under R version 4.4.3roc_model_180 <- timeROC(
T = sad$dias_out,
delta = sad$evento_ampu,
marker = as.numeric(pred_180),
cause = 1,
times = 180,
iid = TRUE
)
plot(roc_model_180, time = 180)
roc_model_90 <- timeROC(
T = sad$dias_out,
delta = sad$evento_ampu,
marker = as.numeric(pred_90),
cause = 1,
times = 90,
iid = TRUE
)
plot(roc_model_90, time = 90)
roc_model_120 <- timeROC(
T = sad$dias_out,
delta = sad$evento_ampu,
marker = as.numeric(pred_120),
cause = 1,
times = 120,
iid = TRUE
)
plot(roc_model_120, time = 120)
MODELO DE CAUSA ESPECÍFICA
#modelo de causa específica
#Opcion con objeto creado
#si evento es 1, es porque se amputo
#si evento es 2, es porque se murió o cicatrizó
table(sad$evento_ampu)##
## 0 1 2
## 102 58 246tiempo_ampu <- Surv(sad$dias_out, sad$evento_ampu==1)
tiempo_compet <- Surv(sad$dias_out, sad$evento_ampu==2)
cox_ampu <- coxph(tiempo_ampu ~
san_elian + edad , data = sad)
summary(cox_ampu)## Call:
## coxph(formula = tiempo_ampu ~ san_elian + edad, data = sad)
##
## n= 406, number of events= 58
##
## coef exp(coef) se(coef) z Pr(>|z|)
## san_elian 0.3424 1.4083 0.0416 8.24 <0.0000000000000002 ***
## edad 0.0361 1.0368 0.0137 2.64 0.0083 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## san_elian 1.41 0.710 1.30 1.53
## edad 1.04 0.965 1.01 1.06
##
## Concordance= 0.845 (se = 0.022 )
## Likelihood ratio test= 93.7 on 2 df, p=<0.0000000000000002
## Wald test = 82.1 on 2 df, p=<0.0000000000000002
## Score (logrank) test = 85.9 on 2 df, p=<0.0000000000000002cox_compet <- coxph(tiempo_compet ~
san_elian + edad
, data = sad)
summary(cox_compet)## Call:
## coxph(formula = tiempo_compet ~ san_elian + edad, data = sad)
##
## n= 406, number of events= 246
##
## coef exp(coef) se(coef) z Pr(>|z|)
## san_elian -0.06442 0.93761 0.01742 -3.70 0.00022 ***
## edad -0.02199 0.97825 0.00656 -3.35 0.00081 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## san_elian 0.938 1.07 0.906 0.970
## edad 0.978 1.02 0.966 0.991
##
## Concordance= 0.606 (se = 0.022 )
## Likelihood ratio test= 28.2 on 2 df, p=0.0000008
## Wald test = 28 on 2 df, p=0.0000008
## Score (logrank) test = 28.2 on 2 df, p=0.0000007#modelo de causa específica
modelo_CSC <- CSC(
formula = Hist(dias_out, evento_ampu) ~ san_elian +edad, data = sad,
cause = 1 # Aclara que CHD es el evento de interés
)
# Ver resultados. Ojo! No funciona si hacés summary(modelo)
# Vamos a tener el modelo de interés y debajo el modelo competitivo
modelo_CSC## CSC(formula = Hist(dias_out, evento_ampu) ~ san_elian + edad,
## data = sad, cause = 1)
##
## Right-censored response of a competing.risks model
##
## No.Observations: 406
##
## Pattern:
##
## Cause event right.censored
## 1 58 0
## 2 246 0
## unknown 0 102
##
##
## ----------> Cause: 1
##
## Call:
## coxph(formula = survival::Surv(time, status) ~ san_elian + edad,
## x = TRUE, y = TRUE)
##
## n= 406, number of events= 58
##
## coef exp(coef) se(coef) z Pr(>|z|)
## san_elian 0.3424 1.4083 0.0416 8.24 <0.0000000000000002 ***
## edad 0.0361 1.0368 0.0137 2.64 0.0083 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## san_elian 1.41 0.710 1.30 1.53
## edad 1.04 0.965 1.01 1.06
##
## Concordance= 0.845 (se = 0.022 )
## Likelihood ratio test= 93.7 on 2 df, p=<0.0000000000000002
## Wald test = 82.1 on 2 df, p=<0.0000000000000002
## Score (logrank) test = 85.9 on 2 df, p=<0.0000000000000002
##
##
##
## ----------> Cause: 2
##
## Call:
## coxph(formula = survival::Surv(time, status) ~ san_elian + edad,
## x = TRUE, y = TRUE)
##
## n= 406, number of events= 246
##
## coef exp(coef) se(coef) z Pr(>|z|)
## san_elian -0.06442 0.93761 0.01742 -3.70 0.00022 ***
## edad -0.02199 0.97825 0.00656 -3.35 0.00081 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## san_elian 0.938 1.07 0.906 0.970
## edad 0.978 1.02 0.966 0.991
##
## Concordance= 0.606 (se = 0.022 )
## Likelihood ratio test= 28.2 on 2 df, p=0.0000008
## Wald test = 28 on 2 df, p=0.0000008
## Score (logrank) test = 28.2 on 2 df, p=0.0000007#score comparado con ipa
score_comp <-Score(
list("FG model" = modelo_fg_san, "CSC model" = modelo_CSC),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sad,
times = c(0:200),
cause = 1,
metrics ="auc",
summary = "IPA"
)#Extraigo las IPAs por tiempo y modelo
ipa <- as.data.frame(score_comp$Brier$score)[, c("model", "times", "IPA")]
ipa## model times IPA
## 1 Null model 0 NaN
## 2 Null model 1 0.00000
## 3 Null model 2 0.00000
## 4 Null model 3 0.00000
## 5 Null model 4 0.00000
## 6 Null model 5 0.00000
## 7 Null model 6 0.00000
## 8 Null model 7 0.00000
## 9 Null model 8 0.00000
## 10 Null model 9 0.00000
## 11 Null model 10 0.00000
## 12 Null model 11 0.00000
## 13 Null model 12 0.00000
## 14 Null model 13 0.00000
## 15 Null model 14 0.00000
## 16 Null model 15 0.00000
## 17 Null model 16 0.00000
## 18 Null model 17 0.00000
## 19 Null model 18 0.00000
## 20 Null model 19 0.00000
## 21 Null model 20 0.00000
## 22 Null model 21 0.00000
## 23 Null model 22 0.00000
## 24 Null model 23 0.00000
## 25 Null model 24 0.00000
## 26 Null model 25 0.00000
## 27 Null model 26 0.00000
## 28 Null model 27 0.00000
## 29 Null model 28 0.00000
## 30 Null model 29 0.00000
## 31 Null model 30 0.00000
## 32 Null model 31 0.00000
## 33 Null model 32 0.00000
## 34 Null model 33 0.00000
## 35 Null model 34 0.00000
## 36 Null model 35 0.00000
## 37 Null model 36 0.00000
## 38 Null model 37 0.00000
## 39 Null model 38 0.00000
## 40 Null model 39 0.00000
## 41 Null model 40 0.00000
## 42 Null model 41 0.00000
## 43 Null model 42 0.00000
## 44 Null model 43 0.00000
## 45 Null model 44 0.00000
## 46 Null model 45 0.00000
## 47 Null model 46 0.00000
## 48 Null model 47 0.00000
## 49 Null model 48 0.00000
## 50 Null model 49 0.00000
## 51 Null model 50 0.00000
## 52 Null model 51 0.00000
## 53 Null model 52 0.00000
## 54 Null model 53 0.00000
## 55 Null model 54 0.00000
## 56 Null model 55 0.00000
## 57 Null model 56 0.00000
## 58 Null model 57 0.00000
## 59 Null model 58 0.00000
## 60 Null model 59 0.00000
## 61 Null model 60 0.00000
## 62 Null model 61 0.00000
## 63 Null model 62 0.00000
## 64 Null model 63 0.00000
## 65 Null model 64 0.00000
## 66 Null model 65 0.00000
## 67 Null model 66 0.00000
## 68 Null model 67 0.00000
## 69 Null model 68 0.00000
## 70 Null model 69 0.00000
## 71 Null model 70 0.00000
## 72 Null model 71 0.00000
## 73 Null model 72 0.00000
## 74 Null model 73 0.00000
## 75 Null model 74 0.00000
## 76 Null model 75 0.00000
## 77 Null model 76 0.00000
## 78 Null model 77 0.00000
## 79 Null model 78 0.00000
## 80 Null model 79 0.00000
## 81 Null model 80 0.00000
## 82 Null model 81 0.00000
## 83 Null model 82 0.00000
## 84 Null model 83 0.00000
## 85 Null model 84 0.00000
## 86 Null model 85 0.00000
## 87 Null model 86 0.00000
## 88 Null model 87 0.00000
## 89 Null model 88 0.00000
## 90 Null model 89 0.00000
## 91 Null model 90 0.00000
## 92 Null model 91 0.00000
## 93 Null model 92 0.00000
## 94 Null model 93 0.00000
## 95 Null model 94 0.00000
## 96 Null model 95 0.00000
## 97 Null model 96 0.00000
## 98 Null model 97 0.00000
## 99 Null model 98 0.00000
## 100 Null model 99 0.00000
## 101 Null model 100 0.00000
## 102 Null model 101 0.00000
## 103 Null model 102 0.00000
## 104 Null model 103 0.00000
## 105 Null model 104 0.00000
## 106 Null model 105 0.00000
## 107 Null model 106 0.00000
## 108 Null model 107 0.00000
## 109 Null model 108 0.00000
## 110 Null model 109 0.00000
## 111 Null model 110 0.00000
## 112 Null model 111 0.00000
## 113 Null model 112 0.00000
## 114 Null model 113 0.00000
## 115 Null model 114 0.00000
## 116 Null model 115 0.00000
## 117 Null model 116 0.00000
## 118 Null model 117 0.00000
## 119 Null model 118 0.00000
## 120 Null model 119 0.00000
## 121 Null model 120 0.00000
## 122 Null model 121 0.00000
## 123 Null model 122 0.00000
## 124 Null model 123 0.00000
## 125 Null model 124 0.00000
## 126 Null model 125 0.00000
## 127 Null model 126 0.00000
## 128 Null model 127 0.00000
## 129 Null model 128 0.00000
## 130 Null model 129 0.00000
## 131 Null model 130 0.00000
## 132 Null model 131 0.00000
## 133 Null model 132 0.00000
## 134 Null model 133 0.00000
## 135 Null model 134 0.00000
## 136 Null model 135 0.00000
## 137 Null model 136 0.00000
## 138 Null model 137 0.00000
## 139 Null model 138 0.00000
## 140 Null model 139 0.00000
## 141 Null model 140 0.00000
## 142 Null model 141 0.00000
## 143 Null model 142 0.00000
## 144 Null model 143 0.00000
## 145 Null model 144 0.00000
## 146 Null model 145 0.00000
## 147 Null model 146 0.00000
## 148 Null model 147 0.00000
## 149 Null model 148 0.00000
## 150 Null model 149 0.00000
## 151 Null model 150 0.00000
## 152 Null model 151 0.00000
## 153 Null model 152 0.00000
## 154 Null model 153 0.00000
## 155 Null model 154 0.00000
## 156 Null model 155 0.00000
## 157 Null model 156 0.00000
## 158 Null model 157 0.00000
## 159 Null model 158 0.00000
## 160 Null model 159 0.00000
## 161 Null model 160 0.00000
## 162 Null model 161 0.00000
## 163 Null model 162 0.00000
## 164 Null model 163 0.00000
## 165 Null model 164 0.00000
## 166 Null model 165 0.00000
## 167 Null model 166 0.00000
## 168 Null model 167 0.00000
## 169 Null model 168 0.00000
## 170 Null model 169 0.00000
## 171 Null model 170 0.00000
## 172 Null model 171 0.00000
## 173 Null model 172 0.00000
## 174 Null model 173 0.00000
## 175 Null model 174 0.00000
## 176 Null model 175 0.00000
## 177 Null model 176 0.00000
## 178 Null model 177 0.00000
## 179 Null model 178 0.00000
## 180 Null model 179 0.00000
## 181 Null model 180 0.00000
## 182 Null model 181 0.00000
## 183 Null model 182 0.00000
## 184 Null model 183 0.00000
## 185 Null model 184 0.00000
## 186 Null model 185 0.00000
## 187 Null model 186 0.00000
## 188 Null model 187 0.00000
## 189 Null model 188 0.00000
## 190 Null model 189 0.00000
## 191 Null model 190 0.00000
## 192 Null model 191 0.00000
## 193 Null model 192 0.00000
## 194 Null model 193 0.00000
## 195 Null model 194 0.00000
## 196 Null model 195 0.00000
## 197 Null model 196 0.00000
## 198 Null model 197 0.00000
## 199 Null model 198 0.00000
## 200 Null model 199 0.00000
## 201 Null model 200 0.00000
## 202 FG model 0 NaN
## 203 FG model 1 0.02497
## 204 FG model 2 0.02497
## 205 FG model 3 0.02497
## 206 FG model 4 0.02497
## 207 FG model 5 0.02497
## 208 FG model 6 0.02491
## 209 FG model 7 0.02348
## 210 FG model 8 0.02339
## 211 FG model 9 0.02568
## 212 FG model 10 0.01370
## 213 FG model 11 0.02123
## 214 FG model 12 0.05061
## 215 FG model 13 0.09874
## 216 FG model 14 0.08498
## 217 FG model 15 0.08275
## 218 FG model 16 0.08238
## 219 FG model 17 0.08531
## 220 FG model 18 0.15003
## 221 FG model 19 0.22796
## 222 FG model 20 0.22796
## 223 FG model 21 0.24986
## 224 FG model 22 0.25465
## 225 FG model 23 0.26802
## 226 FG model 24 0.25917
## 227 FG model 25 0.26927
## 228 FG model 26 0.26927
## 229 FG model 27 0.30532
## 230 FG model 28 0.30470
## 231 FG model 29 0.31251
## 232 FG model 30 0.31251
## 233 FG model 31 0.31251
## 234 FG model 32 0.31251
## 235 FG model 33 0.31251
## 236 FG model 34 0.31219
## 237 FG model 35 0.31219
## 238 FG model 36 0.31219
## 239 FG model 37 0.31211
## 240 FG model 38 0.31211
## 241 FG model 39 0.31211
## 242 FG model 40 0.30443
## 243 FG model 41 0.29462
## 244 FG model 42 0.28899
## 245 FG model 43 0.28899
## 246 FG model 44 0.28899
## 247 FG model 45 0.28899
## 248 FG model 46 0.28505
## 249 FG model 47 0.28505
## 250 FG model 48 0.26565
## 251 FG model 49 0.26520
## 252 FG model 50 0.26475
## 253 FG model 51 0.26475
## 254 FG model 52 0.26475
## 255 FG model 53 0.26428
## 256 FG model 54 0.26428
## 257 FG model 55 0.26428
## 258 FG model 56 0.26428
## 259 FG model 57 0.27185
## 260 FG model 58 0.27185
## 261 FG model 59 0.28075
## 262 FG model 60 0.28075
## 263 FG model 61 0.28075
## 264 FG model 62 0.28075
## 265 FG model 63 0.28258
## 266 FG model 64 0.27434
## 267 FG model 65 0.27434
## 268 FG model 66 0.28180
## 269 FG model 67 0.28127
## 270 FG model 68 0.28216
## 271 FG model 69 0.28216
## 272 FG model 70 0.27655
## 273 FG model 71 0.27655
## 274 FG model 72 0.27594
## 275 FG model 73 0.28217
## 276 FG model 74 0.28217
## 277 FG model 75 0.28217
## 278 FG model 76 0.28217
## 279 FG model 77 0.28217
## 280 FG model 78 0.28160
## 281 FG model 79 0.28160
## 282 FG model 80 0.28160
## 283 FG model 81 0.28160
## 284 FG model 82 0.28073
## 285 FG model 83 0.28689
## 286 FG model 84 0.28689
## 287 FG model 85 0.28689
## 288 FG model 86 0.28689
## 289 FG model 87 0.28689
## 290 FG model 88 0.28689
## 291 FG model 89 0.28606
## 292 FG model 90 0.28883
## 293 FG model 91 0.28814
## 294 FG model 92 0.28749
## 295 FG model 93 0.28749
## 296 FG model 94 0.28749
## 297 FG model 95 0.28681
## 298 FG model 96 0.28681
## 299 FG model 97 0.28681
## 300 FG model 98 0.29269
## 301 FG model 99 0.29269
## 302 FG model 100 0.29269
## 303 FG model 101 0.29198
## 304 FG model 102 0.29178
## 305 FG model 103 0.29104
## 306 FG model 104 0.29034
## 307 FG model 105 0.29034
## 308 FG model 106 0.28979
## 309 FG model 107 0.28166
## 310 FG model 108 0.28166
## 311 FG model 109 0.28166
## 312 FG model 110 0.28166
## 313 FG model 111 0.28166
## 314 FG model 112 0.28080
## 315 FG model 113 0.28534
## 316 FG model 114 0.29374
## 317 FG model 115 0.29290
## 318 FG model 116 0.29290
## 319 FG model 117 0.29119
## 320 FG model 118 0.29119
## 321 FG model 119 0.28863
## 322 FG model 120 0.28845
## 323 FG model 121 0.27765
## 324 FG model 122 0.27765
## 325 FG model 123 0.27765
## 326 FG model 124 0.27765
## 327 FG model 125 0.27765
## 328 FG model 126 0.27765
## 329 FG model 127 0.27765
## 330 FG model 128 0.27661
## 331 FG model 129 0.27661
## 332 FG model 130 0.27661
## 333 FG model 131 0.27661
## 334 FG model 132 0.27661
## 335 FG model 133 0.27661
## 336 FG model 134 0.27661
## 337 FG model 135 0.27557
## 338 FG model 136 0.26828
## 339 FG model 137 0.26828
## 340 FG model 138 0.26828
## 341 FG model 139 0.26828
## 342 FG model 140 0.26603
## 343 FG model 141 0.26488
## 344 FG model 142 0.26488
## 345 FG model 143 0.26488
## 346 FG model 144 0.26399
## 347 FG model 145 0.26324
## 348 FG model 146 0.26324
## 349 FG model 147 0.26182
## 350 FG model 148 0.26182
## 351 FG model 149 0.26182
## 352 FG model 150 0.26182
## 353 FG model 151 0.26259
## 354 FG model 152 0.26259
## 355 FG model 153 0.26259
## 356 FG model 154 0.26259
## 357 FG model 155 0.26259
## 358 FG model 156 0.26259
## 359 FG model 157 0.26259
## 360 FG model 158 0.26259
## 361 FG model 159 0.26259
## 362 FG model 160 0.26259
## 363 FG model 161 0.26259
## 364 FG model 162 0.26259
## 365 FG model 163 0.26259
## 366 FG model 164 0.26259
## 367 FG model 165 0.26259
## 368 FG model 166 0.26259
## 369 FG model 167 0.26407
## 370 FG model 168 0.26407
## 371 FG model 169 0.26407
## 372 FG model 170 0.26170
## 373 FG model 171 0.26170
## 374 FG model 172 0.26170
## 375 FG model 173 0.26170
## 376 FG model 174 0.26170
## 377 FG model 175 0.26098
## 378 FG model 176 0.26098
## 379 FG model 177 0.26098
## 380 FG model 178 0.26098
## 381 FG model 179 0.26098
## 382 FG model 180 0.26098
## 383 FG model 181 0.25947
## 384 FG model 182 0.25947
## 385 FG model 183 0.25947
## 386 FG model 184 0.25947
## 387 FG model 185 0.25947
## 388 FG model 186 0.25947
## 389 FG model 187 0.25947
## 390 FG model 188 0.25947
## 391 FG model 189 0.25791
## 392 FG model 190 0.25618
## 393 FG model 191 0.25618
## 394 FG model 192 0.25516
## 395 FG model 193 0.25516
## 396 FG model 194 0.25516
## 397 FG model 195 0.25516
## 398 FG model 196 0.25401
## 399 FG model 197 0.25222
## 400 FG model 198 0.25082
## 401 FG model 199 0.25082
## 402 FG model 200 0.24877
## 403 CSC model 0 NaN
## 404 CSC model 1 0.02281
## 405 CSC model 2 0.02281
## 406 CSC model 3 0.02281
## 407 CSC model 4 0.02281
## 408 CSC model 5 0.02281
## 409 CSC model 6 0.02275
## 410 CSC model 7 0.02033
## 411 CSC model 8 0.02023
## 412 CSC model 9 0.02206
## 413 CSC model 10 0.00827
## 414 CSC model 11 0.01281
## 415 CSC model 12 0.04269
## 416 CSC model 13 0.09666
## 417 CSC model 14 0.07839
## 418 CSC model 15 0.07569
## 419 CSC model 16 0.07528
## 420 CSC model 17 0.07866
## 421 CSC model 18 0.14860
## 422 CSC model 19 0.23071
## 423 CSC model 20 0.23071
## 424 CSC model 21 0.25294
## 425 CSC model 22 0.25745
## 426 CSC model 23 0.27062
## 427 CSC model 24 0.26089
## 428 CSC model 25 0.27170
## 429 CSC model 26 0.27170
## 430 CSC model 27 0.30914
## 431 CSC model 28 0.30845
## 432 CSC model 29 0.31620
## 433 CSC model 30 0.31620
## 434 CSC model 31 0.31620
## 435 CSC model 32 0.31620
## 436 CSC model 33 0.31620
## 437 CSC model 34 0.31585
## 438 CSC model 35 0.31585
## 439 CSC model 36 0.31585
## 440 CSC model 37 0.31521
## 441 CSC model 38 0.31521
## 442 CSC model 39 0.31521
## 443 CSC model 40 0.30681
## 444 CSC model 41 0.29629
## 445 CSC model 42 0.29003
## 446 CSC model 43 0.29003
## 447 CSC model 44 0.29003
## 448 CSC model 45 0.29003
## 449 CSC model 46 0.28522
## 450 CSC model 47 0.28522
## 451 CSC model 48 0.26458
## 452 CSC model 49 0.26409
## 453 CSC model 50 0.26359
## 454 CSC model 51 0.26359
## 455 CSC model 52 0.26359
## 456 CSC model 53 0.26309
## 457 CSC model 54 0.26309
## 458 CSC model 55 0.26309
## 459 CSC model 56 0.26309
## 460 CSC model 57 0.27191
## 461 CSC model 58 0.27191
## 462 CSC model 59 0.27961
## 463 CSC model 60 0.27961
## 464 CSC model 61 0.27961
## 465 CSC model 62 0.27961
## 466 CSC model 63 0.28166
## 467 CSC model 64 0.27255
## 468 CSC model 65 0.27255
## 469 CSC model 66 0.28014
## 470 CSC model 67 0.27951
## 471 CSC model 68 0.28053
## 472 CSC model 69 0.28053
## 473 CSC model 70 0.27445
## 474 CSC model 71 0.27445
## 475 CSC model 72 0.27380
## 476 CSC model 73 0.28091
## 477 CSC model 74 0.28091
## 478 CSC model 75 0.28091
## 479 CSC model 76 0.28091
## 480 CSC model 77 0.28091
## 481 CSC model 78 0.28029
## 482 CSC model 79 0.28029
## 483 CSC model 80 0.28029
## 484 CSC model 81 0.28029
## 485 CSC model 82 0.27937
## 486 CSC model 83 0.28594
## 487 CSC model 84 0.28594
## 488 CSC model 85 0.28594
## 489 CSC model 86 0.28594
## 490 CSC model 87 0.28594
## 491 CSC model 88 0.28594
## 492 CSC model 89 0.28518
## 493 CSC model 90 0.28815
## 494 CSC model 91 0.28741
## 495 CSC model 92 0.28670
## 496 CSC model 93 0.28670
## 497 CSC model 94 0.28670
## 498 CSC model 95 0.28597
## 499 CSC model 96 0.28597
## 500 CSC model 97 0.28597
## 501 CSC model 98 0.29150
## 502 CSC model 99 0.29150
## 503 CSC model 100 0.29150
## 504 CSC model 101 0.29073
## 505 CSC model 102 0.29059
## 506 CSC model 103 0.28979
## 507 CSC model 104 0.28925
## 508 CSC model 105 0.28925
## 509 CSC model 106 0.28864
## 510 CSC model 107 0.27982
## 511 CSC model 108 0.27982
## 512 CSC model 109 0.27982
## 513 CSC model 110 0.27982
## 514 CSC model 111 0.27982
## 515 CSC model 112 0.27890
## 516 CSC model 113 0.28345
## 517 CSC model 114 0.29263
## 518 CSC model 115 0.29174
## 519 CSC model 116 0.29174
## 520 CSC model 117 0.29058
## 521 CSC model 118 0.29058
## 522 CSC model 119 0.28786
## 523 CSC model 120 0.28785
## 524 CSC model 121 0.27642
## 525 CSC model 122 0.27642
## 526 CSC model 123 0.27642
## 527 CSC model 124 0.27642
## 528 CSC model 125 0.27642
## 529 CSC model 126 0.27642
## 530 CSC model 127 0.27642
## 531 CSC model 128 0.27532
## 532 CSC model 129 0.27532
## 533 CSC model 130 0.27532
## 534 CSC model 131 0.27532
## 535 CSC model 132 0.27532
## 536 CSC model 133 0.27532
## 537 CSC model 134 0.27532
## 538 CSC model 135 0.27422
## 539 CSC model 136 0.26688
## 540 CSC model 137 0.26688
## 541 CSC model 138 0.26688
## 542 CSC model 139 0.26688
## 543 CSC model 140 0.26451
## 544 CSC model 141 0.26329
## 545 CSC model 142 0.26329
## 546 CSC model 143 0.26329
## 547 CSC model 144 0.26231
## 548 CSC model 145 0.26157
## 549 CSC model 146 0.26157
## 550 CSC model 147 0.25997
## 551 CSC model 148 0.25997
## 552 CSC model 149 0.25997
## 553 CSC model 150 0.25997
## 554 CSC model 151 0.26059
## 555 CSC model 152 0.26059
## 556 CSC model 153 0.26059
## 557 CSC model 154 0.26059
## 558 CSC model 155 0.26059
## 559 CSC model 156 0.26059
## 560 CSC model 157 0.26059
## 561 CSC model 158 0.26059
## 562 CSC model 159 0.26059
## 563 CSC model 160 0.26059
## 564 CSC model 161 0.26059
## 565 CSC model 162 0.26059
## 566 CSC model 163 0.26059
## 567 CSC model 164 0.26059
## 568 CSC model 165 0.26059
## 569 CSC model 166 0.26059
## 570 CSC model 167 0.26206
## 571 CSC model 168 0.26206
## 572 CSC model 169 0.26206
## 573 CSC model 170 0.25955
## 574 CSC model 171 0.25955
## 575 CSC model 172 0.25955
## 576 CSC model 173 0.25955
## 577 CSC model 174 0.25955
## 578 CSC model 175 0.25880
## 579 CSC model 176 0.25880
## 580 CSC model 177 0.25880
## 581 CSC model 178 0.25880
## 582 CSC model 179 0.25880
## 583 CSC model 180 0.25880
## 584 CSC model 181 0.25720
## 585 CSC model 182 0.25720
## 586 CSC model 183 0.25720
## 587 CSC model 184 0.25720
## 588 CSC model 185 0.25720
## 589 CSC model 186 0.25720
## 590 CSC model 187 0.25720
## 591 CSC model 188 0.25720
## 592 CSC model 189 0.25554
## 593 CSC model 190 0.25372
## 594 CSC model 191 0.25372
## 595 CSC model 192 0.25256
## 596 CSC model 193 0.25256
## 597 CSC model 194 0.25256
## 598 CSC model 195 0.25256
## 599 CSC model 196 0.25140
## 600 CSC model 197 0.24951
## 601 CSC model 198 0.24805
## 602 CSC model 199 0.24805
## 603 CSC model 200 0.24587# Gráfico
library(ggplot2)
ggplot(ipa, aes(x = times, y = IPA, color = model)) +
geom_line(linewidth = 1.2) +
geom_point(size = 3) +
geom_text(aes(label = percent(IPA, accuracy = 0.1)),
vjust = -0.8, size = 3, show.legend = FALSE) +
scale_y_continuous(labels = percent_format(accuracy = 0.1),
limits = c(0, NA)) +
scale_color_manual(
values = c(
"Null model" = "gray40",
"FG model" = "#2ECC71",
"CSC model" = "#3498DB"
),
drop = FALSE # <- muy importante para que no elimine niveles
) +
labs(
title = "(IPA) en el tiempo",
x = "Dias de seguimiento", y = "IPA", color = "Modelo"
) +
theme_minimal() +
theme(legend.position = "bottom", plot.title.position = "plot")+ scale_x_continuous(
limits = c(30, 200),
breaks = seq(30, 200, by = 10)
)## Warning: Removed 90 rows containing missing values or values outside the scale range
## (`geom_line()`).## Warning: Removed 90 rows containing missing values or values outside the scale range
## (`geom_point()`).## Warning: Removed 90 rows containing missing values or values outside the scale range
## (`geom_text()`).
#Extraigo los AUC por tiempo y modelo
score <- as.data.frame(score_comp$AUC$score)
# Gráfico
ggplot(score, aes(x = times, y = AUC, color = model)) +
geom_line(linewidth = 1.2) +
geom_point(size = 3) +
geom_hline(yintercept = 0.5, linetype = "dashed", alpha = 0.2) +
geom_hline(yintercept = 0.7, linetype = "dotted", alpha = 0.2) +
scale_color_manual(values = c("FG model" = "#E74C3C", "CSC model" = "#3498DB")) +
scale_y_continuous(limits = c(0.45, 0.95), breaks = seq(0.5, 0.8, 0.1)) +
labs(
title = "Discriminacion de Modelos en el Tiempo",
subtitle = "Linea punteada = Azar (0.5), Linea punteada = Adecuada (0.7)",
x = "Dias de seguimiento",
y = "AUC (Area bajo la curva ROC)",
color = "Modelo"
) +
theme_minimal() +
theme(legend.position = "bottom")## Warning: Removed 14 rows containing missing values or values outside the scale range
## (`geom_line()`).## Warning: Removed 14 rows containing missing values or values outside the scale range
## (`geom_point()`).
#Interpretamos que año a año la capacidad discriminativa de los modelos se mantiene alta, tanto en CSC como en FGR#análisis de sensibilidad
#Considerando a todos los que se perdieron antes de los 30 días como amputados
sad_sens<- sad
table(sad_sens$estado_unif)##
## amputacion cerro_herida fallecio persiste
## 58 232 14 102sad_sens <- sad_sens %>%
mutate(
evento_ampu = ifelse(
sad_sens$estado_unif == "persiste" & dias_out < 31,
1,
evento_ampu
)
)
table(sad_sens$evento_ampu)##
## 0 1 2
## 86 74 246mod_fyg_sens <- FGR(
formula = Hist(dias_out, evento_ampu) ~ san_elian + edad
, data = sad_sens,
cause = 1 )
mod_fyg_sens##
## Right-censored response of a competing.risks model
##
## No.Observations: 406
##
## Pattern:
##
## Cause event right.censored
## 1 74 0
## 2 246 0
## unknown 0 86
##
##
## Fine-Gray model: analysis of cause 1
##
## Competing Risks Regression
##
## Call:
## FGR(formula = Hist(dias_out, evento_ampu) ~ san_elian + edad,
## data = sad_sens, cause = 1)
##
## coef exp(coef) se(coef) z p-value
## san_elian 0.243 1.27 0.03190 7.61 0.000000000000027
## edad 0.021 1.02 0.00948 2.21 0.027000000000000
##
## exp(coef) exp(-coef) 2.5% 97.5%
## san_elian 1.27 0.784 1.2 1.36
## edad 1.02 0.979 1.0 1.04
##
## Num. cases = 406
## Pseudo Log-likelihood = -402
## Pseudo likelihood ratio test = 66.6 on 2 df,
##
## Convergence: TRUEMultiplica la isquemia para ver si mejora la calibración, pero no
#voy a multiplicar por 2 isquemia a ver si mejora la calibración del modelo
sadx4= sad
sadx4$isq_num <- as.numeric(as.character(sadx4$isq))
sadx4$zon_num <- as.numeric(as.character(sadx4$zon))
sadx4$inf_num <- as.numeric(as.character(sadx4$inf))
sadx4$ede_num <- as.numeric(as.character(sadx4$ede))
sadx4$local_num <- as.numeric(as.character(sadx4$local))
sadx4$topo_num <- as.numeric(as.character(sadx4$topo))
sadx4$prof_num <- as.numeric(as.character(sadx4$prof))
sadx4$area_num <- as.numeric(as.character(sadx4$area))
sadx4$fase_num <- as.numeric(as.character(sadx4$fase))
sadx4$neu_num <- as.numeric(as.character(sadx4$neu))
sadx4$isq2 <- sadx4$isq * 7## Warning in Ops.factor(sadx4$isq, 7): '*' no es significativo para factoressadx4$san_elian <- rowSums(
sadx4[, c("isq2", "local_num", "topo_num", "zon_num", "inf_num", "ede_num", "neu_num", "prof_num", "area_num", "fase_num")],
na.rm = TRUE
)
mod_fgx4 <- FGR(
Hist(dias_out, evento_ampu) ~ san_elian + edad,
data = sadx4,
cause = 1
)
describe(sadx4$san_elian)## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 406 17.1 3.71 17 17.1 4.45 8 26 18 0.01 -0.46 0.18describe(sad$san_elian)## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 406 18.3 4.09 18 18.3 4.45 8 29 21 -0.1 -0.38 0.2mod_fgx4##
## Right-censored response of a competing.risks model
##
## No.Observations: 406
##
## Pattern:
##
## Cause event right.censored
## 1 58 0
## 2 246 0
## unknown 0 102
##
##
## Fine-Gray model: analysis of cause 1
##
## Competing Risks Regression
##
## Call:
## FGR(formula = Hist(dias_out, evento_ampu) ~ san_elian + edad,
## data = sadx4, cause = 1)
##
## coef exp(coef) se(coef) z p-value
## san_elian 0.2973 1.35 0.0369 8.05 0.00000000000000089
## edad 0.0535 1.05 0.0103 5.17 0.00000024000000000
##
## exp(coef) exp(-coef) 2.5% 97.5%
## san_elian 1.35 0.743 1.25 1.45
## edad 1.05 0.948 1.03 1.08
##
## Num. cases = 406
## Pseudo Log-likelihood = -303
## Pseudo likelihood ratio test = 74.6 on 2 df,
##
## Convergence: TRUEset.seed(123)
sc_boot <- Score(
object = list("Modelo" = mod_fgx4),
formula = Hist(dias_out, evento_ampu) ~ 1,
data = sadx4,
times = 90,
metrics = c("auc", "brier"),
summary = "risks",
cause = 1,
split.method = "bootcv", # ← bootstrap
B = 200, # ← número de muestras bootstrap (200–500 ideal)
conf.int = TRUE,
plots = "calibration"
)## Running crossvalidation algorithm## Fitting the models in 200 learning datasets, then predicting the risks in validation datasets## | | | 0% | | | 1% | | | 2% | |= | 3% | |= | 4% | |= | 5% | |= | 6% | |= | 7% | |== | 8% | |== | 9% | |== | 10% | |== | 11% | |== | 12% | |=== | 13% | |=== | 14% | |=== | 15% | |=== | 16% | |=== | 17% | |==== | 18% | |==== | 19% | |==== | 20% | |==== | 21% | |==== | 22% | |===== | 23% | |===== | 24% | |===== | 25% | |===== | 26% | |===== | 27% | |====== | 28% | |====== | 29% | |====== | 30% | |====== | 31% | |====== | 32% | |======= | 33% | |======= | 34% | |======= | 35% | |======= | 36% | |======= | 37% | |======== | 38% | |======== | 39% | |======== | 40% | |======== | 41% | |======== | 42% | |========= | 43% | |========= | 44% | |========= | 45% | |========= | 46% | |========= | 47% | |========== | 48% | |========== | 49% | |========== | 50% | |========== | 51% | |========== | 52% | |=========== | 53% | |=========== | 54% | |=========== | 55% | |=========== | 56% | |=========== | 57% | |============ | 58% | |============ | 59% | |============ | 60% | |============ | 61% | |============ | 62% | |============= | 63% | |============= | 64% | |============= | 65% | |============= | 66% | |============= | 67% | |============== | 68% | |============== | 69% | |============== | 70% | |============== | 71% | |============== | 72% | |=============== | 73% | |=============== | 74% | |=============== | 75% | |=============== | 76% | |=============== | 77% | |================ | 78% | |================ | 79% | |================ | 80% | |================ | 81% | |================ | 82% | |================= | 83% | |================= | 84% | |================= | 85% | |================= | 86% | |================= | 87% | |================== | 88% | |================== | 89% | |================== | 90% | |================== | 91% | |================== | 92% | |=================== | 93% | |=================== | 94% | |=================== | 95% | |=================== | 96% | |=================== | 97% | |====================| 98% | |====================| 99% | |====================| 100%
## | | | 0%## Calculating the performance metrics in 200 validation data sets## | | | 1% | | | 2% | |= | 3% | |= | 4% | |= | 5% | |= | 6% | |= | 7% | |== | 8% | |== | 9% | |== | 10% | |== | 11% | |== | 12% | |=== | 13% | |=== | 14% | |=== | 15% | |=== | 16% | |=== | 17% | |==== | 18% | |==== | 19% | |==== | 20% | |==== | 21% | |==== | 22% | |===== | 23% | |===== | 24% | |===== | 25% | |===== | 26% | |===== | 27% | |====== | 28% | |====== | 29% | |====== | 30% | |====== | 31% | |====== | 32% | |======= | 33% | |======= | 34% | |======= | 35% | |======= | 36% | |======= | 37% | |======== | 38% | |======== | 39% | |======== | 40% | |======== | 41% | |======== | 42% | |========= | 43% | |========= | 44% | |========= | 45% | |========= | 46% | |========= | 47% | |========== | 48% | |========== | 49% | |========== | 50% | |========== | 51% | |========== | 52% | |=========== | 53% | |=========== | 54% | |=========== | 55% | |=========== | 56% | |=========== | 57% | |============ | 58% | |============ | 59% | |============ | 60% | |============ | 61% | |============ | 62% | |============= | 63% | |============= | 64% | |============= | 65% | |============= | 66% | |============= | 67% | |============== | 68% | |============== | 69% | |============== | 70% | |============== | 71% | |============== | 72% | |=============== | 73% | |=============== | 74% | |=============== | 75% | |=============== | 76% | |=============== | 77% | |================ | 78% | |================ | 79% | |================ | 80% | |================ | 81% | |================ | 82% | |================= | 83% | |================= | 84% | |================= | 85% | |================= | 86% | |================= | 87% | |================== | 88% | |================== | 89% | |================== | 90% | |================== | 91% | |================== | 92% | |=================== | 93% | |=================== | 94% | |=================== | 95% | |=================== | 96% | |=================== | 97% | |====================| 98% | |====================| 99% | |====================| 100%sc_boot$AUC##
## Results by model:
##
## model times AUC lower upper
## <fctr> <num> <char> <char> <char>
## 1: Modelo 90 82.178 75.211 89.281##
## NOTE: Values are multiplied by 100 and given in %.## NOTE: The higher AUC the better.sc_boot$Brier##
## Results by model:
##
## model times Brier lower upper
## <fctr> <num> <char> <char> <char>
## 1: Null model 90 11.293 8.470 14.260
## 2: Modelo 90 8.912 6.620 11.578
##
## Results of model comparisons:
##
## model times reference delta.Brier lower upper
## <fctr> <num> <fctr> <char> <char> <char>
## 1: Modelo 90 Null model -2.381 -3.831 -0.729##
## NOTE: Values are multiplied by 100 and given in %.## NOTE: The lower Brier the better.plotCalibration(sc_boot)
#calculadora para excel
library(tidyr)
#install.packages("writexl")
library(writexl)## Warning: package 'writexl' was built under R version 4.4.3#--------------------------------------------------
# 1. Definir rangos y tiempos
#--------------------------------------------------
times <- c(30, 60, 90, 120, 180)
rango_san_elian <- 6:30
rango_edad <- 30:85
#--------------------------------------------------
# 2. Crear todas las combinaciones posibles
#--------------------------------------------------
tabla_pred <- expand.grid(
san_elian = rango_san_elian,
edad = rango_edad
)
#--------------------------------------------------
# 3. Predecir riesgos con el modelo Fine-Gray
#--------------------------------------------------
pred <- predictRisk(
object = mod_fg,
newdata = tabla_pred,
times = times
)
#--------------------------------------------------
# 4. Pasar predicciones a data frame y nombrar columnas
#--------------------------------------------------
pred_df <- as.data.frame(pred)
colnames(pred_df) <- paste0("riesgo_", times)
#--------------------------------------------------
# 5. Unir combinaciones + predicciones
#--------------------------------------------------
tabla_final <- bind_cols(tabla_pred, pred_df)
#--------------------------------------------------
# 6. (Opcional) convertir riesgos a porcentaje
#--------------------------------------------------
tabla_final_pct <- tabla_final %>%
mutate(across(starts_with("riesgo_"), ~ round(.x , 1)))
#--------------------------------------------------
# 7. Ver primeras filas
#--------------------------------------------------
head(tabla_final_pct)## san_elian edad riesgo_30 riesgo_60 riesgo_90 riesgo_120 riesgo_180
## 1 6 30 0 0 0 0 0
## 2 7 30 0 0 0 0 0
## 3 8 30 0 0 0 0 0
## 4 9 30 0 0 0 0 0
## 5 10 30 0 0 0 0 0
## 6 11 30 0 0 0 0 0#--------------------------------------------------
# 8. Guardar en Excel
#--------------------------------------------------
write_xlsx(tabla_final_pct, "tabla_riesgo_san_elian_edad.xlsx")