trajectory: Trajectory of a deterministic model
In pomp: Statistical Inference for Partially Observed Markov Processes

trajectory

R Documentation

Trajectory of a deterministic model

Description

Compute trajectories of the deterministic skeleton of a Markov process.

Usage

## S4 method for signature 'missing'
trajectory(
  t0,
  times,
  params,
  skeleton,
  rinit,
  ...,
  ode_control = list(),
  format = c("pomps", "array", "data.frame"),
  verbose = getOption("verbose", FALSE)
)

## S4 method for signature 'data.frame'
trajectory(
  object,
  ...,
  t0,
  times,
  params,
  skeleton,
  rinit,
  ode_control = list(),
  format = c("pomps", "array", "data.frame"),
  verbose = getOption("verbose", FALSE)
)

## S4 method for signature 'pomp'
trajectory(
  object,
  params,
  ...,
  skeleton,
  rinit,
  ode_control = list(),
  format = c("pomps", "array", "data.frame"),
  verbose = getOption("verbose", FALSE)
)

## S4 method for signature 'traj_match_objfun'
trajectory(object, ..., verbose = getOption("verbose", FALSE))

Arguments

`t0`	The zero-time, i.e., the time of the initial state. This must be no later than the time of the first observation, i.e., `t0 <= times[1]`.
`times`	the sequence of observation times. `times` must indicate the column of observation times by name or index. The time vector must be numeric and non-decreasing.
`params`	a named numeric vector or a matrix with rownames containing the parameters at which the simulations are to be performed.
`skeleton`	optional; the deterministic skeleton of the unobserved state process. Depending on whether the model operates in continuous or discrete time, this is either a vectorfield or a map. Accordingly, this is supplied using either the `vectorfield` or `map` fnctions. For more information, see skeleton specification. Setting `skeleton=NULL` removes the deterministic skeleton.
`rinit`	simulator of the initial-state distribution. This can be furnished either as a C snippet, an R function, or the name of a pre-compiled native routine available in a dynamically loaded library. Setting `rinit=NULL` sets the initial-state simulator to its default. For more information, see rinit specification.
`...`	additional arguments are passed to `pomp`.
`ode_control`	optional list; the elements of this list will be passed to `ode` if the skeleton is a vectorfield, and ignored if it is a map.
`format`	the format in which to return the results. `format = "pomps"` causes the trajectories to be returned as a single ‘pomp’ object (if a single parameter vector has been furnished to `trajectory`) or as a ‘pompList’ object (if a matrix of parameters have been furnished). In each of these, the `states` slot will have been replaced by the computed trajectory. Use `states` to view these. `format = "array"` causes the trajectories to be returned in a rank-3 array with dimensions `nvar` x `ncol(params)` x `ntimes`. Here, `nvar` is the number of state variables and `ntimes` the length of the argument `times`. Thus if `x` is the returned array, `x[i,j,k]` is the i-th component of the state vector at time `times[k]` given parameters `params[,j]`. `format = "data.frame"` causes the results to be returned as a single data frame containing the time and states. An ordered factor variable, ‘.id’, distinguishes the trajectories from one another.
`verbose`	logical; if `TRUE`, diagnostic messages will be printed to the console.
`object`	optional; if present, it should be a data frame or a ‘pomp’ object.

Details

In the case of a discrete-time system, the deterministic skeleton is a map and a trajectory is obtained by iterating the map. In the case of a continuous-time system, the deterministic skeleton is a vector-field; trajectory uses the numerical solvers in deSolve to integrate the vectorfield.

Value

The format option controls the nature of the return value of trajectory. See above for details.

Examples


  ## The basic components needed to compute trajectories
  ## of a deterministic dynamical system are
  ## rinit and skeleton.

  ## The following specifies these for a simple continuous-time
  ## model: dx/dt = r (1+e cos(t)) x

  trajectory(
    t0 = 0, times = seq(1,30,by=0.1),
    rinit = function (x0, ...) {
      c(x = x0)
    },
    skeleton = vectorfield(
      function (r, e, t, x, ...) {
        c(x=r*(1+e*cos(t))*x)
      }
    ),
    params = c(r=1,e=3,x0=1)
  ) -> po

  plot(po,log='y')

  ## In the case of a discrete-time skeleton,
  ## we use the 'map' function.  For example,
  ## the following computes a trajectory from
  ## the dynamical system with skeleton
  ## x -> x exp(r sin(omega t)).

  trajectory(
    t0 = 0, times=seq(1,100),
    rinit = function (x0, ...) {
      c(x = x0)
    },
    skeleton = map(
      function (r, t, x, omega, ...) {
        c(x=x*exp(r*sin(omega*t)))
      },
      delta.t=1
    ),
    params = c(r=1,x0=1,omega=4)
  ) -> po

  plot(po)


 # takes too long for R CMD check
  ## generate a bifurcation diagram for the Ricker map
  p <- parmat(coef(ricker()),nrep=500)
  p["r",] <- exp(seq(from=1.5,to=4,length=500))
  trajectory(
    ricker(),
    times=seq(from=1000,to=2000,by=1),
    params=p,
    format="array"
  ) -> x
  matplot(p["r",],x["N",,],pch='.',col='black',
    xlab=expression(log(r)),ylab="N",log='x')

pomp documentation built on Sept. 13, 2024, 1:08 a.m.