2.3 統計モデル構築のいくつかの原理

2.3.1探索的データ解析

基本的で重要なこととして、以下のことを考える。
1. 測定データは、量的 or 質的。質的ならば名義 or 順序。
2. ヒストグラムや散布図から分布の形を知る。
3. 他の変数との関係を散布図や箱ひげ図などから知る。
#### 2.3.2 モデルの定式化 １つの応答変数Yといくつかの説明変数を含むモデルでは、Yの確率分布とYと説明変数の関係式という２つの成分をもつ。
\[Y 〜 N (μ, σ^2)\] 一般化線形モデルにおいて、確率分布はすべて指数型分布族（正規分布や二項分布、ポアソン分布 etc.）に属する。
\[E(Y)=α+β_1x_1+β_2x_2+\cdots+β_mx_m\] このときの右辺を線形成分（linear component）と呼ぶ。詳しくは2.4で。

2.3.4 残差とモデルのチェック

残差はモデルの仮定をチェックするための手段である。 それらは互いに独立で、近似的に平均0、分散一定の正規分布に従うべきである。また、説明変数とは無関係である必要がある。

<正規性の確認>

そこで、標準化残差（残差を標準偏差で割ったもの）を正規分布と比較し、度数分布が正規分布に適合しているか、一定の範囲に内側になっているかチェックする。5％なら±1.96。

まず、正規分布を仮定したモデルを考える。応答変数Yiが \[E(Y_i)=μ_i；Y_i〜N（μ_i, σ^2）\] に従うとき、近似的な標準化残差は \[r_i=（y_i-\hat{μ}）/\hat{σ}\] と定義される。
ポアソン分布を仮定したモデルの場合、応答変数Yiは \[E(Y_i)=θ_i；Y_i〜Poisson(θ_i)\] に従う。このときの標準化残差は \[r_i=\frac{y_i-\hat{θ_i}}{\sqrt{\hat{θ}}}\] と表される。
正規性を評価するより感度のよい方法は正規確率プロットを用いることである。
これは各観測値（今回であれば残差）を大きい順に並べ、正規仮定のもとで得られた期待値に対してプロットしたもの。

set.seed(123)
N<-100
d<-rnorm(N, mean=2.4, sd=5)
qqnorm(d)

対角線上あるいはその近く散らばる点は正規分布に従うと考えられ、系統的な偏りや外れ値は正規分布からの外れを表す。

＜説明変数との関係の確認＞

標準化残差の各説明変数やモデルに含まれない潜在変数に対してもプロットする。
変数の効果を適切に表していれば標準化残差は説明変数の値とは無関係に一様分布になるし、曲線的傾向や他の系統パターンがみられれば、他の説明変数に変えたり、追加したりしたほうがよい。

＜分散が一定であることの確認＞

標準化残差を当てはめ値Yiに対してプロットすることで等分散性を確かめる。

＜残差の独立性の確認＞

Yiを測定した順番（時間とか空間）にプロットする。 残差が独立であれば、ランダムにプロットされるはずである。ここでもし何らかの関係性がある場合には特別なモデル化の方法（時間的・空間的自己相関を考慮したモデル）が必要となる。詳しくは11章にて。

2.3.5 推測と解釈

・データのばらつきほとんどすべて説明しきるような複雑なモデルより、データを適度に記述する簡単なモデルが好ましい。
・モデルは自然な解釈がしやすい変数にするとよい。

パラメータの推定値や標準誤差・信頼区間といった推定値の信頼性は、有意水準やp値よりも情報量が多く、推定の精度や“意味のある”といえるほど大きな効果かを考える上で非常に重要である。

2.4 説明変数についての表記法とコーディング

応答変数と各説明変数からなるモデルは \[g[E(Y)]=β_0+β_1x_1+\cdots+β_mx_m・・・・・（2.1）\] と表され、行列表記にすると \[g[E(y)]=Xβ\] ここで、 \[y=\begin{pmatrix}Y_1\\・\\・\\Y_N\end{pmatrix}\] yは反応のベクトルであり \[β=\begin{pmatrix}β_1\\・\\・\\β_p\end{pmatrix}\] βはパラメータのベクトルである。
Xは行列であり、その要素はカテゴリカルな説明変数に対応した定数か量的な説明変数の測定値である。

2.4.1 例：2群における平均

