György Terdik

A new covariance function for
spatio-temporal data analysis with
application to atmospheric pollution and
sensor networking
György Terdik and Subba Rao Tata
UofD, HU & UofM, UK
January 30, 2015
Laboratoire des Signaux et Systemes, Supelec,
Paris, Gif-sur-Yvette
Gy. Terdik & S. Rao T. (UofD, HU & UofM, UK ) Spatio-Temporal Fields
January 30, 2015 Laboratoire des Signaux et S
/ 19

View Slide

Spatio-temporal observations
New York City air pollution data, PM2.5 measure is one
of six primary air pollutants and is a mixture of …ne
particles and gaseous compounds such as sulphur dioxide
(SO2) and nitrogen oxides (NOx)
/ 19

View Slide

In 2002, every 3 days and during the …rst 9 months, 91
equally spaced days, observed at 15 monitoring stations,
X (s, t), see [SM05],
fZ (s, t) = ∆t
X (s, t) : (s, t) : s 2 R2, t 2 Zg
0 10 20 30 40 50 60 70 80 90 100
0
10
20
30
40
50
60
70
80
90
Days
PM2.5
/ 19

View Slide

Missing data: mean by locations and …xed in time, but
3/4 of data at Location #11 is missing, namely from
days 24th to 91st, location s0
, fZ(s0
, t); t = 1, 2, 3...ng.
Sample fZ (si
, t) ; si
= 1, 2, ...m; t = 1, . . . , ng. Bronx,
Brooklyn, Manhattan, Queens, Staten Island
-74.15 -74.1 -74.05 -74 -73.95 -73.9 -73.85 -73.8
40.55
40.6
40.65
40.7
40.75
40.8
40.85
40.9
40.95
X = -73.96
Y = 40.75
Latitudes
Longitudes
/ 19

View Slide

-74.2 -74.15 -74.1 -74.05 -74 -73.95 -73.9 -73.85 -73.8 -73.75 -73.7
40.5
40.55
40.6
40.65
40.7
40.75
40.8
40.85
40.9
40.95
41
X = -73.84
Y = 40.77
/ 19

View Slide

We assume that the random process is spatially and
temporally second order stationary, homogeneous and
isotropic, i.e. E [Z (s, t)] = µ, Var [Z (s, t)] = σ2
Z
< ∞
Cov [Z (s, t) , Z (s + h, t + u)] = c (h, u) = c (khk , u) ,
We note that c (h, 0) and c (0, u) correspond to the
purely spatial and purely temporal covariances.
Spatio-temporal variogram for fZ (h, t)g
2γ (h, u) = Var [Z (s + h, t + u) Z (s, t)] .
2γ (h, u) = 2 [c (0, 0) c (h, u)] ,
for an isotropic process, γ (h, u) = γ (khk , u).
/ 19

View Slide

DFT
Frequency domain in time and Spatial domain in
space, fZ (si
, t) ; t = 1, . . . , ng
DFT at ωk
= 2πk
n
, k = 0, 1, . . . , n
2
.
Jsi
(ωk
) =
1
p
2πn
n
∑
t=1
Z (si
, t) e itωk ,
periodogram
Isi
(ωk
) = jJsi
(ωk
)j2 .
/ 19

View Slide

Jsi
(ωk
) asymptotically independent and Gaussian
E (Jsi
(ωk
)) = 0
Var (Jsi
(ωk
)) = E (Isi
(ωk
)) ' gsi
(ωk
) .
Cov Jsi
(ωk
) , Jsi
ωk0
' 0, k 6= k0.
Isotropy
gsi
(ωk
) = g (ωk
) , for all i
/ 19

View Slide

Fourier T and Spectral Repr
Introduce white noise
Js,e
(ω) =
Z
eis λ
r
n
2π
dZe
(λ, ω) .
Laplacian
∂2
∂s2
1
+ ∂2
∂s2
2
jc (ω)j2
ν
Js
(ω) = Js,e
(ω) .
λ2
1 λ2
2
jc (ω)j2 ν
dZz
(λ, ω) = dZe
(λ, ω) ,
/ 19

View Slide

Spectral density
fz
(λ, ω) = σ2
e
(2π)2
λ2
1
+ λ2
2
+ jc (ω)j2
2ν
Covariance
Cov (Js
(ω) , Js+h
(ω))
= σ2
e
2π
khk
2 jc (ω)j
2ν 1 K2ν 1
(jc (ω)j khk)
Γ(2ν)
K2ν 1
: modi…ed Bessel function of the second kind of
order 2ν 1.
/ 19

