A Review of Nonparametric Hypothesis Tests of Isotropy Properties in Spatial Data
SSubmitted to Statistical Science
A Review of NonparametricHypothesis Tests of IsotropyProperties in Spatial Data
Zachary D. Weller and Jennifer A. Hoeting
Department of Statistics, Colorado State University
Abstract.
An important aspect of modeling spatially-referenced data isappropriately specifying the covariance function of the random field. Apractitioner working with spatial data is presented a number of choicesregarding the structure of the dependence between observations. One ofthese choices is determining whether or not an isotropic covariance func-tion is appropriate. Isotropy implies that spatial dependence does not de-pend on the direction of the spatial separation between sampling locations.Misspecification of isotropy properties (directional dependence) can leadto misleading inferences, e.g., inaccurate predictions and parameter es-timates. A researcher may use graphical diagnostics, such as directionalsample variograms, to decide whether the assumption of isotropy is rea-sonable. These graphical techniques can be difficult to assess, open to sub-jective interpretations, and misleading. Hypothesis tests of the assumptionof isotropy may be more desirable. To this end, a number of tests of direc-tional dependence have been developed using both the spatial and spectralrepresentations of random fields. We provide an overview of nonparametricmethods available to test the hypotheses of isotropy and symmetry in spa-tial data. We summarize test properties, discuss important considerationsand recommendations in choosing and implementing a test, compare sev-eral of the methods via a simulation study, and propose a number of openresearch questions. Several of the reviewed methods can be implementedin R using our package spTest , available on CRAN . MSC 2010 subject classifications:
Primary 62M30, ; secondary 62G10.
Key words and phrases: isotropy, symmetry, nonparametric spatial covari-ance.
1. INTRODUCTION
Early spatial models relied on the simplifying assumptions that the covariancefunction is stationary and isotropic. With the emergence of new sources of spatial (e-mail: [email protected]) (e-mail: [email protected]) ∗ Weller’s work was supported by the National Science Foundation Research Network onStatistics in the Atmospheric and Ocean Sciences (STATMOS) (DMS-1106862). Hoeting’s re-search was supported by the National Science Foundation (AGS-1419558). † The authors would like to thank Peter Guttorp and Alexandra Schmidt for organizing thePan-American Advanced Study Institute (PASI) on spatiotemporal statistics in June 2014 whichinspired this work. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 th November, 2015 a r X i v : . [ s t a t . M E ] N ov ZACHARY D. WELLER data, for instance, remote sensing via satellite, climate model output, or environ-mental monitoring, a variety of methods and models have been developed thatrelax these assumptions. In the case of anisotropy, there are a number of methodsfor modeling both zonal anisotropy (Journel and Huijbregts, 1978, pg. 179-184;Ecker and Gelfand, 2003; Schabenberger and Gotway, 2004, pg. 152; Banerjeeet al., 2014, pg. 31) and geometric anisotropy (Borgman and Chao, 1994; Eckerand Gelfand, 1999). Rapid growth of computing power has allowed the imple-mentation and estimation of these models.When modeling a spatial process, the specification of the covariance functionwill have an effect on kriging and parameter estimates and the associated uncer-tainty (Cressie, 1993, pg. 127-135). Sherman (2011, pg. 87-90) and Guan et al.(2004) use numerical examples to demonstrate the adverse implications of incor-rectly specifying isotropy properties on kriging estimates. Given the variety ofchoices available regarding the properties of the covariance function (e.g., para-metric forms, isotropy, stationarity) and the effect these choices can have on in-ference, practitioners may seek methods to inform the selection of an appropriatecovariance model.A number of graphical diagnostics have been proposed to determine isotropyproperties. Perhaps the most commonly used methods are directional semivari-ograms and rose diagrams (Matheron, 1961; Isaaks and Srivastava, 1989, pg. 149-154). Banerjee et al. (2014, pg. 38-40) suggest using an empirical semivariogramcontour plot to assess isotropy as a more informative method than directionalsample semivariograms. Another technique involves comparing empirical esti-mates of the covariance at different directional lags to assess symmetry for dataon gridded sampling locations (Modjeska and Rawlings, 1983). One caveat of theaforementioned methods is that they can be challenging to assess, are open tosubjective interpretations, and can be misleading (Guan et al., 2004) becausethey typically do not include a measure of uncertainty. Experienced statisticiansmay have intuition about the interpretation and reliability of these diagnostics,but a novice user may wish to evaluate assumptions via a hypothesis test.Statistical hypothesis tests of second order properties can be used to supple-ment and reinforce intuition about graphical diagnostics and can be more objec-tive. Like the graphical techniques, hypothesis tests have their own caveats; forexample, a parametric test of isotropy demands specification of the covariancefunction. A nonparametric method for testing isotropy avoids the potential prob-lems of misspecification of the covariance function and the requirement of modelestimation under both the null and alternative hypothesis, which can be compu-tationally expensive for large datasets. Furthermore, nonparametric methods donot require the common assumption of geometric anisotropy. However, in aban-doning the parametric assumptions about the covariance function, implementinga test of isotropy presents several challenges (see Section 5). A nonparametric testof isotropy or symmetry can serve as another form of exploratory data analysisthat supplements graphical techniques and informs the choice of an appropriatenonparametric or parametric model. Figure 1 illustrates the process for assessingand modeling isotropy properties.In this article we review nonparametric hypothesis tests developed to test theassumptions of symmetry and isotropy in spatial processes. We summarize testsin both the spatial and spectral domain and provide tables that enable convenient imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 th November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY InvestigatingIsotropyPropertiesGraphical Techniques
Hypothesis Tests
EmpiricalSemivariogramContour Plot Rose Diagram DirectionalSemivariograms Parametric
Nonparametric
More Assessment ofIsotropy Needed? Y e s Choice ofNonparametric TestSee Figure 2
Model forCovarianceFunctionNonparametricCovarianceFunction ParametricCovarianceFunction D e t e r m i n e I s o t r o p y P r o p e r t i e s N o ParametricCovarianceFunctionIsotropic AnisotropicIsotropic Anisotropic GeometricAnisotropy ZonalAnisotropyModelEstimationAssessment ofCandidate Models S p e c t r a l D o m a i n S p a t i a l D o m a i n Type ofAnisotropyModel EstimationUnder H and H G e o m e t r i c Z o n a l Model ComparisonUnder H and H e.g., LRT S p e c t r a l D o m a i n S p a t i a l D o m a i n ChooseFinalModel
Fig 1: A flow chart illustrating the process of determining and modeling isotropyin spatial data. The gray boxes indicate the focus of this paper. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Fig 1: A flow chart illustrating the process of determining and modeling isotropyin spatial data. The gray boxes indicate the focus of this paper. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ZACHARY D. WELLER comparisons of test properties. A simulation study evaluates the empirical sizeand power of several of the methods and enables a direct comparison of methodperformance. The simulations also lead to new insights into test performanceand implementation beyond those given in the original works. Finally, we includegraphics that demonstrate considerations for choosing a nonparametric test andillustrate the process of determining isotropy properties.The remainder of this article is organized as follows: Section 2 establishes nota-tion and definitions; Section 3 details the various nonparametric hypothesis testsof isotropy and symmetry and includes Tables 1-3 which facilitate comparisonbetween tests as well as test selection for users; Section 4 describes the simula-tion study comparing the various methods; Sections 5 and 6 provide discussionand conclusions. Additional details on the simulation study are specified in theAppendix.
2. NOTATION AND DEFINITIONS
Here we briefly review key definitions required for tests of isotropy. For addi-tional background, see Schabenberger and Gotway (2004). Let { Y ( s ) : s ∈ D ⊆ R d , d > } be a second order stationary random field (RF). Below we will assumethat d = 2, although many of the results hold for the more general case of d > h = ( h , h ), the semivariogram function describes dependencebetween observatons, Y , at spatial locations separated by lag h and is defined as(2.1) γ ( h ) = 12 Var( Y ( s + h ) − Y ( s )) . The classical, moment-based estimator of the semivariogram (Matheron, 1962) is(2.2) (cid:98) γ ( h ) = 12 |D ( h ) | (cid:88) [ Y ( s ) − Y ( s + h )] , where the sum is over D ( h ) = { s : s , s + h ∈ D} and |D ( h ) | is the number ofelements in D ( h ). The set D ( h ) is the set of sampling location pairs that areseparated by spatial lag h . The covariance function, C ( h ), is an alternative tothe semivariogram for describing spatial dependence and is given by C ( h ) =lim || v ||→∞ γ ( v ) − γ ( h ) if the limit exists.Let { s , . . . , s n } ⊂ D be the finite set of locations at which the random pro-cess is observed, providing the random vector ( Y ( s ) , . . . , Y ( s n )) T . The samplinglocations may follow one of several spatial sampling designs: gridded locations,randomly and uniformly distributed locations, or a cluster design. Some authorsmake the distinction between a lattice process and a geostatistical process ob-served on a grid (Fuentes and Reich, 2010; Schabenberger and Gotway, 2004,pg. 6-10). We do not explore this distinction further and will use the term gridthroughout.It is often of interest to infer the effect of covariates on the process, deduce de-pendence structure, and/or predict Y with quantifiable uncertainty at new loca-tions. To achieve these goals, the practitioner must specify the distributional prop-erties of the spatial process. A common assumption is that the finite-dimensionaljoint distribution is multivariate normal (MVN), in which case we call the RF imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Here we briefly review key definitions required for tests of isotropy. For addi-tional background, see Schabenberger and Gotway (2004). Let { Y ( s ) : s ∈ D ⊆ R d , d > } be a second order stationary random field (RF). Below we will assumethat d = 2, although many of the results hold for the more general case of d > h = ( h , h ), the semivariogram function describes dependencebetween observatons, Y , at spatial locations separated by lag h and is defined as(2.1) γ ( h ) = 12 Var( Y ( s + h ) − Y ( s )) . The classical, moment-based estimator of the semivariogram (Matheron, 1962) is(2.2) (cid:98) γ ( h ) = 12 |D ( h ) | (cid:88) [ Y ( s ) − Y ( s + h )] , where the sum is over D ( h ) = { s : s , s + h ∈ D} and |D ( h ) | is the number ofelements in D ( h ). The set D ( h ) is the set of sampling location pairs that areseparated by spatial lag h . The covariance function, C ( h ), is an alternative tothe semivariogram for describing spatial dependence and is given by C ( h ) =lim || v ||→∞ γ ( v ) − γ ( h ) if the limit exists.Let { s , . . . , s n } ⊂ D be the finite set of locations at which the random pro-cess is observed, providing the random vector ( Y ( s ) , . . . , Y ( s n )) T . The samplinglocations may follow one of several spatial sampling designs: gridded locations,randomly and uniformly distributed locations, or a cluster design. Some authorsmake the distinction between a lattice process and a geostatistical process ob-served on a grid (Fuentes and Reich, 2010; Schabenberger and Gotway, 2004,pg. 6-10). We do not explore this distinction further and will use the term gridthroughout.It is often of interest to infer the effect of covariates on the process, deduce de-pendence structure, and/or predict Y with quantifiable uncertainty at new loca-tions. To achieve these goals, the practitioner must specify the distributional prop-erties of the spatial process. A common assumption is that the finite-dimensionaljoint distribution is multivariate normal (MVN), in which case we call the RF imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY a Gaussian random field (GRF). We are interested in second order properties;thus, hereafter we assume that the mean of the RF is zero.A common simplifying assumption on the spatial dependence structure is thatit is isotropic. Definition . A second order stationary spatial process is isotropic if thesemivariogram, γ ( h ) , [or equivalently, the covariance function C ( h ) ] of the spatialprocess depends on the lag vector h only through its Euclidean length, h = || h || ,i.e., γ ( h ) = γ ( h ) for some function γ ( · ) of a univariate argument. Isotropy implies that the dependence between any two observations depends onlyon the distance between their sampling locations and not on their relative orien-tation. A spatial process that is not isotropic is called anisotropic. Anisotropy isbroadly categorized as either geometric or zonal (Zimmerman, 1993). In practice,if a process is assumed to be anisotropic, it is usually assumed to be geometri-cally anisotropic due to its precise formal and functional definition (Ecker andGelfand, 1999). Geometric anisotropy is governed by a scaling parameter, R , androtation parameter, θ , and implies the anisotropy can be corrected via a lineartransformation of the lag vector or, equivalently, the sampling locations (Cressie,1993, pg. 64).Symmetry is another directional property of the covariance (semivariogram)function, which is often used to describe the spatial variation of processes on agrid. We discuss symmetry properties here as they are a subset of isotropy, andmethods for testing isotropy can often be used to test symmetry. The followingdefinitions come from Lu and Zimmerman (2005) and Scaccia and Martin (2005)where the notation C ( h , h ) denotes the covariance between random variableslocated h columns and h rows apart on a rectangular grid, denoted L . Definition . A second order stationary spatial process on a grid is reflec-tion or axially symmetric if C ( h , h ) = C ( − h , h ) for all ( h , h ) ∈ L . Definition . A second order stationary spatial process on a grid is diag-onally or laterally symmetric if C ( h , h ) = C ( h , h ) for all ( h , h ) ∈ L . Definition . A second order stationary spatial process on a grid is com-pletely symmetric if it is both reflection and laterally symmetric, i.e., C ( h , h ) = C ( − h , h ) = C ( h , h ) = C ( − h , h ) for all ( h , h ) ∈ L . Complete symmetry is a weaker property than isotropy. Isotropy requires that C ( h , h ) depends only on (cid:112) h + h for all h , h . The relationship between theseproperties is given by:(2.3) isotropy = ⇒ complete symmetry = ⇒ reflection symmetrydiagonal symmetry . Thus, rejecting a null hypothesis of reflection symmetry implies evidence againstthe assumptions of reflection symmetry, complete symmetry, and isotropy. How-ever, failure to reject a null hypothesis of reflection symmetry does not imply anassumption of complete symmetry or isotropy is appropriate. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY a Gaussian random field (GRF). We are interested in second order properties;thus, hereafter we assume that the mean of the RF is zero.A common simplifying assumption on the spatial dependence structure is thatit is isotropic. Definition . A second order stationary spatial process is isotropic if thesemivariogram, γ ( h ) , [or equivalently, the covariance function C ( h ) ] of the spatialprocess depends on the lag vector h only through its Euclidean length, h = || h || ,i.e., γ ( h ) = γ ( h ) for some function γ ( · ) of a univariate argument. Isotropy implies that the dependence between any two observations depends onlyon the distance between their sampling locations and not on their relative orien-tation. A spatial process that is not isotropic is called anisotropic. Anisotropy isbroadly categorized as either geometric or zonal (Zimmerman, 1993). In practice,if a process is assumed to be anisotropic, it is usually assumed to be geometri-cally anisotropic due to its precise formal and functional definition (Ecker andGelfand, 1999). Geometric anisotropy is governed by a scaling parameter, R , androtation parameter, θ , and implies the anisotropy can be corrected via a lineartransformation of the lag vector or, equivalently, the sampling locations (Cressie,1993, pg. 64).Symmetry is another directional property of the covariance (semivariogram)function, which is often used to describe the spatial variation of processes on agrid. We discuss symmetry properties here as they are a subset of isotropy, andmethods for testing isotropy can often be used to test symmetry. The followingdefinitions come from Lu and Zimmerman (2005) and Scaccia and Martin (2005)where the notation C ( h , h ) denotes the covariance between random variableslocated h columns and h rows apart on a rectangular grid, denoted L . Definition . A second order stationary spatial process on a grid is reflec-tion or axially symmetric if C ( h , h ) = C ( − h , h ) for all ( h , h ) ∈ L . Definition . A second order stationary spatial process on a grid is diag-onally or laterally symmetric if C ( h , h ) = C ( h , h ) for all ( h , h ) ∈ L . Definition . A second order stationary spatial process on a grid is com-pletely symmetric if it is both reflection and laterally symmetric, i.e., C ( h , h ) = C ( − h , h ) = C ( h , h ) = C ( − h , h ) for all ( h , h ) ∈ L . Complete symmetry is a weaker property than isotropy. Isotropy requires that C ( h , h ) depends only on (cid:112) h + h for all h , h . The relationship between theseproperties is given by:(2.3) isotropy = ⇒ complete symmetry = ⇒ reflection symmetrydiagonal symmetry . Thus, rejecting a null hypothesis of reflection symmetry implies evidence againstthe assumptions of reflection symmetry, complete symmetry, and isotropy. How-ever, failure to reject a null hypothesis of reflection symmetry does not imply anassumption of complete symmetry or isotropy is appropriate. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ZACHARY D. WELLER
The aforementioned properties of isotropy and symmetry were defined in termsof examining the spatial random process in the spatial domain, where secondorder properties depend on the spatial separation, h . Alternatively, a spatialrandom process can be represented in the spectral domain using Fourier analysis.For the purposes of investigating second order properties, we are interested inthe spectral representation of the covariance function, called the spectral densityfunction and denoted f ( ω ), where ω = ( ω , ω ). Under certain conditions andassumptions (Fuentes and Reich, 2010, pg. 62), the spectral density function isgiven by(2.4) f ( ω ) = 1(2 π ) (cid:90) R exp( − i ω T h ) C ( h ) d h , so that the covariance function, C ( h ), and the spectral density function, f ( ω ),form a Fourier transform pair.Properties of the covariance function imply properties of the spectral density.For example, if the covariance function is isotropic (2.1), then the spectral den-sity (2.4) depends on ω only through its length, ω = || ω || , and we can write f ( ω ) = f ( ω ), where f ( · ) is called the isotropic spectral density (Fuentes, 2013).Consequently, second order properties of a second order stationary RF can beexplored via either the covariance function or the spectral density function. Teststatistics quantifying second order properties can be constructed using the pe-riodogram, an estimator of (2.4) and denoted by I ( · ). For a real-valued spatialprocess on a rectangular grid Z ⊂ R , a moment-based periodogram used toestimate (2.4) is(2.5) I ( ω , ω ) = 1(2 π ) n − (cid:88) h = − n +1 n − (cid:88) h = − n +1 (cid:98) C ( h , h ) cos( h ω + h ω ) , where n and n denote the number of rows and columns of the grid and (cid:98) C ( h , h )is a non-parametric estimator of the covariance function. In practice, the peri-odogram (2.5) is used to estimate the spectral density at the Fourier or harmonicfrequencies. The frequency ω = ( ω , ω ) is a Fourier or harmonic frequency if ω j is a multiple of 2 π/n j , j = 1 ,
2. Furthermore, the set of frequencies is limited to { ω j = 2 πk j /n j , k j = 1 , , . . . , n ∗ j } , where n ∗ j is ( n j − / n j is odd and n j / − n j is even.
3. TESTS OF ISOTROPY AND SYMMETRY3.1 Brief History
Matheron (1961) developed one of the earliest hypothesis test of isotropy whenhe used normality of sample variogram estimators to construct a χ test foranisotropy in mineral deposit data. Cabana (1987) tested for geometric anisotropyusing level curves of random fields. Vecchia (1988) and Baczkowski and Mardia(1990) developed tests for isotropy assuming a parametric covariance function.Baczkowski (1990) also proposed a randomization test for isotropy. Despite theseearly works, little work on testing isotropy was published during the 1990s, al-though the PhD dissertation work of Lu (1994) would eventually have an notewor-thy impact on the literature. Then, in the 2000s, a number of nonparametric tests imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Matheron (1961) developed one of the earliest hypothesis test of isotropy whenhe used normality of sample variogram estimators to construct a χ test foranisotropy in mineral deposit data. Cabana (1987) tested for geometric anisotropyusing level curves of random fields. Vecchia (1988) and Baczkowski and Mardia(1990) developed tests for isotropy assuming a parametric covariance function.Baczkowski (1990) also proposed a randomization test for isotropy. Despite theseearly works, little work on testing isotropy was published during the 1990s, al-though the PhD dissertation work of Lu (1994) would eventually have an notewor-thy impact on the literature. Then, in the 2000s, a number of nonparametric tests imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY of second-order properties emerged. Some of the developments used estimates ofthe variogram or covariogram to test symmetry and isotropy properties (Lu andZimmerman, 2001; Guan, 2003; Guan et al., 2004, 2007; Maity and Sherman,2012). These works generally borrowed ideas from two bodies of literature: (a)theory on the distributional and asymptotic properties of semivariogram estima-tors (e.g., Baczkowski and Mardia, 1987; Cressie, 1993, pg. 69-47; Hall and Patil,1994) and (b) subsampling techniques to estimate the variance of statistics de-rived from spatial data (e.g., Possolo, 1991; Politis and Sherman, 2001; Sherman,1996; Lahiri, 2003; Lahiri and Zhu, 2006). Other nonparametric methods usedthe spectral domain to test isotropy and symmetry (Scaccia and Martin, 2002,2005; Lu and Zimmerman, 2005; Fuentes, 2005). These works generally extendedideas used in the time series literature (e.g., Priestley and Rao, 1969; Priestley,1981) to the spatial case. Methods for testing isotropy and symmetry in both thespatial and spectral domains, under the assumption of a parametric covariancefunction, have also been developed recently (Stein et al., 2004; Haskard, 2007;Fuentes, 2007; Matsuda and Yajima, 2009; Scaccia and Martin, 2011). A popular approach to testing second order properties was pioneered in theworks of Lu (1994) and Lu and Zimmerman (2001) who leveraged the asymptoticjoint normality of the sample variogram computed at different spatial lags. Thesubsequent works of Guan et al. (2004, 2007) and Maity and Sherman (2012)built upon these ideas and are the primary focus of this subsection. Lu (1994)details methods for testing axial symmetry. While Lu and Zimmerman (2001),Guan et al. (2004), and Maity and Sherman (2012) focus on testing isotropy,these methods can also be used to test symmetry. Finally, Bowman and Cru-jeiras (2013) detail a more computational approach for testing isotropy. BothLi et al. (2007, 2008b) and Jun and Genton (2012) use an approach analogousto the methods from Lu and Zimmerman (2001), Guan et al. (2004, 2007), andMaity and Sherman (2012) but for spatiotemporal data. Table 1 summarizes testproperties discussed in this section and Section 3.3.Nonparametric tests for anisotropy in the spatial domain are based on a nullhypothesis that is an approximation to isotropy. Under the null hypothesis thatthe RF is isotropic, it follows that the values of γ ( · ) evaluated at any two spatiallags that have the same norm are equal, regardless of the direction of the lags. Tofully specify the most general null hypothesis of isotropy, theoretically, one wouldneed to compare variogram values for an infinite set of lags. In practice a smallnumber of lags are specified. Then it is possible to test a hypothesis consisting ofa set of linear contrasts of the form(3.1) H : A γ ( · ) = as a proxy for the null hypothesis of isotropy, where A is a full row rank matrix(Lu and Zimmerman, 2001). For example, a set of lags, denoted Λ , commonlyused in practice for gridded sampling locations with unit spacing is(3.2) Λ = { h = (1 , , h = (0 , , h = (1 , , h = ( − , } , and the corresponding A matrix under H : A γ ( Λ ) = is(3.3) A = (cid:20) − − (cid:21) . imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 th November, 2015
ZACHARY D. WELLER
One of the first steps in detecting potential anisotropy is the choice of lags, as thetest results will only hold for the particular set of lags considered (Guan et al.,2004). While this choice is subjective, there are several considerations and usefulguidelines for determining the set of lags (see Section 5).For nonparametric tests of symmetry, the null hypothesis of symmetry using(3.1) can be expressed by a countable set of contrasts for a process observedon a grid. Tests of symmetry will be subject to similar practical considers astests of isotropy, and practitioners testing symmetry properties will need tochoose a small set of lags and form a hypothesis that is a surrogate for sym-metry. For example, testing reflection symmetry of a process observed on theinteger grid would require comparing estimates of C ( · ) evaluated at the lag pairs { (1 , , ( − , } , { (2 , , ( − , } , { (1 , , ( − , } , etc.The tests in Lu and Zimmerman (2001), Guan et al. (2004, 2007), and Maityand Sherman (2012) involve estimating either the semivariogram, γ ( · ), or co-variance, C ( · ), function at the set of chosen lags, Λ . Denoting the set of pointestimates of the semivariogram/covariance function at the given lags as (cid:98) G n , thetrue values as G , and normalizing constant a n , a central result for all three meth-ods is that(3.4) a n ( (cid:98) G n − G ) d −−−→ n →∞ M V N ( , Σ ) , under increasing domain asymptotics and mild moment and mixing conditionson the RF. Using the A matrix, an estimate of the variance covariance matrix, (cid:98) Σ , and (cid:98) G n , a quadratic form is constructed, and a p-value can be obtained froman asymptotic χ distribution with degrees of freedom given by the row rank of A . The primary differences between these works are the assumed distribution ofsampling locations, the shape of the sampling domain, and the estimation of G and Σ . These differences are important when choosing a test that is appropriatefor a particular set of data (see Tables 1 and 2 and Figure 4 for more informationabout these differences).Maity and Sherman (2012) develop a test with the fewest restrictions on theshape of the sampling region and distribution of sampling locations. Their testcan be used when the sampling region is any convex subset in R d and the distri-bution of sampling locations in the region follows any general spatial samplingdesign. The test in Guan et al. (2004) also allows convex subsets in R d and is de-veloped for gridded and non-gridded sampling locations but requires non-griddedsampling locations to be uniformly distributed on the domain, i.e., generated bya homogenous Poisson process. The Poisson assumption is dropped in Guan et al.(2007). Lu and Zimmerman (2001) require the domain to be rectangular and theobservations to lie on a grid.Another difference between methods is the form of the nonparametric estima-tor of G . In Lu and Zimmerman (2001), (cid:98) G n is computed using the log of theclassical sample semivariogram estimator (2.2). Guan et al. (2004, 2007) also usethe estimator in (2.2) for gridded sampling locations, but use a kernel estimator of γ ( h ) for non-gridded locations. Maity and Sherman (2012) use a kernel estimatorof the covariance function. When smoothing over spatial lags in R , the kernelis typically given as Nadaraya-Watson (Nadaraya, 1964; Watson, 1964) prod-uct kernel, independently smoothing over horizontal and vertical lags. Commonchoices for the kernel are the Epanechnikov or truncated Gaussian kernels. The imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
One of the first steps in detecting potential anisotropy is the choice of lags, as thetest results will only hold for the particular set of lags considered (Guan et al.,2004). While this choice is subjective, there are several considerations and usefulguidelines for determining the set of lags (see Section 5).For nonparametric tests of symmetry, the null hypothesis of symmetry using(3.1) can be expressed by a countable set of contrasts for a process observedon a grid. Tests of symmetry will be subject to similar practical considers astests of isotropy, and practitioners testing symmetry properties will need tochoose a small set of lags and form a hypothesis that is a surrogate for sym-metry. For example, testing reflection symmetry of a process observed on theinteger grid would require comparing estimates of C ( · ) evaluated at the lag pairs { (1 , , ( − , } , { (2 , , ( − , } , { (1 , , ( − , } , etc.The tests in Lu and Zimmerman (2001), Guan et al. (2004, 2007), and Maityand Sherman (2012) involve estimating either the semivariogram, γ ( · ), or co-variance, C ( · ), function at the set of chosen lags, Λ . Denoting the set of pointestimates of the semivariogram/covariance function at the given lags as (cid:98) G n , thetrue values as G , and normalizing constant a n , a central result for all three meth-ods is that(3.4) a n ( (cid:98) G n − G ) d −−−→ n →∞ M V N ( , Σ ) , under increasing domain asymptotics and mild moment and mixing conditionson the RF. Using the A matrix, an estimate of the variance covariance matrix, (cid:98) Σ , and (cid:98) G n , a quadratic form is constructed, and a p-value can be obtained froman asymptotic χ distribution with degrees of freedom given by the row rank of A . The primary differences between these works are the assumed distribution ofsampling locations, the shape of the sampling domain, and the estimation of G and Σ . These differences are important when choosing a test that is appropriatefor a particular set of data (see Tables 1 and 2 and Figure 4 for more informationabout these differences).Maity and Sherman (2012) develop a test with the fewest restrictions on theshape of the sampling region and distribution of sampling locations. Their testcan be used when the sampling region is any convex subset in R d and the distri-bution of sampling locations in the region follows any general spatial samplingdesign. The test in Guan et al. (2004) also allows convex subsets in R d and is de-veloped for gridded and non-gridded sampling locations but requires non-griddedsampling locations to be uniformly distributed on the domain, i.e., generated bya homogenous Poisson process. The Poisson assumption is dropped in Guan et al.(2007). Lu and Zimmerman (2001) require the domain to be rectangular and theobservations to lie on a grid.Another difference between methods is the form of the nonparametric estima-tor of G . In Lu and Zimmerman (2001), (cid:98) G n is computed using the log of theclassical sample semivariogram estimator (2.2). Guan et al. (2004, 2007) also usethe estimator in (2.2) for gridded sampling locations, but use a kernel estimator of γ ( h ) for non-gridded locations. Maity and Sherman (2012) use a kernel estimatorof the covariance function. When smoothing over spatial lags in R , the kernelis typically given as Nadaraya-Watson (Nadaraya, 1964; Watson, 1964) prod-uct kernel, independently smoothing over horizontal and vertical lags. Commonchoices for the kernel are the Epanechnikov or truncated Gaussian kernels. The imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY T a b l e : P r o p e r t i e s o f n o np a r a m e t r i c t e s t s o f i s o t r o p y . “ D o m a i n ” r e f e r s t o t h e d o m a i nu s e d t o r e p r e s e n tt h e R F ( s p a t i a l o r s p ec t r a l ) , “ T e s t S t a t B a s e d O n ” li s t s t h e n o np a r a m e t r i ce s t i m a t o r u s e d t o c o n s t r u c tt h e t e s t s t a t i s t i c “ D i s t b ’ n ”g i v e s t h e li m i t i n ga s y m p t o t i c d i s t r i bu t i o n o f t h e t e s t s t a t i s t i c , a nd “ G P ” d e n o t e s w h e t h e r t h e t e s t r e q u i r e s d a t a t o b e g e n e r a t e d f r o m a G a u ss i a np r o ce ss . H y p o t h e s i s T e s t P r o p e r t i e s T e s t M e t h o d I s o t r o p y S y mm e t r y D o m a i n T e s t S t a t B a s e d O n A s y m p t o t i c s D i s t b ’ n G P L u a nd Z i mm e r m a n ( ) y e s y e ss p a t i a l s e m i v a r i og r a m i n c d o m a i n χ y e s G u a n e t a l. ( , ) y e s y e ss p a t i a l ( k e r n e l ) a v a r i og r a m i n c d o m a i n χ b n o S c a cc i aa nd M a r t i n ( , ) p a r t i a l y e ss p e c t r a l p e r i o d og r a m i n c d o m a i n Z , t n o L u a nd Z i mm e r m a n ( ) p a r t i a l y e ss p e c t r a l p e r i o d og r a m i n c d o m a i n χ , F n o F u e n t e s ( ) p a r t i a l n o s p e c t r a l s p a t i a l p e r i o d og r a m s h r i n k i n g ( m i x e d ) χ y e s M a i t y a ndSh e r m a n ( ) y e s y e ss p a t i a l k e r n e l c o v a r i og r a m i n c d o m a i n χ n o B o w m a n a nd C r u j e i r a s ( ) y e s n o s p a t i a l v a r i og r a m i n c d o m a i n a pp r o x χ y e s V a n H a l a e t a l. ( ) y e s y e ss p e c t r a l e m p i r i c a lli k e li h oo d s h r i n k i n g ( m i x e d ) χ n o a f o r g r i dd e d s a m p li n g l o c a t i o n s , t h ee s t i m a t o r i n ( . ) i s u s e d w h il e a k e r n e l v a r i og r a m e s t i m a t o r i s u s e d f o r n o n - g r i dd e d s a m p li n g l o c a t i o n s b p - v a l u e s m a y n ee d t o b e a pp r o x i m a t e du s i n g fin i t e s a m p l e a d j u s t m e n t s T a b l e : T e s t i m p l e m e n t a t i o n , p a r t . “ Sub s a m p ” d e fin e s w h e t h e r s p a t i a l s ub s a m p li n g p r o ce du r e s a r e n ee d e d t o p e r f o r m t h e t e s t , “ S & P s i m ” d e n o t e s w h e t h e r o r n o tt h e a u t h o r ( s ) o f t h e m e t h o dp r o v i d e a s i m u l a t i o n o f t e s t s i ze a ndp o w e r ( S ee a l s o T a b l e ) . H y p o t h e s i s T e s t I m p l e m e n t a t i o n T e s t M e t h o d S a m p li n g D o m a i nSh a p e S a m p li n g D e s i g nSub s a m pS & P S i m S o f t w a r e L u a nd Z i mm e r m a n ( ) r e c t a n g u l a r i n R g r i dn o y e s a n o G u a n e t a l. ( , ) c o n v e x s ub s e t s i n R d g r i d / un i f b / n o n - un i f c y e s y e s a R p a c k ag e spTest S c a cc i aa nd M a r t i n ( , ) r e c t a n g u l a r i n R g r i dn o y e s a n o L u a nd Z i mm e r m a n ( ) r e c t a n g u l a r i n R g r i dn o y e s R p a c k ag e spTest F u e n t e s ( ) r e c t a n g u l a r i n R d g r i dn o y e s a n o M a i t y a ndSh e r m a n ( ) c o n v e x s ub s e t s i n R d n o n - un i f c y e s y e s a R p a c k ag e spTest B o w m a n a nd C r u j e i r a s ( ) c o n v e x s ub s e t s i n R d un i f b n o y e s a R p a c k ag e sm V a n H a l a e t a l. ( ) s ub s e t s i n R d n o n - un i f c n o y e s a n o a s i m u l a t e dd a t aa r e G a u ss i a n o n l y b s a m p li n g l o c a t i o n s m u s t b e g e n e r a t e db y h o m og e n e o u s P o i ss o np r o c e ss ,i. e . un i f o r m l y d i s t r i bu t e d o n t h e d o m a i n c s a m p li n g l o c a t i o n s c a nb e g e n e r a t e db y a n y g e n e r a l s a m p li n g d e s i g n imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY T a b l e : P r o p e r t i e s o f n o np a r a m e t r i c t e s t s o f i s o t r o p y . “ D o m a i n ” r e f e r s t o t h e d o m a i nu s e d t o r e p r e s e n tt h e R F ( s p a t i a l o r s p ec t r a l ) , “ T e s t S t a t B a s e d O n ” li s t s t h e n o np a r a m e t r i ce s t i m a t o r u s e d t o c o n s t r u c tt h e t e s t s t a t i s t i c “ D i s t b ’ n ”g i v e s t h e li m i t i n ga s y m p t o t i c d i s t r i bu t i o n o f t h e t e s t s t a t i s t i c , a nd “ G P ” d e n o t e s w h e t h e r t h e t e s t r e q u i r e s d a t a t o b e g e n e r a t e d f r o m a G a u ss i a np r o ce ss . H y p o t h e s i s T e s t P r o p e r t i e s T e s t M e t h o d I s o t r o p y S y mm e t r y D o m a i n T e s t S t a t B a s e d O n A s y m p t o t i c s D i s t b ’ n G P L u a nd Z i mm e r m a n ( ) y e s y e ss p a t i a l s e m i v a r i og r a m i n c d o m a i n χ y e s G u a n e t a l. ( , ) y e s y e ss p a t i a l ( k e r n e l ) a v a r i og r a m i n c d o m a i n χ b n o S c a cc i aa nd M a r t i n ( , ) p a r t i a l y e ss p e c t r a l p e r i o d og r a m i n c d o m a i n Z , t n o L u a nd Z i mm e r m a n ( ) p a r t i a l y e ss p e c t r a l p e r i o d og r a m i n c d o m a i n χ , F n o F u e n t e s ( ) p a r t i a l n o s p e c t r a l s p a t i a l p e r i o d og r a m s h r i n k i n g ( m i x e d ) χ y e s M a i t y a ndSh e r m a n ( ) y e s y e ss p a t i a l k e r n e l c o v a r i og r a m i n c d o m a i n χ n o B o w m a n a nd C r u j e i r a s ( ) y e s n o s p a t i a l v a r i og r a m i n c d o m a i n a pp r o x χ y e s V a n H a l a e t a l. ( ) y e s y e ss p e c t r a l e m p i r i c a lli k e li h oo d s h r i n k i n g ( m i x e d ) χ n o a f o r g r i dd e d s a m p li n g l o c a t i o n s , t h ee s t i m a t o r i n ( . ) i s u s e d w h il e a k e r n e l v a r i og r a m e s t i m a t o r i s u s e d f o r n o n - g r i dd e d s a m p li n g l o c a t i o n s b p - v a l u e s m a y n ee d t o b e a pp r o x i m a t e du s i n g fin i t e s a m p l e a d j u s t m e n t s T a b l e : T e s t i m p l e m e n t a t i o n , p a r t . “ Sub s a m p ” d e fin e s w h e t h e r s p a t i a l s ub s a m p li n g p r o ce du r e s a r e n ee d e d t o p e r f o r m t h e t e s t , “ S & P s i m ” d e n o t e s w h e t h e r o r n o tt h e a u t h o r ( s ) o f t h e m e t h o dp r o v i d e a s i m u l a t i o n o f t e s t s i ze a ndp o w e r ( S ee a l s o T a b l e ) . H y p o t h e s i s T e s t I m p l e m e n t a t i o n T e s t M e t h o d S a m p li n g D o m a i nSh a p e S a m p li n g D e s i g nSub s a m pS & P S i m S o f t w a r e L u a nd Z i mm e r m a n ( ) r e c t a n g u l a r i n R g r i dn o y e s a n o G u a n e t a l. ( , ) c o n v e x s ub s e t s i n R d g r i d / un i f b / n o n - un i f c y e s y e s a R p a c k ag e spTest S c a cc i aa nd M a r t i n ( , ) r e c t a n g u l a r i n R g r i dn o y e s a n o L u a nd Z i mm e r m a n ( ) r e c t a n g u l a r i n R g r i dn o y e s R p a c k ag e spTest F u e n t e s ( ) r e c t a n g u l a r i n R d g r i dn o y e s a n o M a i t y a ndSh e r m a n ( ) c o n v e x s ub s e t s i n R d n o n - un i f c y e s y e s a R p a c k ag e spTest B o w m a n a nd C r u j e i r a s ( ) c o n v e x s ub s e t s i n R d un i f b n o y e s a R p a c k ag e sm V a n H a l a e t a l. ( ) s ub s e t s i n R d n o n - un i f c n o y e s a n o a s i m u l a t e dd a t aa r e G a u ss i a n o n l y b s a m p li n g l o c a t i o n s m u s t b e g e n e r a t e db y h o m og e n e o u s P o i ss o np r o c e ss ,i. e . un i f o r m l y d i s t r i bu t e d o n t h e d o m a i n c s a m p li n g l o c a t i o n s c a nb e g e n e r a t e db y a n y g e n e r a l s a m p li n g d e s i g n imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER kernel estimators also require the choice of a bandwidth parameter, w . Choosingan appropriate bandwidth is one of the challenges of implementing the tests fornon-gridded sampling locations, and the conclusion of the test has the poten-tial to be sensitive to the choice of the bandwidth parameter (see Section 5 forrecommendations on bandwidth selection).Nonparametric tests in the spatial domain also vary in the estimation of Σ , theasymptotic variance-covariance of (cid:98) G n in (3.4). Lu and Zimmerman (2001) use aplug-in estimator, which requires the choice of a parameter, m , that truncates thesum used in estimation. Spatial resampling methods are another approach used toestimate Σ. The method used for spatial resampling and properties of estimatorscomputed from spatial resampling will depend on the underlying spatial samplingdesign (Lahiri, 2003, pg. 281). Guan et al. (2004, 2007) use a moving window ap-proach, creating overlapping subblocks that cover the sampling region. Maity andSherman (2012) employ the grid based block bootstrap (GBBB) (Lahiri and Zhu,2006). The GBBB approach divides the spatial domain into regions, then replaceseach region by sampling (with replacement) a block of the sampling domain hav-ing the same shape and volume as the region, creating a spatial permutation ofblocks of sampling locations. When using the resampling methods, the user mustchose the window or block size and the conclusion of the test has the potential tochange based on the chosen size. Irregularly shaped sampling domains can posea challenge in using the subsampling methods. For example, for an irregularlyshaped sampling domain, many incomplete blocks may complicate the subsam-pling procedure. We summarize guidelines for choosing the window/block size inSection 5.Another approach to testing isotropy in the spatial domain is given by Bowmanand Crujeiras (2013) who take a more empirical and computationally-intensiveapproach. Their methods are available in the R software (R Core Team, 2014)package sm (Bowman and Azzalini, 2014). One caveat of using the sm package isthat the methods are computationally expensive, even for moderate sample sizes.For example, a test of isotropy on 300 uniformly distributed sampling locations ona 10 ×
16 sampling domain took approximately 9.5 minutes where the methodsfrom Guan et al. (2004) took 1.6 seconds using a laptop with 8 GB of memoryand a 2 GHz Intel Core i7 processor. Because of the computational costs, we donot consider this method further.
For gridded sampling locations, nonparametric spectral methods have beendeveloped for testing symmetry (Scaccia and Martin, 2002, 2005; Lu and Zim-merman, 2005) and stationarity (Fuentes, 2005), but none are designed with aprimary goal of testing isotropy. Due to the difficulties presented by non-griddedsampling locations, historically there have been fewer developments using spec-tral methods for non-gridded sampling locations than for gridded (or lattice)data, but this is an area that has received more attention recently (see, e.g.,Fuentes, 2007; Matsuda and Yajima, 2009; Bandyopadhyay et al., 2015). Despitethe challenges, Van Hala et al. (2014) have proposed a nonparametric, empiricallikelihood approach to test isotropy and separability for non-gridded samplinglocations.The primary motivation for using the spectral domain over the spatial domain imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 th November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY are simpler asymptotics in the spectral domain. Unlike estimates of the variogramor covariogram at different spatial lags, estimates of the spectral density at differ-ent frequencies via the periodogram are asymptotically independent under certainconditions (Pagano, 1971; Schabenberger and Gotway, 2004, pg. 78,194). Addi-tionally, in practice, tests of symmetry in the spectral domain are generally notsubject to as many choices (e.g., spatial lag set, bandwidth, block size) as thosein the spatial domain.Analogous to testing isotropy in the spatial domain by using a finite set of spa-tial lags, tests of symmetry in the spectral domain typically involve estimatingand comparing the spectral density (2.4) at a finite set of the Fourier frequencies.For example, axial symmetry (2.2) of the covariance function implies axial sym-metry of the spectral density, f ( ω , ω ) = f ( − ω , ω ), which can be evaluatedby comparing I ( ω , ω ) to I ( − ω , ω ) at a finite set of frequencies. Similarly, thenull hypothesis of isotropy can be approximated by comparing periodogram (2.5)estimates at a set of distinct frequencies with the same norm (Fuentes, 2005).Although most of the current spectral methods are not directly designed to testisotropy, the hypothesis of complete symmetry can be used to reject the assump-tion of isotropy due to (2.3). However, certain types of anisotropy may not bedetected by these tests. For example, a geometrically anisotropic process havingthe major axis of the ellipses of equicorrelation parallel to the x -axis is axiallysymmetric, and the anisotropy wouldn’t be detected by a test of axial symmetry.Scaccia and Martin (2002, 2005) use the periodogram (2.5) to develop a testfor axial symmetry. They propose three test statistics that are a function of theperiodogram values. The first uses the average of the difference in the log of theperiodogram values, log[ I ( ω , ω )] − log[ I ( ω , − ω )]. The second and third teststatistics use the average of standardized periodogram differences, [ I ( ω , ω ) − I ( ω , − ω )] / [ I ( ω , ω ) + I ( ω , − ω )]. These test statistics are shown to asymptot-ically follow a standard normal or t distribution via the Central Limit Theorem,and the corresponding distributions are used to obtain a p-value.Lu and Zimmerman (2005) also use the periodogram as an estimator of thespectral density to test properties of axial and complete symmetry of processeson the integer grid, Z . They use the asymptotic distribution of the periodogramto construct two potential test statistics. Both test statistics leverage the factthat, under certain conditions and at certain frequencies,(3.5) 2 I ( ω , ω ) f ( ω , ω ) iid −−−−−−→ n ,n →∞ χ . Under the null hypothesis of axial or complete symmetry, (3.5) implies that ratiosof periodogram values at different frequencies follow an F (2 ,
2) distribution. Thepreferred test statistic produces a p-value via a Cram´er-von Mises (CvM) good-ness of fit test using the appropriate set of periodogram ratios. Because rejectinga hypothesis of axial symmetry implies rejecting a hypothesis of complete sym-metry, Lu and Zimmerman (2005) recommend a two-stage procedure for testingcomplete symmetry. At the first stage, they test the hypothesis of axial symmetry,and if the null hypothesis is not rejected, they test the hypothesis of completesymmetry. To control the overall type-I error rate at α , the tests at each stagecan be performed using a significance level of α/ imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
2) distribution. Thepreferred test statistic produces a p-value via a Cram´er-von Mises (CvM) good-ness of fit test using the appropriate set of periodogram ratios. Because rejectinga hypothesis of axial symmetry implies rejecting a hypothesis of complete sym-metry, Lu and Zimmerman (2005) recommend a two-stage procedure for testingcomplete symmetry. At the first stage, they test the hypothesis of axial symmetry,and if the null hypothesis is not rejected, they test the hypothesis of completesymmetry. To control the overall type-I error rate at α , the tests at each stagecan be performed using a significance level of α/ imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER likelihood (SFDEL) approach that can be used for inference about spatial co-variance structure. One of the applications of this method is testing isotropy.An advantage of this method over other frequency domain approaches is thatit can be used for non-gridded sampling locations. To implement the test, theuser must select the set of lags and, because the sampling locations are not grid-ded, the number and spacing of frequencies. Van Hala et al. (2014) offer someguidelines for these choices based on the simulations and theoretical considera-tions (e.g., the frequencies need to be asymptotically distant). Once these choicesare made, Van Hala et al. (2014) maximize an empirical likelihood under a mo-ment constraint corresponding to isotropy and show that the log-ratio of the con-strained and unconstrained maximizer asymptotically follows a χ distribution.The SFDEL method relies on the asymptotic independence of the periodogramvalues, and the smallest sample size used in simulations was n = 600. Thus, it isnot clear how the method will perform for smaller sample sizes.Fuentes (2005) introduces a nonparametric, spatially varying spectral densityto represent nonstationary spatial processes. While the method can be used totest the assumption of isotropy, the test requires a large sample size on a finegrid. For this reason and also because the test was primarily designed to test theassumption of stationarity, we do not consider it further.
4. SIMULATION STUDY
We designed a simulation study to compare the empirical size, power, andcomputational costs for four of the methods. For gridded sampling locations, wecompare Lu and Zimmerman (2005)[hereafter, LZ] to Guan et al. (2004)[hereafterdenoted as GSC or GSC-g when we are specifically referring to the test when ap-plied to gridded sampling locations]. For uniformly distributed sampling locationswe compare Maity and Sherman (2012)[MS] to Guan et al. (2004, 2007)[GSC-ufor the method used for uniformly distributed sampling locations].We performed the tests on the same realizations of the RF. Data are sim-ulated on rectangular grids or rectangular sampling domains because they aremore realistic than square domains and simulations on rectangular domains werenot previously demonstrated. We simulate Gaussian data with mean zero andexponential covariance functions with no nugget, a sill of one, and effective rangevalues corresponding to short, medium, and long range dependence. We introducevarying degrees of geometric anisotropy via coordinate transformations governedby a rotation parameter θ and scaling parameter R that define the ellipses ofequicorrelation (see Figure 5 in the Appendix). The parameter θ quantifies theangle between the major axis of the ellipse and the x -axis (counter-clockwise ro-tation) while R gives the ratio of the major and minor axes of the ellipse. Wealso performed simulations that investigate the effect of the lag set, block size,and bandwidth. Although some simulations are given in the original works, oursimulations serve to provide a direct comparison of the effects of changing thesevalues and provide further insight into how to choose them. See the Appendixfor additional simulation details and results.Figures 2 and 3 illustrate a subset of the simulation results compring empiricalsize, power, and computational time (full results in Appendix, Tables 5 and 6).These simulations indicate that nonparametric tests for anisotropy have higherpower for gridded (Figure 5) than for non-gridded (Figure 6) sampling designs. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
We designed a simulation study to compare the empirical size, power, andcomputational costs for four of the methods. For gridded sampling locations, wecompare Lu and Zimmerman (2005)[hereafter, LZ] to Guan et al. (2004)[hereafterdenoted as GSC or GSC-g when we are specifically referring to the test when ap-plied to gridded sampling locations]. For uniformly distributed sampling locationswe compare Maity and Sherman (2012)[MS] to Guan et al. (2004, 2007)[GSC-ufor the method used for uniformly distributed sampling locations].We performed the tests on the same realizations of the RF. Data are sim-ulated on rectangular grids or rectangular sampling domains because they aremore realistic than square domains and simulations on rectangular domains werenot previously demonstrated. We simulate Gaussian data with mean zero andexponential covariance functions with no nugget, a sill of one, and effective rangevalues corresponding to short, medium, and long range dependence. We introducevarying degrees of geometric anisotropy via coordinate transformations governedby a rotation parameter θ and scaling parameter R that define the ellipses ofequicorrelation (see Figure 5 in the Appendix). The parameter θ quantifies theangle between the major axis of the ellipse and the x -axis (counter-clockwise ro-tation) while R gives the ratio of the major and minor axes of the ellipse. Wealso performed simulations that investigate the effect of the lag set, block size,and bandwidth. Although some simulations are given in the original works, oursimulations serve to provide a direct comparison of the effects of changing thesevalues and provide further insight into how to choose them. See the Appendixfor additional simulation details and results.Figures 2 and 3 illustrate a subset of the simulation results compring empiricalsize, power, and computational time (full results in Appendix, Tables 5 and 6).These simulations indicate that nonparametric tests for anisotropy have higherpower for gridded (Figure 5) than for non-gridded (Figure 6) sampling designs. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY T a b l e : T e s t I m p l e m e n t a t i o n , p a r t . T h i s t a b l ec o n t i nu e s t h e li s t o f c h o i ce s a nd c o n s i d e r a t i o n s f o r i m p l e m e n t i n gag i v e n t e s t . “ S a m pS i ze ( S / A ) ” i nd i c a t e s t h e m i n i m u m s a m p l e s i ze s u s e d i n s i m u l a t i o n s ( S ) a nd a pp li c a t i o n s ( A ) p r o v i d e db y t h e a u t h o r ( s ) o f t h e m e t h o d . H y p o t h e s i s T e s t I m p l e m e n t a t i o n T e s t M e t h o d C h o i c e s O t h e r C o n s i d e r a t i o n s S a m pS i z e ( S / A ) L u a nd Z i mm e r m a n ( ) s p a t i a ll ag s e t , t r un c a t i o np a r a m e t e r o p t i m a l t r un c a t i o np a r a m e t e r G u a n e t a l. ( ) g r i dd e dd e s i g n s p a t i a ll ag s e t , w i nd o w s i z e o p t i m a l w i nd o w s i z e , e d g ee ff e c t s , fin i t e s a m p l e a d j u s t m e n t G u a n e t a l. ( ) un i f o r m d e s i g n G u a n e t a l. ( ) n o n - un i f o r m d e s i g n s p a t i a ll ag s e t , k e r n e l f un c t i o n , b a nd w i d t hp a r a m e t e r , w i nd o w s i z e o p t i m a l b a nd w i d t h & w i nd o w s i z e , e d g ee ff e c t s , fin i t e s a m p l e a d j u s t m e n t S c a cc i aa nd M a r t i n ( , )t e s t s t a t i s t i c r e q u i r e s g r i dd e d s a m p li n g l o c a t i o n s ; d e s i g n e d t o t e s t s y mm e t r y L u a nd Z i mm e r m a n ( )t e s t s t a t i s t i c r e q u i r e s g r i dd e d s a m p li n g l o c a t i o n s ; t w o - s t ag e t e s t i n g p r o c e du r e , d e s i g n e d t o t e s t s y mm e t r y ; r e li e s o n a s y m p t o t i c i nd e p e nd e n c e F u e n t e s ( ) k e r n e l f un c t i o n , b a nd w i d t h p a r a m e t e r s , f r e q u e n c y s e t , s p a t i a l k n o t s r e q u i r e s fin e g r i d ; d e s i g n e d t o t e s t s t a t i o n a r i t y M a i t y a ndSh e r m a n ( ) l ag s e t Λ , k e r n e l f un c t i o n , b a nd w i d t h p a r a m e t e r , s ubb l o c k s i z e , nu m b e r o f b oo t s t r a p s a m p l e s o p t i m a l b a nd w i d t h & b l o c k s i z e B o w m a n a nd C r u j e i r a s ( ) b a nd w i d t hp a r a m e t e r c o m pu t a t i o n a ll y i n t e n s i v e V a n H a l a e t a l. ( ) l ag s e t , nu m b e r a nd s p a c i n go f f r e q u e n c i e s o p t i m a l nu m b e r a nd s p a c i n go f f r e q u e n c i e s , r e li e s o n a s y m p t o t i c i nd e p e nd e n c e imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY T a b l e : T e s t I m p l e m e n t a t i o n , p a r t . T h i s t a b l ec o n t i nu e s t h e li s t o f c h o i ce s a nd c o n s i d e r a t i o n s f o r i m p l e m e n t i n gag i v e n t e s t . “ S a m pS i ze ( S / A ) ” i nd i c a t e s t h e m i n i m u m s a m p l e s i ze s u s e d i n s i m u l a t i o n s ( S ) a nd a pp li c a t i o n s ( A ) p r o v i d e db y t h e a u t h o r ( s ) o f t h e m e t h o d . H y p o t h e s i s T e s t I m p l e m e n t a t i o n T e s t M e t h o d C h o i c e s O t h e r C o n s i d e r a t i o n s S a m pS i z e ( S / A ) L u a nd Z i mm e r m a n ( ) s p a t i a ll ag s e t , t r un c a t i o np a r a m e t e r o p t i m a l t r un c a t i o np a r a m e t e r G u a n e t a l. ( ) g r i dd e dd e s i g n s p a t i a ll ag s e t , w i nd o w s i z e o p t i m a l w i nd o w s i z e , e d g ee ff e c t s , fin i t e s a m p l e a d j u s t m e n t G u a n e t a l. ( ) un i f o r m d e s i g n G u a n e t a l. ( ) n o n - un i f o r m d e s i g n s p a t i a ll ag s e t , k e r n e l f un c t i o n , b a nd w i d t hp a r a m e t e r , w i nd o w s i z e o p t i m a l b a nd w i d t h & w i nd o w s i z e , e d g ee ff e c t s , fin i t e s a m p l e a d j u s t m e n t S c a cc i aa nd M a r t i n ( , )t e s t s t a t i s t i c r e q u i r e s g r i dd e d s a m p li n g l o c a t i o n s ; d e s i g n e d t o t e s t s y mm e t r y L u a nd Z i mm e r m a n ( )t e s t s t a t i s t i c r e q u i r e s g r i dd e d s a m p li n g l o c a t i o n s ; t w o - s t ag e t e s t i n g p r o c e du r e , d e s i g n e d t o t e s t s y mm e t r y ; r e li e s o n a s y m p t o t i c i nd e p e nd e n c e F u e n t e s ( ) k e r n e l f un c t i o n , b a nd w i d t h p a r a m e t e r s , f r e q u e n c y s e t , s p a t i a l k n o t s r e q u i r e s fin e g r i d ; d e s i g n e d t o t e s t s t a t i o n a r i t y M a i t y a ndSh e r m a n ( ) l ag s e t Λ , k e r n e l f un c t i o n , b a nd w i d t h p a r a m e t e r , s ubb l o c k s i z e , nu m b e r o f b oo t s t r a p s a m p l e s o p t i m a l b a nd w i d t h & b l o c k s i z e B o w m a n a nd C r u j e i r a s ( ) b a nd w i d t hp a r a m e t e r c o m pu t a t i o n a ll y i n t e n s i v e V a n H a l a e t a l. ( ) l ag s e t , nu m b e r a nd s p a c i n go f f r e q u e n c i e s o p t i m a l nu m b e r a nd s p a c i n go f f r e q u e n c i e s , r e li e s o n a s y m p t o t i c i nd e p e nd e n c e imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER
In both comparisons the methods from GSC have favorable empirical power overthe competitor with a comparable empirical size. As the effective range increases,both empirical size and power tend to increase for the methods from GSC, butthey tend to decrease for MS. GSC-g and LZ have similar computation time,while MS is much more computationally expensive than GSC-u. This differenceis due to the bootstrapping required by MS.Unsurprisingly, as the strength of anisotropy increases (measured by R ), powerincreases for all the methods. For a geometrically anisotropic process, the majorand minor axes of anisotropy are orthogonal. In comparing the effect of the orien-tation of isotropy ( θ ) on the methods, it is important to note that when θ = 0, themajor axis of the ellipse defining the geometric anisotropy is parallel to the x -axisand corresponds to a spatial process that is axially symmetric but not completelysymmetric. When θ = 3 π/ x -axis, and the spatial process is neither axially nor completely sym-metric (see Figure 5 in the Appendix for contours of equal correlation used in thesimulation). The original works generally only simulate data from a geometricallyanisotropic process with the major axis of anisotropy forming a 45-degree anglewith the x -axis; hence, our simulation study more carefully explores the effect ofchanging the orientation of geometric anisotropy. The methods from GSC exhibithigher power when θ = 0 than when θ = 3 π/
8. This is due to the fact that thelag set, Λ , from (3.2) used for the tests contains a pair of spatial lags that areparallel to the major and minor axes of anisotropy when θ = 0, indicating thatan informed choice of spatial lags improves the test’s ability to detect anisotropy.This same result does not hold for MS. It is unclear whether this behavior is ob-served because the method uses the covariogram rather than the semivariogram,the GBBB rather than moving window approach for estimating Σ , or perhapsboth. The simulation results indicate that the LZ test has low empirical power;however, this method was developed to test symmetry properties on square grids,and the choice of a rectangular grid for our simulation study does not allow for alarge number of periodogram ordinates for the second stage of the procedure fortesting the complete symmetry hypothesis.Results from simulations that investigate the effects of changing the lag set, theblock size, and the bandwidth for non-gridded sampling locations are displayedin Tables 7-9, respectively, in the Appendix. For both GSC-u and MS, the lagset in (3.2) provided an empirical size close to the nominal level. Using morelags or longer lags decreased the size and power for GSC-u. This may be dueto the additional uncertainty induced by estimating the covariance between thesemivariance at more lags and the larger variance of semivariance estimates atlonger lags. For MS the longer lags lead to an inflated size and more lags decreasedthe power. In this case, the GBBB may not be capturing the uncertainty incovariance estimates at longer lags with the chosen block size. The MS test wasnot overly sensitive to block size with larger blocks leading to slightly higherpower. MS found that an overly large block size was adverse for test size. ForGSC-u the small and normal sized windows performed at nominal size levels withcomparable power while larger windows were detrimental to test size and power.For GSC-u, we find that choosing a large window tends to lead to overestimationof the asymptotic variance-covariance matrix due to fewer blocks being usedto re-estimate the semivariance. Finally, the results investigating the bandwidth imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 th November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY selection for GSC-u indicate that choosing an overly large bandwidth inflates testsize while choosing too small a bandwidth deflates test size and power. However,the results also indicate that, for the small sample size, test size and power areless negatively affected when approximating the p-value via the finite sampleadjustment.Weller (2015b) demonstrates applications of several of these methods on tworeal data sets. The R package spTest (Weller, 2015a) implements the tests in LZ,GSC, and MS for rectangular grids and sampling regions and is available on theComprehensive R Archive Network ( CRAN ). P r opo r t i on o f R e j e c t i on s . . . . . . . . . GSC-g
GSC-g LZ LZ N = 216
N = 375
N = 216
N = 375 (1.11s) (7.29s) (1.45s) (4.99s)
Method (Time: 1 Test) R = R = R = θ = θ = π Fig 2: Empirical size and power for Guan et al. (2004) [GSC-g] and Lu and Zim-merman (2005) [LZ] for 500 realizations of a mean 0 GRF with gridded samplinglocations using a nominal level of α = 0 .
05. Colors and shapes indicate the typeof anisotropy. Gray points correspond to the isotropic case. The results corre-spond to a “medium” effective range. Computational time for each method isalso displayed.
5. RECOMMENDATIONS
Based on the simulation results we offer recommendations for implementationof nonparametric tests of isotropy. The flow chart in Figure 1 along with Figure4 summarize the steps in the process. Tables 1-3 compare the tests. Table 4summarizes the recommendations provided below.In choosing a nonparametric test for isotropy, the distribution of samplinglocations on the sampling domain is perhaps the most important consideration.Data on a grid simplifies estimation because the semivariogram or covariogram imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Based on the simulation results we offer recommendations for implementationof nonparametric tests of isotropy. The flow chart in Figure 1 along with Figure4 summarize the steps in the process. Tables 1-3 compare the tests. Table 4summarizes the recommendations provided below.In choosing a nonparametric test for isotropy, the distribution of samplinglocations on the sampling domain is perhaps the most important consideration.Data on a grid simplifies estimation because the semivariogram or covariogram imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER T a b l e : G e n e r a l R ec o mm e nd a t i o n s f o r T e s t I m p l e m e n t a t i o n . T h i s t a b l ec o n t a i n s a li s t o f g e n e r a l r ec o mm e nd a t i o n s f o r t e s t i m p l e m e n t a t i o n . T h e s e g u i d e li n e s w ill n o t a pp l y i n a ll s i t u a t i o n s a nd w ill v a r y b a s e d o n a v a r i e t y o ff a c t o r s i n c l ud i n g , bu t n o t li m i t e d t o , t h e s a m p l e s i ze , d e n s i t y o f s a m p li n g l o c a t i o n s , a nd s c a l e o f t h e p r o b l e m . S ee a dd i t i o n a l d i s c u ss i o n i nS ec t i o n . H y p o t h e s i s T e s t C h o i c e s T e s t M e t h o d L ag S e t a B l o c k S i z e B a nd w i d t h P - v a l u e m i n . n G u a n e t a l. ( ) g r i dd e dd e s i g n L e n g t h : s h o r t e r p r e f e rr e d n b < n / n /a fin i t e s a m p l e a d j u s t m e n t G u a n e t a l. ( , ) un i f o r m d e s i g n O r i e n t a t i o n : E q n ( . ) n b (cid:46) n / . < w < . b fin i t e s a m p l e a d j u s t m e n t w h e n n < , a s y m p t o t i c χ w h e n n ≥ M a i t y a ndSh e r m a n ( ) N u m b e r : ( p a i r s ) n b (cid:38) n / e m p i r i c a l b a nd w i d t h a s y m p t o t i c χ a P r i o r k n o w l e d g e ,i f a v a il a b l e , s h o u l db e u s e d t o i n f o r m t h e c h o i c e o f l ag s . b O u r s i m u l a t i o n ss u gg e s tt h e s e b a nd w i d t h v a l u e s a r e r e a s o n a b l e w h e nu s i n ga G a u ss i a n k e r n e l w i t h t r un c a t i o np a r a m e t e r o f . . imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Based on the simulation results we offer recommendations for implementationof nonparametric tests of isotropy. The flow chart in Figure 1 along with Figure4 summarize the steps in the process. Tables 1-3 compare the tests. Table 4summarizes the recommendations provided below.In choosing a nonparametric test for isotropy, the distribution of samplinglocations on the sampling domain is perhaps the most important consideration.Data on a grid simplifies estimation because the semivariogram or covariogram imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER T a b l e : G e n e r a l R ec o mm e nd a t i o n s f o r T e s t I m p l e m e n t a t i o n . T h i s t a b l ec o n t a i n s a li s t o f g e n e r a l r ec o mm e nd a t i o n s f o r t e s t i m p l e m e n t a t i o n . T h e s e g u i d e li n e s w ill n o t a pp l y i n a ll s i t u a t i o n s a nd w ill v a r y b a s e d o n a v a r i e t y o ff a c t o r s i n c l ud i n g , bu t n o t li m i t e d t o , t h e s a m p l e s i ze , d e n s i t y o f s a m p li n g l o c a t i o n s , a nd s c a l e o f t h e p r o b l e m . S ee a dd i t i o n a l d i s c u ss i o n i nS ec t i o n . H y p o t h e s i s T e s t C h o i c e s T e s t M e t h o d L ag S e t a B l o c k S i z e B a nd w i d t h P - v a l u e m i n . n G u a n e t a l. ( ) g r i dd e dd e s i g n L e n g t h : s h o r t e r p r e f e rr e d n b < n / n /a fin i t e s a m p l e a d j u s t m e n t G u a n e t a l. ( , ) un i f o r m d e s i g n O r i e n t a t i o n : E q n ( . ) n b (cid:46) n / . < w < . b fin i t e s a m p l e a d j u s t m e n t w h e n n < , a s y m p t o t i c χ w h e n n ≥ M a i t y a ndSh e r m a n ( ) N u m b e r : ( p a i r s ) n b (cid:38) n / e m p i r i c a l b a nd w i d t h a s y m p t o t i c χ a P r i o r k n o w l e d g e ,i f a v a il a b l e , s h o u l db e u s e d t o i n f o r m t h e c h o i c e o f l ag s . b O u r s i m u l a t i o n ss u gg e s tt h e s e b a nd w i d t h v a l u e s a r e r e a s o n a b l e w h e nu s i n ga G a u ss i a n k e r n e l w i t h t r un c a t i o np a r a m e t e r o f . . imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY P r opo r t i on o f R e j e c t i on s . . . . . . . . . GSC-u
GSC-u MS MS N = 300
N = 450
N = 300
N = 450 (2.17s) (4.44s) (83.40s) (271.22s)
Method (Time: 1 Test) R = R = R = θ = θ = π Fig 3: Empirical size and power for Guan et al. (2004) [denoted GU] and Maityand Sherman (2012) [denoted MS] for 200 realizations of a mean 0 GRF withuniformly distributed sampling locations using a nominal level of α = 0 .
05. Col-ors and shapes indicate the type of anisotropy. Gray points correspond to theisotropic case. The results correspond to a “medium” effective range. Computa-tional time for each method is also displayed.can be estimated at spatial lags that are exactly observed separating pairs ofsampling locations. A grid also allows the option of using easily implementedtests in the spectral domain.Sample size requirements for the asymptotic properties of tests using the spatialdomain to approximately hold will depend on the dependence structure of therandom field. GSC note that convergence of their test statistic is slow in the caseof gridded sampling locations and obtain an approximate p-value via subsamplingrather than the asymptotic χ distribution. Tests using the spectral domain relyon the asymptotic independence of periodogram values, and correlation in finitesamples can lead to an inflated test size (LZ). Based on our simulations, werecommend the sample size be at least 150 for gridded sampling locations and atleast 300 for non-gridded sampling locations. However, power tends be low whenthe sample size is small and/or the anisotropy is weak (Figures 2 and 3).We focus on implementation of the methods that use the spatial domain for theremainder of this section. We discuss the choice of lags, block size, and bandwidthfor the tests in GSC and MS. Due to the large number of choices required toimplement the tests (e.g., block size, bandwidth, kernel function, subsamplingmethod), features of the random field (e.g., sill, range), and properties of the imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
05. Col-ors and shapes indicate the type of anisotropy. Gray points correspond to theisotropic case. The results correspond to a “medium” effective range. Computa-tional time for each method is also displayed.can be estimated at spatial lags that are exactly observed separating pairs ofsampling locations. A grid also allows the option of using easily implementedtests in the spectral domain.Sample size requirements for the asymptotic properties of tests using the spatialdomain to approximately hold will depend on the dependence structure of therandom field. GSC note that convergence of their test statistic is slow in the caseof gridded sampling locations and obtain an approximate p-value via subsamplingrather than the asymptotic χ distribution. Tests using the spectral domain relyon the asymptotic independence of periodogram values, and correlation in finitesamples can lead to an inflated test size (LZ). Based on our simulations, werecommend the sample size be at least 150 for gridded sampling locations and atleast 300 for non-gridded sampling locations. However, power tends be low whenthe sample size is small and/or the anisotropy is weak (Figures 2 and 3).We focus on implementation of the methods that use the spatial domain for theremainder of this section. We discuss the choice of lags, block size, and bandwidthfor the tests in GSC and MS. Due to the large number of choices required toimplement the tests (e.g., block size, bandwidth, kernel function, subsamplingmethod), features of the random field (e.g., sill, range), and properties of the imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER
Nonparametric Hypothesis TestSpatial Sampling DesignRF representation GU MSBowman MS G r i d d e d D e s i g n U n i f o r m D e s i g n G e n e r a l D e s i g n LZ SMFuentes GG Bowman S p e c t r a l D o m a i n S p a t i a l D o m a i n Model forCovarianceFunction(See Figure 1)
Fig 4: Spatial sampling design considerations for choosing a nonparametric hy-pothesis test of isotropy. The method we recommended for testing isotropy ineach situation is given in bold including LZ = Lu and Zimmerman (2005); SM =Scaccia and Martin (2005); GSC-g = Guan et al. (2004) for gridded sampling lo-cations; GSC-u = Guan et al. (2004) for uniformly distributed sampling locations;MS = Maity and Sherman (2012); GSC-n = Guan et al. (2007) for non-uniformsampling locations.sampling design (e.g., density of sampling locations, shape of sampling domain),the recommendations we offer will not apply in all situations. The numerousmoving parts make it challenging to develop general recommendations, especiallywhen choosing a bandwidth.When determining the lag set, Λ , for use in (3.1), the user needs to select(a) the norm of the lags (e.g., Euclidean distance),(b) the orientation (direction) of the lags, and(c) the number of lags.Regarding (a), short lags are preferred. Estimates of the spatial dependence atlarge lags may be less reliable than estimates at shorter lags because they arebased on a smaller number of pairs of observations and hence more variable.Additionally, empirical and theoretical evidence (Lu and Zimmerman, 2001) in-dicates that values of γ ( · ) in two different directions generally exhibit the largestdifference at a lag less than the effective range, the distance beyond which pairsof observations can be assumed to be independent. Finally, there is mathemat-ical support that correctly specifying the covariance function at short lags is imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Fig 4: Spatial sampling design considerations for choosing a nonparametric hy-pothesis test of isotropy. The method we recommended for testing isotropy ineach situation is given in bold including LZ = Lu and Zimmerman (2005); SM =Scaccia and Martin (2005); GSC-g = Guan et al. (2004) for gridded sampling lo-cations; GSC-u = Guan et al. (2004) for uniformly distributed sampling locations;MS = Maity and Sherman (2012); GSC-n = Guan et al. (2007) for non-uniformsampling locations.sampling design (e.g., density of sampling locations, shape of sampling domain),the recommendations we offer will not apply in all situations. The numerousmoving parts make it challenging to develop general recommendations, especiallywhen choosing a bandwidth.When determining the lag set, Λ , for use in (3.1), the user needs to select(a) the norm of the lags (e.g., Euclidean distance),(b) the orientation (direction) of the lags, and(c) the number of lags.Regarding (a), short lags are preferred. Estimates of the spatial dependence atlarge lags may be less reliable than estimates at shorter lags because they arebased on a smaller number of pairs of observations and hence more variable.Additionally, empirical and theoretical evidence (Lu and Zimmerman, 2001) in-dicates that values of γ ( · ) in two different directions generally exhibit the largestdifference at a lag less than the effective range, the distance beyond which pairsof observations can be assumed to be independent. Finally, there is mathemat-ical support that correctly specifying the covariance function at short lags is imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY important for spatial prediction (Stein, 1988). Considering (b), if the process isanisotropic, the ideal choice of Λ and A contrasts lags with the same norm butoriented in the direction of weakest and strongest spatial correlation. Typically,the directions of weakest and strongest spatial correlation will be orthogonal andthus, lags contrasted using the A matrix should also be orthogonal. Prior infor-mation, if available, about the underlying physical/biological process giving riseto the data can also be used to inform the orientation of the lags (Guan et al.,2004). If no prior information about potential anisotropy is available, lags ori-ented in the same directions as those in (3.2) are a good starting set. In regardsto (c), detecting certain types of anisotropy requires a sufficient number of lagsbut using a large number of lags requires a large number of observations (Guanet al., 2004). As a general guideline, we suggest using four lags to construct twocontrasts.Several tests require selection of a window or block size to estimate the variance-covariance matrix. The moving window from GSC creates overlapping subblocksof data by sliding the window over a grid placed on the region. Each of these sub-blocks are used to re-estimate the semivariance. The block size from MS definesthe size of resampled blocks when implementing the GBBB. The GBBB permutesresampled blocks to create a new realization of the process over the entire domain.Choosing the window size in GSC requires balancing two competing goals. First,the moving window should be large enough to create subblocks that are represen-tative of the dependence structure for the entire RF. Second, the window shouldbe small enough to allow for a sufficient number of subblocks to re-estimate thesemivariance, as these values are used to obtain an estimate of the asymptoticvariance-covariance. A window that is too large or too small can potentially leadto under- or over-estimation of the asymptotic variance-covariance. For GSC-u,the windows must be large enough to ensure enough pairs of sampling locationsare in each subblock to compute an estimate of the semivariance without havingto over-smooth. For gridded sampling locations, GSC demonstrate good empiricalsize and power by using moving windows with size of only 2 ×
2. However, theyfind slow convergence to the asymptotic χ distribution, and a p-value is insteadcomputed by approximating the distribution of the test statistic by computing itsvalue for each of the subblocks. Hence, approximating the p-value to two decimalplaces will require at least 100 subblocks over the sampling region. This may notbe possible in practice. For example, a 12 ×
12 grid of sampling locations withmoving windows of size 2 × n b , Sherman (1996) proposes choosing the block sizesuch that n b ≈ cn / for a constant, c , when the spatial dependence does not ex-hibit a large range. In a number of different applications of spatial subsampling, c is typically chosen to be between 0.5 and 2 (Politis and Sherman, 2001; Guan imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
12 grid of sampling locations withmoving windows of size 2 × n b , Sherman (1996) proposes choosing the block sizesuch that n b ≈ cn / for a constant, c , when the spatial dependence does not ex-hibit a large range. In a number of different applications of spatial subsampling, c is typically chosen to be between 0.5 and 2 (Politis and Sherman, 2001; Guan imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLER et al., 2004, 2006). Based on our simulations, we find acceptable empirical sizeand power for GSC-g using small windows and approximating the p-value withthe finite sample adjustment. Thus, we recommend setting n b < n / for GSC-g.For example, we used windows with size 3 × × ×
12 and 25 ×
15, respectively. In the case of uniformly distributed samplinglocations (see Table 8 in the Appendix), the empirical size and power from GSC-uwas negatively affected by a large moving window size; hence, we recommend set-ting c = 1 and choosing n b (cid:46) n / . For the MS test, a small block size negativelyaffected the empirical size and power; thus, we recommend choosing n b (cid:38) n / for this test.Between the choices of a lag set, block size, and bandwidth, choosing an appro-priate bandwidth to smooth over observed spatial lags for non-gridded samplinglocations is the most challenging. For GSC-u the user needs to choose the formof the smoothing kernel as well as the bandwidth for both the entire grid andthe subblocks while MS use an Epanechnikov kernel and empirical bandwidthbased on a user-specified tuning parameter. If the selected bandwidth is too largethen over-smoothing occurs. In oversmoothing, there is very little filtering of thelag distance and direction. The lack of filtering produces similar estimates of thespatial dependence at lags with different directions and distances. If the selectedbandwidth is too small, then there is very little smoothing and estimates of thespatial dependence are based on a small number of pairs of sampling locationsand thus highly variable. Considering the aforementioned effects of the band-width, the bandwidth should decrease as n increases under the usual increasingdomain asymptotics. For example, simulations (not included) indicated a band-width of w = 0 .
65 maintains nominal size when n = 950, but leads to deflated testsize and power when n = 400 on a smaller domain. Garc´ıa-Soid´an et al. (2004),Garc´ıa-Soid´an (2007), and Kim and Park (2012) develop theoretically optimalbandwidths for nonparametric semivariogram estimation, but these works arenot applicable here because they focus on the isotropic case and require an esti-mate of the second derivative of the variogram. We have found that the empiricalbandwidth used by MS tends to produce nominal size (see Table 6). For GSC-u wefind the most consistent results with a bandwidth in the range of 0 . < w < . n < χ dis-tribution. Thus, the user should consider using the finite sample adjustment fornon-gridded sampling locations when the sample size is small there are at least100 subblocks. While it is challenging to choose a bandwidth for GSC-u and thep-value of the test is sensitive to this parameter, the method exhibits nominalsize and substantially higher power than MS when chosen appropriately.
6. CONCLUSIONS
There are several important avenues of future research. Methods to more for-mally characterize the optimal block size and bandwidth parameters for the testsin the spatial domain would enhance the applicability of the tests. The perfor-mance of the tests for non-gridded data in Guan et al. (2004) and Maity and imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
There are several important avenues of future research. Methods to more for-mally characterize the optimal block size and bandwidth parameters for the testsin the spatial domain would enhance the applicability of the tests. The perfor-mance of the tests for non-gridded data in Guan et al. (2004) and Maity and imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY Sherman (2012) are sensitive to these choices and their optimality remains anopen question. Zhang et al. (2014) develop a nonparametric method for estimat-ing the asymptotic variance-covariance matrix of statistics derived from spatialdata that avoids choosing tuning parameters which could simplify test implemen-tation. A second area of future work is further development of nonparametric testsof isotropy for gridded and non-gridded data in the spectral domain. A third areaof further investigation is to compare nonparametric to parametric methods fortesting isotropy, e.g., Scaccia and Martin (2011). A final area of future work isdevelopment of a formal definition and more careful quantification of power ofthe tests. For example, the degree of geometric anisotropy could be quantifiedusing different characteristics of the covariance function, including the ratio ofthe major and minor axes of the ellipse, degree of rotation of the ellipse fromthe coordinate axes, and range of the process. Furthermore, it is important toconsider the effects of density and design of sampling locations, sample size, andthe amount of noise (nugget and sill) in the observations on a test’s ability todetect anisotropy.There is a volume of work on tests for isotropy in other areas of spatial statis-tics. Methods for detecting anisotropy in spatial point process data have beendeveloped, e.g., Schabenberger and Gotway (2004, pg. 200-205), Guan (2003),Guan et al. (2006), and Nicolis et al. (2010). For multivariate spatial data, Jona-Lasinio (2001) proposes a test for isotropy. Gneiting et al. (2007) provide a reviewof potential second order assumptions and models for spatiotemporal geostatis-tical data, and a number of tests for second order properties of spatiotemporaldata have been developed, e.g., Fuentes (2006), Li et al. (2007), Park and Fuentes(2008), Shao and Li (2009), Jun and Genton (2012). Li et al. (2008a) constructa test of the covariance structure for multivariate spatiotemporal data. Testsfor isotropy have also been developed in the computer science literature (e.g.,Molina and Feito, 2002; Chorti and Hristopulos, 2008; Spiliopoulos et al., 2011;Thon et al., 2015).Appropriately specifying the second order properties of the random field isan important step in modeling spatial data, and a number of models have beendeveloped to capture anisotropy in spatial processes. Graphical tools, such asdirectional sample semivariograms, are commonly used to evaluate the assump-tion of isotropy, but these diagnostics can be misleading and open to subjectiveinterpretation. We have presented and reviewed a number of procedures thatcan be used to more objectively test hypotheses of isotropy and symmetry with-out assuming a parametric form for the covariance function. These tests may behelpful for a novice user deciding on an appropriate spatial model. In abandoningparametric assumptions, these hypothesis testing procedures are subject and sen-sitive to choices regarding smoothing parameters, subsampling procedures, andfinite sample adjustments. The test that is most appropriate for a set of data willlargely depend on the sampling design. Additionally, there are trade-offs betweenthe empirical power demonstrated by the tests and the number of choices usermust make to implement the tests (e.g., between Guan et al. (2004) and Maityand Sherman (2012)). We have offered recommendations regarding the variouschoices of method and their implementation and have made the tests available inthe spTest software. Because of the sensitivity of the tests to the various choices,we believe that graphical techniques and nonparametric hypothesis tests should imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 th November, 2015 ZACHARY D. WELLER be used in a complementary role. Graphical techniques can provide an initial in-dication of isotropy properties and inform sensible choices for a hypothesis test,e.g., in choosing the spatial lag set, while hypothesis tests can affirm intuitionabout graphical techniques.
REFERENCES
Baczkowski, A. (1990). A test of spatial isotropy. In
Compstat , pages 277–282. Springer.Baczkowski, A. and Mardia, K. (1987). Approximate lognormality of the sample semi-variogramunder a gaussian process.
Communications in Statistics-Simulation and Computation ,16(2):571–585.Baczkowski, A. and Mardia, K. (1990). A test of spatial symmetry with general application.
Communications in Statistics-Theory and Methods , 19(2):555–572.Bandyopadhyay, S., Lahiri, S. N., Nordman, D. J., et al. (2015). A frequency domain empiricallikelihood method for irregularly spaced spatial data.
The Annals of Statistics , 43(2):519–545.Banerjee, S., Carlin, B. P., and Gelfand, A. E. (2014).
Hierarchical modeling and analysis forspatial data . CRC Press.Borgman, L. and Chao, L. (1994). Estimation of a multidimensional covariance function in caseof anisotropy.
Mathematical geology , 26(2):161–179.Bowman, A. W. and Azzalini, A. (2014).
R package sm : nonparametric smoothing methods(version 2.2-5.4) . University of Glasgow, UK and Universit`a di Padova, Italia.Bowman, A. W. and Crujeiras, R. M. (2013). Inference for variograms. Computational Statistics& Data Analysis , 66:19–31.Cabana, E. (1987). Affine processes: a test of isotropy based on level sets.
SIAM Journal onApplied Mathematics , 47(4):886–891.Chorti, A. and Hristopulos, D. T. (2008). Nonparametric identification of anisotropic (elliptic)correlations in spatially distributed data sets.
Signal Processing, IEEE Transactions on ,56(10):4738–4751.Cressie, N. (1993).
Statistics for Spatial Data: Wiley Series in Probability and Statistics . Wiley-Interscience New York.Ecker, M. D. and Gelfand, A. E. (1999). Bayesian modeling and inference for geometricallyanisotropic spatial data.
Mathematical Geology , 31(1):67–83.Ecker, M. D. and Gelfand, A. E. (2003). Spatial modeling and prediction under stationarynon-geometric range anisotropy.
Environmental and Ecological Statistics , 10(2):165–178.Fuentes, M. (2005). A formal test for nonstationarity of spatial stochastic processes.
Journal ofMultivariate Analysis , 96(1):30–54.Fuentes, M. (2006). Testing for separability of spatial–temporal covariance functions.
Journalof statistical planning and inference , 136(2):447–466.Fuentes, M. (2007). Approximate likelihood for large irregularly spaced spatial data.
Journalof the American Statistical Association , 102(477):321–331.Fuentes, M. (2013). Spectral methods.
Wiley StatsRef: Statistics Reference Online .Fuentes, M. and Reich, B. (2010). Spectral domain.
Handbook of Spatial Statistics , pages 57–77.Garc´ıa-Soid´an, P. (2007). Asymptotic normality of the Nadaraya–Watson semivariogram esti-mators.
Test , 16(3):479–503.Garc´ıa-Soid´an, P. H., Febrero-Bande, M., and Gonz´alez-Manteiga, W. (2004). Nonparametrickernel estimation of an isotropic variogram.
Journal of Statistical Planning and Inference ,121(1):65–92.Gneiting, T., Genton, M., and Guttorp, P. (2007). Geostatistical space-time models, stationarity,separability and full symmetry.
Statistical Methods for Spatio-Temporal Systems , pages 151–175.Guan, Y., Sherman, M., and Calvin, J. A. (2004). A nonparametric test for spatial isotropyusing subsampling.
Journal of the American Statistical Association , 99(467):810–821.Guan, Y., Sherman, M., and Calvin, J. A. (2006). Assessing isotropy for spatial point processes.
Biometrics , 62(1):119–125.Guan, Y., Sherman, M., and Calvin, J. A. (2007). On asymptotic properties of the markvariogram estimator of a marked point process.
Journal of statistical planning and inference ,137(1):148–161.Guan, Y. T. (2003).
Nonparametric methods of assessing spatial isotropy . PhD thesis, TexasA&M University. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Nonparametric methods of assessing spatial isotropy . PhD thesis, TexasA&M University. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY Hall, P. and Patil, P. (1994). Properties of nonparametric estimators of autocovariance forstationary random fields.
Probability Theory and Related Fields , 99(3):399–424.Haskard, K. A. (2007).
An anisotropic Mat´ern spatial covariance model: REML estimation andproperties.
PhD thesis, University of Adelaide.Irvine, K. M., Gitelman, A. I., and Hoeting, J. A. (2007). Spatial designs and properties ofspatial correlation: effects on covariance estimation.
Journal of agricultural, biological, andenvironmental statistics , 12(4):450–469.Isaaks, E. H. and Srivastava, R. M. (1989).
Applied geostatistics , volume 2. Oxford UniversityPress New York.Jona-Lasinio, G. (2001). Modeling and exploring multivariate spatial variation: A test procedurefor isotropy of multivariate spatial data.
Journal of Multivariate Analysis , 77(2):295–317.Journel, A. G. and Huijbregts, C. J. (1978).
Mining geostatistics . Academic press.Jun, M. and Genton, M. G. (2012). A test for stationarity of spatio-temporal random fields onplanar and spherical domains.
Statistica Sinica , 22(4):1737.Kim, T. Y. and Park, J.-S. (2012). On nonparametric variogram estimation.
Journal of theKorean Statistical Society , 41(3):399–413.Lahiri, S. and Zhu, J. (2006). Resampling methods for spatial regression models under a classof stochastic designs.
The Annals of Statistics , 34(4):1774–1813.Lahiri, S. N. (2003).
Resampling methods for dependent data . Springer Science & BusinessMedia.Li, B., Genton, M. G., and Sherman, M. (2007). A nonparametric assessment of proper-ties of space–time covariance functions.
Journal of the American Statistical Association ,102(478):736–744.Li, B., Genton, M. G., and Sherman, M. (2008a). Testing the covariance structure of multivariaterandom fields.
Biometrika , 95(4):813–829.Li, B., Genton, M. G., Sherman, M., et al. (2008b). On the asymptotic joint distribution ofsample space–time covariance estimators.
Bernoulli , 14(1):228–248.Lu, H. and Zimmerman, D. (2001). Testing for isotropy and other directional symmetry prop-erties of spatial correlation. preprint .Lu, H.-C. (1994).
On the distributions of the sample covariogram and semivariogram and theiruse in testing for isotropy . PhD thesis, University of Iowa.Lu, N. and Zimmerman, D. L. (2005). Testing for directional symmetry in spatial dependenceusing the periodogram.
Journal of Statistical Planning and Inference , 129(1):369–385.Maity, A. and Sherman, M. (2012). Testing for spatial isotropy under general designs.
Journalof statistical planning and inference , 142(5):1081–1091.Matheron, G. (1961). Precision of exploring a stratified formation by boreholes with rigidspacing-application to a bauxite deposit. In
International Symposium of Mining Research,University of Missouri , volume 1, pages 407–22.Matheron, G. (1962).
Trait´e de g´eostatistique appliqu´ee . Editions Technip.Matsuda, Y. and Yajima, Y. (2009). Fourier analysis of irregularly spaced data on rd.
Journalof the Royal Statistical Society: Series B (Statistical Methodology) , 71(1):191–217.Modjeska, J. S. and Rawlings, J. (1983). Spatial correlation analysis of uniformity data.
Bio-metrics , pages 373–384.Molina, A. and Feito, F. R. (2002). A method for testing anisotropy and quantifying its directionin digital images.
Computers & Graphics , 26(5):771–784.Nadaraya, E. A. (1964). On estimating regression.
Theory of Probability & Its Applications ,9(1):141–142.Nicolis, O., Mateu, J., and DErcole, R. (2010). Testing for anisotropy in spatial point pro-cesses. In
Proceedings of the Fifth International Workshop on Spatio-Temporal Modelling ,pages 1990–2010.Pagano, M. (1971). Some asymptotic properties of a two-dimensional periodogram.
Journal ofApplied Probability , 8(4):841–847.Park, M. S. and Fuentes, M. (2008). Testing lack of symmetry in spatial–temporal processes.
Journal of statistical planning and inference , 138(10):2847–2866.Politis, D. N. and Sherman, M. (2001). Moment estimation for statistics from marked point pro-cesses.
Journal of the Royal Statistical Society: Series B (Statistical Methodology) , 63(2):261–275.Possolo, A. (1991). Subsampling a random field.
Lecture Notes-Monograph Series , pages 286–294. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Lecture Notes-Monograph Series , pages 286–294. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015 ZACHARY D. WELLERPriestley, M. and Rao, T. S. (1969). A test for non-stationarity of time-series.
Journal of theRoyal Statistical Society. Series B (Methodological) , pages 140–149.Priestley, M. B. (1981).
Spectral analysis and time series . Academic press.R Core Team (2014).
R: A Language and Environment for Statistical Computing . R Foundationfor Statistical Computing, Vienna, Austria.Scaccia, L. and Martin, R. (2002). Testing for simplification in spatial models. In
Compstat ,pages 581–586. Springer.Scaccia, L. and Martin, R. (2005). Testing axial symmetry and separability of lattice processes.
Journal of Statistical Planning and Inference , 131(1):19–39.Scaccia, L. and Martin, R. (2011). Model-based tests for simplification of lattice processes.
Journal of Statistical Computation and Simulation , 81(1):89–107.Schabenberger, O. and Gotway, C. A. (2004).
Statistical methods for spatial data analysis . CRCPress.Shao, X. and Li, B. (2009). A tuning parameter free test for properties of space–time covariancefunctions.
Journal of Statistical Planning and Inference , 139(12):4031–4038.Sherman, M. (1996). Variance estimation for statistics computed from spatial lattice data.
Journal of the Royal Statistical Society. Series B (Methodological) , pages 509–523.Sherman, M. (2011).
Spatial statistics and spatio-temporal data: covariance functions and direc-tional properties . John Wiley & Sons.Spiliopoulos, I., Hristopulos, D. T., Petrakis, M., and Chorti, A. (2011). A multigrid method forthe estimation of geometric anisotropy in environmental data from sensor networks.
Com-puters & Geosciences , 37(3):320–330.Stein, M. L. (1988). Asymptotically efficient prediction of a random field with a misspecifiedcovariance function.
The Annals of Statistics , 16(1):55–63.Stein, M. L., Chi, Z., and Welty, L. J. (2004). Approximating likelihoods for large spatial datasets.
Journal of the Royal Statistical Society: Series B (Statistical Methodology) , 66(2):275–296.Thon, K., Geilhufe, M., and Percival, D. B. (2015). A multiscale wavelet-based test for isotropyof random fields on a regular lattice.
Image Processing, IEEE Transactions on , 24(2):694–708.Van Hala, M., Bandyopadhyay, S., Lahiri, S. N., and Nordman, D. J. (2014). A frequency domainempirical likelihood for estimation and testing of spatial covariance structure. preprint .Vecchia, A. V. (1988). Estimation and model identification for continuous spatial processes.
Journal of the Royal Statistical Society. Series B (Methodological) , pages 297–312.Watson, G. S. (1964). Smooth regression analysis.
Sankhy¯a: The Indian Journal of Statistics,Series A , pages 359–372.Weller, Z. (2015a). spTest: Nonparametric Hypothesis Tests of Isotropy and Symmetry . Rpackage version 0.2.2.Weller, Z. D. (2015b). spTest: an R package implementing nonparametric tests of isotropy. submitted to Journal of Statistical Software, available on arXiv .Zhang, X., Li, B., and Shao, X. (2014). Self-normalization for spatial data.
Scandinavian Journalof Statistics , 41(2):311–324.Zimmerman, D. L. (1993). Another look at anisotropy in geostatistics.
Mathematical Geology ,25(4):453–470. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Mathematical Geology ,25(4):453–470. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY APPENDIXSIMULATION STUDY DETAILS AND FURTHER RESULTS
We define the isotropic exponential covariance function as(6.1) C ( h ) = (cid:26) σ exp( − φh ) if h > ,τ + σ otherwisewhere h = || s i − s j || is the distance between sites s i and s j (Irvine et al., 2007).The corresponding semivariogram is γ ( h ) = ( τ + σ ) − σ exp( − φh ), where τ is the nugget, τ + σ is the sill, and the effective range, ξ , the distance beyondwhich the correlation between observations is less than 0.05, is ξ = − φ log (cid:18) . τ + σ σ (cid:19) . Simulations in Section 4 were performed using the exponential covariance function(6.1) with a partial sill, σ , of 1 and no nugget, τ = 0. We also performedsimulations using different nugget values (results not included). As expected,introducing a nugget had an adverse effect on empirical test size and power. Forthe no nugget simulations, effective ranges, ξ , for isotropic processes were chosento be 3, 6, and 12 corresponding to short, medium, and long range dependence.Geometric anisotropy was introduced by transforming the sampling locationsaccording to a scaling parameter, R , and a rotation parameter, θ . Given an ( R, θ )pair, the coordinates ( x, y ) are transformed to the “anisotropic” coordinates,( x a , y a ) via ( x a , y a ) = ( x, y ) (cid:20) cos θ sin θ − sin θ cos θ (cid:21) (cid:20) R (cid:21) . A realization from the anisotropic process is then created by simulating using thedistance matrix from the transformed coordinates and placing the observed valuesat their corresponding untransformed sampling locations. Figure 5 shows theisotropic exponential correlogram corresponding to τ = 1 and ξ = 6 and contoursof equicorrelation corresponding to the ( R, θ ) values used in the simulation study.Note that a larger value of R corresponds to a more anisotropic process.For the simulations comparing the GSC-g and LZ (Lu and Zimmerman, 2005)tests in Table 5, data were simulated on a subset of the integer grid, Z . Thep-values for the GG test were approximated using a finite sample statistic (Guanet al., 2004), and we used the lag set in (3.2) and A matrix in (3.3). For theresults involving the LZ test, a test of complete symmetry was performed as anapproximation to the null hypothesis of isotropy. The p-values for the LZ test wereobtained using the CvM* statistic. A nominal level of α = 0 .
05 was maintainedby first testing reflection symmetry at α = 0 .
025 then testing complete symmetryat α = 0 .
025 if the hypothesis of reflection symmetry was not rejected. For theGG test, the moving window dimensions were 3 × × ×
12 and 25 ×
15, respectively.For the simulations in Table 6 comparing the GU (Guan et al., 2004) andMS (Maity and Sherman, 2012) tests, data were simulated at random, uniformlydistributed sampling locations on 10 ×
16 and 10 ×
20 sampling domains. Thelag set, Λ , used for both tests is given in (3.2) with A matrix (3.3), and the p-values for both methods were obtained using the asymptotic χ distribution. For imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
20 sampling domains. Thelag set, Λ , used for both tests is given in (3.2) with A matrix (3.3), and the p-values for both methods were obtained using the asymptotic χ distribution. For imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ZACHARY D. WELLER . . . Correlogram
Euclidean distance C o rr e l a t i on -2 -1 - - R = , θ = 0 x lag y l ag -2 -1 - - R = , θ = 0 x lag y l ag -2 -1 - - R = , θ = 0 x lag y l ag -2 -1 - - R = , θ = 3 π /8 x lag y l ag -2 -1 - - R = , θ = 3 π /8 x lag y l ag Contours of Equal Correlation (0.65)
Effective Range = 6, Geometric Anisotropy
Fig 5: Correlogram and contours of equal correlation for the covariance modelsused in the simulation study. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Fig 5: Correlogram and contours of equal correlation for the covariance modelsused in the simulation study. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY semivariogram estimates in GU, we use independent (product) Gaussian (normal)kernels with a truncation parameter of 1.5. The bandwidth for the Gaussian kernelfor smoothing over lags on the entire field and on moving windows was chosenas w = 0 .
75. We used the empirical bandwidth and the product Epanechnikovkernel given in Maity and Sherman (2012) to implement the MS test. For bothtests, a grid with spacing of 1 was laid on the sampling region. Using this grid,the moving window dimensions for the GU test were 4 × ×
2. For the MS test, B = 100 resamples using the GBBBwere used to estimate the asymptotic variance-covariance matrix.For the results in Tables 7 - 9, we simulated mean 0, Gaussian RFs with ex-ponential covariance function with no nugget, a sill of one, and medium effectiverange ( ξ = 6). Sampling locations were generated randomly and uniformly over a16 ×
10 sampling domain. We use the lag set and A matrix from 3.2 and 3.3, re-spectively, unless otherwise noted. All tests were performed using a nominal levelof α = 0 .
05. For the GU tests, we use product Gaussian kernels with a truncationparameter of 1.5. For the MS tests, we use the default product Epanechnikovkernels with empirical bandwidth specified in Maity and Sherman (2012).The simulation results in Table 7 demonstrate the effects of changing the set oflags for the GU and MS tests. For these simulations, the lag set labeled “normal”corresponds to the lag set given in (3.2). The lag set labeled “long” representsthe lags in (3.2) multiplied by 2.5. Finally, the lag set labeled “more” standsfor the lags in (3.2) with the additional pair of lags { h = (1 . , . , h =( − . , . } . The lags h and h are a pair of lags the create approximate22 . ◦ and 112 . ◦ angles, respectively, with the x -axis (counter-clock wise rotation)and have Euclidean length of approximately 1.22. These were chosen to supple-ment the lag pairs ( h , h ) which have unit length and create 0 ◦ and 90 ◦ angleswith the x -axis and ( h , h ) which have length √ ≈ .
41 and create 45 ◦ and 135 ◦ angles with the x -axis. The lag sets are plotted in Figure 6. The A matrix for the“more” lagset was constructed as in (3.3), where orthogonal lags are contrasted.The p-values were calculated using the asymptotic χ distribution with degreesof freedom based on the number of pairs of lags contrasted. For the GU method,we used a bandwidth of 0.75. The moving window dimensions were 4 ×
2. Forthe MS method, we chose block dimensions of 4 × B = 75 resamplesusing the GBBB to estimate the asymptotic variance-covariance matrix. Table8 demonstrates the effects of changing the block size for the GU and MS tests.For these simulations, the labels “small”, “normal”, and “large” correspond tomoving windows/blocks of size 3 ×
2, 4 ×
2, and 5 ×
3, respectively. Because wesimulated n = 300 uniformly distributed sampling locations on a 16 ×
10 domain,we expect 1.875 sampling locations per unit area. Thus, we expect n b = 11.3, 15,and 28.1 points per block for the small, normal, and large block sizes, respectively.We find that the methods tend to have nominal size when n b ≈ n / = 17 .
3. Forboth tests, we used the lags in (3.2), and the blocks are defined by a grid withspacing 0.5 placed on the sampling region (i.e., a 4 × × spTest software). We performedthe tests using a nominal level of α = 0 .
05, and the p-values were calculatedusing the asymptotic χ distribution. For the GU method, we used a bandwidthof 0.75. For the MS method, we used B = 100 resamples using the GBBB to es-timate the asymptotic variance-covariance matrix. Finally, Table 9 demonstrates imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
05, and the p-values were calculatedusing the asymptotic χ distribution. For the GU method, we used a bandwidthof 0.75. For the MS method, we used B = 100 resamples using the GBBB to es-timate the asymptotic variance-covariance matrix. Finally, Table 9 demonstrates imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ZACHARY D. WELLER -2 -1 . . . Normal x lag y l ag -2 -1 . . . Long x lag y l ag -2 -1 . . . More x lag y l ag Different Lag Sets
Fig 6: The lag sets used for the simulations in Table 7.the effects of changing the bandwidth for the GU test. We use bandwidths of w = 0.65, 0.75, and 0.85. The p-values are calculated using both the asymptotic χ distribution and using a finite sample adjustment similar to the one used byGuan et al. (2004) for gridded sampling locations. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Fig 6: The lag sets used for the simulations in Table 7.the effects of changing the bandwidth for the GU test. We use bandwidths of w = 0.65, 0.75, and 0.85. The p-values are calculated using both the asymptotic χ distribution and using a finite sample adjustment similar to the one used byGuan et al. (2004) for gridded sampling locations. imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY Table 5
Empirical size and power for Guan et al. (2004) [denoted GG] and Lu and Zimmerman (2005)[denoted LZ] for 500 realizations of a mean 0 GRF with gridded sampling locations using anominal level of α = 0 . . Computational time for each method is also included. (a) Sample size of n = 216 gridded samplinglocations.
18 cols ×
12 rows grideffective range
R θ
Method 3 6 120 0 GG 0.05 0.07 0.05LZ 0.04 0.04 0.08 √ √ π GG 0.27 0.31 0.34LZ 0.14 0.12 0.132 π GG 0.77 0.85 0.86LZ 0.29 0.33 0.33
Computational Time for 1 Test
GG 1.11 secondsLZ 1.45 seconds (b) Sample size of n = 375 gridded samplinglocations.
25 cols ×
15 rows grideffective range
R θ
Method 3 6 120 0 GG 0.05 0.06 0.07LZ 0.06 0.07 0.07 √ √ π GG 0.47 0.55 0.55LZ 0.16 0.19 0.182 π GG 0.97 0.99 0.98LZ 0.37 0.43 0.45
Computational Time for 1 Test
GG 7.29 secondsLZ 4.99 seconds imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
GG 7.29 secondsLZ 4.99 seconds imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ZACHARY D. WELLER
Table 6
Empirical size and power for Guan et al. (2004) [denoted GU] and Maity and Sherman (2012)[denoted MS] for 200 realizations of a mean 0 GRF with uniformly distributed samplinglocations using a nominal level of α = 0 . . Computational time for each method is alsoincluded. (a) Sample size of n = 300 uniformly dis-tributed sampling locations.
10 height ×
16 width domaineffective range
R θ
Method 3 6 120 0 GU 0.02 0.04 0.05MS 0.04 0.05 0.04 √ √ π GU 0.12 0.13 0.16MS 0.08 0.07 0.042 π GU 0.37 0.51 0.51MS 0.27 0.23 0.21
Computational Time for 1 Test
GU 2.17 secondsMS 83.40 seconds (b) Sample size of n = 450 uniformly dis-tributed sampling locations.
10 height ×
20 width domaineffective range
R θ
Method 3 6 120 0 GU 0.00 0.04 0.05MS 0.05 0.07 0.03 √ √ π GU 0.09 0.18 0.21MS 0.12 0.06 0.082 π GU 0.55 0.58 0.65MS 0.37 0.23 0.21
Computational Time for 1 Test
GU 4.44 secondsMS 162.35 seconds
Table 7
Effects of changing the lag set. Empirical size and power for Guan et al. (2004) [denoted GU]and Maity and Sherman (2012) [MS] for 100 realizations of a mean 0 GRF with n = 400 uniformly distributed sampling locations. The label “normal” corresponds to the lag set in (3.2) ,while “long” represents using longer lags, and “more” denotes using more lags (see Figure 6).
16 width ×
10 height domainLag Set
R θ
Method normal long more0 0 GU 0.02 0.00 0.01MS 0.03 0.14 0.03 √ π GU 0.19 0.07 0.16MS 0.11 0.24 0.072 π GU 0.56 0.17 0.40MS 0.27 0.33 0.21 imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth
Method normal long more0 0 GU 0.02 0.00 0.01MS 0.03 0.14 0.03 √ π GU 0.19 0.07 0.16MS 0.11 0.24 0.072 π GU 0.56 0.17 0.40MS 0.27 0.33 0.21 imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth November, 2015
ONPARAMETRIC HYPOTHESIS TESTS OF ISOTROPY Table 8
Effects of changing the window/block size. Empirical size and power for Guan et al. (2004)[denoted GU] and Maity and Sherman (2012) [MS] for 200 realizations of a mean 0 GRF with n = 300 uniformly distributed sampling locations. The label “normal” corresponds to thewindow/block size of × , while “small” represents using a smaller window, and “large”denotes using a larger window.
16 width ×
10 height domainWindow/Block Size
R θ
Method small normal large0 0 GU 0.06 0.04 0.01MS 0.03 0.04 0.02 √ Table 9
Effects of changing bandwidth. Empirical size and power for Guan et al. (2004) [denoted GU]for 100 realizations of a mean 0 GRF with n = 400 uniformly distributed sampling locationsusing a nominal level of α = 0 . . (a) P-value: asymptotic χ distribution
16 width ×
10 height domainEffective Range
R θ
Bandwidth 3 6 120 0 0.65 0.00 0.00 0.000.75 0.03 0.06 0.040.85 0.06 0.11 0.16 √ π π (b) P-value: finite sample
16 width ×
10 height domainEffective Range
R θ
Bandwidth 3 6 120 0 0.65 0.02 0.03 0.020.75 0.03 0.06 0.060.85 0.07 0.10 0.09 √ π π imsart-sts ver. 2014/10/16 file: sts-isotropy6.tex date: Friday 6 thth