この章はなにか

多くのDataScienceの問題に回帰モデルが用いられているからわかるように、回帰モデルはDataScienceの基本となるので紹介されている。

この節はなにか

回帰モデルで推定された係数の不確か性を考える

係数推定の不確か性

以下のairqualityデータで考える。

data("airquality")
# A tibble: 153 x 1
# airquality$Ozone $Solar.R $Wind $Temp $Month  $Day
# <int>    <int> <dbl> <int>  <int> <int>
#   1               41      190   7.4    67      5     1
# 2               36      118   8      72      5     2
# 3               12      149  12.6    74      5     3
# 4               18      313  11.5    62      5     4
# 5               NA       NA  14.3    56      5     5
# 6               28       NA  14.9    66      5     6
# 7               23      299   8.6    65      5     7
# 8               19       99  13.8    59      5     8
# 9                8       19  20.1    61      5     9
# 10               NA      194   8.6    69      5    10
# # … with 143 more rows

全変数を用いた重回帰モデルにおける、Windの係数などは以下となる。

fit = glm(Ozone ~ ., data = airquality)
summary(fit)$coef['Wind',]
# Estimate    Std. Error       t value      Pr(>|t|) 
# -3.318444e+00  6.445095e-01 -5.148789e+00  1.231276e-06

この結果より、この係数の期待値は-3.318444で、Std. Errorは0.6445095となる。

分析において、そのモデルの仮定が正しい場合は当たり前だがRのsummaryで出力されるStd. Error、つまり標準的なStd. Errorは正しいStd. Error推定となる。しかし、標準的なStd. Errorは仮定が正しくない場合は誤った推定結果となる。
そうなると、不確か性を表す指標Std. Errorが正確ではないので不確か性を考えることができなくなる。

実ビジネスにおいて、仮定が本当に正しいのか？ということは検証することができないので、1章で紹介したノンパラメトリックブートストラップのように仮定を必要としないロバストな手法をできるだけ使った方が良い。

標準的なRのsummaryで出力されるStd. Errorは均一分散を仮定した推定となっている。これは過去記事で紹介したように、均一分散の仮定を含む、ガウスマルコフ仮定の際にOLSが最も良い推定量となることからそのようになっている。しかし、均一分散の仮定は多くの場合満たされない。

knknkn.hatenablog.com

HC standard error

OLSにおいてのみ使えるHeteroskedastic Consistent(HC) standard errorというロバストな手法が存在する。

HC standard errorでは均一分散を仮定として置かないよりロバストなstandard errorを計算する。

HC standard errorはAER libraryのsandwich を用いることで計算できる。まず、sandwich::vcovHCを用いると各変数のHC共分散が出るので、その中で自分同士のHC共分散が通常のHC分散になるので、HC standard errorはその値のルートを取ることで計算できる。

library(AER)

bvar = sandwich::vcovHC(fit)
round(bvar, 1)
# (Intercept) Solar.R  Wind Temp Month  Day
# (Intercept)       432.9     0.1 -13.3 -3.6  -3.2  0.3
# Solar.R             0.1     0.0   0.0  0.0   0.0  0.0
# Wind              -13.3     0.0   0.8  0.1  -0.2 -0.1
# Temp               -3.6     0.0   0.1  0.1  -0.1  0.0
# Month              -3.2     0.0  -0.2 -0.1   1.8  0.0
# Day                 0.3     0.0  -0.1  0.0   0.0  0.1

sqrt(bvar['Wind','Wind'])
# 0.9128877

通常のstandard errorの値は0.6445095となったが、HC standard errorはそれよりやや大きい0.9128877となった。これは前述のように、共分散を仮定として置いていないことが理由となる。

ちなみに、書籍ではHC standard errorと書かれているがWhite standard errors、日本語だとWhiteの頑強な（ロバスト）標準誤差の方がよく呼ばれている名前な気がします。

また、先程1章のようにノンパラメトリックブートストラップのような手法がロバストと書いたが、実際にノンパラメトリックブートストラップで求めた分散は以下より、0.8676731となりHC standard errorと近い値となる。

B = 10000
beta = vector(length=B)
n = nrow(airquality)
for (b in 1:B){
  bs = sample.int(n, n, replace = TRUE)
  bsfit = lm(Ozone ~ ., data = airquality, subset = bs)
  beta[b] = coef(bsfit)['Wind']
}

sd(beta) # 0.8676731

独立の仮定

先程の均一分散の仮定と同様に、OLSの標準的な出力にはデータが独立という仮定も含まれている。このとき、standard errorを不均一分散を考慮したHC standard errorを紹介したが、データ同士に相関がある場合を考慮したclustered robust standard errorというものがある。

これは、同じクラスタ（グループ）のデータには相関があるが、クラスターが別だと関係性が発生しない場合に機能する。例えば、いろいろな地域のなにかしらの値が複合されたデータがにおいて、同じ地域で値は相関するが、別の地域と比べると相関がないという状況はよくある。その場合はstandard errorを計算するときはクラスタをその地域にしたclustered robust standard errorとする。