# General Conditions for Valid Inference in Multi-Way Clustering

Luther Yap Department of Economics, Princeton University. Email:
<EMAIL_ADDRESS>
###### Abstract

This paper proves a new central limit theorem for a sample that exhibits
multi-way dependence and heterogeneity across clusters. Statistical inference
for situations where there is both multi-way dependence and cluster
heterogeneity has thus far been an open issue. Existing theory for multi-way
clustering inference requires identical distributions across clusters (implied
by the so-called separate exchangeability assumption). Yet no such homogeneity
requirement is needed in the existing theory for one-way clustering. The new
result therefore theoretically justifies the view that multi-way clustering is
a more robust version of one-way clustering, consistent with applied practice.
The result is applied to linear regression, where it is shown that a standard
plug-in variance estimator is valid for inference.

##  1 Introduction

Clustering standard errors on multiple dimensions is common and attractive in
applied econometrics because it allows observations to be dependent whenever
they share a cluster on any dimension.111E.g., Dube et al. (2010) clustered on
state and border segment when studying the effect of minimum wages on
employment; Nunn and Wantchekon (2011) clustered on ethnic groups and district
when studying the effect of slave trade on trust; Michalopoulos and
Papaioannou (2013) clustered on country and ethnolinguistic family when
studying the effect of pre-colonial institutions on development. The variance
estimator proposed by Cameron et al. (2011) (henceforth CGM) has thus been
widely applied to contexts with multi-way dependence. Existing justification
for the asymptotic validity of the CGM estimator and other inference
procedures in multi-way clustering relies on separate exchangeability, which
implies the homogeneity of clusters. This paper provides general conditions
such that the plug-in mean estimator is asymptotically normal, and the CGM
variance estimator is consistent, even when clusters are heterogeneous. These
conditions do not include separate exchangeability, and they mimic the
conditions in one-way clustering: the only substantive assumption is that two
observations are independent when they do not share any cluster. Since
asymptotic normality and consistent variance estimation are sufficient for
valid inference, the results in this paper provide sufficient general
conditions for valid inference in multi-way clustering.

An environment with multi-way clustering permits dependence whenever
observations share at least one cluster. To fix ideas, suppose observations
can be partitioned on two different dimensions — state and industry.
Observations in the same state or in the same industry are plausibly
correlated, but two observations in different states and different industries
are assumed to be independent.222This setting permits more general dependence
structures than one-way clustering. If there is one-way clustering by state,
then two observations from different states are automatically independent. In
two-way clustering, two observations from different states are not necessarily
independent because they may share the same industry. The CGM variance
estimator accommodates such dependence, and subsequent literature provided a
theoretical basis for its validity (e.g., Davezies et al. (2021); MacKinnon et
al. (2021)). Menzel (2021) also showed the validity of a bootstrap procedure
for multi-way clustering that is robust to asymptotic non-
normalities.333Menzel (2021) pointed out that a purely interactive data-
generating processes unique to multi-way dependence has an asymptotic
distribution that is not normal. Section 2 will consider this process and show
how the assumptions of this paper rules it out. The theoretical basis for
inference thus far relies on separate exchangeability, the assumption that
random variables are exchangeable on either clustering dimension, though not
necessarily both.

However, as noted by MacKinnon et al. (2021), separate exchangeability implies
identical marginal distributions. Since exchangeability implies identical
distribution, separate exchangeability in the state-industry example implies
that random variable in Alaska and California must be drawn from the same
distribution. In contrast, existing asymptotic theory on one-way clustering
(e.g., Hansen and Lee (2019); Djogbenou et al. (2019)) allows the distribution
of the random variable to be heterogeneous over clusters. The only substantive
assumption is that observations that do not share any cluster are independent.
Since the only available conditions for the validity of multi-way clustering
require separate exchangeability, the literature lacks general conditions for
multi-way clustering that generalize one-way clustering and permit
heterogeneity over clusters. This paper fills the gap, and thus justifies
multi-way clustering as a more robust version of one-way clustering.

######  Example 1.

To illustrate separate exchangeability, consider an additive random effects
model. Individual $i$ who belongs to cluster $g(i)$ on the $G$ dimension and
cluster $h(i)$ on the $H$ dimension has random variable $W_{i}$ generated from
$W_{i}=\alpha_{g(i)}+\gamma_{h(i)}+\varepsilon_{i}$, where cluster-specific
$\alpha_{g},\gamma_{h}$ and individual-specific $\varepsilon_{i}$ are
independent of each other. If we assume separate exchangeability, then
$\alpha_{g}$, $\gamma_{h}$, and $\varepsilon_{i}$ are iid.444To see this, for
individuals $i$ and $j$ where $g(i)\neq g(j)$, $h(i)=h(j)=h$, separate
exchangeability implies
$\alpha_{g(i)}+\gamma_{h}+\varepsilon_{i}\stackrel{{\scriptstyle
d}}{{=}}\alpha_{g(j)}+\gamma_{h}+\varepsilon_{j}$. Since
$\alpha_{g},\gamma_{h}$ and $\varepsilon_{i}$ are independent,
$\varepsilon_{i}\stackrel{{\scriptstyle d}}{{=}}\varepsilon_{j}$ and
$\alpha_{g}\stackrel{{\scriptstyle d}}{{=}}\alpha_{g^{\prime}}$. In contrast,
under one-way cluster asymptotics, the cluster-specific error $\alpha_{g}$ is
allowed to be heteroskedastic. General conditions provided in this paper
permits valid inference even when $\alpha_{g},\gamma_{h},\varepsilon_{i}$ are
heteroskedastic in this model.

The main result is a central limit theorem for multi-way clustering with
heterogeneous cluster sizes and distributions. I apply the theorem to a simple
setting of a linear regression, but it is more broadly applicable to many
other econometric procedures that exhibit a similar clustering structure.

##  2 Setting and Main Result

Consider a setup with two-way clustering on dimensions $G$ and $H$ for random
vectors $\\{W_{i}\\}_{i=1}^{n}$, where
$W_{i}:=(W_{i1},W_{i2},\cdots,W_{iK})^{\prime}\in\mathbb{R}^{K}$ and $i$ is
the unit of observation, for a sequence of populations of size
$n$.555Clustering in more than two dimensions is possible, and derivations are
entirely analogous. For example, $G$ could denote states and $H$ denote
industries. This section establishes a central limit theorem (CLT) for a
weighted sum of the random vector i.e., $\sum_{i}\omega_{i}W_{i}$, where
$\omega_{i}$ are nonstochastic scalar weights, as $n\rightarrow\infty$. For
$C\in\\{G,H\\}$, let $\mathcal{N}^{C}_{c}$ denote the set of observations in
cluster $c$ on dimension $C$ — this partitions the population on the $C$
dimension.

Let $g(i)$ and $h(i)$ denote the cluster that observation $i$ belongs to on
the $G$ and $H$ dimensions respectively. These cluster identities are
nonstochastic and observed. Let $N^{C}_{c}=|\mathcal{N}^{C}_{c}|$ denote the
cluster size for $C\in\\{G,H\\}$ and
$N_{gh}:=|\mathcal{N}^{G}_{g}\cap\mathcal{N}^{H}_{h}|$. These cluster sizes
are allowed to be heterogeneous in a way that will be formalized in the
assumptions below. $W_{i}$ is assumed to be independent of any $W_{j}$ when
$j\notin\mathcal{N}^{G}_{g(i)}\cup\mathcal{N}^{H}_{h(i)}=:\mathcal{N}_{i}$,
i.e., when $i$ and $j$ do not share a cluster on either dimension. Hence,
$\mathcal{N}_{i}$ is the set of observations plausibly dependent with $i$.
This environment is stated as Assumption 1, the main substantive assumption.

######  Assumption 1.

$W_{i}\\!\perp\\!\\!\\!\perp W_{j}$ if $g(i)\neq g(j)$ and $h(i)\neq h(j)$.

Assumption 1 is agnostic about the dependence structure when $W_{i}$ and
$W_{j}$ share at least one cluster. It also allows the data generating process
to be arbitrarily heterogeneous across different clusters, mimicking the
heterogeneity permitted in one-way clustering (e.g., Hansen and Lee (2019)).
Since one-way clustering is a special case of two-way clustering where
everyone is in their own $H$ cluster, the result here generalizes existing
results in one-way clustering. In contrast, existing literature in multi-way
clustering assumes separate exchangeability that additionally imposes
identical distribution over clusters, so they do not immediately generalize
one-way clustering. $\\{W_{i}\\}_{i=1}^{n}$ being separately exchangeable
implies Assumption 1 but the converse is not true.666To illustrate this, let
$N_{gh}=1$ and $W_{gh}$ denote the observation in cluster $g$ and $h$ on the
respective dimensions. Due to Kallenberg (2005), $\\{W_{gh}\\}_{g\geq 1,h\geq
1}$ is separately exchangeable if and only if there exists a representation
$W_{gh}=f(\alpha_{g},\gamma_{h},\varepsilon_{gh})$, where
$(\alpha_{g},\gamma_{h},\varepsilon_{gh})\stackrel{{\scriptstyle
iid}}{{\sim}}U[0,1]$. Then, it is obvious that $W_{gh}\\!\perp\\!\\!\\!\perp
W_{g^{\prime}h^{\prime}}$ for $g\neq g^{\prime}$, $h\neq h^{\prime}$. A
counterexample for the converse is some $W_{gh}=-W_{gh^{\prime}}$. These
random variables are allowed to be perfectly correlated since they share a
cluster under Assumption 1. However, we cannot find a representation $f(.)$,
because that representation implies
$E[W_{gh}|\alpha_{g}]\\!\perp\\!\\!\\!\perp E[W_{gh^{\prime}}|\alpha_{g}]$.

Observations that share a cluster are allowed to be dependent, but they need
not be. Hence, let $A_{ij}:=1[W_{i}\mathchoice{\mathrel{\hbox to0.0pt{\kern
8.33331pt\kern-5.27776pt$\displaystyle\not$\hss}{\\!\perp\\!\\!\\!\perp}}}{\mathrel{\hbox
to0.0pt{\kern
8.33331pt\kern-5.27776pt$\textstyle\not$\hss}{\\!\perp\\!\\!\\!\perp}}}{\mathrel{\hbox
to0.0pt{\kern
1.91675pt\kern-4.11108pt$\scriptstyle\not$\hss}{\\!\perp\\!\\!\\!\perp}}}{\mathrel{\hbox
to0.0pt{\kern
0.41675pt\kern-3.3333pt$\scriptscriptstyle\not$\hss}{\\!\perp\\!\\!\\!\perp}}}W_{j}]$
be a 0-1 indicator for whether $W_{i}$ and $W_{j}$ are actually dependent, so
$A_{ij}=A_{ji}$, and $A_{ii}=1$.777It is insufficient to define the indicator
as $A_{ij}:=1[Cov(W_{i},W_{j})\neq 0]$, since the proof contains third and
fourth moments. For $K=1$, zero covariance between a pair of observations is
insufficient to ensure objects such as $E[W_{i}W_{j}W_{k}]$ and
$E[W_{i}W_{j}W_{k}W_{l}]-E[W_{i}W_{k}]E[W_{j}W_{l}]$ are zero. This notation
allows a particular form of misspecification where the researcher is
conservative and clusters on dimension $G$ when it is not required. Every
observation $W_{i}$ is weighted by nonstochastic scalar $\omega_{i}$. For
positive definite matrix $Q$, let $\lambda_{\operatornamewithlimits{min}}(Q)$
denote the smallest eigenvalue of $Q$. Then, let
$Q_{n}:=Var\left(\sum_{i=1}^{n}\omega_{i}W_{i}\right)$ denote the variance of
the sum and $\lambda_{n}:=\lambda_{\operatornamewithlimits{min}}(Q_{n})$
denote its smallest eigenvalue. For example, when $K=1$ and equal weights are
placed on all observations, $W_{i}$ is a scalar and
$\lambda_{n}=Q_{n}=Var(\sum_{i}\omega_{i}W_{i})$. $K_{0}$ is used throughout
the paper to denote an arbitrary constant.

######  Assumption 2.

For $C\in\\{G,H\\}$, and $k\in\\{1,2,\cdots,K\\}$, there exists $K_{0}<\infty$
such that:

  1. 1.

$E[W_{ik}^{4}]\leq K_{0}$ for all $i$.

  2. 2.

$\frac{1}{\lambda_{n}}\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\omega_{i}|\right)^{2}\rightarrow
0$.

  3. 3.

$\frac{1}{\lambda_{n}}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}A_{ij}|\omega_{i}\omega_{j}|\leq
K_{0}$.

