EconBase
← Back to paper

Experimental Design for Matching

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.

125,279 characters · 20 sections · 71 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Experimental Design for Matching

\RUNAUTHOR{Wang}

\RUNTITLE{Test of Matching}

\TITLE{Experimental Design for Matching}

\ARTICLEAUTHORS{ \AUTHOR{Chonghuan Wang} \AFF{Naveen Jindal School of Management, University of Texas at Dallas, Richardson, TX 75080, \EMAIL{[email removed]}} }

\ABSTRACT{Matching mechanisms play a central role in operations management across diverse fields including education, healthcare, and online platforms. However, experimentally comparing a new matching algorithm against a status quo presents some fundamental challenges due to matching interference, where assigning a unit in one matching may preclude its assignment in the other. In this work, we take a design-based perspective to study the design of randomized experiments to compare two predetermined matching plans on a finite population, without imposing outcome or behavioral models. We introduce the notation of a disagreement set, which captures the difference between the two matching plans, and show that it admits a unique decomposition into disjoint alternating paths and cycles with useful structural properties. Based on these properties, we propose the Alternating Path Randomized Design, which sequentially randomizes along these paths and cycles to effectively manage interference. Within a minimax framework, we optimize the conditional randomization probability and show that, for long paths, the optimal choice converges to $\sqrt{2}-1$, minimizing worst-case variance. We establish the unbiasedness of the Horvitz-Thompson estimator and derive a finite-population Central Limit Theorem that accommodates complex and unstable path and cycle structures as the population grows. Furthermore, we extend the design to many-to-one matchings, where capacity constraints fundamentally alter the structure of the disagreement set. Using graph-theoretic tools, including finding augmenting paths and Euler-tour decomposition on an auxiliary unbalanced directed graph, we construct feasible alternating path and cycle decompositions that allow the design and inference results to carry over. }

\KEYWORDS{Matching, A/B Testing, Randomized Experiment, Causal Inference}

Introduction

Matching is a foundational topic in operations management, with deep theoretical roots and wide-ranging real-world deployment. Across domains such as school choices and admission (e.g., chen2021information,ignacio2025stable,gong2025dynamic), organ exchange and transplantation (e.g., roth2004kidney,ashlagi2021kidney, tang2025multi), refugee settlement (e.g., delacretaz2023matching,ahani2021placement), ridesharing platforms (e.g., freund2024two), dating markets (e.g., hitsch2010matching, rios2023improving) and video game (e.g., chen2021matchmaking), matching mechanisms serve as critical infrastructure for allocating scarce resources, forming partnerships, and ultimately increasing social welfare. The resulting economic and social impact has made matching a central topic of interest to both researchers and practitioners.

A standard approach in both research and practice is to model agent behavior or forecast the matching utilities, and then design mechanisms that perform well under those models. While this modeling approach has yielded important insights, practitioners may still worry about the reliability of such assumptions when deploying mechanisms in real-world systems. This naturally raises the question of how one can empirically evaluate the performance of a newly proposed matching mechanism. Moreoever, as digitization increases the granularity of observable features used in matching, it becomes common to ask whether incorporating additional features actually improves matching quality, and whether such improvements justify the associated engineering and financial costs. These considerations motivate a complementary perspective focused on empirical testing rather than structural modeling. In this work, we study a basic but important question: How can we experimentally compare two given matching plans on the same finite population without imposing outcome or behavioral models?

More concretely, we consider a setting in which two predetermined matching plans can be applied to the same finite population: a treatment matching, generated by a new mechanism or using additional features, and a control matching, representing the status quo. We define the treatment effect as the difference in average outcomes under the treatment and control matchings. Our objective is to design valid randomized experiments that allow us to estimate and infer this treatment effect. We adopt a design-based perspective, meaning that inference relies on the known randomization scheme rather than a model of matching outcomes or agent behavior. This enables us to draw conclusions without imposing outcome or behavioral models, which is essential when matching assumptions themselves are under investigation.

The core technical challenge arises because a unit may be paired with different counterparts under the treatment and control matchings. As a result, in an experiment we can only observe the outcome of one realized match, but never the counterfactual outcome of the alternative match. This lack of individual-level counterfactuals is a familiar obstacle in standard causal inference. Beyond this, the structure of matching introduces hard feasibility constraints that couple the randomization decisions across units. In classical experimental designs, assigning one unit to treatment or control is typically independent of assignments for other units. In contrast, when deciding whether to adopt a match from the treatment or control plan, that decision may preclude other feasible matches due to capacity or compatibility constraints, and the resulting impact can propagate across the matching. We refer to this phenomenon as matching interference. Because matching pairs cannot be treated as independent experimental units, a central task in our setting is to design randomization schemes that carefully manage and quantify such interference to enable valid estimation and inference.

Main Results

We start with arguably the fundamental case of designing randomized experiments for one-to-one matchings. We introduce the notion of the disagreement set, defined as the set of matching pairs that appear in exactly one of the two matchings (treatment or control). This set fully characterizes the difference between the two matchings and therefore contains all pairs that require randomization. A key structural observation is that the disagreement set admits a unique decomposition into disjoint paths and cycles, where edges correspond to matching pairs and alternate between the treatment matching and the control matching. Along any such path or cycle, adjacent edges cannot be simultaneously realized in the experiment due to capacity constraints. This motivates our Alternating Path Randomized Design (AP design), in which randomization proceeds sequentially along each path or cycle to ensure feasibility. Intuitively, the decision to include an edge is randomized conditional on the realized selection of its predecessor along the path or cycle.

We next optimize the conditional randomization probabilities of our design using a minimax framework. We show that, in the worst case and as paths and cycles grow long, the optimal conditional probability converges to $\sqrt{2}-1$, which is strictly less than $0.5$. Intuitively, this reflects a trade-off: reducing the variance contribution of the current edge while accounting for downstream edges along the path or cycle. For inference, we establish a finite-population central limit theorem (CLT) for the Horvitz–Thompson estimator, enabling valid asymptotic inference under the design-based perspective. To our knowledge, this result is not a direct consequence of existing CLTs in the literature. The technical novelty arises from the fact that the path/cycle decomposition can exhibit highly heterogeneous and potentially unstable patterns: the number of components may be large with short lengths, small with extremely long lengths, or fluctuate between these regimes without converging to a stable structure as the population grows. Existing tools typically target either many short components or a few long ones, but not (unstable) mixtures that change with the size of the population. We address this by combining classical Lindeberg–Feller arguments with $\alpha$-mixing techniques through Bolzeno-Weierstrass theorem and the subsequence principle, yielding a CLT that accommodates such instability.

We further consider the experimental design problem for many-to-one matchings. In this setting, a new and significant challenge arises: the path/cycle decomposition of the disagreement set is generally not unique, and many possible decompositions violate capacity constraints and are therefore not admissible for randomization. A key part of our contribution is to identify valid decompositions that preserve feasibility so that our AP design and its inference guarantees remain applicable. We provide a set of sufficient conditions under which a decomposition is admissible, and show that finding such a decomposition can be reduced to two graph-theoretic problems, finding augmenting paths and Euler-tour decompositions, via the construction of an auxiliary unbalanced directed graph and a flow network. This extension from one-to-one to many-to-one matchings is structurally non-trivial, as it requires additional combinatorial feasibility considerations beyond those present in the one-to-one case.

Related Work

Matching, as one of the most common practice in today's platform operations, has attracted much attention in running experiments on platforms. The main focus of the current literature is that the matching introduces interference among experimental units (see, e.g., johari2022experimental). More specifically, individuals in the market compete with each other through a matching process. li2022interference and dhaouadi2023price develop insightful market models to study bias and variance in two-sided markets with interference. weng2024experimental carefully investigate the experimentation on one-side markets, such as the live-interaction platforms and the video games. To mitigate the bias, bright2025reducing consider a matching mechanism based on linear programming and debias the estimators via shadow prices. bajari2023experimental and masoero2024multiple propose the two-sided randomization to reduce the bias, which are also known as multiple randomization designs. wager2021experimenting take an equilibrium model into account and study the local experimentation of the two-sided market. pouget2019variance, doudchenko2020causal, harshaw2023design introduce the framework of bipartite experiments. zhan2024estimating and zhang2025debiasing investigate seller-side experiments by using the MNL model to capture the interference introduced by the two-sided market. If the market has some Markovian properties, farias2022markovian propose a novel model-free difference-in-Q estimator and chen2024experimenting showcase how the treatment locality can help to further reduce the variance. All the mentioned works try to mitigate the bias issue introduced by a fixed matching mechanism. However, they do not discuss how to test or compare two given matching mechanisms, which is the scope of this paper.

More broadly, the paper is related to the growing literature of designing experiments in complex systems with interference. One of the most well-studied interferences is network effect (see, e.g., eckles2017design,athey2018exact,li2022random,ugander2023randomized,viviano2023causal,cortez2023exploiting, candogan2024correlated,holtz2025reducing). A large stream of remedy in the recent literature is to run switchback experiment (see, bojinov2023design, hu2022switchback, xiong2024data, ni2025enhancing, chen2025efficient, ni2025reliable). In a switchback experiment, the entire platform alternates between treatment and control usually on an hourly basis. In addition, recent academic focus has shifted to link the decision-making to experimental design, e.g., dhaouadi2023price, johari2024does, li2025choosing and ni2025decision. johari2024does reveal the scenarios where the interference does not matter in terms of decision-making, although it biases the estimation. In contrast, li2025choosing identify that when testing recommendation algorithms, the interference may alter the sign of the treatment effects leading to wrong decision of deploying the new algorithms or not. ni2025decision adopts a robust optimization approach to suggest the optimal decision rule on adopting the new treatment based on the observations from switchback experiments. For experimentation in operational systems, li2023experimenting investigate the experimentation under queuing congestion, and chen2025bias study the bias of experiment in complex inventory problem. xiong2025automated recently introduces the simulation tool to tackle the challenge that the objective is typically a black-box function of the design. Also, much work centers on identifying individualized treatment rules to maximize the welfare, i.e., policy learning (kitagawa2018should, kallus2021minimax, zhang2024minimax). We refer to wager2024causal for a detailed discussion on policy learning. In addition, causal message passing by shirani2024causal can be powerful to model very general network interference in experimentation.

Preliminaries

Model

We begin with the one-to-one matching setting for clarity. Consider a population of $2N$ agents that need to be matched in pairs. Our formulation accommodates both one-sided environments (e.g., pairing two players in a one-versus-one game) and two-sided environments (e.g., matching $N$ workers to $N$ jobs). In the two-sided case, we assume for simplicity that the population consists of $N$ agents on each side. A matching is represented by a function

equation*[equation* omitted — 55 chars of source]

where $[2N]$ represents the set $\{1,2,\cdots,2N\}$, $m(i)$ denotes the agent matched to $i$, and $m(i)=0$ indicates that $i$ is unmatched. A one-to-one matching is feasible if $m(i)=j>0$ guarantees $m(j)=i$ for any agent $i$, so that matches are symmetric and an agent cannot be matched to multiple partners. We let ${\mathbb{M}}$ denote the set of all the realized match pairs:

equation[equation omitted — 76 chars of source]

where the condition $i<m(i)$ ensures that each pair appears exactly once and precludes the unmatched agents, i.e., the agent with $m(i)=0$.

figure[figure omitted — 273 chars of source]

We are given two matchings on the same set of $2N$ agents: a treatment matching ${\mathbb{M}}^t$ and a control matching ${\mathbb{M}}^c$. The treatment matching ${\mathbb{M}}^t$ can be proposed by a new matching algorithm considering new features and the control matching ${\mathbb{M}}^c$ is from the status quo. For an agent $i$, we denote its partner in ${\mathbb{M}}^t$ by $m^t(i)$ and in ${\mathbb{M}}^c$ by $m^c(i)$. Because each agent can be matched with at most one partner, an agent cannot simultaneously realize both $m^t(i)$ and $m^c(i)$, which motivates the need to compare the two matchings through an experiment rather than direct observation. To formalize the difference between ${\mathbb{M}}^t$ and ${\mathbb{M}}^c$, we define the disagreement set

equation[equation omitted — 158 chars of source]

which contains all match pairs that appear in exactly one of the two matchings. We further partition ${\triangle}{\mathbb{M}}^{(t,c)}$ into

equation*[equation* omitted — 210 chars of source]

which we refer to respectively as t-matches and c-matches. Thus, ${\triangle}{\mathbb{M}}^{(t,c)}_t$ collects the match pairs that would be realized only under the treatment matching, and ${\triangle}{\mathbb{M}}^{(t,c)}_c$ collects those realized only under the control matching. Comparing these pairs is the core object of the experimental design.

exampleFigure (ref) illustrates a toy instance with $2N=10$ agents and two one-to-one matchings. The pair $(5,10)$ appears in both ${\mathbb{M}}^t$ and ${\mathbb{M}}^c$ and is therefore excluded from the disagreement set. The resulting disagreement set is \[ {\triangle}{\mathbb{M}}^{(t,c)} = \{(1,6),(1,7),(2,6),(2,7),(3,8),(3,9),(4,8),(4,9)\}. \] Among these, the t-matches are ${\triangle}{\mathbb{M}}^{(t,c)}_t = \{(1,6),(2,7),(3,8),(4,9)\}$, and the c-matches are ${\triangle}{\mathbb{M}}^{(t,c)}_c = \{(1,7),(2,6),(3,9),(4,8)\}$.

We associate each ordered pair $(i,j)$ with a potential outcome $Y_{i,j}$, representing the outcome that would be realized if $i$ were matched with $j$. We make no structural assumptions or behavior models on $Y_{i,j}$. Under the design-based perspective (e.g., imbens2015causal), these values are treated as fixed but unknown, and randomness arises solely from the experimental design. The outcome $Y_{i,j}$ is revealed only if the match $(i,j)$ is realized in the experiment. Our estimand of interest is the average treatment effect (ATE) of the treatment matching ${\mathbb{M}}^t$ relative to the control matching ${\mathbb{M}}^c$:

equation*[equation* omitted — 179 chars of source]

If $m^t(i)\not= m^c(j)$, then at most one of $Y_{i,m^t(i)}$ and $Y_{i,m^c(i)}$ can be observed, so $\tau$ cannot be directly observed and must be estimated from an experiment. Since match pairs appearing in both ${\mathbb{M}}^t$ and ${\mathbb{M}}^c$ contribute equally to $\bar{Y}_t$ and $\bar{Y}_c$ and cancels out in $\tau$, the disagreement set ${\triangle}{\mathbb{M}}^{(t,c)}$ contains exactly the match pairs relevant for identification.

Given these definitions, an experimental design for matching consists of two elements: (i) a randomization scheme that selects a feasible set of realized matches ${\mathbb{M}}^{\text{EXP}}$ in the presence of matching interference, and (ii) an estimator $\hat{\tau}$ that uses the realized outcomes and the randomization process to estimate $\tau$.

Before introducing our proposed design, it is instructive to consider a simple baseline that highlights the limitations of directly switching between the two matchings.

remark[The Naive Experimental Design] A natural baseline design selects either the treatment matching ${\mathbb{M}}^t$ or the control matching ${\mathbb{M}}^c$ with equal probability $1/2$. Let ${\mathbb{M}}^{\mathrm{Naive}}$ denote the realized matching and define \[ \bar{Y}^{\mathrm{Naive}} = \frac{1}{N} \sum_{(i,j)\in {\mathbb{M}}^{\mathrm{Naive}}} Y_{i,j}. \] Under this design, $\bar{Y}^{\mathrm{Naive}} = \bar{Y}_t$ with probability $1/2$ and $\bar{Y}^{\mathrm{Naive}} = \bar{Y}_c$ with probability $1/2$. A natural inverse-propensity estimator of $\tau$ is \[ \hat{\tau}_{\mathrm{Naive}} = 2\,\bar{Y}^{\mathrm{Naive}}\left( \mathbb{I}\{{\mathbb{M}}^{\mathrm{Naive}}={\mathbb{M}}^t\} - \mathbb{I}\{{\mathbb{M}}^{\mathrm{Naive}}={\mathbb{M}}^c\} \right). \] \begin{proposition}[Bias and Variance of the Naive Design] Under the naive design, the estimator is unbiased $\mathbb{E}(\hat{\tau}_{\text{Naive}})=\tau$ with the variance $\mathrm{Var}(\hat{\tau}_{\text{Naive}})=(\bar{Y}_t+\bar{Y}_c)^2$. \end{proposition} Although unbiased, the estimator exhibits two drawbacks. First, its variance does not vanish as $N\to\infty$ in general. For example, in a two-sided market with $N$ workers and $N$ jobs, suppose each realized match yields an outcome $Y_{i,j}\ge c_0>0$. Then for any two complete matchings, $\mathrm{Var}\!\left(\hat{\tau}_{\mathrm{Naive}}\right)\ge c_0^2$ for any $N>0$, so increasing the experimental population does not improve precision. Second, $\widehat{\tau}^{\mathrm{Naive}}$ takes only two possible values and therefore does not satisfy a CLT that would enable asymptotic inference.

Alternating Path Randomized Design

We now introduce our proposed Alternating Path Randomized Design (AP Design) for comparing two predetermined matchings in the presence of matching interference. The key structural insight underlying AP design is that the disagreement set ${\triangle}{\mathbb{M}}^{(t,c)}$ can be decomposed into several alternating paths or cycles, which are defined as follows.

definition[Alternating Path/Cycle] An alternating path (or alternating cycle) in ${\triangle}{\mathbb{M}}^{(t,c)}$ is an ordered sequence of agents $\{v_1,v_2\cdots,v_{k+1}\}$ such that: \begin{itemize} • For every $i \in \{1,\dots,k\}$, the edge $(v_i, v_{i+1})$ belongs to ${\triangle}{\mathbb{M}}^{(t,c)}$. • The edges alternate between ${\triangle}{\mathbb{M}}^{(t,c)}_t$ and ${\triangle}{\mathbb{M}}^{(t,c)}_c$, i.e., \[(v_i, v_{i+1}) \in {\triangle}{\mathbb{M}}^{(t,c)}_t \;\Rightarrow\; (v_{i+1}, v_{i+2}) \in {\triangle}{\mathbb{M}}^{(t,c)}_c,\] and \[(v_i, v_{i+1}) \in {\triangle}{\mathbb{M}}^{(t,c)}_c \;\Rightarrow\; (v_{i+1}, v_{i+2}) \in {\triangle}{\mathbb{M}}^{(t,c)}_t,\] for all $i \in \{1,\dots,k-1\}$. • If $v_1 \neq v_{k+1}$, then both $v_1$ and \(v_{k+1}\) have degree 1 in ${\triangle}{\mathbb{M}}^{(t,c)}$, and the sequence forms an open path. If $v_1 = v_{k+1}$, the sequence forms a cycle. \end{itemize}

\addtocounter{example}{-1}

example[Continued] In Figure (ref) (third panel), the disagreement decomposes into two alternating cycles, $\{1,6,2,7,1\}$ and $\{3,8,4,9,3\}$. One can also visualize alternation by noting that adjacent edges in each cycle are with different colors, which correspond to different sets , ${\triangle}{\mathbb{M}}^{(t,c)}_t$ and ${\triangle}{\mathbb{M}}^{(t,c)}_c$. In this instance, the decomposition into these two cycles seems to be unique.

The endpoint condition in Definition (ref) ensures that alternating paths cannot be extended further within ${\triangle}{\mathbb{M}}^{(t,c)}$. Without this condition, any connected subsequence of vertices along a maximal alternating path or cycle would also qualify as an alternating path, which would prevent a unique decomposition. With the condition in place, we obtain the following uniqueness result.

proposition[Uniqueness of Decomposition] Let ${\mathbb{M}}^t$ and ${\mathbb{M}}^c$ be two one-to-one matchings. Then the disagreement set ${\triangle}{\mathbb{M}}^{(t,c)}$ admits a unique decomposition into a disjoint collection of alternating paths and alternating cycles. Furthermore, every connected component of the graph induced by ${\triangle}{\mathbb{M}}^{(t,c)}$ is either an alternating path or an alternating cycle.

The proposition follows from two simple observations. First, the connected components of any edge set are uniquely determined. Second, since ${\mathbb{M}}^t$ and ${\mathbb{M}}^c$ are one-to-one matchings, each agent is incident to at most one edge in ${\triangle}{\mathbb{M}}^{(t,c)}_t$ and at most one edge in ${\triangle}{\mathbb{M}}^{(t,c)}_c$. Consequently, every agent in any connected component of the graph induced by ${\triangle}{\mathbb{M}}^{(t,c)}$ has maximum degree two, and its edges necessarily alternate between ${\triangle}{\mathbb{M}}^{(t,c)}_t$ and ${\triangle}{\mathbb{M}}^{(t,c)}_c$, forming either an alternating path or an alternating cycle. We let ${\triangle\mathcal{P}^{(t,c)}}:=\{{\mathcal{P}}_1,\cdots, {\mathcal{P}}_m\}$ denote the unique collection of alternating paths and cycles obtained from decomposing ${\triangle}{\mathbb{M}}^{(t,c)}$. With mild abuse of notation, we may write $(i,j)\in{\triangle\mathcal{P}^{(t,c)}}$ to indicate that $i$ and $j$ are adjacent vertices within the one component pf ${\triangle\mathcal{P}^{(t,c)}}$. Additionally, a useful property of alternating cycles we want to highlight is that they always have even length.

proposition[Even-Length Cycles] Let ${\mathbb{M}}^t$ and ${\mathbb{M}}^c$ be two one-to-one matchings. For any ${\mathcal{P}}_i=\{v_{i,1},\cdots,v_{i,k(i)+1}\}\in {\triangle\mathcal{P}^{(t,c)}}$ that is a cycle, $k(i)$ must be an even integer with $k(i)\ge 4$.

If $k(i)$ were odd, then by the definition of alternation, the first and last edges $(v_{i,1}, v_{i,2})$ and $(v_{i,k(i)}, v_{i,k(i)+1})$ would both belong to ${\triangle}{\mathbb{M}}^{(t,c)}_t$ or both belong to ${\triangle}{\mathbb{M}}^{(t,c)}_c$. Since $v_{i,1} = v_{i,k(i)+1}$ in a cycle, this would imply that $v_{i,1}$ is incident to two distinct edges from ${\triangle}{\mathbb{M}}^{(t,c)}_t$ (or from ${\triangle}{\mathbb{M}}^{(t,c)}_c$), contradicting the one-to-one matching property that each agent can be matched with at most one partner in ${\mathbb{M}}^t$ and in ${\mathbb{M}}^c$. Hence, $k(i)$ must be even. The case $k(i)=2$ is impossible because it would require two distinct agents simultaneously matched to each other in both matchings, contradicting the definition of a disagreement edge. Therefore, $k(i)\ge 4$.

We first note that matching interference is fully captured within each alternating component. Adjacent edges on an alternating path cannot be realized simultaneously, and no such interference arises between different components in ${\triangle\mathcal{P}^{(t,c)}}$. Hence, randomization can be performed independently across cycles and paths, and we describe AP design on a generic component ${\mathcal{P}}_i=\{v_{i,1},\cdots,v_{i,k(i)+1}\} \in {\triangle\mathcal{P}^{(t,c)}}$. Foe each edge $(v_{i,j},v_{i,j+1})$, we define a binary indicator $W_{i,j}\in\{0,1\}$, where $W_{i,j}=1$ indicates that the match is realized in the experiment. The randomization proceeds sequentially along the path/cycle. More specifically, the randomization of $W_{i,j}$ is conditional on $W_{i,j-1}$. If $W_{i,j-1}=1$, $W_{i,j}$ has to be zero, reflecting that adjacent matches cannot both be realized due to capacity constraints. If $W_{i,j-1}=0$, $W_{i,j}$ can be randomized as 1 with probability $p$. More formally, we define our AP design in the following definition.

definition[Alternating Path Randomized Design (AP Design)] Let ${\triangle\mathcal{P}^{(t,c)}} = \{{\mathcal{P}}_1,\cdots,{\mathcal{P}}_i,\cdots,{\mathcal{P}}_m\}$ be the alternating path/cycle decomposition of ${\triangle}{\mathbb{M}}^{(t,c)}$. AP design generates a binary selection vector ${\bm{W}} = \{W_{i,j}\}_{i\in [m],j\in[k(i)]}$ with the following rules: (i) Independence across paths and cycles. Randomization is performed independently across paths/cycles ${\mathcal{P}}_i$. (ii) Path case. If ${\mathcal{P}}_i = (v_{i,1},\dots,v_{i,k(i)+1})$ is an open path, then \begin{equation} \mathbb{P}(W_{i,1}=1)=\frac{p}{1+p}. \end{equation} and for $j=2,\dots,k(i)$, \begin{equation} \mathbb{P}(W_{i,j}=1 \mid W_{i,j-1}=0) = p, \qquad \mathbb{P}(W_{i,j}=1 \mid W_{i,j-1}=1) = 0. \end{equation} (iii) Cycle case. If ${\mathcal{P}}_i = (v_{i,1},\dots,v_{i,k(i)+1})$ is a cycle, for $1\le j\le k(i)-1$, $W_{i,j}$ follows the same rules in \Eqref{eq:first-edge} and \Eqref{eq:conditional-randomization}. The final edge (i.e., $j=k(i)$) is set deterministically: \[ \mathbb{P}(W_{i,k(i)} = 1 \mid W_{i,1}=0, \; W_{i,k(i)-1}=0) = 1, \] and $W_{i,k(i)}=0$ otherwise.

There are several remarks regarding this randomization procedure. First, although $W_{i,k(i)}$ is set deterministically conditional on $W_{i,1}$ and $W_{i,k(i)-1}$ in the cycle case, it remains a random variable because $W_{i,1}$ and $W_{i,k(i)-1}$ are themselves random. Second, we do not introduce an additional randomization step for $W_{i,k(i)}$ because this final decision does not affect the feasibility of any other edges, and setting it deterministically allows us to include as many matches as possible. Similar conditional randomization ideas have appeared, for example, in ni2025enhancing, though in a different context for improving switchback experiments. For notational convenience, we use ${\bm{W}}_i$ to denote the vector of indicators on component ${\mathcal{P}}_i$ and ${\bm{W}}$ for the full assignment on the disagreement set ${\triangle}{\mathbb{M}}^{(t,c)}$. The realized selection ${\bm{W}}$ induces the experimental matching ${\mathbb{M}}^{\text{AP}}$. Match pairs outside ${\triangle}{\mathbb{M}}^{(t,c)}$ can be implemented arbitrarily, as they do not affect the estimation of the average treatment effect; however, implementing them can be useful if one is interested in estimating $\bar{Y}_t$ or $\bar{Y}_c$ themselves.

For the estimation, given any realization of the assignment ${\bm{w}}$, we use the Horvitz-Thompson estimator for the average treatment effect:

equation[equation omitted — 518 chars of source]

Note that if ${\mathbb{I}}\{w_{i,j}=1\}=1$, $Y_{v_{i,j},v_{i,j+1}}$ is observed and if ${\mathbb{I}}\{w_{i,j}=1\}=0$, though $Y_{v_{i,j},v_{i,j+1}}$ is unobservable, ${\mathbb{I}}\{w_{i,j}=1\} \cdot Y_{v_{i,j},v_{i,j+1}}$ equals to zero anyway. Therefore, $\hat{\tau}_{\text{AP}}$ is well-defined under the AP design. For notational convenience, we write

equation*[equation* omitted — 399 chars of source]

where $\Gamma_i$ represents the corresponding fixed non-random path-level treatment effect and $\hat{\Gamma}_i$ is its estimation. In the following lemma, we compute the non-conditioned probability of each $W_{i,j}=1$, $\mathbb{P}(W_{i,j}=1)$ for each edge undet the AP design.

