Dynamically evolving Gaussian spatial fields

Extremes DOI 10.1007/s10687-010-0120-8

Dynamically evolving Gaussian spatial fields Anastassia Baxevani · Krzysztof Podgórski · Igor Rychlik

Received: 10 June 2010 / Revised: 28 July 2010 / Accepted: 28 July 2010 c Springer Science+Business Media, LLC 2010

Abstract We discuss general non-stationary spatio-temporal surfaces that involve dynamics governed by velocity fields. The approach formalizes and expands previously used models in analysis of satellite data of significant wave heights. We start with homogeneous spatial fields. By applying an extension of the standard moving average construction we obtain models which are stationary in time. The resulting surface changes with time but is dynamically inactive since its velocities, when sampled across the field, have distributions centered at zero. We introduce a dynamical evolution to such a field by composing it with a dynamical flow governed by a given velocity field. This leads to non-stationary models. The models are extensions of the earlier discretized autoregressive models which account for a local velocity of traveling surface. We demonstrate that for such a surface its dynamics is a combination of dynamics introduced by the flow and the dynamics resulting from the covariance structure of the underlying stochastic field. We extend this approach to fields that are only locally stationary and have their parameters varying over a larger spatio-temporal horizon.

Research by K. Podgórski was partially supported by the Swedish Research Council Grant 2008-5382. A. Baxevani (B) · I. Rychlik Department of Mathematical Sciences, Chalmers University of Technology and University of Gothenburg, 41296, Gothenburg Sweden e-mail: [email protected] I. Rychlik e-mail: [email protected] K. Podgórski Centre for Mathematical Sciences, Mathematical Statistics, Lund University, Box 118, 22100 Lund, Sweden e-mail: [email protected]

A. Baxevani et al.

Keywords Spectral density · Covariance function · Stationary second order processes · Velocity field AMS 2000 Subject Classifications 60G10 · 60G15 · 60G60

1 Introduction 1.1 Motivation Recent technological advances such as aerial laser and satellite scanning produce increasingly complex environmental data that covers large regions in space and relatively long periods of time. Examples of such data come from marine climate, air quality, and vegetation surveys. To account for all aspects of such spatio-temporal data is a challenging task and a proper stochastic framework has to be carefully designed to capture important features of the data. Typically Gaussian distribution is assumed and thus the main focus is on describing correlation structures. Consequently, non-stationary covariances are subject of intense studies, see for example, Cox and Isham (1988), Gupta and Waymire (1987), Porcu et al. (2009), Stein (2005) and references therein. Explicit formulas are typically limited to separable spatiotemporal covariances for which no interesting dynamics is observed (see Section 2). Thus, for dynamically evolving fields there is a need for non-separable dependence. Dynamics can be most conveniently described through velocity fields. For ocean water waves, where the dynamics is brought by the dispersion relation, statistical distributions of velocities have been studied, see Longuet-Higgins (1957), Baxevani et al. (2003), and Podgórski et al. (2000), and was proved that the extremal crest velocities involve partial derivatives of correlations in space and in time that are non-vanishing in the non-separable case. Dynamics is also present in non-stationary phenomena such as space-time variability of atmospheric pressure that manifests as movement of storms, cyclones, etc. Similarly as in the stationary case, the movement of extremes is associated with (nonstationary) correlation structures. To describe and understand dynamics of extreme episodes one has to propose physically interpretable stochastic models. In this work, this is achieved by using velocity fields to transport stochastic spatial fields generated at various time points. In this sense our approach can be viewed as a hybrid method that combines deterministic dynamics with random fields. The stochastic models with embedded velocity got attention in the turbulence theory, where the stationary covariance function satisfies the so-called Taylor’s hypothesis r (tv, 0) = r (0, t), where v is a constant velocity, see Taylor (1938). The simplest model of this type is the so-called frozen field movement. A more general class satisfying Taylor’s hypothesis was discussed in Gupta and Waymire (1987). Locally stationary

Dynamically evolving Gaussian spatial fields

autoregressive model (with slowly changing velocity v) was used in Baxevani et al. (2009), to model variation in time and space of significant wave height data from satellite records and a stationary buoy or a system of buoys. Our approach includes the above models as special cases but we do not assume a constant velocity field. The proposed general framework allows non-separable spatio-temporal covariances. The motion of stochastic surface can be observed by means of local velocities and we demonstrate that their median coincide with the deterministic velocity. Our work has been applied to describe the correlation structure of significant wave height in an area of the North Atlantic Ocean, see Ailliot et al. (2010), where the velocity field has been modelled as a hidden Markov process. In another application, the proposed framework has been used to compute the safety index for fatigue failure in a vessel, see Mao et al. (2008). More specifically, the moving spatial field of significant wave heights (standard deviation of the sea surface) is observed from a moving vessel and the temporal covariance of encountered significant wave heights is evaluated. This temporal covariance depends both on the ship velocities and on the spatio-temporal moving field and is needed to evaluate the variance of accumulated fatigue damage during a voyage. 1.2 Terminology and notation In environmental applications our main object of studies, a stochastic field X (p, t), has two arguments: a space variable p which represents position and a time variable t. We only study the case p = (x, y), however extension to higher dimensions is straightforward. We use the following terminology. A field is stationary if it is invariant with respect to shifts in time and space, i.e. if for each fixed p0 and t0 , we d

d

have that X (p + p0 , t + t0 ) = X (p, t), where = stands for the equality of underlying probability distributions of stochastic processes. Invariance with respect to shift in space (time) will be referred to as spatial (temporal) stationarity. If the field is d

invariant with respect to the space rotation, i.e. X (Rφ p, t) = X (p, t), where Rφ is the rotation by an angle φ, then we call X isotropic. Finally, a field that is isotropic and stationary in space is referred to as homogeneous. 1.3 The approach Gaussian stationary fields constitute a convenient class of models that have found many applications. In this work they serve as building blocks for more general, non-stationary models. The need to reach beyond stationarity is due to two aspects observed in environmental records: dynamics and spatio-temporal variability due to different properties at different locations (or/and time instants). The nature of these deviations from stationarity is different and they have to be treated accordingly. In our approach, dynamics are introduced through a deterministic flow that transports independently and locally generated stochastic fields. Further, long scale variability is represented by location and time dependent spectra of the underlying locally stationary fields. We start with a spatial-only covariance and introduce temporal