Assumption 2.1 requires the fourth moment to be bounded, which is stronger
than the moment condition in one-way clustering.888See equation (7) of Hansen
and Lee (2019) for the condition in one-way clustering. The proof in one-way
clustering usually verifies a Lindeberg condition because blocks of
observations are independent of each other. With multi-way dependence, we no
longer have independent blocks because each cluster can have observations that
are dependent with observations from a different cluster when these
observations share a cluster on a different dimension. Hence, a different
proof strategy is required. The proof in this paper uses Stein’s method, which
requires stronger moment restrictions, but provides a non-asymptotic bound on
the approximation error — details are in Subsection 2.1.

Assumption 2.2 requires the contribution of the cluster with the largest
weight to be small relative to the total variance. In the special case where
everyone is equally weighted with $\omega_{i}=1$, the condition is simply
$(1/\lambda_{n})\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}\rightarrow 0$.
Intuitively, this condition is required so that the removal of a cluster does
not change the variance substantively. This assumption allows the ratio of any
two cluster sizes to diverge to infinity. It is identical to equation (12) of
Hansen and Lee (2019) when $C=G=H$. Assumption 2.2 also rules out having
components that are perfectly negatively correlated: if the components of the
vector were perfectly negatively correlated, $\lambda_{n}=0$.

Assumption 2.3 is fairly unrestrictive about the convergence rate. To aid
exposition, suppose $\omega_{i}=1$ $\forall i$, $K=1$, and $C$ is taken to be
the clustering dimension that
$\lambda_{n}\asymp\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}A_{ij}$.999To be
clear about the notation, $a\asymp b$ if and only if there exists
$K_{0}<\infty$ such that $a/b,b/a\in[-K_{0},K_{0}]$. Since $E[W_{i}^{2}]$ is
bounded,
$\lambda_{n}\asymp\operatornamewithlimits{max}_{C\in\\{G,H\\}}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}A_{ij}$.
With strong dependence, $A_{ij}=1$ for all $i,j\in\mathcal{N}^{C}_{c}$, so
$\lambda_{n}\asymp\sum_{c}(N^{C}_{c})^{2}$. However, if the researcher were
conservative and clustered on $C$ when the data is indeed iid, then $A_{ij}=1$
if and only if $i=j$, so $\lambda_{n}\asymp n$. Assumption 2.3 has
implications on $\lambda_{n}$, which then determines how strong Assumption 2.2
is. Namely, when $\lambda_{n}\asymp n$, Assumption 2.2 requires
$\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}/n\rightarrow 0$. When
$\lambda_{n}\asymp\sum_{c}(N^{C}_{c})^{2}$, then Assumption 2.2 only requires
$\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}/(\sum_{c^{\prime}}(N^{C}_{c^{\prime}})^{2})\rightarrow
0$. The weaker version of Assumption 2.2 permits balanced panels where the
unit and time dimensions increase at the same rate, while the stronger version
does not.101010To see this, let $M$ denote the number of units and time
periods, so there are $M^{2}$ observations.
$\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}/(\sum_{c^{\prime}}(N^{C}_{c^{\prime}})^{2})=M^{2}/M^{3}=1/M\rightarrow
0$, but $\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}/n=M^{2}/M^{2}=1\neq
o(1)$. The assumption that $(1/\lambda_{n})\sum_{c}(N_{c}^{C})^{2}\leq K_{0}$
matches equation (11) of Hansen and Lee (2019).

######  Remark 1.

Assumption 2.3 rules out the following purely interactive model. As pointed
out by Menzel (2021), this model has an asymptotic distribution that is non-
normal, and there is no analog in one-way clustering. For
$g\in\\{1,\cdots,M\\}$, $h\in\\{1,\cdots,M\\}$ and $N_{gh}=1$, we observe
$W_{gh}=\alpha_{g}\gamma_{h}$, where $\alpha_{g},\gamma_{h}$ are iid with mean
zero and variances $\sigma_{\alpha}^{2}$ and $\sigma_{\gamma}^{2}$
respectively, so there are $M^{2}$ observations. Then,
$\sum_{g,h}W_{gh}/M=\left(\sum_{g}\alpha_{g}/\sqrt{M}\right)\left(\sum_{h}\gamma_{h}/\sqrt{M}\right)\xrightarrow{d}Z_{1}Z_{2}$,
where $Z_{1}$ and $Z_{2}$ are independent standard normal distributions. This
limiting distribution is also known as Gaussian chaos.
$\sum_{g}(N^{G}_{g})^{2}/\lambda_{n}=M^{3}/(M^{2}\sigma_{\alpha}^{2}\sigma_{\gamma}^{2})=M/\sigma_{\alpha}^{2}\sigma_{\gamma}^{2}\rightarrow\infty$
violates Assumption 2.3.

######  Theorem 1.

Under Assumption 1 and 2,
$Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}N(0,I_{K})$.
Further,

  1. 1.

If $E[W_{i}]=0$ $\forall i$, then $Q_{n}^{-1}\hat{Q}_{n}\xrightarrow{p}I_{K}$,
where
$\hat{Q}_{n}:=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}W_{i}W_{j}^{\prime}$.

  2. 2.

If $E[W_{i}]=\mu$ $\forall i$ and
$\frac{1}{\lambda_{n}}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}|\omega_{i}\omega_{j}|\leq
K_{0}$ for some $K_{0}<\infty$, then, for
$\bar{W}=(\sum_{i}\omega_{i}W_{i})/(\sum_{j}\omega_{j})$ and
$\hat{Q}_{n}:=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}(W_{i}-\bar{W})(W_{j}-\bar{W})^{\prime}$,
$\bar{W}\xrightarrow{p}\mu$ and $Q_{n}^{-1}\hat{Q}_{n}\xrightarrow{p}I_{K}$.

The theorem tells us that, under the aforementioned conditions,
$Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])$ is asymptotically
standard normal and the plug-in variance estimator proposed by CGM is
consistent for multi-way clustering. One-way clustering is a special case of
this theorem when one dimension is weakly nested within the other: examples
include $G=H$ so both dimensions are identical, or if we cluster by county and
state (as counties are nested in states), or if everyone is in their own $H$
cluster. A sufficient condition for consistent variance estimation is
$E[W_{i}]=0$, similiar to theorem 3 of Hansen and Lee (2019). This assumption
is sufficient in many applications: for example, linear regressions considered
in Section 3 are identified by requiring the expectation of the residual term
to be zero. Additionally, the condition $E[W_{i}]=\mu$ matches theorem 4 of
Hansen and Lee (2019) for consistent variance estimation. Theorem 1.2 uses a
stronger form of Assumption 2.3 where $A_{ij}=1$ for all
$i,j\in\mathcal{N}^{C}_{c}$.

######  Remark 2.

If $E[W_{i}]\neq 0$, then the variance estimator need not be consistent.
Unlike one-way clustering, it may not even be conservative. Suppose
$E[W_{i}]\neq 0$ for some $i$, and define $\tilde{W}_{i}:=W_{i}-E[W_{i}]$.
Then,
$Q_{n}^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}W_{i}W_{j}^{\prime}=Q_{n}^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{W}_{i}\tilde{W}_{j}^{\prime}\right)+Q_{n}^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}E[W_{i}]E[W_{j}]^{\prime}\right)$.
Since
$Q_{n}^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{W}_{i}\tilde{W}_{j}^{\prime}\right)=o_{P}(1)$
by Theorem 1.1, and $Q_{n}$ is positive semidefinite, whether the asymptotic
variance is over or under estimated depends on whether
$\sum_{i}\sum_{j\in\mathcal{N}_{i}}E[W_{i}]E[W_{j}]^{\prime}$ is positive
semidefinite. Let $K=1$ for exposition. In one-way clustering, the variance is
weakly over-estimated, so inference is conservative. To see this, let
$W^{G}_{g}$ denote the vector of $W_{i}$ such that $g(i)=g$.
$\sum_{i}\sum_{j\in\mathcal{N}_{i}}E[W_{i}]E[W_{j}]=\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}E[W_{i}]E[W_{j}]=\sum_{g}1^{\prime}E[W^{G}_{g}]E[W^{G}_{g}]^{\prime}1\geq
0$. In two-way clustering,
$\sum_{i}\sum_{j\in\mathcal{N}_{i}}E[W_{i}]E[W_{j}]$ can be negative. An
example is where $n=3$: $cov(W_{1},W_{3})=0$ but $cov(W_{1},W_{2})\neq 0$ and
$cov(W_{2},W_{3})\neq 0$, so $W_{1}$ and $W_{2}$ share a cluster in one
dimension and $W_{2}$ and $W_{3}$ share a cluster on a different dimension.
Further, $E[W_{2}]=-1$ and $E[W_{1}]=E[W_{3}]=1$. Then,
$\sum_{i}\sum_{j\in\mathcal{N}_{i}}E[W_{i}]E[W_{j}]=-1$.

###  2.1 Proof Sketch

The proof of Theorem 1 proceeds by first proving a CLT for a scalar random
variable, then applying the Cramer-Wold device to obtain the multivariate CLT.
The scalar CLT is proven using Stein’s method. I adapt the proof strategy from
Ross (2011) to obtain an upper bound on the Wasserstein distance between a
pivotal statistic and the standard normal random variable. By exploiting the
multi-way clustering structure, the upper bound on the distance can be shown
to converge to zero. All details are in Appendix A.

For ease of exposition, consider a simpler environment where $K=1$,
$\omega_{i}=1$ for all $i$, and $A_{ij}=1$ whenever $c(i)=c(j)$ for some $c$,
and $E[W_{i}]=0$. Lemma 4 in Appendix A provides an explicit bound on the
Wasserstein distance. With $d_{W}(.)$ denoting the Wasserstein distance,
$\sigma_{n}^{2}:=Q_{n}$ and $R=\sum_{i}X_{i}/\sigma_{n}$,

$\displaystyle
d_{W}(R,Z)\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[W_{i}W_{j}W_{k}]\right|+\frac{\sqrt{2}}{\sqrt{\pi}\sigma_{n}^{2}}\sqrt{Var\left(\sum_{i=1}^{n}\sum_{j\in\mathcal{N}_{i}}W_{i}W_{j}\right)}$

At this point, my proof departs from the proofs in existing statistical
literature that employ Stein’s method (e.g., Chen and Shao (2004)). Let
$N_{i}:=|\mathcal{N}_{i}|$. Holder’s inequality is employed on objects such as
$\sum_{i}|\sum_{j,k\in\mathcal{N}_{i}}E[W_{i}W_{j}W_{k}]|$. Existing
literature uses the $L^{1}$ norm of moments $E[W_{i}^{3}]$ and $L^{\infty}$
norm of $N_{i}$, resulting in
$\left(\operatornamewithlimits{max}_{m}N_{m}\right)^{2}\sum_{i}E[W_{i}^{3}]$.
In contrast, my proof uses the $L^{\infty}$ norm of $E[W_{i}^{3}]$ and $L^{1}$
norm of $N_{i}$, resulting in
$\operatornamewithlimits{max}_{m}E[W_{m}^{3}]\sum_{i}N_{i}^{2}$. Hence,

$\displaystyle\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[W_{i}W_{j}W_{k}]\right|\leq\frac{1}{\sigma_{n}^{3}}\operatornamewithlimits{max}_{m}E[W_{m}^{3}]\sum_{i}N_{i}^{2}$

Since $\operatornamewithlimits{max}_{m}E[W_{m}^{3}]$ is bounded by Assumption
2.1, it suffices to show $\sum_{i}N_{i}^{2}/\sigma_{n}^{3}\rightarrow 0$. Due
to Assumption 1, $N_{i}\leq N^{G}_{g(i)}+N^{H}_{h(i)}$, so

$\displaystyle\frac{1}{\sigma_{n}^{3}}\sum_{i}N_{i}^{2}\leq\frac{1}{\sigma_{n}^{3}}\sum_{i}(N^{G}_{g(i)}+N^{H}_{h(i)})^{2}$
$\displaystyle\leq\frac{1}{\sigma_{n}^{3}}\operatornamewithlimits{max}_{g,h}(N^{G}_{g}+N^{H}_{h})\sum_{i}(N_{g(i)}+N_{h(i)})$
$\displaystyle\leq\left[\frac{1}{\sigma_{n}}\operatornamewithlimits{max}_{g,h}(N^{G}_{g}+N^{H}_{h})\right]\frac{1}{\sigma_{n}^{2}}\left(\sum_{g}(N^{G}_{g})^{2}+\sum_{h}(N^{H}_{h})^{2}\right)$

