Skip to article frontmatterSkip to article content

Point Processes

CSIC
UCM
IGEO

In this section, we look at point processes (PP). In particular, we will use point processes to outline the framework for how we can jointly model extreme event occurrences and magnitudes. We will start with temporal point processes which will address extreme event occurrences. Then we will go into marked temporal point processes which will incorporate magnitudes. finally we will add the spatial component.


Temporal Point Process

These are processes that are concerned with modeling sequences of random events in continuous time. Let’s say we have an ordered sequence of events at time tt. We denote this as

H={tn}n=1NtnTR+\begin{aligned} \mathcal{H} &= \left\{ t_n \right\}_{n=1}^{N} && && t_n\in\mathcal{T}\subseteq \mathbb{R}^+ \end{aligned}

Typically, T=(0,T]\mathcal{T}=(0,T], but it can be between any two arbitrary time endpoints, e.g., T=[t0,t1]\mathcal{T}=[t_0,t_1]. We will also use the notation of the historical events predating our time event of interest, tt. We denote this as

Ht={tntn<t,tnH}\mathcal{H}_t = \left\{t_n|t_n < t,t_n\in \mathcal{H} \right\}

Lastly, we will define the conditional intensity function (aka the hazard function) as

λ(tHt)=limΔt0p(tn[t,t+Δt]Ht)Δt=E[N(Δt)Ht]Δt\boldsymbol{\lambda} (t|\mathcal{H}_t) = \underset{\Delta t\downarrow 0}{\lim} \frac{p(t_n\in [t,t+\Delta t]|\mathcal{H}_t)} {\Delta t} = \frac{\mathbb{E}\left[ N(\Delta t) |\mathcal{H}_t\right]}{\Delta t}

where Δt\Delta t is the infinitesimal time interval containing tt.

We will use the common shorthand to denote the conditional dependence on the historical dataset Ht\mathcal{H}_t.

λ(t)=λ(tHt)\lambda^*(t) = \boldsymbol{\lambda}(t|\mathcal{H}_t)

We can write out the conditional likelihood function as the probability that we observe an event of interest given all of the history as:

p(H)=...p^*(\mathcal{H}) = ...

This can be decomposed as

p(H)=(n=1Nλ(t))exp(0Tλ(τ)dτ) p^*(\mathcal{H}) = \left(\prod_{n=1}^N\lambda^*(t)\right) \exp\left(-\int_0^T\lambda^*(\tau)d\tau\right)

Learning

In general, we are interested in finding the best parameters of our model given access to potentially many sequences of events. So naturally, we can simply maximize the log-likelihood.

L(θ)=arg maxθnDlogp(H;θ)\boldsymbol{L}(\boldsymbol{\theta}) = \underset{\boldsymbol{\theta}}{\argmax} \hspace{2mm} \sum_{n\in\mathcal{D}}\log p(\mathcal{H};\boldsymbol{\theta})

We can write out the joint log-likelihood of observing H\mathcal{H} within a time interval T=(0,T]\mathcal{T} = (0,T] which is given by

logp(H)=n=1Nlogλ(tn)0Tλ(τ)dτ\log p(\mathcal{H}) = \sum_{n=1}^N\log \lambda^*(t_n) - \int_0^T \lambda^*(\tau)d\tau

The first term is the log-likelihood of the specific events at the specific location, tnt_n, that we observe them. The second term is the probability that we do not observe them anywhere else within the time interval of interest. We can also shorten the notation by introducing the cumulative hazard function as

Λ(T)=0Tλ(τ)dτ\Lambda^*(\mathcal{T}) = \int_{0}^T\lambda^*(\tau)d\tau

This will leave us with

logp(H)=n=1Nlogλ(tn)Λ(T)\log p^*(\mathcal{H}) = \sum_{n=1}^N\log \lambda^*(t_n) - \Lambda^*(\mathcal{T})

There are other alternatives to maximizing the log-likelihood. In general, the loss function can look like

