setwd("~/ONLINEDS")
#Cargar los paquetes necesarios
library(pacman)
## Warning: package 'pacman' was built under R version 3.6.3
p_load(markdown,knitr,dplyr,tidyr, tidyverse, hashmap, lubridate,
summarytools,ggpubr, kableExtra, reshape2,
sf, tmap, readr, devtools, plotly, gganimate, gifski)
#Primera parte datos diarios y acumulados
#Leer datos
datos<- read.csv("Casos_Diarios_Estado_Nacional_Confirmados.csv")
#Casos diarios Sonora
Sonora<- t(datos[datos$nombre=="SONORA" ,])
Sonora <- as.vector(Sonora)
Sonora<-Sonora[4:121]
Sonora <- as.numeric(Sonora)
Sonora <- as.vector(Sonora)
csonora <- cumsum(Sonora)
#Casos diarios Sinaloa
Sinaloa<- t(datos[datos$nombre=="SINALOA" ,])
Sinaloa <- as.vector(Sinaloa)
Sinaloa<-Sinaloa[4:121]
Sinaloa <- as.numeric(Sinaloa)
Sinaloa <- as.vector(Sinaloa)
csinaloa <- cumsum(Sinaloa)
#Fecha
Fecha = seq(from = as.Date("2020-01-05"), to = as.Date("2020-05-01"), by = 'day')
SonSin <-data.frame(Fecha, Sonora, Sinaloa)
cSonSin <-data.frame(Fecha, csonora, csinaloa)
#Graficos de casos diarios confirmados
plot(Sonora)

#Grafica Sonora
ggplot(data = SonSin) +
ggtitle("Casos diarios COVID-19 en Sonora")+
geom_line(mapping = aes(x = Fecha, y = Sonora ))

# Sinaloa y Sonora
ggplot(data=cSonSin) +
geom_line(aes(Fecha, csinaloa, colour = 'Sinaloa')) +
geom_line(aes(Fecha, csonora, colour ='Sonora') )+
xlab('Fecha')+
ylab('Casos Diarios') +
labs(colour = "Estados")+
transition_reveal(Fecha)

