[R] How to apply a system of ordinary differential equations to a cell grid?
Marine Regis
marine.regis at hotmail.fr
Tue Jul 18 00:55:00 CEST 2017
It seems that the method `iteration` in the function `ode` can be useful in my case:
> Method "iteration" is special in that here the function func should
> return the new value of the state variables rather than the rate of
> change. This can be used for individual based models, for difference
> equations, or in those cases where the integration is performed within
> func).
I have tested the method but I don't understand why it doesn't work with one time step:
solve_sir_model <- function (times, parameters) {
sir_model <- function (times, states, parameters) {
with(as.list(c(states, parameters)), {
dSdt <- -beta*S*I
dIdt <- beta*S*I-gamma*I
dRdt <- gamma*I
dNdt <- dSdt + dIdt + dRdt
return(list(c(dSdt, dIdt, dRdt, dNdt)))
})
}
states <- c(S = 99, I = 1, R = 0, N = 100)
return(ode(y = states, times = times, func = sir_model, parms = parameters, method = "iteration"))
}
require(deSolve)
output <- as.data.frame(solve_sir_model(times = 0, parameters = c(beta = 400, gamma = 28)))
Error in iteration(y, times, func, parms, ...) :
times should be equally spaced In addition: Warning messages:
1: In min(x) : no non-missing arguments to min; returning Inf
2: In max(x) : no non-missing arguments to max; returning -Inf
Thank you very much for your help.
Marine
________________________________
De : David Winsemius <dwinsemius at comcast.net>
Envoyé : mercredi 21 juin 2017 22:30
À : Marine Regis
Cc : r-help at r-project.org
Objet : Re: [R] How to apply a system of ordinary differential equations to a cell grid?
> On Jun 21, 2017, at 12:48 PM, Marine Regis <marine.regis at hotmail.fr> wrote:
>
> Hello,
>
> I am developing an agent-based model to simulate the spread of infectious diseases in heterogeneous landscapes composed of habitat polygons (or clumps of connected cells). To simplify the model, I consider a habitat grid (or raster) containing the polygon ID of each cell. In addition, I have epidemiological parameters associated with each polygon ID. At each time step, the parameter values change in the polygon. Thus, the data frame �landscape� (see below) is updated at each time step. Here is an example at t = 0:
>
> landscape <- data.frame(polygon_ID = seq(1, 10, by = 1),
> beta = sample(c(100, 200, 400, 600), 10, replace = TRUE),
> gamma = sample(c(25, 26, 27, 28), 10, replace = TRUE))
>
> To study the disease dynamics, I also am developing a compartmental model based on a system of ordinary differential equations (ODEs). Here is an example to represent the system of ODEs:
>
> solve_sir_model <- function (times, parameters) {
>
> sir_model <- function (times, states, parameters) {
>
> with(as.list(c(states, parameters)), {
>
> dSdt <- -beta*S*I
> dIdt <- beta*S*I-gamma*I
> dRdt <- gamma*I
> dNdt <- dSdt + dIdt + dRdt
>
> return(list(c(dSdt, dIdt, dRdt, dNdt)))
>
> })
> }
>
> states <- c(S = 99, I = 1, R = 0, N = 100)
> return(ode(y = states, times = times, func = sir_model, parms = parameters))
> }
>
> require(deSolve)
> output <- as.data.frame(solve_sir_model(times = seq(0, 5, by = 1), parameters = c(beta = 400, gamma = 28)))
>
> Here is my question: at each time step, is it possible to apply the system of ODEs to each habitat polygon (thus each row) in the data frame �landscape�? I am using lsoda as an ODE solver. Do I need to use another solver to apply the ODEs at each time step?
>
> Thank you very much for your advice.
> Have a nice day
> Marine
>
There's also ode.2D in the same package {deSolve} and it's help page has a 2-d diffusion example that might be cognate.
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
R-help -- Main R Mailing List: Primary help - Homepage - SfS<https://stat.ethz.ch/mailman/listinfo/r-help>
stat.ethz.ch
The main R mailing list, for announcements about the development of R and the availability of new code, questions and answers about problems and solutions using R ...
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
David Winsemius
Alameda, CA, USA
[[alternative HTML version deleted]]
More information about the R-help
mailing list