RegLog
KarlaRamirez
2025-05-01
Regresión Logística
Los datos se obtuvieron de los expedientes del Instituto Nacional de Cardiología Dr. Ignacio Chávez en la Ciudad de México, incluyendo notas médicas, historias clínicas y resultados de exámenes de laboratorio y gabinete, de pacientes adultos con enfermedades cardiacas, atendidos entre enero y agosto de 2006 (Pale Carrión, 2006).
Muestra: 65 pacientes hospitalizados que cumplieron con los criterios de inclusión:
Criterios de inclusión utilizados:
Pacientes mayores de 18 años.
Diagnóstico confirmado de sepsis asociado a infecciones como neumonía, IVUS (infección de vías urinarias), celulitis, mediastinitis o bacteriemia.
Criterios de exclusión utilizados:
Pacientes que no cumplan con las definiciones clínicas de las infecciones mencionadas.
Pacientes embarazadas.
Aquellos bajo tratamiento inmunosupresor, con neutropenia no asociada a sepsis o antecedentes de neoplasia.
Pasos:
- Carga de los datos
library(readxl)
pacientes <- read_excel("D:/Estancia/Base_Proyecto.xlsx")
View(pacientes)Se importa la base de datos junto con la libreria necesaria.
- Conversion de la variable dependiente
pacientes$y_binaria <- ifelse(pacientes$ESTADO_PAC == "sobrevivientes", 1, 0)Se convierte la variable ESTADO_PAC a una versión binaria para utilizarla como variable respuesta en el modelo. Aquí, se codifica:
1 = paciente sobrevivió
2 = paciente NO sobrevivió
- Selección de variables predictoras
Se selccionaron 16 variables clínicas como predictoras y se convierte en una matriz.
X_vars <- c("EDAD", "HTO", "LEUCOS", "PLAQ", "GLUC",
"BUN", "CREAT", "BILIS", "ALB", "COL T", "TRIG",
"COL HDL", "COL LDL", "PCR ING", "PCR 24HRS", "PCR EGRE")
X <- as.matrix(pacientes[, X_vars])
X## EDAD HTO LEUCOS PLAQ GLUC BUN CREAT BILIS ALB COL T TRIG COL HDL
## [1,] 59 32.0 12.0 122 122 37.0 1.70 1.60 3.10 114 136 56
## [2,] 29 26.0 22.0 148 81 22.0 0.90 1.20 2.20 102 207 23
## [3,] 49 41.0 16.7 208 113 24.0 1.20 0.70 4.39 210 196 51
## [4,] 22 36.0 20.1 529 123 70.0 2.20 0.80 2.70 110 124 33
## [5,] 35 35.0 14.5 313 103 11.6 0.90 0.90 3.90 153 129 38
## [6,] 41 34.0 9.0 331 90 13.0 1.20 0.40 3.60 144 109 23
## [7,] 57 34.0 8.3 193 75 30.0 1.64 0.60 3.50 149 82 39
## [8,] 75 37.0 12.9 66 75 85.0 2.70 1.20 3.20 109 93 30
## [9,] 53 37.0 11.4 264 102 11.0 0.90 0.50 3.47 106 95 20
## [10,] 49 37.0 29.0 253 241 33.0 2.30 3.10 2.56 83 86 16
## [11,] 85 43.0 15.5 303 186 21.0 1.00 0.67 4.10 203 50 70
## [12,] 51 28.0 14.5 186 91 36.0 2.10 3.00 2.38 98 128 14
## [13,] 75 29.0 16.0 111 120 48.0 1.70 0.60 3.20 106 101 39
## [14,] 52 52.0 18.5 261 171 18.0 1.80 1.50 3.70 150 77 37
## [15,] 56 44.0 18.0 477 163 30.0 1.50 1.20 3.85 142 73 37
## [16,] 53 42.0 15.4 73 293 23.0 1.20 2.20 3.30 98 138 11
## [17,] 50 29.0 24.2 328 143 28.0 1.40 1.11 3.10 142 223 35
## [18,] 51 28.0 15.1 198 264 21.0 0.90 0.50 2.70 72 126 46
## [19,] 58 29.0 18.0 131 194 17.0 0.90 0.50 3.60 152 155 35
## [20,] 68 26.0 13.1 81 176 56.0 1.70 0.80 3.70 141 137 30
## [21,] 65 31.0 5.0 108 124 25.0 1.00 1.80 2.80 147 63 47
## [22,] 59 40.0 17.6 170 163 31.0 1.40 3.10 2.60 115 96 24
## [23,] 21 34.0 26.0 38 238 13.0 1.10 0.90 4.20 162 124 57
## [24,] 53 36.0 16.3 152 123 13.0 0.80 2.90 3.10 82 59 30
## [25,] 64 37.0 13.4 160 120 10.0 0.90 1.90 3.40 84 80 34
## [26,] 64 41.0 12.0 153 99 34.0 1.14 3.20 2.60 69 95 24
## [27,] 54 35.0 11.1 366 143 13.0 0.80 0.40 3.00 190 179 34
## [28,] 51 33.0 19.5 115 141 27.0 1.30 1.00 3.50 174 192 38
## [29,] 65 43.0 13.0 180 201 59.0 2.80 1.00 3.00 158 172 30
## [30,] 51 28.0 12.1 140 157 13.0 1.10 0.30 3.70 113 168 32
## [31,] 53 30.0 10.0 393 96 13.0 0.90 0.80 3.00 177 130 57
## [32,] 50 26.0 24.0 156 152 20.0 1.20 0.60 4.10 142 312 28
## [33,] 63 30.0 13.4 151 138 19.0 0.40 2.10 2.20 56 170 8
## [34,] 22 32.0 17.0 206 210 12.0 0.80 0.60 2.80 156 198 27
## [35,] 25 39.0 13.0 353 71 12.5 1.00 1.30 3.10 122 124 39
## [36,] 56 28.9 13.5 399 73 34.0 2.20 0.90 3.60 152 70 39
## [37,] 56 30.0 20.5 161 129 25.0 1.00 0.80 2.20 162 94 36
## [38,] 49 36.0 25.0 257 260 54.0 1.80 2.50 2.60 120 70 42
## [39,] 59 44.0 15.7 291 155 125.0 5.80 2.00 3.00 129 180 21
## [40,] 79 41.0 12.1 281 231 37.0 1.60 0.57 3.30 130 122 31
## [41,] 60 38.0 11.4 187 173 36.0 2.10 1.00 2.94 144 101 32
## [42,] 73 26.0 7.5 128 112 16.0 1.20 0.60 2.10 105 89 19
## [43,] 47 49.0 1.8 49 85 133.0 4.70 3.30 2.40 51 154 6
## [44,] 67 32.0 14.4 201 114 30.0 1.10 1.70 3.60 103 55 49
## [45,] 52 42.0 20.0 172 176 38.0 1.90 1.09 3.95 224 330 31
## [46,] 53 44.0 15.0 419 316 12.0 0.90 0.50 3.23 200 217 25
## [47,] 50 31.0 19.0 825 106 72.0 2.10 0.60 2.80 114 72 41
## [48,] 66 29.0 14.9 185 64 64.0 1.60 0.40 3.00 125 73 25
## [49,] 37 21.0 15.0 50 124 45.0 1.60 8.00 1.90 60 76 8
## [50,] 55 45.0 12.5 156 146 33.0 1.40 1.00 3.10 80 80 24
## [51,] 66 26.0 9.7 86 189 63.0 5.80 1.20 2.10 76 108 21
## [52,] 59 46.0 12.9 261 114 22.0 1.60 1.40 3.30 115 202 31
## [53,] 78 30.0 15.0 119 274 43.0 1.80 4.50 2.35 52 173 8
## [54,] 66 30.0 15.1 294 140 50.0 4.40 0.50 2.00 76 108 21
## [55,] 67 28.0 11.0 121 141 47.0 0.90 0.40 3.00 120 151 28
## [56,] 56 30.0 14.8 84 380 118.0 4.00 1.40 4.50 115 202 31
## [57,] 61 31.0 11.7 160 203 17.0 1.70 1.20 2.70 83 60 31
## [58,] 64 43.0 26.0 158 135 31.0 2.30 0.70 2.40 114 100 25
## [59,] 66 37.0 18.7 208 77 39.0 2.40 1.00 2.50 84 80 34
## [60,] 56 30.0 15.0 97 230 70.0 1.35 1.70 2.80 61 91 24
## [61,] 50 48.0 11.4 76 128 42.0 1.20 5.00 2.54 86 96 7
## [62,] 75 32.0 12.3 247 198 22.0 0.60 1.00 2.90 115 160 32
## [63,] 44 47.0 17.7 55 93 30.0 1.40 1.00 3.60 280 313 30
## [64,] 64 28.0 3.9 171 309 36.0 1.60 0.50 3.40 149 48 47
## [65,] 51 29.0 11.0 108 108 39.0 2.40 1.30 2.25 99 218 9
## COL LDL PCR ING PCR 24HRS PCR EGRE
## [1,] 30 14 6 9.4
## [2,] 73 195 422 422.0
## [3,] 119 112 88 17.0
## [4,] 67 136 147 10.5
## [5,] 89 65 228 7.0
## [6,] 99 25 26 27.0
## [7,] 93 71 86 20.0
## [8,] 60 38 52 10.0
## [9,] 67 53 30 8.0
## [10,] 49 215 83 109.0
## [11,] 123 6 32 32.0
## [12,] 60 227 230 26.0
## [13,] 46 42 75 75.0
## [14,] 97 88 110 31.0
## [15,] 90 131 106 69.0
## [16,] 59 74 117 19.0
## [17,] 62 336 201 15.0
## [18,] 95 84 46 33.0
## [19,] 86 77 110 70.0
## [20,] 83 100 89 43.0
## [21,] 87 21 14 7.0
## [22,] 71 151 136 54.0
## [23,] 52 88 113 36.0
## [24,] 40 79 76 31.0
## [25,] 34 98 90 34.0
## [26,] 71 136 111 52.0
## [27,] 100 26 14 14.0
## [28,] 96 184 172 95.0
## [29,] 110 132 134 17.0
## [30,] 100 254 212 11.0
## [31,] 46 114 99 47.0
## [32,] 42 88 72 11.0
## [33,] 15 92 77 15.0
## [34,] 89 89 96 31.0
## [35,] 72 275 326 24.0
## [36,] 99 92 88 33.0
## [37,] 107 59 47 29.0
## [38,] 65 17 18 8.0
## [39,] 72 70 73 40.0
## [40,] 74 195 215 113.0
## [41,] 91 35 30 30.0
## [42,] 68 136 164 164.0
## [43,] 14 284 256 284.0
## [44,] 82 179 67 82.0
## [45,] 127 45 263 263.0
## [46,] 131 68 99 152.0
## [47,] 58 129 135 38.0
## [48,] 88 19 29 16.0
## [49,] 36 57 130 130.0
## [50,] 40 117 115 42.0
## [51,] 33 251 246 84.0
## [52,] 43 164 160 160.0
## [53,] 15 335 269 269.0
## [54,] 33 166 76 76.0
## [55,] 92 98 82 115.0
## [56,] 43 151 183 183.0
## [57,] 40 131 112 112.0
## [58,] 73 302 324 119.0
## [59,] 34 117 117 117.0
## [60,] 70 232 157 157.0
## [61,] 59 66 66 66.0
## [62,] 100 101 82 31.0
## [63,] 187 89 85 85.0
## [64,] 92 115 105 40.0
## [65,] 46 99 198 154.0Se añadió una columna de (1s) para representar e intercepto de β0 del modelo.
# Agregamos el intercepto (columna de 1s)
X <- cbind(1, X) # tamaño: n x (p + 1)
X## EDAD HTO LEUCOS PLAQ GLUC BUN CREAT BILIS ALB COL T TRIG COL HDL
## [1,] 1 59 32.0 12.0 122 122 37.0 1.70 1.60 3.10 114 136 56
## [2,] 1 29 26.0 22.0 148 81 22.0 0.90 1.20 2.20 102 207 23
## [3,] 1 49 41.0 16.7 208 113 24.0 1.20 0.70 4.39 210 196 51
## [4,] 1 22 36.0 20.1 529 123 70.0 2.20 0.80 2.70 110 124 33
## [5,] 1 35 35.0 14.5 313 103 11.6 0.90 0.90 3.90 153 129 38
## [6,] 1 41 34.0 9.0 331 90 13.0 1.20 0.40 3.60 144 109 23
## [7,] 1 57 34.0 8.3 193 75 30.0 1.64 0.60 3.50 149 82 39
## [8,] 1 75 37.0 12.9 66 75 85.0 2.70 1.20 3.20 109 93 30
## [9,] 1 53 37.0 11.4 264 102 11.0 0.90 0.50 3.47 106 95 20
## [10,] 1 49 37.0 29.0 253 241 33.0 2.30 3.10 2.56 83 86 16
## [11,] 1 85 43.0 15.5 303 186 21.0 1.00 0.67 4.10 203 50 70
## [12,] 1 51 28.0 14.5 186 91 36.0 2.10 3.00 2.38 98 128 14
## [13,] 1 75 29.0 16.0 111 120 48.0 1.70 0.60 3.20 106 101 39
## [14,] 1 52 52.0 18.5 261 171 18.0 1.80 1.50 3.70 150 77 37
## [15,] 1 56 44.0 18.0 477 163 30.0 1.50 1.20 3.85 142 73 37
## [16,] 1 53 42.0 15.4 73 293 23.0 1.20 2.20 3.30 98 138 11
## [17,] 1 50 29.0 24.2 328 143 28.0 1.40 1.11 3.10 142 223 35
## [18,] 1 51 28.0 15.1 198 264 21.0 0.90 0.50 2.70 72 126 46
## [19,] 1 58 29.0 18.0 131 194 17.0 0.90 0.50 3.60 152 155 35
## [20,] 1 68 26.0 13.1 81 176 56.0 1.70 0.80 3.70 141 137 30
## [21,] 1 65 31.0 5.0 108 124 25.0 1.00 1.80 2.80 147 63 47
## [22,] 1 59 40.0 17.6 170 163 31.0 1.40 3.10 2.60 115 96 24
## [23,] 1 21 34.0 26.0 38 238 13.0 1.10 0.90 4.20 162 124 57
## [24,] 1 53 36.0 16.3 152 123 13.0 0.80 2.90 3.10 82 59 30
## [25,] 1 64 37.0 13.4 160 120 10.0 0.90 1.90 3.40 84 80 34
## [26,] 1 64 41.0 12.0 153 99 34.0 1.14 3.20 2.60 69 95 24
## [27,] 1 54 35.0 11.1 366 143 13.0 0.80 0.40 3.00 190 179 34
## [28,] 1 51 33.0 19.5 115 141 27.0 1.30 1.00 3.50 174 192 38
## [29,] 1 65 43.0 13.0 180 201 59.0 2.80 1.00 3.00 158 172 30
## [30,] 1 51 28.0 12.1 140 157 13.0 1.10 0.30 3.70 113 168 32
## [31,] 1 53 30.0 10.0 393 96 13.0 0.90 0.80 3.00 177 130 57
## [32,] 1 50 26.0 24.0 156 152 20.0 1.20 0.60 4.10 142 312 28
## [33,] 1 63 30.0 13.4 151 138 19.0 0.40 2.10 2.20 56 170 8
## [34,] 1 22 32.0 17.0 206 210 12.0 0.80 0.60 2.80 156 198 27
## [35,] 1 25 39.0 13.0 353 71 12.5 1.00 1.30 3.10 122 124 39
## [36,] 1 56 28.9 13.5 399 73 34.0 2.20 0.90 3.60 152 70 39
## [37,] 1 56 30.0 20.5 161 129 25.0 1.00 0.80 2.20 162 94 36
## [38,] 1 49 36.0 25.0 257 260 54.0 1.80 2.50 2.60 120 70 42
## [39,] 1 59 44.0 15.7 291 155 125.0 5.80 2.00 3.00 129 180 21
## [40,] 1 79 41.0 12.1 281 231 37.0 1.60 0.57 3.30 130 122 31
## [41,] 1 60 38.0 11.4 187 173 36.0 2.10 1.00 2.94 144 101 32
## [42,] 1 73 26.0 7.5 128 112 16.0 1.20 0.60 2.10 105 89 19
## [43,] 1 47 49.0 1.8 49 85 133.0 4.70 3.30 2.40 51 154 6
## [44,] 1 67 32.0 14.4 201 114 30.0 1.10 1.70 3.60 103 55 49
## [45,] 1 52 42.0 20.0 172 176 38.0 1.90 1.09 3.95 224 330 31
## [46,] 1 53 44.0 15.0 419 316 12.0 0.90 0.50 3.23 200 217 25
## [47,] 1 50 31.0 19.0 825 106 72.0 2.10 0.60 2.80 114 72 41
## [48,] 1 66 29.0 14.9 185 64 64.0 1.60 0.40 3.00 125 73 25
## [49,] 1 37 21.0 15.0 50 124 45.0 1.60 8.00 1.90 60 76 8
## [50,] 1 55 45.0 12.5 156 146 33.0 1.40 1.00 3.10 80 80 24
## [51,] 1 66 26.0 9.7 86 189 63.0 5.80 1.20 2.10 76 108 21
## [52,] 1 59 46.0 12.9 261 114 22.0 1.60 1.40 3.30 115 202 31
## [53,] 1 78 30.0 15.0 119 274 43.0 1.80 4.50 2.35 52 173 8
## [54,] 1 66 30.0 15.1 294 140 50.0 4.40 0.50 2.00 76 108 21
## [55,] 1 67 28.0 11.0 121 141 47.0 0.90 0.40 3.00 120 151 28
## [56,] 1 56 30.0 14.8 84 380 118.0 4.00 1.40 4.50 115 202 31
## [57,] 1 61 31.0 11.7 160 203 17.0 1.70 1.20 2.70 83 60 31
## [58,] 1 64 43.0 26.0 158 135 31.0 2.30 0.70 2.40 114 100 25
## [59,] 1 66 37.0 18.7 208 77 39.0 2.40 1.00 2.50 84 80 34
## [60,] 1 56 30.0 15.0 97 230 70.0 1.35 1.70 2.80 61 91 24
## [61,] 1 50 48.0 11.4 76 128 42.0 1.20 5.00 2.54 86 96 7
## [62,] 1 75 32.0 12.3 247 198 22.0 0.60 1.00 2.90 115 160 32
## [63,] 1 44 47.0 17.7 55 93 30.0 1.40 1.00 3.60 280 313 30
## [64,] 1 64 28.0 3.9 171 309 36.0 1.60 0.50 3.40 149 48 47
## [65,] 1 51 29.0 11.0 108 108 39.0 2.40 1.30 2.25 99 218 9
## COL LDL PCR ING PCR 24HRS PCR EGRE
## [1,] 30 14 6 9.4
## [2,] 73 195 422 422.0
## [3,] 119 112 88 17.0
## [4,] 67 136 147 10.5
## [5,] 89 65 228 7.0
## [6,] 99 25 26 27.0
## [7,] 93 71 86 20.0
## [8,] 60 38 52 10.0
## [9,] 67 53 30 8.0
## [10,] 49 215 83 109.0
## [11,] 123 6 32 32.0
## [12,] 60 227 230 26.0
## [13,] 46 42 75 75.0
## [14,] 97 88 110 31.0
## [15,] 90 131 106 69.0
## [16,] 59 74 117 19.0
## [17,] 62 336 201 15.0
## [18,] 95 84 46 33.0
## [19,] 86 77 110 70.0
## [20,] 83 100 89 43.0
## [21,] 87 21 14 7.0
## [22,] 71 151 136 54.0
## [23,] 52 88 113 36.0
## [24,] 40 79 76 31.0
## [25,] 34 98 90 34.0
## [26,] 71 136 111 52.0
## [27,] 100 26 14 14.0
## [28,] 96 184 172 95.0
## [29,] 110 132 134 17.0
## [30,] 100 254 212 11.0
## [31,] 46 114 99 47.0
## [32,] 42 88 72 11.0
## [33,] 15 92 77 15.0
## [34,] 89 89 96 31.0
## [35,] 72 275 326 24.0
## [36,] 99 92 88 33.0
## [37,] 107 59 47 29.0
## [38,] 65 17 18 8.0
## [39,] 72 70 73 40.0
## [40,] 74 195 215 113.0
## [41,] 91 35 30 30.0
## [42,] 68 136 164 164.0
## [43,] 14 284 256 284.0
## [44,] 82 179 67 82.0
## [45,] 127 45 263 263.0
## [46,] 131 68 99 152.0
## [47,] 58 129 135 38.0
## [48,] 88 19 29 16.0
## [49,] 36 57 130 130.0
## [50,] 40 117 115 42.0
## [51,] 33 251 246 84.0
## [52,] 43 164 160 160.0
## [53,] 15 335 269 269.0
## [54,] 33 166 76 76.0
## [55,] 92 98 82 115.0
## [56,] 43 151 183 183.0
## [57,] 40 131 112 112.0
## [58,] 73 302 324 119.0
## [59,] 34 117 117 117.0
## [60,] 70 232 157 157.0
## [61,] 59 66 66 66.0
## [62,] 100 101 82 31.0
## [63,] 187 89 85 85.0
## [64,] 92 115 105 40.0
## [65,] 46 99 198 154.0La variable respuesta es el vector y, que contiene 0s y 1s según la condición clínica.
# Variable respuesta
y <- pacientes$y_binaria
y## [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [39] 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0- Función de regresión logística IRLS
IRLS (Mínimos cuadrados reponderados iterativamente) es un algoritmo para encontrar los coeficientes óptimos β en modelos como la regresión logística. Se basa en el método de Newton-Raphson, que es un método iterativo para encontrar el máximo de la función log-verosimilitud.
En el caso de la regresión logística
\[ \text{logit}(\pi_i) = \log\left(\frac{\pi_i}{1 - \pi_i}\right) = x_i^{T} \beta \]
Dentro de esta función
eta: calcula el predictor lineal \(\eta_i = x_i^{T} \beta\)pi: aplica la función logística \(\pi_i = \frac{1}{1 + e^{-\eta_i}}\)W: es la matriz diagonal de pesos \(\pi_i (1 - \pi_i)\)g: gradiente de la log-verosimilitudH: matriz Hessiana (segunda derivada de la log-verosimilitud)beta: se actualiza iterativamente hasta converger
logistic_IRLS <- function(X, y, tol = 1e-3, max_iter = 100) {
n <- nrow(X) #número de observaciones
k <- ncol(X) #el número de variables predictoras
beta <- rep(0, k) # Inicializamos beta en 0
iter <- 0 #Establece un contador de iteraciones
gnorm <- 1 #Inicializa el valor de la norma del gradiente
while (gnorm > tol && iter < max_iter) {
eta <- X %*% beta # X * beta; Calcula el predictor lineal
pi <- 1 / (1 + exp(-eta)) # Función logística; Calcula la probabilidad estimada pi
W <- diag(as.vector(pi * (1 - pi))) # Calcula la matriz de pesos;
#W es una matriz diagonal con los pesos derivados de la varianza binomial.
g <- t(X) %*% (y - pi) # Gradiente; Este gradiente nos dice en qué dirección mejorar la log-verosimilitud.
H <- t(X) %*% W %*% X # Hessiana; Es una aproximación de la curvatura de la log-verosimilitud.
beta <- beta + solve(H, g) # Actualización de beta
gnorm <- sqrt(sum(g^2)) # Norma del gradiente;
iter <- iter + 1
}
list(beta = beta, iter = iter, converged = (gnorm <= tol))
}
#Este procedimiento se repite hasta que los coeficientes dejen de cambiar significativamente.Estimación del modelo
Se ejecuta la funcion antes mencionada, para estimar los coeficientes del modelo. El resultado es el factor
betaque contiene los coeficientes asociados a cada variable predictora.
resultado <- logistic_IRLS(X, y) # Muestra el vector de coeficientes ajustados del modelo
coeficientes <- resultado$beta # el vector de coeficientes estimados del modelo desde el objeto resultado
names(coeficientes) <- c("Intercepto", X_vars)
print(coeficientes)## [,1]
## 1.314776e+00
## EDAD -5.417827e-02
## HTO -1.530043e-02
## LEUCOS 1.991209e-01
## PLAQ -2.547096e-03
## GLUC -4.975364e-05
## BUN -3.669661e-02
## CREAT 3.678109e-01
## BILIS 1.641311e-02
## ALB 1.056151e+00
## COL T -1.846474e-02
## TRIG 5.914577e-03
## COL HDL 6.831901e-02
## COL LDL -3.305965e-03
## PCR ING -1.764314e-02
## PCR 24HRS 1.718628e-02
## PCR EGRE -2.994335e-02
## attr(,"names")
## [1] "Intercepto" "EDAD" "HTO" "LEUCOS" "PLAQ"
## [6] "GLUC" "BUN" "CREAT" "BILIS" "ALB"
## [11] "COL T" "TRIG" "COL HDL" "COL LDL" "PCR ING"
## [16] "PCR 24HRS" "PCR EGRE"Interpretación de los coeficientes
Los coeficientes estimados del modelo fueron transformados a odds ratios (\(e^{\beta_j}\)) para facilitar la interpretación. A continuación, se presentan algunos ejemplos destacados:
- EDAD (\(OR = 0.947\)): Por cada año adicional, las probabilidades de supervivencia disminuyen un 5.41%.
- LEUCOS (\(OR = 1.220\)): Cada aumento en la cuenta de leucocitos incrementa las probabilidades de sobrevivir en un 1.99%.
odds_ratios <- exp(coeficientes)
round(odds_ratios, 3)## [,1]
## 3.724
## EDAD 0.947
## HTO 0.985
## LEUCOS 1.220
## PLAQ 0.997
## GLUC 1.000
## BUN 0.964
## CREAT 1.445
## BILIS 1.017
## ALB 2.875
## COL T 0.982
## TRIG 1.006
## COL HDL 1.071
## COL LDL 0.997
## PCR ING 0.983
## PCR 24HRS 1.017
## PCR EGRE 0.971
## attr(,"names")
## [1] "Intercepto" "EDAD" "HTO" "LEUCOS" "PLAQ"
## [6] "GLUC" "BUN" "CREAT" "BILIS" "ALB"
## [11] "COL T" "TRIG" "COL HDL" "COL LDL" "PCR ING"
## [16] "PCR 24HRS" "PCR EGRE"#Los coeficientes 𝛽i representan los log-odds.
#Al aplicar la exponencial, obtenemos odds ratios:
#Un OR > 1 indica mayor riesgo; OR < 1, menor riesgo.Evaluación del modelo
Se generó las predicciones del mdelo (probabilidades de supervivencia) y se clasifican usando un umbral de 0.5. Se obtiene una matriz de contigencia (matriz de confución) para comparar las predicciones con datos reales.
Evaluación del modelo: Matriz de confusión
Se compararon las predicciones del modelo con los valores reales utilizando una matriz de confusión. Los resultados fueron:
| Real: No sobrevivió (0) | Real: Sobrevivió (1) | |
|---|---|---|
| Predijo: 0 | 15 (Verdaderos negativos) | 5 (Falsos negativos) |
| Predijo: 1 | 9 (Falsos positivos) | 36 (Verdaderos positivos) |
pi_hat <- 1 / (1 + exp(-X %*% coeficientes))
pred_clase <- ifelse(pi_hat > 0.5, 1, 0)
table(Predicho = pred_clase, Real = y)## Real
## Predicho 0 1
## 0 15 5
## 1 9 36#Calcula las probabilidades predichas.
#Clasifica según umbral de 0.5.
#Se construye matriz de confusión para comparar valores reales vs. predichos.- Métricas de desempeño del modelo
A partir de esta matriz se calcularon las siguientes métricas de desempeño:
- Precisión total: 78.46%
- Sensibilidad (recall, sobrevivientes): 87.8%
- Especificidad (no sobrevivientes): 62.5%
- Precisión positiva: 80.0%
El modelo presenta buen rendimiento para predecir pacientes que sobrevivieron, aunque con menor capacidad para detectar correctamente a quienes no lo hicieron.
# Definir los valores manualmente según la matriz
VP <- 15 # Verdaderos negativos (predijo 0 y era 0)
FN <- 5 # Falsos negativos (predijo 0 pero era 1)
FP <- 9 # Falsos positivos (predijo 1 pero era 0)
VN <- 36 # Verdaderos positivos (predijo 1 y era 1)
# Exactitud (Accuracy)
accuracy <- (VP + VN) / (VP + VN + FP + FN)
# Sensibilidad (Recall para clase 1)
sensibilidad <- VN / (VN + FN)
# Especificidad (para clase 0)
especificidad <- VP / (VP + FP)
# Precisión (Precision para clase 1)
precision <- VN / (VN + FP)
# Mostrar resultados
cat("Exactitud: ", round(accuracy * 100, 2), "%\n")## Exactitud: 78.46 %cat("Sensibilidad (Recall): ", round(sensibilidad * 100, 2), "%\n")## Sensibilidad (Recall): 87.8 %cat("Especificidad: ", round(especificidad * 100, 2), "%\n")## Especificidad: 62.5 %cat("Precisión: ", round(precision * 100, 2), "%\n")## Precisión: 80 %- Gráfico del valor p
Este gráfico permite evaluar la significancia estadística del poder discriminante de tu modelo (es decir, qué tan buena es la exactitud del modelo para clasificar) sin depender directamente del valor-p de cada coeficiente del modelo.
Prueba de hipótesis sobre el poder discriminante del modelo
Se plantea la siguiente prueba de hipótesis:
- \(H_0\): el modelo clasifica al azar, es decir, \(\gamma = 0.5\).
- \(H_1\): el modelo clasifica mejor que al azar, es decir, \(\gamma > 0.5\).
Se utiliza la estadística Q de Press:
\[ Q = n(2\gamma - 1)^2 \sim \chi^2_1 \]
Con una exactitud observada de 78.46% y \(n
= 65\), se obtiene un valor-p de aproximadamente \(4.5 \times 10^{-6}\).
Dado que \(p < 0.05\), se rechaza la
hipótesis nula y se concluye que el modelo tiene un poder de
clasificación significativamente superior al azar.
# Parámetros
n <- 65 # número total de observaciones
powers <- seq(0.5, 1, length.out = 1000) # rango de poder discriminante
q_vals <- n * (2 * powers - 1)^2 # estadística Q
p_vals <- 1 - pchisq(q_vals, df = 1) # valor-p asociado (distribución chi-cuadrado con 1 g.l.)
# Exactitud observada en tu modelo (de la matriz de confusión)
exactitud_observada <- (15 + 36) / (15 + 5 + 9 + 36) # = 0.7846 aprox.
# Cargar ggplot2
library(ggplot2)
# Crear dataframe
df <- data.frame(power = powers, pvalue = p_vals)
# Graficar
ggplot(df, aes(x = power, y = pvalue)) +
geom_line(color = "steelblue", size = 1) +
geom_vline(xintercept = exactitud_observada, color = "red", linetype = "dashed") +
annotate("text", x = exactitud_observada + 0.02, y = 0.2,
label = paste0("Exactitud observada = ", round(exactitud_observada, 4)),
color = "red", angle = 90, vjust = -0.5, hjust = 0) +
labs(
title = "Valor-p de la estadística Q de Press",
subtitle = "Evaluación de la significancia del clasificador",
x = expression(Poder~discriminante~(gamma)),
y = "Valor-p"
) +
theme_minimal()## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
En el gráfico se observa cómo el valor-p disminuye a medida que la exactitud del modelo aumenta. La línea roja vertical representa la exactitud observada del modelo (78.46%). Como esta se encuentra en una zona donde el valor-p es extremadamente bajo (aproximadamente 0.0000045), se concluye que la capacidad de clasificación del modelo es estadísticamente significativa.
Por lo tanto, se rechaza la hipótesis nula de clasificación aleatoria y se confirma que el modelo tiene un poder discriminante real y útil para la detección de sepsis.
- Comparación con la funcion de R
modelo_glm <- glm(y_binaria ~ ., data = pacientes[, c("y_binaria", X_vars)], family = binomial)
summary(modelo_glm)$coefficients## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.314775e+00 4.900538725 0.268291911 0.788474632
## EDAD -5.417826e-02 0.055446850 -0.977120594 0.328509453
## HTO -1.530044e-02 0.062137322 -0.246235851 0.805499664
## LEUCOS 1.991209e-01 0.097813781 2.035714422 0.041779040
## PLAQ -2.547096e-03 0.003227373 -0.789216310 0.429985591
## GLUC -4.975415e-05 0.006170826 -0.008062802 0.993566884
## BUN -3.669661e-02 0.029248251 -1.254660078 0.209602184
## CREAT 3.678109e-01 0.618058082 0.595107290 0.551771756
## BILIS 1.641313e-02 0.374977057 0.043771037 0.965086914
## ALB 1.056151e+00 1.011104985 1.044551551 0.296230282
## `COL T` -1.846474e-02 0.026700752 -0.691543654 0.489223957
## TRIG 5.914576e-03 0.010315426 0.573371943 0.566392883
## `COL HDL` 6.831900e-02 0.062367059 1.095434008 0.273326542
## `COL LDL` -3.305967e-03 0.025069312 -0.131873072 0.895084688
## `PCR ING` -1.764314e-02 0.010309912 -1.711279444 0.087029539
## `PCR 24HRS` 1.718628e-02 0.009940808 1.728861074 0.083833962
## `PCR EGRE` -2.994334e-02 0.010302170 -2.906508379 0.003654871