A Nonlinear Transform for the Diagonalization of the Bernoulli-Laplace Diffusion Model and Orthogonal Polynomials

TL;DR

This paper introduces a novel method for diagonalizing the transition matrix of the Bernoulli-Laplace diffusion model, enabling explicit calculation of eigenvalues and eigenvectors through a mapping to a differential operator, simplifying analysis of such Markov processes.

Contribution

A new approach mapping the eigenproblem to a differential operator, making the diagonalization of Bernoulli-Laplace and similar models more straightforward and broadly applicable.

Findings

Explicit eigenvalues and eigenvectors derived for the Bernoulli-Laplace model.

Method applicable to other two-urn diffusion models.

Simplifies analysis compared to previous orthogonal polynomial methods.

Abstract

The Bernoulli-Laplace model describes a diffusion process of two types of particles between two urns. To analyze the finite-size dynamics of this process, and for other constructive results we diagonalize the corresponding transition matrix and calculate explicitly closed-form expressions for all eigenvalues and eigenvectors of the Markov transition matrix $T_{BL}$. This is done by a new method based on mapping the eigenproblem for $T_{BL}$ to the associated problem for a linear partial differential operator $L_{BL}$ acting on the vector space of homogeneous polynomials in three indeterminates. The method is applicable to other Two Urns models and is relatively easy to use compared to previous methods based on orthogonal polynomials or group representations.

Equations

$$a_j^{(m+1)}$$

$$p_j$$

$$(x-y{)}^2{G}_{xy}=N(N-1)(\lambda -1)G,$$

$$N_w+N_b=N=N_1+N_2$$

$$p_i$$

$$q_i$$

$$r_i$$

$$L^{\prime }[H(u,v)]=\lambda (N)H(u,v)$$

$$ST{S}^{-1}=D$$

$$H(u,v,w)=H(f(x,y,z),g(x,y,z),h(x,y,z))=G(x,y,z)$$

$$L^{\prime }[H(u,v,w)]=\lambda (N)H(u,v,w)$$

$$G^{(k)}(x,y,z)=\sum _{i=0}^{N_w}{c}_k(i)x^i{y}^{N_1-i}{z}^{N_w-i}.$$

$$L_p=$$

$$L_q=$$

$$L_{BL}[G^{(k)}]$$

$$u=x-yz,y=y,z=z.$$

$$H^{(k)}(u,y,z)=G^{(k)}(x(u,yz),y,z)=\sum _i{b}_i^k{u}^i{y}^{N_1-i}{z}^{N_w-i}$$

$${\partial }_x$$

$${\partial }_y$$

$${\partial }_z$$

$${\partial }_{xy}$$

$${\partial }_{yz}$$

$$N_1{N}_2({\lambda }_k-1)H^{(k)}$$

$$L_{BL}^{\prime }[H^{(k)}]$$

$$=$$

$$N_1{N}_2({\lambda }_k-1)b_i^k$$

$$=$$

$$b_i^k=\frac{\left(N_1-i+1\right)\left(N_w-i+1\right)b_{i-1}^k}{N_1{N}_2({\lambda }_k-1)+i\left(N_w-i\right)+(N_b+1)i}$$

$${\lambda }_k$$

$$=$$

$${\lambda }_0$$

$${\lambda }_1$$

$$b_i^k$$

$$=$$

$$H^{(k)}=\sum _{i=k}^{N^{\prime }}(-1{)}^{i-k}\frac{(k-N_1{)}_{i-k}(k-N_w{)}_{i-k}}{(i-k)!(2k-N{)}_{i-k}}{u}^i{y}^{N_1-i}{z}^{N_w-i}$$

Full text

1 Introduction

The Bernoulli-Laplace (BL) model arise from diffusion theory and is related to the shuffling of cards [8]. Symmetries of the permutation group $S_{N}$ appear naturally in this model and other random walks on groups. Previous solutions of this model have appeared in Diaconis and Shashahani [7] and in the works of Karlin and MacGregor [10]. Group representations are used explicitly in the first; the derivation of a non-standard inner product or equivalently a measure for orthogonal polynomials which are related to the eigensolutions of the BL model appears in the second.

