[R-sig-dyn-mod] Help with function

Hi Austin,

On Mi, 2015-03-11 at 10:03 -0600, Austin Mullings wrote:

I had one more question, which may be a novice question, but how do I

need

to create my dataset in a way that this function will run off of it?
Specifically, I am trying to run it on some mock data (until I learn more
about simulation), and do not know how R and this package expects my
function or dataset column/variable names to be for it to run. Here is a
picture of my mock data, if you wouldn't mind letting me know whether I
have the variable names ordered incorrectly or need to rename them
accordingly.

I think you cannot send pictures via the mailing list. Unfortunately, I
didn't get your question. What do you mean by "data to run your
simulation on"?

Concerning the sorting of variables: In the example below, I
initialized

yini <- runif(6)

resulting in a vector yini of length 6. The i'th element of yini is y[i]
in myfn (at time point 0). This means that the i'th element is the
initial value of the i'th state, which in turn will be the i'th column
next to the time column of the ode() output.

Cheers,
Daniel

On Tue, Feb 17, 2015 at 2:04 AM, Daniel Kaschek <
daniel.kaschek at physik.uni-freiburg.de> wrote:

Hi Austin,

one possibility to do this is the code below. In this example code, I
randomly initialize your parameters and the initial state. Instead of
using the indexes to access parameter values or states you can use the
names as well (see ?ode, first example). However, in my experience, the
function with() slows down your code considerably.

Best wishes,
Daniel


library(deSolve)

myfn <- function(t, y, p) {

  dy <- numeric(6)
  dy[1] <- p[1] - p[2]*y[4] - p[3]*y[1]
  dy[2] <- y[5]*dy[1]*(1+y[6]) - p[5]*y[2] - (p[6]*y[3]*y[4]*y[1] -

p[7])

  dy[3] <- p[8] - p[9]*y[2]
  dy[4] <- y[6]*dy[1] - p[10]*y[4] + p[11]*(p[6]*y[3]*y[4]*y[1] -

p[7]) +

p[4]
  dy[5] <- 1.0 - y[5]*(y[1] + y[4] + 1)
  dy[6] <- 1.0 - y[6]*(y[1] + 1)

  return(list(dy))

}

pars <- runif(11)
yini <- runif(6)
times <- seq(0, 10, .01)


out <- ode(yini, times, myfn, pars)
plot(out)


On Mo, 2015-02-16 at 13:15 -0700, Austin Mullings wrote:

I am trying to get this equation to work with deSolve, and can't

seem to

get it to work. Does someone know how one would write the function

for

this?

dY1(t)/dt = a - bY3(t) Y4(t) - cY1(t),
dY2(t)/dt = Y5(t) [dY1(t)/dt] [1 + Y6(t)] - eY2(t) - {f Y3(t) Y4(t)

Y1(t) -

g},
dY3(t)/dt = h - iY2(t),
dY4(t)/dt = Y6(t) [dY1(t)/dt] - jY4(t) + k {f Y3(t) Y4(t) Y1(t) - g}

+ d,

dY5(t)/dt = 1.0 - Y5(t) [Y1(t) + Y4(t) + 1],
dY6(t)/dt = 1.0 - Y6(t) [Y1(t) + 1].

      [[alternative HTML version deleted]]

--
Daniel Kaschek
Institute of Physics
Freiburg University

Room 210
Phone: +49 761 2038531

--
Daniel Kaschek
Institute of Physics
Freiburg University

Room 210
Phone: +49 761 2038531