Learn R Programming

RHawkes (version 1.0)

mllRH: Minus loglikelihood of a RHawkes model

Description

Calculates the minus loglikelihood of a RHawkes model with given immigration hazard function \(\mu\), offspring density function \(h\) and branching ratio \(\eta\) for event times tms on interval \([0,cens]\).

Usage

mllRH(tms, cens, par, 
      h.fn = function(x, p) dexp(x, rate = 1 / p), 
      mu.fn = function(x, p) {
        exp(dweibull(x, shape = p[1], scale = p[2], log = TRUE) - 
        pweibull(x, shape = p[1], scale = p[2], lower.tail = FALSE, log.p = TRUE))
      }, 
      H.fn = function(x, p) pexp(x, rate = 1 / p), 
      Mu.fn = function(x, p) {
        -pweibull(x, shape = p[1], scale = p[2], lower.tail = FALSE, log.p = TRUE)
      })

Arguments

tms

A numeric vector, with values sorted in ascending order. Event times to fit the RHawkes point process model.

cens

A scalar. The censoring time.

par

A numeric vector containing the parameters of the model, in order of the immigration parameters \(\mu(.)\), offspring parameters \(h(.)\) and lastly the branching ratio \(\eta(.)\).

h.fn

A (vectorized) function. The offspring density function.

mu.fn

A (vectorized) function. The immigration hazard function.

H.fn

A (vectorized) function. Its value at t gives the integral of the offspring density function from 0 to t.

Mu.fn

A (vectorized) function. Its value at t gives the integral of the immigrant hazard function from 0 to t.

Value

The value of the negative log-likelihood.

Examples

Run this code
# NOT RUN {
## earthquake times over 96 years
data(quake);
tms <- sort(quake$time);
# add some random noise to the simultaneous occurring event times
tms[213:214] <- tms[213:214] + 
                    sort(c(runif(1, -1, 1)/(24*60), runif(1, -1, 1)/(24*60)))


## calculate the minus loglikelihood of an RHawkes with some parameters 
## the default hazard function and density functions are Weibull and 
## exponential respectively
mllRH(tms, cens = 96*365.25 , par = c(0.5, 20, 1000, 0.5))
## calculate the MLE for the parameter assuming known parametric forms
## of the immigrant hazard function and offspring density functions.  
system.time(est <- optim(c(0.5, 20, 1000, 0.5), 
                        mllRH, tms = tms, cens = 96*365.25,
                        control = list(maxit = 5000, trace = TRUE),
                        hessian = TRUE)
            )
## point estimate by MLE
est$par
## standard error estimates:
diag(solve(est$hessian))^0.5
# }

Run the code above in your browser using DataLab