Graphical modelling and partial characteristics for multitype and multivariate-marked spatio-temporal point processes

This paper contributes to the multivariate analysis of marked spatio-temporal point process data by introducing different partial point characteristics and extending the spatial dependence graph model formalism. Our approach yields a unified framework for different types of spatio-temporal data including both, purely qualitatively (multivariate) cases and multivariate cases with additional quantitative marks. The proposed graphical model is defined through partial spectral density characteristics, it is highly computationally efficient and reflects the conditional similarity among sets of spatio-temporal sub-processes of either points or marked points with identical discrete marks. The paper considers three applications, two on crime data and a third one on forestry.


Introduction
Spatio-temporal point patterns, where a finite set of pairs of {(s i , t i )} n i=1 with s i ∈ W ⊆ R 2 and t i ∈ T ⊆ R + are the point location and the time of occurrence of the i-th event, respectively, have become ubiquitous in various scientific areas arising in a number of scientific fields, such as infectious disease epidemiology (Gabriel and Diggle, 2009), the study of tornado events (González et al., 2020), fire dynamics (Møller and Díaz-Avalos, 2010) or seismography (Ogata, 1988;Choi and Hall, 1999). In turn, an ever-increasing demand for efficient statistical techniques, which not only account for the spatio-temporal specificity of the data but also facilitate an easy-to-read interpretation, is continuously emerging. Although some progress has been made in the development of spatio-temporal characteristics and models, and different classical point process statistics such as the J-function (Cronie and Van Lieshout, 2015), Ripley's K and pair correlation function (Diggle et al., 1995;Gabriel and Diggle, 2009;Møller and Ghorbani, 2012;Gabriel, 2014), or local indicators of spatio-temporal association functions (Siino et al., 2018) have been extended to the spatio-temporal case, marked spatio-temporal point patterns, where additional qualitative (yielding so-called multitype or multivariate point processes) or quantitative information is available for each pair {(s i , t i )} n i=1 , have not been covered much in the literature. Thus there is an increasing need of efficient statistical techniques allowing for the investigation and analysis of such type of spatio-temporal point processes. For a general review on different spatiotemporal point process statistics and models commonly used at present we refer the interested reader to González et al. (2016). To approach this limitation, the aim of this paper is to contribute to the analysis of multivariate spatio-temporal point processes, where locations and times for a set of different types of points, such as a collection of distinct tree species, is under study. In particular, we consider multivariate marked spatio-temporal point processes where both qualitative and quantitative marks are available for each single pair {(s i , t i )} n i=1 . At present, two different approaches can be identified in the literature focussing on quantitatively marked spatio-temporal point processes where one, potentially timevarying, real-valued mark is attached to each single point location. One strand of the literature, mainly applied in the field of spatio-temporal earthquake research, was covered by Rathbun (1993), Ogata et al. (2003), Choi and Hall (1999) and Marsan and Lengliné (2008)) amongst others. They use mechanistic models (Diggle, 2013) which are defined through a parametric conditional intensity framework and express the infinitesimal expected rate of events at a particular time t i at point location s i conditional on the history of the complete spatio-temporal process up to time t i (see Vere-Jones (2009)). While the conditional intensity formalism provides a complete description of the process, the conditional intensity itself may be intractable or can not be evaluated exactly without numerical methods. Different from the mechanistic modelling approach, Särkkä and Renshaw (2006), Renshaw and Comas (2008), Comas et al. (2009), Renshaw et al. (2009), Cronie and Särkkä (2011), Cronie et al. (2012), Comas et al. (2013), and Redenbach and Särkkä (2013) discussed so-called growthinteraction processes to model the evolution of a quantitatively marked spatial point pattern over equidistant steps in time such as the diameter at breast height (DBH) value for a set of trees recorded over consecutive times. While the above specifications consider the analysis of quantitatively marked spatiotemporal point patterns, the investigation of cross-characteristics through marked ver-sions of the spatio-temporal reduced second-order moment measure and Ripley's Kfunction has just recently been discussed by Iftimi et al. (2019). Unlike the classical second-order summary characteristics such as the spatio-temporal K-function, the corresponding marked version allows to investigate the pairwise interrelation between subsets of points within one quantitatively marked spatio-temporal point pattern, e.g. the pairwise distance between juvenile and adult trees classified subject to a given threshold computed from the quantitative mark itself. That is, the marked spatio-temporal version of Ripley's K-function describes the expected number of further space-time points of type j from an arbitrary space-time point of type i of the process, given that the points in question have space and time separation r ≥ 0 and t ≥ 0, respectively. Besides these marked spatio-temporal point process characteristics, these authors also shortly pointed to extensions of classical cross-and dot-type point process characteristics for the multivariate case over equidistant steps in time in the supplementary material of their paper. Originating in the purely spatial case, these two types of point process characteristics for qualitatively marked processes investigate the pairwise distances between the point locations of two distinct component patterns or between the point locations of one component and those of any alternative patterns. Different from the above approaches, the interest of the present paper are pairwise as well as global structural interrelations between different spatio-temporal components conditional on all remaining components defined in terms of partial spatio-temporal point process characteristics and the spatio-temporal dependence graph model, respectively. While spatio-temporal point process characteristics quantify the conditional interrelation between two distinct components given all remaining components of the spatio-temporal process, the spatio-temporal dependence graph model simultaneously elicits the potential structural interrelations between all spatio-temporal components in form of an undirected graph, and thus allows to detect directed and also induced spatiotemporal interdependencies. That is, our focus is the extension of the more classical concepts of partial correlation into the field of multivariate and multivariate-marked spatio-temporal point processes. Building upon the results of Eckardt (2016), Eckardt and Mateu (2019a) and Eckardt and Mateu (2019b), this paper introduces different partial point process characteristics in the frequency domain. In addition, adopting the ideas of classical multitype point process characteristics, a new dot-type spectra, the dot-spectra, is introduced which reflects the linear interrelation of one component and any alternative patterns included. To the best of our knowledge, the treatment of a combination of discrete with quantitative marks in a context of spatio-temporal point processes is new. If, in addition, we consider partial characteristics we go a step further with respect to the existing literature. Finally, the extension of a spatial dependence graph model to the spatio-temporal context is also new. The remainder of the paper is structured as follows. Section 2 provides some background on the main characteristics of point processes in the spatio-temporal domain. Section 3 develops the main results of the spectral analysis for spatio-temporal point processes. Then, Section 4 presents the spatio-temporal dependence graph model. Applications to crime data and forestry are developed in Sections 5 and 6. The paper ends with some final conclusions.