lemmaFor any ${\mathcal{P}}_i=(v_{i,1},\cdots,v_{i,k(i)+1})\in{\triangle\mathcal{P}^{(t,c)}}$, the unconditional probability of $W_{i,j}=1$ under AP design satisfies, \begin{equation} \mathbb{P}(W_{i,j}=1) = \begin{cases} \frac{1-(-p)^{j-1}}{(1+p)^2}, & if {\mathcal{P}}_i is a cycle and j=k(i), \\ \frac{p}{p+1}, & otherwise. \end{cases} \end{equation}

The lemma shows that only the final edge of an alternating cycle has a different unconditional selection probability. This reflects the fact that, in a cycle, feasibility constraints require the last edge to depend on both its predecessor and the first edge. Moreover, the choice of the initial probability $\mathbb{P}(W_{i,1}=1)={p}/{(1+p)}$ in \Eqref{eq:first-edge} yields the clean closed-form expression in \Eqref{eq:uncondition-prob}.

Since $\mathbb{P}(W_{i,j}=1)$ is strictly positive whenever $p>0$, every edge in ${\triangle}{\mathbb{M}}^{(t,c)}$ has positive selection probability. Under the AP design, the Horvitz-Thompson estimator is therefore unbiased in the standard design-based sense.

proposition[Unbiasedness of the Horvitz-Thompson estimator] For any $p\in(0,1)$, the Horvitz-Thompson estimator $\hat{\tau}_{\text{AP}}$ is unbiased for ATE, i.e., $\mathbb{E}[\hat{\tau}_{\text{AP}}]=\tau$.

Variance Analysis and Optimal Design

We now analyze the variance of the Horvitz--Thompson estimator under the AP design and study how to choose the design parameter $p$ (i.e., the conditional probability) to improve efficiency. The first key observation is that randomization is independent across alternating paths and cycles, which yields a path-wise variance decomposition.

proposition[Variance Decomposition] Under the AP design, the variance of $\hat{\tau}_{\text{AP}}$ can be decomposed into the sum of the variance of the path-wise estimator $\hat{\Gamma}_i$, i.e., \begin{equation*} \mathrm{Var} (\hat{\tau}_{AP})=\frac{1}{N^2} \sum_{{\mathcal{P}}_i\in {\triangle\mathcal{P}^{(t,c)}}} \mathrm{Var} (\hat{\Gamma}_i). \end{equation*}

With Proposition (ref), we can zoom in to the variance of one single alternating path/cycle.

theoremFor ${\mathcal{P}}_i=\{v_{i,1},\cdots,v_{i,k(i)+1}\}\in{\triangle\mathcal{P}^{(t,c)}}$, if ${\mathcal{P}}_i$ is a path rather than a cycle (i.e., $v_{i,k(i)+1}\not= v_{i,1}$), we have \begin{equation} \mathrm{Var}(\hat{\Gamma}_i)= \frac{1}{p} \sum_{j=1}^{k(i)} Y_{v_{i,j},v_{i,j+1}}^2 + 2 \sum_{1\le j<q \le k(i)} p^{q-j-1} Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}. \end{equation} If ${\mathcal{P}}_i$ is a cycle rather than a path, we have \begin{align} &\mathrm{Var}(\hat{\Gamma}_i) = \frac{1}{p} \sum_{j=1}^{k(i)-1} Y_{v_{i,j},v_{i,j+1}}^2 + 2 \sum_{1\le j<q \le k(i)-1} p^{q-j-1} Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}+ \frac{p^2+2p-p^{k(i)-1}}{1+p^{k(i)-1}} Y_{v_{i,k(i)},v_{i,k(i)+1}}^2\notag \\ &+ 2 \sum_{1\le j \le k(i)-1} (-1)^j \left(\frac{(1-(-p)^{k(i)-1-j})(1-(-p)^{j-1})}{(1+p^{k(i)-1})}-1\right)Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,k(i)},v_{i,k(i)+1}}. \end{align}

The cycle expression differs from \Eqref{eq:variance-of-path} only through the contribution of the final edge. Along an alternating path, correlations between edge selections decay geometrically with distance. In contrast, in the cycle case the last edge is also constrained by the first edge, closing the loop and inducing additional dependence between the final edge and all preceding edges. This additional dependence gives rise to the extra terms in \Eqref{eq:variance-of-cycle}.

Note that the variance expressions depend on the potential outcomes $Y_{i,j}$, for which we have not imposed any modeling assumptions. To study the order of magnitude of the variance and to facilitate comparison across designs, we introduce a boundedness assumption on potential outcomes below.

assumption[Bounded Potential Outcomes] There exists a constant $B>0$ such that for all $(i,j)\in{\mathbb{M}}$, the potential outcome satisfies $Y_{i,j}\in[0,B]$.

Assumption (ref) is mild in practice, as the upper bound $B$ may be large and need not be known to the experimenter. Similar boundedness conditions are standard in the design-based causal inference literature (see, e.g., bojinov2023design, chattopadhyay2023design, chen2025efficient). Under this assumption, we can study the worst-case variance of the AP design, defined as the largest variance attainable over all potential outcomes satisfying the bound. Formally, for a fixed design parameter $p$ and alternating path/cycle decomposition ${\triangle\mathcal{P}^{(t,c)}}$, we consider

equation[equation omitted — 500 chars of source]

where the second equality follows from Proposition (ref), since randomization is independent across components. Fortunately, The maximization problem in \Eqref{eq:max-variance} admits a simple solution that the worst-case variance is achieved when all relevant potential outcomes attain their upper bound $B$.

propositionThe solution to the maximization problem of \Eqref{eq:max-variance} is $Y_{i,j}=B$ for all $(i,j)\in {\mathbb{M}}$. Particularly, under Assumption (ref), all the correlation terms in \Eqref{eq:variance-of-path} and \Eqref{eq:variance-of-cycle} are positive.

Proposition (ref) implies that all coefficients multiplying the potential outcomes in Theorem (ref) are nonnegative, including those associated with cross-edge correlations. Consequently, the largest variance is achieved when all potential outcomes attain their maximal value $B$. The intuition is as follows. Along an alternating component ${\mathcal{P}}_i=\{v_{i,1},\cdots,v_{i,k(i)+1}\}$, the adjunct two indicators $W_{i,j}$ and $W_{i,j+1}$ are negatively dependent due to feasibility constraints that if one edge is realized, the adjacent edge must be skipped. However, adjacent edges belong to different matchings (treatment v.s. control), so their contributions enter the estimator with opposite signs. These two effects combine to produce positive covariance terms in the variance expression. Related phenomena have been observed in switchback experiments with temporal structure (e.g., bojinov2023design, chen2025efficient), although the underlying dependence in our setting arises from matching interference rather than time. Building on Proposition (ref), we now characterize the worst-case variance contribution of a single alternating path or cycle under the AP design.

theoremSuppose Assumption (ref) holds, for ${\mathcal{P}}_s=\{v_{s,1},\cdots,v_{s,k(s)+1}\}\in{\triangle\mathcal{P}^{(t,c)}}$, if ${\mathcal{P}}_s$ is a path, we have \begin{equation} \max_{\substack{ Y_{i,j} \in [0,B] \\ \forall (i,j)\in{\mathcal{P}}_s }} \mathrm{Var} (\hat{\Gamma}_s)= \left(\frac{1}{p}+\frac{2}{1-p}\right) B^2 k(s) + \frac{2B^2(p^{k(s)}-1)}{(1-p)^2}. \end{equation} If ${\mathcal{P}}_s$ is a cycle, we have \begin{align} \max_{\substack{ Y_{i,j} \in [0,B] \\ \forall (i,j)\in{\mathcal{P}}_s }} \mathrm{Var} (\hat{\Gamma}_s)&= \left(\frac{1}{p}+\frac{2}{1-p}\right) B^2 (k(s)-1) +\frac{2B^2(p^{k(s)-1}-1)}{(1-p)^2}\notag\\ &\qquad\qquad+\frac{B^2(4+2p-p^2-p^3-p^{k(s)-2}(2+p+p^2))}{(1-p)(1+p^{k(s)-1})}. \end{align} Moreover, for both altering paths and cycles, we can have \begin{equation*} \lim_{k(s)\rightarrow \infty} \frac{\max_{\substack{ Y_{i,j}\le [0,B] \\ \forall (i,j)\in{\mathcal{P}}_s }} \mathrm{Var} (\hat{\Gamma}_s)}{k(s)} = \frac{B^2(1+p)}{p(1-p)}. \end{equation*} Equivalently, under the worst case, $\mathrm{Var} (\hat{\Gamma}_s)=\frac{B^2(1+p)}{p(1-p)} k(s) +o(k(s))$.

Theorem (ref) shows that for each path/cycle, the variance of the estimator grows only linearly in its length. Therefore, aggregating over components yields

equation*[equation* omitted — 252 chars of source]

where we have used $\sum_{s}k(s)=N$ for one-to-one matching. Thus, the variance of the AP estimator decreases at rate $1/N$, despite the dependence induced by matching interference. For comparison, under the naive design in Proposition (ref), the worst-case variance $\mathrm{Var}(\hat{\tau}_{\text{Naive}})$ is $4B^2$, which does not vanish as $N\rightarrow\infty$. Finally, Theorem (ref) implies that although the indicators $W_{i,j}$ are dependent along each component, the induced correlations decay sufficiently fast with distance that asymptotic efficiency is preserved.

Optimal AP Design under Minimax Rule

With Theorem (ref) in hand, we now turn to optimizing the randomization probability $p$. In the previous sections, for the ease of presentation, we set all the paths and cycles to share a common value of $p$. However, since randomization across alternating paths is independent, each alternating path/cycle ${\mathcal{P}}_s$ can in principle adopt its own design parameter $p_s$. In this section, we therefore concentrate on optimizing the probability for one path/cycle. Following a common practice in decision theory (berger2013statistical) and experimental design (bojinov2023design), we formulate the design problem into a minimax problem, under which we seek to minimize the worst case risk, i.e.,

equation[equation omitted — 284 chars of source]

The equation holds since our $\hat{\Gamma}_s$ is unbiased. Thus, we are actually trying to minimize the RHS of \Eqref{eq:path-variance-upper-bound} or \Eqref{eq:cycle-variance-upper-bound} with respect to $p$. For a path/cycle of finite length $k(s)$, a closed-form minimizer of \Eqref{eq:path-variance-upper-bound} or \Eqref{eq:cycle-variance-upper-bound} is generally not available due to the mixed polynomial–exponential structure of the variance expression. However, when the path is long, the normalized worst-case variance is asymptotically linear in $k(s)$ with leading coefficient $(1+p)/(p(1-p))$. This yields a limiting optimal probability.

corollaryLet $p^*_{k(s)}$ denote the optimal solution $p$ to \Eqref{eq:minimax} of a path with length $k(s)$. Then, we have that $\lim_{k(s)\rightarrow\infty}p^*_{k(s)}=\sqrt{2}-1\approx 0.4142$.

Thus, for long alternating paths, the AP design asymptotically favors a conditional randomization probability strictly below $1/2$. Intuitively, choosing $p<1/2$ reflects the fact that we are not only minimizing the variance contribution of the current edge, but must also account for downstream edges and control the propagation of interference along the path.

Table (ref)reports numerical solutions for the optimal $p^*_{k(s)}$ under both the path case and the cycle case for various values $k(s)$. When $k(s)=1$, there is no interference since only one match exists, so we can set $p=1$, and the match is always included. For this reason, the table begins at $k(s)=2$. In the cycle case, as we have discussed in Proposition (ref), the minimum feasible length of any cycle is $k(s)=4$ and $k(s)$ must be even. For both paths and cycles, the optimal value satisfies $p_{k(s)}^*=1$ whenever $k(s)\le 4$. In these short cases, the minimax-optimal strategy coincides with applying the naive design on the component that with probability $1/2$ one adopts all treatment edges and with probability $1/2$ one adopts all control edges. As the length increases, the AP design gradually yields lower worst-case variance compared to the naive design. The case $k(s)=4$ therefore serves as the critical breakpoint beyond which a carefully chosen $p<1$ is beneficial. Consistent with the asymptotic analysis, when $k(s)\ge 50$ the optimal probability $p_{k(s)}^*$ has effectively stabilized near $\sqrt{2}-1$.

table[table omitted — 1,021 chars of source]

Asymptotic Inference

In this section, we present the asymptotic inference result of our estimator $\hat{\tau}_{\text{AP}}$.

theorem[Central Limit Theorem] For AP design, if Assumption (ref) holds, the size of the disagreement set goes to infinity i.e., $|{\triangle}{\mathbb{M}}^{(t,c)}|\rightarrow\infty$ and $N^2\mathrm{Var}(\hat{\tau}_{\text{AP}})\rightarrow \infty$, then, as $N\rightarrow\infty$, \begin{equation*} \frac{\hat{\tau}_{AP}-\tau}{\sqrt{\mathrm{Var}(\hat{\tau}_{AP})}} \xrightarrow{d} \mathcal{N}(0,1), \end{equation*} where $\xrightarrow{d}$ denotes convergence in distribution.

Although Theorem (ref) resembles a standard central limit result, it does not appear to follow from a direct application of well-known CLTs. The difficulty arises from the potentially complex structure of the alternating path/cycle decomposition ${\triangle\mathcal{P}^{(t,c)}}$. Figure (ref) illustrates three qualitatively distinct configurations that may occur. Configuration 1 corresponds to the case where all components in ${\triangle\mathcal{P}^{(t,c)}}$ are uniformly small. In this regime, the classical Lindeberg–Feller CLT (e.g., feller1991introduction) applies directly, as dependence is localized and vanishes in aggregate. Configuration 2 corresponds to the opposite extreme, where there exists a single (or a small fixed number of) components of diverging length. In this regime, we cannot appeal to i.i.d. or Lindeberg–Feller arguments, since we impose no distributional model on $Y_{i,j}$. However, on an infinitely long alternating path, the sequence $W_{i,j}$ behaves (up to boundary effects) as a Markov chain with transition matrix: \[ {\bm{P}} =

pmatrix[pmatrix omitted — 33 chars of source]

, \] so that nonstationary $\alpha$-mixing CLT results (e.g., bradley2017central) may be invoked.