L(θ)=arg maxθExp(xnθ)[f(x)]\boldsymbol{L}(\boldsymbol{\theta}) = \underset{\boldsymbol{\theta}}{\argmax} \hspace{2mm} \mathbb{E}_{x\sim p(x_n|{\boldsymbol{\theta}})} \left[ \boldsymbol{f}(x)\right]

where xx is described by some parametric distribution, p(xnθ)p(x_n|{\boldsymbol{\theta}}), and ff is some criteria. In the above example, our criteria is simply the log-likelihood function. We can use some other generative modeling methods like:

  1. GANs use a parametric model and the loss is some sample quality metric
  2. reinforcement learning uses some policy and the loss is some reward function
  3. variational inference uses some approximate posterior and the loss is the evidence lower bound.

Usages

Prediction. The first obvious use case is prediction. In this case, we have some observed data over a period, T(0,T]\mathcal{T}\in(0,T], and we would like to know what will happen in a forecast period, Tτ(T,T+τ]\mathcal{T}_\tau \in (T,T+\tau].

  1. How much time, τn\tau_n, until the next event?
  2. What type of mark, yny_n, of the next event?
  3. How many events of the type, yny_n, will happen?

We can even sample some potential trajectories of events of the future which helps to answer how many events could possibly happen.

100-Year Events. If the marks distribution is parametric, we can do some post-analysis about the occurrence of events. This can be done through the use of return periods (RPs) or average recurrence intervals (ARI). See sections Annual Exceedence Probability and Average Recurrence Interval for more details.


Example I: Homogeneous Poisson Process (HPP)

In this case, we have a dataset of number of exceedances along a timeline as seen in equation (1). Essentially, we have a vector, t\mathbf{t}, which has the counts per unit time. According to the traditional PP, we can define our assumptions as:

  1. The number of events in any two disjoint intervals are independent
  2. The number of events in any interval [a,b][a,b] for 0t0<t1T0\leq t_0 < t_1 \leq T follows a Poisson distribution with rate λ(t1t0)\lambda(t_1-t_0).
  3. The inter-event times are iid rv that follow the exponential distribution with a rate parameter, λ.

Following these assumptions, we say that our intensity function λ(t)\lambda^*(t) be a constant parameter with no dependence on time.

λ(t)=λ\lambda^*(t) = \lambda

This means that our cumulative Hazard function, Λ(T)\Lambda^*(\mathcal{T}), will also not depend on any of the historical events and it will be constant with time. Plugging our terms into the cumulative hazard function in equation (9) results in

Λ(T)=0Tλdτ=(T0)λ=λT\begin{aligned} \Lambda(\mathcal{T}) &= \int_0^T \lambda d\tau = (T-0) \lambda = \lambda T \end{aligned}

where TT is the interval of interest, e.g., number of years, TyearsT_\text{years}. So, we can plug these two quantities into our log likelihood function in equation (10) as

logp(H)=n=1NlogλλT=NlogλλT\begin{aligned} \log p(\mathcal{H}) &= \sum_{n=1}^N\log \lambda - \lambda T \\ &= N\log\lambda - \lambda T \end{aligned}

As mentioned above, the inter-arrival time is an exponential distribution. Please see section ... for more details.


Example II: Inhomogeneous Poisson Process

Also in this case, in this case, we have a dataset of number of exceedances along a timeline as seen in equation (1). However, we let our intensity function λ(t)\lambda^*(t) be a function parameter with dependence on time but no dependence on any historical events.

λ(t)=λ(t)\lambda^*(t) = \lambda(t)

This means that our cumulative Hazard function, Λ(T)\Lambda^*(\mathcal{T}), will also not depend on any of the historical events but it will depend on time.

Λ(T)=0Tλ(τ)dτ=0Tλ(τ)dτ\begin{aligned} \Lambda(\mathcal{T}) &= \int_{0}^T\lambda^*(\tau)d\tau = \int_0^T \lambda(\tau) d\tau \end{aligned}

So, we can plug these two quantities into our log likelihood function into the equation (10)

