ヘキット・モデルとは？

2023.03.07

ヘキット・モデル

（1）サンプル・セレクション・バイアス

　標本から素直に推定すると、分析結果が偏ってしまう場合があります。これをサンプル・セレクション・バイアスといいます。

　サンプル・セレクション・バイアスの例には、下図があります。

　このとき、母集団では「Xが大きいとき、Yが大きい」という関係性があるのに、標本では「Xが大きいとき、Yが小さい」という関係性が生まれてしまっています。この状態で単回帰分析をすると、深刻なサンプル・セレクション・バイアスが発生します。

（2）ヘキット・モデル

　ヘキット・モデルは、目的変数Yだけが欠落してしまったというサンプル・セレクション・バイアスに対処できるモデルです。「見えないものを見ようとする」モデルであり、個人的にとても感動したモデルです。ヘキット・モデルは、母集団では

$$Y=\beta_0+\beta_1X+u$$

でも、観測された標本では

$$Y=\beta_0+\beta_1X+\gamma 逆ミルズ比+\epsilon$$

であるモデルです。このとき、逆ミルズ比を含めて回帰しないと欠落変数バイアスが発生します。

　説明不足の点は後述します。

（3）バイアスが起きる条件

　サンプル・セレクションは即座にバイアスを引き起こしません。例えば、母集団が1万個、無作為に1000個抽出した場合、バイアスは発生しません。いわゆる「無作為抽出」だからです。無作為抽出は、バイアスを産まないサンプル・セレクションの最も有名な例です。

　問題は「どのような場合にサンプル・セレクションが推定結果を歪める（＝バイアスを引き起こす）のか？」になります。

　ヘキット・モデルによれば、γ≠0のとき、サンプルのみで分析すると、バイアスが発生します。逆ミルズ比が欠落変数となるからです。

$$Y_{観測有}=\beta_0+\beta_1X+\gamma 逆ミルズ比+\epsilon$$

　したがって、「サンプル・セレクション・バイアスがある」とは、「γについてt検定し、帰無仮説γ=0が棄却され、対立仮説γ≠0が採択された場合」と言えます。

（4）逆ミルズ比

　逆ミルズ比とは、確率密度÷累積分布で表されるものです。ただし、ヘキット・モデルでは、標準正規分布に従うvがc以上のときのvの期待値を、逆ミルズ比を捉えてもらって構いません。

$$v　〜　N(0,1)$$

$$逆ミルズ比=E( v|v≧c)$$

　例えば、c=-1のとき、-1以下の値を取らないvの分布は赤になります。計算すると、このときの逆ミルズ比は約0.288になります。

　c=0.5のとき、0.5以下の値を取らないvの分布は赤になります。計算すると、このときの逆ミルズ比は約1.141になります。

　なお、ヘキット・モデルにおける逆ミルズ比は

$$逆ミルズ比=E( v|v≧c)=\frac{\phi(c)}{1-\Phi(c)}$$

$$\phi(c)：v=cにおける標準正規分布の確率密度$$

$$\Phi(c)：v=cまでの標準正規分布の累積分布$$

です。

潜在変数モデル

　ヘキット・モデルは、線形回帰モデルとプロビット・モデルが合わさったモデルです。したがって、一部を潜在変数モデルで解釈するとわかりやすいです。

（1）母集団でのモデル

　実際の現象は

$$Y=\beta_0+\beta_1X+u$$

です。ただし

$$Y：目的変数、X：説明変数、u：誤差項$$

$$\beta_0，\beta_1：構造式の回帰係数$$

$$E(u|X)=0$$

です。これを構造式といいます。

（2）サンプル・セレクション

　観測される場合をS=1、されない場合をS=0と考えます。そして、

$$S^*≧0　ならば　S=1$$

$$S^*<0　ならば　S=0$$

とします。ヘキット・モデルでは、サンプル・セレクションに関わる潜在変数S*は次に従います。

$$S^*= \delta_0+\delta_1X +\delta_2 Z+v$$

　ただし

$$ Z：セレクションにのみ関わる変数，v：誤差項$$

$$\delta_0，\delta_1，\delta_2：セレクション式の回帰係数$$

$$E(v|X)=0$$

です。これをセレクション式といいます。

　また、ヘキット・モデルでは、誤差項vが標準正規分布が成り立っていると仮定します。

$$v　〜　N(0,1)$$

　つまり、セレクション式はプロビット・モデルと同じです。

（3）標本でのモデル

　結論を言えば、観測されるYは

$$Y_{観測有}=\beta_0+\beta_1X+\gamma E( v|v≧c)+\epsilon$$

になります。ただし、誤差項uとvが相関すること

$$E(u_i|v_i)=\gamma v_i$$

を追加的に仮定しており、cは

$$c=-(\delta_0+\delta_1X +\delta_2 Z)$$

です。また、E( v|v≧c)は逆ミルズ比といい、標準正規分布に従うvがc以上のときのvの期待値です。

$$逆ミルズ比=E( v|v≧c)=\frac{\phi(c)}{1-\Phi(c)}$$

（4）なぜ逆ミルズ比が出てくる？

　観測されるYの期待値が

$$E(Y|X,S=1)　=　\beta_0+\beta_1X+\gamma E( v|v≧c)$$