Since $\lambda_{n}=\sigma_{n}$ when $K=1$,
$\operatornamewithlimits{max}_{g,h}(N^{G}_{g}+N^{H}_{h})/\sigma_{n}\rightarrow
0$ by Assumption 2.2 and
$\left(\sum_{g}(N^{G}_{g})^{2}+\sum_{h}(N^{H}_{h})^{2}\right)/\sigma_{n}^{2}$
is bounded by Assumption 2.3. Hence, the term is $o(1)$.

A similar argument is made for the fourth moment that features in
$Var\left(\sum_{i=1}^{n}\sum_{j\in\mathcal{N}_{i}}W_{i}W_{j}\right)$. To
complete the proof for variance estimation, observe that since the fourth
moments exist, the consistency of the plug-in variance estimator can be proven
by using Chebyshev’s inequality and existing intermediate results.

######  Remark 3.

Due to the proof strategy, the intermediate results are informative about the
quality of the normal approximation. With $d_{K}(.)$ denoting the Kolmogorov
distance, proposition 1.2 from Ross (2011) implies that
$d_{K}(R,Z)\leq(2/\pi)^{1/4}\sqrt{d_{W}(R,Z)}$. Since $Z$ is standard normal
in the proof of CLT, the bound on $d_{W}(.)$ also places a bound on the
Kolmogorov distance $d_{K}(.)$. This is then informative of the maximum
distance between the pivotal statistic and the standard normal.

##  3 Application

This section applies Theorem 1 to linear regressions, showing that using the
normal approximation with the CGM estimator is valid. Consider a linear model
where scalar outcome $Y_{i}$ is generated by

$\displaystyle
Y_{i}=D_{i}\theta+W_{i}^{\prime}\gamma+u_{i}=:X_{i}^{\prime}\beta+u_{i}$

$D_{i}\in\mathbb{R}$ is the regressor of interest, $W_{i}\in\mathbb{R}^{K-1}$
is a vector of controls that may include the intercept, and let
$X_{i}=(X_{i1},X_{i2},\cdots,X_{iK})^{\prime}:=(D_{i},W_{i}^{\prime})^{\prime}\in\mathbb{R}^{K}$.
We are interested in estimating $\theta$. The coefficient vector
$\beta:=(\theta,\gamma^{\prime})^{\prime}\in\mathbb{R}^{K}$ is the same for
all individuals. The stochastic residual term $u_{i}$ satisfies
$E[u_{i}|X_{i}]=0$ for all $i$, and is allowed to be multi-way clustered. The
standard OLS estimator is

$\displaystyle\hat{\beta}=\left(\sum_{i=1}^{n}X_{i}X_{i}^{\prime}\right)^{-1}\left(\sum_{i=1}^{n}X_{i}Y_{i}\right)=\beta+\left(\sum_{i=1}^{n}X_{i}X_{i}^{\prime}\right)^{-1}\left(\sum_{i=1}^{n}X_{i}u_{i}\right)$

This object is assumed to be well-defined in that
$\sum_{i=1}^{n}X_{i}X_{i}^{\prime}$ is invertible. Using an equivalent
representation with data matrices, the model is
$Y=D\theta+W\gamma+u=X\beta+u$. Let $M_{W}=I-W(W^{\prime}W)^{-1}W^{\prime}$
denote the annihilator matrix. Let $\tilde{D}:=M_{W}D$ be the $D$ with $W$’s
partialled out, and define $\tilde{Y},\tilde{u}$ in a similar manner. By the
Frisch-Waugh-Lovell theorem (FWL),

$\displaystyle\hat{\theta}$
$\displaystyle=(\tilde{D}^{\prime}\tilde{D})^{-1}\tilde{D}^{\prime}\tilde{Y}=\theta+(\tilde{D}^{\prime}\tilde{D})^{-1}\tilde{D}^{\prime}\tilde{u}=\theta+\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{-1}\left(\sum_{i}\tilde{D}_{i}\tilde{u}_{i}\right)=\hat{\beta}_{1}$

where $\tilde{D}_{i}$ is the $i$th component of $\tilde{D}$, so
$\sum_{i}\tilde{D}_{i}\tilde{u}_{i}=\tilde{D}^{\prime}\tilde{u}=D^{\prime}M_{W}u=\sum_{i}\tilde{D}_{i}u_{i}$.
Let
$\sigma_{n}^{2}:=Var(\hat{\theta})=Var\left(\sum_{i}\tilde{D}_{i}u_{i}/(\sum_{i^{\prime}}\tilde{D}_{i^{\prime}}^{2})\right)$
and
$\hat{\sigma}_{n}^{2}:=\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\hat{u}_{i}\hat{u}_{j}\tilde{D}_{i}\tilde{D}_{j}\right)/\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}$.
Estimated residuals are
$\hat{u}_{i}:=Y_{i}-X_{i}\hat{\beta}=u_{i}-X_{i}(\hat{\beta}-\beta)$. Due to
FWL,
$\hat{u}_{i}=\tilde{Y}_{i}-\tilde{D}_{i}\hat{\theta}=u_{i}-\tilde{D}_{i}(\hat{\theta}-\theta)$.

Inference for $\hat{\theta}$, depends on whether we are conditioning on $X$:
the conditions for asymptotic normality differ slightly between random and
nonrandom $X$. I consider each of them in turn.

###  3.1 Fixed Regressors

First, consider regressions where the $X$’s are nonrandom. An example might be
when the object of interest is the difference between male and female wages.
Their unobserved error may be correlated by state and industry conditional on
$X$, but the gender status $D_{i}$ is fixed. This can be viewed as inference
on a descriptive object.

With $u_{i}$’s having a multi-way clustered structure, we can apply Theorem 1
on $\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{-1}\sum_{i}\tilde{D}_{i}u_{i}$,
where scalar weights are given by
$\omega_{i}=\tilde{D}_{i}/(\sum_{i^{\prime}}\tilde{D}_{i^{\prime}}^{2})$.

######  Assumption 3.

For $C\in\\{G,H\\}$ and nonstochastic $\tilde{D}_{i}$, there exists
$K_{0}<\infty$ such that:

  1. 1.

$E[u_{i}^{4}]\leq K_{0}$, $E[u_{i}]=0$.

  2. 2.

$\frac{\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\tilde{D}_{i}|\right)^{2}}{\sum_{c^{\prime}}\left(\sum_{j\in\mathcal{N}^{C}_{c^{\prime}}}|\tilde{D}_{j}|\right)^{2}}\rightarrow
0$.

  3. 3.

$\frac{\sum_{c^{\prime}}\sum_{i,j\in\mathcal{N}^{C}_{c}}|\tilde{D}_{i}\tilde{D}_{j}|}{Var\left(\sum_{i}\tilde{D}_{i}u_{i}\right)}\leq
K_{0}$.

  4. 4.

$u_{i}\\!\perp\\!\\!\\!\perp u_{j}$ if $g(i)\neq g(j)$ and $h(i)\neq h(j)$.

######  Proposition 1.

Under Assumption 3, $(\hat{\theta}-\theta)/\sigma_{n}\xrightarrow{d}N(0,1)$,
and $\hat{\sigma}_{n}^{2}/\sigma_{n}^{2}\xrightarrow{p}1$.

Assumption 3 works in the environment where there is no misspecification, so
$A_{ij}=1$ whenever $i,j$ share at least one cluster. Hence,
$\sigma_{n}^{2}\asymp\operatornamewithlimits{max}_{C\in\\{G,H\\}}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}|\omega_{i}\omega_{j}|$,
satisfying the conditions of Theorem 1. Consequently, instead of making an
assumption on the contribution of the cluster with the largest weight on the
total variance, a leverage condition in the form of Assumption 3.2 can be
obtained. This condition is also empirically verifiable: the researcher can
calculate
$L_{C}:=\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\tilde{D}_{i}|\right)^{2}/\left(\sum_{c^{\prime}}\left(\sum_{j\in\mathcal{N}^{C}_{c^{\prime}}}|\tilde{D}_{j}|\right)^{2}\right)$,
and check if it is small. As a benchmark, when observations are not clustered
and all weights $\tilde{D}_{i}$ are the same, $L_{C}=1/n$. Hence, if we
believe that $n=30$ is sufficiently large for asymptotics in the iid case,
then $L_{C}<1/30$ may be acceptable.

Proposition 1 implies that the usual inference procedure is still valid even
when the unobserved component is arbitrarily heterogeneous across different
clusters. In contrast, the separate exchangeability of $u_{i}$ requires
$u_{i}$ to be identically distributed across different clusters (e.g., the
unobserved component of wages for women is identically distributed across
states) — it is a strong assumption that is no longer required here. If there
are fixed effects in the model, the vector of indicators can be collected in
$W$ and the argument proceeds as usual.111111Fixed effects account for a shift
in the unobserved component, so separate exchangeability still makes a
restriction on the distribution of the remaining unobserved component.

###  3.2 Stochastic Regressors

Next, consider stochastic $X$. This is the relevant case when considering
causal regressions. For example, we may be interested in the effect of a
randomly assigned opportunity to participate in a job training program $D_{i}$
on wages $Y_{i}$. Both $X_{i}$ and $u_{i}$ are plausibly correlated within
state and within industry. Although $\hat{\theta}=\hat{\beta}_{1}$, we can no
longer apply Theorem 1 to $\sum_{i}\tilde{D}_{i}u_{i}$ because the multi-way
dependence structure breaks once $X_{i}$’s are random.

Define $S_{n}:=\sum_{i=1}^{n}E[X_{i}X_{i}^{\prime}]$ and
$Q_{n}:=Var\left(\sum_{i=1}^{n}X_{i}u_{i}\right)$, and denote their sample
analogs as $\hat{S}_{n}=\sum_{i}X_{i}X_{i}^{\prime}$ and
$\hat{Q}_{n}:=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\hat{u}_{i}\hat{u}_{j}X_{i}X_{j}^{\prime}$.
Let the smallest eigenvalue of $Q_{n}$ be
$\lambda_{n}:=\lambda_{\operatornamewithlimits{min}}(Q_{n})$. The asymptotic
variance of $\hat{\beta}$ and its sample analog are
$V(\hat{\beta}):=S_{n}^{-1}Q_{n}S_{n}^{-1}$ and
$\hat{V}(\hat{\beta}):=\hat{S}_{n}^{-1}\hat{Q}_{n}\hat{S}_{n}^{-1}$
respectively.

Assumption 4 provides sufficient conditions for asymptotic normality of the
estimator $\hat{\beta}$ and consistency of the CGM variance estimator. The
conditions mimic Assumption 2 so that Theorem 1 is applicable to the random
vector $X_{i}u_{i}$. The new condition is a weak regularity condition that
$\lambda_{\operatornamewithlimits{min}}\left(S_{n}/n\right)\geq K_{1}>0$,
mimicking to the rank condition in OLS.

######  Assumption 4.

For $C\in\\{G,H\\}$, and $k\in\\{1,2,\cdots,K\\}$, there exists $K_{0}<\infty$
and $K_{1}>0$:

  1. 1.

$E[u_{i}^{4}|X_{i}]\leq K_{0}$, $E[X_{ik}^{4}]\leq K_{0}$, $E[u_{i}|X_{i}]=0$
for all $i$.

  2. 2.

$\frac{1}{\lambda_{n}}\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}\rightarrow
0$.

  3. 3.

$\frac{1}{\lambda_{n}}\sum_{c}(N^{C}_{c})^{2}\leq K_{0}$.

  4. 4.

$(X_{i}^{\prime},u_{i})^{\prime}\\!\perp\\!\\!\\!\perp(X_{j}^{\prime},u_{j})^{\prime}$
if $g(i)\neq g(j)$ and $h(i)\neq h(j)$.

  5. 5.

$\lambda_{\operatornamewithlimits{min}}\left(\frac{1}{n}S_{n}\right)\geq
K_{1}$.

######  Proposition 2.

Under Assumption 4,
$Q_{n}^{-1/2}S_{n}(\hat{\beta}-\beta)\xrightarrow{d}N(0,I_{K})$, and
$[S_{n}^{-1}Q_{n}S_{n}^{-1}]^{-1}[\hat{S}_{n}^{-1}\hat{Q}_{n}\hat{S}_{n}^{-1}]\xrightarrow{p}I_{K}$.

Proposition 2 is useful for doing F tests on a subvector of $\beta$. The proof
of Proposition 2 proceeds by applying Theorem 1 to $\sum_{i}X_{i}u_{i}$, and
showing that $S_{n}^{-1}\hat{S}_{n}\xrightarrow{p}I_{K}$. The latter requires
the rank condition of Assumption 4.5. It then remains to show that the
remainder terms are asymptotically negligible. Nonetheless, if we are only
interested in $\theta$, using the residualized objects $\hat{\theta}$ and
variance estimator for the residualized object $\hat{\sigma}_{n}^{2}$ is still
valid. This follows from FWL, and the refinement of FWL for variance
estimators in Ding (2021).