logp(;θH)=n=1Nlogλ(tn)Λ(T)=n=1Nlogλ(tn)0Tλ(τ)dτ\begin{aligned} \log p(;\boldsymbol{\theta}|\mathcal{H}) &= \sum_{n=1}^N\log \lambda(t_n) - \Lambda^*(\mathcal{T}) \\ &= \sum_{n=1}^N\log \lambda(t_n) - \int_0^T \lambda(\tau) d\tau \end{aligned}

The difficult part for this equation is the 2nd term which is an integral; however, there are many ways to deal with this. For example, we can use a parametric form for the intensity

λ(t)λθ(t)\lambda(t) \approx \lambda_{\boldsymbol{\theta}}(t)

which would result in a closed-form integral. For example, we could use a log-linear model, a cox process or a Hawkes process to name a few. The game is to:

  1. use a simple parametric function that has a closed form integral form
  2. use a more complex parametric function and approximate the integral using quadrature or discretization strategies.

See the section Temporal Dependencies for more ideas of temporal parameterizations.


Marked Temporal Point Process

These are processes that are concerned with modeling sequences of random events in continuous time along with some additional meta-data, i.e., marks. Marks can be whatever type of meta-data we have available. For example, we could have some magnitude, e.g., temperature, Earthquake magnitude. We could also have some spatial information, i.e., latitude, longitude, and/or altitude.

Firstly, we will have some underlying process which is dependent upon time

yn=y(tn)y:R+RDytnTR+\begin{aligned} y_n &= y(t_n) && && y:\mathbb{R}^+\rightarrow \mathbb{R}^{D_y} && && t_n\in\mathcal{T}\subseteq \mathbb{R}^+ \end{aligned}

Let’s say we have a sequence of time stamps, tnt_n, and their associated marks, yy.

This is given as a sequence of events

H={tn,yn}n=1N\begin{aligned} \mathcal{H} &= \left\{ t_n, y_n \right\}_{n=1}^{N} \end{aligned}

We will also use the notation of the historical events predating time, tt.

Ht={(tn,yn)tn<t,tnH}\mathcal{H}_t = \left\{(t_n,y_n)|t_n < t,t_n\in \mathcal{H} \right\}

From a PP perspective, we can model this as a 2D PP which results in

A={[t0,t1]×[y,)}tTR+y0YRA = \left\{ [t_0, t_1]\times[y,\infty) \right\} \hspace{10mm} t\in\mathcal{T}\subseteq\mathbb{R}^+ \hspace{10mm} y_0\in\mathcal{Y}\subseteq\mathbb{R}

Lastly, we will define the conditional intensity function

λ(tHt)=limΔt0p((tn,yn)[t,t+Δt]Ht)Δt\lambda (t|\mathcal{H}_t) = \underset{\Delta t\downarrow 0}{\lim} \frac{p\left((t_n,y_n)\in [t,t+\Delta t]|\mathcal{H}_t\right)} {\Delta t}

We will use the common shorthand to denote the conditional dependence on the historical dataset Ht\mathcal{H}_t.

λ(t,y)=λ(t,yHt)\boldsymbol{\lambda}^*(t,y) = \boldsymbol{\lambda}(t,y|\mathcal{H}_t)

We can write out the joint density as an autoregressive probability where the arrival time, tnt_n, and the mark, yny_n, is conditioned upon the history.

pn(t,y)=p(t,y(t1,y1),(t2,y2),,(tn1,yn1))p_n^*(t,y) = p(t,y|(t_1, y_1), (t_2, y_2), \ldots, (t_{n-1}, y_{n-1}))

We can write out the joint log-likelihood of observing H\mathcal{H} within a time interval T=[0,T]\mathcal{T} = [0,T] which is given by

logp(H)=n=1Nlogλ(tn,yn)+0Tλ(τ)dτ\log p^*(\mathcal{H}) = \sum_{n=1}^N\log \lambda^*(t_n,y_n) + \int_0^T \lambda^*(\tau)d\tau

Marks Parameterizations

Now, let’s dive a bit into the marks distribution. In general, we can model the marks in three ways: 1) conditionally independent, 2) conditioned on time, and 3) time conditioned on the marks. These can be seen in these equations

