Random Transpositions on Contingency Tables

This section describes the relationship between contingency tables and double cosets of $S_{n}$. Let $n$ be a positive integer and suppose $\lambda=(\lambda_{1},\lambda_{2},\dots,\lambda_{I})$ is a partition of $n$. That is, $\lambda_{i}$ are all positive integers, $\lambda_{1}\geq\lambda_{2}\geq\ldots\geq\lambda_{I}>0$, and $\lambda_{1}+\lambda_{2}+\ldots+\lambda_{I}=n$.

If $G$ is an arbitrary finite group and $H,K$ are two sub-groups of $G$, then the $H-K$ double-cosets are the equivalence classes defined by the equivalence relation

$$s\sim t\iff s=htk\,\,\,\,\,\,\,\,\text{for}\,\,\,\,\,\,\,\,s,t\in G,\,\,\,h\in H,\,\,\,k\in K.$$

Let $\mu=(\mu_{1},\dots,\mu_{J})$ be a second partition of $n$, let $\mathcal{T}_{\lambda,\mu}$ denote the space of contingency tables with row sums given by $\lambda$ and column sums $\mu$. That is, an element of $\mathcal{T}_{\lambda,\mu}$ is an $I\times J$ table with non-negative integer entries, row sums $\lambda_{1},\dots,\lambda_{I}$ and column sums $\mu_{1},\dots,\mu_{J}$. The following bijection is described in further detail in [29], and the induced distribution, along with many more references, in [16].

The mapping from $S_{n}$ to tables is easy to describe: Fix $\sigma\in S_{n}$. Inspect the first $\lambda_{1}$ positions in $\sigma$. Let $T_{11}$ be the number of elements from $\{1,2,\dots,\mu_{1}\}$ occurring in these positions, $T_{12}$ the number of elements from $\{\mu_{1}+1,\dots,\mu_{1}+\mu_{2}\}$, …and $T_{1J}$ the number of elements from $\{n-\mu_{J}+1,\dots,n\}$. In general, $T_{ij}$ is the number of elements from $\{\mu_{1}+\ldots+\mu_{i-1}+1,\dots,\mu_{1}+\ldots+\mu_{j}\}$ which occur in the positions $\lambda_{1}+\lambda_{2}+\ldots+\lambda_{i-1}+1$ up to $\lambda_{1}+\ldots+\lambda_{i}$.

Thanks to Zhihan Li, another way to describe the mapping from $S_{n}$ to contingency tables is to consider the $n\times n$ permutation matrix defined by $\sigma$. The partition $\mu$ divides the columns into $J$ blocks and $\lambda$ divides the rows into $I$ blocks. The table $T$ corresponding to $\sigma$ is defined by $T_{ij}$ as the number of $1$s in the $i$ block of rows and $j$ block of columns. For example, $\sigma=12543$ gives the permutation matrix

which defines the table $\begin{pmatrix}2&0&1\cr 0&2&0\end{pmatrix}$.

The random transpositions chain on $S_{n}$ is easy to describe. If a permutation $\sigma\in S_{n}$ is viewed as an arrangement of a deck of $n$ unique cards, the random transpositions chain can be described: Pick one card with your left hand and one card with your right hand (allowing for the possibility of picking the same card), and swap the two cards. The transition probabilities are generated by the probability measure $Q$ on $S_{n}$ defined

$$Q(\sigma)=\begin{cases}\frac{2}{n^{2}}&\text{if}\,\,\,\ \sigma=(i,j),i<j\\ \frac{1}{n}&\text{if}\,\,\,\ \sigma=id\\ 0&\text{otherwise},\end{cases}$$

(3)

where $(i,j)$ denotes the transposition of elements $i$ and $j$. The transition matrix for the random transposition Markov chain is then

$$P(x,y)=Q(y\cdot x^{-1}),\,\,\,\,\,\,\,\,x,y\in S_{n}.$$

Note that $Q$ is constant on conjugacy classes: If $i<j,k<\ell$ and $x=(ij),y=(k\ell)$ are two transpositions, then