In this paper, we give a third way for deriving exact solutions of all the eigenvectors of the BL model, through a nonlinear transform that triangularizes and then diagonalizes the transition matrix $T_{BL}$. In brief, our method associates a specific linear partial differential operator (LPDO) $L_{BL}$ that acts on the vector space of homogeneous polynomials, $G,$ to the matrix $T_{BL}$. The $L_{BL}$ inherits the symmetries of the BL model; it encodes the tri-diagonal singly-stochastic (column sums are equal to 1 ) and anti-symmetric structure of $T_{BL}$. The components of the (right) eigenvectors (in view of the equal column sums of $T_{BL})$ of $T_{BL}$ is encoded in the coefficients of the homogeneous polynomial $G.$ A classical theory for the symmetries of such LPDOs have been formulated in terms of the Lie algebra of symmetry operators $K$ that commutes with $L_{BL}$ (cf. [13]).

It turns out and we exploit in our method, that the symmetries of $L_{BL}$ appear iu the form of suitable linear and nonlinear transformations $P$ on the independent variables $x,y,$ etc. or indeterminates of $G.$ The ease of use of this method resides in the transparent or explicit way to find these transformations $P.$ Our algorithm is completed by associating the transformed LPDO, $L_{BL}^{\prime},$back to what turns out to be a triangular matrix $T_{BL}^{\prime}$; in other words, the transformation $P$ for $L_{BL}$ encodes a similarity transformation that triangularizes $T_{BL},$ i.e., $PT_{BL}P^{-1}=T_{BL}^{\prime}$, which is then solved directly for its eigenvalues and right eigenvectors.

Here, we give a summary of the Urn models to which the BL model is related as an extension. The Ehrenfest model and the Polya Urn models are two of the early solvable models in the literature [9]. They appear as two of the exactly solved cases in Friedman’s formulation of Urn models where precisely one urn and balls of two colors are drawn and replaced with additions [9]. A dual formulation of Friedman’s Urn models was introduced in a series of recent papers [16, 14]: instead of balls of two colors and one urn, the dual formulation uses two urns and one-colored balls. The latter is more convenient for modelling of certain network science models [3], such as the Voter model where two balls are drawn and returned to the two urns with prescribed probabilities that depend on the order in which they are drawn. This is because many of the network science models are irreversible Markov chain models [4] which have absorbing states. Their transition matrices $T$, unlike $T_{BL}$ for the BL model, are not symmetric, in an essential sense, that is, there are no non-standard inner products for $R^{N}$ in which these matrices $T$ have a symmetric form.

Using a new method based on diagonalization of transition matrices [15], we solved exactly the eigenvectors of several well-known models, including the Ehrenfest model, the Voter model [18, 12, 11, 5, 17, 2], the Moran model for genetic drift, and the Naming game models [19], [20]. Most of these models are irreversible Markov models with absorbing states, and have essentially non-symmetric transition matrices in the sense just mentioned. The BL model however, is based on two urns and balls of two colors. Thus, it is not strictly in the class of Two Urns models to which we recently applied our method. In modifying this method so that it applies to the BL model, we will have shown that the new method is not only easy to use but also flexible in extension to new problems.

One of the main points here is the technical simplicity of uncovering the symmetries of the above LPDOs within our method, through the explicit appearance of the expressions $u=f(x,y)$ in the coefficients of $L.$ We give here the LPDO $L_{V}$ for the Voter model [16, 14] to indicate what we mean: first the propagation equation for the transition matrix $T_{V}$ is given by

[TABLE]

the eigen-problem for the associated LPDO $L_{V}$ $=(x-y)^{2}G_{xy}$ acting on the homogeneous polynomial $G(x,y)=\sum_{j}c_{j}x^{j}y^{N-j}$ (which encodes the components $c_{j}$ of the right-eigenvector of $T_{V})$ is given

[TABLE]

which clearly suggests the transformation $u=x-y,$ $v=y.$ Indeed this triangularized and diagonalized the Voter model and led to its complete solution.

Contrast this ease of use with the fact that triangularization and diagonalization of a given transition matrix of size $N$ has computational complexity $O(N^{3}).$ In other words, exact integration of the Two Urns models via diagonalization of transition matrices are nontrivial problems, that are difficult to solve but once known, the solutions are easy to verify. [16, 14] provides a simple method to find such explicit diagonalization and hence all eigenvectors for a class of transition matrices from the Two Urns models, even when their transition matrices are essentially non-symmetric. Note that the eigenproblem and diagonalization of symmetric matrices have a lower computational complexity.