Conditionally Independent:λ(tn,yn)=λg(tn)p(yn)Temporal Conditioned Marks:λ(tn,yn)=λg(tn)p(yntn)Marks Conditioned Time:λ(tn,yn)=λ(tnyn)\begin{aligned} \text{Conditionally Independent}: && && \lambda^*(t_n,y_n) &= \lambda_g^*(t_n)p^*(y_n) \\ \text{Temporal Conditioned Marks}: && && \lambda^*(t_n,y_n) &= \lambda_g^*(t_n)p^*(y_n|t_n) \\ \text{Marks Conditioned Time}: && && \lambda^*(t_n,y_n) &= \lambda^*(t_n|y_n) \\ \end{aligned}

The first case, we say there is no dependence between the marks. The second case gives us a more flexible parameterization of the marks which influences how the marks behave wrt time. The third cases is the most flexible parameterization and frankly the most correct because we state that the occurrence of events is also conditioned on the marks.

There are some known special cases of these marks. These include:

  1. A compound Poisson process if λg(t)=λ(t)\lambda_g^*(t)=\lambda(t) and f(yt)=f(yt)f^*(y|t)=f(y|t) for deterministic functions λ(t)\lambda(t) and f(yt)f(y|t).
  2. A process with independent marks if λg(t)\lambda_g^*(t) and Hg\mathcal{H}_g-intensity and f(yt)=f(yt)f^*(y|t)=f(y|t)
  3. A process with unpredictable marks if f(yt)=f(yt)f^*(y|t)=f(y|t).

Example III: MPP 4 Extremes

We have a marked point process where we represent it as a 2D point process as shown in equation (22) where we write the joint intensity function for the temporal plane and the mark plane. However, we simplify it to be a parametric form.

λ(t,y)=f(y;θ)\lambda^*(t,y) = \boldsymbol{f}(y;\boldsymbol{\theta})

where f(y;θ)f(y;\boldsymbol{\theta}) is some parametric function in terms of the marks. Now, it’s easier to reason about the cumulative hazard function because it’s an integral of some parametric PDF which has a closed-form double integral.

Λ(A)=0Tyλ(τ,y)dydτ=0Tyf(y;θ)dydτ\Lambda(A) =\int_0^T\int_y^\infty\lambda(\tau,y)dyd\tau =\int_0^T\int_y^\infty \boldsymbol{f}(y;\boldsymbol{\theta})dyd\tau

We recognize that the inner integral for the mark domain is simply the survival function, S\boldsymbol{S}, of the parametric PDF, f\boldsymbol{f}.

y0f(y;θ)dy=10yf(y;θ)dy=1F(y0;θ):=S(y0;θ)\int_{y_0}^\infty\boldsymbol{f}(y;\boldsymbol{\theta})dy = 1 - \int_0^y\boldsymbol{f}(y;\boldsymbol{\theta})dy = 1 - \boldsymbol{F}(y_0;\boldsymbol{\theta}) := \boldsymbol{S}(y_0;\boldsymbol{\theta})

We take the threshold of interest, y0y_0, to be the lower bound of the mark space. So, the remaining outer integral on the temporal domain is a simple homogeneous Poisson process that was done previously in equation (9). Plugging our expression into this equation leaves us with

Λ(A)=0TS(y0;θ)dτ\Lambda(A) = \int_0^T\boldsymbol{S}(y_0;\boldsymbol{\theta})d\tau

Our final log-likelihood expression will be

logp(A)=n=1N(A)logf(y;θ)0TS(y0;θ)dτ\log p^*(A) = \sum_{n=1}^{N(A)} \log \boldsymbol{f}(y;\boldsymbol{\theta}) - \int_0^T \boldsymbol{S}(y_0;\boldsymbol{\theta})d\tau

We can use whatever PDF we want for the marks, e.g., Normal, LogNormal, or T-Student. However, in the literature for extreme values, we typically use the GPD or even the GEVD in some cases.

Annual Exceedence Probabilities. We can also do return periods where try to find the annual exceedence probability or the average recurrence interval. See sections Annual Exceedence Probability and Average Recurrence Interval for more details. For our case, as shown in equations (8) and (14) we equate these to our cumulative distribution function. After solving for the respective TpT_p and TaT_a, we arrive at