Spatio-temporal point processes
To introduce spatio-temporal point processes, we follow González et al. (2016) and references therein as well as Iftimi et al. (2019). A spatio-temporal point process X is, rigorously speaking, a random countable measure defined in a subset +∞). A realisation of a point process is called a point pattern and it is a finite set of pairs where the components are intended to state the spatial location s i ∈ W , and the time associated with that spatial location t i ∈ T . Let N (A × B) denote the number of points of the set (A × B) ∩ X. Stationarity and isotropy for spatio-temporal point processes can be defined as follows. X is called (spatio-temporally) stationary when the process (s, t) + X keeps the distribution of the original process X. On the other hand, X is (spatially) isotropic if, for any rotation r around the origin, the rotated point process rX = {(rs, t) : (s, t) ∈ X} keeps the distribution of X. Given a finite point process, it is frequently convenient to project it onto the spatial and temporal windows, so we can treat separately space and time (Møller and Ghorbani, 2012), For the marked case consider the process {(s i , t i ), m i } n i=1 where m i is a mark in a suitable mark space M. It is then called a spatio-temporal marked point pattern. When M = {1, 2, . . . , k}, k ≥ 2, the process is called a multitype spatio-temporal point process. The associated spatio-temporal point process is called the ground process, and it is denoted by X g .

Spatio-temporal point process descriptors
In analogy with the classical theory of random variables, we would like to deal with the distribution of the points of X in W × T × M where M ⊂ M. The product densities λ (k) , k ≥ 1 describe the probability that there is a point of the process in each of the pairwise disjoints balls with centres in k given points ξ 1 , . . . , ξ k and infinitesimal spatiotemporal marked volumes dξ 1 , . . . , dξ k ⊆ W ×T ×M with dξ i = ds i ×dt i ×ν(dm i ), where ν() is a bounded reference measure on the mark space, and size |dξ i | = ds i dt i ν(dm i ), i = 1, . . . , k. They can be defined by using Campbell's formula, which states that given a marked spatio-temporal point process X, for any non-negative function h on where = indicates that the summation is over distinct k-tuples of marked spatiotemporal events.