$\displaystyle yxy^{-1}=(k\ell)\cdot(ij)\cdot(k\ell)=\begin{cases}(ij)&\text{if}\,\,i\neq k,j\neq\ell\,\,\text{or}\,\,i=k,j=\ell\\ (j\ell)&\text{if}\,\,i=k,j\neq\ell\\ (ki)&\text{if}\,\,j=\ell,i\neq k\end{cases}.$

Since transpositions generate $S_{n}$, this proves that if $x$ is a transposition, then $yxy^{-1}$ is a transposition for any $y\in S_{n}$, and thus $Q(yxy^{-1})=Q(x)$.

The random transpositions chain on $S_{n}$ induces a Markov chain on $\mathcal{T}_{\lambda,\mu}$ as a ‘lumping’. From a table $T\in\mathcal{T}_{\lambda,\mu}$, choose a permutation $x\in S_{n}$ that maps to the double coset corresponding to $T$. From $x$, sample $y$ from $P(x,\cdot)$ and move to the contingency table corresponding to $y$. More details of this lumping are discussed in Section 2.1, and this perspective will be useful for relating the well-studied eigenvalues of the chain on $S_{n}$ to the chain on $\mathcal{T}_{\lambda,\mu}$. It is also possible to write down the transition probabilities independently of the chain on $S_{n}$.

Note if $T_{i_{1},j_{1}}=0$ or $T_{i_{2},j_{2}}=0$, then $F_{(i_{1},j_{1}),(i_{2},j_{2})}(T)$ has negative values and is not an element of $\mathcal{T}_{\lambda,\mu}$. Definition 1.4 does not allow selecting this possibility, since only pairs with $T_{i_{1},j_{1}}T_{i_{2},j_{2}}\geq 1$ will be chosen.

There are other ways of thinking about the Markov chain. Suppose the contingency table is created from $n$ data pairs of the type $(i,j)$, where $i$ indicates the row feature and $j$ indicates the column feature. One step of the chain is to pick two pairs uniformly (possibly picking the same pair) and swapping the column values (or equivalently, swapping the row values). This is exactly the swap Markov chain on bi-partite graphs, which has been studied for various special cases, e.g. [2].

The interpretation of the chain as a lumping of the random transpositions chain on $S_{n}$ is advantageous because the chain on contingency tables inherits its eigenvalues from the chain on $S_{n}$, and these eigenvalues are well-known. Representation theory also gives an expression for the multiplicities of the eigenvalues in terms of Kostka numbers corresponding to partitions. Section 2.1 reviews this background in more detail, and thoroughly explains part (b) of the following theorem, which are the results for eigenvalues and multiplicities.

Theorem 1.6 (a) and (b) arise from the general theory for Markov chains on double cosets developed in [12] and [18] and summarized in Section 2.1. Section 3 contains a proof of Theorem 1.6(c) as well as a description of the multiplicities. It seems unlikely that there would be a simple formula for the number and multiplicity of the eigenvalues, as this would give a way of enumerating the size of the state space $\mathcal{T}_{\lambda,\mu}$, which is a #- P complete problem [11]. The orthogonal polynomials of the Fisher-Yates distribution have not been explicitly computed. Section 3.2 computes linear and quadratic eigenfunctions of $P$, for any size table, which can be used as a starting point for finding orthogonal polynomials.

For $2\times J$ tables, the stationary distribution is multi-variate hypergeometric, and the computation of the eigenvalues simplifies to give a more complete picture.

The eigenvalues and eigenvectors allow an analysis of the convergence rates of this Markov chain. The measure that we study for convergence rate is the mixing time. That is, for a Markov chain on a space $\Omega$ with transition probability $P$ and stationary distribution $\pi$, the mixing time is defined as

$$t_{mix}(\epsilon)=\sup_{x_{0}\in\Omega}\inf\{t>0:\|P^{t}(x_{0},\cdot)-\pi(\cdot)\|_{TV}<\epsilon\},\quad\epsilon>0,$$

where $\|\mu-\pi\|_{TV}=\sup_{A\subset\Omega}|\mu(A)-\pi(A)|$ is the total-variation distance between probability measures.

It is challenging to analyze the mixing time in full generality, but understanding the linear and quadratic eigenfunctions of the chain give the following results:

• •

  For $\lambda=(n-k,k),\mu=(n-\ell,\ell)$, the time to stationarity, averaged over starting positions sampled from $\pi_{\lambda,\mu}$, is bounded above by $(n/4)\log(\min(k,\ell))$.