A. Baxevani et al.

dependence following a classical time series approach in which independent spatial innovations have the assumed spatial covariance structure. We then argue that properly defined velocities when sampled randomly from the surface are centered at zero and that this indicates that the fields are dynamically inactive. We introduce dynamics by means of a velocity field which represents the motion of the surface. This velocity field transports independent stochastic innovations which are created at each time point and weights them by a time dependent scaling that diminishes over time. When these scaled innovations are added up, the result is a field that is non-stationary in both space and time. The spatio-temporal variability is due only to different velocities at different locations and times and thus is dictated by the underlying flow. Analysis of the velocity distributions of these fields can be performed by a method of Baxevani et al. (2003). In this paper it is observed that the velocities are centered at a value which is the sum of the flow velocities and an additional term which represents the non-isotropic character of the underlying stochastic field. It is also important to account for spatial variability that is due to geographic structural variation. The approach is extended to account for this type of non-stationarity in space. This is achieved by taking spectral representations corresponding to stationary fields but making the spectra depend on parameters that vary from one location to another. 1.4 Relation to previous work Our main inspiration is the model that was described in a discretized version in Baxevani et al. (2009). It starts with a spatial stationary Gaussian covariance function r S (p) = σ 2 exp(−|p|2 /(2L 2 )). Then, temporal dependence is introduced by considering the recursive autoregression X (p, t) = ρ X (p − vdt, t − dt) +

1 − ρ 2 t (p),

(1)

with independent (in t) innovations t (p) which have the covariance r S . The model has simple interpretation: at each time step the past surface is moving forward to a new location with velocity v and is modified by an independent innovation with prescribed (fixed) spatial covariance structure. The resulting covariance is of the form r (p, t) = ρ t r S (p − vt), and is stationary. In Baxevani et al. (2009) the model was also extended to account for a lack of stationarity in space. This was accomplished at three different levels. The first was to ensure that the innovations t (p) had general non-stationary covariance r S (p , p). In particular, the dependence of L on location in the Gaussian covariances was considered. Secondly, the auto-regression coefficient ρ was made dependent on the location. Finally, it was also assumed that the velocity v depended on both location p and time t. In this paper, the special case of Eq. 1 where v = 0 will be referred to as the underlying static f ield and will be denoted by X 0 . We note that the first two sources of non-stationarity discussed above are due to the non-dynamical field while the last one is a result of the dynamical flow given by the velocities.


In this paper we extend the model in Baxevani et al. (2009) in several respects. First, we provide with a fully continuous set-up by using time moving averages of independent spatial fields, i.e. we construct the random fields X (p, t) =

∞ −∞

f (t − s) (p; ds).

(2)

The method of building spatio-temporal moving averages through integration of independent stochastic processes has been also used in Schlather (2009), where independent time processes were averaged over space. Here we found it more natural to introduce dynamics by averaging independent spatial fields over time. Secondly, and more importantly, we note that the construction is independent of the form of spatial covariance r S and can lead to fairly general time dependence with correlations ρ(t) ∼ f ∗ f˜(t), where ∗ stands for the operation of convolution and f˜(u) = f (−u). Essentially, for each spatial covariance r S (p) and a general class of temporal correlations ρ(t), we give explicit representations of the Gaussian fields which have covariance structure given by the product r S (p) · ρ(t). The model in Eq. 2 coincides with the one that in a discretized version was given by Eq. 1. The dynamics are expressed by an arbitrary time varying velocity field that generates a flow given by ψ t,h (p) which is the location at time t + h of a point that at time t is at p. Such a flow is incorporated into a stochastic framework by means of the stochastic integral (3) Y (p, t) = f (t, t − s; p) (ψ t,s−t (p); ds). Here, for fixed p and t the value f (t, t − s; p) represents the weight with which the innovation (ψ t,s−t (p), ds) contributes to the value Y (p, t). In this way we obtain a large class of Gaussian spatio-temporal fields that incorporate dynamical evolution of a random surface, see also Schlather (2009). A general scheme of fitting to the actual spatio-temporal data can be obtained by extension of the approach presented in Baxevani et al. (2009). The concept of integration that is used above, is based on the general methods of defining integrals with respect to spectral measures of orthogonal projections in Hilbert spaces (see, for example, Mlak 1991) and is standard in mathematical literature, so here we only sketch fundamentals in Appendix A.1. The generality of the approach allows a natural extension to second order models that goes beyond Gaussian distribution. However, this is not explored here. Steps in this direction have been undertaken for the fields driven by Laplace motion, see Åberg and Podgórski (2008), and will be continued in future research.

2 Spatio-temporal static fields Models for proper spatio-temporal dependence structure of stochastic fields have been intensively studied in the context of environmental and geo-sciences (Cox and

A. Baxevani et al.

Isham 1988; Gupta and Waymire 1987; Porcu et al. 2009; Stein 2005; Schlather 2009). The goal of this section is to present a unified approach to construction of spatio-temporal dependence models. Several explicit formulas for separable and non-separable covariances are given. 2.1 Locally stationary spatial fields Before we turn to the building of spatio-temporal structures let us briefly discuss a way to obtain a fairly general class of spatial non-stationary fields. The starting point is the following spectral representation of a stationary process d X (p) = exp(ip · ω) S(ω) d B(ω), (4) Rn

where the symmetric function S(ω) ≥ 0 is a spectral density and d B is a centered random Gaussian measure whose variance coincides with the Lebesgue measure in Rn . A natural extension of Eq. 4 to non-stationary fields is by considering spectra that depend on location. More precisely, for a family of symmetric spectral densities Sp (ω) ≥ 0 parameterized by p, we define d X (p) = exp(ip · ω) Sp (ω) d B(ω). Rn

The covariance of X is given by

r S (p, p ) = Cov(X (p), X (p )) =

Rn

exp(i(p − p ) · ω) Sp (ω)Sp (ω) dω,

and is non-stationary. If Sp (·) ≈ Sp0 (·) in some neighborhood of p0 , then X (·) in this neighborhood can approximately be viewed as a realization of X (p) = exp(ip · ω) Sp0 (ω) d B(ω). Rn

