fitPP.fun: Fit a non homogeneous Poisson Process In NHPoisson: Modelling and Validation of Non Homogeneous Poisson Processes

Description

This function fits by maximum likelihood a NHPP where the intensity λ(t) is formulated as a function of covariates. It also calculates and plots approximate confidence intervals for λ(t).

Usage

 1 2 3 4 5 fitPP.fun(covariates = NULL, start, fixed=list(), posE = NULL, inddat = NULL, POTob = NULL, nobs = NULL, tind = TRUE, tim = NULL, minfun="nlminb", modCI = "TRUE", CIty = "Transf", clevel = 0.95, tit = "", modSim = "FALSE", dplot = TRUE, xlegend = "topleft", lambdaxlim=NULL,lambdaylim=NULL,...) 

Arguments

 covariates Matrix of the covariates to be included in the linear predictor of the PP intensity (each column is a covariate). It is advisable to give names to the columns of this matrix (using dimnames), since they will be used in the output. Otherwise the default names 'Covariate i' are used. The offset covariates must be included in this matrix. A maximum of 50 covariates are allowed. start Named list of the initial values for the estimation of the β parameters (including fixed parameters). The names of the list must be (compulsory): b0 (for the intercept), b1 (for the first column in covariates), b2 (for the second column), b3 (for the third column), etc. fixed Named list of the fixed β parameters. The elements of this list must be elements of the list start. posE Optional (see Details section). Numeric vector of the position of the PP occurrence points. inddat Optional (see Details section). Index vector equal to 1 for the observations used in the estimation process By default, all the observations are considered. POTob Optional (see Details section). List with elements T and thres that defines the PP resulting from a POT approach; see POTevents.fun for more details. nobs Optional. Number of observations in the observation period; it is only neccessary if POTob, inddat and covariates are NULL. tind Logical flag. If it is TRUE, an independent term is fitted in the linear predictor. It cannot be a character string, so TRUE and not'TRUE' should be used. tim Optional. Time vector of the observation period. By default, a vector 1,...n is considered. minfun Label indicating the function to minimize the negative of the loglikelihood function. There are two possible values: "nlminb" (the default option) and "optim". In the last case, the method of optimization can be chosen with an additional method argument. modCI Logical flag. If it is TRUE, confidence intervals for λ(t) values are calculated. CIty Label indicating the method to calculate the approximate confidence intervals for λ(t). It can be "Transf" for transformed asymptotic intervals (default) or "Delta" for the delta method; see CItran.fun and CIdelta.fun for details. clevel Confidence level of the confidence intervals. tit Character string. A title for the plot. modSim Logical flag. If it is FALSE, information on the estimation process is shown on the screen. For simulation process, the option TRUE should be preferred. dplot Logical flag. If it is TRUE, the fitted intensity is plotted. xlegend Label indicating the position where the legend on the graph will be located. lambdaxlim Optional. Numeric vector of length 2, giving the lowest and highest values which determine the x range. lambdaylim Optional. Numeric vector of length 2, giving the lowest and highest values which determine the y range. ... Further arguments to pass to optim or to nlminb (depending on the value of the minfun argument).

Details

A Poisson process (PP) is usually specified by a vector containing the occurrence points of the process (t_i)_{i=1}^k, (argument posE). Since PP are often used in the framework of POT models, fitPP.fun also provides the possibility of using as input the series of the observed values in a POT model (x_i)_{i=1}^n and the threshold used to define the extreme events (argument POTob).

In the case of PP defined by a POT approach, the observations of the extreme events which are not defined as the occurrence point are not considered in the estimation. This is done through the argument inddat, see POTevents.fun. If the input is provided via argument POTob, index inddat is calculated automatically. See Coles (2001) for more details on the POT approach.

The maximization of the loglikelihood function can be done using two different optimization routines, optim or nlminb, selected in the argument minfun. Depending on the covariates included in the function, one routine can succeed to converge when the other fails.

This function allows us to keep fixed some β parameters (offset terms). This can be used to specify an a priori known component to be included in the linear predictor during fitting. The fixed parameters must be specified in the fixed argument (and also in start); the fixed covariates must be included as columns of covariates.