• •

  For $\lambda=(n-k,k)$ and any $\mu$, the multivariate Hahn polynomials are used to compute an upper bound for the time to stationarity starting from extreme states in which the second row has a single non-zero entry.

• •

  For any $\lambda,\mu$, the mixing time is bounded below by $((n/4)-1)\log(C_{\lambda,\mu}n)$, where $C_{\lambda,\mu}$ is a computable constant.

The mixing time of the random transpositions chain on $S_{n}$ is $(1/2)n\log(n)+cn$, which is an upper bound for the mixing time of the lumped chain on $\mathcal{T}_{\lambda,\mu}$. The results listed above suggest that there is not much of a speed-up.

Note that for two way contingency tables it is easy to sample directly from the Fisher-Yates distribution (by generating a random permutation and using the argument of Lemma 1.2). One reason for interest in the Markov chain mixing time is that, for three and higher way tables, the Markov chain method is the only available approach, so it is important to see how it works in the two way case.

The uniform swap Markov chain on contingency tables was first proposed in [11], with the motivation of studying the uniform distribution [10]. The mixing time of this chain has been analyzed in special cases, e.g.  if $n=\sum\lambda_{i}=\sum\mu_{j}$ and the number of rows and columns is fixed, then the mixing time is of order $n^{2}$. When the table has two rows, the chain mixes in time polynomial in the number of columns and $n$ [25]. Chung, Graham, and Yau [7] indicate that a modified version of the chain converges in time polynomial in the number of rows and columns.

A contingency table can be thought of as a bi-partite graph with prescribed degree sequences, and a similar swap Markov chain defined for any graph with fixed degree sequences to sample from the uniform distribution, e.g.  [2]. There is a large literature on methods for sampling graphs with given degree sequences (not necessarily bi-partite); see [20] for a recent review. Recently the chain was studied for contingency tables with entries in $\mathbb{Z}/q\mathbb{Z}$ and mixing time with cut-off proven at $c_{q}n^{2}\log(n)$ [41].

A Markov chain on $\mathcal{T}_{\lambda,\mu}$ for which the stationary distribution is Fisher-Yates was studied in [17], motivated by studying conditional distributions. The Markov chain is created by taking the Gibbs sampler using the uniform swap chain defined in [11]. This gives the chain: Pick a pair of rows $i_{1},i_{2}$ and columns $j_{1},j_{2}$ uniformly at random. This determines a $2\times 2$ submatrix. Replace it with a $2\times 2$ sub-table with the same margins, chosen from the hypergeometric distribution. As noted in Section 2 of [17], while it is straightforward to sample directly from Fisher-Yates distribution for $2$-way tables, there is not a simple algorithm for sampling from the analagous distribution on $3$-way or higher-dimensional tables. Thus, Markov chains with Fisher-Yates distribution as stationary are valuable for this application; multiway tables are discussed further in Section 5.2.

Section 2 contains an overview of basic definitions and facts for Markov chains, double cosets, and orthogonal polynomials. Section 3 establishes the eigenvalues and multiplicities, as well as the eigenfunctions of the Markov chain as polynomials, and gives formulas for the linear and quadratic polynomial eigenfunctions. These are then used in Section 4 to find explicit upper and lower bounds on the mixing time of the chain, in specific cases. Section 5 discusses some future directions.

The results of this section use basic tools of representation theory; in particular, induced representations and Frobenius reciprocity. For background on these topics, see [27], [30], [47]. In the specific case in hand – the symmetric group – additional tools are available; Young’s rule and the hook-length formula. These are clearly developed in [31]. See also, [29].

Let $H,K$ be subgroups of the finite group $G$, $Q$ a probability on $G$. If $\text{supp}(Q)$ is not contained in a coset of a subgroup, then the random walk on $G$ induced by $Q$ is ergodic with a uniform stationary distribution. This random walk on $G$ is defined by picking a random element from $Q$ and multiplying the current state. That is, transitions are

$$P(x,y)=Q(yx^{-1}).$$

See [8] for an introduction to random walks on groups.

