[R-sig-dyn-mod] ReacTran

[Karline Soetaert]

Titus, 

I do not understand your first question; your grid goes from 0 to 1., so i do not see any upsurge at 10 km distance 


The plotting can be done as follows:

par(mfrow = c(2,2))
image(out, which = "X", mfrow = NULL, grid = Grid$x.mid)
image(out, which = "Y", mfrow = NULL, grid = Grid$x.mid)
matplot.1D(out, subset = time %in% tt[seq(0,100, by = 10)], type = "l", mfrow = NULL)



or even:


par(mfrow = c(2,2))
image(out, mfrow = NULL, grid = Grid$x.mid)
matplot.1D(out, subset = time %in% tt[seq(0,100, by = 10)], type = "l", mfrow = NULL)


karline

________________________________________
From: r-sig-dynamic-models-bounces at r-project.org <r-sig-dynamic-models-bounces at r-project.org> on behalf of Titus Kassem <gtkassem at gmail.com>
Sent: 24 February 2014 17:20
To: r-sig-dynamic-models at r-project.org
Subject: [R-sig-dyn-mod] ReacTran

I have a pde model whose output I am not  sure is correct.

library(ReacTran)
#library(lattice)
library(rootSolve)
library(shape)
library(deSolve)
N <-100
Grid <- setup.grid.1D(x.up = 0, x.down = 1, N = N)

b = 0.003; delta=0.9; a=0.008; mu=0.02; gamma=0.8;

uini <-rep(x = 0.98, times = N)
vini <- rep(x = 0.02, times = N)
yini <- c(uini, vini)
xx<-seq(0,25,length=N)
Victor1D <- function(t, y, parms) {
  u <- y[1:N]
  v <- y[(N+1):(2*N)]
  du<-b-mu*u-(u*v)/(a+delta*v) +tran.1D (C = u, C.up = 0.98, C.down =
0.98,D = 0.01, dx = Grid)$dC
  dv<-(u*v)/(a+delta*v)-(mu+gamma)*v+tran.1D (C = v, C.up = 0.02, C.down =
0.02,
                                              D = 0.01, dx = Grid)$dC

  list(c(du, dv))
}

tt <- seq(from = 0, to = 25, by =0.1)
print(system.time(
  out <- ode.1D(y = yini, func =Victor1D,
                times = tt, parms = NULL, nspec = NULL,
                names = c("X", "Y"), dimens = N)
))

u=out[,2:(N+1)];
v=out[,(N+1):(2*N)];

plot(xx,v[100,] ,type="l",col='black',lwd=2,xlab="Distance,
Km",ylab="Prevalence",
main="b",cex.main=0.6,cex.lab=0.5,cex.axis=0.5,frame=FALSE)
#text(c(0,0.600),labels=c("t=100 yrs"),cex=0.5)
lines(xx,v[50,] ,type="l",col='black',lwd=2, text(c(0,0.33),labels=c("t=50
yrs"),cex=0.5)
)
lines(xx,v[30,] ,type="l",col='black',lwd=2,text(c(0,0.172),labels=c("t=30
yrs"),cex=0.5))
lines(xx,v[20,] ,type="l",col='black',lwd=2,text(c(0,0.11),labels=c("t=20
yrs"),cex=0.5))
lines(xx,v[10,] ,type="l",col='black',lwd=2,text(c(0,0.057),labels=c("t=10
yrs"),cex=0.5))

# image(out, mfrow = NULL, grid = Grid$x.mi

Two questions really: why is this upsurge at -10km distance? How can I do
three-dimensional plot that incorporates time. distance and prevalence?

With kind regards,

Titus

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-dynamic-models mailing list
R-sig-dynamic-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models