#numero semilla
set.seed(42)
# Declaramos variable 'temp' como un archivo temporal
temp <- tempfile()
ArchivoZip = "http://187.191.75.115/gobmx/salud/datos_abiertos/datos_abiertos_covid19.zip"
download.file(ArchivoZip, temp)
# Utilizamos la función 'unz' para extraer el archivo CSV y lo asignamos a la variable 'temp'
ArchivoCsv = unz(temp, "200503COVID19MEXICO.csv")
# Introducimos los datos del CSV en la tabla
datos <- read.csv(file=ArchivoCsv, header = TRUE, encoding = "UTF-8")
# Eliminamos la referencia al archivo temporal y eliminamos los archivos
unlink(temp)
remove(ArchivoCsv)
remove(temp)
remove(ArchivoZip)
#Asignamos otra variable a datos
datos2 <- datos
datos2$Num <- 1
#Diccionarios
d_origen <- hashmap(c(1,2,99), c("USMER","Fuera USMER","No especificado"))
d_sector <- hashmap(c(1,2,3,4,5,6,7,8,9,10,11,12,13,99),
c("CRUZ ROJA","DIF","ESTATAL","IMSS","IMSS-BIENESTAR","ISSSTE","MUNICIPAL","PEMEX",
"PRIVADA","SEDENA","SEMAR","SSA","UNIVERSITARIO","NO ESPECIFICADO"))
d_sexo <- hashmap(c(1,2,99), c("MUJER","HOMBRE","NO ESPECIFICADO"))
d_tipo <- hashmap(c(1,2,99), c("AMBULATORIO","HOSPITALIZADO","NO ESPECIFICADO"))
d_si_no <- hashmap(c(1,2,97,98,99), c("SI","NO","NO APLICA","SE IGNORA","NO ESPECIFICADO"))
d_nacion <- hashmap(c(1,2,99), c("MEXICANA","EXTRANJERA","NO ESPECIFICADO"))
d_resultado <- hashmap(c("1","2","3"),c("Positivo SARS-CoV-2","No positivo SARS-CoV-2","Resultado pendiente"))
d_entidad <- hashmap(c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,
36,97,98,99),
c("AGU","BCN","BCS","CAM","COA","COL","CHP","CHH","CMX","DUR","GUA","GRO","HID","JAL",
"MEX","MIC","MOR","NAY","NLE","OAX","PUE","QUE","ROO","SLP","SIN","SON","TAB","TAM",
"TLA","VER","YUC","ZAC","EUM","N_A","S_I","N_E"))
#municipios
datos2$MUNICIPIO_RES <- paste(sprintf("%02d", datos2$ENTIDAD_RES), sprintf("%03d", datos2$MUNICIPIO_RES))
datos2$MUNICIPIO_RES <- as.factor(gsub("[[:space:]]", "", datos2$MUNICIPIO_RES))
#decodificacion
datos2$ORIGEN <- as.factor(d_origen[[datos2$ORIGEN]])
datos2$SECTOR <- as.factor(d_sector[[datos2$SECTOR]])
datos2$ENTIDAD_UM <- as.factor(d_entidad[[datos2$ENTIDAD_UM]])
datos2$SEXO <- as.factor(d_sexo[[datos2$SEXO]])
datos2$ENTIDAD_NAC <- as.factor(d_entidad[[datos2$ENTIDAD_NAC]])
datos2$ENTIDAD_RES <- as.factor(d_entidad[[datos2$ENTIDAD_RES]])
datos2$TIPO_PACIENTE <- as.factor(d_tipo[[datos2$TIPO_PACIENTE]])
datos2$INTUBADO <- as.factor(d_si_no[[datos2$INTUBADO]])
datos2$NEUMONIA <- as.factor(d_si_no[[datos2$NEUMONIA]])
datos2$NACIONALIDAD <- as.factor(d_nacion[[datos2$NACIONALIDAD]])
datos2$EMBARAZO <- as.factor(d_si_no[[datos2$EMBARAZO]])
datos2$HABLA_LENGUA_INDIG <- as.factor(d_si_no[[datos2$HABLA_LENGUA_INDIG]])
datos2$DIABETES <- as.factor(d_si_no[[datos2$DIABETES]])
datos2$EPOC <- as.factor(d_si_no[[datos2$EPOC]])
datos2$ASMA <- as.factor(d_si_no[[datos2$ASMA]])
datos2$INMUSUPR <- as.factor(d_si_no[[datos2$INMUSUPR]])
datos2$HIPERTENSION <- as.factor(d_si_no[[datos2$HIPERTENSION]])
datos2$OTRA_COM <- as.factor(d_si_no[[datos2$OTRA_COM]])
datos2$CARDIOVASCULAR <- as.factor(d_si_no[[datos2$CARDIOVASCULAR]])
datos2$OBESIDAD <- as.factor(d_si_no[[datos2$OBESIDAD]])
datos2$RENAL_CRONICA <- as.factor(d_si_no[[datos2$RENAL_CRONICA]])
datos2$TABAQUISMO <- as.factor(d_si_no[[datos2$TABAQUISMO]])
datos2$OTRO_CASO <- as.factor(d_si_no[[datos2$OTRO_CASO]])
datos2$RESULTADO <- as.factor(d_resultado[[datos2$RESULTADO]])
datos2$MIGRANTE <- as.factor(d_si_no[[datos2$MIGRANTE]])
datos2$UCI <- as.factor(d_si_no[[datos2$UCI]])
levels(datos2$PAIS_NACIONALIDAD) <- c(levels(datos2$PAIS_NACIONALIDAD), "Se ignora")
datos2$PAIS_NACIONALIDAD[datos2$PAIS_NACIONALIDAD == 99] <- "Se ignora"
levels(datos2$PAIS_ORIGEN) <- c(levels(datos2$PAIS_ORIGEN), "No aplica")
datos2$PAIS_ORIGEN[datos2$PAIS_ORIGEN == 97] <- "No aplica"
#Formato de Fecha
datos$FECHA_ACTUALIZACION <- ymd(datos$FECHA_ACTUALIZACION)
datos2$FECHA_INGRESO <- ymd(datos2$FECHA_INGRESO)
datos2$FECHA_SINTOMAS <- ymd(datos2$FECHA_SINTOMAS)
datos2$FECHA_DEF <- na_if(datos2$FECHA_DEF,"9999-99-99")
datos2$FECHA_DEF <- ymd(datos2$FECHA_DEF)
write.csv(datos2, file = "tabla_procesada.csv")
#Estadistica descriptiva
#view(dfSummary(datos2))
print(dfSummary(datos2, graph.magnif = 0.75), method = 'render', headings = FALSE)
#resumen
tab_resultado <- count(datos2, Resultado = datos2$RESULTADO)
tab_resultado["Porcentaje"] <- round(tab_resultado$n*100 / sum(tab_resultado$n),1)
kable(tab_resultado, format.args = list(big.mark = ",")) %>%
kable_styling(bootstrap_options = "striped", full_width = F, position = "center")
Resultado
|
n
|
Porcentaje
|
No positivo SARS-CoV-2
|
59,704
|
62.3
|
Positivo SARS-CoV-2
|
23,471
|
24.5
|
Resultado pendiente
|
12,664
|
13.2
|
#grafica de pie
tab_resultado_lab <- paste0(tab_resultado$Resultado, "\n (", tab_resultado$Porcentaje, "%)")
ggpie(tab_resultado, "Porcentaje", label = tab_resultado_lab,
lab.pos = "in", lab.font = "black", lab.adjust = 0,
fill = "Resultado", color = "white",
palette = c("#C1D40E", "#FF9999", "#EAEAEA"))