Intensity function
As a first particular case of Eq.
(1), we focus on the intensity measure and intensity function. Usually, the analysis of a spatio-temporal point pattern starts with the estimation and modelling of the intensity function as it rules the univariate distribution of X in W × T × M . Considering the so-called intensity measure given by and we call λ(s, t, m) the first-order intensity function of X. Consider the projection of the process X to only its spatio-temporal coordinates, the resulting process is called the ground process and it is denoted by X g . It can be shown that the intensity satisfies where λ g (s, t) is the intensity of the ground process and f (m) is a conditional density on M in the spatio-temporal location (s, t). In case that X g is stationary, or equivalently homogeneous, then λ g (s, t) ≡ λ > 0. This constant is called the intensity of the ground process and X is said spatio-temporally homogeneous.
The first-order intensity of the ground process can be defined as well as When the first-order intensity function of the ground process λ g (s, t) can be factorised as λ g (s, t) = λ 1 (s)λ 2 (t), where λ 1 (·) and λ 2 (·) are non-negative functions, then the process is called first-order spatio-temporal separable. This separability is often taken as a working assumption in the literature in order to facilitate the estimations. In that case, the effects that are nonseparable could be interpreted as second-order effects. Note that a stationary spatiotemporal point process X is trivially first-order separable as its intensity is constant. Once the sets X space and X time have been defined, it is naturally possible to define the marginal spatial and temporal intensity functions λ space and λ time as so that λ g (s, t) ∝ λ space (s)λ time (t), with λ g , λ space , λ time all being constant when X is homogeneous.
To estimate the spatio-temporal first-order intensity function of the ground process X g , the estimation of the marginal spatial and temporal intensity functions is first presented. For the spatial intensity function, a non-parametric specification which is most frequently used at present is defined in the form of a kernel estimator where k (s) = 1 2 k s , k(·) is a bivariate kernel and > 0 is the bandwidth, and is an edge-correction intended to stabilise the mass of the estimator so that its integral is roughly the number of points n. The marginal temporal intensity function λ time (t) can be estimated in the very same non-parametric fashion. The bandwidth is a sensitive parameter extremely delicate to be chosen; however, there are several methods to approach to a proper value. We note that there are some alternatives to estimate the spatial and temporal intensity components by using parametric or semi-parametric methods. The suitability of these approaches depends on how well we know the data if there are helpful covariates. We note that under separability, given two unbiased estimatorsλ space (·) andλ time (·), an unbiased estimator of the spatio-temporal first-order intensity function of the ground process X g isλ When the spatio-temporal separability is not fulfilled there are some options to properly estimate the intensity, for instance, a non-separable estimator is given by González et al. (2020) where k 1 δ is a one-dimensional Gaussian kernel with bandwidth δ and c δ (v; T ) is the analogous to Diggle's edge-correction for the temporal component.
An important summary statistic for marked spatio-temporal point process is the pair correlation function. This can be defined as the standardised version (and far more useful) of the product density function, where g g is the pair correlation function of the ground process and it is given by g g ((s 1 , t 1 ), (s 2 , t 2 )) = λ g ((s 1 , t 1 ), (s 2 , t 2 )) λ g (s 1 , t 1 )λ g (s 2 , t 2 ) , (s 1 , t 1 ), (s 2 , t 2 ) ∈ W × T.
The advantage of having a standardised product density function is that this function takes the constant value 1 for a spatio-temporal complete random process in the presence of independent marking. So values above or below this benchmark will be easily interpreted towards clustering or regularity. One of the most important working assumptions when dealing with marked spatiotemporal point processes is the concept of second-order intensity-reweighted stationarity defined as follows. A marked spatio-temporal point process X is second-order intensityreweighted stationary (SOIRS) (Gabriel and Diggle, 2009) if g((s 1 , t 1 , m 1 ), (s 2 , t 2 , m 2 )) =ḡ(s 1 − s 2 , t 1 − t 2 , m 1 , m 2 ), for any (s 1 , t 1 ), (s 2 , t 2 ) ∈ W × T , whereḡ is a non-negative function.
If the process is also isotropic,ḡ(s 1 − s 2 , t 1 − t 2 , m 1 , m 2 ) = g 0 (r, t, m 1 , m 2 ), meaning that the pair correlation depends only on the distances r = s 1 − s 2 and t = |t 1 − t 2 |, where g 0 is a non-negative function. One of the most common methods for the estimation of the pair correlation function is the non-parametric kernel approach since such estimator is easy to interpret and implement. Assuming that the spatio-temporal point pattern is given by a sequence of , the estimator is given bŷ where k 1 and k 2δ are kernel functions with spatial and temporal bandwidths and δ, and w ij is an edge-correction factor for correcting the lack of information occurring between points close to the edge of W × T and the unobserved outsider points (see e.g, Gabriel, 2014). A spatio-temporal adaptation of the mean product of marks sited a distance r apart (see e.g, Renshaw, 2002) can be thought of as a natural extension by including the temporal dimension, i.e, for a stationary and isotropic process U (r, t) = λ 2 g g (r, t)S(r, t)ds 1 dt 1 ds 2 dt 2 , where ds 1 and ds 2 are two infinitesimal spatial areas separated by a distance r and dt 1 and dt 2 are two infinitesimal temporal lengths separated by a distance t. S(r, t) represents a spatio-temporal mark correlation function that has not yet been examined in the current literature and that deserves especial attention given its extremely usefulness for analysing complex point patterns. Finally, the marked spatio-temporal K-function was defined in Iftimi et al. (2019) in its general version. We can take advantage of the pair correlation function in the case of SOIRS processes. Let C, D ⊂ M, so the K-function is given by For Poisson processes, the K-function is 2πr 2 t. This statistic can be estimated through the following expression where 1(·) is the indicator function, and w ij is a suitable edge-correction.

Multitype spatio-temporal point patterns
In analogy with the classical theory of multitype point patterns, we can define some useful descriptors. Consider a multitype spatio-temporal point process composed by q types of points, so that X = q i=1 X (i) . It is straightforward to consider a multitype point process as a marked point process where M is the set of indices {1, . . . , q}; thus, the measure ν(·) is the counting measure. Let N i (A × B) denote the number of points of type i of the set (A × B) ∩ X. We can naturally define the spatio-temporal cross-product density function as Note that λ (2) ii ≡ λ (2) in Eq. (4). The spatio-temporal cross-covariance density function for, say, type i-points can be defined by Having the first-and second-order characteristics at hand, two alternative statistics can be defined: (a) the spatio-temporal correlation function Cor((s, t), (s , t )) = (dsdtds dt ) 1/2 ζ ij ((s, t), (s , t )) (ζ i ((s, t), (s , t ))ζ j ((s, t), (s , t ))) 1/2 and (b) the scaled cross-covariance density function All these functions are well defined under some regularity conditions, indeed multiple coincident events are precluded. This assumption implies that which leads to a natural definition of Bartlett's complete covariance density function as where δ(·) is a multivariate Dirac delta function with We note that this function simplifies under second-order stationarity to κ ii (c, h) = λ . Generalising Bartlett's complete covariance density function to the spatio-temporal case, we assume that the spatio-temporal cross-covariance and the complete crosscovariance density functions coincide, such that

