TH completo Carro Vinicki
2025-07-09
##Variables como factor y cambio de nombres
train <- train %>% mutate_at(vars(Died, gender, hypertension, a_fib, CHD, diabetes, anemia, depression, Hyperlipemia, CKD, COPD ),~as.factor(.))
train$heart=train$`heart rate`
train$syst=train$`Systolic blood pressure`
train$dias=train$`Diastolic blood pressure`
train$sat=train$`SP O2`
train$magn=train$`Magnesium ion`
train$probnp=train$`NT-proBNP`
train$resp=train$`Respiratory rate`
train$ureanit=train$`Urea nitrogen`
train$urine=train$`Urine output`
train$sodio=train$`Blood sodium`
train$lactico=train$`Lactic acid`
train$calcio=train$`Blood calcium`
train$potasio=train$`Blood potassium`
train <- train %>% select(-`heart rate`)
train <- train %>% select(-`Systolic blood pressure`)
train <- train %>% select(-`Diastolic blood pressure`,- `Respiratory rate`, -`Magnesium ion`)
train <- train %>% select(-`Urea nitrogen`, -`Urine output`, -`Blood sodium`)
train <- train %>% select(-`Blood calcium`, -`Blood potassium`, -`Lactic acid`, -`NT-proBNP`)##Hago base nueva DIED0, donde la variable muerte=1 es muertos y muerte=0 es vivos al revĆ©s de Died segĆŗn lo hablado con Diego
DIED0=train
DIED0<- DIED0 %>%
mutate(muerte = as.factor(case_when(
Died== 1 ~ "0", Died == 0 ~ "1")))
levels(DIED0$muerte)## [1] "0" "1"describe(DIED0$muerte)## DIED0$muerte
## n missing distinct
## 824 0 2
##
## Value 0 1
## Frequency 108 716
## Proportion 0.131 0.869describe(train$Died)## train$Died
## n missing distinct
## 824 0 2
##
## Value 0 1
## Frequency 716 108
## Proportion 0.869 0.131##Tabla 1a Antecedentes Tabla 1b Continuas
#Creo variable muerte 1 si/no para que la tabla se vea mejor
DIED0<- DIED0 %>% select(-`Died`)
DIED0$muerte1=DIED0$muerte
DIED0$muerte1 <- factor(DIED0$muerte1,
levels = c(0, 1),
labels = c("No muerte", "Muerte"))
catvars = c( "gender", "hypertension", "a_fib", "CHD", "diabetes", "anemia", "depression", "Hyperlipemia", "CKD", "COPD")
vars = c("gender", "hypertension", "a_fib", "CHD", "diabetes", "anemia", "depression", "Hyperlipemia", "CKD", "COPD")
tabla_1a <- CreateTableOne(vars = vars, strata = "muerte1", factorVars = catvars, data = DIED0)
print(tabla_1a)## Stratified by muerte1
## No muerte Muerte p test
## n 108 716
## gender = 2 (%) 52 (48.1) 382 (53.4) 0.365
## hypertension = 1 (%) 67 (62.0) 521 (72.8) 0.029
## a_fib = 1 (%) 64 (59.3) 313 (43.7) 0.004
## CHD = 1 (%) 10 ( 9.3) 66 ( 9.2) 1.000
## diabetes = 1 (%) 42 (38.9) 310 (43.3) 0.448
## anemia = 1 (%) 22 (20.4) 266 (37.2) 0.001
## depression = 1 (%) 8 ( 7.4) 90 (12.6) 0.166
## Hyperlipemia = 1 (%) 39 (36.1) 278 (38.8) 0.664
## CKD = 1 (%) 24 (22.2) 276 (38.5) 0.001
## COPD = 1 (%) 3 ( 2.8) 64 ( 8.9) 0.046library(dplyr)
kableone(tabla_1a)| No muerte | Muerte | p | test | |
|---|---|---|---|---|
| n | 108 | 716 | ||
| gender = 2 (%) | 52 (48.1) | 382 (53.4) | 0.365 | |
| hypertension = 1 (%) | 67 (62.0) | 521 (72.8) | 0.029 | |
| a_fib = 1 (%) | 64 (59.3) | 313 (43.7) | 0.004 | |
| CHD = 1 (%) | 10 ( 9.3) | 66 ( 9.2) | 1.000 | |
| diabetes = 1 (%) | 42 (38.9) | 310 (43.3) | 0.448 | |
| anemia = 1 (%) | 22 (20.4) | 266 (37.2) | 0.001 | |
| depression = 1 (%) | 8 ( 7.4) | 90 (12.6) | 0.166 | |
| Hyperlipemia = 1 (%) | 39 (36.1) | 278 (38.8) | 0.664 | |
| CKD = 1 (%) | 24 (22.2) | 276 (38.5) | 0.001 | |
| COPD = 1 (%) | 3 ( 2.8) | 64 ( 8.9) | 0.046 |
catvars = c( "gender", "hypertension", "a_fib", "CHD", "diabetes", "anemia", "depression", "Hyperlipemia", "CKD", "COPD")
#edad y enojo son continuas
vars = c("age","temperature", "hematocrit","Leucocyte","Platelets","INR","Creatinine","glucose","Chloride","PH","Bicarbonate","PCO2","heart","syst","dias","sat","magn","probnp","resp","ureanit","urine","sodio","lactico", "calcio","potasio", "SP O2")
tabla_2a <- CreateTableOne(vars = vars, strata = "muerte1", factorVars = catvars, data = DIED0)
print(tabla_2a)## Stratified by muerte1
## No muerte Muerte p test
## n 108 716
## age (mean (SD)) 75.98 (14.01) 73.76 (13.01) 0.102
## temperature (mean (SD)) 36.48 (0.73) 36.70 (0.57) <0.001
## hematocrit (mean (SD)) 31.97 (5.43) 31.91 (5.11) 0.913
## Leucocyte (mean (SD)) 13.84 (7.96) 10.37 (4.52) <0.001
## Platelets (mean (SD)) 214.03 (133.44) 248.72 (108.03) 0.003
## INR (mean (SD)) 2.02 (1.40) 1.58 (0.75) <0.001
## Creatinine (mean (SD)) 1.84 (1.20) 1.60 (1.21) 0.046
## glucose (mean (SD)) 151.15 (52.98) 148.52 (50.57) 0.617
## Chloride (mean (SD)) 103.36 (6.04) 102.12 (5.28) 0.026
## PH (mean (SD)) 7.36 (0.08) 7.38 (0.06) <0.001
## Bicarbonate (mean (SD)) 23.48 (5.32) 27.41 (5.01) <0.001
## PCO2 (mean (SD)) 43.16 (12.22) 45.49 (11.54) 0.053
## heart (mean (SD)) 88.24 (14.68) 83.80 (15.83) 0.006
## syst (mean (SD)) 110.39 (16.75) 119.21 (17.05) <0.001
## dias (mean (SD)) 57.67 (9.21) 60.00 (10.84) 0.034
## sat (mean (SD)) 95.90 (3.01) 96.35 (2.14) 0.054
## magn (mean (SD)) 2.18 (0.31) 2.10 (0.24) 0.004
## probnp (mean (SD)) 15554.25 (14867.09) 10814.14 (12818.97) <0.001
## resp (mean (SD)) 21.68 (4.27) 20.51 (3.93) 0.005
## ureanit (mean (SD)) 48.71 (28.71) 34.83 (20.81) <0.001
## urine (mean (SD)) 1251.19 (1042.46) 1999.65 (1271.35) <0.001
## sodio (mean (SD)) 138.18 (5.35) 139.02 (3.97) 0.052
## lactico (mean (SD)) 2.40 (1.36) 1.71 (0.75) <0.001
## calcio (mean (SD)) 8.23 (0.63) 8.54 (0.55) <0.001
## potasio (mean (SD)) 4.34 (0.58) 4.15 (0.38) <0.001
## SP O2 (mean (SD)) 95.90 (3.01) 96.35 (2.14) 0.054kableone(tabla_2a)| No muerte | Muerte | p | test | |
|---|---|---|---|---|
| n | 108 | 716 | ||
| age (mean (SD)) | 75.98 (14.01) | 73.76 (13.01) | 0.102 | |
| temperature (mean (SD)) | 36.48 (0.73) | 36.70 (0.57) | <0.001 | |
| hematocrit (mean (SD)) | 31.97 (5.43) | 31.91 (5.11) | 0.913 | |
| Leucocyte (mean (SD)) | 13.84 (7.96) | 10.37 (4.52) | <0.001 | |
| Platelets (mean (SD)) | 214.03 (133.44) | 248.72 (108.03) | 0.003 | |
| INR (mean (SD)) | 2.02 (1.40) | 1.58 (0.75) | <0.001 | |
| Creatinine (mean (SD)) | 1.84 (1.20) | 1.60 (1.21) | 0.046 | |
| glucose (mean (SD)) | 151.15 (52.98) | 148.52 (50.57) | 0.617 | |
| Chloride (mean (SD)) | 103.36 (6.04) | 102.12 (5.28) | 0.026 | |
| PH (mean (SD)) | 7.36 (0.08) | 7.38 (0.06) | <0.001 | |
| Bicarbonate (mean (SD)) | 23.48 (5.32) | 27.41 (5.01) | <0.001 | |
| PCO2 (mean (SD)) | 43.16 (12.22) | 45.49 (11.54) | 0.053 | |
| heart (mean (SD)) | 88.24 (14.68) | 83.80 (15.83) | 0.006 | |
| syst (mean (SD)) | 110.39 (16.75) | 119.21 (17.05) | <0.001 | |
| dias (mean (SD)) | 57.67 (9.21) | 60.00 (10.84) | 0.034 | |
| sat (mean (SD)) | 95.90 (3.01) | 96.35 (2.14) | 0.054 | |
| magn (mean (SD)) | 2.18 (0.31) | 2.10 (0.24) | 0.004 | |
| probnp (mean (SD)) | 15554.25 (14867.09) | 10814.14 (12818.97) | <0.001 | |
| resp (mean (SD)) | 21.68 (4.27) | 20.51 (3.93) | 0.005 | |
| ureanit (mean (SD)) | 48.71 (28.71) | 34.83 (20.81) | <0.001 | |
| urine (mean (SD)) | 1251.19 (1042.46) | 1999.65 (1271.35) | <0.001 | |
| sodio (mean (SD)) | 138.18 (5.35) | 139.02 (3.97) | 0.052 | |
| lactico (mean (SD)) | 2.40 (1.36) | 1.71 (0.75) | <0.001 | |
| calcio (mean (SD)) | 8.23 (0.63) | 8.54 (0.55) | <0.001 | |
| potasio (mean (SD)) | 4.34 (0.58) | 4.15 (0.38) | <0.001 | |
| SP O2 (mean (SD)) | 95.90 (3.01) | 96.35 (2.14) | 0.054 |
##Modelo univariado
DIED0 %>%
tbl_uvregression(
method = glm,
y = muerte1,
method.args = list(family = binomial),
exponentiate = TRUE,
pvalue_fun = ~style_pvalue(.x, digits = 2)) %>% bold_p() %>% bold_labels()## There was a warning constructing the model for variable "muerte". See message
## below.
## ! glm.fit: algorithm did not converge
## There was a warning running `tbl_regression()` for variable "muerte". See
## message below.
## ! glm.fit: fitted probabilities numerically 0 or 1 occurred, glm.fit: fitted
## probabilities numerically 0 or 1 occurred, glm.fit: fitted probabilities
## numerically 0 or 1 occurred, glm.fit: fitted probabilities numerically 0 or 1
## occurred, glm.fit: fitted probabilities numerically 0 or 1 occurred, glm.fit:
## fitted probabilities numerically 0 or 1 occurred, glm.fit: fitted
## probabilities numerically 0 or 1 occurred, glm.fit: fitted probabilities
## numerically 0 or 1 occurred, glm.fit: fitted probabilities numerically 0 or 1
## occurred, glm.fit: fitted probabilities numerically 0 or 1 occurred, glm.fit:
## fitted probabilities numerically 0 or 1 occurred, glm.fit: fitted
## probabilities numerically 0 or 1 occurred, glm.fit: fitted probabilities
## numerically 0 or 1 occurred, glm.fit: fitted probabilities numerically 0 or 1
## occurred, glm.fit: fitted probabilities numerically 0 or 1 occurred, glm.fit:
## fitted probabilities numerically 0 or 1 occurred, glm.fit: fitted
## probabilities numerically 0 or 1 occurred, glm.fit: fitted probabilities
## numerically 0 or 1 occurred, glm.fit: fitted probabilities numerically 0 or 1
## occurred, and glm.fit: fitted probabilities numerically 0 or 1 occurred| Characteristic | N | OR | 95% CI | p-value |
|---|---|---|---|---|
| age | 824 | 0.99 | 0.97, 1.00 | 0.10 |
| gender | 824 | |||
| Ā Ā Ā Ā 1 | ā | ā | ||
| Ā Ā Ā Ā 2 | 1.23 | 0.82, 1.85 | 0.31 | |
| hypertension | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 1.63 | 1.07, 2.48 | 0.022 | |
| a_fib | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 0.53 | 0.35, 0.80 | 0.003 | |
| CHD | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 1.00 | 0.52, 2.12 | 0.99 | |
| diabetes | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 1.20 | 0.80, 1.83 | 0.39 | |
| anemia | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 2.31 | 1.44, 3.86 | <0.001 | |
| depression | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 1.80 | 0.90, 4.13 | 0.13 | |
| Hyperlipemia | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 1.12 | 0.74, 1.72 | 0.59 | |
| CKD | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 2.20 | 1.38, 3.61 | 0.001 | |
| COPD | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā 1 | 3.44 | 1.25, 14.2 | 0.040 | |
| temperature | 824 | 1.82 | 1.30, 2.57 | <0.001 |
| SP O2 | 824 | 1.08 | 1.00, 1.18 | 0.057 |
| hematocrit | 824 | 1.00 | 0.96, 1.04 | 0.91 |
| Leucocyte | 824 | 0.90 | 0.87, 0.93 | <0.001 |
| Platelets | 824 | 1.00 | 1.00, 1.01 | 0.003 |
| INR | 824 | 0.65 | 0.54, 0.78 | <0.001 |
| Creatinine | 824 | 0.87 | 0.76, 1.01 | 0.051 |
| glucose | 824 | 1.00 | 1.00, 1.00 | 0.62 |
| Chloride | 824 | 0.96 | 0.92, 1.0 | 0.026 |
| PH | 824 | 893 | 36.6, 22,910 | <0.001 |
| Bicarbonate | 824 | 1.19 | 1.14, 1.25 | <0.001 |
| PCO2 | 824 | 1.02 | 1.00, 1.04 | 0.053 |
| heart | 824 | 0.98 | 0.97, 1.00 | 0.007 |
| syst | 824 | 1.04 | 1.02, 1.05 | <0.001 |
| dias | 824 | 1.02 | 1.00, 1.04 | 0.034 |
| sat | 824 | 1.08 | 1.00, 1.18 | 0.057 |
| magn | 824 | 0.34 | 0.16, 0.73 | 0.005 |
| probnp | 824 | 1.00 | 1.00, 1.00 | <0.001 |
| resp | 824 | 0.93 | 0.89, 0.98 | 0.005 |
| ureanit | 824 | 0.98 | 0.97, 0.99 | <0.001 |
| urine | 824 | 1.00 | 1.00, 1.00 | <0.001 |
| sodio | 824 | 1.05 | 1.00, 1.09 | 0.053 |
| lactico | 824 | 0.51 | 0.41, 0.62 | <0.001 |
| calcio | 824 | 2.75 | 1.90, 4.04 | <0.001 |
| potasio | 824 | 0.38 | 0.24, 0.59 | <0.001 |
| muerte | 824 | |||
| Ā Ā Ā Ā 0 | ā | ā | ||
| Ā Ā Ā Ā muerte1 | 118,848,047,226,233,318,016,086 | 0.00, Inf | >0.99 | |
| Abbreviations: CI = Confidence Interval, OR = Odds Ratio | ||||
##Variables lineales. Comportamiento
#resp NO ES LINEAL. baja en el Ćŗltimo quintilo. Puede usarse igualmente como lineal p=0.82
DIED0<- DIED0 %>% mutate(resp_q = ntile(resp, 5))
DIED0%>%
dplyr::group_by(resp_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x = resp_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(resp_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = resp_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles resp", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(resp ~ rcs(urine, 3), data = DIED0)
anova(fit)## Wald Statistics Response: resp
##
## Factor Chi-Square d.f. P
## urine 1.46 2 0.483
## Nonlinear 0.05 1 0.828
## TOTAL 1.46 2 0.483#dicotomizar resp
quantile(DIED0$resp,0.80)## 80%
## 23.9describe(DIED0$resp)## DIED0$resp
## n missing distinct Info Mean pMedian Gmd .05
## 824 0 725 1 20.66 20.49 4.451 14.66
## .10 .25 .50 .75 .90 .95
## 15.97 17.86 20.20 23.25 25.76 27.55
##
## lowest : 11.1379 11.8182 12 12.2121 12.3333
## highest: 33.08 33.1154 33.85 34.6934 40.9DIED0 <- DIED0 %>% mutate(respfactor = as.factor(case_when(resp< 23.9 ~ "0",resp>= 23.9 ~"1")))
describe(DIED0$respfactor)## DIED0$respfactor
## n missing distinct
## 824 0 2
##
## Value 0 1
## Frequency 661 163
## Proportion 0.802 0.198levels(DIED0$respfactor)## [1] "0" "1"modresp1<-glm(muerte~respfactor, family = "binomial", data=DIED0)
tab_model(modresp1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 8.18 | 6.45Ā āĀ 10.53 | <0.001 |
| respfactor [1] | 0.43 | 0.28Ā āĀ 0.68 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.017 | ||
#temperatura. El primer quintilo es mƔs bajo que los demƔs. Se puede usar igualmente en forma lineal
DIED0<- DIED0 %>% mutate(temperature_q = ntile(temperature, 5))
DIED0%>%
dplyr::group_by(temperature_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x = temperature_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(temperature_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = temperature_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles temperature", y = "ProporciĆ³n de muertes") +
theme_minimal()
#dicotomizar temp
quantile(DIED0$temperature,0.20)## 20%
## 36.2describe(DIED0$temperature)## DIED0$temperature
## n missing distinct Info Mean pMedian Gmd .05
## 824 0 583 1 36.67 36.66 0.6547 35.80
## .10 .25 .50 .75 .90 .95
## 35.97 36.30 36.66 37.02 37.44 37.67
##
## lowest : 33.25 34.1507 34.3241 34.6056 34.6722
## highest: 38.2698 38.2778 38.3651 38.4136 38.6435DIED0 <- DIED0 %>% mutate(tempfactor = as.factor(case_when(temperature< 36.2 ~ "0",temperature>= 36.2 ~"1")))
describe(DIED0$tempfactor)## DIED0$tempfactor
## n missing distinct
## 824 0 2
##
## Value 0 1
## Frequency 164 660
## Proportion 0.199 0.801levels(DIED0$tempfactor)## [1] "0" "1"modtemp1<-glm(muerte~tempfactor, family = "binomial", data=DIED0)
tab_model(modtemp1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 4.12 | 2.84Ā āĀ 6.17 | <0.001 |
| tempfactor [1] | 1.86 | 1.17Ā āĀ 2.91 | 0.007 |
| Observations | 824 | ||
| R2 Tjur | 0.009 | ||
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(temperature, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## temperature 15.72 2 0.0004
## Nonlinear 3.21 1 0.0733
## TOTAL 15.72 2 0.0004#leucocitos: podrĆa tomarse como lineal
DIED0<- DIED0 %>% mutate(Leucocyte_q = ntile(Leucocyte, 5))
DIED0%>%
dplyr::group_by(Leucocyte_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x = Leucocyte_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(Leucocyte_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = Leucocyte_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles leucocitos", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(temperature, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## temperature 15.72 2 0.0004
## Nonlinear 3.21 1 0.0733
## TOTAL 15.72 2 0.0004#plaquetas
#plaquetas no es lineal
DIED0<- DIED0 %>% mutate(Platelets_q = ntile(Platelets, 5))
DIED0%>%
dplyr::group_by(Platelets_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x = Platelets_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(Platelets_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = Platelets_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles plaquetas", y = "ProporciĆ³n de muertes") +
theme_minimal()
library(rms)
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(Platelets, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## Platelets 27.3 2 <.0001
## Nonlinear 18.2 1 <.0001
## TOTAL 27.3 2 <.0001#INR puede incluirse como lineal
DIED0<- DIED0 %>% mutate(INR_q = ntile(INR, 5))
DIED0%>%
dplyr::group_by(INR_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x = INR_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(INR_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = INR_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles rin", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(INR, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## INR 21.24 2 <.0001
## Nonlinear 0.73 1 0.393
## TOTAL 21.24 2 <.0001#lactico. puede usarse lineal o dicotomizar el ultimo quintil
DIED0<- DIED0 %>% mutate(lactico_q = ntile(lactico, 5))
DIED0%>%
dplyr::group_by(lactico_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =lactico_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(lactico_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = lactico_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles lactico", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(lactico, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## lactico 43.91 2 <.0001
## Nonlinear 1.07 1 0.302
## TOTAL 43.91 2 <.0001#dicotomizar lactico
quantile(DIED0$lactico,0.80)## 80%
## 2.3describe(DIED0$lactico)## DIED0$lactico
## n missing distinct Info Mean pMedian Gmd .05
## 824 0 246 1 1.801 1.675 0.8847 0.825
## .10 .25 .50 .75 .90 .95
## 1.000 1.200 1.600 2.125 2.895 3.402
##
## lowest : 0.5 0.533333 0.6 0.625 0.65
## highest: 6.16667 6.45 6.5 6.6 7.25DIED0 <- DIED0 %>% mutate(lacticofactor = as.factor(case_when(lactico< 2.3 ~ "0",lactico>= 2.3 ~"1")))
describe(DIED0$lacticofactor)## DIED0$lacticofactor
## n missing distinct
## 824 0 2
##
## Value 0 1
## Frequency 657 167
## Proportion 0.797 0.203levels(DIED0$lacticofactor)## [1] "0" "1"modlac1<-glm(muerte~lacticofactor, family = "binomial", data=DIED0)
tab_model(modlac1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 9.77 | 7.58Ā āĀ 12.84 | <0.001 |
| lacticofactor [1] | 0.26 | 0.17Ā āĀ 0.40 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.050 | ||
#PH. No lineal. Dicotomiza en primer quintilo
DIED0<- DIED0 %>% mutate(PH_q = ntile(PH, 5))
DIED0%>%
dplyr::group_by(PH_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =PH_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(PH_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x = PH_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles PH", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(PH, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## PH 23.20 2 <.0001
## Nonlinear 4.64 1 0.0313
## TOTAL 23.20 2 <.0001#dicotomizar PH
quantile(DIED0$PH,0.20)## 20%
## 7.33describe(DIED0$PH)## DIED0$PH
## n missing distinct Info Mean pMedian Gmd .05
## 824 0 331 1 7.379 7.38 0.06996 7.273
## .10 .25 .50 .75 .90 .95
## 7.301 7.340 7.380 7.421 7.459 7.480
##
## lowest : 7.09 7.15833 7.16923 7.182 7.1975
## highest: 7.53 7.534 7.55 7.56 7.58DIED0 <- DIED0 %>% mutate(PHfactor = as.factor(case_when(PH< 7.33 ~ "0",PH>= 7.33 ~"1")))
describe(DIED0$PHfactor)## DIED0$PHfactor
## n missing distinct
## 824 0 2
##
## Value 0 1
## Frequency 153 671
## Proportion 0.186 0.814levels(DIED0$PHofactor)## Warning: Unknown or uninitialised column: `PHofactor`.## NULLmodPH1<-glm(muerte~PHfactor, family = "binomial", data=DIED0)
tab_model(modPH1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 3.37 | 2.34Ā āĀ 4.99 | <0.001 |
| PHfactor [1] | 2.43 | 1.54Ā āĀ 3.78 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.019 | ||
#bicarbonato: nolineal dicotomizar primer quintilo
DIED0<- DIED0 %>% mutate(Bicarbonate_q = ntile(Bicarbonate, 5))
DIED0%>%
dplyr::group_by(Bicarbonate_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =Bicarbonate_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(Bicarbonate_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =Bicarbonate_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles Bicarbonate", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(Bicarbonate, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## Bicarbonate 66.24 2 <.0001
## Nonlinear 9.41 1 0.0022
## TOTAL 66.24 2 <.0001#dicotomizar bicarbonato
quantile(DIED0$Bicarbonate,0.20)## 20%
## 22.7DIED0 <- DIED0 %>% mutate(bicarbonatefactor = as.factor(case_when(Bicarbonate< 22.7 ~ "0",Bicarbonate>= 22.7 ~"1")))
modbic1<-glm(muerte~bicarbonatefactor, family = "binomial", data=DIED0)
tab_model(modbic1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 2.17 | 1.58Ā āĀ 3.03 | <0.001 |
| bicarbonatefactor [1] | 5.04 | 3.29Ā āĀ 7.73 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.076 | ||
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(Bicarbonate, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## Bicarbonate 66.24 2 <.0001
## Nonlinear 9.41 1 0.0022
## TOTAL 66.24 2 <.0001#ureanit. lineal
DIED0<- DIED0 %>% mutate(ureanit_q = ntile(ureanit, 5))
DIED0%>%
dplyr::group_by(ureanit_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =ureanit_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(ureanit_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =ureanit_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles ureanit", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(ureanit, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## ureanit 31.66 2 <.0001
## Nonlinear 1.01 1 0.316
## TOTAL 31.66 2 <.0001#urine. No lineal. parto en el 3er quintilo
DIED0<- DIED0 %>% mutate(urine_q = ntile(urine, 5))
DIED0%>%
dplyr::group_by(urine_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =urine_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(urine_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =urine_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles urine", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(urine, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## urine 51.2 2 <.0001
## Nonlinear 12.8 1 0.0004
## TOTAL 51.2 2 <.0001#dicotomizar urine
quantile(DIED0$urine,0.40)## 40%
## 1372DIED0 <- DIED0 %>% mutate(urinefactor = as.factor(case_when(urine< 1372 ~ "0",urine>= 1372 ~"1")))
modur1<-glm(muerte~urinefactor, family = "binomial", data=DIED0)
tab_model(modur1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 3.52 | 2.73Ā āĀ 4.60 | <0.001 |
| urinefactor [1] | 3.73 | 2.44Ā āĀ 5.79 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.048 | ||
#Creatinine. Lineal.
DIED0<- DIED0 %>% mutate(Creatinine_q = ntile(urine, 5))
DIED0%>%
dplyr::group_by(Creatinine_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =Creatinine_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(Creatinine_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =Creatinine_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles urine", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(Creatinine, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## Creatinine 4.16 2 0.125
## Nonlinear 0.54 1 0.461
## TOTAL 4.16 2 0.125#Cloro.No Lineal. Al dicotomizarlo se vuelve NS
DIED0<- DIED0 %>% mutate(Chloride_q = ntile(urine, 5))
DIED0%>%
dplyr::group_by(Chloride_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =Chloride_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(Chloride_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =Chloride_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles urine", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(Chloride, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## Chloride 11.01 2 0.0041
## Nonlinear 5.45 1 0.0196
## TOTAL 11.01 2 0.0041#dicotomizar cloro
quantile(DIED0$Chloride,0.40)## 40%
## 101DIED0 <- DIED0 %>% mutate(clorofactor = as.factor(case_when(Chloride< 101 ~ "0",Chloride>= 101 ~"1")))
modur1<-glm(muerte~clorofactor, family = "binomial", data=DIED0)
tab_model(modur1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 8.47 | 6.03Ā āĀ 12.30 | <0.001 |
| clorofactor [1] | 0.68 | 0.44Ā āĀ 1.04 | 0.084 |
| Observations | 824 | ||
| R2 Tjur | 0.004 | ||
#syst. Baja primer quintil. no lineal. Se dicotomiza
DIED0<- DIED0 %>% mutate(syst_q = ntile(syst, 5))
DIED0%>%
dplyr::group_by(syst_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =syst_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(syst_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =syst_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles syst", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(syst, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## syst 31.98 2 <.0001
## Nonlinear 4.92 1 0.0266
## TOTAL 31.98 2 <.0001#dicotomizar syst
quantile(DIED0$syst,0.20)## 20%
## 103DIED0 <- DIED0 %>% mutate(systfactor = as.factor(case_when(syst< 103 ~ "0",syst>= 103 ~"1")))
modur1<-glm(muerte~systfactor, family = "binomial", data=DIED0)
tab_model(modur1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 3.25 | 2.30Ā āĀ 4.69 | <0.001 |
| systfactor [1] | 2.65 | 1.71Ā āĀ 4.08 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.025 | ||
#mAGNESIO. Lineal
DIED0<- DIED0 %>% mutate(magn_q = ntile(magn, 5))
DIED0%>%
dplyr::group_by(magn_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =magn_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(magn_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =magn_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles magn", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(magn, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## magn 11.36 2 0.0034
## Nonlinear 3.31 1 0.0688
## TOTAL 11.36 2 0.0034#probnp lineal
DIED0<- DIED0 %>% mutate(probnp_q = ntile(probnp, 5))
DIED0%>%
dplyr::group_by(probnp_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =probnp_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(probnp_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =probnp_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles probnp", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(probnp, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## probnp 13.44 2 0.0012
## Nonlinear 2.73 1 0.0988
## TOTAL 13.44 2 0.0012#heart lineal
DIED0<- DIED0 %>% mutate(heart_q = ntile(heart, 5))
DIED0%>%
dplyr::group_by(heart_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =heart_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(heart_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =heart_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles heart", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(heart, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## heart 8.32 2 0.0156
## Nonlinear 1.98 1 0.1592
## TOTAL 8.32 2 0.0156#SP O2 sat no impresiona lineal
DIED0<- DIED0 %>% mutate(sat_q = ntile(sat, 5))
DIED0%>%
dplyr::group_by(sat_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =sat_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(sat_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =sat_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles heart", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(sat, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## sat 4.84 2 0.0891
## Nonlinear 1.15 1 0.2836
## TOTAL 4.84 2 0.0891#calcio. no lineal. dicotomiza en el 2do quintil
DIED0<- DIED0 %>% mutate(calcio_q = ntile(calcio, 5))
DIED0%>%
dplyr::group_by(calcio_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =calcio_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(calcio_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =calcio_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles calcio", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(calcio, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## calcio 37.77 2 <.0001
## Nonlinear 7.36 1 0.0067
## TOTAL 37.77 2 <.0001#dicotomizar calcio
quantile(DIED0$calcio,0.40)## 40%
## 8.37DIED0 <- DIED0 %>% mutate(calciofactor = as.factor(case_when(calcio< 8.37 ~ "0",calcio>= 8.37 ~"1")))
modur1<-glm(muerte~calciofactor, family = "binomial", data=DIED0)
tab_model(modur1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 3.62 | 2.81Ā āĀ 4.74 | <0.001 |
| calciofactor [1] | 3.49 | 2.29Ā āĀ 5.40 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.043 | ||
#potasio. No lineal. Baja en el ultimo quintilo
DIED0<- DIED0 %>% mutate(potasio_q = ntile(potasio, 5))
DIED0%>%
dplyr::group_by(potasio_q) %>%
dplyr::summarise(tar = mean (as.numeric(muerte)-1)) %>%
ggplot(aes(x =potasio_q, y = tar)) +
geom_point() +
geom_smooth(method = "lm", col = "blue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.5)## `geom_smooth()` using formula = 'y ~ x'
DIED0 %>%
dplyr::group_by(potasio_q) %>%
dplyr::summarise(tar = mean(as.numeric(muerte)-1)) %>%
ggplot(aes(x =potasio_q, y = tar)) +
geom_col(fill = "steelblue") +
geom_text(aes(label = round(tar, 2)), vjust = -0.3) +
labs(x = "Cuartiles potasio", y = "ProporciĆ³n de muertes") +
theme_minimal()
dd <- datadist(DIED0); options(datadist = 'dd')
fit <- lrm(muerte ~ rcs(potasio, 3), data = DIED0)
anova(fit)## Wald Statistics Response: muerte
##
## Factor Chi-Square d.f. P
## potasio 24.25 2 <.0001
## Nonlinear 5.18 1 0.0229
## TOTAL 24.25 2 <.0001#dicotomizar potasio
quantile(DIED0$potasio,0.80)## 80%
## 4.45DIED0 <- DIED0 %>% mutate(potasiofactor = as.factor(case_when(potasio< 4.45 ~ "0",potasio>= 4.45 ~"1")))
modur1<-glm(muerte~potasiofactor, family = "binomial", data=DIED0)
tab_model(modur1)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 8.08 | 6.38Ā āĀ 10.41 | <0.001 |
| potasiofactor [1] | 0.46 | 0.30Ā āĀ 0.72 | 0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.015 | ||
##ingreso de variables antecedentes al modelo
#hipertension, FA, anemia, CKD, COPD
modeloa<-glm(muerte~hypertension, family = "binomial", data=DIED0)
modelob<-glm(muerte~hypertension+ a_fib, family = "binomial", data=DIED0)
modeloc<-glm(muerte~hypertension+a_fib+anemia, family = "binomial", data=DIED0)
modelod<-glm(muerte~hypertension+a_fib+anemia+ CKD, family = "binomial", data=DIED0)
modelof<-glm(muerte~hypertension+ a_fib+ anemia+ CKD+COPD, family = "binomial", data=DIED0)
#comparaciĆ³n de modelos
ta<-tab_model(modeloa,modelob,modeloc,modelod,modelof, transform=NULL,show.se = T, show.dev = T, show.loglik = T, show.intercept = F,show.stat = F, show.reflvl = T, digits = 2, digits.p = 3, show.r2 = F,collapse.ci = F)
knitr::asis_output(ta$knitr)| Ā | muerte | muerte | muerte | muerte | muerte | |||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Predictors | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p |
| a_fib1 | -0.64 | 0.21 | -1.06Ā āĀ -0.23 | 0.002 | -0.60 | 0.21 | -1.02Ā āĀ -0.19 | 0.005 | -0.64 | 0.21 | -1.07Ā āĀ -0.23 | 0.003 | -0.64 | 0.21 | -1.06Ā āĀ -0.22 | 0.003 | ||||
| anemia1 | 0.81 | 0.25 | 0.33Ā āĀ 1.33 | 0.001 | 0.73 | 0.26 | 0.24Ā āĀ 1.25 | 0.004 | 0.74 | 0.26 | 0.26Ā āĀ 1.27 | 0.004 | ||||||||
| CKD1 | 0.68 | 0.25 | 0.20Ā āĀ 1.19 | 0.007 | 0.72 | 0.25 | 0.24Ā āĀ 1.23 | 0.004 | ||||||||||||
| COPD1 | 1.32 | 0.61 | 0.29Ā āĀ 2.75 | 0.029 | ||||||||||||||||
| hypertension1 | 0.49 | 0.22 | 0.06Ā āĀ 0.91 | 0.022 | 0.51 | 0.22 | 0.08Ā āĀ 0.93 | 0.019 | 0.52 | 0.22 | 0.09Ā āĀ 0.95 | 0.017 | 0.41 | 0.22 | -0.03Ā āĀ 0.84 | 0.067 | 0.39 | 0.22 | -0.06Ā āĀ 0.82 | 0.083 |
| Observations | 824 | 824 | 824 | 824 | 824 | |||||||||||||||
| Deviance | 635.056 | 625.651 | 614.108 | 606.316 | 599.479 | |||||||||||||||
| log-Likelihood | -317.528 | -312.826 | -307.054 | -303.158 | -299.739 | |||||||||||||||
#sin hipertensiĆ³n. Se decide elimiar hipertensiĆ³n
modelog<-glm(muerte~ a_fib+ anemia+ CKD+COPD, family = "binomial", data=DIED0)
anova(modelof, modelog, test="LRT")## Analysis of Deviance Table
##
## Model 1: muerte ~ hypertension + a_fib + anemia + CKD + COPD
## Model 2: muerte ~ a_fib + anemia + CKD + COPD
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 818 599
## 2 819 602 -1 -2.94 0.087 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##ingreso de variables continuas al modelo
modeloh<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor, family = "binomial", data=DIED0)
modeloi<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+ PHfactor, family = "binomial", data=DIED0)
modeloj<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+ PHfactor+bicarbonatefactor, family = "binomial", data=DIED0)
modelok<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor, family = "binomial", data=DIED0)
modelol<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor+respfactor, family = "binomial", data=DIED0)
modelom<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor+respfactor, family = "binomial", data=DIED0)
modelon<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor+respfactor+urinefactor, family = "binomial", data=DIED0)
modeloo<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor, family = "binomial", data=DIED0)
modelop<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor, family = "binomial", data=DIED0)
modeloq<-glm(muerte~ a_fib+ anemia+ CKD+COPD+tempfactor+PHfactor+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor, family = "binomial", data=DIED0)
#comparaciĆ³n de modelos
ta<-tab_model(modeloh,modeloi,modeloj,modelok,modelol,modeloo,modelop, modeloq,transform=NULL,show.se = T, show.dev = T, show.loglik = T, show.intercept = F,show.stat = F, show.reflvl = T, digits = 2, digits.p = 3, show.r2 = F,collapse.ci = F)
knitr::asis_output(ta$knitr)| Ā | muerte | muerte | muerte | muerte | muerte | muerte | muerte | muerte | ||||||||||||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Predictors | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p |
| a_fib1 | -0.56 | 0.22 | -0.99Ā āĀ -0.13 | 0.011 | -0.62 | 0.22 | -1.06Ā āĀ -0.19 | 0.005 | -0.56 | 0.23 | -1.01Ā āĀ -0.11 | 0.015 | -0.58 | 0.23 | -1.04Ā āĀ -0.13 | 0.012 | -0.56 | 0.23 | -1.02Ā āĀ -0.10 | 0.017 | -0.40 | 0.24 | -0.88Ā āĀ 0.07 | 0.099 | -0.45 | 0.25 | -0.93Ā āĀ 0.04 | 0.071 | -0.44 | 0.25 | -0.93Ā āĀ 0.05 | 0.077 |
| anemia1 | 0.73 | 0.26 | 0.24Ā āĀ 1.25 | 0.004 | 0.76 | 0.26 | 0.27Ā āĀ 1.29 | 0.004 | 0.75 | 0.27 | 0.24Ā āĀ 1.29 | 0.005 | 0.70 | 0.27 | 0.18Ā āĀ 1.26 | 0.010 | 0.70 | 0.28 | 0.17Ā āĀ 1.25 | 0.012 | 0.73 | 0.28 | 0.20Ā āĀ 1.31 | 0.009 | 0.79 | 0.28 | 0.24Ā āĀ 1.36 | 0.006 | 0.77 | 0.28 | 0.23Ā āĀ 1.34 | 0.007 |
| bicarbonatefactor1 | 1.49 | 0.24 | 1.02Ā āĀ 1.96 | <0.001 | 1.32 | 0.24 | 0.84Ā āĀ 1.80 | <0.001 | 1.34 | 0.25 | 0.85Ā āĀ 1.83 | <0.001 | 1.14 | 0.26 | 0.64Ā āĀ 1.65 | <0.001 | 0.94 | 0.27 | 0.42Ā āĀ 1.46 | <0.001 | 0.86 | 0.27 | 0.32Ā āĀ 1.38 | 0.002 | ||||||||
| calciofactor1 | 0.84 | 0.25 | 0.36Ā āĀ 1.34 | 0.001 | 0.92 | 0.26 | 0.42Ā āĀ 1.43 | <0.001 | ||||||||||||||||||||||||
| CKD1 | 0.86 | 0.25 | 0.38Ā āĀ 1.37 | 0.001 | 0.92 | 0.26 | 0.43Ā āĀ 1.44 | <0.001 | 1.05 | 0.26 | 0.55Ā āĀ 1.59 | <0.001 | 0.92 | 0.27 | 0.40Ā āĀ 1.47 | 0.001 | 0.86 | 0.27 | 0.35Ā āĀ 1.41 | 0.001 | 0.92 | 0.28 | 0.40Ā āĀ 1.48 | 0.001 | 0.82 | 0.28 | 0.29Ā āĀ 1.39 | 0.003 | 0.95 | 0.29 | 0.40Ā āĀ 1.53 | 0.001 |
| COPD1 | 1.35 | 0.61 | 0.32Ā āĀ 2.78 | 0.026 | 1.59 | 0.61 | 0.54Ā āĀ 3.03 | 0.010 | 1.28 | 0.62 | 0.20Ā āĀ 2.73 | 0.041 | 1.27 | 0.64 | 0.17Ā āĀ 2.75 | 0.046 | 1.28 | 0.65 | 0.16Ā āĀ 2.78 | 0.049 | 1.36 | 0.66 | 0.21Ā āĀ 2.88 | 0.040 | 1.37 | 0.66 | 0.23Ā āĀ 2.88 | 0.036 | 1.44 | 0.66 | 0.30Ā āĀ 2.95 | 0.029 |
| lacticofactor1 | -0.96 | 0.24 | -1.43Ā āĀ -0.48 | <0.001 | -0.91 | 0.24 | -1.39Ā āĀ -0.43 | <0.001 | -0.94 | 0.25 | -1.43Ā āĀ -0.44 | <0.001 | -0.91 | 0.26 | -1.41Ā āĀ -0.40 | <0.001 | -0.91 | 0.26 | -1.42Ā āĀ -0.40 | <0.001 | ||||||||||||
| PHfactor1 | 1.10 | 0.24 | 0.61Ā āĀ 1.57 | <0.001 | 0.69 | 0.26 | 0.16Ā āĀ 1.20 | 0.009 | 0.75 | 0.27 | 0.22Ā āĀ 1.28 | 0.005 | 0.72 | 0.27 | 0.18Ā āĀ 1.25 | 0.008 | 0.55 | 0.28 | -0.01Ā āĀ 1.09 | 0.052 | 0.51 | 0.29 | -0.06Ā āĀ 1.06 | 0.075 | 0.28 | 0.30 | -0.32Ā āĀ 0.87 | 0.350 | ||||
| potasiofactor1 | -0.73 | 0.30 | -1.32Ā āĀ -0.13 | 0.016 | ||||||||||||||||||||||||||||
| respfactor1 | -0.72 | 0.26 | -1.21Ā āĀ -0.21 | 0.005 | -0.72 | 0.26 | -1.23Ā āĀ -0.20 | 0.006 | -0.76 | 0.27 | -1.28Ā āĀ -0.23 | 0.005 | -0.76 | 0.27 | -1.28Ā āĀ -0.22 | 0.005 | ||||||||||||||||
| systfactor1 | 0.71 | 0.25 | 0.20Ā āĀ 1.20 | 0.006 | 0.66 | 0.26 | 0.14Ā āĀ 1.16 | 0.011 | 0.64 | 0.26 | 0.12Ā āĀ 1.15 | 0.015 | ||||||||||||||||||||
| tempfactor1 | 0.62 | 0.24 | 0.14Ā āĀ 1.09 | 0.010 | 0.56 | 0.25 | 0.06Ā āĀ 1.04 | 0.024 | 0.47 | 0.26 | -0.05Ā āĀ 0.97 | 0.073 | 0.33 | 0.27 | -0.20Ā āĀ 0.85 | 0.207 | 0.35 | 0.27 | -0.18Ā āĀ 0.87 | 0.185 | 0.19 | 0.28 | -0.37Ā āĀ 0.72 | 0.497 | 0.24 | 0.28 | -0.32Ā āĀ 0.78 | 0.385 | 0.22 | 0.28 | -0.34Ā āĀ 0.77 | 0.428 |
| urinefactor1 | 1.02 | 0.25 | 0.53Ā āĀ 1.52 | <0.001 | 0.95 | 0.26 | 0.46Ā āĀ 1.46 | <0.001 | 0.96 | 0.26 | 0.46Ā āĀ 1.47 | <0.001 | ||||||||||||||||||||
| Observations | 824 | 824 | 824 | 824 | 824 | 824 | 824 | 824 | ||||||||||||||||||||||||
| Deviance | 596.150 | 577.297 | 539.753 | 524.617 | 517.008 | 490.324 | 478.836 | 473.222 | ||||||||||||||||||||||||
| log-Likelihood | -298.075 | -288.649 | -269.876 | -262.309 | -258.504 | -245.162 | -239.418 | -236.611 | ||||||||||||||||||||||||
##Agregado de varibles continuas al modelo sin temperatura ni ph ni fa
#modelo sin temperatura ni ph
modeloq1<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor, family = "binomial", data=DIED0)
modeloq2<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor+ureanit, family = "binomial", data=DIED0)
modeloq3<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor+ureanit+INR, family = "binomial", data=DIED0)
modeloq4<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor++ureanit+INR+Creatinine, family = "binomial", data=DIED0)
modeloq5<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor++ureanit+INR+Creatinine+magn, family = "binomial", data=DIED0)
modeloq6<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor++ureanit+INR+Creatinine+magn+probnp, family = "binomial", data=DIED0)
modeloq7<-glm(muerte~ anemia+CKD+COPD+bicarbonatefactor+lacticofactor+respfactor+urinefactor+systfactor+calciofactor+potasiofactor++ureanit+INR+Creatinine+magn+probnp+heart, family = "binomial", data=DIED0)
#comparaciĆ³n de modelos
ta<-tab_model(modeloq1,modeloq2,modeloq3,modeloq4,modeloq5,modeloq6,modeloq7,transform=NULL,show.se = T, show.dev = T, show.loglik = T, show.intercept = F,show.stat = F, show.reflvl = T, digits = 2, digits.p = 3, show.r2 = F,collapse.ci = F)
knitr::asis_output(ta$knitr)| Ā | muerte | muerte | muerte | muerte | muerte | muerte | muerte | |||||||||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Predictors | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p | Log-Odds | std. Error | CI | p |
| anemia1 | 0.78 | 0.28 | 0.25Ā āĀ 1.35 | 0.005 | 0.89 | 0.29 | 0.34Ā āĀ 1.49 | 0.002 | 0.89 | 0.29 | 0.33Ā āĀ 1.49 | 0.002 | 0.91 | 0.30 | 0.35Ā āĀ 1.52 | 0.002 | 0.92 | 0.30 | 0.35Ā āĀ 1.52 | 0.002 | 0.92 | 0.30 | 0.35Ā āĀ 1.52 | 0.002 | 0.92 | 0.30 | 0.36Ā āĀ 1.52 | 0.002 |
| bicarbonatefactor1 | 0.92 | 0.27 | 0.40Ā āĀ 1.44 | 0.001 | 0.63 | 0.28 | 0.07Ā āĀ 1.17 | 0.025 | 0.64 | 0.28 | 0.08Ā āĀ 1.18 | 0.023 | 0.72 | 0.29 | 0.15Ā āĀ 1.27 | 0.012 | 0.79 | 0.29 | 0.22Ā āĀ 1.36 | 0.006 | 0.80 | 0.29 | 0.22Ā āĀ 1.37 | 0.006 | 0.76 | 0.29 | 0.18Ā āĀ 1.34 | 0.010 |
| calciofactor1 | 0.90 | 0.25 | 0.40Ā āĀ 1.40 | <0.001 | 1.01 | 0.26 | 0.50Ā āĀ 1.53 | <0.001 | 1.00 | 0.26 | 0.49Ā āĀ 1.53 | <0.001 | 1.02 | 0.26 | 0.50Ā āĀ 1.54 | <0.001 | 1.13 | 0.27 | 0.60Ā āĀ 1.67 | <0.001 | 1.13 | 0.27 | 0.60Ā āĀ 1.67 | <0.001 | 1.14 | 0.27 | 0.61Ā āĀ 1.69 | <0.001 |
| CKD1 | 0.92 | 0.29 | 0.37Ā āĀ 1.50 | 0.001 | 1.62 | 0.34 | 0.98Ā āĀ 2.31 | <0.001 | 1.60 | 0.34 | 0.96Ā āĀ 2.29 | <0.001 | 1.42 | 0.35 | 0.75Ā āĀ 2.13 | <0.001 | 1.43 | 0.35 | 0.77Ā āĀ 2.14 | <0.001 | 1.42 | 0.35 | 0.76Ā āĀ 2.14 | <0.001 | 1.39 | 0.35 | 0.72Ā āĀ 2.11 | <0.001 |
| COPD1 | 1.40 | 0.65 | 0.27Ā āĀ 2.89 | 0.032 | 1.49 | 0.70 | 0.29Ā āĀ 3.09 | 0.033 | 1.45 | 0.70 | 0.25Ā āĀ 3.05 | 0.038 | 1.46 | 0.69 | 0.27Ā āĀ 3.04 | 0.034 | 1.55 | 0.70 | 0.35Ā āĀ 3.16 | 0.027 | 1.55 | 0.70 | 0.35Ā āĀ 3.15 | 0.026 | 1.54 | 0.69 | 0.35Ā āĀ 3.14 | 0.026 |
| Creatinine | 0.31 | 0.19 | -0.03Ā āĀ 0.70 | 0.104 | 0.29 | 0.19 | -0.05Ā āĀ 0.68 | 0.120 | 0.29 | 0.19 | -0.05Ā āĀ 0.68 | 0.121 | 0.29 | 0.19 | -0.05Ā āĀ 0.69 | 0.118 | ||||||||||||
| heart | -0.01 | 0.01 | -0.03Ā āĀ 0.01 | 0.193 | ||||||||||||||||||||||||
| INR | -0.15 | 0.12 | -0.40Ā āĀ 0.09 | 0.224 | -0.16 | 0.12 | -0.40Ā āĀ 0.08 | 0.202 | -0.15 | 0.12 | -0.40Ā āĀ 0.09 | 0.224 | -0.15 | 0.12 | -0.40Ā āĀ 0.09 | 0.224 | -0.18 | 0.13 | -0.43Ā āĀ 0.07 | 0.163 | ||||||||
| lacticofactor1 | -0.93 | 0.26 | -1.43Ā āĀ -0.42 | <0.001 | -0.93 | 0.26 | -1.44Ā āĀ -0.41 | <0.001 | -0.86 | 0.27 | -1.39Ā āĀ -0.32 | 0.002 | -0.86 | 0.27 | -1.39Ā āĀ -0.32 | 0.002 | -0.85 | 0.27 | -1.38Ā āĀ -0.30 | 0.002 | -0.85 | 0.27 | -1.38Ā āĀ -0.31 | 0.002 | -0.83 | 0.28 | -1.36Ā āĀ -0.28 | 0.003 |
| magn | -1.06 | 0.52 | -2.06Ā āĀ -0.03 | 0.041 | -1.05 | 0.52 | -2.06Ā āĀ -0.02 | 0.042 | -1.05 | 0.52 | -2.05Ā āĀ -0.02 | 0.042 | ||||||||||||||||
| potasiofactor1 | -0.82 | 0.29 | -1.38Ā āĀ -0.26 | 0.004 | -0.49 | 0.30 | -1.07Ā āĀ 0.12 | 0.109 | -0.42 | 0.31 | -1.02Ā āĀ 0.19 | 0.171 | -0.44 | 0.31 | -1.04Ā āĀ 0.18 | 0.158 | -0.44 | 0.31 | -1.05Ā āĀ 0.18 | 0.153 | -0.45 | 0.31 | -1.05Ā āĀ 0.17 | 0.152 | -0.43 | 0.31 | -1.04Ā āĀ 0.19 | 0.162 |
| probnp | 0.00 | 0.00 | -0.00Ā āĀ 0.00 | 0.876 | 0.00 | 0.00 | -0.00Ā āĀ 0.00 | 0.752 | ||||||||||||||||||||
| respfactor1 | -0.76 | 0.27 | -1.28Ā āĀ -0.23 | 0.004 | -0.64 | 0.28 | -1.19Ā āĀ -0.09 | 0.022 | -0.62 | 0.28 | -1.17Ā āĀ -0.06 | 0.027 | -0.61 | 0.28 | -1.15Ā āĀ -0.05 | 0.031 | -0.61 | 0.28 | -1.16Ā āĀ -0.05 | 0.030 | -0.61 | 0.28 | -1.16Ā āĀ -0.05 | 0.031 | -0.50 | 0.29 | -1.07Ā āĀ 0.09 | 0.090 |
| systfactor1 | 0.66 | 0.26 | 0.14Ā āĀ 1.17 | 0.011 | 0.69 | 0.27 | 0.15Ā āĀ 1.21 | 0.011 | 0.67 | 0.27 | 0.13Ā āĀ 1.19 | 0.014 | 0.64 | 0.27 | 0.10Ā āĀ 1.17 | 0.018 | 0.72 | 0.28 | 0.17Ā āĀ 1.26 | 0.009 | 0.72 | 0.28 | 0.17Ā āĀ 1.26 | 0.009 | 0.65 | 0.28 | 0.09Ā āĀ 1.20 | 0.021 |
| ureanit | -0.03 | 0.01 | -0.04Ā āĀ -0.02 | <0.001 | -0.03 | 0.01 | -0.04Ā āĀ -0.02 | <0.001 | -0.04 | 0.01 | -0.05Ā āĀ -0.02 | <0.001 | -0.03 | 0.01 | -0.05Ā āĀ -0.02 | <0.001 | -0.03 | 0.01 | -0.05Ā āĀ -0.02 | <0.001 | -0.04 | 0.01 | -0.05Ā āĀ -0.02 | <0.001 | ||||
| urinefactor1 | 1.08 | 0.25 | 0.59Ā āĀ 1.57 | <0.001 | 0.92 | 0.26 | 0.42Ā āĀ 1.43 | <0.001 | 0.91 | 0.26 | 0.41Ā āĀ 1.43 | <0.001 | 0.96 | 0.26 | 0.45Ā āĀ 1.48 | <0.001 | 0.90 | 0.26 | 0.39Ā āĀ 1.42 | 0.001 | 0.90 | 0.26 | 0.39Ā āĀ 1.43 | 0.001 | 0.92 | 0.26 | 0.41Ā āĀ 1.45 | <0.001 |
| Observations | 824 | 824 | 824 | 824 | 824 | 824 | 824 | |||||||||||||||||||||
| Deviance | 478.073 | 449.584 | 448.082 | 445.040 | 441.010 | 440.985 | 439.292 | |||||||||||||||||||||
| log-Likelihood | -239.037 | -224.792 | -224.041 | -222.520 | -220.505 | -220.492 | -219.646 | |||||||||||||||||||||
##Efectos principales
modEFPRINC<- glm(muerte ~ a_fib+COPD + CKD + lacticofactor + bicarbonatefactor+ ureanit+urinefactor+ calcio+ respfactor +anemia, family = "binomial", data = DIED0)
summary(modEFPRINC)##
## Call:
## glm(formula = muerte ~ a_fib + COPD + CKD + lacticofactor + bicarbonatefactor +
## ureanit + urinefactor + calcio + respfactor + anemia, family = "binomial",
## data = DIED0)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.8647 1.8703 -2.60 0.00930 **
## a_fib1 -0.3640 0.2490 -1.46 0.14377
## COPD1 1.4549 0.7120 2.04 0.04100 *
## CKD1 1.5426 0.3308 4.66 0.000003108 ***
## lacticofactor1 -0.9415 0.2603 -3.62 0.00030 ***
## bicarbonatefactor1 0.7288 0.2733 2.67 0.00766 **
## ureanit -0.0325 0.0058 -5.60 0.000000022 ***
## urinefactor1 0.9257 0.2571 3.60 0.00032 ***
## calcio 0.8262 0.2325 3.55 0.00038 ***
## respfactor1 -0.6461 0.2760 -2.34 0.01923 *
## anemia1 0.9072 0.2936 3.09 0.00200 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 640.10 on 823 degrees of freedom
## Residual deviance: 459.33 on 813 degrees of freedom
## AIC: 481.3
##
## Number of Fisher Scoring iterations: 6tab_model(modEFPRINC)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 0.01 | 0.00Ā āĀ 0.29 | 0.009 |
| a fib [1] | 0.69 | 0.42Ā āĀ 1.13 | 0.144 |
| COPD [1] | 4.28 | 1.27Ā āĀ 22.23 | 0.041 |
| CKD [1] | 4.68 | 2.50Ā āĀ 9.18 | <0.001 |
| lacticofactor [1] | 0.39 | 0.23Ā āĀ 0.65 | <0.001 |
| bicarbonatefactor [1] | 2.07 | 1.21Ā āĀ 3.53 | 0.008 |
| ureanit | 0.97 | 0.96Ā āĀ 0.98 | <0.001 |
| urinefactor [1] | 2.52 | 1.53Ā āĀ 4.20 | <0.001 |
| calcio | 2.28 | 1.46Ā āĀ 3.64 | <0.001 |
| respfactor [1] | 0.52 | 0.31Ā āĀ 0.91 | 0.019 |
| anemia [1] | 2.48 | 1.42Ā āĀ 4.50 | 0.002 |
| Observations | 824 | ||
| R2 Tjur | 0.279 | ||
anova(modEFPRINC)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: muerte
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 823 640
## a_fib 1 9.1 822 631 0.00254 **
## COPD 1 5.5 821 626 0.01948 *
## CKD 1 14.0 820 612 0.00018 ***
## lacticofactor 1 29.2 819 582 0.00000006411 ***
## bicarbonatefactor 1 41.5 818 541 0.00000000012 ***
## ureanit 1 37.3 817 503 0.00000000103 ***
## urinefactor 1 16.8 816 487 0.00004123349 ***
## calcio 1 11.3 815 475 0.00079 ***
## respfactor 1 5.5 814 470 0.01856 *
## anemia 1 10.5 813 459 0.00119 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1##CalibraciĆ³n modelo efectos principales
modEFPRINC<- glm(muerte ~ a_fib+COPD + CKD + lacticofactor + bicarbonatefactor+ ureanit+urinefactor+ calcio+ respfactor +anemia, family = "binomial", data = DIED0)
#HOSMER & LEMERSHOW. LA P ME DA . NO RECHAZA LA HIPĆTESIS NULA.ES COMO QUIERO QUE ME DE
pred_modelo <- as.data.frame(predict(modEFPRINC, DIED0, type = "response")) %>% rename(pred = `predict(modEFPRINC, DIED0, type = "response")`)
pred_modelo$verdad <- DIED0$muerte
H_L <- HLtest(modEFPRINC,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ a_fib + COPD + CKD + lacticofactor + bicarbonatefactor +
## ureanit + urinefactor + calcio + respfactor + anemia, family = "binomial",
## data = DIED0)
## ChiSquare df P_value
## 13.8 8 0.0881HL_DF[,1:4]## cut total obs exp
## 1 [0.0269,0.647] 83 47 48.382
## 2 (0.647,0.807] 82 23 21.873
## 3 (0.807,0.88] 82 13 12.816
## 4 (0.88,0.917] 83 5 8.279
## 5 (0.917,0.942] 82 8 5.731
## 6 (0.942,0.958] 82 5 4.031
## 7 (0.958,0.971] 83 1 2.869
## 8 (0.971,0.98] 82 6 2.065
## 9 (0.98,0.988] 82 0 1.320
## 10 (0.988,0.999] 83 0 0.634HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.0269,0.647] | 83 | 47 | 48.382 | -0.199 | 1 |
| (0.647,0.807] | 82 | 23 | 21.873 | 0.241 | 2 |
| (0.807,0.88] | 82 | 13 | 12.816 | 0.051 | 3 |
| (0.88,0.917] | 83 | 5 | 8.279 | -1.140 | 4 |
| (0.917,0.942] | 82 | 8 | 5.731 | 0.948 | 5 |
| (0.942,0.958] | 82 | 5 | 4.031 | 0.482 | 6 |
| (0.958,0.971] | 83 | 1 | 2.869 | -1.103 | 7 |
| (0.971,0.98] | 82 | 6 | 2.065 | 2.739 | 8 |
| (0.98,0.988] | 82 | 0 | 1.320 | -1.149 | 9 |
| (0.988,0.999] | 83 | 0 | 0.634 | -0.796 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 108 716str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
#curva roc. Area 0.83 sin temperatura.vamos con este modelo prueba i3
roc<-roc(pred_modelo$verdad, pred_modelo$pred)## Setting levels: control = 0, case = 1## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
##MODELO DE EFECTOS PRINCIPALES SIN A_FIB. SON MENOS VARIABLES
modEFPRINC1<- glm(muerte ~ COPD + CKD + lacticofactor + bicarbonatefactor+ ureanit+urinefactor+ calcio+ respfactor +anemia, family = "binomial", data = DIED0)
summary(modEFPRINC1)##
## Call:
## glm(formula = muerte ~ COPD + CKD + lacticofactor + bicarbonatefactor +
## ureanit + urinefactor + calcio + respfactor + anemia, family = "binomial",
## data = DIED0)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.87256 1.86854 -2.61 0.00912 **
## COPD1 1.44229 0.70565 2.04 0.04096 *
## CKD1 1.55474 0.33179 4.69 0.000002788 ***
## lacticofactor1 -0.95172 0.26026 -3.66 0.00026 ***
## bicarbonatefactor1 0.71961 0.27418 2.62 0.00867 **
## ureanit -0.03320 0.00581 -5.71 0.000000011 ***
## urinefactor1 0.97211 0.25497 3.81 0.00014 ***
## calcio 0.80482 0.23171 3.47 0.00051 ***
## respfactor1 -0.64954 0.27565 -2.36 0.01845 *
## anemia1 0.94086 0.29229 3.22 0.00129 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 640.10 on 823 degrees of freedom
## Residual deviance: 461.48 on 814 degrees of freedom
## AIC: 481.5
##
## Number of Fisher Scoring iterations: 6anova(modEFPRINC1)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: muerte
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 823 640
## COPD 1 6.1 822 634 0.01386 *
## CKD 1 12.6 821 621 0.00038 ***
## lacticofactor 1 30.3 820 591 0.000000036271 ***
## bicarbonatefactor 1 42.9 819 548 0.000000000058 ***
## ureanit 1 40.5 818 508 0.000000000192 ***
## urinefactor 1 18.8 817 489 0.000014227547 ***
## calcio 1 10.3 816 479 0.00134 **
## respfactor 1 5.6 815 473 0.01826 *
## anemia 1 11.5 814 461 0.00071 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#HOSMER & LEMERSHOW. LA P ME DA . NO RECHAZA LA HIPĆTESIS NULA.ES COMO QUIERO QUE ME DE
pred_modelo <- as.data.frame(predict(modEFPRINC1, DIED0, type = "response")) %>% rename(pred = `predict(modEFPRINC1, DIED0, type = "response")`)
pred_modelo$verdad <- DIED0$muerte
H_L <- HLtest(modEFPRINC1,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ COPD + CKD + lacticofactor + bicarbonatefactor +
## ureanit + urinefactor + calcio + respfactor + anemia, family = "binomial",
## data = DIED0)
## ChiSquare df P_value
## 4.75 8 0.784HL_DF[,1:4]## cut total obs exp
## 1 [0.0297,0.662] 83 45 48.187
## 2 (0.662,0.808] 82 22 21.732
## 3 (0.808,0.879] 82 16 12.646
## 4 (0.879,0.916] 83 7 8.264
## 5 (0.916,0.941] 82 7 5.912
## 6 (0.941,0.959] 82 3 4.085
## 7 (0.959,0.969] 83 3 2.998
## 8 (0.969,0.979] 82 4 2.129
## 9 (0.979,0.988] 82 1 1.403
## 10 (0.988,0.999] 83 0 0.644HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.0297,0.662] | 83 | 45 | 48.187 | -0.459 | 1 |
| (0.662,0.808] | 82 | 22 | 21.732 | 0.057 | 2 |
| (0.808,0.879] | 82 | 16 | 12.646 | 0.943 | 3 |
| (0.879,0.916] | 83 | 7 | 8.264 | -0.440 | 4 |
| (0.916,0.941] | 82 | 7 | 5.912 | 0.448 | 5 |
| (0.941,0.959] | 82 | 3 | 4.085 | -0.537 | 6 |
| (0.959,0.969] | 83 | 3 | 2.998 | 0.001 | 7 |
| (0.969,0.979] | 82 | 4 | 2.129 | 1.282 | 8 |
| (0.979,0.988] | 82 | 1 | 1.403 | -0.340 | 9 |
| (0.988,0.999] | 83 | 0 | 0.644 | -0.803 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 108 716str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
#curva roc. Area 0.83 sin temperatura.vamos con este modelo prueba i3
roc<-roc(pred_modelo$verdad, pred_modelo$pred)## Setting levels: control = 0, case = 1## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
##AGREGA DIABETES : MODELO DE EFECTOS PRINCIPALES. DIABETES AGREGA
VALIDEZ EXTERNA
modEFPRINC1diab<- glm(muerte ~ COPD + diabetes+CKD + lacticofactor + bicarbonatefactor+ ureanit+urinefactor+ calcio+ respfactor +anemia, family = "binomial", data = DIED0)
summary(modEFPRINC1diab)##
## Call:
## glm(formula = muerte ~ COPD + diabetes + CKD + lacticofactor +
## bicarbonatefactor + ureanit + urinefactor + calcio + respfactor +
## anemia, family = "binomial", data = DIED0)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.90681 1.87882 -2.61 0.00901 **
## COPD1 1.50434 0.71412 2.11 0.03515 *
## diabetes1 0.26903 0.26733 1.01 0.31426
## CKD1 1.52540 0.33317 4.58 0.0000046848 ***
## lacticofactor1 -0.95906 0.26057 -3.68 0.00023 ***
## bicarbonatefactor1 0.70080 0.27451 2.55 0.01068 *
## ureanit -0.03471 0.00603 -5.76 0.0000000086 ***
## urinefactor1 0.94203 0.25695 3.67 0.00025 ***
## calcio 0.80753 0.23311 3.46 0.00053 ***
## respfactor1 -0.64876 0.27608 -2.35 0.01878 *
## anemia1 0.94009 0.29296 3.21 0.00133 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 640.10 on 823 degrees of freedom
## Residual deviance: 460.46 on 813 degrees of freedom
## AIC: 482.5
##
## Number of Fisher Scoring iterations: 6anova(modEFPRINC1diab)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: muerte
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 823 640
## COPD 1 6.1 822 634 0.01386 *
## diabetes 1 1.2 821 633 0.27214
## CKD 1 11.6 820 621 0.00065 ***
## lacticofactor 1 30.5 819 591 0.000000033002 ***
## bicarbonatefactor 1 42.6 818 548 0.000000000066 ***
## ureanit 1 43.3 817 505 0.000000000047 ***
## urinefactor 1 17.4 816 487 0.000030689495 ***
## calcio 1 10.0 815 477 0.00154 **
## respfactor 1 5.5 814 472 0.01914 *
## anemia 1 11.4 813 460 0.00074 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova(modEFPRINC1,modEFPRINC1diab, test="LRT")## Analysis of Deviance Table
##
## Model 1: muerte ~ COPD + CKD + lacticofactor + bicarbonatefactor + ureanit +
## urinefactor + calcio + respfactor + anemia
## Model 2: muerte ~ COPD + diabetes + CKD + lacticofactor + bicarbonatefactor +
## ureanit + urinefactor + calcio + respfactor + anemia
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 814 461
## 2 813 460 1 1.02 0.31#HOSMER & LEMERSHOW. LA P ME DA . NO RECHAZA LA HIPĆTESIS NULA.ES COMO QUIERO QUE ME DE
pred_modelo <- as.data.frame(predict(modEFPRINC1diab, DIED0, type = "response")) %>% rename(pred = `predict(modEFPRINC1diab, DIED0, type = "response")`)
pred_modelo$verdad <- DIED0$muerte
H_L <- HLtest(modEFPRINC1diab,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ COPD + diabetes + CKD + lacticofactor +
## bicarbonatefactor + ureanit + urinefactor + calcio + respfactor +
## anemia, family = "binomial", data = DIED0)
## ChiSquare df P_value
## 5.31 8 0.724HL_DF[,1:4]## cut total obs exp
## 1 [0.0245,0.654] 83 47 48.154
## 2 (0.654,0.805] 82 22 21.981
## 3 (0.805,0.881] 82 16 12.587
## 4 (0.881,0.915] 83 4 8.330
## 5 (0.915,0.94] 82 7 5.869
## 6 (0.94,0.958] 82 5 4.046
## 7 (0.958,0.97] 83 3 2.977
## 8 (0.97,0.979] 82 3 2.065
## 9 (0.979,0.988] 82 1 1.350
## 10 (0.988,0.999] 83 0 0.641HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.0245,0.654] | 83 | 47 | 48.154 | -0.166 | 1 |
| (0.654,0.805] | 82 | 22 | 21.981 | 0.004 | 2 |
| (0.805,0.881] | 82 | 16 | 12.587 | 0.962 | 3 |
| (0.881,0.915] | 83 | 4 | 8.330 | -1.500 | 4 |
| (0.915,0.94] | 82 | 7 | 5.869 | 0.467 | 5 |
| (0.94,0.958] | 82 | 5 | 4.046 | 0.475 | 6 |
| (0.958,0.97] | 83 | 3 | 2.977 | 0.014 | 7 |
| (0.97,0.979] | 82 | 3 | 2.065 | 0.651 | 8 |
| (0.979,0.988] | 82 | 1 | 1.350 | -0.301 | 9 |
| (0.988,0.999] | 83 | 0 | 0.641 | -0.801 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 108 716str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
#curva roc. Area 0.83 sin temperatura.vamos con este modelo prueba i3
roc<-roc(pred_modelo$verdad, pred_modelo$pred)## Setting levels: control = 0, case = 1## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
##Modelo efectos principales + diabetes: interacciĆ³n con ckd
#EL TERMINO DE INTERACCIĆN NO ES SIGNIFICATIVO. DEJO DIABETES PARA AUMENTAR VALIDEZ EXTERNA PERO SACO EL TĆRMINO DE INTERACCIĆN
modEFPRINC1diabCKD<- glm(muerte ~ COPD + diabetes+CKD + lacticofactor + bicarbonatefactor+ ureanit+urinefactor+ calcio+ respfactor +anemia + diabetes:CKD, family = "binomial", data = DIED0)
summary(modEFPRINC1diabCKD)##
## Call:
## glm(formula = muerte ~ COPD + diabetes + CKD + lacticofactor +
## bicarbonatefactor + ureanit + urinefactor + calcio + respfactor +
## anemia + diabetes:CKD, family = "binomial", data = DIED0)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.91274 1.87952 -2.61 0.00895 **
## COPD1 1.49550 0.71266 2.10 0.03586 *
## diabetes1 0.15439 0.30854 0.50 0.61680
## CKD1 1.32565 0.42319 3.13 0.00173 **
## lacticofactor1 -0.96386 0.26110 -3.69 0.00022 ***
## bicarbonatefactor1 0.70983 0.27491 2.58 0.00982 **
## ureanit -0.03487 0.00604 -5.78 0.0000000077 ***
## urinefactor1 0.94276 0.25725 3.66 0.00025 ***
## calcio 0.81326 0.23328 3.49 0.00049 ***
## respfactor1 -0.65393 0.27660 -2.36 0.01807 *
## anemia1 0.93287 0.29286 3.19 0.00145 **
## diabetes1:CKD1 0.41363 0.56833 0.73 0.46674
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 640.10 on 823 degrees of freedom
## Residual deviance: 459.93 on 812 degrees of freedom
## AIC: 483.9
##
## Number of Fisher Scoring iterations: 6anova(modEFPRINC1diabCKD)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: muerte
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 823 640
## COPD 1 6.1 822 634 0.01386 *
## diabetes 1 1.2 821 633 0.27214
## CKD 1 11.6 820 621 0.00065 ***
## lacticofactor 1 30.5 819 591 0.000000033002 ***
## bicarbonatefactor 1 42.6 818 548 0.000000000066 ***
## ureanit 1 43.3 817 505 0.000000000047 ***
## urinefactor 1 17.4 816 487 0.000030689495 ***
## calcio 1 10.0 815 477 0.00154 **
## respfactor 1 5.5 814 472 0.01914 *
## anemia 1 11.4 813 460 0.00074 ***
## diabetes:CKD 1 0.5 812 460 0.46760
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova(modEFPRINC1,modEFPRINC1diabCKD, test="LRT")## Analysis of Deviance Table
##
## Model 1: muerte ~ COPD + CKD + lacticofactor + bicarbonatefactor + ureanit +
## urinefactor + calcio + respfactor + anemia
## Model 2: muerte ~ COPD + diabetes + CKD + lacticofactor + bicarbonatefactor +
## ureanit + urinefactor + calcio + respfactor + anemia + diabetes:CKD
## Resid. Df Resid. Dev Df Deviance Pr(>Chi)
## 1 814 461
## 2 812 460 2 1.55 0.46#HOSMER & LEMERSHOW. LA P ME DA . NO RECHAZA LA HIPĆTESIS NULA.ES COMO QUIERO QUE ME DE
pred_modelo <- as.data.frame(predict(modEFPRINC1diabCKD, DIED0, type = "response")) %>% rename(pred = `predict(modEFPRINC1diabCKD, DIED0, type = "response")`)
pred_modelo$verdad <- DIED0$muerte
H_L <- HLtest(modEFPRINC1diabCKD,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ COPD + diabetes + CKD + lacticofactor +
## bicarbonatefactor + ureanit + urinefactor + calcio + respfactor +
## anemia + diabetes:CKD, family = "binomial", data = DIED0)
## ChiSquare df P_value
## 5.11 8 0.746HL_DF[,1:4]## cut total obs exp
## 1 [0.0248,0.661] 83 47 48.356
## 2 (0.661,0.803] 82 21 21.726
## 3 (0.803,0.878] 82 16 12.593
## 4 (0.878,0.915] 83 6 8.423
## 5 (0.915,0.942] 82 6 5.841
## 6 (0.942,0.959] 82 5 4.074
## 7 (0.959,0.97] 83 2 2.921
## 8 (0.97,0.98] 82 4 2.069
## 9 (0.98,0.988] 82 1 1.356
## 10 (0.988,0.999] 83 0 0.641HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.0248,0.661] | 83 | 47 | 48.356 | -0.195 | 1 |
| (0.661,0.803] | 82 | 21 | 21.726 | -0.156 | 2 |
| (0.803,0.878] | 82 | 16 | 12.593 | 0.960 | 3 |
| (0.878,0.915] | 83 | 6 | 8.423 | -0.835 | 4 |
| (0.915,0.942] | 82 | 6 | 5.841 | 0.066 | 5 |
| (0.942,0.959] | 82 | 5 | 4.074 | 0.459 | 6 |
| (0.959,0.97] | 83 | 2 | 2.921 | -0.539 | 7 |
| (0.97,0.98] | 82 | 4 | 2.069 | 1.343 | 8 |
| (0.98,0.988] | 82 | 1 | 1.356 | -0.305 | 9 |
| (0.988,0.999] | 83 | 0 | 0.641 | -0.801 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 108 716str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
#curva roc. Area 0.83 sin temperatura.vamos con este modelo prueba i3
roc<-roc(pred_modelo$verdad, pred_modelo$pred)## Setting levels: control = 0, case = 1## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
##CREO QUE EL MODELO DE EFECTOS PRINCIPALES ES MODEFECPRIN
modEFPRINC1diab<- glm(muerte ~ COPD + diabetes+CKD + lacticofactor + bicarbonatefactor+ ureanit+urinefactor+ calcio+ respfactor +anemia, family = "binomial", data = DIED0)
summary(modEFPRINC1diab)##
## Call:
## glm(formula = muerte ~ COPD + diabetes + CKD + lacticofactor +
## bicarbonatefactor + ureanit + urinefactor + calcio + respfactor +
## anemia, family = "binomial", data = DIED0)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.90681 1.87882 -2.61 0.00901 **
## COPD1 1.50434 0.71412 2.11 0.03515 *
## diabetes1 0.26903 0.26733 1.01 0.31426
## CKD1 1.52540 0.33317 4.58 0.0000046848 ***
## lacticofactor1 -0.95906 0.26057 -3.68 0.00023 ***
## bicarbonatefactor1 0.70080 0.27451 2.55 0.01068 *
## ureanit -0.03471 0.00603 -5.76 0.0000000086 ***
## urinefactor1 0.94203 0.25695 3.67 0.00025 ***
## calcio 0.80753 0.23311 3.46 0.00053 ***
## respfactor1 -0.64876 0.27608 -2.35 0.01878 *
## anemia1 0.94009 0.29296 3.21 0.00133 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 640.10 on 823 degrees of freedom
## Residual deviance: 460.46 on 813 degrees of freedom
## AIC: 482.5
##
## Number of Fisher Scoring iterations: 6anova(modEFPRINC1diab)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: muerte
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 823 640
## COPD 1 6.1 822 634 0.01386 *
## diabetes 1 1.2 821 633 0.27214
## CKD 1 11.6 820 621 0.00065 ***
## lacticofactor 1 30.5 819 591 0.000000033002 ***
## bicarbonatefactor 1 42.6 818 548 0.000000000066 ***
## ureanit 1 43.3 817 505 0.000000000047 ***
## urinefactor 1 17.4 816 487 0.000030689495 ***
## calcio 1 10.0 815 477 0.00154 **
## respfactor 1 5.5 814 472 0.01914 *
## anemia 1 11.4 813 460 0.00074 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tab_model(modEFPRINC1diab)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 0.01 | 0.00Ā āĀ 0.28 | 0.009 |
| COPD [1] | 4.50 | 1.33Ā āĀ 23.54 | 0.035 |
| diabetes [1] | 1.31 | 0.78Ā āĀ 2.23 | 0.314 |
| CKD [1] | 4.60 | 2.45Ā āĀ 9.07 | <0.001 |
| lacticofactor [1] | 0.38 | 0.23Ā āĀ 0.64 | <0.001 |
| bicarbonatefactor [1] | 2.02 | 1.17Ā āĀ 3.44 | 0.011 |
| ureanit | 0.97 | 0.95Ā āĀ 0.98 | <0.001 |
| urinefactor [1] | 2.57 | 1.56Ā āĀ 4.27 | <0.001 |
| calcio | 2.24 | 1.43Ā āĀ 3.58 | 0.001 |
| respfactor [1] | 0.52 | 0.31Ā āĀ 0.91 | 0.019 |
| anemia [1] | 2.56 | 1.47Ā āĀ 4.65 | 0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.277 | ||
MODEFECPRIN=modEFPRINC1diab##Base test
library(readr)
test_TH <- read_csv("C:/Users/Administrador/Downloads/test_TH.csv")## Rows: 352 Columns: 36
## āā Column specification āāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāā
## Delimiter: ","
## dbl (36): Died, age, gender, hypertension, a_fib, CHD, diabetes, anemia, dep...
##
## ā¹ Use `spec()` to retrieve the full column specification for this data.
## ā¹ Specify the column types or set `show_col_types = FALSE` to quiet this message.View(test_TH)
test_TH1=test_TH
#preparaciĆ³n de la base test
test_TH1$heart=test_TH$`heart rate`
test_TH1$syst=test_TH$`Systolic blood pressure`
test_TH1$dias=test_TH$`Diastolic blood pressure`
test_TH1$magn=test_TH$`Magnesium ion`
test_TH1$probnp=test_TH$`NT-proBNP`
test_TH1$resp=test_TH$`Respiratory rate`
test_TH1$ureanit=test_TH$`Urea nitrogen`
test_TH1$urine=test_TH$`Urine output`
test_TH1$sodio=test_TH$`Blood sodium`
test_TH1$lactico=test_TH$`Lactic acid`
test_TH1$calcio=test_TH$`Blood calcium`
test_TH1$potasio=test_TH$`Blood potassium`
summary(test_TH1$lactico)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.50 1.23 1.60 1.85 2.15 8.33quantile(test_TH1$lactico,0.60)## 60%
## 1.79test_TH1 <- test_TH1 %>% mutate_at(vars(Died, gender, hypertension, a_fib, CHD, diabetes, anemia, depression, Hyperlipemia, CKD, COPD ),~as.factor(.))
#dicotomizar resp
quantile(test_TH1$resp,0.80)## 80%
## 24.4describe(test_TH1$resp)## test_TH1$resp
## n missing distinct Info Mean pMedian Gmd .05
## 352 0 337 1 21.05 20.88 4.496 15.00
## .10 .25 .50 .75 .90 .95
## 16.07 18.05 20.80 23.68 26.00 28.46
##
## lowest : 12.5 12.6 12.697 13 13.5758
## highest: 30.8667 32 33.1034 35.2667 35.75test_TH1 <- test_TH1 %>% mutate(respfactor = as.factor(case_when(resp< 23.9 ~ "0",resp>= 23.9 ~"1")))
#dicotomizar urine
quantile(test_TH1$urine,0.40)## 40%
## 1301test_TH1 <- test_TH1 %>% mutate(urinefactor = as.factor(case_when(urine< 1372 ~ "0",urine>= 1372 ~"1")))
#dicotomizar bicarbonato
quantile(test_TH1$Bicarbonate,0.20)## 20%
## 22.8test_TH1 <- test_TH1 %>% mutate(bicarbonatefactor = as.factor(case_when(Bicarbonate< 22.7 ~ "0",Bicarbonate>= 22.7 ~"1")))
#dicotomizar lactico
quantile(test_TH1$lactico,0.80)## 80%
## 2.3describe(test_TH1$lactico)## test_TH1$lactico
## n missing distinct Info Mean pMedian Gmd .05
## 352 0 154 1 1.85 1.692 0.95 0.8386
## .10 .25 .50 .75 .90 .95
## 1.0000 1.2312 1.6000 2.1479 2.8129 3.6858
##
## lowest : 0.5 0.6 0.7 0.74 0.79
## highest: 6.12222 6.725 6.77778 7.9 8.33333test_TH1 <- test_TH1 %>% mutate(lacticofactor = as.factor(case_when(lactico< 2.3 ~ "0",lactico>= 2.3 ~"1")))
#crear muerte en TH1
test_TH1 <- test_TH1 %>% mutate(muerte = as.factor(case_when(
Died == "1" ~ "0",
Died == "0" ~ "1"
)))##predict
#Predict
test_TH1$prediccion_prob <- predict(MODEFECPRIN, newdata = test_TH1, type = "response")
summary(test_TH1$prediccion_prob)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.085 0.854 0.935 0.862 0.972 0.998roc_obj <- roc(test_TH1$muerte, test_TH1$prediccion_prob)## Setting levels: control = 0, case = 1## Setting direction: controls < casesauc(roc_obj)## Area under the curve: 0.755roc<-roc(test_TH1$muerte, test_TH1$prediccion_prob)## Setting levels: control = 0, case = 1
## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
#slope e intercept
#HOSMER & LEMERSHOW.
pred_modelo <- as.data.frame(predict(MODEFECPRIN, test_TH1, type = "response")) %>% rename(pred = `predict(MODEFECPRIN, test_TH1, type = "response")`)
pred_modelo$verdad <- test_TH1$muerte
H_L <- HLtest(MODEFECPRIN,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ COPD + diabetes + CKD + lacticofactor +
## bicarbonatefactor + ureanit + urinefactor + calcio + respfactor +
## anemia, family = "binomial", data = DIED0)
## ChiSquare df P_value
## 5.31 8 0.724HL_DF[,1:4]## cut total obs exp
## 1 [0.0245,0.654] 83 47 48.154
## 2 (0.654,0.805] 82 22 21.981
## 3 (0.805,0.881] 82 16 12.587
## 4 (0.881,0.915] 83 4 8.330
## 5 (0.915,0.94] 82 7 5.869
## 6 (0.94,0.958] 82 5 4.046
## 7 (0.958,0.97] 83 3 2.977
## 8 (0.97,0.979] 82 3 2.065
## 9 (0.979,0.988] 82 1 1.350
## 10 (0.988,0.999] 83 0 0.641HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.0245,0.654] | 83 | 47 | 48.154 | -0.166 | 1 |
| (0.654,0.805] | 82 | 22 | 21.981 | 0.004 | 2 |
| (0.805,0.881] | 82 | 16 | 12.587 | 0.962 | 3 |
| (0.881,0.915] | 83 | 4 | 8.330 | -1.500 | 4 |
| (0.915,0.94] | 82 | 7 | 5.869 | 0.467 | 5 |
| (0.94,0.958] | 82 | 5 | 4.046 | 0.475 | 6 |
| (0.958,0.97] | 83 | 3 | 2.977 | 0.014 | 7 |
| (0.97,0.979] | 82 | 3 | 2.065 | 0.651 | 8 |
| (0.979,0.988] | 82 | 1 | 1.350 | -0.301 | 9 |
| (0.988,0.999] | 83 | 0 | 0.641 | -0.801 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 51 301str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
##comparacion coeficiente
modelo_test <- glm(formula(MODEFECPRIN), data = test_TH1, family = binomial)
coef_train <- coef(MODEFECPRIN)
ic_train <- confint(MODEFECPRIN)## Waiting for profiling to be done...coef_test <- coef(modelo_test)
comparacion_coef <- data.frame(
Coeficiente = names(coef_train),
Beta_train = round(coef_train, 3),
IC95_inf = round(ic_train[, 1], 3),
IC95_sup = round(ic_train[, 2], 3),
Beta_test = round(coef_test, 3),
Test_dentro_del_IC95 = coef_test >= ic_train[, 1] & coef_test <= ic_train[, 2]
)
comparacion_coef## Coeficiente Beta_train IC95_inf IC95_sup Beta_test
## (Intercept) (Intercept) -4.907 -8.652 -1.276 -2.499
## COPD1 COPD1 1.504 0.287 3.159 -0.401
## diabetes1 diabetes1 0.269 -0.250 0.801 0.543
## CKD1 CKD1 1.525 0.895 2.205 0.816
## lacticofactor1 lacticofactor1 -0.959 -1.468 -0.444 -0.837
## bicarbonatefactor1 bicarbonatefactor1 0.701 0.157 1.235 0.126
## ureanit ureanit -0.035 -0.047 -0.023 -0.027
## urinefactor1 urinefactor1 0.942 0.442 1.452 1.140
## calcio calcio 0.808 0.360 1.275 0.541
## respfactor1 respfactor1 -0.649 -1.185 -0.099 -0.601
## anemia1 anemia1 0.940 0.384 1.537 0.401
## Test_dentro_del_IC95
## (Intercept) TRUE
## COPD1 FALSE
## diabetes1 TRUE
## CKD1 FALSE
## lacticofactor1 TRUE
## bicarbonatefactor1 FALSE
## ureanit TRUE
## urinefactor1 TRUE
## calcio TRUE
## respfactor1 TRUE
## anemia1 TRUEknitr::kable(comparacion_coef, caption = "ComparaciĆ³n de coeficientes entre Train y Test")| Coeficiente | Beta_train | IC95_inf | IC95_sup | Beta_test | Test_dentro_del_IC95 | |
|---|---|---|---|---|---|---|
| (Intercept) | (Intercept) | -4.907 | -8.652 | -1.276 | -2.499 | TRUE |
| COPD1 | COPD1 | 1.504 | 0.287 | 3.159 | -0.401 | FALSE |
| diabetes1 | diabetes1 | 0.269 | -0.250 | 0.801 | 0.543 | TRUE |
| CKD1 | CKD1 | 1.525 | 0.895 | 2.205 | 0.816 | FALSE |
| lacticofactor1 | lacticofactor1 | -0.959 | -1.468 | -0.444 | -0.837 | TRUE |
| bicarbonatefactor1 | bicarbonatefactor1 | 0.701 | 0.157 | 1.235 | 0.126 | FALSE |
| ureanit | ureanit | -0.035 | -0.047 | -0.023 | -0.027 | TRUE |
| urinefactor1 | urinefactor1 | 0.942 | 0.442 | 1.452 | 1.140 | TRUE |
| calcio | calcio | 0.808 | 0.360 | 1.275 | 0.541 | TRUE |
| respfactor1 | respfactor1 | -0.649 | -1.185 | -0.099 | -0.601 | TRUE |
| anemia1 | anemia1 | 0.940 | 0.384 | 1.537 | 0.401 | TRUE |
##modelo sin false
modjao<- glm(muerte ~ diabetes + lacticofactor + ureanit+urinefactor+ calcio+ respfactor +anemia, family = "binomial", data = DIED0)
tab_model(modjao)| Ā | muerte | ||
|---|---|---|---|
| Predictors | Odds Ratios | CI | p |
| (Intercept) | 0.00 | 0.00Ā āĀ 0.08 | 0.001 |
| diabetes [1] | 1.42 | 0.87Ā āĀ 2.34 | 0.166 |
| lacticofactor [1] | 0.29 | 0.18Ā āĀ 0.48 | <0.001 |
| ureanit | 0.97 | 0.96Ā āĀ 0.98 | <0.001 |
| urinefactor [1] | 2.62 | 1.63Ā āĀ 4.25 | <0.001 |
| calcio | 2.79 | 1.84Ā āĀ 4.30 | <0.001 |
| respfactor [1] | 0.47 | 0.28Ā āĀ 0.79 | 0.004 |
| anemia [1] | 2.89 | 1.69Ā āĀ 5.14 | <0.001 |
| Observations | 824 | ||
| R2 Tjur | 0.211 | ||
summary(modjao)##
## Call:
## glm(formula = muerte ~ diabetes + lacticofactor + ureanit + urinefactor +
## calcio + respfactor + anemia, family = "binomial", data = DIED0)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -5.99948 1.77773 -3.37 0.00074 ***
## diabetes1 0.34819 0.25123 1.39 0.16575
## lacticofactor1 -1.22293 0.24518 -4.99 0.000000610 ***
## ureanit -0.02614 0.00488 -5.35 0.000000086 ***
## urinefactor1 0.96140 0.24437 3.93 0.000083475 ***
## calcio 1.02495 0.21643 4.74 0.000002183 ***
## respfactor1 -0.76008 0.26091 -2.91 0.00358 **
## anemia1 1.06076 0.28235 3.76 0.00017 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 640.10 on 823 degrees of freedom
## Residual deviance: 498.38 on 816 degrees of freedom
## AIC: 514.4
##
## Number of Fisher Scoring iterations: 6#HOSMER & LEMERSHOW. LA P ME DA . NO RECHAZA LA HIPĆTESIS NULA.ES COMO QUIERO QUE ME DE
pred_modelo <- as.data.frame(predict(modjao, DIED0, type = "response")) %>% rename(pred = `predict(modjao, DIED0, type = "response")`)
pred_modelo$verdad <- DIED0$muerte
H_L <- HLtest(modjao,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ diabetes + lacticofactor + ureanit + urinefactor +
## calcio + respfactor + anemia, family = "binomial", data = DIED0)
## ChiSquare df P_value
## 3.89 8 0.867HL_DF[,1:4]## cut total obs exp
## 1 [0.101,0.672] 83 42 42.98
## 2 (0.672,0.791] 82 20 21.12
## 3 (0.791,0.866] 82 12 13.57
## 4 (0.866,0.905] 83 11 9.42
## 5 (0.905,0.93] 82 9 6.79
## 6 (0.93,0.949] 82 7 4.91
## 7 (0.949,0.962] 83 4 3.69
## 8 (0.962,0.973] 82 1 2.69
## 9 (0.973,0.982] 82 1 1.84
## 10 (0.982,0.999] 83 1 1.00HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.101,0.672] | 83 | 42 | 42.98 | -0.150 | 1 |
| (0.672,0.791] | 82 | 20 | 21.12 | -0.243 | 2 |
| (0.791,0.866] | 82 | 12 | 13.57 | -0.425 | 3 |
| (0.866,0.905] | 83 | 11 | 9.42 | 0.515 | 4 |
| (0.905,0.93] | 82 | 9 | 6.79 | 0.846 | 5 |
| (0.93,0.949] | 82 | 7 | 4.91 | 0.946 | 6 |
| (0.949,0.962] | 83 | 4 | 3.69 | 0.163 | 7 |
| (0.962,0.973] | 82 | 1 | 2.69 | -1.029 | 8 |
| (0.973,0.982] | 82 | 1 | 1.84 | -0.618 | 9 |
| (0.982,0.999] | 83 | 1 | 1.00 | 0.000 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 108 716str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
#curva roc. Area 0.83 sin temperatura.vamos con este modelo prueba i3
roc<-roc(pred_modelo$verdad, pred_modelo$pred)## Setting levels: control = 0, case = 1## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
##predict modjao
#Predict
test_TH1$prediccion_prob <- predict(modjao, newdata = test_TH1, type = "response")
summary(test_TH1$prediccion_prob)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.129 0.839 0.929 0.864 0.964 0.994roc_obj <- roc(test_TH1$muerte, test_TH1$prediccion_prob)## Setting levels: control = 0, case = 1## Setting direction: controls < casesauc(roc_obj)## Area under the curve: 0.76roc<-roc(test_TH1$muerte, test_TH1$prediccion_prob)## Setting levels: control = 0, case = 1
## Setting direction: controls < casesplot(roc, main = "Curva ROC", col = "red", lwd = 3)
text(0.8, 0.2, paste("AUC ROC =", round(roc$auc, 3)), adj = c(0, 1))
#slope e intercept
#HOSMER & LEMERSHOW.
pred_modelo <- as.data.frame(predict(modjao, test_TH1, type = "response")) %>% rename(pred = `predict(modjao, test_TH1, type = "response")`)
pred_modelo$verdad <- test_TH1$muerte
H_L <- HLtest(modjao,10)
HL_DF <- cbind(as.data.frame(H_L$table), indice = as.factor(c(1,2,3,4,5,6,7,8,9,10)))
H_L## Hosmer and Lemeshow Goodness-of-Fit Test
##
## Call:
## glm(formula = muerte ~ diabetes + lacticofactor + ureanit + urinefactor +
## calcio + respfactor + anemia, family = "binomial", data = DIED0)
## ChiSquare df P_value
## 3.89 8 0.867HL_DF[,1:4]## cut total obs exp
## 1 [0.101,0.672] 83 42 42.98
## 2 (0.672,0.791] 82 20 21.12
## 3 (0.791,0.866] 82 12 13.57
## 4 (0.866,0.905] 83 11 9.42
## 5 (0.905,0.93] 82 9 6.79
## 6 (0.93,0.949] 82 7 4.91
## 7 (0.949,0.962] 83 4 3.69
## 8 (0.962,0.973] 82 1 2.69
## 9 (0.973,0.982] 82 1 1.84
## 10 (0.982,0.999] 83 1 1.00HL_DF %>% knitr::kable()| cut | total | obs | exp | chi | indice |
|---|---|---|---|---|---|
| [0.101,0.672] | 83 | 42 | 42.98 | -0.150 | 1 |
| (0.672,0.791] | 82 | 20 | 21.12 | -0.243 | 2 |
| (0.791,0.866] | 82 | 12 | 13.57 | -0.425 | 3 |
| (0.866,0.905] | 83 | 11 | 9.42 | 0.515 | 4 |
| (0.905,0.93] | 82 | 9 | 6.79 | 0.846 | 5 |
| (0.93,0.949] | 82 | 7 | 4.91 | 0.946 | 6 |
| (0.949,0.962] | 83 | 4 | 3.69 | 0.163 | 7 |
| (0.962,0.973] | 82 | 1 | 2.69 | -1.029 | 8 |
| (0.973,0.982] | 82 | 1 | 1.84 | -0.618 | 9 |
| (0.982,0.999] | 83 | 1 | 1.00 | 0.000 | 10 |
ggplot(HL_DF, aes(x=indice))+geom_point(aes(y=exp, color = "blue"))+geom_point(aes(y=obs, color = "red"))+scale_color_manual(labels = c("Esperado", "Observado"), values = c("blue", "red"))
#slope e intercept
table(pred_modelo$verdad)##
## 0 1
## 51 301str(pred_modelo$verdad)## Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...pred_modelo$verdad <- as.numeric(as.character(pred_modelo$verdad))
pred_modelo <- pred_modelo# esta parece ser una observacion duplicada
a <- valProbggplot(pred_modelo$pred, pred_modelo$verdad)
a$ggPlot
modelo_test <- glm(formula(modjao), data = test_TH1, family = binomial)
coef_train <- coef(modjao)
ic_train <- confint(modjao)## Waiting for profiling to be done...coef_test <- coef(modelo_test)
comparacion_coef <- data.frame(
Coeficiente = names(coef_train),
Beta_train = round(coef_train, 3),
IC95_inf = round(ic_train[, 1], 3),
IC95_sup = round(ic_train[, 2], 3),
Beta_test = round(coef_test, 3),
Test_dentro_del_IC95 = coef_test >= ic_train[, 1] & coef_test <= ic_train[, 2]
)
comparacion_coef## Coeficiente Beta_train IC95_inf IC95_sup Beta_test
## (Intercept) (Intercept) -5.999 -9.544 -2.563 -3.286
## diabetes1 diabetes1 0.348 -0.139 0.849 0.658
## lacticofactor1 lacticofactor1 -1.223 -1.704 -0.740 -0.813
## ureanit ureanit -0.026 -0.036 -0.017 -0.023
## urinefactor1 urinefactor1 0.961 0.487 1.448 1.154
## calcio calcio 1.025 0.609 1.459 0.645
## respfactor1 respfactor1 -0.760 -1.267 -0.242 -0.661
## anemia1 anemia1 1.061 0.526 1.637 0.533
## Test_dentro_del_IC95
## (Intercept) TRUE
## diabetes1 TRUE
## lacticofactor1 TRUE
## ureanit TRUE
## urinefactor1 TRUE
## calcio TRUE
## respfactor1 TRUE
## anemia1 TRUEknitr::kable(comparacion_coef, caption = "ComparaciĆ³n de coeficientes entre Train y Test")| Coeficiente | Beta_train | IC95_inf | IC95_sup | Beta_test | Test_dentro_del_IC95 | |
|---|---|---|---|---|---|---|
| (Intercept) | (Intercept) | -5.999 | -9.544 | -2.563 | -3.286 | TRUE |
| diabetes1 | diabetes1 | 0.348 | -0.139 | 0.849 | 0.658 | TRUE |
| lacticofactor1 | lacticofactor1 | -1.223 | -1.704 | -0.740 | -0.813 | TRUE |
| ureanit | ureanit | -0.026 | -0.036 | -0.017 | -0.023 | TRUE |
| urinefactor1 | urinefactor1 | 0.961 | 0.487 | 1.448 | 1.154 | TRUE |
| calcio | calcio | 1.025 | 0.609 | 1.459 | 0.645 | TRUE |
| respfactor1 | respfactor1 | -0.760 | -1.267 | -0.242 | -0.661 | TRUE |
| anemia1 | anemia1 | 1.061 | 0.526 | 1.637 | 0.533 | TRUE |