The substantial challenge arises in Configuration 3, where ${\triangle\mathcal{P}^{(t,c)}}$ contains a mixture of finitely long and infinitely long components. In this case, no single CLT applies uniformly. Moreover, as $N$ grows, the configuration itself may fluctuate among the three types, and even within Configuration 3 the relative proportions of short and long components need not converge. At a high level, our proof partitions ${\triangle\mathcal{P}^{(t,c)}}$ into two groups: a long-path group ${\triangle\mathcal{P}^{(t,c)}}_L$ and a short-path group ${\triangle\mathcal{P}^{(t,c)}}_S$. We define the path/cycle whose variance $\mathrm{Var}(\hat{\Gamma}_i)$ is larger than $(\sum_i \mathrm{Var}(\hat{\Gamma}_i))^{2/3}$ as the long path, otherwise as a short path. Both $|{\triangle\mathcal{P}^{(t,c)}}_L|$ and $|{\triangle\mathcal{P}^{(t,c)}}_S|$ may diverge, and neither group is assumed to be stabilized as $N$ grows. sing the Bolzano–Weierstrass theorem together with the subsequence principle (Theorem 2.6 in billingsley2013convergence), we show that every subsequence of contributions from ${\triangle\mathcal{P}^{(t,c)}}_L$ and ${\triangle\mathcal{P}^{(t,c)}}_S$ admits a further subsequence along which a normal limit holds. On a side note, for infinite components, the only deviation from the Markovian structure comes from the final edge of a cycle, whose boundary effect vanishes relative to the component length. For finite components, the interior dependence does not affect the limit as long as the component size is uniformly bounded. Together, these arguments establish the asymptotic normality of $\hat{\tau}_{\text{AP}}$ under the stated conditions.

figure[figure omitted — 204 chars of source]

A practical caveat here is that the closed-form expressions in \Eqref{eq:variance-of-path} and \Eqref{eq:variance-of-cycle} are not identifiable from observed data. In particular, terms of the form $Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}$ with $q=j+1$ cannot be observed, since matching interference guarantees that two consecutive edges along an alternating component are never simultaneously realized. As a consequence, the true variance cannot be estimated directly. To construct confidence intervals and yield implementable inference, we therefore derive an estimable upper bound on the true variance. Using the Cauchy–Schwarz inequality, we obtain an upper bound for the path case that depends only on observable quantities and can be estimated unbiasedly from the realized data. The cycle case follows the same principle, and the explicit expression is deferred to Appendix (ref).

propositionFor ${\mathcal{P}}_i=\{v_{i,1},\cdots,v_{i,k(i)+1}\}\in{\triangle\mathcal{P}^{(t,c)}}$, if ${\mathcal{P}}_i$ is a path, we have the following upper bound $\tilde{\sigma}^2_i$ of $\mathrm{Var}(\hat{\Gamma}_i)$, \begin{equation} \tilde{\sigma}^2_i := \left(\frac{1}{p}+1\right) \sum_{j=1}^{k(i)} Y_{v_{i,j},v_{i,j+1}}^2+\sum_{j=2}^{k(i)-1} Y_{v_{i,j},v_{i,j+1}}^2 + 2 \sum_{\substack{1\le j \le k(i)-2 \\ j+2\le q \le k(i)}} p^{q-j-1} Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}. \end{equation} Moreover, $\tilde{\sigma}^2_i$ can be estimated by \begin{align*} \hat{\sigma}^2_i&:= \left(\frac{1}{p}+1\right) \sum_{j=1}^{k(i)} \frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y^2_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}+\sum_{j=2}^{k(i)-1} \frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y^2_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}\\ &\qquad\qquad\qquad + 2 \sum_{\substack{1\le j \le k(i)-2 \\ j+2\le q \le k(i)}} p^{q-j-1} \frac{Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}{\mathbb{I}}\{w_{i,j}=1,w_{i,q}=1\}}{\mathbb{P}(W_{i,j}=1,W_{i,q}=1)}, \end{align*} where $\mathbb{P}(W_{i,j}=1)=p/(p+1)$ and $\mathbb{P}(W_{i,j}=1,W_{i,q}=1)=(p^2-(-p)^{q-j+1})/(p+1)^2$. The estimator $\hat{\sigma}^2_i$ is an unbiased estimator of $\tilde{\sigma}^2_i$, i.e., $\mathbb{E}[\hat{\sigma}^2_i]=\tilde{\sigma}^2_i$.

There are alternative ways to obtain conservative variance bounds. For example, we can bound all the cross term $Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}$ by $1/2(Y_{v_{i,j},v_{i,j+1}}^2+Y_{v_{i,q},v_{i,q+1}}^2)$, regardless of whether they are observable or not. In addition, we note that Proposition (ref) does not rely on the boundedness assumption (Assumption (ref)). If we do wish to exploit boundedness, the unobserved outcomes can always be bounded by the constant $B$. Using Proposition (ref), $\sqrt{\mathrm{Var}(\hat{\tau}_{\text{AP}})}$ can be estimated slightly conservatively via $\hat{\sigma}^2:={1}/{N^2} \sum_{{\mathcal{P}}_i\in {\triangle\mathcal{P}^{(t,c)}}} \hat{\sigma}^2_i$. Then, we can form an approximate asymptotic $\alpha$-level confidence interval for the AP estimator by \[\left[\hat{\tau}_{\text{AP}}-c_{\frac{1+\alpha}{2}}\hat{\sigma},\hat{\tau}_{\text{AP}}+c_{\frac{1+\alpha}{2}}\hat{\sigma}\right],\] where $c_{(1+\alpha)/2}$ is the $(1+\alpha)/2$ quantile of the standard normal distribution.

Many-to-One Matching

We now extend the AP design to two-sided many-to-one environments, in which agents on one side of the market may match with multiple agents on the other side. Such settings arise frequently in applications, including student–school assignments (chen2016school), split-liver transplantation (tang2025multi, tang2025split), ridesharing and carpooling platforms (castillo2025matching), and refugee resettlement (bansak2024dynamic). In each of these examples, suppliers possess a capacity that constrains the number of demands they may serve, while each demand may be served by at most one supplier.

Formally, let the set of suppliers be ${\mathbb{S}}:=\{s_1,\cdots,s_N\}$ and the set of demands be ${\mathbb{D}}:=\{d_1,\cdots,d_M\}$. Each supplier $s_i\in{\mathbb{S}}$ has capacity $C_i$, meaning it may be matched with at most $C_i$ demands, while each demand may be matched with at most one supplier. For simplicity of exposition we assume homogeneous capacities $C_i=C_0$ for all $s_i$, though our results extend to heterogeneous capacities without conceptual difficulty. Let $m_{i,j}\in\{0,1\}$ indicate whether supplier $i$ is matched with demand $j$. A feasible matching ${\mathbb{M}}$ in many-to-one matching with capacity constraint is then defined by

equation[equation omitted — 191 chars of source]

and we let ${\mathcal{M}}$ denote the set of all feasible many-to-one matchings.

Our inferential goal remains the same: given two feasible matchings ${\mathbb{M}}^t \in {\mathcal{M}}$ and ${\mathbb{M}}^c \in {\mathcal{M}}$, we seek to experimentally compare their average performance without imposing outcome or behavioral models. For a potential outcome $Y_{i,j}$ defined for each feasible pair $(i,j)$, we define the ATE analogously to the one-to-one case as

equation[equation omitted — 137 chars of source]

The normalization by $C_0 N$ reflects the maximum number of possible matches. Normalizing by the number of demand $M$ is also reasonable and would yield an equivalent estimand up to scale. Finally, following the philosophy of (ref), we define the disagreement set as

equation[equation omitted — 187 chars of source]

If we can identify a collection of alternating paths and cycles whose randomized realizations form a feasible matching in ${\mathcal{M}}$, then the AP design and the associated inference results from the one-to-one setting continue to apply without modification. Thus, the new challenges introduced by the many-to-one setting arise not from inference, but from ensuring that the decomposition of the disagreement set admits a feasible randomized implementation that respects all capacity constraints. In the one-to-one case this issue does not arise. Proposition (ref) guarantees that the alternating path/cycle decomposition is unique, and every alternating component is automatically feasible under randomization. In contrast, in the many-to-one case, (i) the disagreement set may admit multiple distinct alternating decompositions, and (ii) not all decompositions are feasible, since randomization along certain paths may violate supplier capacities.

We therefore begin by introducing a set of sufficient conditions under which an alternating decomposition of the disagreement set yields a feasible randomized implementation under the AP design for the many-to-one matching.

theoremLet ${\triangle\mathcal{P}^{(t,c)}}:=\{{\mathcal{P}}_1,\cdots, {\mathcal{P}}_p\}$ be a collection of alternating paths and cycles. If the following conditions hold, then the AP design yields a feasible randomized matching in ${\mathbb{M}}$ and the central limit theorem in Theorem (ref) continues to hold: \begin{enumerate} • Each match in ${\triangle}{\mathbb{M}}^{(t,c)}$ is covered exactly once by ${\triangle\mathcal{P}^{(t,c)}}$; • Each supplier appears in at most $C_0$ alternating paths and cycles in ${\triangle\mathcal{P}^{(t,c)}}$; • Each demand appears in at most one alternating path or cycle in ${\triangle\mathcal{P}^{(t,c)}}$; • Within each ${\mathcal{P}}_i\in{\triangle\mathcal{P}^{(t,c)}}$, no supplier or demand is traversed or visited more than once. \end{enumerate}

Note that one implicit condition is that each path or cycle must be an alternating path or cycle. The first condition ensures that every disagreement edge is accounted for exactly once. If an edge is omitted or duplicated, the resulting estimator would be biased. Conditions 2–4 jointly guarantee feasibility under the AP design. In particular, Condition 4 implies that any supplier (or demand) appears on a given path at most twice, once as an incoming edge and once as an outgoing edge, and under the AP design at most one of these two edges can be realized due to alternating interference. Combined with Condition 2, this ensures that no supplier exceeds its capacity across all paths. Condition 3 analogously guarantees demand feasibility. The necessity of Condition 4 may not be immediately obvious. Figure (ref) illustrates a decomposition in which $C_0=2$ and Conditions 1–3 are satisfied, but Condition 4 fails because supplier $s_2$ appears three times along the same path. Under the AP design, it is then possible for edges $(d_3,s_2)$ and $(d_4,s_2)$ from Path 1 and $(d_2,s_2)$ from Path 2 to be realized simultaneously, violating the capacity constraint of supplier $s_2$. Condition 4 rules out precisely such cases.

figure[figure omitted — 226 chars of source]

We next describe an efficient procedure for constructing a collection of alternating paths and cycles that satisfies the conditions in Theorem (ref). The key idea is to represent the disagreement set as an auxiliary directed graph and reduce the decomposition problem to two classical graph-theoretic tasks: finding augmenting paths and performing Euler–tour decompositions.

definition[Auxiliary Unbalanced Directed Graph] Given a disagreement set ${\triangle}{\mathbb{M}}^{(t,c)}$, define its auxiliary directed graph ${\mathcal{G}}=(\mathcal{V},\mathcal{E})$ as follows. Vertices. The vertex set $\mathcal{V}$ consists of: \begin{enumerate} • Supplier vertices $v_{s_i}$ for each $s_i \in {\mathbb{S}}$; • Degree-1 demand vertices $v_{d_j}$ for each $d_j \in {\mathbb{D}}$ that is matched with exactly one supplier in ${\triangle}{\mathbb{M}}^{(t,c)}$. \end{enumerate} Directed edges. The edge set $\mathcal{E}$ consists of two types: \begin{enumerate} • Supplier-to-supplier edges. For each demand $d_k$ such that both $(s_i,d_k)$ and $(s_j,d_k)$ lie in ${\triangle}{\mathbb{M}}^{(t,c)}$, if $(s_i,d_k)\in{\mathbb{M}}^t$ and $(s_j,d_k)\in{\mathbb{M}}^c$, then include a directed edge from $s_i$ to $s_j$, $(v_{s_i},v_{s_j})$ in $\mathcal{E}$ and label it by $d_k$. • Supplier–demand edges. For each degree-1 demand vertex $v_{d_j}$, if $(s_k,d_j)\in{\mathbb{M}}^t$ then include a directed edge $(v_{s_k},v_{d_j})$. If $(s_k,d_j)\in{\mathbb{M}}^c$ then include $(v_{d_j},v_{s_k})$ instead. \end{enumerate}

In Figure (ref), we showcase the auxiliary unbalanced directed graph of the disagreement set in Figure (ref). There are several useful observations and properties of the graph that we want to highlight. Each directed edge represents one or two matches in the disagreement set. Demands of degree two in ${\triangle}{\mathbb{M}}^{(t,c)}$ (i.e., the demands that are matched to two different suppliers in the treatment matching and the control matching) induce supplier-to-supplier edges, while degree-one demands induce supplier–demand edges. Multiple directed edges may exist between the same pair of suppliers, and they may point in the same or opposite directions. For each supplier node in ${\mathcal{V}}$, the number of the outgoing edges (i.e., the out-degree, $\text{deg}^+(v)$) is equal to the number of matches associated with the supplier in the treatment matching. Correspondingly, the in-degree of each supplier node, $\text{deg}^-(v)$, is the number of matches in the control matching associated with node. Thus, for each supplier node $v_{s_i}$, neither $\text{deg}^-(v_{s_i})$ nor $\text{deg}^+(v_{s_i})$ can be larger than $C_0$. Each degree-1 demand node $v_{d_i}$ has either $\text{deg}^-(v_{d_i})=1$ or $\text{deg}^+(v_{d_i})=1$ depending on whether the associated match is in the treatment matching or the control matching. Moreover, the graph is generally unbalanced, meaning that for at least some vertices $v$ such that $\text{deg}^-(d)\not = \text{deg}^+(d)$. Every disagreement set ${\triangle}{\mathbb{M}}^{(t,c)}$ induces a unique auxiliary unbalanced directed graph ${\mathcal{G}}$. Conversely, every simple directed path or simple directed cycle in ${\mathcal{G}}$ corresponds uniquely to an alternating path or cycle in ${\triangle\mathcal{P}^{(t,c)}}$. We denote this inverse mapping by ${\mathcal{A}}^{-1}$, so that for any directed path $\pi$ in ${\mathcal{G}}$, the object ${\mathcal{A}}^{-1}(\pi)$ is an alternating path or cycle in ${\triangle\mathcal{P}^{(t,c)}}$.