######  Corollary 1.

Under Assumption 4, $(\hat{\theta}-\theta)/\sigma_{n}\xrightarrow{d}N(0,1)$,
and $\hat{\sigma}_{n}/\sigma_{n}\xrightarrow{p}1$.

The practitioner’s takeaway from Proposition 2 is that the existing CGM
variance estimator can be used for valid inference with multi-way clustering.
With Corollary 1, $\hat{\theta}$ and $\hat{\sigma}_{n}^{2}$ can be used as the
mean and variance estimators respectively. These results provide the formal
theoretical guarantee for using the estimator, under weaker conditions that
permits heterogeneity across clusters.

Besides the application mentioned, Theorem 1 also has implications on the
conditions required for valid inference when the random variable is multi-way
clustered in many other econometric models, including design-based settings
and instrument variables models. Inference for estimators based on moment
conditions can be done by straightforward application of Theorem 1 as in the
linear regression case.

##  Appendix A Proof of Theorem 1

The proof strategy is as follows. I first prove Lemma 1, which is a central
limit theorem (CLT) for scalars that permits weights on the random variable.
The proof of Lemma 1 relies on Lemmas 2 to 7. Lemmas 2 to 4 derive an upper
bound on the Wasserstein distance between a pivotal statistic and standard
normal $Z$. Lemmas 5 to 7 then show that the derived upper bound is $o(1)$.
With Lemma 1, the multivariate CLT of Theorem 1 is obtained by using the
Cramer-Wold device. The remainder of the proof proceeds in the following
order: (i) introduce definitions and notation, (ii) state Lemma 1, (iii) state
and prove Lemmas 2 to 7, (iv) prove Lemma 1, (v) state and prove Lemma 8 that
is required for consistent variance estimation, then (vi) complete the proof
of Theorem 1.

The following definitions and notations are used throughout the proof. Let
$d_{W}(X,Y)$ denote the Wasserstein distance between random variables $X$ and
$Y$, so $d_{W}(X,Y)=0$ if and only if the distributions of $X$ and $Y$ are
identical. The norms of functions are defined as the sup norm i.e.,
$||f||=\sup_{x\in D}|f(x)|$. For vector $a$, $||a||=(a^{\prime}a)^{1/2}$ is
the Euclidean norm, and for positive semi-definite matrix $A$ and
$\lambda_{\operatornamewithlimits{max}}(A)$ denoting the largest eigenvalue,
$||A||=\sqrt{\lambda_{\operatornamewithlimits{max}}(A^{\prime}A)}$ denotes the
spectral norm, and $A^{1/2}$ denotes the symmetric matrix such that
$A^{1/2}A^{1/2}=A$.
$\sum_{i\in\mathcal{N}^{G}_{g}}\sum_{j\in\mathcal{N}^{G}_{g}}$ is abbreviated
as $\sum_{i,j\in\mathcal{N}^{G}_{g}}$. The dependency neighborhood of $i$,
$\mathcal{N}_{i}\subseteq\\{1,\cdots,n\\}$, is defined as the set of
observations where $i\in\mathcal{N}_{i}$ and $X_{i}$ is independent of
$\\{X_{j}\\}_{j\neq\mathcal{N}_{i}}$, and $N_{i}:=|\mathcal{N}_{i}|$ is the
number of observations in $i$’s dependency neighborhood. In the rest of this
proof, $X_{i}$ denotes a scalar random variable while $W_{i}\in\mathbb{R}^{K}$
as stated in the main text is a random vector.

Every scalar random variable $X_{i}$ is weighted by nonstochastic
$\omega_{i}$. Denote the variance of the sum as
$\sigma_{n}^{2}:=Var\left(\sum_{i=1}^{n}\omega_{i}X_{i}\right)$. We are
interested in the asymptotic distribution of
$(1/\sigma_{n})\sum_{i=1}^{n}\omega_{i}X_{i}$. If all observations are equally
weighted, $\omega_{i}=1$ $\forall i$.

######  Assumption 5.

For $C\in\\{G,H\\}$, there exists $K_{0}<\infty$ such that:

  1. 1.

$E[X_{i}]=0$ and $E[X_{i}^{4}]\leq K_{0}<\infty$ for all $i$.

  2. 2.

$\frac{1}{\sigma_{n}^{2}}\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\omega_{i}|\right)^{2}\rightarrow
0$

  3. 3.

$\frac{1}{\sigma_{n}^{2}}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}A_{ij}|\omega_{i}\omega_{j}|\leq
K_{0}<\infty$

  4. 4.

$X_{i}\\!\perp\\!\\!\\!\perp X_{j}$ if $g(i)\neq g(j)$ and $h(i)\neq h(j)$.

######  Lemma 1.

Under Assumption 5,
$(1/\sigma_{n})\sum_{i=1}^{n}\omega_{i}X_{i}\xrightarrow{d}N(0,1)$, where
$\sigma_{n}^{2}:=Var\left(\sum_{i=1}^{n}\omega_{i}X_{i}\right)$. Further,
using feasible estimator
$\hat{\sigma}_{n}^{2}:=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}X_{j}$,
$\hat{\sigma}_{n}^{2}/\sigma_{n}^{2}\xrightarrow{p}1$.

######  Lemma 2.

If $R$ is a random variable and $Z$ has a standard normal distribution, and we
define the family of functions
$\mathcal{F}=\\{f:||f||,||f^{\prime\prime}||\leq
2,||f^{\prime}||\leq\sqrt{2\pi}\\}$, then
$d_{W}(R,Z)\leq\sup_{f\in\mathcal{F}}|E[f^{\prime}(R)-Rf(R)]|$.

###### Proof.

See Ross (2011) theorem 3.1. ∎

######  Lemma 3.

Let $X_{1},\cdots,X_{n}$ be random variables such that
$E[X_{i}]=0,\sigma_{n}^{2}=Var(\sum_{i}X_{i})$, and define
$R=\sum_{i}X_{i}/\sigma_{n}$. If
$R_{i}:=\sum_{j\neq\mathcal{N}_{i}}X_{j}/\sigma_{n}$, then

$\displaystyle E[Rf(R)]$
$\displaystyle=E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(f(R)-f(R_{i})-(R-R_{i})f^{\prime}(R))\right]+E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(R-R_{i})f^{\prime}(R)\right]$

###### Proof.

Start from right hand side:

$\displaystyle E$
$\displaystyle\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(f(R)-f(R_{i})-(R-R_{i})f^{\prime}(R))\right]+E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(R-R_{i})f^{\prime}(R)\right]$
$\displaystyle=E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(f(R)-f(R_{i}))\right]=E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}f(R)\right]+E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}f(R_{i})\right]$
$\displaystyle=E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}f(R)\right]=E[Rf(R)]$

The first equality in the final line comes from the fact that $R_{i}$ is
independent of $X_{i}$ based on how dependency neighborhoods are defined.
Hence, $E[X_{i}f(R_{i})]=0$. ∎

######  Lemma 4.

Let $X_{1},\cdots,X_{n}$ be random variables such that,
$E[X_{i}]=0,\sigma_{n}^{2}=Var(\sum_{i}X_{i})$, and define
$R=\sum_{i}X_{i}/\sigma_{n}$. Let the collection $(X_{1},\cdots,X_{n})$ have
dependency neighborhoods $\mathcal{N}_{i}$, $i=1,\cdots,n$. Then for $Z$ a
standard normal random variable,

$d_{W}(R,Z)\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[X_{i}X_{j}X_{k}]\right|+\frac{\sqrt{2}}{\sqrt{\pi}\sigma_{n}^{2}}\sqrt{Var\left(\sum_{i=1}^{n}\sum_{j\in\mathcal{N}_{i}}X_{i}X_{j}\right)}$
(1)

###### Proof.

Due to Lemma 2, to bound $d_{W}(R,Z)$ from above, it is sufficient to bound
$|E[f^{\prime}(R)-Rf(R)]|$, where $||f||,||f^{\prime\prime}||\leq
2,||f^{\prime}||\leq\sqrt{2/\pi}$. Define $R_{i}:=\sum_{j\neq
N_{i}}X_{j}/\sigma_{n}$, so $X_{i}$ is independent of $R_{i}$.

$\displaystyle|E[$ $\displaystyle
f^{\prime}(R)-Rf(R)]|=|E[f^{\prime}(R)]-E[Rf(R)]|$
$\displaystyle\leq\left|E[f^{\prime}(R)]-E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(f(R)-f(R_{i})-(R-R_{i})f^{\prime}(R))\right]-E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(R-R_{i})f^{\prime}(R)\right]\right|$
$\displaystyle\leq\left|E\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(f(R)-f(R_{i})-(R-R_{i})f^{\prime}(R))\right]\right|+\left|E\left[f^{\prime}(R)\left(1-\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(R-R_{i})\right)\right]\right|$

The first inequality applies Lemma 3, and the second inequality applies the
triangle inequality. Consequently, it is sufficient to show that the first
term is bounded by the corresponding first term of Equation (1), and the
second term is bounded by the corresponding second term.

Consider the first term. By Taylor expansion of $f(R_{i})$ around $f(R)$, and
the triangle inequality, the term that generates the third moment is:

$\displaystyle|E$
$\displaystyle\left[\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(f(R)-f(R_{i})-(R-R_{i})f^{\prime}(R))\right]|\leq\frac{||f^{\prime\prime}||}{2\sigma_{n}}\left|\sum_{i=1}^{n}E[X_{i}(R-R_{i})^{2}]\right|$
$\displaystyle=\frac{1}{\sigma_{n}^{3}}\left|\sum_{i=1}^{n}E\left[X_{i}\left(\sum_{j\in\mathcal{N}_{i}}X_{j}\right)^{2}\right]\right|=\frac{1}{\sigma_{n}^{3}}\left|\sum_{i=1}^{n}\sum_{j,k\in\mathcal{N}_{i}}E[X_{i}X_{j}X_{k}]\right|\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[X_{i}X_{j}X_{k}]\right|$

Turning now to the second term,

$\displaystyle\left|E\left[f^{\prime}(R)\left(1-\frac{1}{\sigma_{n}}\sum_{i=1}^{n}X_{i}(R-R_{i})\right)\right]\right|\leq\frac{||f^{\prime}||}{\sigma_{n}^{2}}\left|E\left[\sigma_{n}^{2}-\sigma_{n}\sum_{i=1}^{n}X_{i}(R-R_{i})\right]\right|$
$\displaystyle\leq\frac{||f^{\prime}||}{\sigma_{n}^{2}}E\left|\sigma_{n}^{2}-\sum_{i=1}^{n}X_{i}\left(\sum_{j\in
N_{i}}X_{j}\right)\right|\leq\frac{||f^{\prime}||}{\sigma_{n}^{2}}E\left[\left(\sigma_{n}^{2}-\sum_{i=1}^{n}X_{i}\left(\sum_{j\in
N_{i}}X_{j}\right)\right)^{2}\right]^{1/2}1^{1/2}$
$\displaystyle\leq\frac{\sqrt{2}}{\sqrt{\pi}\sigma_{n}^{2}}\sqrt{Var\left(\sum_{i=1}^{n}\sum_{j\in
N_{i}}X_{i}X_{j}\right)}$

∎

######  Lemma 5.

$E[|X_{i}X_{j}X_{k}|]\leq\operatornamewithlimits{max}_{m}E[|X_{m}|^{3}]$,
$E[|X_{i}X_{j}X_{k}X_{l}|]\leq\operatornamewithlimits{max}_{m}E[|X_{m}|^{4}]$,
and
$|E[X_{i}X_{k}]E[X_{j}X_{l}]|\leq\operatornamewithlimits{max}_{m}E[|X_{m}|^{4}]$.

###### Proof.

By the arithmetic mean — geometric mean (AM-GM) inequality,

$E|X_{i}X_{j}X_{k}|\leq\frac{1}{3}\left(E|X_{i}|^{3}+E|X_{j}|^{3}+E|X_{k}|^{3}\right)\leq\operatornamewithlimits{max}_{m}E[|X_{m}|^{3}]$

A similar argument yields
$E[|X_{i}X_{j}X_{k}X_{l}|]\leq\operatornamewithlimits{max}_{m}E[|X_{m}|^{4}]$.
For the final result, first observe that $E[X_{i}X_{k}]^{2}\pm
2E[X_{i}X_{k}]E[X_{j}X_{l}]+E[X_{j}X_{l}]^{2}=(E[X_{i}X_{k}]\pm
E[X_{j}X_{l}])^{2}\geq 0$. Hence,