Spectral analysis of spatio-temporal point processes
Extending Dorai-Raj (2001) to the present context, and requiring orderliness and discrete equidistant points in time, we define the auto-and cross-spectral density functions of the i-th and j-th components of a multivariate second-order stationary spatiotemporal point process as the Fourier transform of the complete spatio-temporal autoand cross-covariance density functions at frequencies w st = (w s , w u ). By this, the auto-spectrum of the i-th component is ii (w u ) is the convolution of the auto-spectral density functions for the spatial and temporal components, w s = (w p , w q ) is a two-dimensional array of spatial frequencies and w u is a vector of temporal frequencies. Applying Bochner's theorem, κ ii (c, k) can be recovered by the inverse Fourier transform of (8), Further, under second-order spatio-temporal separability, (8) simplifies and allows for the decomposition into which yields, by inverse Fourier operation, again the complete spatio-temporal autocovariance density function in the form of Likewise, the spatio-temporal cross-spectral density function ij (w u ), which measures the linear interrelation of the spatio-temporal components N i and N j , is defined by As such that it suffices to consider only one cross-spectrum. However, at the same time, as the spatio-temporal cross-covariance density function is not necessarily symmetric is a complex-valued function which can be decomposed into its real and imaginary parts either in terms of Cartesian or polar coordinates yielding the co-spectrum C ij (w s ), the quadrature spectrum Q ij (w s , w u ), the amplitude spectrum a ij (w s , w u ) = mod (f ij (w s , w u )) and the phase spectrum ). The amplitude spectrum represents the relative magnitude of the power attributable to frequencies (w s , w u ) while the phase spectrum indicates the similarity of two patterns up to linear shifts (cf. Chatfield (1989); Priestley (1981)).
Although the spatio-temporal cross-spectrum provides insights into the linear interrelation of two components at frequencies w st , it is preferable to compute the spatiotemporal spectral coherence function, satisfying 0 ≤ |R ij (w st )| 2 ≤ 1. For a subset of J components, a different spectral coherence function which quantifies the extend to which the i-th component is determinable from the J components is the multiple coherence function |R m iJ (w st )| 2 , Defining κ i• (c, h) as the complete dot-type cross-covariance function between the i-th and all N \ {i}-th components and substituting κ i• (c, h) for κ ij (c, h) in (9) yields the dot-type spatio-temporal cross-spectrum from which, in turn, the dot-type spatio-temporal spectral coherence function |R i• (w s , w u )| 2 and multiple coherence function |R m iJ (w st )| 2 can be computed. We note that |R m iJ (w st )| 2 and |R i• (w st )| 2 coincide whenever J equals N \ {i}. Unlike the ordinary cross-spectral characteristics, the above dot-type functions express the linear interrelation between one particular component and the set of all remaining components.
Besides the spectral coherence functions, a different spectral quantity which measures the linear effect of the j-th on i-th (resp. i-th on j-th) component is the spatio-temporal gain spectrum G i|j (w s , w u ) (resp. G j|i (w s , w u )) defined by where G j|i (w s , w u ) can be computed analogous to G i|j (w s , w u ). Defining a dot-type version of the above functions yields the expression where f •• is the auto-spectrum defined over all components except i. This function provides information on the linear effect of all the alternative components on the i-th component.
Analogous to (8), the marked spatio-temporal auto-spectrum for the i-th component of a multivariate-marked point process is defined as the Fourier transform of the auto-type spatio-temporal mean product of marks U ii , and is given by The cross-term expression of (13) is obtained in the same way through the Fourier transform of U ij . Recapitulating Bochner's theorem, we note that again both quantities U ii and U ij can uniquely be recovered by the inverse Fourier operations. As for the multitype case, we note that the marked cross-spectrum could be extended to a dottype version by substituting U ij by the dot-type spatio-temporal mean product of marks U i• . This statistic would then include dot-type versions of the spatio-temporal pair and mark correlation functions which need further rigorous investigations in future research. Adopting the ideas of Renshaw andFord (1983, 1984), two different spectral representations of the spatio-temporal (marked) spectra can be defined. These spectra, the R-and the Θ-spectrum at the temporal frequency w u can be computed directly by converting the spatial frequency w s into polar form w θ with = p 2 + q 2 and θ = tan −1 (p/q). This yields the R t -spectrumf R ( , w u ) = 1 n θf . . , 170 • where n θ is the number of periodogram ordinates for which θ − 5 • < θ ≤ θ + 5 • , respectively. Herê f • denotes the polar form computed from the ordinary (marked) spatio-temporal crossspectra at temporal frequency w u . While the R-spectrum provides useful information on the scales of point patterns under the assumption of isotropy, the Θ-spectrum could be used to investigate directional features of the point pattern.

Estimation of spatio-temporal spectral density functions
Next, the estimation of both functions from empirical data by means of spatio-temporal auto-and cross-periodograms is presented. Assume that we have a multivariate spatiotemporal point pattern in a bounded region W × T ⊂ R 2 × R + where W is required to be a rectangular region with sides of length l 1 and l 2 and independent of t. Let {s i (t)} = {(x ik (t); y ik (t))} with i = 1, . . . , n i (resp {s j (t)}) denote the locations of points of type i (resp. of type j) recorded at time t, t ∈ T . To ease notation, we will be using only the subindex k for the coordinates, whenever no confusion arises. For simplicity, we assume that the locations have been scaled to the unit square prior to analysis and that T is an ordered set of consecutive times recorded for equidistant steps in discrete time. The spatio-temporal auto-and cross-periodograms result from the DFT of the point locations {s i (t)} and {s j (t)} at time t ∈ T . For events of type i, the DFT is defined as where a i (p, q, u) and b i (p, q, u) are the real and the imaginary parts of F i (p, q, u), p = 0, 1, 2, . . . , q = ±1, ±2, . . . and u = − T −1 2 , . . . , T 2 . In general, p and q are assumed to be independent of u. Under second-order separability, (14) factorises to where F (t) i (p, q) is the Fourier transform of the spatial frequencies (w p , w q ) for events of type i at time t. From this expression, the spatio-temporal auto-periodogram for frequencies w s = (2πp/n i , 2πq/n i ) and w u = 2πu itself follows as where h = t − t is the time lag and F i (·). The computation of the spatio-temporal cross-periodogram follows analogously to (16) leading to where p, q and u are defined as above.
Likewise, in the presence of both one qualitative and one quantitative marks for each point location, both the marked spatio-temporal auto-and the marked spatio-temporal cross-periodograms result from the discrete Fourier transforms of the marked locations {s i (t), m i (s i (t))} and {s j (t), m j (s j (t))} at time t ∈ T . Assuming that the marked locations have been scaled to the unit square prior to the analysis, the spatial component where µ(m k (x k (t), y k (t)) is the mean computed over all quantitative marks for the i-th component. Plugging-in this expression for F , y k (t))) × exp(−2πı(px k (t) + qy k (t))) (18) where p, q and u are defined as above.