Thus such random fields provide a convenient way of modelling spatial phenomena that are locally stationary but globally non-stationary. Example 1 (Non-stationary locally isotropic covariance) In this example, n is an arbitrary natural number although we are mainly interested in n = 1, 2. For alternative derivations of the covariances, see also Stein (2005), Schlather (2009), and Porcu et al. (2009). In Baxevani et al. (2009), we considered isotropic (invariant on rotation) spectra that locally had the so-called Gaussian form Sp (ω) =

s 2 (p)L n (p) exp − L 2 (p)|ω|2 /2 , n/2 2π


where s 2 (p) is the variance at a location p. The covariance of processes with such spectra is given by s(p)s(p ) L(p)L(p ) n/2 r S (p, p ) = 2 π Rn 2 (p) + L 2 (p ) L |ω|2 dω × exp i(p − p ) · ω − 4 n/2 s(p)s(p ) L(p)L(p ) = exp(i(p − p ) · ω) n 2 π R

× exp −ω T −1 ω/2 dω, where = (p, p ) = 2/ L 2 (p) + L 2 (p ) · I, with I being the identity matrix in Rn . Using the formula for the characteristic function of n-dimensional Gaussian vector we obtain −n/2

2L(p)L(p ) T r S (p, p ) = s(p)s(p ) exp −(p − p ) (p − p )/2 . L 2 (p) + L 2 (p ) (5) Consequently, if s(p) and L(p) are approximately constant in some region, then the correlation is approximately invariant with respect to isometries of Rn . In such a case, we refer to the field as locally isotropic. Obviously, by taking an arbitrary positive definite matrix we obtain an anisotropic extension of the model. Here and in the rest of this paper we use the convention that vectors are column matrices, and for a matrix A its transpose is denoted by AT . 2.2 Building spatio-temporal dependence In Appendix A.1, a notion of integral which gives meaning to the general spatiotemporal field X (p, t) = f (t, s; p) (p; ds), (6) is introduced for a deterministic kernel f , where (·; ds) is a Gaussian field-valued measure that is uniquely characterized by the time dependent spatial covariances r S (p, p ; s). Here, they will be referred to as spatial covariances governing X . As an example of r S (p, p ; s), one can consider the covariances of the previous subsection. Dependence on time can be introduced quite arbitrarily by making the spectra Sp to also depend on time t. The model above is the most general form of static f ields that are discussed in this work. Since the fields we are interested in are centered Gaussian fields, to compare the different models is enough to compare their covariance functions, f (t, s; p) f (t , s; p ) · r S (p, p ; s) ds. r (p, p ; t, t ) = Cov(X (p, t), X (p , t )) = (7)

A. Baxevani et al.

One important simplification of the model is obtained by taking r S (p, p ; s) = r S (p−p ; s). Several other specifications of the above model will be used throughout the paper. They are listed below in order of increasing complexity. We commence the presentation with some notation. We write r (p, p ; t, t ) for the covariance between X (p, t) and X (p , t ). In the presence of spatial stationarity we write r (p; t, t ) while if the field is stationary in time we write r (p, p ; t). Finally if the field X is stationary both in space and time we write r (p; t).

Separable stationary moving average This case is defined by taking f (t, s; p) = f (t −s), and thus the kernel f does not depend on p, and stationary spatial covariance r S (p − p ) that is independent of time. This case is the completely stationary case, with covariance given by r (p; t) = r S (p) · r T (t), where r T (t) = f (t − s) f (−s) ds. If additionally we assume that r S (p) is isotropic we obtain the subcase of a homogeneous moving average field.

Separable time stationary moving average A generalization of the previous case is when the spatial stationarity is dropped and the covariance is given by r (p, p ; t) = r S (p, p ) · r T (t),

(8)

where r S (p, p ) is a non-stationary spatial covariance. Observe the temporal stationarity of the model.

Separable covariance model This is the case when the kernel f is independent of the space variable p, and the spatial covariance is independent of the time variable t. In this case the covariance can be still presented as a product of the spatial and temporal covariances r (p, p ; t, t ) = r S (p, p ) · r T (t, t ), ∞ where r T (t, t ) = −∞ f (t, s) · f (t , s) ds. This is sometimes referred to as multiplicative separability of the model.

Heteroscedastic moving average This is the case of time dependent spatial covariance structure with space independent kernel f (t − t ) for which ∞ r (p, p ; t, t ) = f (t − s) · f (t − s) · r S (p, p ; s) ds. (9) −∞

We note that typically this model is non-stationary both in time and space. The terminology is borrowed from the general theory of time series as a spatial analog of the non-constant variance innovation case. We also consider a heterodscedastic,


space-stationary moving average which is defined by stationary covariances r S (p − p ; s). Clearly in this case, there is stationarity in space.

Temporal stationary moving average Here we assume a time independent (homoscedastic) spatial covariance structure with space dependent kernel f (t, s; p) = f (t − s; p) for which r (p, p ; t) = r S (p, p ) · f p ∗ f˜p (t),

(10)

where f p ∗ f˜p is the convolution of f p (s) = f (s; p) with f˜p (s) = f (−s; p ). We note that there is stationarity in the time direction at any fixed spatial position p while temporal models differ at various locations. Example 2 An interesting example of the temporal stationary moving average can be obtained by taking for f p a density of gamma random variable with the scale parameter λ(p) and the shape parameter τ , i.e. f p (t) =

1 t τ −1 e−t/λ(p) . (τ )λ(p)τ

Then, in each fixed space point p the temporal dependence is represented by Mattérn covariance while r S (p, p ; t) is given by r S (p, p ) multiplied by a density of generalized asymmetric Laplace distribution, see Kotz et al. (2001) for explicit formulas in terms of the Bessel functions. Remark 1 (Temporal moving averages) The temporal moving average field is introduced by Eq. 8 for f (t, s) = f (t − s) and by assuming that the spatial covariance function r S is independent of time. The relation to moving averages appearing in time series analysis can be more explicitly seen through the approximation of the field by a sum. Let s = k t, k = −M, . . . , M for some large M and t = n t. Then, X (p, t) ≈

M

f ((n − k) t) · k (p) ·

√

t,

(11)

k=−M

where k (p) are independent (in time) Gaussian fields with Cov(k (p), k (p )) = r S (p, p ) that are related to the fields by k (p) =

(p; [k t, (k + 1) t]) . √

t

The relation (11) can be rewritten as X n (p) =

lim

M→∞, t→0

M k=−M

αk n−k (p),

(12)