figure[figure omitted — 212 chars of source]

With the auxiliary graph in place, our task reduces to finding a collection of simple directed paths and cycles that cover all edges of ${\mathcal{G}}$ while respecting the capacity constraints. Define the source and sink vertex sets as

equation*[equation* omitted — 184 chars of source]

and call vertices with $\text{deg}^+(v)=\text{deg}^-(v)$ balanced. If every vertex in the graph is balanced, then ${\mathcal{G}}$ is an Eulerian digraph, meaning that its edges can be decomposed entirely into directed cycles. Our approach is therefore to first find a set of simple directed paths that originate in ${\mathcal{S}}$ and terminate in ${\mathcal{D}}$. Removing these paths from ${\mathcal{G}}$ yields a residual subgraph in which all vertices are balanced, and hence the residual is Eulerian.

To extract a collection of simple directed paths, we interpret the directed edges of $G$ as arcs of a unit-capacity flow network. Construct a flow network ${\mathcal{G}}'=({\mathcal{V}}',{\mathcal{E}}')$ as follows: retain all original vertices ${\mathcal{V}}$ and edges ${\mathcal{E}}$, and assign to each original edge $e\in {\mathcal{E}}$ a capacity $c(e)=1$. Add two additional vertices, a super-source $\mathbf{s}$ and a super-sink $\mathbf{t}$. For each $v\in {\mathcal{S}}$, add an arc $(\mathbf{s},v)$ with capacity $c(\mathbf{s},v)=\text{deg}^+(v)-\text{deg}^-(v)$; for each $v\in {\mathcal{D}}$, add an arc $(v,\mathbf{t})$ with capacity $c(v,\mathbf{t})=\text{deg}^-(v)-\text{deg}^+(v)$. Figure (ref) includes a constructed flow network. Since $\sum_{v\in V} \text{deg}^+(v)-\text{deg}^-(v)=0$, we have \[ \sum_{v\in {\mathcal{S}}} \text{deg}^+(v)-\text{deg}^-(v) = \sum_{v\in {\mathcal{D}}} \text{deg}^-(v)-\text{deg}^+(v):=C', \] so the total supply out of $\mathbf{s}$ equals the total demand into $\mathbf{t}$. Therefore, the max flow on ${\mathcal{G}}'$ from $\mathbf{s}$ to $\mathbf{t}$ is with capacity $C'$. Because all original edges have unit capacity, no original edge is used more than once when carrying flow. Furthermore, standard augmenting-path algorithms (e.g., Edmonds--Karp algorithm edmonds1972theoretical, see also Algorithm (ref) in the appendix) ensure that every augmenting path encountered during the construction is simple in ${\mathcal{G}}'$, and therefore the projected path in ${\mathcal{G}}$ is also simple. For simplicity, we define the mapping ${\mathcal{F}}^{-1}$ which projects the path ${\mathcal{P}}'$ on ${\mathcal{G}}'$ back to the original path ${\mathcal{P}}$ on ${\mathcal{G}}$. Thus the support of the max flow on original edges decomposes into edge-disjoint simple directed paths starting at vertices in ${\mathcal{S}}$ and ending at vertices in ${\mathcal{D}}$. Let $\widetilde{{\mathcal{P}}}_p$ denote this family of paths on the auxiliary unbalanced directed graph.

Remove from ${\mathcal{G}}$ all edges used by the paths in $\mathcal{P}$ along with the nodes that subsequently have both in-degree and out-degree equal to 0. Let the resultingdirected graph be denoted by ${\mathcal{G}}_{\mathrm{rem}}=({\mathcal{V}}_{\mathrm{rem}},{\mathcal{E}}_{\mathrm{rem}})$. For each $v\in {\mathcal{S}}$, exactly $\text{deg}^+(v)-\text{deg}^-(v)$ outgoing edges are removed and no incoming edges are removed. For each $v\in {\mathcal{D}}$, exactly $\text{deg}^-(v)-\text{deg}^+(v)$ incoming edges are removed and no outgoing edges are removed. Moreover, for each originally balanced vertex $v$, the removal of a path that passes through $v$ eliminates one incoming and one outgoing edge leaving the vertex still balanced. Hence for all $v\in {\mathcal{V}}_{\mathrm{rem}}$, the remaining in-degree and out-degree satisfy $\deg^+_{{\mathcal{G}}_{\mathrm{rem}}}(v) = \deg^-_{{\mathcal{G}}_{\mathrm{rem}}}(v)$. Thus, ${\mathcal{G}}_{\mathrm{rem}}$ is an Eulerian directed graph (possibly disconnected). A classical graph-theoretic result (see, e.g., bang2008digraphs) states that every Eulerian directed graph decomposes into edge-disjoint simple directed cycles. For completeness, we present an explicit cycle-decomposition procedure in Algorithm (ref) in the appendix. Applying this result to ${\mathcal{G}}_{\mathrm{rem}}$ yields a family $\widetilde{{\mathcal{P}}}_c$ of simple directed cycles whose edges partition ${\mathcal{E}}_{\mathrm{rem}}$.

Finally, projecting each path in $\widetilde{{\mathcal{P}}}_p$ and each cycle in $\widetilde{{\mathcal{P}}}_c$ back by ${\mathcal{A}}^{-1}$ can offer us a corresponding collection of alternating paths and alternating cycles in the original disagreement set. The resulting set of the alternating paths and cycles satisfies all the conditions in Theorem (ref). The details of the whole procedure can be found in Algorithm (ref). Once such a decomposition is obtained, we can directly apply the AP design and all inference guarantees developed for the one-to-one matching case.

algorithm[algorithm omitted — 1,983 chars of source]

Numerical Results