We aim here to highlight this method’s ease of use, relative to the group representation method and the method of orthogonal polynomials. Moreover, the BL model differs significantly from the original Two Urns subclass of models for which our method was initially formulated. Thus, we also aim to show that, with the specific introduction of a nonlinear change of independent variables, this method can be applied to more complex models than the original class of models. Since the BL transition matrix $T_{BL}$ is from a reversible Markov chain with a stationary distribution [4], it is non-symmetric singly-stochastic only in a trivial sense. In other words, there exists (a difficult to find) non-standard inner product for $R^{N},$ in which $T_{BL}$ is symmetric and hence doubly-stochastic. Karlin and McGregor [10] using their powerful integral representation method of finding an explicit way to symmetrize $T$ by introducing a non-standard inner product (or orthogonal measure) into the problem, have related the right eigenvectors of the transition matrix $T_{BL}$ of the $BL$ model to the orthogonal polynomials called the dual Hahn polynomials. A third aim of this paper is therefore to re-derive from the diagonalization of $T_{BL},$ this non-standard inner product in which the dual Hahn polynomials are an orthogonal polynomial system. Note that this non-standard inner product, once found, yields a symmetric version of $T_{BL}$which is an example of a Jacobi operator that arise in the the classical moments problem [1], and is related to orthogonal polynomials via the Riemann-Hilbert method [6].

The paper is organized as follows: section 3 concerns the calculation of the right eigenvectors and eigenvalues of $T$ in closed form by method introduced in [16, 15, 14]; in view of the fact that these right eigenvectors are not orthogonal in the usual Euclidean inner-product, section 4 concerns the transformations needed to calculate the orthonormal system of left eigenvectors of $T_{BL},$and also the derivation of the non-standard inner product in which the right eigenvectors are now orthogonal; section 5 concerns the elementary proofs, based only on the eigenvectors and eigenvalues of $T,$ for the tight upper and lower bounds for times to stationarity in the BL model [7] ; section 6 concerns a numerically exact evaluation of the expression for the TV norm in these bounds using the eigenvectors and eigenvalues of $T_{BL}$ directly, hence slightly sharper estimates for the mixing times of the BL model.

Beyond the balanced special case of the BL problem treated in detail in this paper, the same generating function method can be used to prove similar tight bounds for mixing times in the other cases.

2 Transition matrix of the BL model

Let the transition matrix $T_{BL}$ be defined so that $(T_{BL\text{ }})_{ij}=Pr\{n_{w}(m+1)=i\;|\;n_{w}(m)=j\}$, so that the sum of each column is $1$. In the general BL model for balls of two colors and two urns, $N_{1},$ $N_{2},$ $N_{w},$ $N_{b}$ are fixed parameters satisfying the constraints

[TABLE]

where $N$ equals total number of balls in the model. For $i=0,...,N_{w}\leq N_{1},$ (where by abuse of notation $i$ stands for both the row label of transposed matrix $T_{BL}^{t}$ and the number of white balls in urn 1, $n_{w}),$ the transition probabilities are explicitly given by

[TABLE]

.

3 Diagonalization - Right eigenvectors of the

general BL model

In [16, 14], we developed an explicit method for exactly integrating or solving a 5-parameters subclass of a class of Two Urns models which is parametrized by six real parameters. Our method is based on a relationship between certain banded stochastic matrices $T$ (such as tridiagonal and pentadiagonal non-symmetric transition matrices of markov chain models) and the LPDOs acting on the vector space of homogeneous polynomials, $G(x,y)$ of finite order in two indeterminates. The symmetries of the LPDO, $L,$ associated with a given non-symmetric singly stochastic matrix from this solvable subclass of the Two Urns models, are identified and used explicitly to transform from the original indeterminates ( independent variables $x,y$ say) to suitable new variables (such as $u=f(x,y),v=g(x,y)$). In the new variables $u,v,$ the transformed LPDO, $L^{\prime},$ acts on the (again homogeneous of same order as $G(x,y))$ polynomial $H(u,v).$ We have shown in [16, 14] that the transformed eigen-problem

[TABLE]

for a well-defined subclass of such Two Urns problems is equivalent (via the inverse of the original relationship between banded matrix and LPDO) to the eigen-problem for a triangular matrix, which can then be solved explicitly for both right and left eigenvectors. In other words, at the end of this brief summary, the symmetries of LPDO $L$ inherited from the original banded stochastic matrix $T,$generate an explicit similarity transformation, $P,$ such that