$\displaystyle|E[X_{i}X_{k}]E[X_{j}X_{l}]|$
$\displaystyle\leq\frac{1}{2}(E[X_{i}X_{k}]^{2}+E[X_{j}X_{l}]^{2})\leq\frac{1}{2}(E[X_{i}^{2}X_{k}^{2}]+E[X_{j}^{2}X_{l}^{2}])$
$\displaystyle\leq\frac{1}{4}(E[X_{i}^{4}]+E[X_{j}^{4}]+E[X_{k}^{4}]+E[X_{l}^{4}])\leq\operatornamewithlimits{max}_{m}E[X_{m}^{4}]$

∎

######  Lemma 6.

Under Assumption 5,
$\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[\omega_{i}\omega_{j}\omega_{k}X_{i}X_{j}X_{k}]\right|\rightarrow
0$.

###### Proof.

Note that $E[X_{i}X_{j}X_{k}]=0$ whenever one of $\\{X_{i},X_{j},X_{k}\\}$ is
independent of the other two, so
$E[\omega_{i}\omega_{j}\omega_{k}X_{i}X_{j}X_{k}]$ is nonzero only if
$A_{ij}$, $A_{ik}$, or $A_{jk}$ is nonzero. Apply the triangle inequality and
push the absolute value into the expectation.

$\displaystyle\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[\omega_{i}\omega_{j}\omega_{k}X_{i}X_{j}X_{k}]\right|$
$\displaystyle\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}(A_{ij}+A_{jk}+A_{ik})E[\omega_{i}\omega_{j}\omega_{k}X_{i}X_{j}X_{k}]\right|$
$\displaystyle\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\sum_{j,k\in\mathcal{N}_{i}}(A_{ij}+A_{jk}+A_{ik})|\omega_{i}\omega_{j}\omega_{k}|E[|X_{i}X_{j}X_{k}|]$
$\displaystyle\leq\frac{\operatornamewithlimits{max}_{m}E[|X_{m}|^{3}]}{\sigma_{n}^{3}}\sum_{i=1}^{n}\sum_{j,k\in\mathcal{N}_{i}}|\omega_{i}\omega_{j}\omega_{k}|(A_{ij}+A_{jk}+A_{ik})$

The last inequality applies Lemma 5. Observe
$\operatornamewithlimits{max}_{m}E[|X_{m}|^{3}]\leq K_{0}$ since the 4th
moment exists, so it remains to show that the remaining terms are $o(1)$.

$\displaystyle\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\sum_{j,k\in\mathcal{N}_{i}}(A_{ij}+A_{jk}+A_{ik})|\omega_{i}\omega_{j}\omega_{k}|$
$\displaystyle\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left(\sum_{j,k\in\mathcal{N}^{G}_{g(i)}}+\sum_{j,k\in\mathcal{N}^{H}_{h(i)}}\right)(A_{ij}+A_{jk}+A_{ik})|\omega_{i}\omega_{j}\omega_{k}|$

It is sufficient to consider the $G$ dimension as the $H$ dimension is
analogous.

$\displaystyle\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\sum_{j,k\in\mathcal{N}^{G}_{g(i)}}(A_{ij}+A_{jk}+A_{ik})|\omega_{i}\omega_{j}\omega_{k}|$
$\displaystyle=\frac{3}{\sigma_{n}^{3}}\sum_{g}\sum_{i,j,k\in\mathcal{N}^{G}_{g}}A_{ij}|\omega_{i}\omega_{j}\omega_{k}|$
$\displaystyle\frac{1}{\sigma_{n}^{3}}\sum_{g}\sum_{i,j,k\in\mathcal{N}^{G}_{g}}A_{ij}|\omega_{i}\omega_{j}||\omega_{k}|$
$\displaystyle\leq\left(\frac{\operatornamewithlimits{max}_{g}\sum_{k\in\mathcal{N}^{G}_{g}}|\omega_{k}|}{\sigma_{n}}\right)\frac{1}{\sigma_{n}^{2}}\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}A_{ij}|\omega_{i}\omega_{j}|=o(1)$

Convergence occurs because
$(1/\sigma_{n}^{2})\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}A_{ij}|\omega_{i}\omega_{j}|<\infty$
by Assumption 5.3 and
$\operatornamewithlimits{max}_{g}\sum_{k\in\mathcal{N}^{G}_{g}}|\omega_{k}|/\sigma_{n}=\left(\operatornamewithlimits{max}_{g}\left(\sum_{k\in\mathcal{N}^{G}_{g}}|\omega_{k}|\right)^{2}/\sigma_{n}^{2}\right)^{1/2}=o(1)$
by Assumption 5.2. ∎

######  Lemma 7.

Under Assumption 5,
$\frac{1}{\sigma_{n}^{4}}Var\left(\sum_{i=1}^{n}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}X_{j}\right)=o(1)$.

###### Proof.

$\displaystyle\frac{1}{\sigma_{n}^{4}}$ $\displaystyle
Var\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}X_{j}\right)=\frac{1}{\sigma_{n}^{4}}E\left[\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}X_{j}\right)^{2}\right]-\frac{1}{\sigma_{n}^{4}}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}E[\omega_{i}\omega_{j}X_{i}X_{j}]\right)^{2}$
$\displaystyle=\frac{1}{\sigma_{n}^{4}}\sum_{i}\sum_{j}\sum_{k\in\mathcal{N}_{i}}\sum_{l\in\mathcal{N}_{j}}(E[\omega_{i}\omega_{j}\omega_{k}\omega_{l}X_{i}X_{j}X_{k}X_{l}]-E[\omega_{i}\omega_{k}X_{i}X_{k}]E[\omega_{j}\omega_{l}X_{j}X_{l}])$
$\displaystyle=\frac{1}{\sigma_{n}^{4}}\sum_{i}\sum_{j}\sum_{k\in\mathcal{N}_{i}}\sum_{l\in\mathcal{N}_{j}}\omega_{i}\omega_{j}\omega_{k}\omega_{l}(E[X_{i}X_{j}X_{k}X_{l}]-E[X_{i}X_{k}]E[X_{j}X_{l}])$

When $(X_{i},X_{k})\\!\perp\\!\\!\\!\perp(X_{j},X_{l})$,
$E[X_{i}X_{j}X_{k}X_{l}]=E[X_{i}X_{j}]E[X_{k}X_{l}]$. Hence, we only have to
consider where there is at least one pair that is correlated i.e., when
$A_{ij}$, $A_{il}$, $A_{kj}$, or $A_{kl}$ is not zero. As before, with finite
4th moment and Lemma 5, it is sufficient to show

$\displaystyle\frac{1}{\sigma_{n}^{4}}\sum_{i}\sum_{j}\sum_{k\in\mathcal{N}_{i}}\sum_{l\in\mathcal{N}_{j}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|(A_{ij}+A_{il}+A_{kj}+A_{kl})=o(1)$

It is sufficient to consider the $A_{ij}$ term because everything else is
analogous.

$\displaystyle\sum_{i}\sum_{j}\sum_{k\in\mathcal{N}_{i}}\sum_{l\in\mathcal{N}_{j}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}\leq\sum_{i}\left(\sum_{j\in\mathcal{N}^{G}_{g(i)}}+\sum_{j\in\mathcal{N}^{H}_{h(i)}}\right)\left(\sum_{k\in\mathcal{N}^{G}_{g(i)}}+\sum_{k\in\mathcal{N}^{H}_{h(i)}}\right)\left(\sum_{l\in\mathcal{N}^{G}_{g(j)}}+\sum_{l\in\mathcal{N}^{H}_{h(j)}}\right)|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}$

The first and last terms of the summation take the form:

$\displaystyle\sum_{i}\sum_{j\in\mathcal{N}^{G}_{g(i)}}\sum_{k\in\mathcal{N}^{G}_{g(i)}}\sum_{l\in\mathcal{N}^{G}_{g(j)}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}$
$\displaystyle=\sum_{g}\sum_{i,j,k,l\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}\leq\left(\operatornamewithlimits{max}_{g}\sum_{k,l\in\mathcal{N}^{G}_{g}}|\omega_{k}||\omega_{l}|\right)\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}|A_{ij}$

Since
$\frac{1}{\sigma_{n}^{2}}\operatornamewithlimits{max}_{h}\sum_{i,k\in\mathcal{N}^{H}_{h}}|\omega_{i}||\omega_{k}|=o(1)$
and
$\frac{1}{\sigma_{n}^{2}}\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}|A_{ij}<\infty$
by Assumption 5, these terms are $o(1)$ when divided by $\sigma_{n}^{4}$.

The interactive terms have the form:

$\displaystyle\sum_{i}$
$\displaystyle\sum_{j\in\mathcal{N}^{G}_{g(i)}}\sum_{k\in\mathcal{N}^{G}_{g(i)}}\sum_{l\in\mathcal{N}^{H}_{h(j)}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}$
$\displaystyle=\sum_{i,j,k}\sum_{g}1[i\in\mathcal{N}^{G}_{g}]1[j\in\mathcal{N}^{G}_{g}]1[k\in\mathcal{N}^{G}_{g}]\sum_{l}\sum_{h}1[j\in\mathcal{N}^{H}_{h}]1[l\in\mathcal{N}^{H}_{h}]|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}$
$\displaystyle=\sum_{j}\sum_{i,k}\sum_{g}1[i\in\mathcal{N}^{G}_{g}]1[j\in\mathcal{N}^{G}_{g}]1[k\in\mathcal{N}^{G}_{g}]A_{ij}|\omega_{i}\omega_{j}\omega_{k}|\sum_{h}\sum_{l}1[j\in\mathcal{N}^{H}_{h}]1[l\in\mathcal{N}^{H}_{h}]|\omega_{l}|$
$\displaystyle\leq\left(\operatornamewithlimits{max}_{j}\sum_{h}\sum_{l}1[j\in\mathcal{N}^{H}_{h}]1[l\in\mathcal{N}^{H}_{h}]|\omega_{l}|\right)\left(\sum_{g}\sum_{i,j,k\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}\omega_{k}|A_{ij}\right)$
$\displaystyle=\left(\operatornamewithlimits{max}_{h}\sum_{l\in\mathcal{N}^{H}_{h}}|\omega_{l}|\right)\left(\sum_{g}\sum_{i,j,k\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}\omega_{k}|A_{ij}\right)$
$\displaystyle=\left(\operatornamewithlimits{max}_{h}\sum_{l\in\mathcal{N}^{H}_{h}}|\omega_{l}|\right)\left(\operatornamewithlimits{max}_{g}\sum_{k\in\mathcal{N}^{G}_{g}}|\omega_{k}|\right)\left(\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}|A_{ij}\right)$

Since
$\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}|A_{ij}/\sigma_{n}^{2}\leq
K_{0}$ and
$\operatornamewithlimits{max}_{g}\sum_{k\in\mathcal{N}^{G}_{g}}|\omega_{k}|/\sigma_{n}=o(1)$,

$\displaystyle\frac{1}{\sigma_{n}^{4}}$
$\displaystyle\sum_{i}\sum_{j\in\mathcal{N}^{G}_{g(i)}}\sum_{k\in\mathcal{N}^{G}_{g(i)}}\sum_{l\in\mathcal{N}^{H}_{h(j)}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|A_{ij}$
$\displaystyle\leq\left(\frac{1}{\sigma_{n}}\operatornamewithlimits{max}_{h}\sum_{l\in\mathcal{N}^{H}_{h}}|\omega_{l}|\right)\left(\frac{1}{\sigma_{n}}\operatornamewithlimits{max}_{g}\sum_{k\in\mathcal{N}^{G}_{g}}|\omega_{k}|\right)\left(\frac{1}{\sigma_{n}^{2}}\sum_{g}\sum_{i,j\in\mathcal{N}^{G}_{g}}|\omega_{i}\omega_{j}|A_{ij}\right)=o(1)$

∎

###### Proof of Lemma 1.

Apply Lemma 4 on random variable $\omega_{i}X_{i}$ to obtain:

$\displaystyle
d_{W}(R,Z)\leq\frac{1}{\sigma_{n}^{3}}\sum_{i=1}^{n}\left|\sum_{j,k\in\mathcal{N}_{i}}E[\omega_{i}\omega_{j}\omega_{k}X_{i}X_{j}X_{k}]\right|+\frac{\sqrt{2}}{\sqrt{\pi}\sigma_{n}^{2}}\sqrt{Var\left(\sum_{i=1}^{n}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}X_{j}\right)}$