We evaluate the real-world performance of our AP design on one of the largest workforce datasets provided by Revelio Lab\footnotemark \footnotetext{See, \url{https://wrds-www.wharton.upenn.edu/pages/about/data-vendors/revelio-labs/}.}. The dataset contains rich job–worker matching information across a broad range of occupations and industries. Following tang2024match, we focus on a subset consisting of mid-level software engineers in the United States during the period 2010–2015, with employment duration as the primary outcome of interest.

To construct a controlled experimental environment, we cluster workers and jobs into 50 groups each. Using the full employment history data, we form a ground-truth matching-value matrix $\mathbf{Y}^*$ where each entry is defined as the average of an indicator for whether an employment spell exceeds six months, taken over all worker–job matches between worker group $i$ and job group $j$. With full access to $\mathbf{Y}^*$, the entry $Y^*_{i,j}$ can be interpreted as the probability that a worker from group $i$ remains employed for longer than six months when matched to a job from group $j$. By construction, $Y^*_{i,j}\in[0,1]$. Throughout the numerical study, we interpret this probability as a proxy for worker satisfaction with the assigned job: a higher value of $Y^*_{i,j}\in[0,1]$ indicates a more stable and hence more desirable match from the worker’s perspective.

Table (ref) reports a variance comparison in a controlled setting with $n$ workers and $n$ jobs, where one worker and one job are drawn from each group and $n$ ranges from 10 to 50. The two matching plans under comparison are constructed as follows. Under the treatment plan, worker $k$ is matched to job $k$. Under the control plan, worker $k$ is matched to job $k+1$ for $1\le k\le n-1$, and worker $n$ is matched to job $1$, forming a cyclic shift. The table compares the AP design (with randomization probability $p=0.5$) against the naive design described in Remark (ref). “True Var” denotes the true variance of the estimator computed using the full matching-value matrix $\mathbf{Y}^*$. “Empirical Var” is the sample variance obtained from 500 independent randomizations. “Var UB” refers to the theoretical variance upper bound, and “Est. Var UB” denotes its data-driven estimate derived from Proposition (ref). Consistent with our theoretical analysis, the naive design exhibits a variance that does not diminish as $n$ increases. In contrast, the variance of the AP design decreases steadily with the size of the problem, illustrating its scalability. Moreover, the estimated variance upper bound closely tracks the true variance across all values of $n$, suggesting that the proposed bound is both informative and practically tight.

table[table omitted — 1,189 chars of source]

To verify the central limit theorem, we plot out the normalized empirical distribution of 5000 independent runs of the randomization in Figure (ref), together with a standard normal as a comparison. Furthermore, Figure (ref) presents a Q–Q plot comparing the empirical quantiles of the sample with those of a standard normal distribution. The points align closely with the 45-degree reference line across the entire range, indicating that the empirical distribution exhibits strong agreement with normality. Only negligible deviations are observed in the extreme tails. This visual evidence is consistent with the results of the Shapiro–Wilk test (test statistic $= 0.9999$, p-value $= 0.2817$), which fails to reject the null hypothesis of normality. Collectively, these findings indicate that the distribution of the data does not significantly deviate from a standard normal distribution.

figure[figure omitted — 483 chars of source]

Figure (ref) illustrates how the variance of the AP estimator and its estimable upper bound from Proposition (ref) vary with the randomization probability $p$ when $n=50$. The true variance and the upper bound nearly coincide across the entire range of $p$, providing further empirical evidence that the proposed variance bound is tight for this dataset. The empirically optimal randomization probability is $p^*\simeq 0.417$, which is close to the worst-case optimal value $\sqrt{2}-1$. An additional observation is that the variance curve is relatively flat in a neighborhood around $p^*$, indicating that the performance of the AP design is robust to moderate misspecification of the optimized randomization probability.

figure[figure omitted — 244 chars of source]

Discussion and Concluding Remarks

In this work, we study a fundamental yet underexplored challenge in operations management: how to experimentally compare matching mechanisms in an assumption-free, design-based manner in the presence of matching interference. While much of the existing literature focuses on mitigating bias or improving efficiency within a fixed matching mechanism, we take a complementary perspective and develop a principled framework for directly testing two predetermined matching plans, such as a newly proposed algorithm against a status quo, without imposing outcome or behavioral models.

Our key methodological contribution is the Alternating Path Randomized Design, which exploits the structure of the disagreement set between two matchings. In the one-to-one setting, we show that this disagreement set admits a unique decomposition into alternating paths and cycles, along which sequential randomization can be conducted while respecting feasibility constraints. Within a minimax framework, we characterize the variance properties of the resulting Horvitz–Thompson estimator and show that, for long paths and cycles, the optimal conditional randomization probability converges to $\sqrt{2}-1$. We further establish a finite-population central limit theorem that remains valid even when the configuration of paths and cycles is complex, heterogeneous, and potentially unstable as the population grows.

We extend the framework to many-to-one matchings, where capacity constraints introduce substantial new challenges. By mapping the disagreement structure to an auxiliary unbalanced directed graph, we reduce the construction of feasible alternating paths and cycles to classical graph-theoretic problems involving augmenting paths and Euler-tour decompositions. This extension highlights that adapting experimental designs from one-to-one to many-to-one settings is far from mechanical and requires careful attention to feasibility and interference.

Overall, this work offers a unified, design-based methodology for testing matching mechanisms under interference, with direct relevance for platforms and policymakers seeking to validate new matching algorithms before full-scale deployment. We hope these results will encourage further research on experimental design in complex allocation systems where interference is intrinsic rather than incidental.

\ECSwitch

center[center omitted — 81 chars of source]

Technical Tools

The first important technical tool is the non-stationary $\alpha$-Mixing CLT by bradley2017central in their Theorem 1.1.

theorem[Non-Stationary $\alpha$-Mixing CLT, bradley2017central] Suppose \(d\) is a positive integer. For each \(n \in \mathbb{N}\), suppose \(L_n := (L_{n1}, L_{n2}, \dots, L_{nd})\) is an element of \(\mathbb{N}^d\), and suppose \(X^{(n)} := (X^{(n)}_k, k \in B(L_n))\) is an array of random variables such that for each \(k \in B(L_n)\), \(\mathbb{E} X^{(n)}_k = 0\) and \(\mathbb{E}(X^{(n)}_k)^2 < \infty\). Suppose the following mixing assumptions hold: \[ \alpha(m) := \sup_n \alpha(X^{(n)}, m) \to 0 \quad \text{as } m \to \infty, \] \[ \rho'(1) := \sup_n \rho'(X^{(n)}, 1) < 1. \] For each \(n \in \mathbb{N}\), define the random sum \(S(X^{(n)}, L_n) = \sum_{k \in B(L_n)} X^{(n)}_k\), define the quantity \(\sigma_n^2 := \mathbb{E}[S(X^{(n)}, L_n)]^2\), and assume that \(\sigma_n^2 > 0\). Suppose also that the Lindeberg condition \[ \forall \varepsilon > 0, \quad \lim_{n \to \infty} \frac{1}{\sigma_n^2} \sum_{k \in B(L_n)} \mathbb{E}[(X^{(n)}_k)^2 {\mathbb{I}}(|X^{(n)}_k| > \varepsilon \sigma_n)] = 0 \] holds. Then \(\sigma_n^{-1} S(X^{(n)}, L_n) \xrightarrow{d} N(0, 1)\) as \(n \to \infty\), where \(\xrightarrow{d}\) denotes convergence in distribution.

The theorem relies on the following concepts defined on a probability space \((\Omega, \mathcal{F}, P)\):

itemize• Strong Mixing Coefficient: For any two \(\sigma\)-fields \(A, B \subseteq \mathcal{F}\), \[ \alpha(A, B) := \sup_{A \in A, B \in B} |P(A \cap B) - P(A)P(B)|. \] • Maximal Correlation Coefficient: For any two \(\sigma\)-fields \(A, B \subseteq \mathcal{F}\), \[ \rho(A, B) := \sup \{ |\mathrm{Corr}(f, g)| : f \in L^2_{\text{real}}(A), g \in L^2_{\text{real}}(B) \}. \] • Random Field: \(X := (X_k, k \in \mathbb{Z}^d)\) is a random field, not necessarily stationary. • Box \(B(L)\): For \(L = (L_1, L_2, \dots, L_d) \in \mathbb{N}^d\), \[ B(L) := \{ k = (k_1, k_2, \dots, k_d) \in \mathbb{N}^d : \forall u \in \{1, 2, \dots, d\}, 1 \leq k_u \leq L_u \}, \] with \(|B(L)| = L_1 \cdot L_2 \cdot \dots \cdot L_d\). • \(\alpha\)-Mixing for Random Fields: The random field \(X\) is \(\alpha\)-mixing if \[ \alpha(X, n) := \sup \alpha(\sigma(X_k, k \in Q), \sigma(X_k, k \in S)), \] where the supremum is over all pairs of nonempty, disjoint sets \(Q, S \subset \mathbb{Z}^d\) such that there exist \(u \in \{1, 2, \dots, d\}\) and \(j \in \mathbb{Z}\) with \(Q \subset \{ k : k_u \leq j \}\) and \(S \subset \{ k : k_u \geq j + n \}\), and \(\alpha(X, n) \to 0\) as \(n \to \infty\). • \(\rho'\) Coefficient: The coefficient \(\rho'(X^{(n)}, 1)\) is the supremum of \(\rho(A, B)\) over pairs of nonempty, disjoint sets \(Q, S \subset B(L_n)\) such that there exist \(u \in \{1, 2, \dots, d\}\) and \(j \in \mathbb{Z}\) with \(Q \subset \{ k : k_u \leq j \}\) and \(S \subset \{ k : k_u \geq j + 1 \}\).
theorem[Lindeberg-Feller Central Limit Theorem, feller1991introduction] Consider a triangular array of random variables \(\{X_{n,k} : n \in \mathbb{N}, k = 1, \dots, n\}\), where for each fixed \(n\), the variables \(X_{n,1}, \dots, X_{n,n}\) are independent, with \(\mathbb{E}[X_{n,k}] = 0\) and \(\mathrm{Var}(X_{n,k}) = \sigma_{n,k}^2 < \infty\). Let \(s_n^2 = \sum_{k=1}^n \sigma_{n,k}^2\), and assume \(s_n^2 > 0\). Suppose the Lindeberg condition holds: for every \(\varepsilon > 0\), \[ \lim_{n \to \infty} \frac{1}{s_n^2} \sum_{k=1}^n \mathbb{E}\left[ X_{n,k}^2 {\mathbb{I}}(\{|X_{n,k}| > \varepsilon s_n\}) \right] = 0. \] Then, as \(n \to \infty\), \[ \frac{1}{s_n} \sum_{k=1}^n X_{n,k} \xrightarrow{d} \mathcal{N}(0,1), \] where \(\xrightarrow{d}\) denotes convergence in distribution.

The next technical theorem we are relying on is Martingale Rosenthal Inequality, which is the Thereom 6.6.7 in de1999decoupling.

theorem[Martingale Rosenthal Inequality, de1999decoupling] Let $p \ge 2$. Let $M_n = \sum_{i=1}^n d_i$ where $\{M_n, \mathcal{F}_n\}$ is a martingale with martingale difference sequence $\{d_i\}$. Set $\|X\|_p = (\mathbb{E}|X|^p)^{1/p}$ for any random variable $X$, and define \[ U_{n,p} = \left( \sum_{j=1}^n \|d_j\|_p^p \right)^{1/p}, \qquad s_n(M) = \left( \sum_{i=1}^n \mathbb{E}(|d_i|^2 \mid \mathcal{F}_{i-1}) \right)^{1/2}. \] Then there exist constants $0 < c_p, C_p < \infty$ such that \[ c_p \big[ \| s_n(M) \|_p^p + U_{n,p}^p \big] \le \| M_n \|_p^p \le C_p \big[ \| s_n(M) \|_p^p + U_{n,p}^p \big]. \]

Another useful form of Bernstein inequality for reversible Markov Chain by Paulin2012ConcentrationIF is as following.

theorem[Bernstein Inequality, Theorem 3.9 by Paulin2012ConcentrationIF] Let $X_1, \ldots, X_n$ be a stationary reversible Markov chain with finite state space $\Omega$, absolute spectral gap $\gamma^*$ and the stationary distribution $\pi$. Let $f_1, f_2, \ldots, f_n \in L^2(\pi)$ be functions satisfying $\lvert f_i(x) - \mathbb{E}_\pi(f_i)\rvert \le C$ for every $x \in \Omega$. Let \[ S' := \sum_{i=1}^n f_i(X_i) \quad\text{and}\quad V_{S'} := \sum_{i=1}^n \mathrm{Var}_\pi(f_i). \] Then for every $n \ge 1$ and $t \ge 0$, we have \begin{equation} \mathbb{P}_\pi\!\left( \bigl|S' - \mathbb{E}_\pi(S')\bigr| \ge t \right) \le 2 \exp\!\left( -\,\frac{t^2 \cdot \bigl(2\gamma^* - (\gamma^*)^2\bigr)}{8V_{S'} + 20tC} \right). \end{equation}

The following theorem is useful to bounded the fourth moment for a strong mixing sequence established by rio2017asymptotic in their Theorem 2.1.

theorem[Fourth–moment bound under strong mixing {rio2017asymptotic}] Let $\{X_k\}_{k\ge1}$ be a sequence of centered real–valued random variables with $\mathbb E|X_k|^4<\infty$ for all $k$, and put $S_n:=\sum_{k=1}^n X_k$. Let $\alpha(m)$ denote the strong mixing coefficients of $\{X_k\}_{k\ge 1}$ and \[ \alpha^{-1}(u)\ :=\ \inf\{\,m\ge 0:\ \alpha(m)\le u\,\},\qquad u\in(0,1]. \] For each $k$, let $Q_k$ be the quantile function of $|X_k|$, i.e. \[ Q_k(u)\ :=\ \inf\{\,t\ge 0:\ \mathbb P(|X_k|>t)\le u\,\},\qquad u\in(0,1]. \] Define \[ M_{4,\alpha,n}\big(Q_k\big)\ :=\ \sum_{k=1}^n \int_0^1 \big[\alpha^{-1}(u)\wedge n\big]^3\, Q_k^4(u)\,du . \] Then, for every $n\ge1$, \[ \mathbb E\big(S_n^4\big)\ \le\ 3\Bigg(\sum_{i=1}^n\sum_{j=1}^n \big|\mathbb E(X_i X_j)\big|\Bigg)^{\!2} \;+\;48\, M_{4,\alpha,n}\big(Q_k\big). \]

Algorithm Details

In this section, Algorithm (ref) shows the details of the classical Edmonds--Karp used in Algorithm (ref). Algorithm (ref) relies on another standard routine of breath-first search, which is presented in Algorithm (ref). Algortihm (ref) is the routine of Eulerian Cycle Decomposition.

algorithm[algorithm omitted — 823 chars of source]
algorithm[algorithm omitted — 959 chars of source]
algorithm[algorithm omitted — 1,946 chars of source]

Proof to Lemma (ref)

proof{Proof to Lemma (ref).} We first consider the normal case where ${\mathcal{P}}_i$ is not a cycle or $j\not= k(i)$. In this case, for $j>1$, we can have \begin{align*} \mathbb{P}(W_{i,j}=1)&=\mathbb{P}(W_{i,j}=1, W_{i,j-1}=0)\\ &=\mathbb{P}(W_{i,j}=1\mid W_{i,j-1}=0)\cdot \mathbb{P}(W_{i,j-1}=0)\\ &=p\cdot(1-\mathbb{P}(W_{i,j-1}=1)), \end{align*} where the first equation due to the fact that $W_{i,j}=1$ implies $W_{i,j-1}=0$, and the last equation holds according to AP design. By simple calculation, we can know that $\mathbb{P}(W_{i,j}=1)-\frac{p}{p+1}=-p\cdot (\mathbb{P}(W_{i,j-1}=1)-\frac{p}{p+1})$ allowing us to use the recursion. With the initial value $\mathbb{P}(W_{i,1}=1)=\frac{p}{1+p}$, we have \begin{equation*} \mathbb{P}(W_{i,j}=1)=\frac{p}{p+1}. \end{equation*} For the last edge of a cycle ${\mathcal{P}}_i$ (i.e., $j=k(i)$), we can have \begin{align} \mathbb{P}(W_{i,j}=1)&=\mathbb{P}(W_{i,j}=1, W_{i,j-1}=0,W_{i,1}=0)\notag\\ &=\mathbb{P}(W_{i,j}=1\mid W_{i,j-1}=0,W_{i,1}=0)\cdot \mathbb{P}(W_{i,j-1}=0\mid W_{i,1}=0)\cdot \mathbb{P}( W_{i,1}=0)\notag \\ &=(1-\frac{p}{p+1})\cdot\mathbb{P}(W_{i,j-1}=0\mid W_{i,1}=0). \end{align} Now we only need to calculate $\mathbb{P}(W_{i,j-1}=0\mid W_{i,1}=0)$. \begin{align*} \mathbb{P}(W_{i,j-1}=0\mid W_{i,1}=0)&=1- \mathbb{P}(W_{i,j-1}=1\mid W_{i,1}=0)\\ &=1- \mathbb{P}(W_{i,j-1}=1,W_{i,j-2}=0\mid W_{i,1}=0)\\ &=1- \mathbb{P}(W_{i,j-1}=1\mid W_{i,j-2}=0, W_{i,1}=0)\mathbb{P}(W_{i,j-2}=0\mid W_{i,1}=0)\\ &=1- p\cdot \mathbb{P}(W_{i,j-2}=0\mid W_{i,1}=0). \end{align*} Again, by recursion, we can have \begin{equation} \mathbb{P}(W_{i,j-1}=0\mid W_{i,1}=0)-\frac{1}{p+1}=(-p)^{j-2} \left(\mathbb{P}(W_{i,1}=0\mid W_{i,1}=0)-\frac{1}{p+1}\right)=-\frac{(-p)^{j-1}}{p+1}, \end{equation} which leads to $\mathbb{P}(W_{i,j-1}=0\mid W_{i,1}=0)=\frac{1}{p+1}-\frac{(-p)^{j-1}}{p+1}$. By plugging this result back to \Eqref{eq:cycle-noncondition-prob}, we finish the prove. A side note here is that when $j=k(i)$ and ${\mathcal{P}}_i$ is a cycle, $j$ will always be larger than 4 and thus $j-2$ in the above equations is well defined. \qed

Proof to Theorem (ref)

proof{Proof to Theorem (ref)} A more compact way to write $\hat{\Gamma}_i$ is \begin{equation} \hat{\Gamma}_i=\sum_{j=1}^{k(i)}(-1)^{j+1}\frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}. \end{equation} For simplicity, we introduce the notation $\tilde{Y}_{v_{i,j},v_{i,j+1}}$ to denote $\frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}$. Based on \Eqref{eq:compact-estimator}, the variance of $\hat{\Gamma}_i$ can be decomposed as follows, \begin{align} \mathrm{Var}(\hat{\Gamma}_i) &= \mathrm{Var}\left(\sum_{j=1}^{k(i)}(-1)^{j+1} \tilde{Y}_{v_{i,j},v_{i,j+1}} \right)\notag\\ &=\sum_{j=1}^{k(i)} \mathrm{Var}\left(\tilde{Y}_{v_{i,j},v_{i,j+1}}\right)+2 \sum_{0\le j<q \le k(i)} (-1)^{j+q+2} \mathrm{Cov}\left(\tilde{Y}_{v_{i,j},v_{i,j+1}},\tilde{Y}_{v_{i,q},v_{i,q+1}}\right). \end{align} We start with the simpler case where ${\mathcal{P}}_i$ is a path. For the first term in \Eqref{eq:variance-decompose}, we can have \begin{align} \mathrm{Var}\left(\tilde{Y}_{v_{i,j},v_{i,j+1}}\right)&=\frac{Y_{v_{i,j},v_{i,j+1}}^2}{(\mathbb{P}(W_{i,j}=1))^2} \mathrm{Var}({\mathbb{I}}\{w_{i,j}=1\})\notag\\ &=\frac{1-\mathbb{P}(W_{i,j}=1)}{\mathbb{P}(W_{i,j}=1)} Y_{v_{i,j},v_{i,j+1}}^2\notag\\ &=\frac{1}{p}Y_{v_{i,j},v_{i,j+1}}^2, \end{align} where the second equation holds based on the variance of a binary distribution, and plugging in \Eqref{eq:uncondition-prob} leads to the last equation. For the second term in \Eqref{eq:variance-decompose}, we have the follows. \begin{align} \mathrm{Cov}\left(\tilde{Y}_{v_{i,j},v_{i,j+1}},\tilde{Y}_{v_{i,q},v_{i,q+1}}\right)&= \mathbb{E}\left[\frac{{\mathbb{I}}\{W_{i,j}=1\}{\mathbb{I}}\{W_{i,q}=1\} }{\mathbb{P}(W_{i,j}=1)\mathbb{P}(W_{i,q}=1)}{Y}_{v_{i,j},v_{i,j+1}}\cdot{Y}_{v_{i,q},v_{i,q+1}}\right]-{Y}_{v_{i,j},v_{i,j+1}}\cdot {Y}_{v_{i,q},v_{i,q+1}}\notag\\ \notag &=\left(\frac{\mathbb{P}(W_{i,j}=1,W_{i,q}=1)}{\mathbb{P}(W_{i,j}=1)\mathbb{P}(W_{i,q}=1)}-1\right) {Y}_{v_{i,j},v_{i,j+1}}\cdot {Y}_{v_{i,q},v_{i,q+1}}\\ &= \left(\frac{\mathbb{P}(W_{i,q}=1 \mid W_{i,j}=1)}{\mathbb{P}(W_{i,q}=1)}-1\right) {Y}_{v_{i,j},v_{i,j+1}}\cdot {Y}_{v_{i,q},v_{i,q+1}}. \end{align} Thus, we now only need to calculate $\mathbb{P}(W_{i,q}=1 \mid W_{i,j}=1)$. \begin{align*} \mathbb{P}(W_{i,q}=1 \mid W_{i,j}=1)&=\mathbb{P}(W_{i,q}=1, W_{i,q-1}=0 \mid W_{i,j}=1)\\ &=\mathbb{P}(W_{i,q}=1\mid W_{i,q-1}=0 ) \mathbb{P}(W_{i,q-1}=0\mid W_{i,j}=1) \\ &=p \cdot (1-\mathbb{P}(W_{i,q-1}=1\mid W_{i,j}=1)). \end{align*} With the above equation, we can resort to recursion that \begin{equation} \mathbb{P}(W_{i,q}=1 \mid W_{i,j}=1)=\frac{p}{p+1}+(-p)^{q-j}\left(\mathbb{P}(W_{i,j}=1 \mid W_{i,j}=1)-\frac{p}{p+1}\right)=\frac{p+(-p)^{q-j}}{p+1}. \end{equation} Putting \Eqref{eq:uncondition-prob} and \Eqref{eq:one-condition-one} back to \Eqref{eq:cov-last-step}, we can have \begin{align} \mathrm{Cov}\left(\tilde{Y}_{v_{i,j},v_{i,j+1}},\tilde{Y}_{v_{i,q},v_{i,q+1}}\right)=-(-p)^{q-j-1}{Y}_{v_{i,j},v_{i,j+1}}\cdot {Y}_{v_{i,q},v_{i,q+1}} \end{align} With plugging \Eqref{eq:variance-single-term} and \Eqref{eq:cov-single-term} to \Eqref{eq:variance-decompose}, we finish the proof for the path case. Under the cycle case, for $j=1,\cdots,k(i)-1$, $\mathrm{Var}(\tilde{Y}_{v_{i,j},v_{i,j+1}})$ still follows \Eqref{eq:variance-single-term}. Also, for any $1\le j<q\le k(i)-1$, the covariance $\mathrm{Cov}(\tilde{Y}_{v_{i,j},v_{i,j+1}},\tilde{Y}_{v_{i,q},v_{i,q+1}})$ are still of the same form in \Eqref{eq:cov-single-term}. What we need to additionally calculate is the term related to the last edge, i.e., $\mathrm{Var}(\tilde{Y}_{v_{i,k(i)},v_{i,k(i)+1}})$ and $\mathrm{Cov}(\tilde{Y}_{v_{i,j},v_{i,j+1}},\tilde{Y}_{v_{i,k(i)},v_{i,k(i)+1}})$ for $1\le j\le k(i)-1$. The variance part is similar as before with a different plug-in of the probability. \begin{align} \mathrm{Var}\left(\tilde{Y}_{v_{i,k(i)},v_{i,k(i)+1}}\right) &=\frac{1-\mathbb{P}(W_{i,k(i)}=1)}{\mathbb{P}(W_{i,k(i)}=1)} Y_{v_{i,k(i)},v_{i,k(i)+1}}^2\notag\\ &=\frac{p^2+2p+(-p)^{k(i)-1}}{1-(-p)^{k(i)-1}} Y_{v_{i,k(i)},v_{i,k(i)+1}}^2. \end{align} For the covariance term, \Eqref{eq:cov-last-step} still holds and the key challenge here is to calculate $\mathbb{P}(W_{i,k(i)}=1 \mid W_{i,j}=1)$ for $1\le j \le k(i)-1$. \begin{align*} \mathbb{P}(W_{i,k(i)}=1 \mid W_{i,j}=1) &=\mathbb{P}(W_{i,k(i)-1}=0, W_{i,1}=0\mid W_{i,j}=1)\\ &=\mathbb{P}(W_{i,k(i)-1}=0\mid W_{i,1}=0, W_{i,j}=1) \mathbb{P}( W_{i,1}=0 \mid W_{i,j}=1)\\ &=\mathbb{P}(W_{i,k(i)-1}=0\mid W_{i,j}=1) \mathbb{P}( W_{i,1}=0 \mid W_{i,j}=1)\\ &=\mathbb{P}(W_{i,k(i)-1}=0\mid W_{i,j}=1) \frac{\mathbb{P}( W_{i,j}=1 \mid W_{i,1}=0)\mathbb{P}( W_{i,1}=0)}{\mathbb{P}( W_{i,j}=1)}, \end{align*} where the first equation holds due to our AP design that $W_{i,k(i)}=1$ if and only if $W_{i,k(i)-1}=W_{i,1}=0$, the third equality is based on the fact that $W_{i,k(i)-1}$ is independent with $W_{i,1}$ conditioned on $W_{i,j}$, and the last equality adopts the Bayes' law. By \Eqref{eq:zero-to-zero} and \Eqref{eq:one-condition-one}, we can have \begin{align*} &\mathbb{P}(W_{i,k(i)-1}=0\mid W_{i,j}=1)=1- \mathbb{P}(W_{i,k(i)-1}=1\mid W_{i,j}=1) = \frac{1-(-p)^{k(i)-1-j}}{p+1},\\ &\mathbb{P}(W_{i,j}=1\mid W_{i,1}=0)=1- \mathbb{P}(W_{i,j}=0\mid W_{i,1}=0) = \frac{p+(-p)^{j}}{p+1}. \end{align*} Together with \Eqref{eq:uncondition-prob}, we can now have, \begin{equation} \mathbb{P}(W_{i,k(i)}=1 \mid W_{i,j}=1)= \frac{(1-(-p)^{j-1})(1-(-p)^{k(i)-1-j})}{(1+p)^2}. \end{equation} By plugging \Eqref{eq:last-edge-one-on-one} into \Eqref{eq:cov-last-step}, we can have \begin{align*} \mathrm{Cov}(\tilde{Y}_{v_{i,j},v_{i,j+1}},\tilde{Y}_{v_{i,k(i)},v_{i,k(i)+1}})=\left(\frac{(1-(-p)^{k(i)-1-j})(1-(-p)^{j-1})}{(1-(-p)^{k(i)-1})}-1\right)Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,k(i)},v_{i,k(i)+1}}. \end{align*} We finish the proof. \qed

Proof to Proposition (ref)

proof{Proof to Proposition (ref)} The core idea in this proof is to verify that all the coefficients in Theorem (ref) on ${\mathcal{P}}_i$ are positive. The only challenging term requiring checking is the covariance with the last edge, \begin{align*} &(-1)^j \left(\frac{(1-(-p)^{k(i)-1-j})(1-(-p)^{j-1})}{(1+p^{k(i)-1})}-1\right)\\ &= (-1)^j \frac{-(-p)^{j-1}-(-p)^{k(i)-1-j}+(-p)^{k(i)-2}-p^{k(i)-1}}{1+p^{k(i)-1}}\\ &= \frac{p^{j-1}+p^{k(i)-1-j}+(-1)^{j}p^{k(i)-2}+(-1)^{j+1}p^{k(i)-1}}{1+p^{k(i)-1}}\\ &\ge 0, \end{align*} where the last step holds since $p^{j-1}\ge -p^{k(i)-1}$ and $p^{k(i)-1-j}\ge -p^{k(i)-2}$. Therefore, with all the positive coefficients, the maximization problem \Eqref{eq:max-variance} is equivalent to maximizing $Y_{v_{i,j},v_{i,j+1}}^2$ and $Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}$ for all $i,j,q$. Thus, the optimal solution is $Y_{i,j}=B$ for all $(i,j)\in{\mathbb{M}}$. We finish the proof. \qed