[TABLE]

where $D$ is diagonal, and $S$ contains the eigenvectors of $T.$

This method can be formalized as an Algorithm as follows: Given the input of a singly stochastic transition matrix $T$ of size $N+1$,

(I) Choose a suitable homogeneous polynomial of finite degree $N$, $G$ that has the components $c_{i}$ of a right eigenvector of $T$ as coefficients of the monomials $x^{i}y^{j}z^{k}$; part of this choice is the number of indeterminates in $G$. For example, the Voter model of size N (number of balls) with a transition matrix $T_{V}(N)$ which is a $N+1$ by $N+1$ real matrix, requires a homogeneous polynomial $G_{V}$ of degree $N$ in the indeterminates $x,y$ because there are two urns.

(II) Associate the recursion inplicit in given Markov matrix $T$ to a LPDO, $L$ which acts on the homogeneous polynomial $G$; the basic elements of this association scheme are the standard linear differential operators for increasing, decreasing and not changing the numbers of balls in each urn (which correspond in the example below to the probabilities $p,q,r$ prescibed by the transition matrix), and a set of multiplication type linear operators that correspond to shifts.

(III) A transformation to new independent variables,(for instance, $u=f(x,y,z),$ $v=g(x,y,z),$ $w=h(x,y,z)$) is chosen to satisfy two conditions:

$(A)$ the transformed polynomial

[TABLE]

is a homogeneous polynomial of the same finite degree as $G;$

$(B)$ $u=f(x,y,z),$ $v=g(x,y,z),$ $w=h(x,y,z)$ is a transformation based on the symmetries of $L$ (cf. [13]), that is, the combinations $f(x,y,z),$ $g(x,y,z),$ $h(x,y,z)$ of the original variables $x,y,z$ appear naturally in the coefficients of the LPDO, $L.$

These conditions $(A)$ and $(B)$ are clearly not sufficient to ensure the transformed LPDO eigenproblem

[TABLE]

is associated with a similar triangular matrix $T^{\prime}$ which explicitly yields all its eigenvectors $b_{i}$. That they are sufficient has to be proved either in each problem to which we apply the Algorithm, or for a class of models as in the case of the Two Urns models.

(IV) Using the transformation in step (III), derive the corresponding transformed LPDO, $L^{\prime}$ that acts on the transformed polynomial $H(u,v,w)$.

(V) Without explicitly calculating the transformed matrix $T^{\prime}$ which is associated with the transformed LPDO, $L^{\prime}$ in step (IV), check that the transformed eigen-problem for $L^{\prime}$ is indeed a recursion system for the transformed eigenvectors $b_{i}$ that can be solved explicitly, i.e., it is equivalent to a triangular linear system of equations. Solve for the eigenvalues and then the eigenvector components $b_{i}$, and if required transform back to the original components $c_{i}$. These are the main outputs of the Algorithm.

(VI) Use the eigenvectors in step (V) to diagonalize the original matrix if necessary. This is the end of the ALgorithm.

Now we apply the Algorithm to the BL model. Given the transition matrix of the BL model (cf. section 2), it will be obvious that three independent variables (instead of the two before) should be used to formulate the BL problem. In step (I) of the Algorithm, we adopt the anzatz that the $LPDO,$ $L_{BL},$ associated with the above $N$ by $N$ matrix $T_{BL}$, now acts on a homogeneous polynomial $G(x,y,z)$ in three indeterminates, $x,y,z.$ We encode the entries $c_{k}(i)$, $i=0,...,N_{w}$ of the $k-th$ eigenvector of the transition matrix for the BL model as follows:

[TABLE]

where $i=$ number of white balls in urn 1 (also denoted $n_{w})$. The choice of three independent variables to encode the components of an eigenvector of $T_{BL}$ in the homogeneous polynomial $G(x,y,x)$ is now made obvious by this explicit expression for $G.$

In step (II) of the Algorithm, we derive from the original eigen-problem for transition matrix $T_{BL},$ an LPDO, $L_{BL}$, that acts on $G^{(k)}$. Towards that aim, we note, in particular, the entries for $p_{i}$ and $q_{i}$ in $T_{BL}$ correspond respectively to the following linear differential operators with coefficients that are monomials in $x,y,$and $z,$

[TABLE]

