In cmjt/stelfi: Hawkes and Log-Gaussian Cox Point Processes Using Template Model Builder

Spatiotemporal self-exciting models {#stelfi}

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)$

The fit_stelfi() function

args(fit_stelfi)

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)

show_hawkes(list(times = times, params = c(mu = 0.32, alpha = 232, beta = 232)))

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)

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()

cmjt/stelfi documentation built on Oct. 25, 2023, 2:53 p.m.