deSolve package can parameters include a matrix?

Question

I'm trying to code an SEIR model that is age-stratified; that is, in my differential equations I have a parameter for mass action that is the sum of beta*(proportion infected)*(number susceptible) over 20 age classes. The transmission coefficient (beta) is calculated from a contact matrix. The contact matrix has 20 columns and rows which represent the age classes (rows=person i, columns=person j), and contains the probability of contact between two people in any age classes. I designed it and read it into R. My problem is I don't know how (or if) I can use a matrix inside my parameters with deSolve. This code below I wrote doesn't work, I believe because of the matrix / I get this error:

Error in beta * S : non-numeric argument to binary operator

Before I fool with it too much, I'd like to know if it is possible to use a matrix like this as a parameter for this model.

mat <-as.matrix(read.csv("H:/IBS 796R/contactmatrix.csv", header=F))

times <- seq(0,20,by=1/52)
parameters <- c(mu=0,v=1/75,N=1,p=0,delta=2.4,beta=mat*0.04,sigma=1/8,gamma=1/15)
xstart <- c(S=0.06,E=0,I=0.001,R=0)

SEIR0 <- function(t,x,parameters){
    S=x[1]
    E=x[2]
    I=x[3]
    R=x[4]
    with(as.list(parameters), {
        dS=v*S -beta*S*I/N -delta*S
        dE=beta*S*1/N -E*(sigma+delta)
        dI=sigma*E -I*(gamma+delta)
        dR=gamma*I-delta*R
        res=c(dS,dE,dI,dR)
        list(res)
    })
}

out <- as.data.frame(lsoda(xstart,times,SEIR0,parameters))

Also, if I print the parameters, this is what beta looks like:

$beta.V1
 [1] 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04
[12] 4e-04 4e-04 8e-03 8e-03 8e-03 8e-03 8e-03 8e-03 8e-03

$beta.V2
 [1] 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04 4e-04
[12] 4e-04 4e-04 8e-03 8e-03 8e-03 8e-03 8e-03 8e-03 8e-03

....through $beta.V20. So I think it's creating 20 vectors, each with 20 arguments...I think each is a row of the original matrix 'mat' multiplied by the constant 0.04? However, when I multiply mat*0.04 outside "parameters" I get the expected matrix. I'm struggling a bit with how to realize these equations using deSolve and would appreciate any advice on whether it's possible. Thanks in advance.

agstudy agstudy · Accepted Answer · 2013-03-11T06:01:17

The error occurs in this line :

dS=v*S -beta*S*I/N -delta*S

non-numeric argument to binary operator means you try to multiply a function for example by a numeric. You can reproduce it by I*1

Error in I * 1 : non-numeric argument to binary operator`

Here, R can' find beta , and beta is interpreted as the Special Functions of Mathematics, so the error. You need to define parameters as

# a list 
list(mu=0,v=1/75,N=1,p=0,delta=2.4,beta=mat*0.04,sigma=1/8,gamma=1/15)

and

 ## you get a vector mu,N,p,delta,beta1,bet2,...  
c(mu=0,v=1/75,N=1,p=0,delta=2.4,beta=mat*0.04,sigma=1/8,gamma=1/15)

I think you can even rewrite your function as:

SEIR0 <- function(t,x,parameters){
  with(as.list(c(parameters, x)), {
    dS = v*S -beta*S*I/N -delta*S    ## matrix
    dE = beta*S*1/N -E*(sigma+delta) ## matrix
    dI = sigma*E -I*(gamma+delta)
    dR = gamma*I-delta*R
    res = c(dS,dE,dI,dR)
    list(res)                        ## different of the structure of xstart
  })
}

This will correct the problem above but the ODE will not work because the number of derivatives returned by SEIR0 must be equal to the length of the initial conditions xstart vector(4 here).

I suggest for example :

  res <- c(dS=mean(dS),dE=mean(dE),dI=dI,dR=dR)
  list(res)

deSolve package can parameters include a matrix?

1 Answers