The double cosets of $H\backslash G/K$ partition the space $G$ and any Markov chain on $G$ defines a random process on the set of double cosets by keeping track of which double coset the process on $G$ is in, giving a ‘lumped’ process. Proposition 2.1 gives a condition for when the random walk induced by $Q$ on $G$ is also a Markov chain on the double cosets. For the statement we fix double coset representatives and write $x$ for the double coset $HxK$. The result is proven in [12].

Since transpositions form a conjugacy class in $S_{n}$, the probability $Q$ from (3) satisfies Proposition 2.1, so lumped to contingency tables gives a Markov process. The following lemma also directly proves that the lumping of the random transpositions chain to contingency tables is Markov, and equivalent to Definition 1.4.

The driving Markov chain on $S_{n}$ for contingency tables is the random transpositions Markov chain studied in [14]. This uses the tools of representation theory and the character theory of the symmetric group. An expository account, aimed at probabilists and statisticians is in Chapter 3D in [8]. One of the main results used below is an explicit determination of the eigenvalues of this chain.

Recall the probability measure which defines the random walk on $S_{n}$:

$$Q(\sigma)=\begin{cases}\frac{2}{n^{2}}&\text{if}\,\,\,\,\sigma=(i,j),i<j\\ \frac{1}{n}&\text{if}\,\,\,\,\sigma=id,\\ 0&\text{otherwise},\end{cases}.$$

This measure $Q$ is not concentrated on a single conjugacy class (the identity is not in the conjugacy class of transpositions). However, we can write

$$Q=\frac{1}{n}I+\frac{n-1}{n}\widetilde{Q},$$

where $I(\sigma)=\mathbf{1}(\sigma=id)$ and $\widetilde{Q}=\frac{1}{\binom{n}{2}}\mathbf{1}(\sigma\in\mathcal{C})$, where $\mathcal{C}\subset S_{n}$ denotes the conjugacy class of transpositions. (This $\widetilde{Q}$ is the same as the one discussed in Remark 3.) Then $\widetilde{Q}$ is concentration on a single conjugacy class $\mathcal{C}$, so the eigenvalues of the random walk are equal to the character ratio,

$$\widetilde{\beta}_{\rho}=\frac{1}{\chi_{\rho}(1)}\sum_{s\in G}Q(s)=\frac{\chi_{\rho}(\mathcal{C})}{\chi_{\rho}(1)},$$

for each partition $\rho$ of $n$.. The formula for this character ratio was determined by Frobenius; see [38] for a modern exposition, and [44] for a proof of this special case:

$$\widetilde{\beta}_{\rho}=\frac{1}{n(n-1)}\sum_{j=1}^{k}\left[\rho_{j}^{2}-(2j-1)\rho_{j}\right].$$

The eigenvalues for the walk driven by $Q$ are then related:

$$\beta_{\rho}=\frac{1}{n}+\frac{n-1}{n}\widetilde{\beta}_{\rho}.$$

This information is summarized in the following lemma; see [14] for more details.

In [14], the chain is proven to mix in total variation distance, with cut-off, after $(n/2)\log(n)$. This gives an initial upper bound on the mixing time of the random transpositions chain on contingency tables.

This section reviews orthogonal polynomials and records the formula for multivariate Hahn polynomials. See [19] for a thorough exposition on multivariate orthogonal polynomials; especially Section 3 for sufficient conditions for an orthonormal basis to exists for a probability distribution.

Let $\pi$ be a probability distribution on a finite space $\Omega\subset\mathbb{N}^{d}$. Let $\ell^{2}(\pi)$ be the space of functions $f:\Omega\to\mathbb{R}$ with the inner product

$$\left\langle f,g\right\rangle_{\pi}=\operatorname{{\mathbb{E}}}_{\pi}[f(X)g(X)]=\sum_{x\in\Omega}f(x)g(x)\pi(x).$$

A set of functions $\{q_{m}\}_{0\leq m<|\Omega|}$ are orthogonal in $\ell^{2}(\pi)$ if

$$\left\langle q_{m},q_{\ell}\right\rangle_{\pi}=d_{m}^{2}\textbf{1}(m=\ell).$$

For the following lemma, $\mathbf{m}=(m_{1},\dots,m_{N})$ denotes an index vector and $|\mathbf{m}|=\sum m_{i}$ is the total degree of the polynomial defined by the vector.