Data collected on a rectangular lattice occur frequently in many areas such as field trials, geostatistics, remotely sensed data, and image analysis. Models for the spatial process often make simplifying assumptions, including axial symmetry and separability. We consider methods for testing these assumptions and compare tests based on sample covariances, tests based on the sample spectrum, and model-based tests

L. SCACCIA

MARTIN R.J

SCACCIA, LUISA

MARTIN R. J.

Archivio istituzionale della ricerca - Università di Macerata

1Testing for Simplification in Spatial ModelsL. Scaccia1 and R. J. Martin21 Dipartimento di Scienze Statistiche, Universita` di Perugia,06100 Perugia, Italy2 Department of Probability and Statistics, University of Sheffield,S3 7RH Sheffield, UKAbstract. Data collected on a rectangular lattice occur frequently in manyareas such as field trials, geostatistics, remotely sensed data, and image analy-sis. Models for the spatial process often make simplifying assumptions, includ-ing axial symmetry and separability. We consider methods for testing theseassumptions and compare tests based on sample covariances, tests based onthe sample spectrum, and model-based tests.Keywords. Autoregressive process, axial symmetry, doubly-geometric pro-cess, lattice process, separability, spatial process.1 IntroductionData collected on a regular two-dimensional lattice arise in many areas. Mostmodels used for such data assume axial symmetry, but can still be difficult touse. Many of the problems of two-dimensional modelling can be overcome byusing separable processes. Under this assumption, time series methods canbe used to analyze and model lattice data (Martin, 1990).However, no tests for separability and axial symmetry are known apart fromGuo and Billard (1998), who suggested the Wald test to test separability. Thismodel-based test compares the fit of an AR(1)·AR(1) process and the moregeneral Pickard process to the data, and actually tests for axial symmetryand separability together.We propose and compare some different model-free and model-based testsfor testing axial symmetry and separability. The paper is structured as fol-lows. In section 2 we give some notation and definitions, illustrate the simu-lation method and present the different models used. The different tests foraxial symmetry and separability are discussed in section 3. In section 4 wecompare the performances of the tests.2 Definitions, models and simulation method2.1 NotationWe assume that data occur on an n1 by n2 rectangular lattice, with rowsindexed by i1 = 1, . . . , n1 and columns by i2 = 1, . . . , n2. Row and columnlags are g1 and g2, with gj = −(nj − 1), . . . , 0, . . . , (nj − 1), for j = 1, 2.The sites are ordered lexicographically, so that (i1, i2) precedes (i1, i2 + 1)for i2 < n2, and (i1, n2) precedes (i1 + 1, 1). Data can be considered as arealization of random variables Yi1,i2 . Let the vector Y contain the Yi1,i2 insite order and assume a constant mean. We assume Y ∼ N(0, V σ2) whereV = V (α) is an n by n positive definite matrix depending on the q-vector ofparameters α.2Assuming second-order stationarity, C(g1, g2) = Cov (Yi1,i2 , Yi1+g1,i2+g2) isthe covariance at lags g1, g2, with C(0, 0) = σ2. The spectrum is f(ω1, ω2) =∑∞g1=−∞∑∞g2=−∞ C(g1, g2) cos(g1ω1 + g2ω2)/(2pi)2 (Priestley, 1981, section9.7). Note that C(g1, g2) = C(−g1,−g2) and f(ω1, ω2) = f(−ω1,−ω2) al-ways.2.2 Axial symmetry and separabilityFor an axially or reflection symmetric process C(g1, g2) = C(g1,−g2),∀g1, g2and, equivalently, f(ω1, ω2) = f(ω1,−ω2),∀ω1, ω2. This means that the co-variances and the spectrum are both symmetric about the axes.For a process to be separable, C(g1, g2) ∝ C(g1, 0) · C(0, g2),∀ g1, g2 andf(ω1, ω2) ∝ f(ω1, 0) · f(0, ω2),∀ ω1, ω2 are also required. Therefore, the co-variances and the spectrum are determined, up to a multiplicative constant,by the margins. Clearly separability implies axial symmetry.It is often desirable in practice that a process should be axially symmetric,and that the covariances should have a simple form. Separable processes haveboth these properties, and are still flexible enough to provide a reasonablerepresentation of many planar structures (Martin, 1990).For a separable process, the V matrix can be expressed as V = Vx ⊗ Vz,where Vx and Vz are two smaller matrices that arise from the two underlyingone-dimensional processes, and this implies that its determinant and inverse(required for exact Gaussian maximum likelihood, generalised least-squaresestimation and exact Gaussian simulation) are easily determinable, which isa major advantage of this subclass of processes (Martin, 1996).2.3 Models simulated and simulation methodWe compared the different tests proposed by simulation: separable (or axiallysymmetric) processes were used to simulate the null distribution of tests forseparability (or axial symmetry), while non-separable (or non-symmetric)processes were used to simulate the distribution of the tests under the al-ternative hypothesis. We indicate tests for axial symmetry by S and testsfor separability by R, using different superscripts for their basis: (c) samplecovariances, (p) sample spectrum, and (m) model-based.The separable process we used is the AR(1)·AR(1) process (Martin, 1979):Yi1,i2 = α1Yi1−1,i2 + α2Yi1,i2−1 − α1α2Yi1−1,i2−1 + i1,i2where the i1,i2 are assumed to be independently distributed as N(0, σ2 ).As an axially symmetric, non-separable model we used a particular case ofa second-order conditional autoregressive process, the CAR(2)SD with sym-metric diagonal term (Balram and Moura, 1993). It can be written as:E(Yi1,i2 |·) = β1(Yi1−1,i2 + Yi1+1,i2) + β2(Yi1,i2−1 + Yi1,i2+1)+ β3(Yi1−1,i2−1 + Yi1+1,i2+1 + Yi1−1,i2+1 + Yi1+1,i2−1)with constant conditional variance. Conditioning is on Yl1,l2 , ∀(l1, l2) 6= (i1, i2).A CAR(2)SD with β3 = −β1β2, β1 = α1/(1 + α21) and β2 = α2/(1 + α22) hasthe same covariance structure as an AR(1)·AR(1) with parameters α1, α2.The non-separable, non-axially-symmetric model used, referred to as thePickard process (Pickard, 1980, Tory and Pickard 1992) can be written as:Yi1,i2 = α1Yi1−1,i2 + α2Yi1,i2−1 + α3Yi1−1,i2−1 + i1,i2 .3The AR(1)·AR(1) is a special case of the Pickard process with α3 = −α1α2.Each process was simulated 1 000 times to estimate the distribution ofthe test statistics under the null and the alternative hypotheses. Choosing amatrix T such that V = TT ′, the observation vector y can be simulated asy = T where  is a random vector of n independent N(0, 1) observations.3 Testing axial symmetry and separability3.1 Tests based on sample covariancesSince axial symmetry and separability are usually defined in terms of thecovariance structure, it seems natural to try and test these two hypothe-ses using the sample covariances. The sample estimator of the covarianceswe used is c(g1, g2) =∑n1−g1i1=1∑n2−g2i2=1(Yi1,i2 − Y )(Yi1+g1,i2+g2 − Y )/n andc(g1,−g2) =∑n1−g1i1=1∑n2−g2i2=1(Yi1,i2+g2 − Y )(Yi1+g1,i2 − Y )/n, for g1, g2 ≥ 0,where Y =∑n1i1=1∑n2i2=1Yi1,i2/n. Writing c(g1, g2) as a quadratic form in Yand using standard results (Cressie, 1993 section 2.4.2), it is easy, given V ,to calculate cov[c(g1, g2), c(g?1 , g?2)], from which it is clear that the c(g1, g2)are highly dependent and have non-constant variance.3.1.1 Testing axial symmetryWe can consider basing tests for axial symmetry on the sample differencesF (g1, g2) = c(g1, g2) − c(g1,−g2) for g1, g2 ≥ 1. Under the null hypothesisE[F (g1, g2)] = 0 ∀ g1, g2. However the F (g1, g2) have non-constant varianceand they are correlated. For 1 ≤ gj ≤ aj , let F denote the vector containingthe F (g1, g2) in lexicographic order, and let cov(F ) = Σ. Since the c(g1, g2)are asymptotically normally distributed, F ∼ N(0, Σ) and a possible test usesS(c) = F ′Σˆ−1F.Under axial symmetry, S(c) ∼ χ2a, asymptotically, with a = a1 ·a2. Notice thatestimating Σ requires estimating V and, therefore, this test is not completelymodel-free. We fitted the AR(1)·AR(1) to estimate V .3.1.2 Testing separabilityModjeska and Rawlings (1983) suggested testing separability by using thesingular value decomposition of the matrix M of sample covariances c(g1, g2)for −(a1 − 1) ≤ g1 ≤ a1 − 1 and 0 ≤ g2 ≤ a2 − 1. Since C(g1, g2) ∝C(g1, 0) ·C(0, g2) under the null hypothesis, the population M has rank onefor a separable process. If λ1 is the largest singular value of M , in absolutevalue, a test for separability could use:R(c) = λ21/∑λ2i .Its distribution needs to be simulated. Note that R(c) ∈ [0, 1], low valuesleading to rejection of separability.3.2 Tests based on sample spectrumThe main advantage of the sample spectrum is that it has better asymptoticproperties than the sample covariances, being asymptotically independent atdifferent frequencies. For ωj a multiple of 2pi/nj , the sample spectrum weused is I(ω1, ω2) =∑n1−1g1=−n1+1∑n2−1g2=−n2+1 c(g1, g2) cos(g1ω1 + g2ω2)/(2pi)2.4Then, for ωj 6= 0, pi, we have I(ω1, ω2)/f(ω1, ω2)→ i.i.d. Exp(1), E[I(ω1, ω2)]→ f(ω1, ω2) and Var[I(ω1, ω2)]→ f2(ω1, ω2), as n1, n2 →∞. Exact formulafor the mean and the variance of the one-dimensional sample spectrum arein Priestley (1981, section 6.1.3).3.2.1 Testing axial symmetryConsider the differences G(ω1, ω2) = I(ω1, ω2) − I(ω1,−ω2) (with ωj =2pikj/nj , kj = 1, 2, . . . , bj , and bj = (nj−1)/2 if nj is odd, or bj = nj/2−1 ifnj is even). Clearly under axial symmetry E[G(ω1, ω2)] = 0 and Var[G(ω1, ω2)]→ 2f2(ω1, ω2) as n1, n2 →∞. Since the G(ω1, ω2) do not have constant vari-ance, we considered two modifications. Firstly, taking the logarithm of theI(ω1, ω2) approximately stabilizes the variance and also reduces the non-normality. For ωj 6= 0, pi we then have, asymptotically log[I(ω1, ω2)] ∼Gumbel (log[f(ω1, ω2)], 1) from which Var {log[I(ω1, ω2)]} → pi2/6 as n1, n2→ ∞. Given b = b1 · b2, a possible test uses:S(p)1 = D√b/√pi2/3,where D =∑ω1,ω2D(ω1, ω2)/b, D(ω1, ω2) = log[I(ω1, ω2)]− log[I(ω1,−ω2)].Assuming b large enough for the central limit theorem to hold, under ax-ial symmetry S(p)1 ∼ N(0, 1). Other tests, for example, the sign test or theWilcoxon test on the differences D(ω1, ω2) or G(ω1, ω2), can also be used.A second modification is to standardize the differences G(ω1, ω2), esti-mating f(ω1, ω2) by [I(ω1, ω2) + I(ω1,−ω2)]/2. Thus, if H(ω1, ω2) denotesG(ω1, ω2)/[I(ω1, ω2) + I(ω1,−ω2)], then another possible test usesS(p)2 =√3b H.If the I(ω1, ω2) are independent and exponentially distributed, H(ω1, ω2) hasa uniform distribution in the interval [−1, 1], and the sum converges rapidlyto a N(0, b/3), so S(p)2 ∼ N(0, 1) asymptotically under axial symmetry.3.2.2 Testing separabilityThe same idea in section 3.1.2 can be applied to a matrix of I(ω1, ω2) to testseparability (we use R(p)1 to denote this test statistic).As an alternative, under separability, log[I(ω1, ω2)] and log[I(ω1,−ω2)] canbe regarded as two sample realizations of the same value of the log spectrumlog[f(ω1, ω2)] = d + log[f(ω1, 0)] + log[f(0, ω2)] where d is a constant. Sotesting for separability reduces to testing for lack of interaction in a two-wayclassification table with two realizations in each cell and, var(log[I(ω1, ω2)])being approximately constant, this can be done using the statistic:R(p)2 = SSrc/SSewhere SSrc is the mean of squares due to the interaction, and SSe is theresidual mean of squares. Under separability, and assuming the log differencesare approximately normal, R(p)2 ∼ F(b1−1)(b2−1),b, asymptotically.53.3 Model-based testsAssuming a specific model behind the data allows the use of model-based testsfor axial symmetry and separability. The idea is to fit a particular unrestrictedmodel (non-separable or non-symmetric) to the data, then to restrict it tobe axially symmetric or separable, imposing constraints on its parametersand to test which model is more appropriate. The comparison can use thegeneralised likelihood ratio test (GLRT) or tests which are asymptoticallyequivalent, such as the Wald test or the Score test.3.3.1 Testing axial symmetry and separability togetherUnder the assumption Y ∼ N(0, σ2), the log-likelihood can be written as`(α, σ2; y) = −(n/2) log(2pi) − (n/2) log(σ2) − log | V | /2 − y′V −1y/(2σ2).The GLRT statistic is thenS(m) = 2[`(αˆ, σˆ2; y)− `(αˆ0, σˆ20 ; y)].Now, to test for axial symmetry and separability together we can considerthe Pickard as the unrestricted model and the AR(1)·AR(1) as the restrictedone, with the constraint α3 = −α1α2 (see section 2.3). Alternatively, theWald (Guo and Billard, 1998) and the Score tests could be used, althoughwe found that GLRT is preferred for smaller lattice sizes. Under the null texthypothesis, the three test statistics are asymptotically distributed as a χ21.3.3.2 Testing separabilitySeparability can be tested comparing the fit of a CAR(2)SD with the fit ofan AR(1)·AR(1), applying the restriction β3 = −β1β2 (see section 2.3) to theparameters of the CAR(2)SD. The GLRT can be expressed this time as:R(m) = 2[`(βˆ, σˆ2; y)− `(βˆ0, σˆ20 ; y)].4 Results and discussionFor the tests on axial symmetry, a good correspondence between the simu-lated and theoretical distributions was found for S(p)1 , S(p)2 , S(m) on 11 by 11or larger lattices. For S(c) with a1 = a2 = 4, a 15 by 15 lattice was required.Tests S(p)1 , S(p)2 are completely general, the main assumption required beingstationarity, while S(c), S(m) require partial or complete specification of amodel. On the other hand, S(m) and S(c) are more powerful (see Table 1).Note that the power partially depends on how far the Pickard process is fromthe AR(1)·AR(1), measured by |α3 + α1α2|. Also note that S(m) is testingfor separability at the same time.Separability proved to be harder to test. Tests R(c) and R(p)1 are not invari-ant, their distributions strongly depending on the parameters of the processconsidered, and also can have values very close to 1 under the alternativehypothesis. Simulations show R(p)2 to be distributed as expected on 11 by11 lattices, although this test does not seem to be very powerful. Also, itassumes axial symmetry, which therefore should be tested first. Test R(m) isvery powerful (see Table 2), although it requires specification of a model.As well as the tests mentioned, we have also considered some other tests,and considered variants of those mentioned - for example using different aj .We have evaluated the tests on other parameter values than those mentioned,6and on some other processes - for example some separable ARMA processes.We have some theoretical results to help explain the power for different pa-rameters under the alternative hypothesis. We have tried the tests on a realdata set - a SAR (synthetic aperture radar) image for which the postulatedform of the point-spread function implies separability. The theoretical distri-bution of the data is Exponential, so further work is needed to see how thetests perform for theoretical SAR data. For further details, see Scaccia (2000)and Scaccia & Martin (2002).Parameters S(p)1 S(p)2 S(c) S(m) |α3 + α1α2|α = (0.1, 0.2, 0.6) 97.4% 97.1% 99.9% 100.0% 0.62α = (0.6, 0.7,−0.8) 87.0% 86.3% 99.5% 100.0% 0.38α = (0.3, 0.4, 0.2) 55.2% 59.0% 78.1% 99.8% 0.32α = (0.1, 0.6, 0.2) 44.6% 45.8% 61.2% 99.3% 0.26α = (0.3, 0.6,−0.15) 5.1% 5.3% 5.8% 6.8% 0.03Table 1. Simulated power of tests for axial symmetry for different Pickard pro-cesses, on a 15 by 15 lattice when the level of the tests is equal to 5%.Parameters R(p)2 R(m) |β3 + β1β2|β = (0, 0, 0.25) 88.8% 100.0% 0.25β = (0.11, 0.07, 0.16) 13.6% 98.7% 0.17β = (0.2, 0.2, 0.05) 9.1% 92.7% 0.09β = (0.25, 0.25, 0) 8.3% 86.8% 0.06β = (0.3, 0.2,−0.03) 5.3% 12.3% 0.03Table 2. Simulated power of tests for separability for different CAR(2)SD processes,on a 15 by 15 lattice when the level of the tests is equal to 5%.ReferencesBalram, N. and Moura, J.M.F. (1993). Noncausal Gauss Markov randomfields: parameter structure and estimation. IEEE Trans. Inform. Theory,39, 1333-1355.Cressie, N.A.C. (1993). Statistics for Spatial Data. New York: Wiley.Guo, J.H. and Billard, L. (1998). Some inference results for causal autore-gressive processes on a plane. J. of Time Ser. An., 19, 681-691.Martin, R.J. (1979). A subclass of lattice processes applied to a problem inplanar sampling. Biometrika, 66, 209-217.Martin, R.J. (1990). The use of time-series models and methods in the analy-sis of agricultural field trials. Commun. Statist. -Theory Meth.. 19, 55-81.Martin, R.J. (1996). Some results on unilateral ARMA lattice processes. J.Statist. Plann. Inf.. 50, 395-411.Modjeska, J.S. and Rawlings, J.O. (1983). Spatial correlation analysis ofuniformity data. Biometrics, 39, 373-384.Pickard, D.K. (1980). Unilateral Markov fields. Adv. Appl. Prob., 12, 655-671.Priestley, M.B. (1981). Spectral Analysis and Time Series. London: Aca-demic Press.Scaccia, L. (2000). Testing for Simplification in Spatial Models. Ph.D. thesis,University of Perugia, Italy.Scaccia, L. and Martin R.J. (2002). Testing axial symmetry and separabilityof lattice processes. Preprint.Tory, E.M. and Pickard, D.K. (1992). Unilateral Gaussian fields. Adv. Appl.Prob., 24, 95-112.

Testing for simplification in spatial models

Abstract

Similar works

Full text

Available Versions

Archivio istituzionale della ricerca - Università di Macerata