A. Baxevani et al.

√ with αk = t · f (k t), which is the well known form of a discrete moving average time series. Example 3 (Temporal Ornstein–Uhlenbeck field) A special case of the separable temporal moving average model mentioned in Example 2 is obtained by taking f (t) = e−λt 1[0,∞) (t). In this case, X (p, t) =

t −∞

e−λ(t−s) (p; ds)

(13)

and since additionally r S (p, p ; s) = r S (p, p ), its covariance is given by r (p, p ; t) = r S (p, p ) ·

1 −λ|t| e . 2λ

(14)

This example corresponds to the case considered in Baxevani et al. (2009), where the autoregression model of order one X (p, t) = ρ X (p, t − t) +

1 − ρ 2 t (p),

has been discussed. It is clear from Remark 1 that the above is a discretization of the Ornstein–Uhlenbeck model with ρ = e−λ t . The space dependent Ornstein–Uhlenbeck model is obtained as a special case of temporal stationary moving average (Eq. 10) by taking a space dependent λ(p) in which case we obtain

r S (p, p ) e−λ(p )·t ; if t > 0, r (p, p ; t) = (15) λ(p) + λ(p ) e−λ(p)·t ; if t < 0. We note stationarity in time as in any other space dependent moving average. Example 4 (Temporal Gaussian dependence) Another case of the temporal stationary moving average model (Eq. 10) can be obtained by taking the Gaussian kernel 2 2 f p (t) = f (t; p) = π −1/4 · e−t /L (p) . By Eq. 10 we have r (p, p ; t) = r S (p, p ) · f p ∗ f˜p (t), and since the convolution of Gaussian kernels is again a Gaussian kernel, the resulting covariance is stationary in t and given by

r (p, p ; t) = r S (p, p ) ·

1 1 + 2 2 L (p) L (p )

−1/2

·e

−

(t−t )2 L 2 (p)+L 2 (p )

.

(16)

Use of space dependent kernel f (t, s; p) is natural for building non-stationary spatio-temporal correlation. By analogy to the approach of Section 2.1, one can


alternatively consider space dependent temporal spectra. More specifically, one can consider the model ei(p,t)·(ω,τ ) Sp (ω)SpT (τ ) d B(ω, τ ), (17) Rn+1

where SpT (τ ) is a location dependent temporal spectrum. This alternative construction is usually equivalent to the one based on space dependent symmetric kernels as stated in the following result. Theorem 1 Let a Gaussian random f ield X 1 (p, t) be def ined through Eq. 17 and a Gaussian space dependent moving average X 2 (p, t) be def ined through Eq. 10. Further assume that R eitτ SpT (τ ) dτ = f p ∗ f˜p (t), where f˜p (t) = f p (−t), so that the covariances in time at a f ixed point p are the same for X 1 and X 2 . If the kernels f p are symmetric and have non-negative Fourier transform, then the spatio-temporal d

covariances for both processes agree and consequently X 1 = X 2 . Proof The spatio-temporal covariance of Eq. 17 is given by Cov(X 1 (p, t), X 1 (p , t )) = ei(p−p )ω Sp (ω) · Sp (ω) dω · Rn · ei(t−t )τ SpT (τ ) · SpT (τ ) dτ R ei(t−t )τ SpT (τ ) · SpT (τ ) dτ, = r S (p, p ) · R

while for space dependent moving averages the covariance is given by Cov(X 2 (p, t), X 2 (p , t )) = r S (p, p ) · f p ∗ f˜p (t − t ). By taking Fourier transform and using the notation F h(τ ) = e−iτ t h(t) dt we obtain that 2 F f p ∗ f˜p = F f p = 2π SpT , so F f p = 2π SpT . On the other hand F

−1

2π

SpT

·

SpT

(t − t ) =

R

ei(t−t )τ

SpT (τ ) · SpT (τ ) dτ.

Consequently, by the inverse Fourier theorem and F f p ∗ f˜p = F f p · F f˜p = 2π SpT SpT

we obtain equality of the covariances.

A. Baxevani et al.

Remark 2 The above result assumes symmetric kernels. The Ornstein–Uhlenbeck process is not covered by it since this process is not represented by a symmetric kernel. The spectra are given by SpT (τ ) =

1 λ2 (p) + τ 2

,

while the kernel approach leads to the correlation given in Eq. 15. Thus the equivalence of the models would mean that e−λ(p )t / λ(p) + λ(p ) is equal to R

e−itτ

1 (λ2 (p) + τ 2 )1/2

·

1 (λ2 (p ) + τ 2 )1/2

dτ

which obviously is not true. We conclude therefore that the symmetry of kernels can not be dropped from the assumptions of Theorem 7.

3 Distributions of apparent velocities To define the motion of a surface is a non-trivial task and a proper definition of velocity emerges as a fundamental issue in describing dynamics of surfaces. There is no unique approach to this problem. For a comprehensive treatment of the issue we refer to Baxevani et al. (2003), Longuet-Higgins (1957), and Podgórski et al. (2000). Below we focus on a conceptually simple surface velocity that was first introduced in the pioneering work of Longuet-Higgins (1957). Once the concept of apparent velocity is decided for, it is possible to find its statistical distribution utilizing, for example, Palm measures on the crossing levels. Generalized Rice formula can be utilized to effectively evaluate such distributions (Rice 1944, 1945; Zähle 1984). There is a number of results dealing with velocities in recent works, see Azais and Wschebor (2009), Azais et al. (2005) and Baxevani et al. (2003). Here we consider only velocities at a fixed (non-random) point of the surface although the crossing level distributions could be obtained using the results of the above mentioned work. Apparent velocity distributions have been used in safety and extreme analysis. For example in Mao et al. (2008) the apparent velocity V = (Vx , Vy ) of the significant wave height has entered the formula for covariance of the encountered field z(t): Cov(z(t1 ), z(t2 )) = σ 2 ρ S ((vx − Vx ) (t2 − t1 ), (v y − Vy ) (t2 − t1 ))ρT (t2 − t1 )

= r z (t2 − t1 ),

(18)

where ρ S denotes the spatial correlation and ρT the temporal one, for a vessel sailing with constant velocity (vx , v y ). Encountered dangerous waves overtaking a vessel is


