Kemal Koçaklı R Notebook

library(deSolve)
library(ggplot2)

parameters <- c(d = 0.5, f = 0.5, r = 0.5, s=13, a=73, B=80)
state <- c(X = 200, I = 20, Y = 1)
#X number of susceptible foxes
#I number of incubating foxes
#Y number of infected foxes
#b birth rate per capita
#d death rate per capita
#r population growth rate
#s rate of incubation
#a death rate of rabid foxes
#b transmission coefficient
#K fox carrying capacity

N=100
Lorenz <- function(t, state, parameters) {
  with(as.list(c(state, parameters )), {
    dX <- r*X - f*X*N - B*X*Y
    dI <- B*X*Y - (s+d+f*N)*I
    dY <- f*t - (a+d+f*N)*Y
    list(c(dX, dI, dY))
  })
}
times <- seq(0, 1, by = 0.01)
out <- ode(y = state, times = times, func = Lorenz, parms = parameters)
head(out)

     time          X         I           Y
[1,] 0.00 200.000000 20.000000 1.000000000
[2,] 0.01  77.009903 51.705353 0.290852015
[3,] 0.02  41.071046 32.795173 0.084635815
[4,] 0.03  24.080212 18.258112 0.024689681
[5,] 0.04  14.512375  9.828561 0.007283975
[6,] 0.05   8.816238  5.236199 0.002250508

outf <- data.frame(times, out)
ggplot(outf, aes(x = X, y = Y)) +
  geom_point() 

X exeni hassas olan tilki sayısı Y ekseni ise enfekte olan tilki sayısını göstermektedir. transmission coefficient (80) ise enfeksiyonun ayayılma hızıdır. Örneğin 80 ni 650 yaptığımızda enfekte olan tilki sayısı dramatik bir şekilde artıyor. Bu değeri 0 girdiğimizde ise daha yavaş olüler gerçekleşiyor.