であることを説明します。まず、X、S=1のときのYの条件付き期待値は

$$E(Y|X,S=1)　=　\beta_0+\beta_X+E(u|X, S=1)$$

　uとXは独立なので

$$E(Y|X,S=1)　=　\beta_0+\beta_X+E(u|S=1)$$

　S=1の条件を潜在変数モデルにしたがって言い換えると

$$=\beta_0+\beta_1X+E(u|　\delta_0+\delta_1X +\delta_2 Z+v≧0)$$

　移項して

$$=\beta_0+\beta_1X+E(u|　v≧-\delta_0-\delta_1X -\delta_2 Z)$$

　誤差項vの条件がついたuの期待値はγvと仮定しており、vにも条件がついていることを考慮すると

$$=\beta_0+\beta_1X+\gamma E( v|　v≧-\delta_0-\delta_1X -\delta_2 Z)$$

です。ここで

$$c=-\delta_0-\delta_1X -\delta_2 Z$$

とおくと

$$=\beta_0+\beta_1X+\gamma E( v|v≧c)$$

になります。この最後の項が逆ミルズ比です。

　最後に逆ミルズ比がでてくるのは、サンプル・セレクションによってXと相関する誤差項が生まれてしまったからです。相関する誤差項を取り除いた新しい誤差項がεです。

$$Y_{観測有}=\beta_0+\beta_1X+\gamma E( v|v≧c)+\epsilon$$

Yの予測

　母集団でのYの期待値は

$$E(Y|X)=\beta_0+\beta_1X$$

です。ただし、標本でのYの期待値は

$$E(Y|X,S=1)=\beta_0+\beta_1X+\gamma E( v|v≧c)$$

$$c=-\delta_0-\delta_1X -\delta_2 Z$$

です。

X→Yの効果

　Xが微小1単位大きいと、Yが平均的に大きい量を、限界効果といいます。

$$限界効果＝\frac{d E(Y|X)}{d X}$$

　母集団にて説明変数Xが目的変数Xに与える限界効果は、β1です。なお、サンプル・セレクション・バイアスがあるので

$$\frac{d E(Y|X,S=1)}{d X}≠\beta_1$$

なことに注意してください。

　なお、この限界効果を因果関係と解釈するには、統計的因果推論の考え方が必要で留保すべきです。

Rコード

（1）ヘキット・モデル

#パッケージのインストールと呼び出し
install.packages("sampleSelection")　#PCにつき1回やればよい
library(sampleSelection)

#ヘックマンの2段階推定
lm0 <- heckit(s ~ x + z, #セレクション式
                y ~ x,  #構造式
                method = "2step", 　#2段階推定
                data=?) 　#データ「？」を使う
summary(lm00)
　
#最尤法
lm00 <- heckit(s ~ x + z,  #セレクション式
                y ~ x,   #構造式
                method = "mls",  　#最尤法
                data=?) 　#データ「？」を使う
summary(lm00)

（2）パッケージ：truncorm

　作図では、打ち切り標準正規分布に関するパッケージを用います。↓

#パッケージのインストール
install.packages("truncnorm")
　
#パッケージの呼び出し
library(truncnorm)

　主な関数はこれです。↓

#母平均0、標準偏差1の正規分布を、aからbまでに限定したときの、xにおける確率密度を求める関数
dtruncnorm(x, a=-Inf, b=Inf, mean = 0, sd = 1)
　
#上の累積分布バージョン
ptruncnorm(x, a=-Inf, b=Inf, mean = 0, sd = 1)

　説明書はこちら→truncnorm: The Truncated Normal Distribution

（3）作図

#打ち切り点を-1に設定
c <- -1
　
#-5から5の間に1000個の連続する数字を作成
x <- seq(-5,5,length=1000)
　
#標準正規分布の確率密度を描画
plot(x, dnorm(x, mean=0, sd=1),  #横軸にx、縦軸に標準正規分布の確率密度
  type="l",　#折れ線グラフ
  lwd=2,　#線の太さを2
  xlim=c(-4,4),  #横軸を-4から4に
  ylim=c(0,0.5),  #縦軸を0から0.5に
  xaxp=c(-4,4,8))  #横軸を-4から4までを8つに分けて目盛りをつける
　
#打ち切り標準正規分布の確率密度を描画
lines(x, dtruncnorm(x, a=c, mean=0, sd=1), type="l", col="red",lwd=2)
　
#逆ミルズ比を求める
M <- dnorm(c,mean=0,sd=1)/(1-pnorm(c,mean=0,sd=1))
　
#逆ミルズ比を描画
abline(v=M,col="blue",lwd=2)

#打ち切り点を0.5に設定
c <- 0.5
x <- seq(-5,5,length=1000)
plot(x, dnorm(x, mean=0, sd=1),
  type="l",
  lwd=2,
  xlim=c(-4,4),
  ylim=c(0,1.5),
  xaxp=c(-4,4,8))
lines(x, dtruncnorm(x, a=c, mean=0, sd=1), type="l", col="red",lwd=2)
M <- dnorm(c,mean=0,sd=1)/(1-pnorm(c,mean=0,sd=1))
abline(v=M,col="blue",lwd=2)