where $G_{yz}^{(k)}=\frac{\partial^{2}}{\partial y\partial z}G^{(k)}$ for example. In addition, it is part of the association scheme that multiplication in the LPDO (cf. [13]) by the coefficient $\frac{x}{yz}$ (resp. $\frac{yz}{x})$ represents down (resp. up) shifts in the index $i$ within the discrete recursion equations of the original eigen-problem for matrix $T_{BL.}$ The $L_{BL}$ associated with the eigen-problem of the tridiagonal Markov matrix $T_{BL}$is given by:

[TABLE]

In step(III) of the Algorithm, we note that the symmetries of $L_{BL}$ with respect to transformations of its independent variables, is expressed in the factor $\left(x-yz\right)$ in its coefficients. This suggests the transformation to the new independent variables

[TABLE]

Since the transformed homogeneous polynomial is now given by

[TABLE]

in terms of the (new) components $b_{i}^{k}$ of the $k-th$ right eigenvector, this transformation clearly satisfies both necessary conditions (A) and (B) in step (III) of the Algorithm.

To prove that it is sufficient for our purpose of obtaining the eigenvalues and eigenvectors exactly and for diagonalization, we proceed by direct calculations.

In step (IV) using the following obvious identities for the transformation of partial derivatives

[TABLE]

the transformed LPDO, $L_{BL}^{\prime}$ in $H^{(k)}$, $k=0,...,$ $N_{w}$ is

[TABLE]

In step (V), by reversing the derivation of the original $L_{BL}$ through the association scheme [13], this $L_{BL}^{\prime}$ in $H$ is shown to be equivalent to the following triangular system for the (right) eigen-problem of the transformed matrix $T_{BL}^{\prime}$:

[TABLE]

We have therefore verified the sufficiency of the transformation where $u=x-yz$ for triangularizing (and later diagonalizing) $T_{BL}.$ This triangular system implies the recursion

[TABLE]

which can be solved directly.

For nontrivial eigensolutions for $k=0,...,N_{w},$ the denominator in $b_{i}^{k}$ must vanish, yielding the following exact expressions for the eigenvalues,

[TABLE]

In the case $N_{1}=N_{w},$

[TABLE]

The eigenvectors (in the transformed variables of $H)$ are given explicitly by:

[TABLE]

Using these coefficients in the definition for $H$ gives

[TABLE]

We summarize the consequences of the above steps of the Algorithm on the BL model in the following theorem:

Theorem 3.1**.**

*In the above Algorithm for the BL model, for any size $N$ of the model, the LPDO, $L_{BL}^{\prime}$, after the transformation (3.9) on the independent variables, is equivalent to a triangular linear system (3.21) which has (right) eigenvectors given by (3.31) and eigenvalues (3.25). The (right) eigenvectors of the original BL matrix $T_{BL}$ are in turn given by (3.38). *

3.1 Hypergeometric functions and dual Hahn polynomials

The next to final step left in this part of the paper is step (VI) in the Algorithm, to invert the above similarity transformation to obtain explicitly the closed-form expressions for the original components of the right-eigenvectors $c_{k}(i)$ of $T_{BL.}$ For this purpose, let $h^{(k)}(u)=H^{(k)}(u,1,1)$. Then,

[TABLE]

Using the hypergeometric identity,

[TABLE]

and the fact that any multiple of an eigenvector remains an eigenvector, we take the polynomial for the right eigenvector components to be

[TABLE]

whose coefficients are the original components $c_{k}(i)$ of the $k-th$ right eigenvector corresponding to $\lambda_{k}$ prior to the transformation above. These expressions are equivalent to the dual Hahn polynomials [10].

For the hypergeometric representation of the eigenvectors to be well defined, we require $N_{2}\geq N_{w}$. There is no loss in generality with this assumption, because we can relabel $N_{1}\leftrightarrow N_{2}$ and $N_{w}\leftrightarrow N_{b}$ so that the assumption holds. From the solution for $g^{(k)}$, we expand in $x^{l},$

[TABLE]

to find the explicit form for the components of the $k-th$ right eigenvectors:.

[TABLE]

Notice that the solution is the $kth$ order backwards difference of the components of the hypergeometric coefficients of $g$.