y=Q(yp;θ)\begin{aligned} y &= \boldsymbol{Q}(y_p;\boldsymbol{\theta}) \end{aligned}

where Q\boldsymbol{Q} is the quantile function for the PDF/CDF and the probability in the yy domain, ypy_p.

Annual Exceedence Probability:yp=11/TaTa[1,)Average Recurrence Interval:yp=exp(1/Tp),Tp[0,)\begin{aligned} \text{Annual Exceedence Probability}: && && y_p &= 1 - 1 / T_a && && T_a\in[1,\infty)\\ \text{Average Recurrence Interval}: && && y_p &= \exp ( - 1 / T_p), && && T_p\in[0,\infty) \end{aligned}

Decoupled Marked Temporal Point Process

We can decompose this joint intensity measure into its conditional dependencies, i.e., the mark depends on the time.

pn(t,y)=pn(ytn=t)pn(t)p_n^*(t,y) = p_n(y|t_n=t) \cdot p_n^*(t)

The term, pn(ytn=t)p_n^*(y|t_n=t) is either a probability density function or a probability mass function depending upon whether the marks are continuous or discrete. Now, we can write the conditional intensity for the marked TPP as

λ(t,y)=λg(t)f(yt)\boldsymbol{\lambda}^*(t,y) = \lambda_g^*(t) \cdot f^*(y|t)

where λg(t)\lambda_g^*(t) is the ground intensity and f(yt)f^*(y|t) is the conditional mark density function. Notice how the arrival times λg(t)\lambda_g^*(t) are similar to the unmarked case except that now this intensity measure may depend on past marks.

Finally, we can write out the joint log-likelihood of observing H\mathcal{H} within a time interval T=[0,T]\mathcal{T} = [0,T] which is given by

logp(H)=n=1Nlogλg(tn)0Tλg(τ)dτ+n=1Nlogf(yntn)\log p(\mathcal{H}) = \sum_{n=1}^N\log \lambda_g^*(t_n) - \int_0^T \lambda_g^*(\tau)d\tau + \sum_{n=1}^N\log f^*(y_n|t_n)

The first two terms are the ground intensity likelihoods for the temporal rate and the third term is the marks likelihood.


Example IV: MHPP 4 Extremes

In this example, we have a DMTPP for extremes. We have a marked point process where we represent it as a 2D point process as shown in equation (22). In this case, we decouple the intensity function as shown above section. For the ground intensity, we have a HPP. For the marks, we have iid parametric distribution

Ground Intensity:λg(t)=λMarks:f(yt)=f(y,θ)\begin{aligned} \text{Ground Intensity}: && && \lambda^*_g(t) &= \lambda \\ \text{Marks}: && && f^*(y|t) &= \boldsymbol{f}(y,\boldsymbol{\theta}) \end{aligned}

We can plug in these terms into the equation (41) to obtain:

logp(A)=n=1N(A)logλ0Tλdτ+n=1N(A)logf(yn,θ)\log p^*(A) = \sum_{n=1}^{N(A)}\log \lambda - \int_0^T\lambda d\tau + \sum_{n=1}^{N(A)}\log \boldsymbol{f}(y_n,\boldsymbol{\theta})

Because we have the homogeneous rate parameter cases, we get the Marked Homogeneous Poisson Process (MHPP). Using the portion of the HPP, we can plug in the terms found in equation (14) into our equation above. This gives us

logp(A)=NlogλλT+n=1N(A)logf(yntn,θ)\log p^*(A) = N\log \lambda - \lambda T + \sum_{n=1}^{N(A)}\log \boldsymbol{f}(y_n|t_n,\boldsymbol{\theta})

We see that both likelihood terms we decoupled as there are no dependencies between parameters of the temporal likelihood (the first two terms) and the marks likelihood (the third term) so they can be solved independently. In the case of extremes, one option is to use the GPD as the marked distribution. We can write out the new log-likelihood as

