Menganalisa Prediksi Penyebaran Covid 19 di wilayah DKI Jakarta 90 hari kedepan dari data Covid 19 DKI Jakarta pada bulan September 2021 dengan jumlah penduduk pada bulan tersebut 10.609.681 Jiwa menggunakan Model SIR.

Model SIR adalah model epidemik yang paling umum digunakan dalam melihat persebaran penyakit. Model SIR yang diperkenalkan oleh Kermack dan McKendrick pada tahun 1927 dan pertama kali digunakan untuk melihat dinamika penyebaran populasi akibat penyakit menular.

library(readxl)
#baca data excel
data.covid <- read_excel("DATACOVIDJKT2.xlsx")
head(data.covid)
## # A tibble: 6 x 4
##   Tanggal             T_Infeksi T_Sembuh T_Meninggal
##   <dttm>                  <dbl>    <dbl>       <dbl>
## 1 2021-09-01 00:00:00    851256   831293       13302
## 2 2021-09-02 00:00:00    851686   832130       13312
## 3 2021-09-03 00:00:00    852029   832952       13322
## 4 2021-09-04 00:00:00    852389   833437       13332
## 5 2021-09-05 00:00:00    852692   833765       13342
## 6 2021-09-06 00:00:00    852909   834489       13359

Membuat data harian dan menggabungkannya dalam data.covid

hari <- 1:30
data.covid <- cbind(data.covid,hari)

Visualisasi Data Covid

library(ggplot2)
plot.covid <- ggplot(data = data.covid, aes(x=hari,y=data.covid$T_Infeksi))+
  geom_point()+geom_smooth(method = "loess", formula = y~x,se=T)+ylab("Jumlah")+
  labs(title = "Total TerInfeksi COVID-19 DKI JKT Sept 2021")+theme(axis.text.x=element_text(
    angle = 45,hjust = 1
  ))
plot.covid
## Warning: Use of `data.covid$T_Infeksi` is discouraged. Use `T_Infeksi` instead.
## Use of `data.covid$T_Infeksi` is discouraged. Use `T_Infeksi` instead.

Membuat pemodelan penyebaran virus corona dengan model SIR berdasarkan jumlah penduduk DKI JKT 2021 Mengubah data ke tipe numerik.

Infected <- as.numeric(data.covid$T_Infeksi)
Day <- as.numeric(data.covid$hari)
N <- 10609681  #Jumlah Penduduk DKI 2021
R <-as.numeric(data.covid$T_Sembuh)

Model SIR

library(deSolve)
## Warning: package 'deSolve' was built under R version 4.1.3
SIR <- function(time, state, parameters){
  par <- as.list(c(state, parameters))
  with(par,{
    dS <- -beta/N*I*S
    dI <- beta/N*I*S-gamma*I
    dR <- gamma*I
    list(c(dS,dI,dR))
  })
}

init <- c(S=N-Infected[1], I=Infected[1], R=0)
RSS <- function(parameters){
  names(parameters) <- c("beta","gamma")
  out <- ode(y=init, times=Day, func=SIR,parms=parameters)
  fit <- out[,3]
  sum(Infected-fit)^2
}

opt <- optim(c(0.5,0.5),RSS,method="L-BFGS-B",lower=c(0,0),upper = c(1,1))
opt$message
## [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

Memanggil nilai beta dan gamma

opt_par <- setNames(opt$par,c("beta","gamma"))
opt_par
##      beta     gamma 
## 0.6578507 0.3380236

Membuat data persebaran dengan waktu persebaran 90 hari

t <- 1:90
fit <- data.frame(ode(y=init, times=t,func=SIR, parms = opt_par))

Visualisasi model SIR

col <- 1:3
matplot(fit$time, fit[,2:4],type="l",xlab = "Day", ylab="Jumlah Penduduk",lwd = 2,lty=1, col=col,log="y")
## Warning in xy.coords(x, y, xlabel, ylabel, log = log, recycle = TRUE): 1 y value
## <= 0 omitted from logarithmic plot
points(Day,Infected)
legend("bottomright",c("Rentan","Terinfeksi","Sembuh"),lty=1,lwd=2,col=col)
title("Model SIR COVID-19 DKI JKT 90 hari mulai Sep 2021",outer=T,line=-2)

Menentukan puncak persebaran dan waktunya

fit[fit$I==max(fit$I),"I",drop=F]
##         I
## 7 1983705

Menentukan jumlah kematian degan tingkat rasio kematian 0.4%

max_infected <- max(fit$I)
max_infected*0.004 
## [1] 7934.82