The above treatment of the eigen-problem by transforming via symmetries, the independent variables of the associated LPDO, $L_{BL},$ is equivalent to a similarity transformation of the transition matrix $T_{BL}$. Let $\mathbf{w}=\mathbf{Pv}$ for some transformation matrix $\mathbf{P}$. Then, the eigen-problem for $\mathbf{w}$ is given by $\mathbf{P}\mathbf{T}_{BL}\mathbf{P}^{-1}\mathbf{w}=\lambda\mathbf{w}$. The above calculations is equivalent to the matrix $\mathbf{P}$ such that the new matrix $\mathbf{T}_{BL}^{\prime}=\mathbf{P}\mathbf{TP}^{-1}$ is lower triangular. The last step in this section is to diagonalize $T_{BL.}.$ We do this by diagonalizing the matrix triangular matrix $\mathbf{T}_{BL}^{\prime}=\mathbf{W\Lambda W}^{-1}$. Here, $\mathbf{\Lambda}=diag(\lambda_{0},\ldots,\lambda_{N})$ and $\mathbf{W}$ are the eigenvectors of $\mathbf{T}_{BL}^{\prime}$. The components of these eigenvectors are $b_{i}$ corresponding to eigenvalue $\lambda_{k}$. Diagonalization of $\mathbf{T}_{BL}^{\prime}$ allows us to explicitly diagonalize the original transition matrix as

[TABLE]

Note that the matrix of eigenvectors is given by $\mathbf{P}^{-1}\mathbf{W}$.

4 Symmetrizing transform, orthogonal measure and dual Hahn polynomials

For transition matrix $T_{BL},$ let $Z$ be given by (where we drop the subscript BL herein, i.e., $T=T_{BL}$)

[TABLE]

Recall the detailed balance of $T$ and its stationary distribution, given by $T_{ij}\pi_{j}=T_{ji}\pi_{i}$ [4] follows from the reversibility and ergodicity of the BL model. Note that $Z$ is the symmetric version of the transition matrix:

[TABLE]

Therefore, $Z$ has an orthonormal set of left eigenvectors, $w_{k}^{T}$. Let $W$ be a matrix whose columns are $w_{k}$. The spectral decomposition of $Z$ by left eigenvectors is given by

[TABLE]

By the definition of $Z_{ij}$, the transformation from $T$ to $Z$ can be expressed as

[TABLE]

where $D$ is a diagonal matrix whose diagonal entries are $\sqrt{\pi_{i}}$. So, arbitrary powers of $T$ is given by

[TABLE]

Defining a new transformation by $S,$

[TABLE]

Since $W$ has a specific normalization, we can equate $S$ with $DW$ after applying the appropriate normalization for $S$. That is, for diagonal matrix $\Delta$, we take

[TABLE]

We can choose any normalization for the right eigenvectors given in $S$, and $\Delta$ will properly renormalize them. Here, we solve for $\Delta$ and $W$ by appealing to the orthogonality of $W$:

[TABLE]

Therefore

[TABLE]

Computing the matrix multiplication on the left side yields the diagonal entries of $\Delta$ denoted by $\Delta_{k}$ given by

[TABLE]

in terms of the right eigenvectors of $T_{BL}.$ Now that we have $\Delta$, we have the representations for both the left-eigenvectors $w_{k}(i)$ and right-eigenvectors $v_{k}(i)$ of $Z$ given by

[TABLE]

in terms of the right-eigenvectors $c_{k}(i)$ of the original $T_{BL}$ that was obtained by our method in section 2.

From Eq. (57), we also have an explicit formula for $S^{-1}$ given by

[TABLE]

So, by Eq. (54), we have

[TABLE]

Computing the matrix multiplication gives the following spectral decomposition

[TABLE]

as the explicit representation of $Pr\{n(m)=i\;|\;n(0)=j\}$ in the BL model.

Since $T^{0}=I$, take $m=0$ in Eq. (63) to find the stationary distribution of the BL model

[TABLE]

Note that this is the orthogonality relation for the right-eigenvectors of $T_{BL}$ with orthogonal measure $\Delta_{k}^{2}$ given in Eq. (63). We have derived the orthogonal measure $\Delta_{k}^{2}$ in which the dual Hahn are an orthogonal polynomial system [10].

We summarize these results on the derivation of a non-standard inner product or orthogonal measure in which the original transition matrix of the BL model becomes a symmetric real matrix and the (right) eigenvectors are the system of orthogonal dual Hahn polynomials:

Theorem 4.1**.**

