# Chapter 6 Spatiotemporal self-exciting models

Although the Hawkes process is traditionally formulated as a temporal point process, it is also possible to formulate a spatiotemporal version of the Hawkes process.

For the spatiotemporal Hawkes processes fitted by this package, temporal self-excitement follows an exponential decay function. The self-excitement over space follows a Gaussian distribution centered at the triggering event. There are two formulations of this model. The default is that the Gaussian function has a fixed covariance matrix, independent of time. Alternatively, covariance can be directly proportional to time, meaning that the self-excitement radiates out from the center over time. This can be appropriate when the mechanism causing self-excitement travels at a finite speed, but is very memory-intensive.

The intensity function used by stelfi is

$\lambda(s,t) = \mu + \alpha \Sigma_{i:\tau_i<t}(\text{exp}(-\beta * (t-\tau_i)) G_i(s-x_i, t - \tau_i)).$

Here $$\mu$$ is the background rate, $$\beta$$ is the rate of temporal decay, $$\alpha$$ is the increase in intensity after an event, $$\tau_i$$ are the event times, $$x_i$$ are the event locations (in 2D Euclidean space) and $$G_i(s-x_i, t - \tau_i)$$ is the spatial self-excitement kernel.

For the version with time-independent spatial excitement: $$G_i(s-x_i, t - \tau_i) = f_X(s - x_i)$$ where $$f_X$$ is the density function of $$X \sim \text{N}(0, \Sigma)$$

For the version with time-dependent spatial excitement: $$G_i(s-x_i, t - \tau_i) = f_X(s - x_i)$$ where $$f_X$$ is the density function of $$X \sim \text{N}(0, (t-\tau_i)\Sigma)$$

## 6.1 The fit_stelfi() function

args(fit_stelfi)
## function (times, locs, sf, smesh, parameters, covariates, GMRF = FALSE,
##     time_independent = TRUE, tmb_silent = TRUE, nlminb_silent = TRUE,
##     ...)
## NULL

### 6.1.1 Spatiotemporal Hawkes process

data(xyt, package = "stelfi")
locs <- data.frame(x = xyt$x, y = xyt$y)
times <- xyt$t domain <- sf::st_as_sf(xyt$window)
ggplot(data = locs, aes(x = x, y = y, size = times)) +
geom_point() +
geom_sf(data = domain, fill = NA, inherit.aes = FALSE) +
theme_void() bnd <- INLA::inla.mesh.segment(as.matrix(sf::st_coordinates(domain)[, 1:2]))
smesh <- INLA::inla.mesh.2d(boundary = bnd,
max.edge = 0.75, cutoff = 0.3)
param <- list(mu = 1, alpha = 3, beta = 6, xsigma = 0.2,
ysigma = 0.2, rho = 0)
fit <- fit_stelfi(times = times, locs = locs, sf = domain, smesh = smesh, parameters = param)
get_coefs(fit)
##            Estimate  Std. Error
## mu       0.31682408  0.06432051
## coefs   -1.14940861  0.20301649
## alpha  232.98697937 14.47443211
## beta   232.98698056 14.47443202
## xsigma   0.55456406  0.02336852
## ysigma   0.28344566  0.01440343
## rho      0.02843629  0.05804435
show_hawkes(list(times = times, params = c(mu = 0.32, alpha = 232, beta = 232))) ### 6.1.2 Including a random field

## Note: just because you can, doesn't mean you should...
param <- list(mu = 1, alpha = 20, beta = 200, kappa = 2, tau = 0.1, xsigma = 0.2, ysigma = 0.2, rho = 0)
fit <- fit_stelfi(times = times, locs = locs, sf = domain, smesh = smesh, parameters = param, GMRF = TRUE)
get_coefs(fit)
##             Estimate   Std. Error
## mu       0.003432009  0.004005054
## coefs   -5.674609405  1.166970622
## alpha  259.064098445 15.626830249
## beta   259.064100335 15.626830029
## xsigma   0.548451292  0.021614775
## ysigma   0.296733579  0.014288649
## rho      0.024276592  0.056616125
## kappa    9.040965870 10.519143632
## tau      0.004644277  0.009560608
show_hawkes(list(times = times, params = c(mu = 0.003, alpha = 259, beta = 259))) get_fields(fit, smesh) |>
show_field(smesh = smesh) +
geom_sf(data = mesh_2_sf(smesh), fill = NA) +
geom_sf(data = domain, fill = NA, col = "black", lwd = 2)  +
theme_void() 