datos3 <- subset(datos2, RESULTADO == "Positivo SARS-CoV-2")
#casos por fecha
tab_fecha_in <- count(datos3, Fecha = datos3$FECHA_INGRESO)
tab_fecha_in["Acumulado"] <- cumsum(tab_fecha_in$n)
g_sum <- ggbarplot(tab_fecha_in, x = "Fecha", y = "n",
fill = "steelblue", color = "black",)
g_acum <- ggline(tab_fecha_in, x = "Fecha", y = "Acumulado",
color = "steelblue", point.color = "black", point.size = 0.2) +
scale_y_continuous(labels=function(x) format(x, big.mark = ",", scientific = FALSE))
ggarrange(g_sum, g_acum, labels = c("A)", "B)"), ncol = 2, nrow = 1)

tab_fecha_in2 <- tab_fecha_in %>% mutate(Dif_dias = difftime(tab_fecha_in$Fecha, lag(tab_fecha_in$Fecha,1)))
tab_fecha_in2$Dif_dias <- replace_na(tab_fecha_in2$Dif_dias,1)
tab_fecha_in2$Dif_dias <- as.integer(tab_fecha_in2$Dif_dias)
tab_fecha_in2["Dias"] <- cumsum(tab_fecha_in2$Dif_dias)
datos_log <- data.frame(Dias = tab_fecha_in2$Dias, Acumulado = tab_fecha_in2$Acumulado)
datos_log <- head(datos_log, nrow(datos_log) - 3)
mod_log <- nls(Acumulado ~ SSlogis(Dias, phi1, phi2, phi3), data = datos_log)
summary(mod_log)
##
## Formula: Acumulado ~ SSlogis(Dias, phi1, phi2, phi3)
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## phi1 3.969e+04 1.123e+03 35.33 <2e-16 ***
## phi2 1.122e+02 5.167e-01 217.13 <2e-16 ***
## phi3 9.368e+00 1.279e-01 73.25 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 164.3 on 61 degrees of freedom
##
## Number of iterations to convergence: 1
## Achieved convergence tolerance: 9.344e-07
p_log1 <- max(datos_log$Dias)
p_log2 <- round(p_log1+5, -1)
p_log3 <- p_log2 + 15
p_log4 <- p_log3 + 15
pred_log <- predict(mod_log,data.frame(Dias = c(p_log1, p_log2, p_log3, p_log4)))
max_pred <-round(pred_log[4]+500,-3)
pred_log
## [1] 22802.28 27671.01 36491.00 38997.84
## attr(,"gradient")
## phi1 phi2 phi3
## [1,] 0.5745591 -1035.49864 -311.1435
## [2,] 0.6972387 -894.24693 -745.9656
## [3,] 0.9194799 -313.63338 -763.7917
## [4,] 0.9826459 -72.23958 -291.5892
dia1_log <- min(datos3$FECHA_INGRESO)
alpha <- coef(mod_log) #Coeficientes
plot(Acumulado ~ Dias, data = datos_log, main = "Modelo de crecimiento logístico para SARS-Cov2 en México",
xlab = "Días", ylab = "Casos positivos acumulados", xlim = c(1, p_log4+15), ylim = c(1, max_pred)) #Casos
curve(alpha[1]/(1 + exp(-(x - alpha[2])/alpha[3])), add = T, col = "blue") #Modelo
points(p_log2, pred_log[2], pch = "x", cex = 1.3)
points(p_log3, pred_log[3], pch = "x", cex = 1.3)
points(p_log4, pred_log[4], pch = "x", cex = 1.3)
text(p_log2, pred_log[2], labels = paste0(dia1_log + p_log2, '\n', "[", round(pred_log[2], -2) ," casos]"),
adj = c(1.1,0))
text(p_log3, pred_log[3], labels = paste0(dia1_log + p_log3, '\n', "[", round(pred_log[3], -2) ," casos]"),
adj = c(1.1,0))
text(p_log4, pred_log[4], labels = paste0(dia1_log + p_log4, '\n', "[", round(pred_log[4], -2) ," casos]"),
pos = 1, offset = 0.8)
