# Finding the MLE for a univariate exponential Hawkes process

The univariate exponential Hawkes process is a self-exciting point process with an event arrival rate of:

$\lambda(t) = \mu + \sum\limits_{t_i<t}{\alpha e^{-\beta(t-t_i)}}$

where $t_1,..t_n$ are the event arrival times.

The log likelihood function is

$- t_n \mu + \frac{\alpha}{\beta} \sum{( e^{-\beta(t_n-t_i)}-1 )} + \sum\limits_{i<j}{\ln(\mu+\alpha e^{-\beta(t_j-t_i)})}$

which can be calculated recursively:

$- t_n \mu + \frac{\alpha}{\beta} \sum{( e^{-\beta(t_n-t_i)}-1 )} + \sum{\ln(\mu+\alpha R(i))}$

$R(i) = e^{-\beta(t_i-t_{i-1})} (1+R(i-1))$

$R(1) = 0$

What numerical methods can I use to find the MLE? What is the simplest practical method to implement?

• I have had success fitting $\mu$ and $\alpha$ by maximizing the MLE the LBFGS implementation in scipy. The log-likelihood is not concave in $\beta$ though, so I simply iterated over a range of $\beta$ values and picked the one with the maximum likelihood. Note that $\alpha < \beta$ is required for stationarity of the process. – Emaad Ahmed Manzoor Jan 5 '17 at 3:56
• curious, what is the correct form of the λ(t) function using the values of R(i) instead of resumming at each step? – crow Sep 17 '17 at 21:53

## 3 Answers

The Nelder-Mead simplex algorithm seems to work well.. It is implemented in Java by the Apache Commons Math library at https://commons.apache.org/math/ . I've also written a paper about the Hawkes processes at Point Process Models for Multivariate High-Frequency Irregularly Spaced Data .

felix, using exp/log transforms seems to ensure positivity of the parameters. As for the small alpha thing, search the arxiv.org for a paper called "limit theorems for nearly unstable hawkes processes"

• Welcome to the site, @StephenCrowley. If you have your own question, please do not post it as ( / as part of) an answer. Click on the gray "ASK QUESTION" button at the top of the page & ask it there. If you have a question for clarification from the OP, you should ask it in a comment to the question post above. (Although frustratingly, you cannot do that until you reach 50 rep.) – gung Jan 14 '13 at 18:58

I solved this problem using the nlopt library. I found a number of the methods converged quite quickly.

• I assume you're familiar with T. Ozaki (1979), Maximum likelihood estimation of Hawkes' self-exciting point processes, Ann. Inst. Statist. Math., vol. 31, no. 1, 145-155. – cardinal Sep 23 '12 at 15:26
• Could you give more details of what you did? It seems there is a problem with setting constraints and also that large beta is indistinguishable from zero alpha (they both look Poisson). – felix Aug 14 '14 at 16:49

You could also do a simple maximization. In R:

neg.loglik <- function(params, data, opt=TRUE) {
mu <- params
alpha <- params
beta <- params
t <- sort(data)
r <- rep(0,length(t))
for(i in 2:length(t)) {
r[i] <- exp(-beta*(t[i]-t[i-1]))*(1+r[i-1])
}
loglik <- -tail(t,1)*mu
loglik <- loglik+alpha/beta*sum(exp(-beta*(tail(t,1)-t))-1)
loglik <- loglik+sum(log(mu+alpha*r))
if(!opt) {
return(list(negloglik=-loglik, mu=mu, alpha=alpha, beta=beta, t=t,
r=r))
}
else {
return(-loglik)
}
}

# insert your values for (mu, alpha, beta) in par
# insert your times for data
opt <- optim(par=c(1,2,3), fn=neg.loglik, data=data)

• How do you ensure mu, alpha and beta are not set to negative values? – felix Aug 14 '14 at 8:50
• You can set the lower and upper parameters in the optim call. – assumednormal Aug 14 '14 at 16:20
• Not for Nelder-Mead you can't can you which is the default? (See stat.ethz.ch/R-manual/R-devel/library/stats/html/optim.html ) . Also, I don't think there is any way to distinguish huge beta from zero alpha is there so a general optimization seems doomed. – felix Aug 14 '14 at 16:44