$\gdef \vec#1{\boldsymbol{#1}} \\ \gdef \rank {\mathrm{rank}} \\ \gdef \det {\mathrm{det}} \\ \gdef \Bern {\mathrm{Bern}} \\ \gdef \Bin {\mathrm{Bin}} \\ \gdef \Mn {\mathrm{Mn}} \\ \gdef \Cov {\mathrm{Cov}} \\ \gdef \Po {\mathrm{Po}} \\ \gdef \HG {\mathrm{HG}} \\ \gdef \Geo {\mathrm{Geo}}\\ \gdef \N {\mathrm{N}} \\ \gdef \LN {\mathrm{LN}} \\ \gdef \U {\mathrm{U}} \\ \gdef \t {\mathrm{t}} \\ \gdef \F {\mathrm{F}} \\ \gdef \Exp {\mathrm{Exp}} \\ \gdef \Ga {\mathrm{Ga}} \\ \gdef \Be {\mathrm{Be}} \\ \gdef \NB {\mathrm{NB}} \\ \gdef \indep {\mathop{\perp\!\!\!\!\perp}} \\ \gdef \tr {\mathrm{tr}}$
0. 導入
＜分割表の統計モデル＞で扱った3パターンとして、
①どこもfixされない場合
②総和がfixされる場合
③列和または行和がfixされる場合
を考えていきます。
具体的には、各パターンについて、
(ⅰ)フルモデル
(ⅱ)条件$H_0$が追加で課されたモデル
を考えていきます。
以下は＜分割表の統計モデル＞で扱った内容の抜粋です。
0-1. どこもfixされない場合
各セルの度数はポアソン分布に従うというモデルを想定することができ、同時確率関数$p(x)$は、
$$\begin{aligned} p(x) &= \prod_{i,j} e^{-\theta_{ij}} \frac{\theta_{ij}^{x_{ij}}}{x_{ij}!} ~~~~~ \mathrm{(A)} \\ &{\scriptsize(X_ij \sim \Po(\theta_{ij})とした)} \end{aligned}$$となり、パラメータの自由度は$IJ$となります。
$\mathrm{(A)}$のみ条件が課されたモデルがフルモデルです。
パラメータ$\theta_{ij}$について、
$$\begin{aligned} \log \theta_{ij} &= \mu + \alpha_i + \beta_j + \gamma_{ij} ~~~~~ \mathrm{(B)} \\ &{\scriptsize(ただし、\textcolor{red}{\sum_i \alpha_i = \sum_j \beta_j = \sum_i \gamma_{ij} = \sum_j \gamma_{ij} = 0})} \end{aligned}$$というパラメータ変換を行い、 帰無仮説$H_0$を、$H_0: \gamma_{ij} = 0 ~~ {\small (\forall i, j)} ~~~~~ \mathrm{(C)}$、とすると、$\mathrm{(A),(B),(C)}$が課されたモデルが『条件$H_0$が追加で課されたモデル』です。
0-2. 総和がfixされる場合
各セルの度数は$IJ$項分布に従うというモデルを想定することができ、同時確率関数$p(x)$は、
$$\begin{aligned} p(x) &= \frac{n!}{\prod_{i,j} x_{ij}!} \prod_{i,j} p_{ij}^{x_{ij}} ~~~~~ \mathrm{(D)} \\ &{\scriptsize(X_ij \sim \Bin(n,p_{ij})とした)} \end{aligned}$$となり、パラメータの自由度は$(IJ-1)$となります。
$\mathrm{(D)}$のみ条件が課されたモデルがフルモデルです。
（ただし、$\sum_{i,j} p_{ij} = 1, \sum_{i,j} x_{ij} = n$、は前提）
帰無仮説$H_0$を、$H_0: p_{ij} = p_{i \bullet} p_{\bullet j} ~~ {\small (\forall i, j)} ~~~~~ \mathrm{(E)}$、とすると、$\mathrm{(D),(E)}$が課されたモデルが『条件$H_0$が追加で課されたモデル』です。
0-3. 列和または行和がfixされる場合
（ここでは列和がfixされる場合を扱います）
列和がfixされる場合、各セルの度数は$I$項分布$(\times J)$に従うというモデルを想定することができ、同時確率関数$p(x)$は、
$$\begin{aligned} p(x) &= \prod_{j} \{ \frac{n_j!}{\prod_i x_{ij}!} \prod_i p_{ij}^{x_{ij}} \} ~~~~~ \mathrm{(F)} \end{aligned}$$となり、パラメータの自由度は$(I-1) \times J$となります。
$\mathrm{(F)}$のみ条件が課されたモデルがフルモデルです。
（ただし、$\sum_{i} p_{ij} = 1 ~~ {\small (\forall j)}, \sum_{i} x_{ij} = n_j ~~ {\small (\forall j)}$、は前提）
帰無仮説$H_0$を、$H_0: p_{i1} = p_{i2} = \cdots = p_{iJ} ~~ (\forall i) ~~~~~ \mathrm{(G)}$、とすると、$\mathrm{(F),(G)}$が課されたモデルが『条件$H_0$が追加で課されたモデル』です。
それでは分割表の推定を各パターンについて順にみていきましょう。
1. どこもfixされない場合
1-1. フルモデル
尤度を$L(\theta)$とすると、$\mathrm{(A)}$より、
$$\begin{aligned} \log L(\theta) &= \log p(x) \\[10px] &= \sum_{i,j} (-\theta_{ij} + x_{ij} \log \theta_{ij}) + {\small (\thetaに無関係の項)} \\[20px] \Rightarrow \frac{\partial}{\partial \theta_{ij}} \log L(\theta) &= -1 + \frac{x_{ij}}{\theta_{ij}} = 0 ~~ {\small (と置く)} \\[10px] \Rightarrow \hat{\theta}_{ij} &= x_{ij} \end{aligned}$$
1-2. 条件$H_0$が追加で課されたモデル
尤度を$L'(\mu,\vec \alpha, \vec \beta)$とすると、$\mathrm{(A),(B),(C)}$より、
$$\begin{aligned} \log L'(\mu,\vec \alpha, \vec \beta) &= \log p(x) \\[10px] &= \sum_{i,j} \{ -\exp(\mu + \alpha_i + \beta_j) + x_{ij} (\mu + \alpha_i + \beta_j) \} + {\small (\mu,\vec \alpha, \vec \betaに無関係の項)} \end{aligned}$$
これを$\mu, \alpha_i, \beta_j$でそれぞれ偏微分して$=0$と置くと以下の様になる。
$$\begin{aligned} \frac{\partial}{\partial \mu} \log L'(\mu,\vec \alpha, \vec \beta) &= \sum_{i,j} \{ -\exp(\mu + \alpha_i + \beta_j) + x_{ij} \} = 0 ~~ {\small (と置く)} \\[10px] \Rightarrow e^{\mu} \cdot \sum_i e^{\alpha_i} \cdot \sum_j e^{\beta_j} &= x_{\bullet \bullet} ~~~~~ \mathrm{(H)} \\ &{\scriptsize (x_{\bullet \bullet} = \sum_{i,j} x_{ij}と置いた)} \\[20px] \frac{\partial}{\partial \alpha_i} \log L'(\mu,\vec \alpha, \vec \beta) &= \sum_{j} \{ -\exp(\mu + \alpha_i + \beta_j) + x_{ij} \} = 0 ~~ {\small (と置く)} \\[10px] \Rightarrow e^{\mu} \cdot e^{\alpha_i} \cdot \sum_j e^{\beta_j} &= x_{i \bullet} ~~~~~ \mathrm{(I)} \\ &{\scriptsize (x_{i \bullet} = \sum_{j} x_{ij}と置いた)} \\[20px] \frac{\partial}{\partial \beta_j} \log L'(\mu,\vec \alpha, \vec \beta) &= \sum_{i} \{ -\exp(\mu + \alpha_i + \beta_j) + x_{ij} \} = 0 ~~ {\small (と置く)} \\[10px] \Rightarrow e^{\mu} \cdot \sum_i e^{\alpha_i} \cdot e^{\beta_j} &= x_{\bullet j} ~~~~~ \mathrm{(J)} \\ &{\scriptsize (x_{\bullet j} = \sum_{i} x_{ij}と置いた)} \\ \end{aligned}$$
$\mathrm{(H),(I)}$から、
$$\begin{aligned} \frac{e^{\alpha_i}}{\sum_k e^{\alpha_k}} &= \frac{x_{i \bullet}}{x_{\bullet \bullet}} \\[10px] \Rightarrow e^{\sum_i \alpha_i} &= (\sum_k e^{\alpha_k})^I \cdot (\frac{1}{x_{\bullet \bullet}})^I \cdot \prod_i x_{i \bullet} \\ &{\scriptsize(e^{\alpha_i} (i=1,\ldots,I)を掛け合わした)} \\[20px] \Rightarrow \sum_k e^{\alpha_k} &= \frac{x_{\bullet \bullet}}{(\prod_i x_{i \bullet})^{\frac{1}{I}}} \\ &{\scriptsize(e^{\sum_i \alpha_i} = e^0 = 1であった)} \end{aligned}$$
これと$\mathrm{(K)}$より、
$$\begin{aligned} e^{\alpha_i} &= \frac{x_{i \bullet}}{x_{\bullet \bullet}} \cdot \frac{x_{\bullet \bullet}}{(\prod_i x_{i \bullet})^{\frac{1}{I}}} \\[10px] \Rightarrow \hat{\alpha}_i &= \log \frac{x_{i \bullet}}{(\prod_i x_{i \bullet})^{\frac{1}{I}}} \end{aligned}$$
同様に$\mathrm{(H),(J)}$から、
$$\begin{aligned} \hat{\beta}_j = \log \frac{x_{\bullet j}}{(\prod_j x_{\bullet j})^{\frac{1}{J}}} \end{aligned}$$が導出され、これらと$\mathrm{(H)}$から、
$$\begin{aligned} \hat{\mu} = \log \frac{ (\prod_i x_{i \bullet})^{\frac{1}{I}} \cdot (\prod_j x_{\bullet j})^{\frac{1}{J}} }{x_{\bullet \bullet}} \end{aligned}$$となる。
これらの結果を用いると、
$$\begin{aligned} \hat{\theta}_{ij} &= e^{\hat{\mu} + \hat{\alpha}_i + \hat{\beta}_j} \\[10px] &= \frac{x_{i \bullet} x_{\bullet j}}{x_{\bullet \bullet}} \end{aligned}$$となる。
（補足1）
計算を解く方法としてはラグランジュの未定乗法を用いたものも可能です。
具体的には、
$$\begin{aligned} \mathcal{L} := \log p(x)-a\sum_i \alpha_i-b\sum_j \beta_j \end{aligned}$$とおいて、$\mu, \alpha_i, \beta_j$でそれぞれ偏微分して$=0$として解き進める形になります。
今回についてはラグランジュの未定乗数法よりも上記計算の方が早そうだったので、そちらを採用しました。
（補足2）
計算過程も重要ですが、結果が重要です。
$(i,j)$セルの期待度数は$\hat{\theta}_{ij}$となりますが、フルモデルの場合には、$\hat{\theta}_{ij}=x_{ij}$、と$x_{ij}$のみから推定値が構成されます。一方で、条件$H_0$が追加で課されたモデルの場合には、$\hat{\theta}_{ij}=\frac{x_{i \bullet} x_{\bullet j}}{x_{\bullet \bullet}}$、と$x_{\bullet \bullet}, x_{i \bullet}, x_{\bullet j}$から推定値が構成されます。
この結果は直感的にも納得のいくものですね。
2. 総和がfixされる場合
2-1. フルモデル
尤度を$L(\vec p)$とすると、$\mathrm{(D)}$より、
$$\begin{aligned} \log L(\vec p) &= \log p(x) \\[10px] &= \sum_{i,j} x_{ij} \log p_{ij} + {\small (p_{ij}に無関係の項)} \end{aligned}$$となる。
上式の前半部分を制約条件$\sum_{i,j} p_{ij} = 1$の下で最大化するために、以下ラグランジュの未定乗数法を用いる。
$$\begin{aligned} \mathcal{L} := \sum_{i,j} x_{ij} \log p_{ij}-a(\sum_{i,j} p_{ij}-1) \end{aligned}$$とおいて、これを$p_{ij}$で偏微分して$=0$とすると、
$$\begin{aligned} \frac{\partial \mathcal{L}}{\partial p_{ij}} &= \frac{x_{ij}}{p_{ij}}-a = 0 \\[10px] \Rightarrow \hat{p}_{ij} &= \frac{x_{ij}}{a} \end{aligned}$$
これと制約条件$\sum_{i,j} p_{ij} = 1$から、$a = x_{\bullet \bullet}$、となる。
以上より、
$$\begin{aligned} \hat{p}_{ij} = \frac{x_{ij}}{x_{\bullet \bullet}} \end{aligned}$$
2-2. 条件$H_0$が追加で課されたモデル
尤度を$L'(\vec p)$*とすると、$\mathrm{(D),(E)}$より、
$$\begin{aligned} \log L'(\vec p) &= \log p(x) \\[10px] &= \sum_{i,j} x_{ij} \log p_{ij} + {\small (p_{ij}に無関係の項)} \\[10px] &= \sum_{i,j} x_{ij} \log (p_{i \bullet} p_{\bullet j}) + {\small (p_{ij}に無関係の項)} \\[10px] &= \sum_{i,j} x_{ij} (\log p_{i \bullet} + \log p_{\bullet j}) + {\small (p_{ij}に無関係の項)} \end{aligned}$$
(*：2-1では$\vec p$は$p_{ij}$が並んだベクトルを意図しており、2-2では$\vec p$は$p_{i \bullet}, p_{\bullet j}$が並んだベクトルを意図しています。そのため本当であれば2-2では$\vec p$とは別の記号として$\vec p’$などとすべきですがわかりやすさを優先しました。)
この前半部分を制約条件$\sum_{i} p_{i \bullet} = \sum_{j} p_{\bullet j} = 1$の下で最大化するために、以下ラグランジュの未定乗数法を用いる。
$$\begin{aligned} \mathcal{L} := \sum_{i,j} x_{ij} (\log p_{i \bullet} + \log p_{\bullet j})-a(\sum_i p_{i \bullet}-1)-b(\sum_j p_{\bullet j}-1) \end{aligned}$$とおいて、これを$p_{i \bullet}, p_{\bullet j}$で偏微分して$=0$とすると、
$$\begin{aligned} \frac{\partial \mathcal{L}}{\partial p_{i \bullet}} &= \sum_j \frac{x_{ij}}{p_{i \bullet}}-a = 0 \\[10px] \Rightarrow \hat{p}_{i \bullet} &= \frac{x_{i \bullet}}{a} \\[20px] \frac{\partial \mathcal{L}}{\partial p_{\bullet j}} &= \sum_i \frac{x_{ij}}{p_{\bullet j}}-b = 0 \\[10px] \Rightarrow \hat{p}_{\bullet j} &= \frac{x_{\bullet j}}{b} \end{aligned}$$
これと制約条件$\sum_{i} p_{i \bullet} = \sum_{j} p_{\bullet j} = 1$から、$a = x_{\bullet \bullet}, b = x_{\bullet \bullet}$、となる。
以上より、
$$\begin{aligned} \hat{p}_{i \bullet} &= \frac{x_{i \bullet}}{x_{\bullet \bullet}} \\[10px] \hat{p}_{\bullet j} &= \frac{x_{\bullet j}}{x_{\bullet \bullet}} \\[10px] \hat{p}_{ij} &= \hat{p}_{i \bullet} \cdot \hat{p}_{\bullet j} \\[10px] &= \frac{x_{i \bullet} x_{\bullet j}}{x_{\bullet \bullet}^2} \end{aligned}$$
（補足）
計算過程も重要ですが、結果が重要です。$(i,j)$セルの期待度数は$x_{\bullet \bullet} \hat{p}_{ij}$となりますが、フルモデルの場合には、$x_{\bullet \bullet} \hat{p}_{ij} = x_{ij}$、と$x_{ij}$のみから推定値が構成されます。一方で、条件$H_0$が追加で課されたモデルの場合には、$x_{\bullet \bullet} \hat{p}_{ij} = \frac{x_{i \bullet} x_{\bullet j}}{x_{\bullet \bullet}}$、と$x_{\bullet \bullet}, x_{i \bullet}, x_{\bullet j}$から推定値が構成されます。
『2. 総和がfixされる場合』におけるフルモデルの場合、条件$H_0$が追加で課されたモデルの場合、の$(i,j)$セルの期待度数はそれぞれ、『1. どこもfixされない場合』におけるそれらと等しいですね！
3. 列和または行和がfixされる場合
（ここでは行和がfixされる場合を扱います）
3-1. フルモデル
尤度を$L(\vec p)$とすると、$\mathrm{(D)}$より、
$$\begin{aligned} \log L(\vec p) &= \log p(x) \\[10px] &= \sum_{j} \sum_{i} x_{ij} \log p_{ij} + {\small (p_{ij}に無関係の項)} \end{aligned}$$となる。
この前半部分を制約条件$\sum_{i} p_{ij} = p_{\bullet j}$の下で最大化するために、以下ラグランジュの未定乗数法を用いる。
$$\begin{aligned} \mathcal{L} := \sum_{j} \sum_{i} x_{ij} \log p_{ij}-a(\sum_{i} p_{ij}-p_{\bullet j}) \end{aligned}$$とおいて、これを$p_{ij}$で偏微分して$=0$とすると、
$$\begin{aligned} \frac{\partial \mathcal{L}}{\partial p_{ij}} &= \frac{x_{ij}}{p_{ij}}-a = 0 \\[10px] \Rightarrow \hat{p}_{ij} &= \frac{x_{ij}}{a} \end{aligned}$$
これと制約条件$\sum_{i} p_{ij} = p_{\bullet j}$から、$a = x_{\bullet j}$、となる。
以上より、
$$\begin{aligned} \hat{p}_{ij} = \frac{x_{ij}}{x_{\bullet j}} \end{aligned}$$
3-2. 条件$H_0$が追加で課されたモデル
尤度を$L'(\vec p)$*とすると、$\mathrm{(F),(G)}$より、
$$\begin{aligned} \log L'(\vec p) &= \log p(x) \\[10px] &= \sum_{j} \sum_{i} x_{ij} \log p_{ij} + {\small (p_{ij}に無関係の項)} \\[10px] &= \sum_{j} \sum_{i} x_{ij} \log p_{i \bullet} + {\small (p_{ij}に無関係の項)} \\ &{\scriptsize (p_{i \bullet} := p_{i1} (=p_{i2}=\cdots=p_{iJ}))} \end{aligned}$$
(*：3-1では$\vec p$は$p_{ij}$が並んだベクトルを意図しており、3-2では$\vec p$は$p_{i \bullet}$が並んだベクトルを意図しています。そのため本当であれば3-2では$\vec p$とは別の記号として$\vec p’$などとすべきですがわかりやすさを優先しました。)
この前半部分を制約条件$\sum_{i} p_{i \bullet} = 1$の下で最大化するために、以下ラグランジュの未定乗数法を用いる。
$$\begin{aligned} \mathcal{L} := \sum_{j} \sum_{i} x_{ij} \log p_{i \bullet}-a(\sum_{i} p_{i \bullet}-1) \end{aligned}$$とおいて、これを$p_{i \bullet}$で偏微分して$=0$とすると、
$$\begin{aligned} \frac{\partial \mathcal{L}}{\partial p_{i \bullet}} &= \sum_j \frac{x_{ij}}{p_{i \bullet}}-a = 0 \\[10px] \Rightarrow \hat{p}_{i \bullet} &= \frac{x_{i \bullet}}{a} \end{aligned}$$
これと制約条件$\sum_{i} p_{i \bullet} = 1$から、$a = x_{\bullet \bullet}$、となる。
以上より、
$$\begin{aligned} (\hat{p}_{i j} =)\hat{p}_{i \bullet} &= \frac{x_{i \bullet}}{x_{\bullet \bullet}} \end{aligned}$$
（補足）
計算過程も重要ですが、結果が重要です。$(i,j)$セルの期待度数は$x_{\bullet j} \hat{p}_{ij}$となりますが、フルモデルの場合には、$x_{\bullet j} \hat{p}_{ij} = x_{ij}$、と$x_{ij}$のみから推定値が構成されます。一方で、条件$H_0$が追加で課されたモデルの場合には、$x_{\bullet j} \hat{p}_{ij} = \frac{x_{i \bullet} x_{\bullet j}}{x_{\bullet \bullet}}$、と$x_{\bullet \bullet}, x_{i \bullet}, x_{\bullet j}$から推定値が構成されます。
『3. 列和または行和がfixされる場合』におけるフルモデルの場合、条件$H_0$が追加で課されたモデルの場合、の$(i,j)$セルの期待度数はそれぞれ、『1. どこもfixされない場合』『2. 総和がfixされる場合』におけるそれらと等しいですね！
そうなのです。最終的には、3パターンのフルモデル・条件$H_0$が追加で課されたモデルの$(i,j)$セルの期待度数は一致するという結果がえられます（重要）。そしてこの期待度数は直感的に納得いくものであり、暗記しておいてもよいかと思います。
まとめ.
- 分割表における推定ではどこがfixされるかで3パターンにわかれ、フルモデルか条件$H_0$が追加で課されたモデルかで2パターンにわかれるため、合計6パターンが考えられる。
- フルモデルにおける$(i,j)$セルの期待度数は、どの3パターンについても、$x_{ij}$、となる。
- 条件$H_0$が追加で課されたモデルにおける$(i,j)$セルの期待度数は、どの3パターンについても、$\frac{x_{i \bullet} x_{\bullet j}}{x_{\bullet \bullet}}$、となる。