another situation in which distribution of encountered velocities is needed. This was studied in Åberg et al. (2008). Velocities have been also used for evaluating motion of extreme episodes. For example, in Baxevani et al. (2003) and Baxevani and Rychlik (2007), rough wave distributions have been obtained under different spatio-temporal scenarios. 3.1 Velocities of a random field Let X (p, t) as before be a Gaussian random field defined through the stochastic integral in Eq. 6, for a sufficiently smooth kernel f so that the process has well-defined partial derivative fields. We introduce velocity in an arbitrary but fixed direction. Since a simple rotation would allow us to obtain velocity in any direction we focus on the velocity along the direction of the x-axis. Indeed, for Rφ =

cos φ − sin φ

sin φ cos φ

,

the matrix representing rotation by angle φ, we can consider X (Rφ p, t) =

f (t, s; Rφ p) (Rφ p; ds),

(19)

which is the field X (p, t) rotated by angle φ. Here (Rφ p; ds) is the Gaussian field φ measure governed by covariance r S (p, p ; s) = r S (Rφ p, Rφ p ; s). Thus to obtain the results in an arbitrary direction φ from the ones given along the x-axis one only has to φ substitute r S for r S and f φ (t, s; p) = f (t, s; Rφ p) for f (t, s; p). Later on, when the dynamical flow is governed by a velocity field v(p, t), one instead has to substitute v(p, t) by the rotated velocity field vφ (p, t) = v(Rφ p, t). The speed, V , in the x-direction is the x-coordinate of the slope of the tangent plane to the up-crossing contour attached to the point (x, y, t). Thus we have V =−

Xt , Xx

(20)

at the points (x, y, t) such that X x = X x (x, y, t) > 0 and X = X (x, y, t). Here X x = X x (x, y, t) and X t = X t (x, y, t) are the first order partial derivatives of X = X (p, t) with respect to x and t, respectively and p = (x, y). Remark 3 If we know the scalar velocities along all directions we can obtain the vector velocity π by integration. For this, at any point (x, y), we define a new velocity V(x, y) = 0 (cos(φ), sin(φ)) · V (x, y, φ) dφ, where V (x, y, φ) for the special case φ = 0 is the velocity defined in Eq. 20 and for other φ is its analog in the φ direction, as described above.

A. Baxevani et al.

3.2 Distributions To obtain the one-dimensional marginal distribution of the velocity V defined in Eq. 20, we use some standard facts from the theory of centered Gaussian random vectors. For a jointly Gaussian vector (X x , X t ), we can write X t = E(X t |X x ) + s X · Z , where E(X t |X x ) = Cov(X x , X t )/Var (X x ) · X x and Z = [X t − E(X t |X x )]/s X is a standard Gaussian random variable independent of X x and s X2 = d

Hence, since X x = independent of Z , we get V =−

Var (X x )Var (X t ) − Cov(X x , X t )2 . Var (X x )

√ Var (X x )Z 1 , where Z 1 is a standard normal variable Xt Cov(X x , X t ) sX = − −√ ·C Xx Var (X x ) Var (X x )

(21)

where C is a random Cauchy variable defined as the ratio Z /Z 1 of two independent normal variables. We use this velocity to describe local dynamics of a stochastic field. We say that a Gaussian stochastic field X (p, t) does not exhibit any organized movement at the point (x, y) and at time t in the direction x, if the median of the distribution of V is equal to zero, i.e. if Cov(X x (p, t), X t (p, t)) = 0.

Theorem 2 A heteroscedastic space-stationary moving average X (p, t) that is def ined by Eq. 9 with r S (p, p ; s) = r S (p − p ; s) does not exhibit any organized movement. Proof Indeed, since the field X (p, t) is governed by a stationary covariance function r S , this follows from Lemma 2 in the Appendix A. Notice that in contrast to the case above, if the kernel is space dependent then the resulting field may exhibit non-trivial dynamics as their velocities are no-longer centered at zero. Theorem 3 For X (p, t) def ined by Eq. 6 with space-stationary innovations def ined through covariance r S (p; s), the center of the velocity in Eq. 20 equals x f (t, s; p) · f t (t, s; p) · r S (0; s) ds, − (22) A where A = | f x (t, s; p)|2r S (0; s) + | f (t, s; p)|2r Sx x (0; s) ds.


Proof This is a direct consequence of Eq. 52 in the Appendix A, Section A.2, since the first order partial derivative with respect to x of the covariance r S (p; s) equals zero when evaluated at p = 0 for the reasons explained in Theorem 1 of the Appendix A. Another way of obtaining dynamics is through space varying scaling. That type of “dynamics” may be not desirable and the following example shows that eliminating the space scaling variability (in terms of variance) should probably precede analysis of the ‘real’ dynamics. Example 5 (Deterministic rescaling) Let X (p, t) be a heteroscedastic spacestationary moving average X (p, t) that is defined by Eq. 9 with r S (p, p ; s) = r S (p − p ; s). By Theorem 2, it does not exhibit any organized motion since r Xxt (p, t) = 0. Consider now a positive deterministic field A(p, t) and define Y (p, t) = A(p, t)X (p, t). Then the covariance of Y (p, t) is given by rY (p, p ; t, t ) = A(p, t)A(p , t )r X (p − p ; t − t ). Consequently by Lemma 2 of the Appendix A, we have Cov(Y x (p, t), Y t (p, t)) = A x (p, t)At (p, t) · r X (0; 0) + A x (p, t)A(p, t) · r Xt (0; 0),

2 Var (Y x (p, t)) = A x (p, t) · r X (0; 0) + A2 (p, t) · r Xx x (0; 0), The center of velocity is thus given by −

A x At r X + A x Ar Xt . (A x )2r X + A2r Xx x

In general, the field Y has non-trivial dynamics (unless A does not depend on space variable) identified by the velocity center as given above. Since the underlying field X has no dynamics we conclude that the organized movement of Y is only due to the deterministic rescaling A. To avoid the dynamics illustrated in the above example, in practice, the variable space rescaling can be eliminated by dividing√the data by local in space standard deviation, i.e. by replacing X (p, t) by X (p, t)/ Var (X (p, t)) (for example by local estimation of the variance). Therefore, typically, we consider a version of the presented models for which r S (p, p; t) = σ 2 (t) or even r S (p, p; t) = 1, i.e. the variance of innovations is space independent. In such case, as long as the kernel f (t, s) is independent of the space variable p, the variance of X (p, t) is only time dependent. For this more general model the result of Theorem 2 remains valid as stated in the following Theorem. The proof is a direct consequence of Lemma 3 of the Appendix A.

