dpaircase: Expected distributions of distances

In reconhub/vimes: VIsualisation and Monitoring of EpidemicS

Description

These functions compute the expected distributions of distances between pairs of cases given a case reporting probability 'pi'. Analytical results are used for some special cases, including:

Usage

dempiric(p, pi, alpha = 0.001)

dgenetic(x, gamma_shape, gamma_rate = 1, gamma_scale = 1/gamma_rate,
  poisson_rate, pi, alpha = 0.001)

dpaircase(x, type = c("temporal", "genetic", "spatial", "empiric"),
  gamma_shape, gamma_rate = 1, gamma_scale = 1/gamma_rate, poisson_rate,
  sd_spatial, p, pi, alpha = 0.001)

dspatial(x, sd, pi, alpha = 0.001)

dtemporal(x, shape, rate = 1, scale = 1/rate, pi, alpha = 0.001)

Arguments

p	A numeric vector providing the probability mass function, or empirical frequencies, of pairwise distances.
pi	The reporting probability, i.e. the proportion of cases of the outbreak that have been reported.
alpha	The probability threshold to be used to determine the maximum value of generations between two successive cases to consider; this value ('max_kappa') will be the smallest k so that p(k > max_kappa) < alpha. Defaults to 0.001.
x	vector of quantiles.
gamma_shape	shape of the gamma distribution used for the serial interval
gamma_rate	an alternative way to specify the scale of the gamma distribution used for the serial interval
gamma_scale	scale of the gamma distribution used for the serial interval
poisson_rate	rate (i.e. mean) of the poisson distribution used for the per time unit genetic mutation rate
type	type of distance to be considered (one of "temporal","genetic", "spatial" or "empiric").
sd_spatial	standard deviation of the Normal spatial kernel.
sd	standard deviation of the Normal spatial kernel.
shape	shape of the gamma distribution used for the serial interval
rate	an alternative way to specify the scale of the gamma distribution used for the serial interval
scale	scale of the gamma distribution used for the serial interval

Details

• temporal distances: the serial interval is assumed to be gamma-distributed; the procedure returns a discretized, weighted convolution of gamma distributions.

• genetic distances: assumed ... TBC

Author(s)

Anne Cori ([e-mail]) and Thibaut Jombart ([e-mail]).

Examples

## COMPARE DEMPIRIC AND DTEMPORAL

## Note in this comparison we are not expecting to get exactly the same
## results since dempiric does the convolution between discretised gamma
## distributions whilst dtemporal does the convolution between gamma
## distributions.

## compute empirical distribution correponding to exponential(mean 50)

mean_exp <- 50
x <- 0:300
reporting_rate <- 0.5
p <- dgamma(x, shape = mean_exp, rate = 1)

## computes pdf of a gamma distr with shape mean_exp and scale = rate = 1
## (i.e. an exponential distr with mean mean_exp)

## use this as an empirical distribution to feed into dempiric

empiric_exp_distr_r50 <- dpaircase(x, type = "empiric",
                                   p = p,
                                   pi = reporting_rate)
temporal_distr_r50 <- dpaircase(x, type = "temporal",
                                gamma_shape = mean_exp,
                                gamma_rate = 1,
                                pi = reporting_rate)

## compare the two

correlation <- cor(empiric_exp_distr_r50,
                   temporal_distr_r50)


## graphical comparison

plot(x, empiric_exp_distr_r50, xlab = "Time", ylab = "pdf",
     main = "Time between linked cases",
     cex.main = 1, pch = 3)
mtext("SI ~ exp(mean=50), pi = 0.5", side = 3)
lines(x, temporal_distr_r50, col = "red")
legend("topright", c("dempiric","dtemporal"),
       pch = c(3, -1), lwd = c(-1, 1),
       col = c("black","red"), bty = "n")



## COMPARE DEMPIRIC AND DGENETIC

## compute empirical distribution correponding to
## an Exponential(mean 50)-Poisson(mean 0.6) mixture

mean_exp <- 50
mutation_rate <- 0.6
x <- 0:300
reporting_rate <- 0.5
prob <- 1-mutation_rate/(mutation_rate+1)


## pmf of a negative binomial distr with parameters size and prob

p <- dnbinom(x, size = mean_exp,prob = prob)


## use this as an empirical distribution to feed into dempiric

empiric_exp_distr_r50 <- dpaircase(x, type = "empiric",
                                   p = p, pi = reporting_rate)
genetic_distr_r50 <- dpaircase(x, type = "genetic",
                               gamma_shape = mean_exp,
                               gamma_rate = 1,
                               poisson_rate = mutation_rate,
                               pi = reporting_rate)

## compare the two

correlation <- cor(empiric_exp_distr_r50,
                   genetic_distr_r50)


## graphical comparison

plot(x, empiric_exp_distr_r50,
     xlab = "Number of mutations", ylab = "pmf",
     main = "Mutations between linked cases",
     cex.main = 1, pch = 3)
mtext("SI ~ exp(mean=50), pi = 0.6", side = 3)
lines(x, genetic_distr_r50, col = "red")
legend("topright", c("dempiric", "dgenetic"),
       pch = c(3, -1), lwd = c(-1, 1),
       col = c("black","red"), bty = "n")

reconhub/vimes documentation built on May 27, 2019, 4:03 a.m.