logp(A)=NlogλλT+n=1N(A)logfGPD(y;θ)\log p^*(A)= N\log\lambda - \lambda T + \sum_{n=1}^{N(A)} \log \boldsymbol{f}_{\text{GPD}}(y;\boldsymbol{\theta})

This is known as the Poisson-GPD algorithm within the EVT literature Chavez-Demoulin * et al., 2005Nemukula & Sigauke, 2020. We can also parameterize the intensity with the parameterization in equation (65) where we have some new free parameters θ={y0,μ,σ,κ}\boldsymbol{\theta} = \left\{ y_0, \mu, \sigma, \kappa\right\}.

Alternatively, we can use the GEVD as the marked distribution.

logp(A)=NlogλλT+n=1N(A)logfGEVD(y;θ)\log p^*(A)= N\log\lambda - \lambda T + \sum_{n=1}^{N(A)} \log \boldsymbol{f}_{\text{GEVD}}(y;\boldsymbol{\theta})

This has no specific name in the EVT literature, however there are a few papers which use this distribution to motivate the GPD via the PP. In this work, we name this the Poisson-GEVD. We can parameterize GEVD parameters in terms of the GPD parameters which are needed for the intensity parameter given in equations (53) where we introduce free parameters θ={μy0,σy0,κy0}\boldsymbol{\theta} = \left\{ \mu_{y_0}, \sigma_{y_0}, \kappa_{y_0}\right\}.

Annual Exceedence Probability. Lastly, something of great interest within the EVT community is to characterize the return periods. Recall section Average Recurrence Interval whereby we showed that the ARI can be related to the conditional CDF. Recall that this is given as

exp(1/Tˉ)=exp(λS(y;θ))\exp(-1/\bar{T}) = \exp(-\lambda\boldsymbol{S}(y;\boldsymbol{\theta}))

where S\boldsymbol{S} is the survival function for the marks distribution and λ is the average rate of occurrences over the threshold, y0y_0. After rearranging this equation and simplifying, we can recognize that this is simply the quantile function of the marks distribution.

y=Q(yp;θ)\begin{aligned} y &= \boldsymbol{Q}(y_p;\boldsymbol{\theta}) \end{aligned}

where Q\boldsymbol{Q} is the quantile function for the PDF/CDF and the probability in the yy domain, ypy_p.

Annual Exceedence Probability:yp=11/(λTR)TR[1,)Average Recurrence Interval:yp=exp(1/(λTˉ)),Tˉ[0,)\begin{aligned} \text{Annual Exceedence Probability}: && && y_p &= 1 - 1 / (\lambda T_R) && && T_R\in[1,\infty)\\ \text{Average Recurrence Interval}: && && y_p &= \exp \left( - 1 / (\lambda \bar{T})\right), && && \bar{T}\in[0,\infty) \end{aligned}

Example VII: MIPP 4 Extremes

In this example, we have a DMTPP for extremes. We have a marked point process where we represent it as a 2D point process as shown in equation (22). In this case, we decouple the intensity function as shown above section. For the ground intensity, we have a HPP or IPP. For the marks, we have iid parametric distribution

Ground Intensity:λg(t)=λ(t):=λtMarks:f(yt)=f(ytn,θ):=f(yθ(tn))\begin{aligned} \text{Ground Intensity}: && && \lambda^*_g(t) &= \lambda(t) := \lambda_t \\ \text{Marks}: && && f^*(y|t) &= \boldsymbol{f}(y|t_n, \boldsymbol{\theta}) := \boldsymbol{f}(y|\boldsymbol{\theta}(t_n)) \end{aligned}

where the parameters of the marks distribution are time dependent

θt=θ(t)\boldsymbol{\theta}_t = \boldsymbol{\theta}(t)

We can plug in these terms into the equation (41) to obtain:

logp(A)=n=1N(A)logλ(t)0Tλ(τ)dτ+n=1N(A)logf(ynθ(tn))\log p^*(A) = \sum_{n=1}^{N(A)}\log \lambda(t) - \int_0^T\lambda(\tau)d\tau + \sum_{n=1}^{N(A)}\log \boldsymbol{f}\left(y_n|\boldsymbol{\theta}(t_n)\right)