Spatio-temporal dependence graph model
Although some progress has been made in the analysis of marked spatio-temporal point patterns and various point process characteristics can be found in the literature, the need for efficient exploratory techniques for multivariate and multivariatemarked spatio-temporal point patterns which allow for the simultaneous investigation of potential conditional interrelations among all component patterns still remains. To overcome this limitation, this section extends the framework of the spatial dependence graph model to introduce a new class of spatio-temporal dependence graph models which allows for the joint analysis of potential direct and indirect interrelations in multivariate and multivariate-marked spatio-temporal point patterns. In addition, using the same ideas underpinning the graphical model, different partial spatio-temporal point process characteristics are introduced which represent the pair interrelation between to (marked) components conditional on all remaining patterns. To put it differently, we are interested in the partial linear interrelations between two component processes which remain conditional on all alternative components as expressed by partial spatio-temporal spectral characteristics, i.e. the partial spatio-temporal spectral density function f ij|V\{i,j} (w st ), which are next presented.

Partial spatio-temporal spectral properties
To formalise the concept of partial spectral properties, let N V denote a d-variate spatiotemporal point patterns indexed by V = 1, . . . , d where d ≥ 3 and N V\{i,j} denote all alternative components of N V except N i and N j . Different from the ordinary spectral properties which do not help to distinguish between direct and induced interrelations, the objective of interest of this section are linear interrelations between any pair of distinct components (N i , N j ) conditional on N V\{i,j} . That is, the pairwise linear interrelation of N i and N j which remains after the linear effect of all alternative components has been removed. In this respect, the partial cross-spectrum f ij|V\{i,j} (w st ) can be regarded as the cross-spectrum of two residual processes i and j computed from N i and N j . Analogously to (10), rescaling of the partial cross-spectral density function yields the partial spectral coherence function |R ij|V\{i,j} (w st )| 2 , which is also bounded between zero and one. However, different from the ordinary spectral coherence function, this function expresses the linear interrelation of two component processes which remains after the linear effect of all remaining component processes has been removed by orthogonal projection. In this sense, the partial spectral coherence can be understood as the partial correlation defined as a function of frequencies w st such that N i and N j are conditionally independent at all spatial and temporal lags given N V\{i,j} (N i ⊥ ⊥ N j | N V\{i,j} ) if |R ij|V\{i,j} (w) st | 2 vanishes at all frequencies w (cf. Brillinger (1981); Rosenberg et al. (1989)). Having only three distinct components i, j and k under study, the calculation can be simplified though Besides, an alternative spectral characteristic called the absolute rescaled inverse spectral density function |d ij (w) st | can be calculated from the negative of the partial spectral coherency function, that is |d ij (w st )| = −R ij|V\{i,j} (w st ) which measures the strength of the linear partial interrelation between N i and N j at frequencies w (cf. Dahlhaus (2000)).

Estimation of partial spectral spatio-temporal density functions
Whilst Section 4.1 formalises the concept of partial point process spectral properties, we now briefly review different computational methods for the calculation of the partial cross-spectrum from empirical data. Applying well-known results from the theory of the multivariate normal and Brillinger (1981, Theorem 8.3.1.), a first approach which requires the inversion of a (d−2)×(d−2) matrix is to compute the partial cross-spectrum using the formula where is a (d − 2) × 1 matrix. As an alternative approach, a step-wise procedure for (20) can be implemented by recursively applying algebraic operations as described by Bendat (1978). However, as the required number of stepwise calculations strictly depends on the number of components under study, this recursive calculus is computationally inefficient in high dimensional settings. Finally, a less computationally intensive approach has been introduced in Dahlhaus (2000) and effectively extended to the spatial domain by e.g. Eckardt (2016) where, under regularity assumptions, the partial spectra can be obtained from the inverse (w st ) of the spectral matrix f (w st ) such that We note that expressions for the multivariate-marked case can analogously be computed from the inverse of the marked cross-spectra m ij (w st ). Likewise, for subsets X i , X J and X K where X J ∪ X K ⊂ X V\{i} , replacing ij(w st ) by iK (w st ), the inverse of f iK (w st ), yields the partial dot-type spectra f iK|K (w st ). Different from ordinary partial spectral characteristics, f iK|J (w st ) describes the partial interrelation between component X i and subset X K conditional on subset X J which is identical to f ik|V\{i,j} if and only if X K reduces exactly to component X k and X J = X V\{i,k} .