例として、2.2.1の慢性病状のモデル式は \[E(Y_{jk})=θ_j；Y_{jk}〜Poisson(θ_j)，j=1,2\] 行列表記にすると \[y=\begin{pmatrix}Y_{1,1}\\Y_{1,2}\\・\\・\\Y_{1,26}\\Y_{2,1}\\・\\・\\Y_{2, 23}\end{pmatrix}，β=\begin{pmatrix}θ_1\\θ_2\end{pmatrix}，X=\begin{pmatrix}1&0\\1&0\\・&・\\・&・\\1&0\\0&1\\・&・\\・&・\\0&1\end{pmatrix}\] Xの上段でE(Y1, k)に対応するθ1を拾い出し、下段でE(Y2,k)に対応するθ2を拾い出す。このモデルから群平均θ1とθ2を推定したり、比較したりできる。

2.4.2 例：2群における単回帰

2.2.2の出生時体重と妊娠期間のモデル式は \[E(Y_{jk})=μ_{jk}=α_j+β_jx_{jk}；Y_{jk}〜N(μ_{jk}, σ^2)\] 行列表記にすると \[y=\begin{pmatrix}Y_{1,1}\\Y_{1,2}\\・\\・\\Y_{1,K}\\Y_{2,1}\\・\\・\\Y_{2, K}\end{pmatrix}，β=\begin{pmatrix}α_1\\α_2\\β_1\\β_2\end{pmatrix}，X=\begin{pmatrix}1&0&x_{11}&0\\1&0&x_{12}&0\\・&・&・&・\\・&・&・&・\\1&0&x_{1K}&0\\0&1&0&x_{21}\\・&・&・&・\\・&・&・&・\\0&1&0&x_{2K}\end{pmatrix}\] Xの上段で男児、下段で女児、3,4列目で妊娠期間の値が入る。

2.4.3 例：2群の平均の比較に対するいくつかの定式化

\[(a)\quad E(Y_{1k})=β_1とE(Y_{2k})=β_2\] この場合 \[β=\begin{pmatrix}β_1\\β_2\end{pmatrix}, \quad 群1：\begin{pmatrix}1\\0\end{pmatrix}, \quad 群2：\begin{pmatrix}0\\1\end{pmatrix}\] .

\[(b)\quad E(Y_{1k})=μ+α_1とE(Y_{2k})=μ+α_2\] このとき、μが総平均を表し、α1とα2がμからの各群の差を表す。 \[β=\begin{pmatrix}μ\\α_1\\α_2\end{pmatrix}, \quad 群1：\begin{pmatrix}1&1&0\end{pmatrix}, \quad 群2：\begin{pmatrix}1&0&1\end{pmatrix}\] ２組の観測値から推定できるのは２つのパラメータまでなので、工夫が必要。
.

\[(c)\quad E(Y_{1k})=μとE(Y_{2k})=μ+α\] これは、群１を参照群とし群２の付加的効果をαとした。端点制約によるパラメータ化という方法。 \[β=\begin{pmatrix}μ\\α\end{pmatrix}, \quad 群1：\begin{pmatrix}1&0\end{pmatrix}, \quad 群2：\begin{pmatrix}1&1\end{pmatrix}\] .

\[(d)\quad E(Y_{1k})=μ+αとE(Y_{2k})=μ-α\] これは、２つの群を対称に扱う零和制約という方法で、αは群間の差を表す。 \[β=\begin{pmatrix}μ\\α\end{pmatrix}, \quad 群1：\begin{pmatrix}1&1\end{pmatrix}, \quad 群2：\begin{pmatrix}1&-1\end{pmatrix}\]

2.4.4 例：順序カテゴリーの説明変数

\[E(Y_{1k})=μ\\E(Y_{2k})=μ+α_1\\E(Y_{3k})=μ+α_1+α_2\] このとき \[β=\begin{pmatrix}μ\\α_1\\α_2\end{pmatrix},\quad群1：\begin{pmatrix}1&0&0\end{pmatrix},\quad群2：\begin{pmatrix}1&1&0\end{pmatrix},\quad群3：\begin{pmatrix}1&1&1\end{pmatrix}\] この表記方法により、α1は群１と比べた群２の相対的効果を表し、α2は群2と比べた群3の相対的効果を表すことができる。
.