A. Baxevani et al.

Theorem 4 A heteroscedastic moving average X (p, t) with innovations which have space independent variance does not exhibit any organized movement. Despite a possibility of introducing dynamics through the space variable kernels as shown in Theorem 3, it is difficult to give a natural interpretation to the obtained center of velocities. In the next section we turn to a more direct method of introducing dynamics in stochastic fields. The method is based on using a deterministic flow generated by velocity fields.

4 Dynamics in the models Contribution of this main part of the paper is two-fold. First, the means of combining spatio-temporal fields discussed in the previous sections with a deterministic dynamical flow are presented. The velocity field defining the flow is then compared with the apparent velocities. It is shown that locally the center of the apparent velocity distribution coincides with the deterministic one. 4.1 Constant velocity dynamics If in Eq. 3, we use constant velocity field v(p, t) = v = (v1 , v2 ) and space independent kernels f , we obtain the field Y (p, t) =

∞ −∞

f (t, s) (p + v(s − t); ds)

(23)

with covariance

r (p, p ; t, t ) =

∞

−∞

f (t, s) · f (t , s) · r S (p + v · (s − t), p + v · (s − t ); s) ds.

A notable special case is given by a spatial stationary innovation covariance r S as r (p − p ; t, t ) =

∞

∞

f (t, s) · f (t , s) · r S (p − p − v(t − t ); s) ds

in which we observe that the dynamic field is equivalent to the static field subordinated to the deterministic dynamics, i.e. Y (p, t) = X (p + v · t, t). Theorem 5 The center of velocities in the x-direction of the f ield Y def ined in Eq. 23, is given by xy | f (t, s)|2r S (0; s) ds v1 + v2 . | f (t, s)|2r Sx x (0; s) ds


If in addition the innovations are homogeneous (isotropic and stationary), then this velocity equals the constant f low velocity component v1 . Proof By Lemma 1 of the Appendix A we have Cov(Y x (p, t), Y t (p, t)) =

∞

f (t, s) · f t (t, s) · r Sx (0; s) ds

−∞

+ v1 + v2 = v1

∞

−∞ ∞

−∞ ∞

−∞

+ v2

| f (t, s)|2 · r Sx x (0; s) ds xy

| f (t, s)|2 · r S (0; s) ds

| f (t, s)|2 · r Sx x (0; s) ds ∞

−∞

xy

| f (t, s)|2 · r S (0; s) ds,

and Var (Y (p, t)) = x

∞ −∞

| f (t, s)|2 · r Sx x (0; s) ds. xy

For the homogeneous field we can use Lemma 2 to conclude that r S (0; s) = 0 and the result follows. Example 6 Consider the Ornstein–Uhlenbeck type time dependence. Then

r (p, p ; t, t ) = e−λ(t+t )

t∧t −∞

e2λs r S (p − p − v(t − t ); s) ds.

If additionally r S does not depend on time, then r (p, p ; t, t ) =

1 r S (p − p − v(t − t )) · e−λ|t−t | . 2λ

This type of spatio-temporal covariance has been used in Gupta and Waymire (1987) and Cox and Isham (1988) for modeling rainfall. 4.2 Spatio-temporal dynamical models In the static scheme described in Section 2, a stochastic field X (p) has been built from independent innovation fields (p; ds) that occurred at time s and were summed up while weighted by f (s). Dynamics can be introduced to this model by assuming that the contribution to a field Y (p, t) from the innovation field (·; ds) that occurred at time s is not evaluated at the point p but at the point ψ t,s−t (p) that corresponds to the

A. Baxevani et al.

location at time s of what at time t is at p. This has been presented in the previous section for the constant velocity dynamics where ψ t,h (p) = p + v · h, and can be generalized as follows. Let us consider a flow ψ t,h (p) obtained from a velocity field v(p, t) satisfying the transport equation ψ t,h (p) = p +

t

t+h

v(ψ t,u−t (p), u) du = p +

h 0

v(ψ t,s (p), t + s) ds,

(24)

i.e. a point with the initial location p at t relocates after h time units to ψ t,h (p). In what follows it will be convenient to use the notation ψ(p, t, h) for ψ t,h (p) and ψ t (p, t, h) = ψ tt,h (p) for the partial derivative of ψ with respect to t. We note the following basic properties ψ(p, t, 0) = p, ˜ = ψ(p, t, h + h), ˜ ψ(ψ(p, t, h), t + h, h) ∂ψ (p, t, h) = v(ψ t,h (p), t + h), ∂h ∂ψ (p, t, 0) = v(p, t). ∂h

(25)

Construction of the stochastic field at a fixed location p and a fixed time t Y (p, t) =

∞

−∞

f (t, s) (ψ t,s−t (p); ds)

(26)

with dynamics driven by ψ is obtained from the following elements • • • •

(p, (s, s + ds])—the field generated at time s. The field is assumed to be independent at different time instants s, f (t, s)—a weight function defining how much the spatial field (p, (s, s + ds]) should contribute to the field Y (p, t) , ψ t,s−t (p)—the location at time s of a flow element that at time t resides at p, (ψ t,s−t (p), (s, s + ds])—the field contributes to Y (p, t) at time s. Consequently, the contribution to Y (p, t) at time s is coming from

s+ds

t (p; ds) := s

(ψ t,s−t (p); ds)

(27)

∞ multiplied by f (t, s) and the integral Y (p, t) = −∞ f (s, t) t (p; ds) in its essence does not differ from the one defined in Eq. 6. One has simply to consider r St (p, p ; s) = r S (ψ t,s−t (p), ψ t,s−t (p ); s)


instead of r S (p, p ; s). Thus if we have two fields X (p, t) and Y (p , t ) with corresponding functions f and g, we obtain the cross-correlation formula ∞ Cov(X (p, t), Y (p , t )) = f (s, t) · g(s, t ) · r S (ψ t,s−t (p), ψ t ,s−t (p ); s) ds. −∞

