A markov chain sampler for contingency table exact inference

Summary

In the inference of contingency table, when the cell counts are not large enough for asymptotic approximation, conditioning exact method is used and often computationally impractical for large tables. Instead, various sampling methods can be used. Based on permutation, the Monte Carlo sampling may become again impractical for large tables. For this, existing the Markov chain method is to sample a few elements of the table at each iteration and is inefficient. Here we consider a Markov chain, in which a sub-table of user specified size is updated at each iteration, and it achieves high sampling efficiency. Some theoretical properties of the chain and its applications to some commonly used tables are discussed. As an illustration, this method is applied to the exact test of the Hardy-Weinberg equilibrium in the population genetics context.

Appendix: Proof of the Theorem

• (i) Irreducibility: We need to show that ∀T, T′ ∈ Γ₀, T′ can be obtained from T through finite number of transitions and P(T′|T ) > 0 (Hastings, 1970). Let (i₁, i₂, j₁, j₂) be the basic move which decrease 1 count at the positions (i₁, j₁) and (i₂, j₂) each, and increase 1 count at the positions (i₁, j₂) and (i₂, j₁) each. Apparently, any of the basic moves keeps the boundary condition unchanged. In case of row plus column constraint with sub-table requirement (b), and other cases with sub-table requirement (a), T′ can be obtained by finite number $S=\frac{1}{4}\sum_{ij}|T_{ij}^{'}-T_{ij}|$ of basic moves e_i’s from T, each basic move involves four count changes in four positions in T, which can be covered by a sub-table D_i with dimension no less than two (in the row plus column constraint case, can be covered by at most two sub-tables with dimension no less than three), since P(D_i) > 0 and P(e_i|D_i) > 0, and the transition from T to T′ can possibly be achieved through some other ways, so

  $$ P(T^{'}|T)\geq\prod_{i=1}^SP(D_i)P(e_i|D_i)>0. $$

  Reversibility: For any T and T′ ∈ Г₀, by the proof above, there are finite number of intermediate states T₁,..., T_m−1, and sub-tables D₁,..., D_m, such that T, the T_i’s and T′ are only differ on a sub-table D_i, and T′ is obtained by successive transition from T via T₁, ..., T_m, let T₀ = T and T_m= T′. It is easy to check that $P(T_{i-1})P(T_i|B_{D_i})=P(T_i)P(T_{i-1}|B_{D_i}),$ by the definition of the transition probability, so we have

  $$ P(T)P(T^{'}|T)=P(T)\prod_{i=1}^mP(D_i)P(T_i|T_{i-1})=\prod_{i=1}^mP(D_i)P(T_{i-1})P(T_i|B_{D_i})\\=\prod_{i=1}^mP(D_i)P(T_{i-1})P(T_i|B_{D_i})=P(T^{'})\prod_{i=1}^mP(D_i)P(T_i|T_{i-1}|B_{D_i})=P(T^{'})P(T|T^{'}). $$

  The reversibility ensures P(·) is the invariant distribution of the chain, the irreducibility ensures P(·) is the unique equilibrium distribution of the chain, and the state space of the chain is finite, thus (6) and (7) are direct results in Chung (1960, p.99).

• (ii) Liu, et al (1995) proved the results for Gibbs sampler, the same method applies here, and the required conditions (a), (b) and (c) in Liu, et al (1995) are satisfied by this sampler since the state space is finite and the chain is irreducible. Then define the forward operator and follow steps similarly there.

• (iii) since P(·) is an (initial) stationary distribution of the chain, E_P|g(T)l < ∞ and irreducibility of the chain is the same as Markov ergodic, Theorem 3.6.7 in Stout (1974) gives the result.