Applying Lemma 6 and 7 on each of the two terms, $d_{W}(R,Z)=o(1)$. Proof for
consistency of the variance estimator is equivalent to proving that
$(\hat{\sigma}_{n}^{2}-\sigma_{n}^{2})/\sigma_{n}^{2}=o_{P}(1)$. By
Chebyshev’s inequality,

$\displaystyle
P\left(\frac{\hat{\sigma}_{n}^{2}-\sigma_{n}^{2}}{\sigma_{n}^{2}}>\epsilon\right)$
$\displaystyle\leq\frac{1}{\epsilon^{2}}\frac{1}{\sigma_{n}^{4}}E[(\hat{\sigma}_{n}^{2}-\sigma_{n}^{2})^{2}]=\frac{Var\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}X_{j}\right)}{\epsilon^{2}\sigma_{n}^{4}}=o_{P}(1)$

The convergence in the last step occurs by Lemma 7. ∎

######  Lemma 8.

Under Assumption 1, 2.1 and 2.2, $\forall i$,
$||(1/(\sum_{i}\omega_{i}))\sum_{i}\omega_{i}(W_{i}-E[W_{i}])||\xrightarrow{p}0$.

###### Proof.

It suffices to show convergence elementwise. Let $X_{i}$ denote a scalar
components of $W_{i}$, i.e., $X_{i}=W_{im}$, where $m\in\\{1,2,\cdots,K\\}$.
By Chebyshev’s inequality, and
$\operatornamewithlimits{max}_{m,k}E[W_{mk}^{2}]<K_{0}$,

$\displaystyle P$
$\displaystyle\left(\frac{1}{\sum_{i}\omega_{i}}\sum_{i}\omega_{i}(X_{i}-E[X_{i}])>\epsilon\right)\leq\frac{1}{\epsilon^{2}}\frac{1}{\left(\sum_{i}\omega_{i}\right)^{2}}E\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}(X_{i}-E[X_{i}])(Y_{i}-E[Y_{i}])\right)$
$\displaystyle\leq\frac{K_{0}}{\epsilon^{2}\left(\sum_{j}\omega_{j}\right)^{2}}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}$

Hence, it suffices to show
$(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j})/\left(\sum_{j}\omega_{j}\right)^{2}=o(1)$.
Observe

$\frac{\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}}{\left(\sum_{j}\omega_{j}\right)^{2}}\leq\frac{\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{j}|}{\left(\sum_{j}\omega_{j}\right)}\frac{\left(\sum_{j}\omega_{j}\right)}{\left(\sum_{j}\omega_{j}\right)}$

so it suffices to show
$\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{j}|/\left(\sum_{j}\omega_{j}\right)=o(1)$.
Since
$\lambda_{n}\leq\sum_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{i}\omega_{j}|\operatornamewithlimits{max}_{m}E[W_{mk}^{2}]\leq\left(\sum_{j}|\omega_{j}|\right)^{2}\operatornamewithlimits{max}_{m}E[W_{mk}^{2}]$,

$\displaystyle\frac{\left(\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{j}|\right)^{2}}{\left(\sum_{j}\omega_{j}\right)^{2}}=\frac{\left(\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{j}|\right)^{2}\operatornamewithlimits{max}_{m}E[W_{mk}^{2}]}{\left(\sum_{j}\omega_{j}\right)^{2}\operatornamewithlimits{max}_{m}E[W_{mk}^{2}]}\leq\operatornamewithlimits{max}_{m}E[W_{mk}^{2}]\frac{\left(\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{j}|\right)^{2}}{\lambda_{n}}=o(1)$

Convergence occurs due to Assumption 2.2 and
$\operatornamewithlimits{max}_{m}E[W_{mk}^{2}]<K_{0}$. ∎

###### Proof of Theorem 1.

To show that
$Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}N(0,I_{K})$,
due to the Cramer-Wold device, it suffices to show that $\forall
l\in\mathbb{R}^{K}$,
$l^{\prime}Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}l^{\prime}N(0,I_{K})$.
If $l$ is a vector of zeroes, then
$l^{\prime}Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}l^{\prime}N(0,I_{K})$
is immediate. For $||l||>0$, it suffices to show
$(1/||l||)l^{\prime}Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}(1/||l||)l^{\prime}N(0,I_{K})=N(0,1)$.
For all nonstochastic $l\in\mathbb{R}^{K}\backslash\\{0\\}$, let
$\sigma_{n}^{2}(l):=Var\left(\sum_{i}(l/||l||)^{\prime}\left(Q_{n}/\lambda_{n}\right)^{-1/2}\omega_{i}(W_{i}-E[W_{i}])\right)$,
so the following hold:

  1. 1.

$E\left[\left(\left(\frac{l}{||l||}\right)^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1/2}(W_{i}-E[W_{i}])\right)\right]=0$
and
$E\left[\left(\left(\frac{l}{||l||}\right)^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1/2}(W_{i}-E[W_{i}])\right)^{4}\right]\leq
K_{0}$ for all $i$.

  2. 2.

$\frac{1}{\sigma_{n}^{2}(l)}\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\omega_{i}|\right)^{2}\rightarrow
0$.

  3. 3.

$\frac{1}{\sigma_{n}^{2}(l)}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}A_{ij}|\omega_{i}\omega_{j}|\leq
K_{0}$.

  4. 4.

$\left(\left(\frac{l}{||l||}\right)^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1/2}(W_{i}-E[W_{i}])\right)\\!\perp\\!\\!\\!\perp\left(\left(\frac{l}{||l||}\right)^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1/2}W_{j}\right)$
if $g(i)\neq g(j)$ and $h(i)\neq h(j)$.

For item 1, since
$\lambda_{n}:=\lambda_{\operatornamewithlimits{min}}(Q_{n})$, all eigenvalues
of $Q_{n}/\lambda_{n}$ must be at least 1. Hence, all eigenvalues of
$(Q_{n}/\lambda_{n})^{-1/2}$ are bounded above by 1. This implies
$|(l/||l||)^{\prime}(Q_{n}/\lambda_{n})^{-1/2}|\leq K_{1}$ for some arbitrary
constant $K_{1}<\infty$. Item 1 then follows from Assumption 2.1. Observe that
$\sigma_{n}^{2}(l)=(l/||l||)^{\prime}(Q_{n}/\lambda_{n})^{-1/2}Q_{n}(Q_{n}/\lambda_{n})^{-1/2}(l/||l||)=1/\lambda_{n}$.
Then, Assumption 2.2 yields item 2 and Assumption 2.3 yields item 3. Item 4 is
immediate from Assumption 1. By applying Lemma 1,
$(1/\sigma_{n}(l))(l/||l||)^{\prime}(Q_{n}/\lambda_{n})^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}N(0,1)$.
By using $\sigma_{n}^{2}(l)=1/\lambda_{n}$, this is equivalent to
$(l/||l||)^{\prime}Q_{n}^{-1/2}\sum_{i=1}^{n}\omega_{i}(W_{i}-E[W_{i}])\xrightarrow{d}N(0,1)$
as required.

Proof of Theorem 1.1

Turning to consistent variance estimation, it suffices to show that for all
$l\in\mathbb{R}^{K}$ such that $||l||=1$,
$P(l^{\prime}Q_{n}^{-1}(\hat{Q}_{n}-Q_{n})l>\epsilon)\rightarrow 0$. Now,
impose the assumption that $E[W_{i}]=0$.

$\displaystyle P(l^{\prime}Q_{n}^{-1}(\hat{Q}_{n}-Q_{n})l>\epsilon)$
$\displaystyle\leq\frac{1}{\epsilon^{2}}E\left[\left(l^{\prime}(Q_{n}^{-1}(\hat{Q}_{n}-Q_{n}))\right)^{2}\right]$
$\displaystyle=\frac{1}{\epsilon^{2}}E\left[\left(l^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\frac{1}{\lambda_{n}}(\hat{Q}_{n}-Q_{n})\right)^{2}\right]\leq\frac{1}{\epsilon^{2}}E\left[\left(l_{0}^{\prime}\frac{1}{\lambda_{n}}(\hat{Q}_{n}-Q_{n})\right)^{2}\right]$

where $l_{0}$ is a vector whose entries are all bounded above by some
arbitrary constant $K_{1}<\infty$ by a similar argument as before. Hence, it
suffices to show that
$(1/\lambda_{n})(\hat{Q}_{n}-Q_{n})\xrightarrow{p}0_{K\times K}$, where
$0_{K\times K}$ is a $K\times K$ matrix of zeroes. Since
$\hat{Q}_{n}-Q_{n}=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}W_{i}W_{j}^{\prime}-E[\omega_{i}\omega_{j}W_{i}W_{j}^{\prime}]$,
it suffices to show convergence elementwise. Let $X_{i}$ and $Y_{i}$ denote
scalar components of $W_{i}$, i.e., $X_{i}=W_{im},Y_{i}=W_{ip}$, where
$m,p\in\\{1,2,\cdots,K\\}$.

$\displaystyle P$
$\displaystyle\left(\frac{1}{\lambda_{n}}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}(X_{i}Y_{j}-E[X_{i}Y_{j}])>\epsilon\right)\leq\frac{1}{\epsilon^{2}}\frac{1}{\lambda_{n}^{2}}Var\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}X_{i}Y_{j}\right)$
$\displaystyle\leq\frac{1}{\epsilon^{2}\lambda_{n}^{2}}\sum_{i}\sum_{j}\sum_{k\in\mathcal{N}_{i}}\sum_{l\in\mathcal{N}_{j}}\left|E[\omega_{i}\omega_{j}\omega_{k}\omega_{l}X_{i}X_{j}Y_{k}Y_{l}]-E[\omega_{i}\omega_{k}X_{i}Y_{k}]E[\omega_{j}\omega_{l}X_{j}Y_{l}]\right|$
$\displaystyle\leq\frac{K_{0}}{\lambda_{n}^{2}}\sum_{i}\sum_{j}\sum_{k\in\mathcal{N}_{i}}\sum_{l\in\mathcal{N}_{j}}|\omega_{i}\omega_{j}\omega_{k}\omega_{l}|(A_{ij}+A_{il}+A_{kj}+A_{kl})=o(1)$

The inequality in the last line is obtained due to Holder’s inequality and
finite moments. An argument similar to that of Lemma 7 yields the $o(1)$
equality.

Proof of Theorem 1.2

Now assume $E[W_{i}]=\mu$. Using Lemma 8, $\bar{W}\xrightarrow{p}\mu$ is
immediate, i.e., $\bar{W}=\mu+o_{P}(1)$. To ease notation, let
$\tilde{W}_{i}:=W_{i}-\mu$. Hence,
$Q_{n}=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}E[\tilde{W}_{i}\tilde{W}_{j}^{\prime}]$.

$\displaystyle\hat{Q}_{n}$
$\displaystyle=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}(W_{i}-\bar{W})(W_{j}-\bar{W})^{\prime}=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}(\tilde{W}_{i}+o_{P}(1))(\tilde{W}_{j}+o_{P}(1))^{\prime}$
$\displaystyle=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}\tilde{W}_{i}\tilde{W}_{j}^{\prime}+2\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}\tilde{W}_{i}1_{K}^{\prime}o_{P}(1)+\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}1_{K}1_{K}^{\prime}o_{P}(1)$

Since
$Q_{n}^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}\tilde{W}_{i}\tilde{W}_{j}^{\prime}=1+o_{P}(1)$
by Theorem 1.1, it then remains to show that each of the two remaining terms
are $o_{P}(1)$ when pre-multiplied by $Q_{n}^{-1}$.

$\displaystyle\left(\frac{1}{\lambda}_{n}Q_{n}\right)^{-1}\frac{1}{\lambda}_{n}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}1_{K}1_{K}^{\prime}$
$\displaystyle\leq
K_{0}\left(\frac{1}{\lambda}_{n}Q_{n}\right)^{-1}1_{K}1_{K}^{\prime}=O(1)1_{K}1_{K}^{\prime}$

The first inequality is due to the assumption that
$(1/\lambda_{n})\sum_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{i}\omega_{j}|\leq
K_{0}$, and the $O(1)$ term occurs due to the eigenvalues of
$(Q_{n}/\lambda_{n})^{-1}$ being bounded above by 1. Take some component
$\tilde{X}_{i}$ of $\tilde{W}_{i}$. For all $\epsilon>0$, there exists
$M_{\epsilon}=K_{0}^{2}/\epsilon<\infty$ such that:

