Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- run_shiny(model = "SIRS (w/out demography)",
- neweqns = mySIRS,
- ics = c(S = 9999, I = 1, R = 0),
- parm0 = c(beta = 5e-5, gamma = 1/7, delta = 0.1),
- parm_names = c("Transmission rate", "Recovery rate", "Loss of immunity"),
- parm_min = c(beta = 1e-5, gamma = 1/21, delta = 1/9999999),
- parm_max = c(beta = 9e-5, gamma = 1 , delta = 1))
- run_shiny(model = "mySIR")
- neweqns = mySIR,
- ics = c(S = 9999, I = 1, R = 0),
- parm0 = c(Transmission rate = 1/3, Recovery rate = 1/7),
- parm_names = c("Transmission rate", "Recovery rate"),
- parm_min = c(Transmission rate = 1/30, Recovery rat = 1/21),
- parm_max = c(Transmission rate = 0, Recovery rate = 1))
- mySIR <- function(t, y, parms) {
- with(as.list(c(y, parms)),{
- # Change in Susceptibles
- dS <- - beta * S * I
- # Change in Infecteds
- dI <- beta * S * I - gamma * I
- # Change in Recovereds
- dR <- gamma * I
- return(list(c(dS, dI, dR)))
- })
- }
- parm0 = c(R0 = 3, Ip = 7)
- parm_names = c("R0", "Infectious period")
- parm_min = c(R0 = 0, Ip = 1)
- parm_max = c(R0 = 20, Ip = 21
- ## Load deSolve package
- library(deSolve)
- ## Create an SIR function
- sir <- function(time, state, parameters) {
- with(as.list(c(state, parameters)), {
- dS <- -beta * S * I
- dI <- beta * S * I - gamma * I
- dR <- gamma * I
- return(list(c(dS, dI, dR)))
- })
- }
- ### Set parameters
- ## Proportion in each compartment: Susceptible 0.999999, Infected 0.000001, Recovered 0
- init <- c(S = 1-1e-6, I = 1e-6, R = 0.0)
- ## beta: infection parameter; gamma: recovery parameter
- parameters <- c(beta = 1.4247, gamma = 0.14286)
- ## Time frame
- times <- seq(0, 200, by = 1)
- ## Solve using ode (General Solver for Ordinary Differential Equations)
- out <- ode(y = init, times = times, func = mySIR, parms = parameters)
- ## change to data frame
- out <- as.data.frame(out)
- ## Delete time variable
- out$time <- NULL
- ## Show data
- head(out, 10)
- mySIR <- function(t, y, parms) {
- +
- + with(as.list(c(y, parms)),{
- +
- + # Change in Susceptibles
- + dS <- - beta * S * I
- +
- + # Change in Infecteds
- + dI <- beta * S * I - gamma * I
- +
- + # Change in Recovereds
- + dR <- gamma * I
- +
- + return(list(c(dS, dI, dR)))
- + })
- + }
- init <- c(S = 0.999999, I = 0.000001, R = 0.0)
- parameters <- c(beta = 1.1, gamma = 0.3)
- times <- seq(0, 60, by = 1)
- out <- ode(y = init, times = times, func = mySIR, parms = parameters)
- out <- as.data.frame(out)
- out$time <- NULL
- head(out, 10)
- #after data output
- matplot(x = times, y = out, type = "l",
- xlab = "Time", ylab = "Susceptible and Recovered", main = "SIR Model",
- lwd = 1, lty = 1, bty = "l", col = 2:4)
- legend(40, 0.7, c("Susceptible", "Infected", "Removed"), pch = 1, col = 2:4, bty = "n")
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement