Chunk con todos los elementos para modelar el caso.
library(deSolve)
InitialConditions <- c(crimenes.urb = 100,
crimenes.suburb = 50,
crimenes.rural = 200
)
times <- seq(0,#initial time
50,#end time
1/4) #time step
politica <- function(t, state, parameters) {
with(as.list(c(state,parameters)), {
#Endogenous auxiliary variables
crimenes.totales = crimenes.urb + crimenes.suburb + crimenes.rural
#Flow variables
incremento.crimenes.urb = poblacion.urb*tasa.incremento
incremento.crimenes.suburb = poblacion.suburb*tasa.incremento
incremento.crimenes.rural = poblacion.rural*tasa.incremento
reduccion.crimenes.urb = ifelse(t>=años.politica.efectiva, tasa.reduccion.urb*crimenes.urb, 0)
reduccion.crimenes.suburb = ifelse(t>=años.politica.efectiva, tasa.reduccion.suburb*crimenes.suburb, 0)
reduccion.crimenes.rural = ifelse(t>=años.politica.efectiva, tasa.reduccion.rural*crimenes.rural, 0)
#State (stock) variables (d)
d.crimenes.urb = incremento.crimenes.urb - reduccion.crimenes.urb
d.crimenes.suburb = incremento.crimenes.suburb - reduccion.crimenes.suburb
d.crimenes.rural = incremento.crimenes.rural - reduccion.crimenes.rural
list(c(d.crimenes.urb, d.crimenes.suburb, d.crimenes.rural), crimenes.totales = crimenes.totales)
})
}
parameters<-c(poblacion.urb = 5000000,
poblacion.suburb = 7000000,
poblacion.rural = 4000000,
tasa.incremento = 10/1000000,
tasa.reduccion.urb = 5/100,
tasa.reduccion.suburb = 7/100,
tasa.reduccion.rural = 10/100,
años.politica.efectiva = 10
)
intg.method<-c("rk4")
out <- ode(y = InitialConditions,
times = times,
func = politica,
parms = parameters,
method =intg.method )
plot(out,
col=c("blue"))
#Mejora de la salida de los gráficos.
library(ggplot2)
library(tidyr)
out_df <- as.data.frame(out)
# 1. Gráfico de crÃmenes por zona
out_long <- pivot_longer(out_df,
cols = c(crimenes.urb, crimenes.suburb, crimenes.rural),
names_to = "Zona",
values_to = "Crimenes")
ggplot(out_long, aes(x = time, y = Crimenes, color = Zona)) +
geom_line() +
labs(title = "CrÃmenes por zona", x = "Tiempo", y = "CrÃmenes")
# 2. Gráfico de crÃmenes totales
ggplot(out_df, aes(x = time, y = crimenes.totales)) +
geom_line(color = "red") +
labs(title = "CrÃmenes totales", x = "Tiempo", y = "CrÃmenes totales")