$\displaystyle
P\left(\left|\frac{1}{\lambda_{n}}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}\tilde{X}_{i}\right|\geq
M_{\epsilon}\right)$
$\displaystyle\leq\frac{1}{\lambda_{n}M_{\epsilon}}E\left[\left|\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}\tilde{X}_{i}\right|\right]$
$\displaystyle\leq\frac{1}{M_{\epsilon}}\operatornamewithlimits{max}_{i}E[|\tilde{X}_{i}|]\frac{1}{\lambda_{n}}\sum_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{i}\omega_{j}|\leq\frac{K_{0}}{K_{0}^{2}/\epsilon}=\epsilon$

Hence,
$Q_{n}^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}\tilde{W}_{i}1_{K}^{\prime}=1_{K}1_{K}^{\prime}O_{P}(1)$.
Since $O_{P}(1)o_{P}(1)=o_{P}(1)$, the result is obtained. ∎

##  Appendix B Proof of Propositions

###### Proof of Proposition 1.

We have
$\hat{\theta}-\theta=\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{-1}\left(\sum_{i}\tilde{D}_{i}u_{i}\right)=\sum_{i}\omega_{i}u_{i}$,
where $\omega_{i}:=\tilde{D}_{i}/(\sum_{j}\tilde{D}_{j}^{2})$. Let
$\sigma_{n}^{2}:=Var(\omega_{i}u_{i})$. Apply Theorem 1 with $K=1$ to
$\sum_{i}\omega_{i}u_{i}$. Assumption 1 and Assumption 2.1 are automatically
satisfied for clustered random variable $u_{i}$ and weight $\omega_{i}$.
Assumption 2.2 is satisfied because

$\displaystyle\frac{1}{\sigma_{n}^{2}}\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\omega_{i}|\right)^{2}\leq\frac{\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\tilde{D}_{i}|\right)^{2}}{K_{0}\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}\sum_{c^{\prime}}\left(\sum_{j\in\mathcal{N}^{C}_{c^{\prime}}}|\tilde{D}_{j}|\right)^{2}}=\frac{\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}\operatornamewithlimits{max}_{c}\left(\sum_{i\in\mathcal{N}^{C}_{c}}|\tilde{D}_{i}|\right)^{2}}{\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}\sum_{c^{\prime}}\left(\sum_{j\in\mathcal{N}^{C}_{c^{\prime}}}|\tilde{D}_{j}|\right)^{2}}\rightarrow
0$

where the first inequality comes from Assumption 3.3 and convergence occurs
due to Assumption 3.2. Assumption 2.3 is satisfied because

$\displaystyle\frac{1}{\sigma_{n}^{2}}\sum_{c}\sum_{i,j\in\mathcal{N}^{C}_{c}}A_{ij}|\omega_{i}\omega_{j}|=\frac{\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}\sum_{c^{\prime}}\sum_{i,j\in\mathcal{N}^{C}_{c}}|\tilde{D}_{i}\tilde{D}_{j}|}{\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}Var\left(\sum_{i}\tilde{D}_{i}u_{i}\right)}<\infty$

Hence, Theorem 1 yields
$(\hat{\theta}-\theta)/\sigma_{n}\xrightarrow{d}N(0,1)$.

To prove consistent variance estimation, it suffices to show
$(\hat{\sigma}_{n}^{2}-\sigma_{n}^{2})/\sigma_{n}^{2}=o_{P}(1)$.

$\displaystyle\hat{\sigma}_{n}^{2}=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}u_{i}\omega_{j}u_{j}-2\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}^{2}\omega_{j}u_{j}\right)(\hat{\theta}-\theta)+\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}^{2}\omega_{j}^{2}\right)(\hat{\theta}-\theta)^{2}$

By Theorem 1,
$\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}u_{i}\omega_{j}u_{j}-\sigma_{n}^{2}\right)/\sigma_{n}^{2}=o_{P}(1)$.
Since
$(\hat{\theta}-\theta)^{2}/\sigma_{n}^{2}\xrightarrow{d}Z^{2}=\chi^{2}_{1}$,

$\displaystyle\frac{\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}^{2}\omega_{j}^{2}\right)(\hat{\theta}-\theta)^{2}}{\sigma_{n}^{2}}$
$\displaystyle=\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}^{2}\omega_{j}^{2}\right)O_{P}(1)$
$\displaystyle\frac{\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{D}_{i}^{2}\tilde{D}_{j}^{2}\right)}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{4}}$
$\displaystyle\leq\frac{\left(\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{D}_{j}^{2}\right)}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{3}}\frac{\sum_{i}\tilde{D}_{i}^{2}}{\sum_{i}\tilde{D}_{i}^{2}}\leq\frac{\left(\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{D}_{j}^{2}\right)}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}O(1)$
$\displaystyle\leq\left(\frac{\operatornamewithlimits{max}_{g}\left(\sum_{j\in\mathcal{N}^{G}_{g}}|\tilde{D}_{j}|\right)^{2}}{\sum_{g^{\prime}}\left(\sum_{j\in\mathcal{N}^{G}_{g^{\prime}}}|\tilde{D}_{j}|\right)^{2}}+\frac{\operatornamewithlimits{max}_{h}\left(\sum_{j\in\mathcal{N}^{H}_{h}}|\tilde{D}_{j}|\right)^{2}}{\sum_{h^{\prime}}\left(\sum_{j\in\mathcal{N}^{H}_{h^{\prime}}}|\tilde{D}_{j}|\right)^{2}}\right)O(1)=o(1)$

Convergence occurs due to Assumption 3.2, so
$\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}^{2}\omega_{j}^{2}\right)(\hat{\theta}-\theta)^{2}/\sigma_{n}^{2}=o_{P}(1)$.
Finally,

$\displaystyle\frac{\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}^{2}\omega_{j}u_{j}\right)(\hat{\theta}-\theta)}{\sigma_{n}^{2}}$
$\displaystyle=\frac{\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{D}_{i}^{2}\tilde{D}_{j}u_{j}\right)}{\sigma_{n}\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{3}}O_{P}(1)$

Applying Markov and Minkowski inequalities,

$\displaystyle P$
$\displaystyle\left(\frac{\left|\sum_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{D}_{i}^{2}\tilde{D}_{j}u_{j}\right|}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{3}\sigma_{n}}>\epsilon\right)\leq\frac{1}{\epsilon}\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{3}\sigma_{n}}E\left[\left|\sum_{i}\sum_{j\in\mathcal{N}_{i}}\tilde{D}_{i}^{2}\tilde{D}_{j}u_{j}\right|\right]$
$\displaystyle\leq\frac{1}{\epsilon}\frac{1}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{3}\sigma_{n}}\sum_{i}\sum_{j\in\mathcal{N}_{i}}E|\tilde{D}_{i}^{2}\tilde{D}_{j}u_{j}|\leq\frac{1}{\epsilon}\frac{\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}E|\tilde{D}_{j}u_{j}|}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}\sigma_{n}}\frac{\sum_{i}\tilde{D}_{i}^{2}}{\sum_{i}\tilde{D}_{i}^{2}}=o(1)$

Convergence occurs because

$\frac{\operatornamewithlimits{max}_{i}\left(\sum_{j\in\mathcal{N}_{i}}E|\tilde{D}_{j}u_{j}|\right)^{2}}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}\leq\frac{\operatornamewithlimits{max}_{j}E|u_{j}|^{2}\operatornamewithlimits{max}_{i}\left(\sum_{j\in\mathcal{N}_{i}}|\tilde{D}_{j}|\right)^{2}}{\left(\sum_{i}\tilde{D}_{i}^{2}\right)^{2}}=o(1)$

∎

For Proposition 2, I first prove a consistency result.

######  Lemma 9.

Under Assumption 1, 2.1 and 2.2, and $E[W_{i}]=0$ $\forall i$,
$||(1/(\sum_{i}\omega_{i}))\sum_{i}\omega_{i}(W_{i}W_{i}^{\prime}-E[W_{i}W_{i}^{\prime}])||\xrightarrow{p}0$.

###### Proof.

It suffices to show convergence elementwise. Let $X_{i}$ and $Y_{i}$ denote
scalar components of $W_{i}$, i.e., $X_{i}=W_{im},Y_{i}=W_{ip}$, where
$m,p\in\\{1,2,\cdots,K\\}$. By Chebyshev’s inequality, and
$\operatornamewithlimits{max}_{m,k}E[W_{mk}^{4}]<K_{0}$,

$\displaystyle P$
$\displaystyle\left(\frac{1}{\sum_{i}\omega_{i}}\sum_{i}\omega_{i}(X_{i}Y_{i}-E[X_{i}Y_{i}])>\epsilon\right)$
$\displaystyle\leq\frac{1}{\epsilon^{2}}\frac{1}{\left(\sum_{i}\omega_{i}\right)^{2}}E\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}(X_{i}Y_{i}-E[X_{i}Y_{i}])(X_{j}Y_{j}-E[X_{j}Y_{j}])\right)\leq\frac{K_{0}}{\epsilon^{2}\left(\sum_{j}\omega_{j}\right)^{2}}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j}$

Hence, it suffices to show
$(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\omega_{i}\omega_{j})/\left(\sum_{j}\omega_{j}\right)^{2}=o(1)$.
This follows from a similar argument as Lemma 8. ∎

###### Proof of Proposition 2.

$E[u_{i}|X_{i}]=0$ implies $E[X_{i}u_{i}]=0$ by law of iterated expectations.
Since $E[u_{i}^{4}|X_{i}]\leq K_{0}$,
$E[u_{i}^{4}X_{ik}^{4}]=E[E[u_{i}^{4}|X_{i}]X_{ik}^{4}]\leq
K_{0}E[X_{ik}^{4}]\leq K_{0}^{2}$ is bounded. By Theorem 1,
$Q^{-1/2}_{n}\sum_{i=1}^{n}X_{i}u_{i}\xrightarrow{d}N(0,I_{K})$.

To complete the normality result, I show that
$S_{n}^{-1}\hat{S}_{n}\xrightarrow{p}I_{K}$, which is the same as showing that
$||S_{n}^{-1}(\hat{S}_{n}-S_{n})||\xrightarrow{p}0$. By applying Lemma 9 with
$\omega_{i}=1$,
$(1/n)(\hat{S}_{n}-S_{n})=(1/n)\sum_{i}(X_{i}X_{i}^{\prime}-E[X_{i}X_{i}^{\prime}])=o_{P}(1)$.
Hence, it suffices that $(S_{n}/n)^{-1}$ has bounded eigenvalues, i.e.,
$\lambda_{\operatornamewithlimits{min}}(S_{n}/n)\geq K_{1}>0$, which is true
by Assumption 4.5. Since
$\hat{\beta}-\beta=\hat{S}_{n}^{-1}\sum_{i}X_{i}u_{i}$, by Slutsky’s lemma,
$Q_{n}^{-1/2}S_{n}(\hat{\beta}-\beta)\xrightarrow{d}N(0,I_{K})$.

Next, proceed to consistent variance estimation. Showing that
$||Q_{n}^{-1}\hat{Q}_{n}-I_{K}||=o_{P}(1)$ is equivalent to showing that,
$\forall l\in\mathbb{R}^{K}$,
$l^{\prime}\left(Q_{n}^{-1}(\hat{Q}_{n}-Q_{n})\right)l=o_{P}(1)$.

$\displaystyle\hat{Q}_{n}$
$\displaystyle:=\sum_{i}\sum_{j\in\mathcal{N}_{i}}\hat{u}_{i}\hat{u}_{j}X_{i}X_{j}^{\prime}=\sum_{i}\sum_{j\in\mathcal{N}_{i}}(u_{i}-X_{i}^{\prime}(\hat{\beta}-\beta))(u_{j}-X_{j}^{\prime}(\hat{\beta}-\beta))X_{i}X_{j}^{\prime}$
$\displaystyle=\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}u_{j}X_{i}X_{j}^{\prime}-2\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}X_{j}^{\prime}(\hat{\beta}-\beta)X_{i}X_{j}^{\prime}\right)+\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}X_{i}^{\prime}(\hat{\beta}-\beta)X_{j}^{\prime}(\hat{\beta}-\beta)X_{i}X_{j}^{\prime}\right)$

By Theorem 1,
$l^{\prime}Q_{n}^{-1}(\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}u_{j}X_{i}X_{j}^{\prime}-Q_{n})l=o_{P}(1)$.
Hence, it remains to show:

$\displaystyle\left|\left|Q_{n}^{-1}\left[-2\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}X_{j}^{\prime}(\hat{\beta}-\beta)X_{i}X_{j}^{\prime}\right)+\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}X_{i}^{\prime}(\hat{\beta}-\beta)X_{j}^{\prime}(\hat{\beta}-\beta)X_{i}X_{j}^{\prime}\right)\right]\right|\right|=o_{P}(1)$

