Global processes Problems such as global warming require modeling of processes that take place on the globe (an oriented sphere). Optimal prediction of quantities such as global mean temperature need models for global covariances. Note: spherical covariances can take values in [-1,1]not just imbedded in R3. Also, stationarity and isotropy are identical concepts on the sphere.
Isotropic covariances on the sphere Isotropic covariances on a sphere are of the form C(p,q) = a iPi (cos pq ) i= 0 where p and q are directions, pq the angle between them, and Pi the Legendre polynomials. Example: ai=(2i+1)i 1 2 C(p,q) =
2 1 1 2cos pq + Global temperature Global Historical Climatology Network 7280 stations with at least 10 years of data. Subset with 839 stations with data 1950-1991 selected. Isotropic correlations The Fourier transform g :R d R G() =F () = (s)e xp(i T s)d s
g(s) =F (G ) = 1 1 d (2p ) T e xp(-i s)G ()d
Properties of Fourier transforms Convolution F (f ) =F (f)F () Scaling 1 F (f(ag)) = F( / a) a Translation F (f(gb)) =e xp(ib)F (f) Parcevals theorem 2
2 f(s) ds = F() d Relates space integration to frequency integration. Decomposes variability. Aliasing d Observe field at lattice Z of spacing . Since T m T
exp(i k) = e xp(i + 2p k) T = e xp(i Tk)e xp(i2pm Tk) the frequencies and =+2pm/ are aliases of each other, and indistinguishable. The highest distinguishable frequency is p/, the Nyquist frequency.
Illustration of aliasing Aliasing applet Spectral representation Stationary processes Z(s) = T exp(is )d Y() Rd
Spectral process Y has stationary increments 2 E dY() =d F() If F has a density f, it is called the spectral density. i(s1-s2 )T Cov(Z(s 1 ),Z(s 2 )) = e f()d R2 Estimating the spectrum For process observed on nxn grid, estimate spectrum by periodogram
2 1 i T j In,n () = z(j)e 2 (2pn) jJ 2 2j = ; J = { (n 1) / 2 ,...,n (n 1) / 2 } n Equivalent to DFT of sample covariance
Properties of the periodogram Periodogram values at Fourier frequencies (j,k)p/ are uncorrelated asymptotically unbiased not consistent To get a consistent estimate of the spectrum, smooth over nearby frequencies Some common isotropic spectra
Squared exponential s2 f()= e xp( 2pa 2 / 4a) 2 C() =s 2 e xp(a ) Matrn
2 f() =j(a 2 + )n1 pj(a )n K n (a ) C() = n1 2 G(n + 1)a 2 n A simulated process js1 ks2 Z(s) = jk cos 2p + + U jk
m n j=0 k=15 15 15 g jk =e xp( j+ 6 k tan(20 ) ) Thetford canopy heights 39-year thinned commercial plantation of Scots pine in Thetford Forest, UK
Density 1000 trees/ha 36m x 120m area surveyed for crown height Focus on 32 x 32 subset Spectrum of canopy heights Whittle likelihood Approximation to Gaussian likelihood using periodogram: IN,N () l () = lo f(; ) +
f(; ) where the sum is over Fourier frequencies, avoiding 0, and f is the spectral density Takes O(N logN) operations to calculate instead of O(N3). Using non-gridded data Consider Y(x) = 2 h(x s)Z(s)d s where
h(x) =1( x i / 2, i =1, 2) Then Y is stationary with spectral density 1 2 fY () = 2 H() fZ () Viewing Y as a lattice process, it has spectral density 2 2p 2p f,Y () = H( +
) fZ ( + ) Z2 Estimation 1 Let Yn2 (x) = nx h(s x)Z(s ) i
i iJx where Jx is the grid square with center x and nx is the number of sites in the square. Define the tapered periodogram 2 1 ix T Ig1Y 2 () = 1 (x)Yn 2 (x)e 2
n 1 (x) where g1 (x) =n x / n. The Whittle likelihood is approximately I1 ,Y 2 (2pj/ n) n L = l o f
(2pj/ n) + ,Y Y (2p )2 f (2pj/ n) j ,Y
n2 A simulated example Estimated variogram Evidence of anisotropy 15o red 60o green 105o blue 150o brown QuickTime and a
TIFF (Uncompressed) decompressor are needed to see this picture. Another view of anisotropy s e2 = 127.1(259) s e2 = 154.6 (134) s2 = 68.8 (255) = 10.7 (45) s2 = 141.0 (127) = 29.5 (35)
Geometric anisotropy Recall that if C(x,y) = C( x y ) we have an isotropic covariance (circular isocorrelation curves). If C(x,y) = C( Ax Ay ) for a linear transformation A, we have geometric anisotropy (elliptical isocorrelation curves). General nonstationary correlation structures are typically locally geometrically anisotropic. Lindgren & Rychlik
transformation x = (2x + y + 109.15) / 2 y = 4(x + 2y 154.5) / 3 QuickTime and a TIFF (Uncompressed) decompressor are needed to see this picture.