pomp (version 3.1)

distributions: Probability distributions

Description

pomp provides a number of probability distributions that have proved useful in modeling partially observed Markov processes. These include the Euler-multinomial family of distributions and the the Gamma white-noise processes.

Usage

reulermultinom(n = 1, size, rate, dt)
deulermultinom(x, size, rate, dt, log = FALSE)
rgammawn(n = 1, sigma, dt)

Arguments

n

integer; number of random variates to generate.

size

scalar integer; number of individuals at risk.

rate

numeric vector of hazard rates.

dt

numeric scalar; duration of Euler step.

x

matrix or vector containing number of individuals that have succumbed to each death process.

log

logical; if TRUE, return logarithm(s) of probabilities.

sigma

numeric scalar; intensity of the Gamma white noise process.

Value

reulermultinom

Returns a length(rate) by n matrix. Each column is a different random draw. Each row contains the numbers of individuals that have succumbed to the corresponding process.

deulermultinom

Returns a vector (of length equal to the number of columns of x) containing the probabilities of observing each column of x given the specified parameters (size, rate, dt).

rgammawn

Returns a vector of length n containing random increments of the integrated Gamma white noise process with intensity sigma.

C API

An interface for C codes using these functions is provided by the package. Visit the package homepage to view the pomp C API document.

Details

If $N$ individuals face constant hazards of death in $k$ ways at rates $r_1,r_2,\dots ,r_k$, then in an interval of duration $\mathrm{\Delta }t$, the number of individuals remaining alive and dying in each way is multinomially distributed: $$(N-\sum _{i=1}^k\mathrm{\Delta }{n}_i,\mathrm{\Delta }{n}_1,\dots ,\mathrm{\Delta }{n}_k)\sim \mathrm{Multinomial}(N;p_0,p_1,\dots ,p_k),$$ where $\mathrm{\Delta }{n}_i$ is the number of individuals dying in way $i$ over the interval, the probability of remaining alive is $p_0=\exp (-\sum _i{r}_i\mathrm{\Delta }t)$, and the probability of dying in way $j$ is $$p_j=\frac{r_j}{\sum _i{r}_i}(1-\exp (-\sum _i{r}_i\mathrm{\Delta }t)).$$ In this case, we say that $$(\mathrm{\Delta }{n}_1,\dots ,\mathrm{\Delta }{n}_k)\sim \mathrm{Eulermultinom}(N,r,\mathrm{\Delta }t),$$ where $r=(r_1,\dots ,r_k)$. Draw $m$ random samples from this distribution by doing

    dn <- reulermultinom(n=m,size=N,rate=r,dt=dt),

where r is the vector of rates. Evaluate the probability that $x=(x_1,\dots ,x_k)$ are the numbers of individuals who have died in each of the $k$ ways over the interval $\mathrm{\Delta }t=$dt, by doing

    deulermultinom(x=x,size=N,rate=r,dt=dt).

Breto & Ionides (2011) discuss how an infinitesimally overdispersed death process can be constructed by compounding a multinomial process with a Gamma white noise process. The Euler approximation of the resulting process can be obtained as follows. Let the increments of the equidispersed process be given by

    reulermultinom(size=N,rate=r,dt=dt).

In this expression, replace the rate $r$ with $r\mathrm{\Delta }W/\mathrm{\Delta }t$, where $\mathrm{\Delta }W\sim \mathrm{Gamma}(\mathrm{\Delta }t/{\sigma }^2,{\sigma }^2)$ is the increment of an integrated Gamma white noise process with intensity $\sigma$. That is, $\mathrm{\Delta }W$ has mean $\mathrm{\Delta }t$ and variance ${\sigma }^2\mathrm{\Delta }t$. The resulting process is overdispersed and converges (as $\mathrm{\Delta }t$ goes to zero) to a well-defined process. The following lines of code accomplish this:

    dW <- rgammawn(sigma=sigma,dt=dt)

    dn <- reulermultinom(size=N,rate=r,dt=dW)

or

    dn <- reulermultinom(size=N,rate=r*dW/dt,dt=dt).

He et al. (2010) use such overdispersed death processes in modeling measles.

For all of the functions described here, access to the underlying C routines is available: see below.

References

C. and E. L. Ionides. Compound Markov counting processe and their applications to modeling infinitesimally over-dispersed systems. Stochastic Processes and their Applications 121, 2571--2591, 2011.

D. He, E.L. Ionides, & A.A. King. Plug-and-play inference for disease dynamics: measles in large and small populations as a case study. Journal of the Royal Society Interface 7, 271--283, 2010.

See Also

Other information on model implementation: Csnippet, accumulators, covariate_table(), dmeasure_spec, dprocess_spec, parameter_trans(), pomp-package, prior_spec, rinit_spec, rmeasure_spec, rprocess_spec, skeleton_spec, transformations, userdata

Examples

# NOT RUN {
print(dn <- reulermultinom(5,size=100,rate=c(a=1,b=2,c=3),dt=0.1))
deulermultinom(x=dn,size=100,rate=c(1,2,3),dt=0.1)
## an Euler-multinomial with overdispersed transitions:
dt <- 0.1
dW <- rgammawn(sigma=0.1,dt=dt)
print(dn <- reulermultinom(5,size=100,rate=c(a=1,b=2,c=3),dt=dW))

# }