Observe that
$X_{i}^{\prime}(\hat{\beta}-\beta)=\left(X_{i}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)\left(Q_{n}^{-1/2}S_{n}(\hat{\beta}-\beta)\right)=\left(X_{i}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)(Z_{K}+1_{K}o_{P}(1))$,
where $1_{K}$ is a $K$-vector of ones. Hence, addressing the second term,

$\displaystyle
X_{i}^{\prime}(\hat{\beta}-\beta)X_{j}^{\prime}(\hat{\beta}-\beta)$
$\displaystyle=\left(X_{i}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)(Z_{K}+1_{K}o_{P}(1))(Z_{K}+1_{K}o_{P}(1))^{\prime}\left(X_{j}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)^{\prime}$
$\displaystyle=\left(X_{i}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)(I_{K}O_{P}(1)+o_{P}(1))\left(X_{j}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)^{\prime}$
$\displaystyle=X_{i}^{\prime}S_{n}^{-1}Q_{n}S_{n}^{-1}X_{j}O_{P}(1)$

This implies

$\displaystyle Q_{n}^{-1}$
$\displaystyle\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}X_{i}^{\prime}(\hat{\beta}-\beta)X_{j}^{\prime}(\hat{\beta}-\beta)X_{i}X_{j}^{\prime}\right)=Q_{n}^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\left(X_{i}^{\prime}S_{n}^{-1}Q_{n}S_{n}^{-1}X_{j}\right)X_{i}X_{j}^{\prime}\right)O_{P}(1)$
$\displaystyle=\frac{1}{n^{2}}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\left(X_{i}^{\prime}\left(\frac{1}{n}S_{n}\right)^{-1}\left(\frac{1}{\lambda_{n}}Q_{n}\right)\left(\frac{1}{n}S_{n}\right)^{-1}X_{j}\right)X_{i}X_{j}^{\prime}\right)O_{P}(1)$

The eigenvalues of $(Q_{n}/\lambda_{n})$ are bounded. To see this, it suffices
to show that there exists $K_{0}<\infty$ such that
$\lambda_{\operatornamewithlimits{max}}(Q_{n})/\lambda_{n}\leq K_{0}$. Due to
finite moments, $Q_{n}:=Var(\sum_{i}X_{i})\leq K_{0}1_{K\times
K}\sum_{c}(N_{c}^{C})^{2}$. Since $(\sum_{c}(N_{c}^{C})^{2})/\lambda_{n}\leq
K_{0}$ by Assumption 4, $\lambda_{n}K_{0}\geq\sum_{c}(N_{c}^{C})^{2}$, which
implies $\lambda_{n}\geq(\sum_{c}(N_{c}^{C})^{2})/K_{0}$. Hence,

$\displaystyle\frac{\lambda_{\operatornamewithlimits{max}}(Q_{n})}{\lambda_{n}}\leq\frac{\sum_{c}(N_{c}^{C})^{2}K_{0}}{\sum_{c}(N_{c}^{C})^{2}\frac{1}{K_{0}}}=K_{0}^{2}$

Recall that $(S_{n}/n)^{-1}$ has bounded eigenvalues. The proof of Theorem 1
also showed that $(Q_{n}/\lambda_{n})^{-1}$ has bounded eigenvalues. By using
Markov and Minkowski inequalities,

$\displaystyle P$
$\displaystyle\left(\frac{1}{n^{2}}\left|l^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\left(X_{i}^{\prime}\left(\frac{1}{n}S_{n}\right)^{-1}\left(\frac{1}{\lambda_{n}}Q_{n}\right)\left(\frac{1}{n}S_{n}\right)^{-1}X_{j}\right)X_{i}X_{j}^{\prime}\right)l\right|>\epsilon\right)$
$\displaystyle\leq\frac{1}{n^{2}\epsilon}E\left[\left|l^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\left(\sum_{i}\sum_{j\in\mathcal{N}_{i}}\left(X_{i}^{\prime}\left(\frac{1}{n}S_{n}\right)^{-1}\left(\frac{1}{\lambda_{n}}Q_{n}\right)\left(\frac{1}{n}S_{n}\right)^{-1}X_{j}\right)X_{i}X_{j}^{\prime}\right)l\right|\right]$
$\displaystyle\leq\frac{1}{n^{2}\epsilon}\sum_{i}N_{i}\operatornamewithlimits{max}_{m,k}E[X_{mk}^{4}]K_{0}\leq\frac{\operatornamewithlimits{max}_{i}N_{i}}{n}\frac{n}{n}K_{0}\rightarrow
0$

where $K_{0}\in\mathbb{R}$ is an arbitrary (finite) constant. Convergence
occurs due to Assumption 4.2, which implies
$\operatornamewithlimits{max}_{i}N_{i}/n\rightarrow 0$. This occurs due to the
result that
$\operatornamewithlimits{max}_{i}\sum_{j\in\mathcal{N}_{i}}|\omega_{j}|/\left(\sum_{j}\omega_{j}\right)=o(1)$
in the proof of Lemma 8, using $\omega_{i}=1$.

Going back to the first term,

$\displaystyle Q_{n}^{-1}$
$\displaystyle\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}X_{j}^{\prime}(\hat{\beta}-\beta)X_{i}X_{j}^{\prime}=Q_{n}^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}\left(X_{i}^{\prime}S_{n}^{-1}Q_{n}^{1/2}\right)(Z_{K}+1_{K}o_{P}(1))X_{i}X_{j}^{\prime}$
$\displaystyle=\frac{1}{n\sqrt{\lambda_{n}}}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}\left(X_{i}^{\prime}\left(\frac{1}{n}S_{n}\right)^{-1}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{1/2}\right)X_{i}X_{j}^{\prime}O_{P}(1)$

By using Markov and Minkowski inequalities,

$\displaystyle P$
$\displaystyle\left(\frac{1}{n\sqrt{\lambda_{n}}}\left|l^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}\left(X_{i}^{\prime}\left(\frac{1}{n}S_{n}\right)^{-1}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{1/2}\right)X_{i}X_{j}^{\prime}l\right|>\epsilon\right)$
$\displaystyle\leq\frac{1}{n\sqrt{\lambda_{n}}\epsilon}E\left[\left|l^{\prime}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{-1}\sum_{i}\sum_{j\in\mathcal{N}_{i}}u_{i}\left(X_{i}^{\prime}\left(\frac{1}{n}S_{n}\right)^{-1}\left(\frac{1}{\lambda_{n}}Q_{n}\right)^{1/2}\right)X_{i}X_{j}^{\prime}l\right|\right]$
$\displaystyle\leq\frac{1}{n\sqrt{\lambda_{n}}\epsilon}\sum_{i}\sum_{j\in\mathcal{N}_{i}}\operatornamewithlimits{max}_{m_{1},m_{2},k}E\left[\left|X_{m_{1}k}u_{m_{1}}X_{m_{2}}^{2}\right|\right]K_{0}$
$\displaystyle\leq\frac{1}{n\sqrt{\lambda_{n}}\epsilon}\sum_{i}N_{i}\operatornamewithlimits{max}_{m_{1},m_{2},k}E\left[|X_{m_{1}k}u_{m_{1}}|^{2}\right]^{1/2}E\left[|X_{m_{2}}^{2}|^{2}\right]^{1/2}K_{0}$
$\displaystyle\leq\frac{\operatornamewithlimits{max}_{i}N_{i}}{\sqrt{\lambda_{n}}}\frac{1}{\epsilon}\operatornamewithlimits{max}_{m_{1},m_{2},k}E[X_{m_{1}k}^{2}u_{m_{1}}^{2}]^{1/2}E[X_{m_{2}}^{4}]^{1/2}K_{0}=o(1)$

The penultimate inequality occurs due to Holder’s inequality. Observe that
$\operatornamewithlimits{max}_{i}N_{i}/\sqrt{\lambda_{n}}=o(1)$ if and only if
$\operatornamewithlimits{max}_{c}(N^{C}_{c})^{2}/\lambda_{n}=o(1)$, which is
given by Assumption 4.2. Convergence in the last step occurs because
$\operatornamewithlimits{max}_{i}N_{i}/\sqrt{\lambda_{n}}=o(1)$, and finite
moments.

Hence, it has been shown that $Q_{n}^{-1}\hat{Q}_{n}\xrightarrow{p}I_{K}$.
Then,
$[S_{n}^{-1}Q_{n}S_{n}^{-1}]^{-1}[\hat{S}_{n}^{-1}\hat{Q}_{n}\hat{S}_{n}^{-1}]\xrightarrow{p}I_{K}$
by the continuous mapping theorem. ∎

###### Proof of Corollary 1.

By Proposition 2,
$(\hat{\beta}_{1}-\beta_{1})/[V(\hat{\beta})]^{1/2}_{11}\xrightarrow{d}N(0,1)$.
Since $\hat{\theta}=\hat{\beta}_{1}$,
$[V(\hat{\beta})]_{11}=V(\hat{\beta}_{1})=V(\hat{\theta})=\sigma_{n}^{2}$.
Hence, $(\hat{\theta}-\theta)/\sigma_{n}\xrightarrow{d}N(0,1)$.

A further implication of Proposition 2 is that
$[\hat{V}(\hat{\beta})]_{11}/[V(\hat{\beta})]_{11}\xrightarrow{p}1$. Using
theorem 3 of Ding (2021), the Liang-Zeger estimators (Liang and Zeger, 1986)
are numerically equivalent regardless of whether the long regression or the
residualized regression were used. Since the CGM estimator is a function of
the Liang-Zeger estimators,
$\hat{\sigma}_{n}^{2}=[\hat{V}(\hat{\beta})]_{11}$. Hence,
$\hat{\sigma}_{n}^{2}/\sigma_{n}^{2}\xrightarrow{p}1$. ∎

## References

  * Cameron et al. (2011) Cameron, A. C., J. B. Gelbach, and D. L. Miller (2011): “Robust inference with multiway clustering,” _Journal of Business & Economic Statistics_, 29, 238–249.
  * Chen and Shao (2004) Chen, L. H. and Q.-M. Shao (2004): “Normal approximation under local dependence,” _The Annals of Probability_ , 32, 1985–2028.
  * Davezies et al. (2021) Davezies, L., X. D’Haultfœuille, and Y. Guyonvarch (2021): “Empirical process results for exchangeable arrays,” _The Annals of Statistics_ , 49, 845–862.
  * Ding (2021) Ding, P. (2021): “The Frisch–Waugh–Lovell theorem for standard errors,” _Statistics & Probability Letters_, 168, 108945.
  * Djogbenou et al. (2019) Djogbenou, A. A., J. G. MacKinnon, and M. Ø. Nielsen (2019): “Asymptotic theory and wild bootstrap inference with clustered errors,” _Journal of Econometrics_ , 212, 393–412.
  * Dube et al. (2010) Dube, A., T. W. Lester, and M. Reich (2010): “Minimum wage effects across state borders: Estimates using contiguous counties,” _The review of economics and statistics_ , 92, 945–964.
  * Hansen and Lee (2019) Hansen, B. E. and S. Lee (2019): “Asymptotic theory for clustered samples,” _Journal of econometrics_ , 210, 268–290.
  * Kallenberg (2005) Kallenberg, O. (2005): _Probabilistic symmetries and invariance principles_ , vol. 9, Springer.
  * Liang and Zeger (1986) Liang, K.-Y. and S. L. Zeger (1986): “Longitudinal data analysis using generalized linear models,” _Biometrika_ , 73, 13–22.
  * MacKinnon et al. (2021) MacKinnon, J. G., M. Ø. Nielsen, and M. D. Webb (2021): “Wild bootstrap and asymptotic inference with multiway clustering,” _Journal of Business & Economic Statistics_, 39, 505–519.
  * Menzel (2021) Menzel, K. (2021): “Bootstrap With Cluster-Dependence in Two or More Dimensions,” _Econometrica_ , 89, 2143–2188.
  * Michalopoulos and Papaioannou (2013) Michalopoulos, S. and E. Papaioannou (2013): “Pre-colonial ethnic institutions and contemporary African development,” _Econometrica_ , 81, 113–152.
  * Nunn and Wantchekon (2011) Nunn, N. and L. Wantchekon (2011): “The slave trade and the origins of mistrust in Africa,” _American Economic Review_ , 101, 3221–52.
  * Ross (2011) Ross, N. (2011): “Fundamentals of Stein’s method,” _Probability Surveys_ , 8, 210–293.