Spatial dependence graph model for spatial point processes
Adopting the ideas of Eckardt (2016) and Eckardt and Mateu (2019a), we now define the spatio-temporal dependence graph model (henceforth STDGM) aiming to relate the structure of an undirected graph to the partial interrelation structure of a multivariate and multivariate-marked spatio-temporal point pattern. To this end, we identify the vertices of an undirected graph with the components of any such spatio-temporal point process such that an edge between the two vertices v i and v j is missing if and only if the component processes N i and N j are conditionally uncorrelated after removal of the linear effect of N V\{i,j} , e.g. if both components are homogeneous spatio-temporal Poisson processes conditional on all remaining components. This assumption is equivalent to observing a vanishing partial cross-spectrum f ij|V\{i,j} (w st ), inverse ij (w st ), partial spectral coherence function R ij|V\{i,j} (w st ) or absolute rescaled inverse spectral density function d ij (w st ) at all frequencies w st in space and time. This leads to the following definition of a STDGM. Let N V be a multivariate or multivariate-marked spatio-temporal point process on Hence, a spatial dependence graph model is an undirected graph in which conditional interrelations can be identified from non-missing edges. Precisely, two components N i and N j are said to be conditionally uncorrelated at all spatial and all temporal lags after removing the linear effect of all remaining components if the unordered pair {v i , v j }, i = j is not in E.

Partial characteristics in the spatio-temporal domain
Adopting Bochner's theorem, different partial spatio-temporal domain characteristics can be computed directly from the partial spectral characteristics through the inverse Fourier transformation. Whence, for the qualitative marks, applying the inverse transformation yields the spatio-temporal partial complete auto-covariance function, and the partial complete cross-covariance function Recapitulating that κ ij = ζ ij , notice that (24) is the partial cross-covariance density function ζ ij|V \{i,j} . Further, the partial correlation Cor ij|V \{i,j} and the partial scaled covariance density function ξ ij|V \{i,j} are then defined as Likewise, similar partial point process characteristics can be computed from the inverse transformation of multivariate-marked as well as dot-type partial cross-spectra expression yielding a rich toolbox of novel numerical summary characteristic for spatialtemporal point data. E.g., for the multivariate-marked case, inverse transformation of the marked auto-and cross-spectral density functions yield the partial auto-type mean product of marks U ii|V \{i,j} (c, t) and the partial cross-type mean product of marks U ij|V \{i,j} (c, t), respectively, which properties and interpretation has not been investigated yet and highly welcomes deeper evaluations in future research.

Multivariate spatio-temporal crime data
This section covers the application of the STDGM to a spatio-temporal crime dataset provided under the Open Government Licence by the British Home Office for London and has been downloaded from http://data.police.uk/data/. This data contains the longitude and latitude for a set of 14 pre-classified crime categories at street-level, either within a one mile radius of a single point or within a custom area of a street recorded by the Metropolitan Police as well as the month of occurrence for each single crime event. For our analysis we preselected all crime events which have been collected within a four-month period from April to July 2016 yielding a sample of 343427 single crime events from which 77139 events have been recorded in April, 86915 in May, 85761 in June and, lastly, 93612 in July. Finally, excluding any duplicated events our data reduces to 127328 events in total.
To give a first impression about the temporal variation of the spatio-temporal point pattern, the numbers of crimes and different numerical summary characteristics per month have been computed. The monthly numbers and spatial first-order intensity functions per crime category are reported in Table 1. Inspecting Table 1, we found that most of the crime events have been categorised as anti-social behaviour and also violence and sexual offences, while possession of weapons appeared least often. Further, an increase in numbers of cases of anti-social behaviour from April to July can be observed contrasting with a decrease in numbers of cases of violence and sexual behaviour in the same period. The intensity reports the mean number of events per unit area of the region considered in London.
We note that this current dataset contains only four temporal instants, and this lack of temporal information basically prevents from running a formal spatio-temporal anal- ysis using second-order characteristics as in Section 2. In addition, so many points in space for each type of crime make computational burden when using for example the spatial K-function. Two more points are in order. The crime data happens on the streets of London, and the network structure is key in the spatial structure of the events. Note that there is a large hole within the spatial region having to do with the network itself. This is not considered into account in the functions shown in Section 2 for a more classical spatio-temporal analysis. A second aspect is that this classical approach does not consider conditional relationships, and can only measure global bivariate cross-relationships when marks are discrete. Due to these number of drawbacks, we considered the bivariate spatial K-function per month for the pairs (Criminal damage and arson vs Violence and sexual offences) and (Robbery vs Vehicle crime). The corresponding K-functions per month are displayed in Figure 1, where the Poisson line is also depicted. In addition, they show a sort of regular, inhibitory structure between the two types of crimes for each pair considered. It is clear that time is not rightly considered here, and the analysis of the spatial structure neglects the remaining information from the other types. These facts motivate our new graphical modelling approach as follows.