Proof to Theorem (ref)

proof{Proof to Theorem (ref)} For the alternating path ${\mathcal{P}}_s=\{v_{s,1},\cdots,v_{s,k(s)+1}\}\in{\triangle\mathcal{P}^{(t,c)}}$, by plugging Proposition (ref) into \Eqref{eq:variance-of-path}, we can get \begin{align} \max_{\substack{ |Y_{i,j}|\le B \\ \forall (i,j)\in{\mathcal{P}}_s }} \mathrm{Var} (\hat{\Gamma}_s) &= \frac{k(s)}{p} B^2 + 2 \sum_{1\le j<q \le k(s)} p^{q-j-1} B^2\notag\\ &=\frac{k(s)}{p} B^2 + 2\sum_{1\le j < k(s)-1}p^{-j-1} B^2 \sum_{j<q\le k(s)} p^q\notag\\ &=\frac{k(s)}{p} B^2 + 2\sum_{1\le j < k(s)-1}p^{-j-1} B^2 \frac{p^{j+1}(1-p^{k(s)-j})}{1-p}\notag\\ &=\frac{k(s)}{p} B^2 + \frac{2B^2}{1-p}\sum_{1\le j < k(s)-1} (1-p^{k(s)-j})\notag\\ &=\frac{k(s)}{p} B^2 + \frac{2B^2}{1-p} \left(k(s)-1-\frac{p^{k(s)}-p}{p-1}\right). \end{align}

Proof to Theorem (ref)

First, we need to formally need the triangular array since for each given $N$, the matching can vary dramatically. Let us denote $Z_N:=N(\hat{\tau}_{\text{AP}}-\tau)$ and we have

equation[equation omitted — 152 chars of source]

where the new subscript $N$ is the index and recall that

align*[align* omitted — 707 chars of source]

We now partition ${\triangle\mathcal{P}^{(t,c)}}_N$ into giant sequences ${\mathcal{G}}_N$ and small sequences ${\mathcal{S}}_N$, where

align*[align* omitted — 329 chars of source]

For simplicity, we introduce the notation that $Z_{i,N}=\hat{\Gamma}_{i,N}- \Gamma_{i,N}$, $Z_{G,N}=\sum_{{\mathcal{P}}_{i,N} \in {\mathcal{G}}_N} Z_{i,N}$ and $Z_{S,N}=\sum_{{\mathcal{P}}_{i,N} \in {\mathcal{S}}_N} Z_{i,N}$. Then, we have $Z_N=Z_{S,N}+Z_{L,N}$, and by independence across all the paths, $\mathrm{Var}(Z_N)=\mathrm{Var}(Z_{S,N})+\mathrm{Var}(Z_{L,N})$. In addition, since $\mathbb{E}[\hat{\tau}_{\text{AP}}]=\tau$, we can have $\mathrm{Var}(Z_N)=\mathrm{Var}(N\hat{\tau}_{\text{AP}})$. Also, we can know that $|{\mathcal{G}}_N |\le (\mathrm{Var}(Z_N))^{\frac{1}{3}}$.

Now, we are going to show that if $\mathrm{Var}(Z_{G,N})/\mathrm{Var}(Z_N)\rightarrow \alpha$ (or $\mathrm{Var}(Z_{S,N})/\mathrm{Var}(Z_N)\rightarrow \alpha$) for any $\alpha>0$, $Z_{G,N}/\sqrt{\mathrm{Var}(Z_{G,N})}$ (correspondingly, $Z_{G,N}/\sqrt{\mathrm{Var}(Z_{G,N})}$) converges to the standard normal distribution. Note that this is still not enough, since $\mathrm{Var}(Z_{G,N})/\mathrm{Var}(Z_N)$ or $\mathrm{Var}(Z_{S,N})/\mathrm{Var}(Z_N)$ may not exist at all.

Case 1: If $\mathrm{Var}(Z_{S,N})/\mathrm{Var}(Z_N)\rightarrow \alpha$ for some $\alpha>0$. Let us consider the ${\mathcal{P}}_{i,N} \in {\mathcal{S}}_N$, and we denote ${\mathcal{P}}_{i,N}:=\{v_{i,N,1},\cdots,v_{i,N,k(i,N)}, v_{i,N,k(i,N)+1}\}$. WLOG, we assume $(v_{i,N,1},v_{i,N,2})\in{\triangle}{\mathbb{M}}^{(t,c)}_{t,N}$. With this notation, we can have

equation*[equation* omitted — 172 chars of source]

From \Eqref{eq:uncondition-prob}, we know that $\mathbb{P} (W_{i,j}=1)\ge \min\{p(1-p),(1-p)^2\}$. Together with our assumption of the bounded $Y_{v_{i,j},v_{i,j+1}}$, we can upper bound the absolute value of $\widetilde{Y}_{v_{i,j},v_{i,j+1}}:=(-1)^{j+1}\frac{Y_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}$ by $\frac{B}{\min\{p(1-p),(1-p)^2\}}$. For simplicity, we introduce the positive constant $C_{p,0}$ to denote that $|\widetilde{Y}_{v_{i,j},v_{i,j+1}}|\le C_{p,0}$ uniformly. We also adopt $\widetilde{W}_{i,j}$ to denote ${{\mathbb{I}}\{w_{i,j}=1\}}-\mathbb{P}(W_{i,j}=1)$, and thus $Z_{i,N}=\sum_{j=1}^{k(i,N)} \widetilde{Y}_{v_{i,j},v_{i,j+1}} \widetilde{W}_{i,j}$. We purposely split out the last edge of ${\mathcal{P}}_{i,N}$ because if ${\mathcal{P}}_i$ is a cycle, the last edge $W_{i,k(i,N)}$ does not satisfy the Markov property, and we define $Z_{i,N}^{(-1)}:=\sum_{j=1}^{k(i,N)-1} \widetilde{Y}_{v_{i,j},v_{i,j+1}} \widetilde{W}_{i,j}$.