View Slide

General, Rd, ν 2 integer, correlation function is given
by
ρ (khk , ω) =
(khk jc (ω)j)2ν d
2
22ν d
2
1Γ 2ν d
2
K
2ν d
2
(khk jc (ω)j) ,
and
C (0, ω) = σ2
e
(2π)d
2 2d
2 jc (ω)j2
2ν d
2
Γ 2ν d
2
Γ (2ν)
= g (ω) .
/ 19

View Slide

Parameters
C (0, ω) = σ2
e
2 (2ν 1) jc (ω)j2
2ν 1
= g (ω) .
jc (ω)j2 and common spectral density g (ω), ARMA,
FARMA etc..
Estimation of parameters
/ 19

View Slide

Correlation
ρ (khk , ω) =
C (khk , ω)
C (0, ω)
=
1
22ν 2Γ (2ν 1)
(khk jc (ω)j)2ν 1 K2ν 1
(jc (ω)j khk) .
Special case ν = 1,
ρ (khk , ω) = khk jc (ω)j K1
(khk jc (ω)j) .
/ 19

View Slide

Spatio-temporal prediction
Js0
(ω) =
1
p
(2πn)
n
∑
t=1
Z (s0
, t) e itω,
and by inversion, we have
Z (s0
, t) =
r
n
2π
π
Z
π
Js0
(ω) eitωdω.
In other words given fJs0
(ω) , for all π ω πg,
we can uniquely recover the sequence
fZ (s0
, t) ; t = 1, . . . , ng.
/ 19

View Slide

Given
J0
m
(ω) = [Js1
(ω) , Js2
(ω) , . . . , Jsm
(ω)] .
We note
E [Jm
(ω)] = 0
E [Jm
(ω) Jm
(ω)] = Fm
(ω) ,
Fm
(ω) = (C (ksi
sj
k , ω) ; i, j = 1, 2, . . . , m), and
each element C (ksi
sj
k , ω) is given above.
/ 19

View Slide

J0
m+1
(ω) = [J0
(ω) , J0
m
(ω)] ,
which has zero mean, and variance covariance matrix
E [Jm+1
(ω) Jm+1
(ω)]
=
C0
(0, ω) E (J0
(ω) J 0
m
(ω))
E (Jm
(ω) J0
(ω)) E (Jm
(ω) Jm
(ω))
=
C0
(0, ω) G0
0
(ω)
G0
(ω) Fm
(ω)
,
/ 19

View Slide

C0
(0, ω) = E (J0
(ω) J0
(ω)) = C (0, ω),
G0
0
(ω) = E [J0
(ω) J 0
m
(ω)]
= [C (ks0
s1
k , ω) , . . . , C (ks0
sm
k , ω)]
E [J0
(ω) jJm
(ω)] = G0
0
(ω) F 1
m
(ω) Jm
(ω)
the minimum mean square error
σ2
m
(ω) = C (0, ω) G0
0
(ω) F 1
m
(ω) G0
(ω) .
ˆ
J0
(ω) = ˆ
G0
0
(ω) ˆ
F 1
m
(ω) Jm
(ω) .
/ 19

View Slide

Air pollution data, common spectral density g (ω),
ARMA(1,1), parameters estimation by Whittle method,
see [ST13] for details
0 20 40 60 80 100
-20
0
20
40
60
80
100
120
Predicted
Measured & means
/ 19

View Slide

-73.9
-73.8
-73.7
-73.6
-73.5
-73.4
40.5
40.6
40.7
40.8
40.9
41
2
4
6
8
10
12
/ 19

View Slide

S. K. Sahu and K. V. Mardia.
A bayesian kriged kalman model for short-term
forecasting of air pollution levels.
Journal of the Royal Statistical Society: Series C
(Applied Statistics) 54(1), 223–244 (2005).
T. Subba Rao and Gy. Terdik.
A space-time covariance function for spatio-temporal
random processes and spatio-temporal prediction
(kriging).
ArXiv e-prints (November 2013).
/ 19

View Slide

György Terdik

György Terdik

S³ Seminar

More Decks by S³ Seminar

Other Decks in Research

Featured

Transcript