Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
106,139 characters · 23 sections · 0 citation commands
Spatial Correlation Robust Inference
\setcounter{page}{1}
Prompted by advances in both data availability and theory in economic geography, international trade, urban economics, development and other fields, empirical work using spatial data has become commonplace in economics. These applications highlight the importance of econometric methods that appropriately account for spatial correlation in real-world settings. While important advances have been made, researchers arguably lack practical methods that allow for reliable inference about parameters estimated from spatial data for the wide-range spatial designs and correlation patterns encountered in applied work.\footnote{\citeasnoun{Ibragimov10} , \citeasnoun{Sun_Kim_2012} and \citeasnoun{Bester_Conley_Hansen_Vogelsang_2016}, for instance, find nontrivial size distortions of modern methods even in arguably fairly benign designs, and \citeasnoun{Kelly2019} reports very large distortions under spatial correlations calibrated to real-world data.} This paper takes a step forward in this regard.
Specifically, we consider the problem of constructing a confidence interval (or test of a hypothesized value) for the mean of a spatially-sampled random variable. We propose a confidence interval constructed in the usual way, i.e., as the sample mean plus and minus an estimate of its standard error multiplied by a critical value. The novelty is that the standard error and critical value are constructed so the resulting confidence interval has the desired large-sample coverage probability (say, $95\%$) for a relatively wide range of correlation patterns and spatial designs. The analysis is described for the mean, but the required modifications for regression coefficients or parameters in GMM settings follow from standard arguments.
To be more precise, suppose that a random variable $y$ is associated with a location $s\in \mathcal{S}$, where $\mathcal{S}\subset \mathbb{R}^{d}$. Figure (ref) shows three one-dimensional ($d=1$ ) spatial designs. Panel (a) shows the familiar case of regularly spaced locations, corresponding to the standard time series setting; panels (b) and (c) show randomly selected locations drawn from a density $g$, where $g$ is uniform in panel (b) and triangular in panel (c). Figure 2 shows two geographic examples, so $d=2$, for the U.S. state of Texas. In panel (a), the locations are randomly selected from a uniform distribution, while in panel (b) they are more likely to be sampled from areas with high economic activity, here measured by light intensity as seen from space.\footnote{ The light data are from \citeasnoun{Henderson_et_al_light}.} In much of our analysis, we will assume that locations are i.i.d. draws from a distribution with density $g$, and so will encompass the irregularly spaced time series and Texas examples.
Adding some notation, suppose
where $y_{l}$ is associated with the spatial location $s_{l}$, $\mu$ is the mean of $y_{l}$, and $u_{l}$ is an unobserved error, assumed to be covariance stationary with mean zero and covariance function $\mathbb{E} [u(r)u(s)]=\sigma_{u}(r-s)$. Let $\overline{y}$ denote the sample mean, and consider the usual t-statistic
where $\hat{\sigma}^{2}$ is an estimator for the variance of $\sqrt{n}( \overline{y}-\mu)$. Tests of the null hypothesis $H_{0}:\mu=\mu_{0}$ reject when $|\tau|>\func{cv}$, where $\func{cv}$ is the critical value, and the corresponding confidence interval for $\mu$ has endpoints $\overline{y}\pm \func{cv}\hat{\sigma}/\sqrt{n}$. Inference methods in this class differ in their choice of $\hat{\sigma}^{2}$ and critical value $\func{cv}$.
The case of regularly-spaced time series observations (panel (a) of Figure (ref)) is the most well-studied version of this problem. Here $\func{Var}(\sqrt{n}(\overline{y}-\mu ))$ is the long-run variance of $y$. Classic choices for $\hat{\sigma}^{2}$ are kernel-based consistent estimators such as those proposed in \citeasnoun{Newey87} and \citeasnoun {Andrews91}, and associated standard normal critical values. A more recent literature initiated by \citeasnoun{Kiefer00} and \citeasnoun{Kiefer05} accounts for the sampling uncertainty of kernel-based $\hat{\sigma}^{2}$ by considering \textquotedblleft fixed-$b$\textquotedblright\ asymptotics where the bandwidth is a fixed fraction of the sample size, which leads to a corresponding upward adjustment of the critical value. Closely related are projection estimators of $\hat{\sigma}^{2}$ where the number of projections is treated as fixed in the asymptotics, as in M\"{u}ller (2004, 2007)\nocite {Muller04}\nocite{Muller07c}, \citeasnoun{Phillips05}, \citeasnoun{Sun13}, and others, leading to Student-t critical values. These newer methods are found to markedly improve size control under moderate serial correlation compared to inference based on standard normal critical values.
In the general spatial case, the variance of $\overline{y}$ depends on the correlation between all of the observations, and this in turn depends on two distinct features of the problem. The first is the correlation between observations at arbitrary locations (say $r$ and $s$); this is given by the covariance function $\sigma_{u}(r-s)$. The second feature is which locations in $\mathcal{S}$ are likely to be sampled; this is given by the spatial density $g$. Only the first of these features is important in the regularly-spaced time series example because the locations do not vary from one application to the next.
Most existing suggestions for spatial inference are derived under the assumption that the locations are (asymptotically) uniformly distributed, corresponding to a constant density $g$: This includes the consistent kernel-based estimator in \citeasnoun{Conley99}, the spatial analogue of the fixed- $b$ kernel approach analyzed in \citeasnoun{Bester_Conley_Hansen_Vogelsang_2016}, as well as the spatial projection-based estimator put forward in \citeasnoun {Sun_Kim_2012}. Exceptions include \citeasnoun{Kelejian2007}, who derive a consistent kernel for $\hat{\sigma}^{2}$ under assumptions that can accommodate arbitrary locations $s_{l}$, and the cluster approach suggested by Ibragimov and M\"{u}ller (2010, 2015)\nocite{Ibragimov10}\nocite {Ibragimov15} and \citeasnoun{Bester11} (also see \citeasnoun{Cao2020}).
This paper makes progress over this literature by developing a method that (i) accounts for sampling uncertainty in $\hat{\sigma}^{2}$ in a spatial context while allowing for nonuniform spatial densities $g$; (ii) is valid under generic weakly correlated $u_{l}$; (iii) also controls size under a restricted but nonparametric form of strongly correlated $u_{l}$. The last property sets it apart from all previously mentioned methods; in a time series setting, \citeasnoun{Robinson05} and \citeasnoun{Muller14} derive inference under parametric forms of strong dependence, and \citeasnoun{Dou_2019} derives optimal inference under a non-parametric form of strong dependence under a simplifying Whittle-type approximation to the implied covariance matrices.
Our method works as follows: First, a benchmark parametric model is specified for the covariance function, say $\sigma_{u}^{0}(\cdot)= \sigma_{u}^{0}(\cdot|c)$, where $c$ is a persistence parameter with larger values indicating less dependence. For a given lower bound on the persistence parameter, say $c_{0}$, a hypothetical covariance matrix for $ (y_{1},...,y_{n})^{\prime}$ is constructed using $\sigma_{u}^{0}( \cdot|c_{0}) $ evaluated at the actual sample locations $(s_{1},...,s_{n})$. The eigenvectors of the demeaned version of this covariance matrix are the (population) principal components of the residuals $\hat{u}_{l}=y_{l}- \overline{y}$ under $\sigma_{u}^{0}(\cdot|c_{0})$, and the sample variance of $q$ of these principal components is the estimator $\hat{\sigma}^{2}$. The critical value is chosen to ensure coverage for all $c\geq c_{0}$. The number of principal components $q$ is chosen to minimize the expect length of the confidence interval in the model where $u_{l}$ is i.i.d. For shorthand, we refer to the method as spatial correlation principal components, abbreviated SCPC.
Intuitively, variance estimators $\hat{\sigma}^{2}$ that are quadratic forms in $\hat{u}$ are sums of squares of weighted averages of $\hat{u}$. Under spatial correlation, most weighted averages are less variable than $ \overline{y}$, leading to a downward biased $\hat{\sigma}^{2}$. SCPC selects the linear combinations of $\hat{u}$ that are most variable, so that the bias is as small as possible in the benchmark model with parameter $c_{0}$.
The remainder of the paper studies this method. Section 2 provides the specifics for SCPC. These specifics raise a variety of issues that are the focus of the remaining sections of the paper. In particular, Section 3 lays out the analytic framework used to study the large-sample and finite-sample Gaussian properties of spatial t-statistics. We use the framework to analyze SCPC, but several of the results in Section 3 encompass other methods, notably \textquotedblleft fixed-$b$\textquotedblright\ kernel-based methods, and general projection estimators with a fixed number of basis functions. We find that in contrast to the regularly spaced time series case, such t-statistics with analogously adjusted critical values are not generically valid under weak correlation as soon as the spatial density function is not uniform. We develop an alternative approach to the construction of critical values that restores validity, and this is used for SCPC inference. Section 4 thus shows that SCPC has the desired large-sample coverage probability under generic weak correlation. Moreover, Section 4 provides a set of (easily verifiable) sufficient conditions that guarantee coverage under arbitrary mixtures of a set of strong correlation patterns in a finite-sample Gaussian setting. Section 4 also investigates the finite-sample coverage probability of SCPC confidence sets when there is heteroskedasticity across locations or measurement errors in locations --- two problems faced in some applications. Section 5 addresses the question of efficiency of SCPC by computing a lower bound on the expected length of confidence intervals for any inference method that controls coverage in a particular class of spatial correlations. Comparing the expected length of SCPC to this lower bound provides a measure of the efficiency of the method. Section 6 compares the properties of SCPC to other methods that have been proposed in the literature, and the results suggest that SCPC dominates these methods over the range of covariance functions and spatial designs considered. Section 7 discusses extensions and implementation issues. First, it discusses how the results developed in the body of the paper for inference about the population mean can be applied to inference problems about regression coefficients or parameters in GMM models. It then discusses two important computational issues involved in computing the critical value and computing the required eigenvectors for the construction of SCPC in very large-$n$ applications. Finally, Section 7 provides a sketch of the generalization of the SCPC method to multivariate (F-test) settings. Proofs are collected in the appendix.
This section provides details for computing the SCPC t-statistic, critical value and associated confidence interval. The construction of SCPC raises a variety of questions about its properties, many of which are posed here and discussed in detail in the remaining sections of the paper.
The construction of the SCPC t-test and confidence interval involves, among other things, various covariance matrices and probability calculations. We stress at the outset that these are used to describe the required calculations, and they are not assumptions about the probability distribution of the data under study. Those assumptions will be listed in Section 3 and, it will turn out, are significantly more general than what would follow from the description in this section.
Let $\mathbf{y}=(y_{1},y_{2},...,y_{n})^{\prime }$ and similarly for $ \mathbf{s}=(s_{1},s_{2},...,s_{n})^{\prime },$ $\mathbf{u} =(u_{1},u_{2},...,u_{n})^{\prime }$ and the vector of residuals $\mathbf{ \hat{u}}=(\hat{u}_{1},\hat{u}_{2},...,\hat{u}_{n})^{\prime }$. Let $\mathbf{l }$ denote an $n\times 1$ vector of $1$s, and $\mathbf{M}=\mathbf{I}-\mathbf{l }(\mathbf{l}^{\prime }\mathbf{l})^{-1}\mathbf{l}^{\prime }$. Consider a benchmark model for $u_{l}$ with a parametric covariance function $\func{Cov} (u(r),u(s))=\sigma _{u}^{0}(r-s|c)$, where smaller values of the scalar parameter $c$ indicate stronger correlations. In the following, we focus on the simple Gaussian exponential (`AR(1)') model where $\sigma _{u}^{0}(r-s|c)=\exp (-c||r-s||)$ for $c>0$. Let $\mathbf{\Sigma }(c)$ denote the $n\times n$ covariance matrix with $\mathbf{\Sigma }(c)_{ij}=\exp (-c||s_{i}-s_{j}||)$, so that $\mathbf{\Sigma }(c)$ is the covariance matrix of $u(s)$ evaluated at the sample locations $\mathbf{s}.$ Let $c_{0}$ denote a pre-determined value of $c$ that is meant to capture an upper bound on the spatial persistence in the data. (The choice of $c_{0}$ is discussed below). Let $\mathbf{r}_{1},\mathbf{r}_{2},...,\mathbf{r}_{n}$ denote the eigenvectors of $\mathbf{M}\mathbf{\Sigma }(c_{0})\mathbf{M}$ corresponding to the eigenvalues ordered from largest to smallest, and normalized so that $ n^{-1}\mathbf{r}_{j}^{\prime }\mathbf{r}_{j}=1$ for all $j$. The scalar variable $n^{-1/2}\mathbf{r}_{j}^{\prime }\mathbf{\hat{u}}$ has the interpretation as the $j$th population principle component of $ \mathbf{\hat{u}|}$$\mathbf{s}$ $\sim \mathcal{N}(\mathbf{0},\mathbf{M} \mathbf{\Sigma }(c_{0})\mathbf{M)}$. The SCPC estimator of $\sigma ^{2}$ based on the first $q$ of these principal components is
and the corresponding SCPC t-statistic is
The critical value $\func{cv}_{\text{SCPC}}(q)$ of the level-$\alpha$ SCPC test is chosen so that size is equal to $\alpha$ under the Gaussian benchmark model with $c\geq c_{0}$. That is, $\func{cv}_{\text{SCPC}}(q)$ satisfies
where $\mathbb{P}_{\mathbf{\Sigma}(c)}^{0}$ means that the probability is computed in the benchmark model $\mathbf{y|s}\sim\mathcal{N}(\mu_{0}\mathbf{l },\mathbf{\Sigma}(c))$.
The final ingredient in the method is the choice of $q$. Let $\mathbb{E} ^{1}[2\hat{\sigma}_{\text{SCPC}}(q)\func{cv}_{\text{SCPC}}(q)|\mathbf{s]}$ denote the expected length of the confidence interval constructed using $ \tau _{\text{SCPC}}(q)$ under the Gaussian i.i.d. model $\mathbf{y|s}\sim \mathcal{N}(\mathbf{l}\mu ,\mathbf{I})$. (The superscript \textquotedblleft 1\textquotedblright\ on $\mathbb{E}$ differentiates this from the benchmark model with covariance matrix $\mathbf{\Sigma }(c)$.) SCPC chooses $q_{\func{ SCPC}_{\text{{}}}}$ to make this length as small as possible, that is $q$ solves
with the equality exploiting that $q\hat{\sigma}_{\text{SCPC}}^{2}(q)\mathbf{ |s}\sim \chi _{q}^{2}$ in the Gaussian i.i.d. model.
U.S. states spatial correlation designs. Before making two additional remarks about the SCPC\ method, we introduce a set of spatial correlation designs that will be used throughout the paper. The idea is to consider a set of real world designs to learn about the usefulness of the SCPC and other methods in practice. In particular, we randomly draw $ n=500$ locations within the boundaries of the 48 contiguous states of the U.S. (we also considered $n=1000$ draws, and found nearly identical results in all exercises). The density of locations $g$ within each state is either uniform ($g_{\text{uniform}}$), or it is proportional to light measured from space ($g_{\text{light}}$) as a proxy for economic activity. We draw five sets of 500 independent locations under each density $g\in\{g_{\text{uniform} },g_{\text{light}}\}$ and $\bar{\rho}_{0}\in\{0.02,0.10\}$ for each state, for a total of 240 (= 48 states $\times$ 5 location draws) sets of locations $\{s_{l}\}_{l=1}^{500}$ and associated covariances under each of the four $ (g,\bar{\rho}_{0})$ pairs.
This section outlines a large-sample framework used to study SCPC and other spatial t-statistics. The first two subsections introduce notation and the asymptotic sampling framework. With these in hand, the remainder of the section summarizes the large-sample distribution of various statistics including the SCPC and kernel-based t-statistics. Proofs are provided in the appendix.
Some of this notation has been introduced earlier, but is repeated here for easy reference.
The sample mean is denoted by $\overline{y}_{n}$, where here and elsewhere we append the subscript $n$ for clarity in the asymptotic analysis. The residual is $\hat{u}_{l}=y_{l}-\overline{y}_{n}$. Let $\mathbf{y} _{n}=(y_{1},...,y_{n})^{\prime}$, and similarly for $\mathbf{u}_{n}$, $\hat{ \mathbf{u}}_{n}$ and $\mathbf{s}_{n}$. The vector $\mathbf{l}_{n}$ is a $ n\times1$ vector of 1s, and $\mathbf{M}_{n}=\mathbf{I}_{n}-\mathbf{l}_{n}( \mathbf{l}_{n}^{\prime}\mathbf{l}_{n})^{-1}\mathbf{l}_{n}$, so that $\mathbf{ \hat{u}}_{n}=\mathbf{M}_{n}\mathbf{u}_{n}$.
Generically, we consider estimators $\hat{\sigma}_{n}^{2}$ that are quadratic forms in $\hat{\mathbf{u}}_{n}$. Let $\mathbf{Q}_{n}$ be a positive semidefinite matrix with $\mathbf{Q}_{n}\mathbf{l}_{n}=\mathbf{0}$. We consider estimators of the form
where the final equality follows from $\mathbf{Q}_{n}\mathbf{l}_{n}=\mathbf{0 }$.
Two leading examples of estimators in this class are kernel-based estimators and orthogonal-projections estimators. For kernel-based estimators, let $ k(r,s)$ denote a positive semi-definite kernel, $k:\mathcal{S}\times \mathcal{S}\mapsto \mathbb{R}$. Let $\mathbf{K}_{n}$ denote an $n\times n$ matrix with $(l,\ell )$ element equal to $k(s_{l},s_{\ell })$ and let $ \mathbf{Q}_{n}=\mathbf{M}_{n}\mathbf{K}_{n}\mathbf{M}_{n}$. Then $\hat{\sigma }_{n}^{2}=n^{-1}\sum_{l}\sum_{\ell }k(s_{l},s_{\ell })\hat{u}_{l}\hat{u} _{\ell }=$ $n^{-1}\hat{\mathbf{u}}_{n}^{\prime }\mathbf{Q}_{n}\hat{\mathbf{u} }_{n}$. For orthogonal-projection estimators, let $\mathbf{\hat{W}}_{n}$ be an $n\times q$ matrix with $j$th column given by $\mathbf{\hat{w}}_{j}$ satisfying $n^{-1}\mathbf{\hat{W}}_{n}^{\prime }\mathbf{\hat{W}}_{n}=q^{-1} \mathbf{I}_{q}$ and $\mathbf{\hat{W}}_{n}^{\prime }\mathbf{l}_{n}=\mathbf{0}$ (the `hat' notation is a reminder that $\mathbf{\hat{W}}$ depends on the locations $\mathbf{s}_{n}$, which are random). With $\mathbf{Q}_{n}=\mathbf{ \hat{W}}_{n}\mathbf{\hat{W}}_{n}^{\prime }$, the orthogonal projection estimator is $\hat{\sigma}_{n}^{2}=\sum_{j=1}^{q}(n^{-1/2}\mathbf{\hat{w}} _{j}^{\prime }\hat{\mathbf{u}}_{n})^{2}=n^{-1}\hat{\mathbf{u}}_{n}^{\prime } \mathbf{Q}_{n}\hat{\mathbf{u}}_{n}$. The SCPC estimator is an orthogonal-projection estimator using the first $q$ eigenvectors of $\mathbf{ M}_{n}\mathbf{\Sigma }(c_{0})\mathbf{M}_{n}$, scaled to have length $1/\sqrt{ q},$ as the columns of $\mathbf{\hat{W}}_{n}$.
For quadratic form estimators $\hat{\sigma}_{n}^{2}(\mathbf{Q}_{n})$, under the null hypothesis the squared t-statistic is a ratio of quadratic forms in $\mathbf{u}_{n}$
The spatial locations $s$ are chosen from $\mathcal{S}$, a compact subset of $\mathbb{R}^{d}$. Sample locations are selected as i.i.d. draws from a distribution $G$ with density $g$, where $g(s)$ is continuous and positive for all $s\in\mathcal{S}$.
The average pairwise correlation of $y$, conditional on the sample locations is $\bar{\rho}_{n}=\frac{1}{n(n-1)}\sum_{l=1}^{n}\sum_{\ell\neq l}\func{Cor} \left(y_{l},y_{\ell}\left\vert \mathbf{s}_{n}\right.\right)$. When $ \overline{\rho}_{n}=0$, $\mathbf{y}_{n}\left\vert \mathbf{s}_{n}\right.$ is white noise. When $\overline{\rho}_{n}=O_{p}(1)$ (and not $o_{p}(1)$), we will say the process exhibits strong correlation. When $ \overline{\rho}_{n}=O_{p}(1/c_{n}^{d})$ where $c_{n}$ is a sequence of constants with $c_{n}\rightarrow\infty$, we follow \citeasnoun{Lahiri_2003} and say the process exhibits weak correlation.
The following asymptotic framework, adapted from \citeasnoun{Lahiri_2003}, is useful for representing weak and strong correlation. Let $B$ be a zero-mean stationary random field on $\mathbb{R}^{d}$ with continuous covariance function $\mathbb{E}[B(s)B(r)]=\sigma_{B}\left(s-r\right)$, and $B$ and $ \{s_{l}\}_{l=1}^{n}$ are independent. To avoid pathological cases, we further assume $\int\sigma_{B}(s)ds>0$ and that $B$ is nonsingular in the sense that $\inf_{||f||=1}\int\int f(r)f(s)\sigma_{B}(s-r)dG(r)dG(s)>0$ with $||f||^{2}=\int f^{2}(s)dG(s)$.
Let $c_{n}$ denote a sequence of constants with either $c_{n}\rightarrow \infty $ or $c_{n}=c>0$. We consider a triangular-array framework with $ u_{l}=B(c_{n}s_{l})$ for $s_{l}\in \mathcal{S}$, so that $\sigma _{u}(s)=\sigma _{B}(c_{n}s)$. The sequence $c_{n}$ determines the `infill' and `outfill' nature of the asymptotics. To see this, note that the volume of the relevant domain for the random field $B$ is $c_{n}^{d}\func{vol}( \mathcal{S)}$, where $\func{vol}(\mathcal{S)}$ is the volume of $\mathcal{S} . $ The average number of sample points per unit of volume is then $ n/(c_{n}^{d}\func{vol}(\mathcal{S)}).$ If $c_{n}^{d}\propto n,$ the volume of the domain is increasing, while the number of points per unit of volume is not; this is the usual outfill asymptotic sampling scheme. On the other hand, when $c_{n}=c$, a constant, the volume of the domain is fixed, and the number of points per unit of volume is proportional to $n$; this is the usual infill sampling. Finally, when $c_{n}\rightarrow \infty $ with $ c_{n}^{d}=o(n)$ the sampling scheme features both infill and outfill asymptotics. A calculation shows that $\overline{\rho } _{n}=O_{p}(1/c_{n}^{d})$, so the sequence $c_{n}$ characterizes weak and strong correlation as described above. With this background, let $ a_{n}=c_{n}^{d}/n$; we will assume that $a_{n}\rightarrow a\in \lbrack 0,\infty )$.
Finally, we specify a set of weighting functions. To simplify the problem, we initially consider weights that are nonrandom. For $j=1,\ldots ,q$, let $ w_{j}:\mathcal{S}\mapsto \mathbb{R} $ denote a set of continuous functions that satisfy $\int w_{j}(s)dG(s)=0$ and $\int w_{j}^{2}(s)dG(s)>0$. We introduce the following notation involving these functions: $\mathbf{w}(s)$ is a $q\times 1$ vector-valued continuous function with $\mathbf{w}(s)=(w_{1}(s),...,w_{q}(s))^{\prime }$; $ \mathbf{w}^{0}(s)=(1,\mathbf{w}(s)^{\prime })^{\prime }$; $\mathbf{W}_{n}$ is a $n\times q$ matrix with $l$th row given by $\mathbf{w}(s_{l})^{\prime }$ , and $\mathbf{W}_{n}^{0}$ is a $n\times (q+1)$ matrix with $l$th row given by $\mathbf{w}^{0}(s_{l})^{\prime }$ so that $\mathbf{W}_{n}^{0}=[\mathbf{l} _{n},\mathbf{W}_{n}]$.
With this background, we now present the large-sample analysis.
As is evident from equation ((ref)) the squared t-statistic is a ratio of squares of weighted average of the elements of $\mathbf{u}_{n}$. This subsection discusses the large-sample distribution of such weighted averages. These results involve weak converge (i.e., convergence in distribution) where our interest lies in these limits conditional on the locations $\mathbf{s}_{n}$. With this in mind, for $\mathbf{X}_{n}$ and $ \mathbf{X}$ $p$-dimensional random vectors, we use the notation $\mathbf{X} _{n}|\mathbf{s}_{n}\Rightarrow _{p}\mathbf{X}$ to denote $\mathbb{E}[h( \mathbf{X}_{n})|\mathbf{s}_{n}]\overset{p}{\rightarrow }\mathbb{E}[h(\mathbf{ X})]$ for any bounded continuous function $h:\mathbb{R}^{p}\mapsto \mathbb{R} $. This notion of weak convergence in probability is slightly weaker than almost sure weak convergence of conditional distributions, but still ensures that the limiting distribution is not induced by the randomness in the locations $\mathbf{s}_{n}$.
Lemma (ref) characterizes the large-sample behavior of sums of the form $\sum_{l=1}^{n}\mathbf{w}^{0}(s_{l})u(s_{l}).$ For the weak correlation result, we invoke the mixing and moment assumptions of \citeasnoun {Lahiri_2003} on $B$ that underlie his Theorem 3.2.
This section presents a useful representation for the limiting distribution of $\tau_{n}^{2}(\mathbf{W}_{n}\mathbf{W}_{n}^{\prime})$ under the assumptions of Lemma (ref).
For SCPC and other estimators, the weights in $\mathbf{w}(s)$ are estimated using the sample locations $\mathbf{s}_{n}$. The conditions under which Lemma (ref) continues to hold for such estimated weights is given in the following theorem.
This subsection discusses how these results can be generalized so they apply to kernel-based variance estimators, $\hat{\sigma}_{n}^{2}(\mathbf{M}_{n} \mathbf{K}_{n}\mathbf{M}_{n})$ and associated t-statistics $\tau _{n}^{2}( \mathbf{M}_{n}\mathbf{K}_{n}\mathbf{M}_{n})$, where the $n\times n$ matrix $ \mathbf{K}_{n}$ has $(l,\ell )$ element equal to $k(s_{l},s_{\ell })$ for a positive semidefinite continuous kernel $k:\mathcal{S\times S}\mapsto \mathbb{R} $. Since in our framework, $s_{l}\in \mathcal{S}$ for a fixed sampling region $\mathcal{S}$, and $k$ does not depend on $n$, these kernel estimators are spatial analogues of fixed-$b$ time series long-run variance estimators considered by \citeasnoun{Kiefer05}, as also investigated by \citeasnoun {Bester_Conley_Hansen_Vogelsang_2016}.
Let $\mathbf{\hat{K}}_{n}=\mathbf{M}_{n}\mathbf{K}_{n}\mathbf{M}_{n}$, and note that the $(l,\ell )$ element of $\mathbf{\hat{K}}_{n}$ is $\hat{k} _{n}(s_{l},s_{\ell })$ with
To begin, consider a simpler problem using a kernel that replaces the sample means in ((ref)) with populations means
By Mercer's Theorem, $\overline{k}(r,s)$ has the representation
where $\{\lambda _{i},\varphi _{i}\}$ are the eigenvalues and eigenfunctions of $\overline{k}$, with eigenvalues ordered from largest to smallest, $ \mathbb{\int }\varphi _{i}(s)dG(s)=0$ and $\mathbb{\int }\varphi _{i}(s)\varphi _{j}(s)dG(s)=\mathbf{1[}i=j]$.
Consider the problem with a truncated version of $\overline{k}$,
We can directly apply Theorem (ref) using $ w_{j}(s)=\lambda_{j}^{1/2}\varphi_{j}(s)$. Specifically, let $\mathbf{\bar{K} }_{n,q}$ be an $n\times n$ matrix with $(l,\ell)$ element equal to $ \overline{k}_{q}(s_{l},s_{\ell})$. Then $\mathbf{u}_{n}^{\prime}\mathbf{\bar{ K}}_{n,q}\mathbf{u}_{n}=\mathbf{u}_{n}^{\prime}\mathbf{W}_{n}\mathbf{W} _{n}^{\prime}\mathbf{u}_{n}$ so that $\tau_{n}^{2}(\mathbf{\bar{K}} _{n,q})=\tau_{n}^{2}(\mathbf{W}_{n}\mathbf{W}_{n}^{\prime})$, and $\mathbb{P} \left(\tau_{n}^{2}(\mathbf{\bar{K}}_{n,q})>\func{cv}^{2}|\mathbf{s} _{n}\right)\overset{p}{\rightarrow}\mathbb{P}\left(Z_{0}^{2}>\sum_{i=1}^{q}(- \frac{\omega_{i}}{\omega_{0}})Z_{i}^{2}\right)$ by Theorem (ref).
To extend this result to the original problem, it is useful to reformulate it in terms of eigenvalues of linear operators. Specifically, denote by $ \mathcal{L}_{G}^{2}$ the Hilbert space of functions $\mathcal{S}\mapsto \mathbb{R} $ with inner product $\langle f_{1},f_{2}\rangle =\int f_{1}(s)f_{2}(s)dG(s)$ . Normalize $\mathbf{\Omega }_{wc}=\kappa \mathbf{V}_{1}+(1-\kappa )\mathbf{V }_{2}$, as in ((ref)). A tedious but straightforward calculation (see ((ref))\ in the appendix) shows that the eigenvalues $\omega _{i}$ of $\mathbf{A}=\mathbf{D}(\func{cv})\mathbf{\Omega }$ with $\mathbf{\Omega =\{\Omega }_{sc}\mathbf{,\Omega }_{wc}\}$ are also the eigenvalues of finite rank self-adjoint linear operators $\mathcal{L} _{G}^{2}\mapsto \mathcal{L}_{G}^{2}$, namely $R_{sc}T_{q}R_{sc}$ and $ R_{wc}T_{q}R_{wc}$ in the strong and weak correlation case, respectively, where
This suggests that the limiting rejection probability for the original non-truncated $\bar{k}$ might be characterized by the (potentially infinite) number of eigenvalues of the operators $RTR:\mathcal{L}_{G}^{2}\mapsto \mathcal{L}_{G}^{2}$ with $R\in \{R_{wc},R_{sc}\}$, where
The following theorem shows this to be the case, and it also includes the generalization to sample demeaned kernels ((ref)) instead of ((ref)).
The proof of Theorem (ref) involves showing that in large samples, the difference between the eigenfunctions of the sample demeaned kernel ((ref)) and the population demeaned kernel ((ref)) becomes small. The following lemma extends and adapts previous results by \citeasnoun{Rosasco2010} to the case of sample demeaned kernels.
Part (a) shows convergence of the eigenspace corresponding to unique eigenvalues, and part (b) shows convergence of the eigenvalues.
Beyond its use in the proof of Theorem (ref), Lemma (ref) can be used to establish the large sample distribution of the SCPC t-statistic for nonrandom $q$ and critical value $\func{cv}$. Note that in this application of Lemma (ref), we are interested in the eigenfunctions of the covariance kernel $k^{0}(r,s)= \sigma_{u}^{0}(r-s|c_{0})$ of the benchmark model, rather than the eigenfunctions of a kernel that defines a kernel-based variance estimator.
Recall from Section (ref) that $\mathbf{r}_{i}$ is the eigenvector of $\mathbf{M}_{n}\mathbf{\Sigma }_{n}(c_{0})\mathbf{M}_{n}$ corresponding to the $i$th largest eigenvalue, normalized to satisfy $n^{-1}\mathbf{r} _{i}^{\prime }\mathbf{r}_{i}=1$. Let $\varphi _{i}^{0}$ be the eigenfunction of the kernel $\overline{k}^{0}(r,s)$ corresponding to the $i$th largest eigenvalue $\lambda _{i}^{0}$, where $k^{0}(r,s)=\sigma _{u}^{0}(r-s|c_{0})$ and $\bar{k}^{0}$ is the demeaned version of $k^{0}$ in analogy to ((ref)). Lemma (ref) and a slightly extended version of Theorem (ref) (see Lemma (ref) in the appendix) then yields the following corollary.
This section presents two results on size control of spatial t-statistics, the first asymptotic and the second a finite-sample result, and applies these to SCPC.
As discussed above (see equation ((ref))), under weak correlation, the asymptotic rejection probability of $\tau _{n}$ for finite $ q$ can be studied via $\mathbf{\Omega }_{wc}(\kappa )=\kappa \mathbf{V} _{1}+(1-\kappa )\mathbf{V}_{2}$, where the covariance function of $u$ and the sequence $c_{n}$ affects the large-sample distribution of $\tau _{n}$ only through the scalar $\kappa \in \lbrack 0,1)$. Thus, if $\overline{\func{ cv}}$ is such that $\sup_{0\leq \kappa <1}\mathbb{P}\left( \sum_{i=0}^{q}\omega _{i}(\kappa ,\overline{\func{cv}})Z_{i}^{2}>0\right) =\alpha $, where $\{\omega _{i}(\kappa ,\overline{\func{cv}})\}_{i=0}^{q}$ are the eigenvalues of $\mathbf{A}(\kappa ,\overline{\func{cv}})=\mathbf{D}( \overline{\func{cv}})\mathbf{\Omega }_{wc}(\kappa )$, then setting $\func{cv} _{n}\geq \overline{\func{cv}}$ for all $n$ yields inference that is asymptotically robust under all forms of weak correlation covered by Theorem (ref) (ii). In the case of a kernel-based variance estimator, the same holds as long as $\overline{\func{cv}}$ satisfies $ \sup_{0\leq \kappa <1}\mathbb{P}\left( \sum_{i=0}^{\infty }\omega _{i}(\kappa ,\overline{\func{cv}})Z_{i}^{2}>0\right) =\alpha $ where $ \{\omega _{i}(\kappa ,\overline{\func{cv}})\}_{i=0}^{\infty }$ are the eigenvalues of the linear operator $L(f)(s)=\int \sqrt{\kappa +(1-\kappa )g(s)}\left( 1-\overline{\func{cv}}^{2}\overline{k}(s,r)\right) \sqrt{\kappa +(1-\kappa )g(r)}f(r)dG(r)$.
The value $\overline{\func{cv}}$ depends on the spatial density $g$, which can be seen directly by inspecting the form of $\mathbf{\Omega}_{wc}$ and the operator $L$. In principle, one could use these expressions to estimate $ \overline{\func{cv}}$ directly. But this would involve estimates of the spatial density $g$, which leads to difficult bandwidth an other choices. We now discuss a simpler approach.
Consider a benchmark model $B^{0}$ that satisfies the assumptions of Theorem (ref) (ii), such as the Gaussian exponential model introduced in Section (ref). Let $\sigma_{B}^{0}$ denote the covariance kernel of $B^{0}$, and suppose $c_{n,0}$, is chosen so that $ a_{n,0}=c_{n,0}^{d}/n\rightarrow a_{0}=0.$ For instance, $c_{n,0}=c_{0}>0$ satisfies this condition, as does $c_{n,0}=n^{1/d}/\log(n)$. Note that for this model $\kappa=0$. Suppose $\func{cv}_{n}=\func{cv}_{n}(\mathbf{s}_{n})$ satisfies
where $\mathbb{P}_{\mathbf{\Sigma}(c)}^{0}$ is computed under the benchmark model, that is under $\mathbf{u}_{n}|\mathbf{s}_{n}\sim\mathcal{N}(0,\mathbf{ \Sigma}(c))$ with $\mathbf{\Sigma}(c)$ the covariance matrix of $ (B^{0}(cs_{1}),...,B^{0}(cs_{n}))^{\prime}$.
The intuition for Theorem (ref) is as follows. The critical value $ \func{cv}_{n}$ in ((ref)) is valid in the benchmark model for all $c\geq c_{n,0}$ and $n$. Thus, it is also valid along arbitrary sequences $c_{n}\geq c_{n,0}$. Since the $c_{n,0}$ model has $\kappa=0$, there exists sequences $c_{n}\geq c_{n,0}$ that induce any $ \kappa\in\lbrack0,1)$ in the benchmark model; thus different sequences $ c_{n} $ in the benchmark model recreate any possible limit distribution under generic weak correlation, so that size control in the benchmark model for all $c\geq c_{n,0}$ translates into size control under generic weak correlation.
For SCPC, the benchmark covariance kernel for $B^{0}$ is exponential $\sigma _{B}^{0}(r,s)=\exp (-||r-s||)$ and (from equation ((ref))) the critical value is chosen to satisfy ((ref)) with equality. Thus, with a fixed value of $c_{0}$, the SCPC t-test $\tau _{\text{SCPC}}(q)$ controls size in large samples under generic weak correlation.\footnote{ Technically, the SCPC choice of $q$ in ((ref)) is also a function of the locations of $\mathbf{s}_{n}$, so $q_{\text{SCPC}}$ is random. However, the argument that establishes Theorem (ref) can be extended under this complication as long as $q_{\text{SCPC}}\leq q_{\max }$ almost surely for some finite and fixed $q_{\max }$. See Theorem (ref) in the appendix for a formal statement.}
In addition and by construction, the SCPC critical value is chosen to satisfy the size constraint for all values of $c\geq c_{0}$ in the benchmark model. Thus, size is controlled by construction also in strong-correlation models with exponential covariance kernels for all $c\geq c_{0}$.
The asymptotic results of the last subsection are comforting, but in finite samples, the robustness of a spatial t-statistic with critical value chosen according to ((ref)) still depends on the choice of $c_{n,0}$ and the benchmark model. This motivates investigating size control in finite samples, which potentially includes `strong' correlation cases.
We restrict attention to Gaussian models where $\mathbf{y}\sim \mathcal{N}( \mathbf{l}\mu ,\mathbf{\Sigma })$ for some $\mathbf{\Sigma }$ and implicitly condition on $\mathbf{s}$, and we also omit the dependence on $n$ to ease notation. In this finite sample conditional framework, the distinction between $\mathbf{W}$ and $\mathbf{\hat{W}}$ is immaterial, so for simplicity, we write $\tau ^{2}(\mathbf{W}\mathbf{W}^{\prime })$ for the t-statistic.\footnote{ This also covers kernel variance estimators by setting $q=T-1$ and using the Choleksy decomposition $\mathbf{MKM}=\mathbf{WW}^{\prime }$.}
Let $\mathcal{V}$ denote a set of covariance matrices. A test using the t-statistic $\tau ^{2}(\mathbf{W}\mathbf{W}^{\prime })$ with critical value $ \func{cv}$ is robust for\ $\mathcal{V}$ if $\sup_{\mathbf{\Sigma }\in \mathcal{V}}\mathbb{P}_{\mathbf{\Sigma }}(\tau ^{2}(\mathbf{W}\mathbf{W} ^{\prime })>\func{cv}^{2})\leq \alpha $. For a finite or parametric set of $ \mathcal{V}$, $\sup_{\mathbf{\Sigma }\in \mathcal{V}}\mathbb{P}_{\mathbf{ \Sigma }}(\tau ^{2}(\mathbf{W}\mathbf{W}^{\prime })>\func{cv}^{2})$ can be established numerically. We therefore focus on an analytical robustness result for a non-parametric class $\mathcal{V}$.
Specifically, we establish a set of readily verifiable sufficient conditions to check robustness for sets $\mathcal{V}$ that are composed of mixtures of parametric covariance matrices $\mathbf{\Sigma }^{p}(\theta )$ for $\theta \in \Theta $. We then apply this result to a set of Mat\'{e}rn covariance matrices with parameter $\theta $ and investigate the robustness of SCPC over arbitrary mixtures of these Mat\'{e}rn models. In addition, we use the result to study the robustness of a popular projection based t-test in a regularly spaced time series setting.
Consider a benchmark model with $\mathbf{\Sigma }=\mathbf{\Sigma }_{0}$, and suppose that $\func{cv}$ has been chosen so that $\mathbb{P}_{\mathbf{\Sigma }_{0}}(\tau ^{2}(\mathbf{W}\mathbf{W}^{\prime })>\func{cv}^{2})=\alpha .$ We are interested in conditions under which
for a probability distribution $F$.
Let $\lambda_{j}(\cdot)$ denote the $j$th largest eigenvalue of some matrix.
The critical value for the SCPC t-test is chosen to control size in exponential models with $c\geq c_{0}$, where $c_{0}$ is calibrated to a value $\overline{\rho}_{0}$. Because $\overline{\rho}$ is monotone in $c$, the resulting SCPC t-test controls size for all $\overline{\rho}\leq \overline{\rho}_{0}$ in the exponential model by construction.
Let $\mathbf{\Sigma}^{p}(\theta)$ denote the covariance matrix associated with a parameter $\theta$, with average pairwise correlation $\overline{\rho} (\theta)$. Let $\Theta_{\overline{\rho}_{L},\overline{\rho}_{U}}=\{$$\theta| \overline{\rho}_{L}\leq\overline{\rho}(\theta)\leq\overline{\rho}_{U}\}$ denote the set of values of $\theta$ that induce correlations between $ \overline{\rho}_{L}$ and $\overline{\rho}_{U}$. If the inequalities in Theorem (ref) are satisfied for all values of $\theta\in\Theta_{ \overline{\rho}_{L},\overline{\rho}_{U}}$, then the SCPC t-test controls size for all mixtures of $\mathbf{\Sigma}^{p}(\theta)$ in this set.
In this section we consider $\mathbf{\Sigma }^{p}(\theta )$ computed from Mat \'{e}rn processes with parameter $\theta =(\nu ,c)$, where $\nu $ and $c$ are positive constants. If $u$ follows a Mat\'{e}rn process, its covariance function $\sigma _{u}(r-s)$ depends on the locations only through $d=||r-s||$ . For $\nu \in \{1/2,3/2,5/2,\infty \}$, the Mat\'{e}rn covariance functions are
For any $\mathbf{\Sigma }(c_{0})$ it is straightforward to compute the bounds $\overline{\rho }_{L}$ and $\overline{\rho }_{U}$ such that the inequalities in Theorem (ref) are satisfied for all values of $ \theta \in \Theta _{\overline{\rho }_{L},\overline{\rho }_{U}}$ with $\nu \in \{1/2,3/2,5/2,\infty \}$ and $c>0$. We carried out this exercise for the U.S. states spatial correlation designs of Section (ref) (the calculations for one set of locations take less than a second). We find $ \overline{\rho }_{L}\leq 0.001$ and $\overline{\rho }_{U}=\overline{\rho } _{0}\in \{0.02,0.10\}$, with very few minor exceptions.
We conclude that SCPC controls size in finite Gaussian samples for a wide range of Mat\'{e}rn process mixtures that imply $\overline{\rho}\leq \overline{\rho}_{0}$, at least for this set of spatial designs.
The spatial design is fixed for regularly-spaced time series, so the theorem can provide general robustness results. Consider, for instance, the equal weighted cosine (EWC) projection estimator of M\"{u}ller (2004, 2007), \citeasnoun {Lazarus_etal_JBES_2018} and Dou (2019) where $\mathbf{w}(s)=\sqrt{2/q}(\cos \pi s,\cos (2\pi s),\ldots ,\cos (q\pi s))$. Suppose the critical value $ \func{cv}_{n}$ is chosen so that size is controlled in a Gaussian AR(1) with coefficient $\exp (-c_{0}/n)$, and $q$ is chosen to minimize expected length in the i.i.d. model. For $c_{0}=10$, $c_{0}=25$ and $c_{0}=50$, we obtain $ q=5,7$ and $10$, respectively, for all $n\in \{50,100,500\}$. Call this test the EWC$(c_{0})$ t-test.
Calculations based on Theorem (ref) for these values of $c_{0}$ and $n$ show that the EWC$(c_{0})$ t-test controls size for arbitrary mixtures of AR(1) processes with coefficients $\exp (-c/n)$, $c\geq c_{0}$. By taking the limit in $n$ and using standard local-to-unity weak convergence results (as in \citeasnoun{Muller14}), one can further apply Theorem 1 to the limiting covariance matrices $\mathbf{\Omega }_{0}$ and $\mathbf{ \Omega }(\theta )$ to study asymptotic robustness of the EWC$(c_{0})$ t-test with an asymptotically justified critical value (which are equal to $\func{cv }=3.53$, $2.71$, $2.40$ for $c_{0}=10$, $25$, $50$, respectively). Another numerical calculation based on Theorem (ref) then shows that these EWC$(c_{0})$ t-tests control asymptotic size for underlying processes that are arbitrary mixtures of local-to-unity models with parameters $c\geq c_{0}$.
Moreover, let $f_{n,0}:[-\pi ,\pi ]\mapsto \lbrack 0,\infty )$ be the spectral density of an AR(1) process with coefficient $\exp (-c_{0}/n)$, so $ f_{n,0}(\omega )\propto (1-2e^{-c_{0}/n}\cos \omega +e^{-2c_{0}/n})^{-1}$. A spectral density $f_{n,1}$ would naturally be considered less persistent than $f_{n,0}$ if $f_{n,1}(\omega )/f_{n,0}(\omega )$ is (weakly) monotonically increasing in $|\omega |$. Denote all such functions by $ \mathcal{F}_{n}$. Define
so $M$ measures by how much $f_{n,1}(\omega )/f_{n,0}(\omega )$ increases over $[0,\pi ]$, and denote by $\mathcal{F}_{n}^{\bar{M}}$ all functions in $ \mathcal{F}_{n}$ with $M\leq \bar{M}$ for some $\bar{M}>1$. Then for any $ f_{n,1}\in \mathcal{F}_{n}^{\bar{M}}$, there exists a CDF $H$ on $[0,\pi ]$ such that
so $f_{n,1}$ has a representation as a scale mixture of $f_{n,0}(\omega )+( \bar{M}-1)\mathbf{1}[|\omega |\geq \theta ]f_{n,0}(\omega )$, $0\leq \theta \leq \pi $. After translating this back into a corresponding mixture of covariance matrices $\mathbf{\Sigma }^{p}(\theta )$, an application of Theorem (ref) shows that the EWC$(c_{0})$ t-test also controls size in this class, for $(c_{0},\bar{M})\in \{(10,10),(25,10),(50,5)\}$ and all $n\in \{50,100,500\}$. These results refine corresponding results in \citeasnoun{Dou_2019} that are based on a Whittle-type diagonal approximation to $ \mathbf{\Sigma }$.
Taking limits as $n\rightarrow \infty $ yields a corresponding asymptotic robustness statement: The function $f_{0}: \mathbb{R} \mapsto \lbrack 0,\infty )$ with $f_{0}(\omega )=(\omega ^{2}+c_{0}^{2})^{-1} $ is proportional to the `local-to-zero' spectral density (cf. M\"{u}ller and Watson (2016, 2017))\nocite{Muller14b}\nocite {Muller15d} of a local-to-unity process with parameter $c_{0}.$\ Consider any process whose local-to-zero spectral density $f_{1}$ is such that $ f_{1}(\omega )/f_{0}(\omega )$ is monotonically increasing in $|\omega |$ with $\lim_{\omega \rightarrow \infty }f_{1}(\omega )/f_{0}(\omega )\leq \bar{M}f_{1}(0)/f_{0}(0)$ and that satisfies the CLT in M\"{u}ller and Watson (2016, 2017). A numerical calculation based on Theorem (ref) then shows that the EWC$(c_{0})$ t-tests for $(c_{0},\bar{M} )\in \{(10,10),(25,10),(50,5)\}$ controls size in large samples under all such processes.
The SCPC t-test is not robust to heteroskedasticity or measurement error in locations by construction. For example, suppose that $u(s)=h(s)\tilde{u}(s)$ , where $\tilde{u}$ is homoskedastic and satisfies the assumptions outlined above for $u$, and $h:\mathcal{S}\mapsto \mathbb{R} $ is a non-random function that induces heteroskedasticity in the $u$ process. The linear combinations of $u$ studied in Lemma (ref) are now $\sum_{l=1}^{n}\mathbf{w}^{0}(s_{l})u(s_{l})=$$ \sum_{l=1}^{n}\mathbf{w}_{h}^{0}(s_{l})\tilde{u}(s_{l})$ where $\mathbf{w} _{h}^{0}(s)=\mathbf{w}^{0}(s)h(s)$. The results of the lemma and subsequent theorems then follow with $\mathbf{w}_{h}^{0}$ replacing $\mathbf{w}^{0}$. But, the test statistic and critical value is computed using $\mathbf{w}^{0}$ , not $\mathbf{w}_{h}^{0}$, so that size control is not guaranteed, even in large samples. An analogous problem arises when the locations $s_{i}$ are measured with error.
In both cases, the particulars of the size distortion depend on the distribution of spatial locations, $g$, the weights $\mathbf{w}^{0}$ (which in turn depend on the value of $\overline{\rho}_{0}$ used to calibrate $ c_{0} $), the function $h$ in the heteroskedastic model and the distribution of the measurement error for the locations.
We summarize two experiments that illustrate and quantify the size distortions in the U.S. states spatial correlation designs. The first experiment is a heteroskedastic model with $\log h$ increasing or decreasing linearly from $\log h(s)=0$ to $\log h(s)=\log 3$ moving from the most westward to the most eastward location, the experiment is repeated with $h$ increasing or decreasing moving north to south, and we record the largest of the four rejection frequencies. Panel (a) of Figure (ref) plots the CDF of rejection frequencies for nominal 5% SCPC tests for each $( \overline{\rho }_{0},g)$ pair. For these designs, the resulting size distortions are not large, except for a few states with $\overline{\rho } _{0}=0.02$ and the light spatial density $g$, where rejection frequencies approach 10%.
The second experiment investigates location measurement error of a form studied in \citeasnoun{Conley_Molinari_2007}. Specifically for each location, $ s_{i}^{\ast }=s_{i}+e_{i}$ where $s_{i}^{\ast }$ is the measured location, $ s_{i}$ is the true location and $e_{i}$ is the measurement error. The error term is $e_{i}=(e_{1,i},e_{2,i})$ with $e_{1,i}$ the north-south and $ e_{2,i} $ the east-west coordinate and $e_{j,i}$ i.i.d.$\mathcal{U}(-\delta ,\delta ) $\ over $j$ and $i$, and $\text{$\delta =0.0375H$}$ with $H$ the length of the smallest square that encompasses all locations, corresponding to \textquotedblleft level 4\textquotedblright\ errors in Conley and Molinari's (2007) classification. The CDFs for the rejection frequencies are shown in panel (b) of Figure (ref). Evidently, measurement error of this sort has little effect on the size of SCPC under uniformly distributed locations, but can have a substantial effect for highly concentrated spatial distributions, especially when $\overline{\rho } _{0}=0.02$.
Figure (ref) showed the expected length of the SCPC confidence interval relative to the length of an oracle confidence interval that uses the true value of $\func{Var}(\sqrt{n}( \overline{y}-\mu))$ conditional on the observed locations $\mathbf{s}$. (As before, in this subsection we keep the conditioning on $\mathbf{s}$ and the dependence on $n$ implicit.) For studying efficiency, a more relevant comparison involves the expected length of the SCPC confidence interval relative to a confidence interval that, like SCPC, does not depend on the true (unknown) value of $\func{Var}(\sqrt{n}(\overline{y}-\mu))$. Ideally, such a comparison would involve SCPC and the most efficient method for constructing a confidence interval. We undertake such a comparison here.
To be specific, let $\func{CS}(\mathbf{y})\subset \mathbb{R} $ denote a confidence set for $\mu$ constructed from $\mathbf{y}$. We restrict attention to location and scale equivariant confidence sets, that is $\func{CS}$ satisfies $\func{CS}(a_{\mu}+a_{\sigma}\mathbf{y} )=\{\mu_{0}:(\mu_{0}-a_{\mu})/a_{\sigma}\in\func{CS}(\mathbf{y})\}$ for all $ \mathbf{y}$, $a_{\mu}\in \mathbb{R} $ and $a_{\sigma}>0$. As in Section (ref), we focus on the Gaussian model $\mathbf{y}\sim\mathcal{N}(\mathbf{l}\mu, \mathbf{\Sigma})$. We want to compare the SCPC interval with a confidence interval that, like SCPC, has good coverage $\mathbb{P}_{\mathbf{\Sigma} }(\mu\in\func{CS}(\mathbf{y}))$ over a range of potential spatial correlation patterns $\mathbf{\Sigma\in}\mathcal{V}$. The metric for measuring efficiency is the expected length $\mathbb{E}^{1}[\int\mathbf{1} [x\in\func{CS}(\mathbf{y})]dx]$ in the i.i.d. model $\mathbf{y}\sim\mathcal{N }(\mathbf{l}\mu,\mathbf{I})$.
Our choice of $\mathcal{V}$ is motivated by the structure of the SCPC benchmark covariance matrix $\mathbf{\Sigma}(c_{0})$. The idea is to include in $\mathcal{V}$ covariance matrices that are weakly less persistent than $ \mathbf{\Sigma}(c_{0})$, and that cannot be easily distinguished from the i.i.d. model. To characterize these covariance matrices, note that $\mathbf{ \Sigma}(c_{0})$ is generated from $u$, an isotropic random field with covariance function $\sigma_{u}(s,r)=\exp(-c_{0}||s-r||)$. Isotropy implies that the spectrum of this random field $F_{0}: \mathbb{R} ^{d}\mapsto\lbrack0,\infty)$ at frequency $\mathbf{\omega}\in \mathbb{R} ^{d}$ can be written as function of the scalar $\omega=||\mathbf{\omega}||$, that is $F_{0}(\mathbf{\omega})=f_{0}(\omega)$ for some $f_{0}: \mathbb{R} \mapsto\lbrack0,\infty)$. As is well known, the exponential covariance model for $d=2$ corresponds to a spectral density function $f_{0}$ proportional to $(c_{0}+\omega^{2})^{-3/2}$. By scale invariance of both $\func{CS}$ and the SCPC interval, it is without loss of generality to set $f_{0}$ equal to
For some $\bar{\omega}>0$, define $f_{\Delta}(\omega)=\mathbf{1}[|\omega|\leq \bar{\omega}](f_{0}(\omega)-f_{\Delta}(\bar{\omega}))$, and let $ f_{R}(\omega)=f_{0}(\omega)-f_{\Delta}(\omega)$, so that
For $0\leq|\omega|\leq\bar{\omega}$, the density $f_{\Delta}$ is equal to $ f_{0}(\omega)-f_{0}(\bar{\omega})$, so that the remainder $f_{R}(\omega)$ is a continuous density that is flat for $|\omega|\leq\bar{\omega}$, and that follows the same decline as $f_{0}$ for $|\omega|>\bar{\omega}$. Since both $ f_{\Delta}(\omega)$ and $f_{R}(\omega)$ are non-negative, we have the corresponding identity in covariance matrices
where $\mathbf{\Sigma}_{\Delta}(\bar{\omega})$ and $\mathbf{\Sigma}_{R}(\bar{ \omega})$ are induced by the isotropic random fields with spectral densities $F_{\Delta}(\mathbf{\omega})=f_{\Delta}(||\mathbf{\omega||)}$ and $F_{R}( \mathbf{\omega})=f_{R}(||\mathbf{\omega||)}$, respectively.
Now consider the covariance matrix
where $\lambda _{1}(\mathbf{\Sigma }_{R}(\bar{\omega}))$ is the largest eigenvalue of $\mathbf{\Sigma }_{R}(\bar{\omega})$. Since $f_{R}(\omega )$ is monotonically decreasing in $|\omega |$, also $\mathbf{\Sigma }_{R}(\bar{ \omega})$ contributes to the persistence of $\mathbf{\Sigma }(c_{0})$ in ( (ref)), so replacing it with white noise of weakly larger variance should make inference about $\mu $ under $\mathbf{\bar{\Sigma}}( \bar{\omega})$ no harder than under $\mathbf{\Sigma }(c_{0})$.\footnote{ In the regularly-spaced time series setting, white noise amounts to a flat spectrum, so $\mathbf{\Sigma }_{0}(\bar{\omega})$ corresponds to an underlying spectral density equal to $f_{\Delta }(\omega )+f_{0}(\bar{\omega} )$, which is the \textquotedblleft kinked\textquotedblright\ spectral density considered by \citeasnoun{Dou_2019}. For arbitrary locations, however, the domain of the spectrum doesn't fold onto the interval $[-\pi ,\pi ]$, so that white noise cannot mathematically be represented by a flat spectrum.} Said differently, a method that is robust under correlation patterns weakly less persistent than $\mathbf{\Sigma }(c_{0})$ should continue to have good coverage after replacing medium and high frequency variation in $\mathbf{y}$ by white noise, that is, under $\mathbf{\bar{\Sigma}}(\bar{\omega})$. This motivates the set $\mathcal{V}=\{\mathbf{\bar{\Sigma}}(\bar{\omega})|\bar{ \omega}>0\}$.
A calculation shows that in the U.S. states spatial correlation designs, the SCPC interval has good coverage properties under this $\mathcal{V}$. With $ \alpha_{\text{SCPC}}(\bar{\omega})=\mathbb{P}_{\mathbf{\bar{\Sigma}}(\bar{ \omega})}(\tau_{\text{SCPC}}^{2}>\func{cv}_{\text{SCPC}}^{2})$ for the nominal 5% level SCPC test, for most designs, $\sup_{\bar{\omega} \geq0}\alpha_{\text{SCPC}}(\bar{\omega})$ is equal or very close to 5%, and it never exceeds 8%. To keep things on an equal footing, we allow $\func{CS} $ the same degree of undercoverage, that is we consider the problem
In words, we seek the invariant confidence set with the shortest expected length in the i.i.d. location model among all confidence sets that are as robust as the SCPC interval under $\mathbf{\bar{\Sigma}}(\bar{\omega})$, $ \bar{\omega}>0$.
Since $\bar{\omega}$ is one-dimensional, one can apply the numerical techniques of \citeasnoun{Elliott15} and \citeasnoun{Muller15c} (also see \citeasnoun {Mueller20}) to obtain an informative lower bound on the objective $\inf_{ \func{CS}}\mathbb{E}^{1}[\int\mathbf{1}[x\in\func{CS}(\mathbf{y})]dx]$ that holds for any equivariant $\func{CS}(\mathbf{y})$ that satisfies the constraint in ((ref)).
We compute such lower bounds in the U.S. states spatial correlation designs. Panel (a) of Figure (ref) shows the CDFs of the length of SCPC confidence intervals relative to the lower bounds for the 240 designs in each $(\overline{\rho }_{0},g)$ pair. The expected lengths of SCPC are within 7% of the efficiency bound for all designs when $\overline{\rho }_{0}=0.02$. When $\overline{\rho }_{0}=0.10$, so that spatial correlation is high, and the spatial locations are highly concentrated as under the light design, the expected length of the SCPC confidence interval can be more that 15% longer than the efficiency bound. In part, this is because the implied efficient confidence sets are complicated and rather uninterpretable functions of $\mathbf{y}$ in this case. We thus repeat the exercise for confidence sets constrained to be symmetric around $\overline{y}$ by imposing $\func{CS}(a_{\mu }+a_{\sigma } \mathbf{y})=\{\mu _{0}:(\mu _{0}-a_{\mu })/a_{\sigma }\in \func{CS}(\mathbf{y })\}$ for all $\mathbf{y}$, $a_{\mu }\in \mathbb{R} $ and $a_{\sigma }\neq 0$. The results are summarized in panel (b), and we can see that SCPC comes closer to the resulting higher bound on confidence interval length.
This section compares SCPC with other methods that have been proposed, focusing on size and expected length of confidence intervals in the benchmark Gaussian model with exponential covariance kernel and parameter $ c_{0}$ (calibrated by $\bar{\rho}_{0}$). We consider two kernel-based methods, two versions of a cluster method, and one projection method. All these methods are t-statistic based tests of the form considered in Section 3.
The kernel based methods use a Bartlett kernel, $k(s,r)=k_{\text{Bartlett} }(||s-r||/b)$. The methods differ in their choice of bandwidth $b$ and critical value. The first method uses a standard normal critical value with $ b$ chosen so the resulting test has size as close as possible to $5\%$. This is a version of the method proposed by \citeasnoun{Conley99}, but with an oracle choice for the bandwidth. The second method sets $b=\max_{l,\ell}||s_{l}-s_{ \ell}||$ and chooses the critical value to obtain exact coverage under $ \mathbf{\Sigma}=\mathbf{I}$. This is the spatial analogue of the method suggested by \citeasnoun{Kiefer00} (KVB) for regularly spaced time series. The cluster methods follow the approach of \citeasnoun{Ibragimov10} (IM) with student-t $_{q}$ critical values and is implemented with $q=4$ and $q=9$ equal-sized clusters.\footnote{ The assignment of locations to clusters is performed sequentially, where at each step, we minimize (across yet unassigned locations) the maximal distance over clusters (among those that have not yet been assigned $n/q$ locations). Cluster distances are computed from the northwest, northeast, southeast and southwest corners of the location circumscribing rectangle, and in the $q=9$ case, also from the mid-points of the four sides of this rectangle, and its center.} The projection method follows \citeasnoun{Sun_Kim_2012} . It uses a student-t$_{q}$ critical value and $q$ low-frequency Fourier weights orthogonalized using the sample locations, where $q$ is chosen as a function of the exponential model parameter $c_{0}$ using the formula in their equation (8). The first and last method are thus tailored to the true value $c_{0}$, just like SCPC.
We analyze these methods in the U.S. states spatial correlation designs, augmented to also include the value $\overline{\rho }_{0}=0.001$ for the average pairwise correlation to investigate performance under `weak' spatial correlations. Figure (ref) summarizes the results for size control and expected lengths by plotting the CDFs for each $(\overline{\rho }_{0},g)$ pair. The first column shows the null rejection frequency for each method; by construction, the rejection frequency for SCPC is at most $5\%$ in all designs. The expected lengths in the second and third column use size-corrected critical values to ensure 95% coverage under $\mathbf{\Sigma }(c_{0})$, and are given in multiples of the expected length of the (non-adjusted) SCPC method. The second column reports these relative expected lengths under $\mathbf{\Sigma } =\mathbf{I}$, and the third column under $\mathbf{\Sigma }(c_{0})$.
Looking at the first column, the kernel and cluster methods have null rejection probabilities close to $5\%$ when $\overline{\rho }_{0}=0.001$, but exhibit significant size distortions for $\overline{\rho }_{0}=0.02$ or $ 0.10$. Evidently, the kernel and cluster methods substantially underestimate the variance of $\overline{y}$ for the latter two values of $\overline{\rho } _{0}$. In contrast, the Fourier projection method has relatively small size distortions under $g=g_{\text{uniform}}$ but can have substantial size distortions under $g=g_{\text{light}}$, even when $\bar{\rho}_{0}=0.001$. This is consistent with the implications of Theorem (ref): the student-t critical value for the projection method is appropriate when $ \mathbf{\Omega }\propto \mathbf{I}$, which it is under weak-correlation with $g$ uniform, but not otherwise, even for large $q$ (cf. Remark (ref)).
The relative lengths shown in the second column are above unity, sometimes by a wide margin, indicating that SCPC is closer to the efficiency bound computed in Section (ref) than these alternative methods, at least for the designs considered here. The third column shows that this continues hold for lengths computed under $\mathbf{\Sigma }(c_{0})$ with a few exceptions. Notably, the expected length of the size-adjusted 9-cluster method is smaller than SCPC when $\overline{\rho }_{0}=0.10$. This apparent good performance comes at the cost of substantially longer confidence intervals in the i.i.d. model.
This section discusses extensions of the method to regression and GMM models, some computational issues, and the multivariate extension of SCPC.
The extension of these results to regression and GMM problems follows from standard arguments. For example, consider the linear regression problem
where $\beta$ is the (scalar) parameter of interest, $\mathbf{z}_{l}$ are additional controls in the regression, and $(w_{l},x_{l},\mathbf{z}_{l})$ are associated with location $s_{l}$. Let $\tilde{x}_{l}=x_{l}-\mathbf{S} _{xz}\mathbf{S}_{zz}^{-1}\mathbf{z}_{l}$ denote the residual from regressing $x_{l}$ on $\mathbf{z}_{l}$, where we use the notation $\mathbf{S} _{ab}=n^{-1}\sum_{l=1}^{n}\mathbf{a}_{l}\mathbf{b}_{l}^{\prime}$ for any vectors $\mathbf{a}_{l}$ and $\mathbf{b}_{l}$. Suppose $\mathbf{S}_{\tilde{x} \tilde{x}}\overset{p}{\rightarrow}\sigma_{\tilde{x}\tilde{x}}^{2}>0$ and
Then
where $\sigma^{2}=\sigma_{\tilde{x}\varepsilon}^{2}/\sigma_{\tilde{x}\tilde{x }}^{4}$. Spatial correlation affects inference in this model through $ \sigma_{\tilde{x}\varepsilon}^{2}$ which incorporates potential correlation between $\tilde{x}_{l}\varepsilon_{l}$ and $\tilde{x}_{\ell}\varepsilon_{ \ell}$ at spatial locations $s_{l}$ and $s_{\ell}$.
Thus, suppose that $\tilde{x}_{l}\varepsilon _{l}$ satisfies the assumptions previously made for $u_{l}$. Then a straightforward calculation shows that setting
in the analysis of the previous sections leads to analogous results with $ \beta $ replacing $\mu $ as the parameter of interest. The extension to GMM inference is analogous; see, for instance, Section 4.4 of \citeasnoun{Mueller2020}.
We highlight two computational issues. The first involves the calculation of the SCPC critical value, and the second involves the problem of computing the eigenvectors $\mathbf{r}_{j}\ $of $\mathbf{M\Sigma }(c_{0})\mathbf{M}$ when $n$ is very large.
The critical value $\func{cv}=\func{cv}_{\text{SCPC}}(q)$ solves $ \sup_{c\geq c_{0}}\mathbb{P}_{\mathbf{\Sigma }(c)}(\tau ^{2}(q^{-1}\sum_{j=1}^{q}\mathbf{r}_{j}\mathbf{r}_{j}^{\prime })>\func{cv} ^{2})=\alpha $ or equivalently (from Theorem (ref)) $ \sup_{c\geq c_{0}}\mathbb{P}\left( Z_{0}^{2}>\sum_{i=1}^{q}\eta _{i}Z_{i}^{2}\right) =\alpha $ where $\eta _{i}=-\omega _{i}/\omega _{0}$, $ \omega _{i}$ are the eigenvalues of $\mathbf{\hat{W}}^{0\prime }\mathbf{ \Sigma }(c)\mathbf{\hat{W}}^{0}\mathbf{D}(\func{cv})$ with $\mathbf{\hat{W}} ^{0}=[\mathbf{l},\mathbf{r}_{1}/\sqrt{q},\ldots ,\mathbf{r}_{q}/\sqrt{q}]$ and $Z_{j}\sim $i.i.d. $\mathcal{N}(0,1)$. \citeasnoun{Bakirov05} show that
which is readily evaluated by numerical quadrature. Thus $\func{cv}_{\text{ SCPC}}(q)$ can be obtained by combining a root-finder with a grid search over $c\geq c_{0}$.
The second problem involves computing the eigenvectors $\mathbf{r} _{j}=(r_{j,1},\ldots ,r_{j,n})^{\prime }$ of the $n\times n$ matrix $\mathbf{ M\Sigma }(c_{0})\mathbf{M}$ when $n$ is very large (say, larger than $n=2000$ ). Here we can leverage the eigenfunction convergence result in Lemma (ref) as discussed in Section (ref): In the notation defined there, we seek to approximate $\mathbf{r}_{j}=(\hat{\varphi} _{j}^{0}(s_{1}),\ldots ,\hat{\varphi}_{j}^{0}(s_{n}))^{\prime }$. Consider a random subset of size $\tilde{n}<n$ of the observed locations $\{\tilde{s} _{l}\}_{l=1}^{\tilde{n}}\subset \{s_{l}\}_{l=1}^{n}$, and let $\mathbf{ \tilde{\Sigma}}(c_{0})$ be the implied $\tilde{n}\times \tilde{n}$ covariance matrix of $(u(\tilde{s}_{1}),\ldots ,u(\tilde{s}_{n}))^{\prime }$ using the benchmark covariance function $\sigma _{u}^{0}(r-s|c_{0})=\exp [-c_{0}||r-s||]$. Let the eigenvector corresponding to the $j$th largest eigenvalue $\tilde{\lambda}_{j}$ of $\mathbf{\tilde{\Sigma}}(c_{0})$ be $ \mathbf{\tilde{r}}_{j}=(\tilde{r}_{1,j},\ldots ,\tilde{r}_{\tilde{n} ,j})^{\prime }$ with $\tilde{n}^{-1}\mathbf{\tilde{r}}_{j}^{\prime }\mathbf{ \tilde{r}}_{j}=1$. As long as $\tilde{n}\rightarrow \infty $ and $\lambda _{q+1}>\lambda _{q}$, Lemma (ref) implies that the span of the $\mathcal{S}\mapsto \mathbb{R} $ functions
converges to the eigenspace spanned by $\varphi _{j}^{0}$, $j=1,\ldots ,q$, just like the full sample estimators $\hat{\varphi}_{j}^{0}$. Thus, it is formally justified to approximate the value of $\hat{\varphi}_{j}^{0}$ at locations $\{s_{l}\}_{l=1}^{n}\ni s_{\ell }\notin \{\tilde{s}_{l}\}_{l=1}^{ \tilde{n}}$ via $r_{j,\ell }=\hat{\varphi}_{j}^{0}(s_{\ell })\approx \tilde{ \varphi}_{j}^{0}(s_{\ell })$---this is a version of the so-called Nystr\"{o} m method (see, for instance, \citeasnoun{Rasmussen05} for discussion and references).
In practice, such approximations can be carried out for several random subsets of $\tilde{n}$ locations, followed by a (sample) principle component analysis to extract the best approximation to the space spanned by the first $q$ eigenvectors. The resulting algorithm has $O(n)$ running time (in contrast to the $O(n^{2})$ running time of a basic implementation of \citeasnoun {Conley99}-type kernel estimators). We provide corresponding STATA and Matlab code in the replication files.
Consider the case where $\mathbf{y}_{l}=\mathbf{\mu }+\mathbf{u}_{l}$ with $ \mathbf{y}_{l}$, $\mathbf{\mu }$ and $\mathbf{u}_{l}$ $m\times 1$ vectors, and we seek to test the hypothesis $H_{0}:\mathbf{\mu }=\mathbf{\mu }_{0}$. Suppose the observations conditional on $\mathbf{s}$ are generated by the model
where $\mathbf{B}(s)$ is an $ \mathbb{R} ^{m}$-valued mean-zero stationary random field on $ \mathbb{R} ^{d}$ with covariance function $\mathbb{E}[\mathbf{B}(s)\mathbf{B} (r)^{\prime }]=\mathbf{\sigma }_{B}(r-s)$. Let $\mathbf{Y}$ and $\mathbf{U}$ be the $n\times m$ matrices of observations and innovations, respectively, and\ $\mathbf{\bar{y}}=n^{-1}\sum_{l=1}^{n}\mathbf{y}_{l}$ the sample mean. The natural analogue to the t-statistic $\tau ^{2}(\mathbf{\hat{ W}\hat{W}}^{\prime })$ is Hotelling's-$T^{2}$ statistic
One would expect that under mixing and moment conditions similar to those of Lemma (ref) (ii)
Note that $T^{2}(\mathbf{\hat{W}\hat{W}}^{\prime })$ is invariant to the transformation $\mathbf{Y\rightarrow YH}$ for nonsingular $\mathbf{H}$. For the purposes of studying the limit distribution of $T^{2}(q)$ under weak correlation, it is thus without loss of generality to normalize $\mathbf{ \sigma }_{B}(\cdot )$ such that the limit covariance matrix in ((ref))\ becomes
where $\mathbf{\kappa }$ is a $m\times 1$ vector with elements in $[0,1)$.
For the extension of the SCPC method, consider a benchmark model indexed by $ \mathbf{c}=(c_{1},\ldots,c_{m})$ where $\limfunc{vec}(\mathbf{Y)|s}\sim \mathcal{N}(\mathbf{\mu}\otimes\mathbf{l}_{n},\mathbf{\Sigma(c))}$ with $ \mathbf{\Sigma(c})=\limfunc{diag}(\mathbf{\Sigma}(c_{1}),\ldots,\mathbf{ \Sigma}(c_{m}))$, and $\mathbf{\Sigma}(c)$ is as in Section (ref). Let $\mathbf{c}_{0}=c_{0}\mathbf{l}_{m}$, a $m\times1$ vector of identical elements $c_{0}$. The SCPC test statistic $T_{\text{SCPC}}^{2}(q)$ is a special case of ((ref))\ with the columns of $\hat{\mathbf{W}} $ equal to the first $q$ eigenvectors of $\mathbf{\Sigma}(c_{0})$, scaled to have length $1/\sqrt{q}$, and with critical value $\func{cv}_{\text{SCPC} }^{T}$ chosen to satisfy
under the null hypothesis, where $\mathbf{c}\geq\mathbf{c}_{0}$ is understood as an elementwise inequality. The value of $q$ that minimizes the expected volume of the confidence ellipsoid under $\limfunc{vec}(\mathbf{Y)|s }\sim\mathcal{N}(\mathbf{\mu}\otimes\mathbf{l},\mathbf{I}_{m}\mathbf{\otimes I}_{n}\mathbf{)}$ is
where $\mathbf{S}_{q}$ is distributed Wishart with $q$ degrees of freedom, and the equality follows from Bartlett's decomposition of a Wishart random matrix, and the formulas for the expectation of a $\chi$ random variable and the volume of an $m$ dimensional ellipsoid.
Since appropriate choices of $c_{j,n}\rightarrow \infty $, $j=1,\ldots ,m$ in the benchmark model can replicate the normalized limit distributions ((ref)) for all $\mathbf{\kappa }$, by the same arguments that lead to Theorem (ref), $T_{\text{SCPC}}^{2}(q)$ controls size under all weak correlation patterns that induce ((ref)). And as in Section (ref), it is straightforward to adapt $T_{\text{SCPC} }^{2}(q)$ to test $m$ restrictions in linear regression and GMM problems. We omit details for brevity. Generalizing the results about the small sample robustness of $\tau _{\text{SCPC}}$ under potentially strong correlations in Theorem (ref) to $T_{\text{SCPC}}^{2}$ is interesting but challenging, and beyond the scope of this paper.