Learn R Programming

surveillance (version 1.5-4)

algo.glrpois: Poisson regression charts

Description

Poisson regression charts for the monitoring of surveillance time series.

Usage

algo.glrpois(disProgObj,control = list(range=range,c.ARL=5, 
         mu0=NULL, Mtilde=1, M=-1, change="intercept",theta=NULL,
         dir=c("inc","dec"),ret=c("cases","value")))

Arguments

disProgObj
object of class disProg to do surveillance for
control
A list controlling the behaviour of the algorithm [object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object],[object Object]

Value

  • survResalgo.glrpois returns a list of class survRes (surveillance result), which includes the alarm value for recognizing an outbreak (1 for alarm, 0 for no alarm), the threshold value for recognizing the alarm and the input object of class disProg. The upperbound slot of the object are filled with the current $GLR(n)$ value or with the number of cases that are necassary to produce an alarm at any timpoint $<=n$. both="" lead="" to="" the="" same="" alarm="" timepoints,="" but="" "cases" has an obvious interpretation.

encoding

latin1

source

Poisson regression charts for the monitoring of surveillance time series (2006), H�hle{Hoehle}, M., SFB386 Discussion Paper 500.

Details

This function implements the seasonal Poisson charts based on generalized likelihood ratio (GLR) as described in the SFB Discussion Paper 500. A moving-window generalized likelihood ratio detector is used, i.e. the detector has the form$$N = \inf\left{ n : \max_{1\leq k \leq n} \left[ \sum_{t=k}^n \log \left{ \frac{f_{\theta_1}(x_t|z_t)}{f_{\theta_0}(x_t|z_t)} \right} \right] \geq c_\gamma \right}$$where instead of $1\leq k \leq n$ the GLR statistic is computed for all $k \in {n-M, \ldots, n-\tilde{M}+1}$. To achieve the typical behaviour from $1\leq k\leq n$ use Mtilde=1 and M=-1.

So $N$ is the time point where the GLR statistic is above the threshold the first time: An alarm is given and the surveillance is resetted starting from time $N+1$. Note that the same c.ARL as before is used, but if mu0 is different at $N+1,N+2,\ldots$ compared to time $1,2,\ldots$ the run length properties differ. Because c.ARL to obtain a specific ARL can only be obtained my Monte Carlo simulation there is no good way to update c.ARL automatically at the moment. Also, FIR GLR-detectors might be worth considering.

At the moment, window limited ``intercept'' charts have not been extensively tested and are at the moment not supported. As speed is not an issue here this doesn't bother too much. Therefore, a value of M=-1 is always used in the intercept charts.

See Also

algo.rkiLatestTimepoint

Examples

Run this code
##Simulate data and apply the algorithm
S <- 1 ; t <- 1:120 ; m <- length(t)
beta <- c(1.5,0.6,0.6)
omega <- 2*pi/52
#log mu_{0,t}
base <- beta[1] + beta[2] * cos(omega*t) + beta[3] * sin(omega*t) 
#Generate example data with changepoint and tau=tau
tau <- 100
kappa <- 0.4
mu0 <- exp(base)
mu1 <- exp(base  + kappa) 

#Generate data
set.seed(42)
x <- rpois(length(t),mu0*(exp(kappa)^(t>=tau)))
s.ts <- create.disProg(week=1:length(t),observed=x,state=(t>=tau))

#Plot the data
plot(s.ts,legend=NULL,xaxis.years=FALSE)

#Run 
cntrl = list(range=t,c.ARL=5, Mtilde=1, mu0=mu0,
             change="intercept",ret="value",dir="inc")
glr.ts <- algo.glrpois(s.ts,control=c(cntrl))
lr.ts  <- algo.glrpois(s.ts,control=c(cntrl,theta=0.4))

plot(glr.ts,xaxis.years=FALSE)
plot(lr.ts,xaxis.years=FALSE)

Run the code above in your browser using DataLab