Here, our free parameters of our model are {λt,θt}\left\{ \lambda_t, \theta_t\right\}. Similar to the HPP case, we can use the GEVD or the GPD.

Return Periods. Lastly, something of great interest within the EVT community is to characterize the return periods. Recall section Average Recurrence Interval whereby we showed that the ARI can be related to the conditional CDF. Recall that this is given as

exp(1/Tˉ)=exp(λtS(y;θt))\exp(-1/\bar{T}) = \exp(-\lambda_t\boldsymbol{S}(y;\boldsymbol{\theta}_t))

where S\boldsymbol{S} is the survival function for the marks distribution and λ is the average rate of occurrences over the threshold, y0y_0. After rearranging this equation and simplifying, we can recognize that this is simply the quantile function of the marks distribution.

y=Q(yp;θt)\begin{aligned} y &= \boldsymbol{Q}(y_p;\boldsymbol{\theta}_t) \end{aligned}

where Q\boldsymbol{Q} is the quantile function for the PDF/CDF and the probability in the yy domain, ypy_p.

Return Period:yp=11/(λtTR)TR[1,)Average Recurrence Interval:yp=exp(1/(λtTˉ)),Tˉ[0,)\begin{aligned} \text{Return Period}: && && y_p &= 1 - 1 / (\lambda_t T_R) && && T_R\in[1,\infty)\\ \text{Average Recurrence Interval}: && && y_p &= \exp \left( - 1 / (\lambda_t \bar{T})\right), && && \bar{T}\in[0,\infty) \end{aligned}

Spatial Point Process

These are processes that are concerned with modeling sequences of random events in continuous space and time. Let’s say we have a sequence

H={(tn,sn)}n=1NtnTR+snΩRDs\begin{aligned} \mathcal{H} &= \left\{ (t_n,\mathbf{s}_n) \right\}_{n=1}^N && && t_n\in\mathcal{T}\subseteq\mathbb{R}^+ && && \mathbf{s}_n\in\mathcal{\Omega}\subseteq\mathbb{R}^{D_s} \end{aligned}

We will also use the notation of the historical events predating time, tt.

Ht={(tn,sn)tn<t,tnH}\mathcal{H}_t = \left\{(t_n,\mathbf{s}_n)|t_n < t,t_n\in \mathcal{H} \right\}

Lastly, we will define the conditional intensity function

λ(t,sHt)=limΔt0,Δs0p(tn[t,t+Δt],snΩ(s,Δs)Ht)Ω(s,Δs)Δt\lambda (t,\mathbf{s}|\mathcal{H}_t) = \underset{\Delta t\downarrow 0, \Delta \mathbf{s} \downarrow 0}{\lim} \frac{p(t_n\in [t,t+\Delta t], \mathbf{s}_n \in \Omega(\mathbf{s},\Delta \mathbf{s})|\mathcal{H}_t)} {|\Omega(\mathbf{s},\Delta \mathbf{s})|\Delta t}

We will use the common shorthand to denote the conditional dependence on the historical dataset Ht\mathcal{H}_t.

λ(t,s)=λ(t,sHt)\boldsymbol{\lambda}^*(t,\mathbf{s}) = \boldsymbol{\lambda}(t,\mathbf{s}|\mathcal{H}_t)

Finally, we can write out the joint log-likelihood of observing H\mathcal{H} within a time interval T=[0,T]\mathcal{T} = [0,T] which is given by

logp(H)=n=1Nlogλ(tn,sn)0TΩλ(τ,s)dsdτ\log p(\mathcal{H}) = \sum_{n=1}^N\log \boldsymbol{\lambda}^*(t_n,\mathbf{s}_n) - \int_0^T \int_\mathcal{\Omega}\boldsymbol{\lambda}^*(\tau,\mathbf{s})d\mathbf{s}d\tau

Marked Spatiotemporal Point Process

A marked spatiotemporal processes that are concerned with modeling sequences of random events in continuous space and time which come with some underlying function for the marks. Firstly, we will have some underlying process which is dependent upon time and space