Now we are going to verify the Lindeberg condition. We start with bounding the fourth moment of $Z_{i,N}$, i.e., $\mathbb{E}[Z_{i,N}^4]$. First, we still separate out the last edge out because $W_{i,1},\cdots,W_{i,k(i,N)-1}$ is always a Markov chain with the transition matrix: \[ {\bm{P}} =

pmatrix[pmatrix omitted — 33 chars of source]

, \] the eigenvalues of ${\bm{P}}$ are $1$ and $-p$. The absolute spectral gap $\gamma^*$ is $1-|-p|=1-p$. We can have the fast mixing properties. To be more specific, we have

align[align omitted — 497 chars of source]

where the inequality is based on the Cauthy-Schwarz inequality. For the first term, we are applying the Rosenthal-type inequality for the strong mixing sequence by rio2017asymptotic (see, also in Theorem (ref)). Define $Q_j(u)$ as the quantile function of $\widetilde{Y}_{v_{i,j},v_{i,j+1}}\widetilde{W}_{i,j}$ for $u\in(0,1)$. Since $|\widetilde{W}_{i,j}|\le 1$, $Q_k(u)\le|\widetilde{Y}_{v_{i,j},v_{i,j+1}}|$. As the chain is geometrically mixing, the strong mixing coefficients $\alpha(m)\le c_\alpha \rho^m$ with $\rho\in(0,1)$. Thus, there exists positive constants $a_{p}$ and $b_{p}$ only dependent on $p$ such that, for all $u\in(0,1)$,

equation*[equation* omitted — 75 chars of source]

where $\alpha^{-1}(u):=\inf\{m\ge 0: \alpha(m)\le u\}$. Following the notation in Theorem (ref),

align[align omitted — 550 chars of source]

For any $1\le \gamma\le\kappa \le k(i,N)-1$, we have

align*[align* omitted — 459 chars of source]

where $\sigma_W^2=p/(1+p)^2$. With this, we can further have,

align[align omitted — 751 chars of source]

Applying \Eqref{eq:M-upper-bound} and \Eqref{eq:cov-upper-bound} to Theorem (ref), we can have

align[align omitted — 1,029 chars of source]

where the last inequality again adopts the fact that all the correlation terms are positive. Putting \Eqref{eq:Z-1-4-bound} back to \Eqref{eq:fourth-moment-decompose}, we get

align*[align* omitted — 893 chars of source]

where $C_{p,2}:=(384 C_\alpha(p))\vee 8$. Therefore, for the set ${\mathcal{S}}_N$ under the condition that $\mathrm{Var}(Z_{S,N})/\mathrm{Var}(Z_N)\rightarrow \alpha$, we have for any $\epsilon>0$,

align*[align* omitted — 1,020 chars of source]

where the first inequality holds due to the fact that $\mathbb{E}[X^2 {\mathbb{I}}\{|X|>s\}]\le \mathbb{E}[X^4]/s^2$, the third inequality is because of the definition of ${\mathcal{S}}_N$, the second equality is based on the indepedence among all the paths, and the last step is because $\mathrm{Var}(Z_{S,N})/\mathrm{Var}(Z_N)\rightarrow \alpha$ and $\mathrm{Var}(Z_N)\rightarrow\infty$. Therefore, we can conclude that if the limit of $\mathrm{Var}(Z_{S,N})/\mathrm{Var}(Z_N)$ exists, we can guarantee that as $N\rightarrow\infty$

equation[equation omitted — 165 chars of source]

Case 2: If $\mathrm{Var}(Z_{G,N})/\mathrm{Var}(Z_N)\rightarrow \beta$ for some $\beta>0$. Let us consider the ${\mathcal{P}}_{i,N} \in {\mathcal{G}}_N$, and we denote ${\mathcal{P}}_{i,N}:=\{v_{i,N,1},\cdots,v_{i,N,k(i,N)}, v_{i,N,k(i,N)+1}\}$. WLOG, we assume $(v_{i,N,1},v_{i,N,2})\in{\triangle}{\mathbb{M}}^{(t,c)}_{t,N}$. Following the same notation in Case 1, we can have

equation*[equation* omitted — 108 chars of source]

Then, we can get

equation[equation omitted — 379 chars of source]

Let us first focus on the second term of \Eqref{eq:giant-variance}. Based on the definition of ${\mathcal{G}}_N$, we can guarantee that $|{\mathcal{G}}_N|\le (\mathrm{Var}(Z_N))^{\frac{1}{3}}$. As the assignment across different paths is independent from each other,

align[align omitted — 495 chars of source]

where the last inequality holds due to the facts that $\mathrm{Var}(\widetilde{W}_{i,j})\le \frac{1}{4}$ and $|{\mathcal{G}}_N|\le (\mathrm{Var}(Z_N))^{\frac{1}{3}}$. Therefore, we can have that the second term of \Eqref{eq:giant-variance} will converge to zero as $\mathrm{Var}(Z_{G,N})/\mathrm{Var}(Z_N)\rightarrow \beta$, i.e.,

equation[equation omitted — 345 chars of source]

For the first term in \Eqref{eq:giant-variance}, we are applying the non-stationary $\alpha$-mixing CLT by bradley2017central (see, also Theorem (ref)). For ${\mathcal{P}}_{i,N}$, we define ${\mathcal{C}}_{i,N}^{(-1)}$ as $W_{i,1},\cdots,W_{i,k(i,N)-1}$, which is a Markov chain. We then concatenate the ${\mathcal{C}}_{i,N}^{(-1)}$ for every ${\mathcal{P}}_{i,N}$ in ${\mathcal{G}}_N$ together and get a long chain ${\mathcal{C}}_N$. For the chain ${\mathcal{C}}_N$, in general, it is not a Markov Chain anymore since between the paths, they are independence. However, the independence still preserves an $\alpha$-mixing sequence. The maximal correlation coefficient just decreases when the two $\sigma$-fields are cross the paths. In order to make the Lindeberg condition hold, we now show that $\mathrm{Var}(\sum_{{\mathcal{P}}_{i,N}\in {\mathcal{G}}_N}Z_{i,N}^{(-1)})$ goes to infinity. Consider that

align[align omitted — 936 chars of source]

where the third equality holds due to \Eqref{eq:tail-variance-bound} and $\mathrm{Var}(Z_{G,N})/\mathrm{Var}(Z_N)\rightarrow \beta$ and the last equality is based on the observation that

align*[align* omitted — 418 chars of source]

and then the limit goes to zero. Thus, \Eqref{eq:variance-variance-wo-last-edge} indicates that $\mathrm{Var}(\sum_{{\mathcal{P}}_{i,N}\in {\mathcal{G}}_N}Z_{i,N}^{(-1)})$ goes to infinity. Moreover, since $\widetilde{Y}_{v_{i,j},v_{i,j+1}}\widetilde{W}_{i,j}$ is bounded by a constant, the Lindeberg condition easily holds. Therefore, we know that

equation*[equation* omitted — 204 chars of source]

By \Eqref{eq:last-edge-small}, \Eqref{eq:variance-variance-wo-last-edge}, and Slutsky's theorem, we can have that as $N\rightarrow\infty$,

equation[equation omitted — 165 chars of source]

Now, we are finalizing the proof by establishing the result for the case when $\alpha$ and $\beta$ may not exist. Define the following notations,

equation*[equation* omitted — 145 chars of source]

It is natural to see $\alpha_N+\beta_N=1$ and both of them are in $[0,1]$. Thus, $(\alpha_N,\beta_N)$ always lies in a compact set. By the Bolzano-Weierstrass theorem, for every subsequence $(\alpha_{N_k},\beta_{N_k})$ of $(\alpha_N,\beta_N)$, there exists a further subsubsequence $(\alpha_{N_{k_j}},\beta_{N_{k_j}})$ converging to some $(\alpha,\beta)$ with $\alpha+\beta=1$. Therefore, by Case 1 (\Eqref{eq:case-1-CLT}) and Case 2 (\Eqref{eq:case-2-CLT}), we can establish, along the subsubsequence $(\alpha_{N_{k_j}},\beta_{N_{k_j}})$,

equation*[equation* omitted — 292 chars of source]

where $Z_1$ and $Z_2$ are i.i.d. normal. By the continuous mapping theorem,

equation*[equation* omitted — 367 chars of source]

and we know that $\sqrt{\alpha_{N_{k_j}}} Z_1+\sqrt{\beta_{N_{k_j}}} Z_2 \sim {\mathcal{N}}(0,\alpha_{N_{k_j}}+\beta_{N_{k_j}})$ which is ${\mathcal{N}}(0,1)$. Thus, we are same to draw the conclusion that for any subsequence of $(\alpha_N,\beta_N)$, we can always have a further subsubsequence, along which our estimator converge to the same distribution ${\mathcal{N}}(0,1)$. By the subsequence principle (see, Theorem 2.6. in billingsley2013convergence), we can guarantee that if $\mathrm{Var}(Z_N)$ goes to infinity,

equation*[equation* omitted — 93 chars of source]

We finish the proof. \qed

Technical Details to Proposition (ref)

propositionFor ${\mathcal{P}}_i=\{v_{i,1},\cdots,v_{i,k(i)+1}\}\in{\triangle\mathcal{P}^{(t,c)}}$, if ${\mathcal{P}}_i$ is a path, we have the following upper bound $\tilde{\sigma}^2_i$ of $\mathrm{Var}(\hat{\Gamma}_i)$, \begin{equation} \tilde{\sigma}^2_i := \left(\frac{1}{p}+1\right) \sum_{j=1}^{k(i)} Y_{v_{i,j},v_{i,j+1}}^2+\sum_{j=2}^{k(i)-1} Y_{v_{i,j},v_{i,j+1}}^2 + 2 \sum_{\substack{1\le j \le k(i)-2 \\ j+2\le q \le k(i)}} p^{q-j-1} Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}. \end{equation} Moreover, $\tilde{\sigma}^2_i$ can be estimated by \begin{align*} \hat{\sigma}^2_i&:= \left(\frac{1}{p}+1\right) \sum_{j=1}^{k(i)} \frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y^2_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}+\sum_{j=2}^{k(i)-1} \frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y^2_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)}\\ &\qquad\qquad\qquad + 2 \sum_{\substack{1\le j \le k(i)-2 \\ j+2\le q \le k(i)}} p^{q-j-1} \frac{Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}{\mathbb{I}}\{w_{i,j}=1,w_{i,q}=1\}}{\mathbb{P}(W_{i,j}=1,W_{i,q}=1)}, \end{align*} where $\mathbb{P}(W_{i,j}=1)=p/(p+1)$ and $\mathbb{P}(W_{i,j}=1,W_{i,q}=1)=(p^2-(-p)^{q-j+1})/(p+1)^2$. Similarly, if ${\mathcal{P}}_i$ is a cycle, the upper bound $\tilde{\sigma}^2_i$ of $\mathrm{Var}(\hat{\Gamma}_i)$ is \begin{align*} \tilde{\sigma}^2_i :=& \left(\frac{1}{p}+2\right) \sum_{j=1}^{k(i)-1} Y_{v_{i,j},v_{i,j+1}}^2+ 2 \sum_{\substack{1\le j \le k(i)-3 \\ j+2\le q \le k(i)-1}} p^{q-j-1} Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}\\ &+ \left(\frac{p^2+2p-p^{k(i)-1}}{1+p^{k(i)-1}}+2\right) Y_{v_{i,k(i)},v_{i,k(i)+1}}^2 \\ &+ 2 \sum_{2\le j \le k(i)-2} (-1)^j \left(\frac{(1-(-p)^{k(i)-1-j})(1-(-p)^{j-1})}{(1+p^{k(i)-1})}-1\right)Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,k(i)},v_{i,k(i)+1}}. \end{align*} We can estimate the upper bound by \begin{align*} \hat{\sigma}^2_i&:= \left(\frac{1}{p}+2\right) \sum_{j=1}^{k(i)-1} \frac{{\mathbb{I}}\{w_{i,j}=1\}\cdot Y^2_{v_{i,j},v_{i,j+1}}}{\mathbb{P}(W_{i,j}=1)} + 2 \sum_{\substack{1\le j \le k(i)-3 \\ j+2\le q \le k(i)-1}} p^{q-j-1} \frac{Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,q},v_{i,q+1}}{\mathbb{I}}\{w_{i,j}=1,w_{i,q}=1\}}{\mathbb{P}(W_{i,j}=1,W_{i,q}=1)}\\ &+\left(\frac{p^2+2p-p^{k(i)-1}}{1+p^{k(i)-1}}+2\right) \frac{{\mathbb{I}}\{w_{i,k(i)}=1\} Y_{v_{i,k(i)},v_{i,k(i)+1}}^2}{\mathbb{P}(W_{i,k(i)}=1)}\\ &+ 2 \sum_{2\le j \le k(i)-2} (-1)^j \left(\frac{(1-(-p)^{k(i)-1-j})(1-(-p)^{j-1})}{(1+p^{k(i)-1})}-1\right)\frac{{\mathbb{I}}\{w_{i,j}=w_{i,k(i)}=1\}Y_{v_{i,j},v_{i,j+1}}Y_{v_{i,k(i)},v_{i,k(i)+1}}}{\mathbb{P}(W_{i,j}=1,W_{i,k(i)}=1)}, \end{align*} where $\mathbb{P}(W_{i,j}=1)=p/(p+1)$ for $1\le j \le k(i)-1$, $\mathbb{P}(W_{i,k(i)}=1)=(1-(-p)^{k(i)-1})/(1+p)^2$, $\mathbb{P}(W_{i,j}=1,W_{i,q}=1)=(p^2-(-p)^{q-j+1})/(p+1)^2$ for $q\not= k(i)$ and $\mathbb{P}(W_{i,j}=1,W_{i,k(i)}=1)=\frac{p(1-(-p)^{j-1})(1-(-p)^{k(i)-1-j})}{(1+p)^3}$. Furthermore, the estimator $\hat{\sigma}^2_i$ is an unbiased estimator of $\tilde{\sigma}^2_i$, i.e., $\mathbb{E}[\hat{\sigma}^2_i]=\tilde{\sigma}^2_i$.