Cross-sectional graphical modelling
To investigate the structural interrelations among the 14 crime categories from a crosssectional perspective and to evaluate the findings of the STDGM, we first discuss the SDGMs computed for each month separately. To this end, we split the data per month into four subsets and computed separate spatial auto-and cross-periodograms. To con- trol for possible variation in strength of the partial interrelations between different pairs of crimes, we considered a threshold level of ξ = 0.6 to discover partial interrelations with an intermediate effect size. That is, for each SDGM an edge is drawn between the nodes i and j if the supremum of the empirical absolute rescaled inverse spectral density function for components i and j equals or exceeds ξ for at least one frequency w for p = 0, 1, . . . , 16 and q = −16, . . . , 15 for the particular month. In this case, the point distributions of the components i and j recorded for a particular time t, t = 1, . . . , 4 are said to be interrelated. The resulting monthly SDGMs are depicted in Figure 2. Looking at these plots, we found 6 isolated nodes in Figures 2a to 2c and 7 isolated nodes in Figure 2d. For Figure 2a and Figure 2b, we observed an identical set of isolated nodes (public order, other theft, robbery, drugs, theft from the person, bicycle theft) while different sets of isolated nodes are shown in Figure 2c (possession of weapons, shoplifting, drugs, other crime, robbery and theft from the person) and Figure 2d (public order, bicycle theft, vehicle crime, criminal damage and arson, theft from the person, other theft and robbery). For all isolated nodes, we concluded that none of these crime patterns are interrelated to any alternative crime pattern included from a cross-sectional perspective. We note that comparing the isolated nodes over the complete period under study, only two crimes remain isolated throughout the four months, namely robbery and theft from the person, while most alternative isolated nodes appeared in, at most, three SDGMs. Further, other crime and shoplifting (resp. criminal damage and arson and vehicle crime) are not interconnected to any alternative crime in June (resp. July). We outline that the isolated nodes would imply that the distributions of point locations of any of these crimes obey complete spatial randomness conditional on all alternative crime patterns included by definition in the spatial dependence graph model.  Besides these isolated nodes, several subgraphs can be identified. Inspecting the upper panel, two subgraphs are shown in Figure 2a (a 2-node and a 6-node subgraph) while all alternative crimes are joined in only one 8-node subgraph in Figure 2b. Turning to the lower panel, we again observed a 2-node and a 6-node subgraph for June ( Figure  2c) whereas all non-isolated nodes form a subgraph in Figure 2d. Further, looking at all four spatial dependence graph models, we found that anti-social behaviour is only directly connected with violence and sexual offences which implies that anti-social behaviour is conditionally independent of all remaining crimes given violence and sexual offences from a purely spatial perspective for all months. Comparing this finding with Table 1, we found that anti-social behaviour, which appeared most often in all four months, is directly connected to the second most often crime category. At the same time, we also found that burglary is linked to violence and sexual offences throughout the complete period from a cross-sectional perspective. Except for June (Figure 2c), a direct interrelation can also be detected for possession of weapons and other crime. Interestingly, unlike anti-social behaviour and violence and sexual offences, we observed an opposite relation between the numbers of possession of weapons and other crime throughout all four months with high numbers for other crime while possession of weapons appeared least often.