*The orthogonal measure in (4.12) symmetrizes $T_{BL}$, is related to the (left) and (right) eigenvectors of $T_{BL}$ by (4.13) and (4.14), and yields the spectral decomposition (4.17). *

5 Bounds of mixing times - elementary proofs

We will discuss first the case $N_{1}=N_{2}=N/2,$ for our method gave the eigenvalues of the BL model to be $\lambda_{k}=1-\frac{4k\left(N-k+1\right)}{N^{2}},$ $\lambda_{1}=1-\frac{4}{N}.$ A heuristic estimate of the number of switches $q$ needed to mix the colors is thus,

[TABLE]

and therefore $q=\frac{1}{4}N\log N$ gives an idea of how many switches or time steps are needed until the variation distance between $\rho_{q}$ and $\pi$ is of order $O(1/N).$ The lower bound is obtained along the lines of Diaconis et al, that is, by an application of the Chebyshev’s inequality. However, all estimates of the relevant mean values and variance needed to apply the Chebyshev’s inequality are constructed explicitly from the properties of the eigenvalues and eigenvectors of the BL model. We will state the theorems below but since the proofs are similar to those in [7], we have provided the details in an appendix.

Theorem 5.1**.**

*For $m=\frac{1}{4}N\log N+(\frac{c}{2}-\frac{\log 2}{8})N,$ for $c>0,$ there is a universal constant $A>0$ such that $E_{\pi}\left[\left\|\rho_{m}(j;.)-\pi\right\|_{V}\right]\leq Ae^{-2c}.$ *

Theorem 5.2**.**

*If $m=\frac{1}{8}N\ln N-\frac{cN}{2},$ then $2||\rho_{m}-\pi||_{V}\geq 1-e^{4c}$ *

6 Exact calculations of mixing times

Given that we can calculate $\rho_{m}(i,j)$ exactly by Eq. (63), and the stationary distribution is given by

[TABLE]

the total variational distance can be exactly computed for all time steps by

[TABLE]

This solution is shown in Figure 1, with the upper bound given in Figure 2.

Appendix A Proofs of Theorems 5.1 and 5.2

The eigenvalues and eigenvectors of the BL problem can be used to construct an upper bound on the variation distance between $the$ probability distribution after $m$. As we proved in section 3, the right-eigenvectors of the original matrix $T_{BL}$ are not orthogonal in the standard inner product of $R^{N}$ but a different inner product weighted by $\Delta_{k}^{2}$ can be used to derive orthogonality of a related system of right-eigenvectors $v_{j}$ of $T=T_{BL}.$

A.1 Upper bound

Let $j=0$ to define

[TABLE]

where

[TABLE]

Then, $\rho_{0}(i;j)=\delta_{ij},$ and hence

[TABLE]

implies that

[TABLE]

Using

[TABLE]

and averaging over initial data $n_{w}(0)=j,$ $j=0,...,N/2,$we have

[TABLE]

for some $A$ independent of $N.$ This proves theorem 5.1.

A.2 Lower bound

In terms of the right-eigenvectors $v_{k},$ $k=0,...,N/2,$ (with col sum =1 )

[TABLE]

and

[TABLE]

Next for $m=\frac{1}{8}N\log N-c\frac{N}{2}$, we get $E[v_{1}]=\frac{v_{1}(0)}{\sqrt{N}}e^{2c}$. A similar calculation gives

[TABLE]

Next we deduce $v_{1}(0)$ and $v_{2}(0)$ from $v_{1}^{2}=Av_{2}+B$, $v_{1}(i)=C(N/4-i)$, and the orthogonality of $v_{i}$:

[TABLE]

Therefore, taking $C=\frac{4}{\sqrt{N}}$, we fine $v_{1}(0)\sim\sqrt{N}$. Furthermore, $Av_{2}(0)=v_{1}(0)-b$. Now, we have

[TABLE]

So with the same normalization as above for $v_{1},$ we deduce $Var\{v_{1}^{\prime}\}$ is uniformly bounded by constant $2b,$ since $B=b+O(\log N/N).$ Now, by Chebyschev’s inequality,

[TABLE]

and

[TABLE]

Thus, if $K\subset\{0,....,N/2\}$ such that $|v_{1}|\leq k$ for $k\in K,$ we deduce

[TABLE]

Choose $k=d\sqrt{N}\lambda_{1}^{m}$ to obtain

[TABLE]

which proves theorem 5.2.