Differential equation modeling with Julia

In the last few examples of epidemiological modeling, I focused on stochastic models. However, most published models are in the form of a set of differential equations. In R, I often use simecol. The ODE package in Julia offers some functionality for solving ordinary differential equations, which should be very familiar to users of Matlab/Octave.

Here's the old favourite, the susceptible-infected-recovered model again.

# Load libraries
using ODE
using DataFrames
using Gadfly

# Define the model
# Has to return a column vector
function SIR(t,x)
    S=x[1]
    I=x[2]
    R=x[3]
    beta=x[4]
    gamma=x[5]
    N=S+I
    dS=-beta*S*I/N
    dI=beta*S*I/N-gamma*I
    dR=gamma*I
    return([dS;dI;dR;0;0])
end

# Initialise model
t = linspace(0,500,101);
inits=[9999,1,0,0.1,0.05];

# Run model
result=ode23(SIR,t,inits);

# Collate results in DataFrame
df=DataFrame();
df["t"]=result[1];
df["S"]=result[2][:,1];
df["I"]=result[2][:,2];
df["R"]=result[2][:,3];

# Plot using Gadfly
p=Gadfly.plot(df,x="t",y="I",Geom.line)

There are a few things to note; firstly, the syntax of the ode23 command expects only the time and a state vector. For parameter values, one can hard-code them in the model description (not good) or pass them as additional variables (with gradient 0, so they don't change). I don't really like either option. Perhaps one can use macros to do much the same job, but I'm not that well acquainted with Julia yet.