Spatio-temporal dependence graph model results
We now discuss the results of the STDGM computed from the crime data over the four-month period. Unlike the cross-sectional analysis, the estimation of the STDGM is related to Fourier transformations of the spatial frequencies over time. We note that, as time is assumed to be recorded in equidistant steps in discrete time, the interval length has an important impact on the estimation of the STDGM. As pointed out by Didelez (2003), large intervals result in a marginalisation over time and information on short-term dependencies between different components might be lost. At the same time, additional correlation could emerge due to common causes which occurred in the meantime.
To control for possible variation in strength of the partial interrelations between different pairs of crimes over time, we consider a threshold level of ξ = 0.6 in order to detect conditional partial interrelations with an intermediate effect size such that an edge is drawn between the nodes i and j if the supremum of the empirical absolute rescaled inverse spectral density function for components i and j equals or exceeds ξ for at least one frequency w st for p = 0, . . . , 16, q = −15, . . . , 16 and u = −2, . . . , 2. That is, edges indicate that the strength of the linear partial interrelation between two component processes is greater than or equal to ξ = 0.6 over all four months. In this particular case, the spatial point distributions of the components i and j are said to be interrelated over time. To state this in a different manner, as the STDGM is defined through the Fourier transform of the spatial frequencies (p, q) at times t, edges represent periodicities of (p, q) over time. The resulting STDGM is depicted in Figure 3. Inspecting the STDGM, eight isolated nodes (public order, other crime, anti-social behaviour, drugs, bicycle theft, shoplifting, theft from the person, possession of weapons) and one 6-node subgraph can be identified. For the isolated nodes, we concluded that none of these crimes are interrelated to any alternative crime included over the fourmonth period. Further, comparing the isolated nodes of the STDGM with the four SDGMs depicted in Figure 2, no link is drawn joining anti-social behaviour and violence and sexual offences for the spatio-temporal case while this interrelation occurred in all cross-sectional plots from April to July 2016. This implies that, although both crimes are interrelated from a cross-sectional perspective, no periodic structures can be found over monthly time intervals. At the same time, as for the analysis of time series, periodicities might be detected for alternative interval lengths.
Turning to the 6-node subgraph, we observed that the spatio-temporal patterns of robbery as well as of vehicle crime are conditionally independent of all remaining crime patterns given the spatio-temporal distribution of other theft. Interestingly, we also observed that burglary is again linked to violence and sexual offences which also holds for the purely spatial dependence structures as depicted in Figures 2a to 2d. This implies that the interrelations of both crimes are also periodic over monthly time intervals. We emphasise that these findings would not have been detected by the classical spatial analysis as conditioning nor partialisation are not able in their case.  To control for possible variation in strength of the partial interrelations between different pairs of marked locations over time, we consider a threshold level of ξ = 0.5 in order to detect conditional partial interrelations with an intermediate effect size such that an edge is drawn between the nodes i and j if the supremum of the empirical absolute rescaled inverse spectral density function for components i and j equals or exceeds ξ for at least one frequency w st for p = 0, . . . , 16, q = −15, . . . , 16 and u = −2, . . . , 2. In this respect, edges indicate that the strength of the linear partial interrelation between two component processes is greater than or equal to ξ = 0.5 over all four years. The resulting STDGM is depicted in Figure 4. Inspecting the STDGM, three isolated nodes (fraxAmer, caryTome, franCaro) as well as one pair and one triplet of connected edges can be identified. For the isolated nodes, we concluded that the quantitative marks and the locations of these particular species are both unrelated to the quantitative marks and the locations of any alternative species included over the complete temporal period. We remind that the present formulation of the STDGM only allows for the detection of linear partial interrelations and potential non-linear interrelation are not captured. Turning to the pair and triplet of interconnected nodes, we found that the marks and locations of (a) caryGlab and carpCaro as well as (b) cercCana, acerRubrand cornFlor are not interrelated to those of any alternative species in the sample except of those species contained in the corresponding 2-node and 3-node subgraphs. At the same time, focussing on (b) and applying basic graph terminology we found that the marked locations of cercCana and cornFlor are conditionally uncorrelated given those of acerRubrand as acerRubrand serves as a separator in this subgraph. For completeness, we consider the only analysis we can run from the more classical spatio-temporal point of view. Two drawbacks are in order here. One is that as we have only four temporal instants, we can not provide a deep spatio-temporal analysis. Second, we can only show marginal analysis of the whole problem, one with bivariate K-functions, and the other with mark weighted K-functions, but not an overall analysis. We then report the cross spatial K-functions for each temporal bin for all connected types in Figure 4 as well as individual mark weighted K-functions (see e.g. Penttinen et al., 1992). K-functions are displayed in Figures 5 and 6. In all the cases the differences betweenK ij (r) and the benchmark clearly suggest that the components are dependent, which goes in the line found in Figure 4. The mark-weighted K-function for each recorded time is subtracted from its theoretical value under independence between locations and marks, i.e, from the classical K-function. We find quite strong indications of a substantial deviation from independence in this case, again reinforcing the results in Figure 4. 6.2 Multivariate-marked spatio-temporal crime data As a second application, the STDGM computed from multivariate-marked crime data is discussed next. This data contains information on a subset of point locations for both property and violent crimes provided as longitude and latitude, the precise date and time of the event, and different attributes of incidents reported in the Analytical Services Application (ASAP) crime report database by the District of Columbia Metropolitan Police Department (MPD). It has been provided under an Open Government Licence and downloaded from https://dcatlas.dcgis.dc.gov/crimecards/. This data is shared via an automated process where addresses are geocoded to the District's Master Address Repository and assigned to the appropriate street block. To compute a STDGM from this source, we initially extracted the month and year from the original time indication and calculated the duration of police investigation at place in seconds from two additional date and time indications in the data. This duration was then considered as quantitative mark yielding the desired multivariate-marked representation.
Restricted to the year 2019 and excluding any cases with either missing date, time or duration information yields a sample of 28999 events. Taking 12 months into account, we preselected a set of 27680 point locations from this sample restricted to a set of five distinct pre-specified crime categories. This data then serves as input for the STDGM where we, paralleling the STDGM of the Duke forest data, considered a threshold level of ξ = 0.5 to control for possible variation in strength of the partial interrelations between different pairs of marked locations over time. This implies that i and j are joint by an edge if the supremum of the empirical absolute rescaled inverse spectral density function for components i and j equals or exceeds ξ for at least one frequency w st for p = 0, . . . , 16, q = −15, . . . , 16 and u = −2, . . . , 2. The resulting STDGM is depicted in Figure 7.
Inspecting the STDGM, we found one isolated node (motor vehicle theft) and two pairs of pairwise interconnected nodes indicating that the marked locations of (a) robbery and burglary, and (b) other theft and car theft meaning theft from car. Interestingly, for (a) the marks and locations of two different types of crimes (one property and one violent crime) are interrelated to each other over the complete 12 month period, whereas both crimes of (b) as well as motor vehicle theft are property related crimes. Neglecting the spatial component and inspecting the number of incident for all five crimes exclusively over the complete period under study, both other theft and car theft show a periodic behaviour with peaks in the summer time while all alternative crimes in the sampled data reflect less fluctuation over time.  so we can conclude that the components are dependent, reinforcing the results in Figure  7. Once again, the mark-weighted K-function for each recorded time is centred by using its the classical K-function. As the lines oscillate around the theoretical value, we can say that there is not enough evidence of substantial deviation from the independence, meaning that the spatio-temporal locations, in this case, neglect somehow the values of the marks. This is different from what we have found in Figure 7, and this is a misleading result due to the marginal analysis done for the marks only in a spatial context. The STDGM highlights relations that with only a classical method could be missing.

Conclusions
The statistical investigation of potential interrelations in marked spatio-temporal point process is an still unresolved and yet highly challenging field of research which has just very recently been started to be explored. When extending classical cross-type characteristics for spatial marked point processes to the spatio-temporal domain, these mark statistics quickly become infeasible and computational burdensome when having large amounts of point locations in time, space or in space-time. To overcome these limitations, the present paper contributes to the multivariate analysis of spatio-temporal point process data by introducing different partial point characteristics and extending the spatial dependence graph model formalism yielding a unified framework for different types of spatio-temporal data including both, purely qualitatively (multivariate) cases and the so-called multivariate-marked spatio-temporal point processes where both qualitative and quantitative information is available for each point location. The proposed graphical model, defined through partial spectral densities characteristics, is highly computationally efficient and reflects in the multivariate-marked case the conditional similarity among sets of spatio-temporal sub-processes of marked points with identical discrete marks. In addition to the definition of the spatio-temporal dependence graph model, different partial spatio-temporal point process characteristics are introduced in the frequency spatio-temporal domain which enhance the classical methodology toolbox in multiple ways providing alternative information besides classical univariate and bivariate crossand dot-type point process statistics as well as traditional multivariate dimensionality reduction techniques. Finally, a new class of spectral characteristics is introduced which mirrors the ideas of classical spatial dot-type point point process characteristics to the frequency domain. These dot-type spectral characteristics reflect the interrelation between a particular pattern and one (resp. two) alternative subset (resp. subsets) of components where the pattern of interest is excluded.