yn=y(sn,tn)y:RDs×R+RDysnΩRDstnTR+\begin{aligned} y_n &= y(\mathbf{s}_n, t_n) && && y:\mathbb{R}^{D_s}\times\mathbb{R}^+\rightarrow \mathbb{R}^{D_y} && && \mathbf{s}_n\in\mathcal{\Omega}\subseteq\mathbb{R}^{D_s} && && t_n\in\mathcal{T}\subseteq \mathbb{R}^+ \end{aligned}

Now, let’s say we have a sequence

H={(tn,sn),yn}n=1N\begin{aligned} \mathcal{H} &= \left\{ (t_n,\mathbf{s}_n), y_n \right\}_{n=1}^N \end{aligned}

We will also use the notation of the historical events predating time, tt.

Ht={(tn,sn,yn)tn<t,tnH}\mathcal{H}_t = \left\{(t_n,\mathbf{s}_n, y_n)|t_n < t,t_n\in \mathcal{H} \right\}

Lastly, we will define the conditional intensity function

λ(t,s,yHt)=limΔt0,Δs0p(tn[t,t+Δt],snΩ(s,Δs),ynYHt)Ω(s,Δs)Δt\lambda (t,\mathbf{s},y|\mathcal{H}_t) = \underset{\Delta t\downarrow 0, \Delta \mathbf{s} \downarrow 0}{\lim} \frac{p(t_n\in [t,t+\Delta t], \mathbf{s}_n \in \Omega(\mathbf{s},\Delta \mathbf{s}), y_n\in \mathcal{Y}|\mathcal{H}_t)} {|\Omega(\mathbf{s},\Delta \mathbf{s})|\Delta t}

We will use the common shorthand to denote the conditional dependence on the historical dataset Ht\mathcal{H}_t.

λ(t,s,y)=λ(t,s,yHt)\boldsymbol{\lambda}^*(t,\mathbf{s},y) = \boldsymbol{\lambda}(t,\mathbf{s},y|\mathcal{H}_t)

Finally, we can write out the joint log-likelihood of observing H\mathcal{H} within a time interval T=[0,T]\mathcal{T} = [0,T] and space interval Ω\mathcal{\Omega} which is given by

logp(H)=n=1Nlogλ(tn,sn)+n=1Nlogf(ynsn,tn)0TΩλ(τ,s)dsdτ\log p(\mathcal{H}) = \sum_{n=1}^N\log \boldsymbol{\lambda}^*(t_n,\mathbf{s}_n) + \sum_{n=1}^N\log \boldsymbol{f}^*(y_n|\mathbf{s}_n, t_n) - \int_0^T \int_\mathcal{\Omega}\boldsymbol{\lambda}^*(\tau,\mathbf{s})d\mathbf{s}d\tau

Literature Review

Applications. Mannshardt-Shamseldin et al., 2010 investigate the differences in precipitation extremes during the 21st century as a trend stemming from global warming. In Cooley & Sain, 2010, they compare the extreme precipitation simulated in a regional climate model over its spatial domain where they apply a Bayesian Hierarchical model for MHPP model.

References
  1. Chavez-Demoulin *, V., Davison, A. C., & McNeil, A. J. (2005). Estimating value-at-risk: a point process approach. Quantitative Finance, 5(2), 227–234. 10.1080/14697680500039613
  2. Nemukula, M. M., & Sigauke, C. (2020). A Point Process Characterisation of Extreme Temperatures: an Application to South African Data. Environmental Modeling & Assessment, 26(2), 163–177. 10.1007/s10666-020-09718-6
  3. Mannshardt-Shamseldin, E. C., Smith, R. L., Sain, S. R., Mearns, L. O., & Cooley, D. (2010). Downscaling extremes: A comparison of extreme value distributions in point-source and gridded precipitation data. The Annals of Applied Statistics, 4(1). 10.1214/09-aoas287
  4. Cooley, D., & Sain, S. R. (2010). Spatial Hierarchical Modeling of Precipitation Extremes From a Regional Climate Model. Journal of Agricultural, Biological, and Environmental Statistics, 15(3), 381–402. 10.1007/s13253-010-0023-9