(28) The technical but standard details of the construction above are omitted. They mainly hinge on the definition (Eq. 27). The median velocity of the obtained field Y (p, t) is given in the next result. Theorem 6 Let the f ield measure (p; ds) be driven by innovations that are stationary in space, so r Sx (0; s) = 0. Then the center of the velocity in the x-direction is given by xx xy rS rS t ∂ψ t,s−t (p) T 2 | f (t, s)| · ψ t,s−t (p) − v(ψ t,s−t (p), s) ds ∂x yx yy rS rS xx V = , xy rS rS T ∂ψ (p) ∂ψ (p) t,t−s | f (t, s)|2 · t,t−s ds ∂x ∂x yx yy rS rS (29) yy

xy

yx

where r Sx x , r S , r S , r S are all evaluated at ψ t,s−t (p). If additionally it is assumed that the innovations are isotropic, then x yy x (p) · αt,s−t (p) · r Sx x + yt,s−t (p) · βt,s−t (p) · r S ds | f (t, s)|2 · xt,s−t , V =

2 x 2 yy 2 x xx | f (t, s)| · xt,s−t (p) · r S + yt,s−t (p) · r S ds (30) where t αt,h (p) = xt,h (p) − v1 (ψ t,h (p), t + h), t βt,h (p) = yt,h (p) − v2 (ψ t,h (p), t + h). t (p), y t (p)), and v(p, t) = Here ψ t,h (p) = (xt,h (p), yt,h (p)), ψ tt,h (p) = (xt,h t,h (v1 (p, t), v2 (p, t).

Proof The proof is a direct consequence of Eq. 21 after applying Lemma 2 and the formulas for covariances given in Eq. 52, both in the Appendix A. Example 7 (Temporal dynamic Ornstein–Uhlenbeck) A dynamic modification of the Ornstein–Uhlenbeck model discussed in Example 3 is obtained by taking Y (p, t) =

t

−∞

e−λ(t−s) t (p; ds)

(31)

A. Baxevani et al.

with covariance

Cov(Y (p, t), Y (p , t )) = e−λ(t+t )

t∧t

−∞

e−2λs r S (ψ t,s−t (p), ψ t ,s−t (p ); s) ds. (32)

We note that the covariance is no longer separable even if r S (p, p ; s) = r S (p, p ). Example 8 (Dynamic autoregression field of order one) A discretized way of introducing dynamics represented by a flow ψ for some suitably chosen time lag dt is through the following recurrence Y (p, t) = ρY (ψ t,−dt (p), t − dt) +

1 − ρ 2 t (p),

(33)

where ρ = ρ(dt) = e−λdt for some λ > 0 and t (p) represent independent in time fields with the spatial covariance Cov(t (p), t (p )) = r S (p, p ; t). This example has been discussed in detail in Baxevani et al. (2009). In the next result, it is shown that the last two examples describe in fact the same model. Theorem 7 Consider a spatio-temporal centered Gaussian f ield, Y (p, t) def ined by the recursive formula (Eq. 33). This f ield has a covariance function that converges with time increment dt decreasing to zero to the covariance function (Eq. 32) scaled by 2λ. Proof Let t = kdt for some k ∈ Z. Using the flow properties given in Eq. 25, the recursive formula in Eq. 33 can be rewritten in a non-recursive way assuming that the series below is convergent Y (p, kdt) =

∞

ρ j 1 − ρ 2 (k− j)dt (ψ kdt,− jdt (p))

j=0

=

k 1 − ρ2 ρ k−l ldt (ψ kdt,(l−k)dt (p)).

(34)

l=−∞

For t ≤ t = k dt, the covariance function is given by Cov(Y (p, t), Y (p , t )) = (1 − ρ 2 )

k

ρ k+k −2l

l=−∞

· r S (ψ kdt,(l−k)dt (p), ψ k dt,(l−k )dt (p ); ldt).

(35)


For small values of dt, 1 − ρ 2 ≈ 2λ · dt, and therefore letting dt → 0 we have lim Cov(Y (p, t), Y (p , t )) = 2λ

dt→0

t

−∞

e−λ(t+t −2u)r S (ψ t,u−t (p), ψ t ,u−t (p ); u) du. (36)

The most complicated dynamics appear when the temporal dependence of the field Y (p, t) varies in space and additionally there are dynamics introduced by a deterministic flow ψ. A model of this kind can be written as Y (p, t) =

∞

−∞

f (t, s; p) t (p; ds),

and has the covariance Cov(Y (p, t), Y (p , t ))=

∞ −∞

f (t, s; p)· f (t , s, p )·r S (ψ t,s−t (p)ψ t ,s−t (p ); s) ds.

Example 9 (Space varying Ornstein–Uhlenbeck field driven by a dynamical flow) This model is the extension of Eq. 31 that is obtained by letting the parameter λ in Eq. 13 depend on p, i.e. by considering Y (p, t) =

t −∞

e−λ(p)(t−s) t (p; ds),

(37)

for some function λ(p). Then Cov(Y (p, t), Y (p , t ))

= e−λ(p )t −λ(p)t

t∧t

−∞

e(λ(p)+λ(p ))s r S (ψ t,s−t (p), ψ t ,s−t (p ); s) ds. (38)

In particular, when the flow is generated by a constant velocity v (the case that is important in local approximations of more general fields) and r S (p, p ; s) = r S (p − p ), we obtain for t < t Cov(Y (p, t), Y (p , t )) = r S p − p − v(t − t ) ·

1 · e−λ(p )(t −t) . (39) λ(p) + λ(p )

Acknowledgments The authors would like to thank Jörg Wegener for his comments on surface velocities and the anonymous referees for numerous suggestions. All these greatly improved the presentation of the results. They would also like to acknowledge financial support from Gothenburg Mathematical Modelling Center.

A. Baxevani et al.

Appendix A In this appendix, a formal presentation of the following method of constructing spatio-temporal stochastic fields is detailed. First, independent spatial fields are generated at a time point, and then they are convoluted to build in a temporal dependence. A standard Hilbert space construction of integral is used, see for example Mlak (1991).

A.1 Integration with respect to field valued random measure For each t ∈ R, let r S (p, p ; t) be a spatial covariance in p and p (non-negative definite function). We interpret it as a spatial covariance of independent innovations created at time t. We assume that for each a < b, the following integral is well defined as a function of p and p : b r(a,b] (p, p ) = r S (p, p ; s) ds, a

and thus corresponds to a certain spatial covariance function. It follows from the additivity of the covariance function with respect to independent fields and its correspondence to the additivity of the integral that there exists a family of Gaussian spatial fields (p; (a, b]) centered at zero such that (i) For each a < b, c < d ∈ R we have r(a,b]∩(c,d] (p, p ) = Cov((p; (a, b]), (p ; (c, d])). ∞ (ai , bi ], where (ai , bi ] are disjoint intervals, we have with (ii) For (a, b] = i=1 probability one (p; (a, b]) =

∞

(p; (ai , bi ]).

i=1

Thus is a σ -additive measure having as values Gaussian random fields and (ii) is evocative of Lebesgue integration. Consequently, for a step function f (t) =

n

αi 1(ai ,bi ] (t),

i=1

where (ai , bi ] are disjoint we define X (p) := X (p) =

n

(40)

f (s)(p; ds) as

αi (p; (ai , bi ]).

(41)

i=1

The function f (t) can be viewed as the weights with which the independent fields (p; ds) are weighted and add up to build the field X (p). It follows immediately (i)–(ii), that the integral is a Gaussian centered field with covariance ∞ from n 2 2 (s)r (p, p ; s) ds = f S i=1 αi r(ai ,bi ] (p, p ). −∞


The remainder of the construction of X (p) = f (s)(p; ds) extends it to any complex-valued function f that satisfies ∞ | f |2 (s) · r S (p, p ; s) ds < ∞, −∞

p .

for each p, This can be done using standard measure theoretic arguments that are skipped here. In particular, it follows that for the fields X and Y with corresponding f and g satisfying the condition above, we have ∞ f (s) · g(s) · r S (p, p ; s) ds. (42) r X,Y (p, p ) = Cov(X (p), Y (p )) = −∞

A.2 Partial derivative fields In this section we derive the partial mean square derivatives of the field X (p, t) defined as the stochastic integral (Eq. 6). For simplicity of presentation, we consider only the case p ∈ R2 . Generalization to higher dimensions is straightforward. In all the results it assumed without further mention that all required derivatives of stochastic fields exist in the mean square sense. Let p and p be two points in R2 with first coordinate x and x respectively. The covariance functions that are considered here can in general depend on six variables: x, y, x , y , t, t . Partial derivatives of these covariances with respect to these vari ables will be indicated by super-scripts. So for example, r x x (p, p ; s) stands for the 2 ∂ second order partial derivative ∂ x∂ x r (p, p ; s). For a field X (p, t) that depends only on three variables with generic names x, y, t, the derivatives are written in a similar manner. For example X x (p , t ) = ∂∂ Xx (x , y , t ). The second order derivative of the covariance function with respect to the first spatial coordinates x and x is given by

xx r(a,b] (p, p ) =

∂2 r(a,b] (p, p ) = ∂ x∂ x

a

b

r Sx x (p, p ; s)ds,

(43)

which is again a covariance function. It follows from the additivity of the covariance function with respect to independent fields and its correspondence to the additivity of the integral that there exists a family of Gaussian spatial fields x (p; (a, b]) centered at zero with properties analogous to those of the field (p; (a, b]) with the governing covariance r Sx x (p, p ; s). By standard arguments and under suitable regularity conditions it can be shown that X x (p) = f (s)x (p; ds). (44) To extend this to calculus of partial derivatives for a process X (p, t), basic facts about mean square derivatives can be employed. A detailed treatment of this can be found in Adler and Taylor (2007). Here we just present some basic principles and

A. Baxevani et al.

resulting formulas. For a field X (t), t = (t1 , . . . , tn ) ∈ Rn the mean square partial derivatives are defined as X ti =

X (t + h · ei ) − X (t) ∂X (t) = lim , h→0 ∂ti h

(45)

where ei is the vector with the ith element 1 and all others zero while convergence is in the mean square sense. Derivatives of higher order are defined in an analogous way. It is straightforward to see that the covariance function of such partial derivatives must be given by ∂ k X (s) ∂ k X (t) ∂ 2k r (s, t) Cov , , (46) = ∂si1 . . . ∂sik ∂ti1 . . . ∂tik ∂si1 ∂ti1 . . . ∂sik ∂tik where r (s, t) = Cov(X (s), X (t)). Lemma 1 For a stationary f ield X (t), the following relations between the f ield and its mean square derivatives (when they exist) are true. (i) The f ield and its f irst order partial derivatives, i.e., X (t) and X t j (t), when evaluated at the same point are uncorrelated. (ii) The f irst and second order partial derivatives of the f ield evaluated at t, i.e., X ti (t) and X t j tk (t), are uncorrelated for any ti , t j and tk . (iii) If additionally the f ield is isotropic, i.e. the covariance function depends only on the Euclidean length |t| of the vector t so we can write r (t) = r (|t|), then the f irst order derivatives of the f ield X ti and X t j , i = j, are uncorrelated. Proof For a stationary field X (t) with (mean square) derivatives of orders α + β and γ + δ for α, β, γ , δ ∈ {0, 1, 2, . . . }, Formula 46 takes the equivalent form α+β α+β+γ +δ r (t) ∂ X (t) ∂ γ +δ X (t) α+β ∂ Cov , |t=0 . (47) = (−1) ∂ α ti ∂ β t j ∂ γ tk ∂ δ tl ∂ α ti ∂ β t j ∂ γ tk ∂ δ tl Remember that stationarity implies r (t) = r (−t), (with some abuse in notation) which in turn means that all odd ordered derivatives of the covariance r are identically zero. Hence it follows from Eq. 47 for β = γ = δ = 0 and α = 1 that X (t) and X t j (t) are uncorrelated for every j and all t, since the first order derivative of the covariance at zero equals zero. Also for β = γ = δ = 1 and α = 0, we obtain X ti (t) and X t j tk (t) are uncorrelated for all i, j, k and every t. This proves statements (i) and (ii). To establish (iii) note that it is enough to consider only the two dimensional case t = (t1 , t2 ). If a field is stationary and isotropic then the spectral measure S is also isotropic, i.e. if Rφ is the rotation by an angle φ, then S = S ◦ Rφ . The covariance between X t1 and X t2 is given by t1 t2 r (t) = ω1 ω2 eiω·t d S(ω). R2


Thus for ω = Rπ/2 ω, we have r t1 t2 (0) = =

ω1 ω2 d S(ω)

R2

ω1 ω2 >0

ω1 ω2 d S(ω) +

=−

ω1 ω2