The estimation of the \hat β covariance matrix is based on the asymptotic distribution of the MLE \hat β, and calculated as the inverse of the negative of the hessian matrix. Confidence intervals for λ(t) can be calculated using two approaches specified in the argument CIty. See Casella (2002) for more details on ML theory and delta method.

Value

An object of class mlePP, which is a subclass of mle. Consequently, many of the generic functions with mle methods, such as logLik or summary, can be applied to the output of this function. Some other generic functions related to fitted models, such as AIC or BIC, can also be applied to mlePP objects.

Note

A homogeneous Poisson process (HPP) can be fitted as a particular case, using an intensity defined by only an intercept and no covariate.

References

Cebrian, A.C., Abaurrea, J. and Asin, J. (2015). NHPoisson: An R Package for Fitting and Validating Nonhomogeneous Poisson Processes. Journal of Statistical Software, 64(6), 1-24.

Coles, S. (2001). An introduction to statistical modelling of extreme values. Springer.

Casella, G. and Berger, R.L., (2002). Statistical inference. Brooks/Cole.

Kutoyants Y.A. (1998).Statistical inference for spatial Poisson processes. Lecture notes in Statistics 134. Springer.

POTevents.fun, globalval.fun, VARbeta.fun, CItran.fun, CIdelta.fun

Examples

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 #model fitted using as input posE and inddat and no confidence intervals data(BarTxTn) covB<-cbind(cos(2*pi*BarTxTn$dia/365), sin(2*pi*BarTxTn$dia/365), BarTxTn$TTx,BarTxTn$Txm31,BarTxTn$Txm31**2) BarEv<-POTevents.fun(T=BarTxTn$Tx,thres=318, date=cbind(BarTxTn$ano,BarTxTn$mes,BarTxTn$dia)) mod1B<-fitPP.fun(covariates=covB, posE=BarEv$Px, inddat=BarEv$inddat, tit="BAR Tx; cos, sin, TTx, Txm31, Txm31**2", start=list(b0=-100,b1=1,b2=-1,b3=0,b4=0,b5=0)) #model fitted using as input a list from POTevents.fun and with confidence intervals tiempoB<-BarTxTn$ano+rep(c(0:152)/153,55) mod2B<-fitPP.fun(covariates=covB, POTob=list(T=BarTxTn$Tx, thres=318), tim=tiempoB, tit="BAR Tx; cos, sin, TTx, Txm31, Txm31**2", start=list(b0=-100,b1=1,b2=-1,b3=0,b4=0,b5=0),CIty="Delta",modCI=TRUE, modSim=TRUE) #model with a fixed parameter (b0) mod1BF<-fitPP.fun(covariates=covB, posE=BarEv$Px, inddat=BarEv\$inddat, tit="BAR Tx; cos, sin, TTx, Txm31, Txm31**2", start=list(b0=-89,b1=1,b2=10,b3=0,b4=0,b5=0), fixed=list(b0=-100)) 

Example output

Loading required package: stats4
Number of events:  137
Number of excesses over threshold 318 : 253

Number of observations  not used in the estimation process:  116
Total number of time observations:  8415
Number of events:  137

Convergence code:  0
Convergence attained
Loglikelihood:  -522.727

Estimated coefficients:
b0      b1      b2      b3      b4      b5
-89.289   2.534   1.425  -0.006   0.557  -0.001
Full coefficients:
b0      b1      b2      b3      b4      b5
-89.289   2.534   1.425  -0.006   0.557  -0.001
attr(,"TypeCoeff")
[1] "Fixed: No  fixed parameters"

Number of observations  not used in the estimation process:  116
Total number of time observations:  8415
Number of events:  137

Convergence code:  0
Convergence attained
Loglikelihood:  -522.889

Estimated coefficients:
b1     b2     b3     b4     b5
2.674  1.488 -0.004  0.630 -0.001
Full coefficients:
b0       b1       b2       b3       b4       b5
-100.000    2.674    1.488   -0.004    0.630   -0.001
attr(,"TypeCoeff")
[1] "Fixed: b0"


NHPoisson documentation built on Feb. 19, 2020, 5:07 p.m.