井出草平の研究ノート

高次元スパース線形モデルにおけるパラメータの信頼区間を構築するためのブートストラップLasso+部分リッジ法[R]

stat.paperswithcode.com

HDCI(High-Dimensional Confidence Intervals)はは、Lasso回帰とRidge回帰の欠点を補うという点にある。

  • Lasso回帰は、変数選択を行う際に多くの変数の係数をゼロにするため、少ししか関連性のない変数が推計から排除される可能性がある。
  • Ridge回帰は、すべての変数を選択するため、多くの無関係な変数がモデルに含まれることになる。

HDCIは、まずLasso回帰で重要な変数を選択し、その後Partial Ridge回帰を用いてこれらの変数の係数を再推定することで、Lassoの変数選択の強みとRidgeのバイアス補正の強みを組み合わせてこの欠点を補正する方法といえる。

特に高次元スパース線形モデルの係数に対する信頼区間を構築するためという点が大きいようだ。高次元スパース線形モデルとは、説明変数の数が観測データの数よりも多い(高次元)場合において、多くの説明変数の係数がゼロまたはゼロに近い(スパース)状況での線形回帰モデルである。こうしたモデルでは、説明変数の多さから過学習のリスクが高まり、モデルの解釈や計算が難しくなる。

HDCIこの方法は次のステップで構成される。

  1. 変数選択: Lasso回帰を用いて重要な変数を選択する。
  2. 偏り補正: Partial Ridge回帰を使用して選択された変数の係数を再推定し、偏りを補正する。
  3. 信頼区間の構築: 推定された係数に基づいて信頼区間を構築する。

「小さな非ゼロ係数が存在する場合」とは、回帰モデルで予測変数の係数が非常に小さいがゼロではない状況を指す。これは、その変数がわずかにではあるがモデルに影響を与えていることを示す。例えば、Lasso回帰では、モデルは重要でない特徴量の係数をゼロに近づける一方で、わずかに関連性のある特徴量(feature)には小さな非ゼロ係数を持たせることがある。

Lasso回帰では、少ししか影響を持っていない独立変数も残ることがある。これは、ペナルティ項を追加することで多くの特徴量の係数をゼロに近づけるが、完全にゼロにならない係数も存在するためである。非常に小さいながらも統計的に有意な変数が残ることがあり、Lasso回帰は変数選択と正則化の効果を同時に持つ手法とされる。

HDCI(High-Dimensional Confidence Intervals)法は、高次元データにおける信頼区間を推定する方法であり、小さな非ゼロ係数の存在やサンプルサイズが少ない場合でも有効な推定を提供する。変数選択バイアスを軽減し、信頼区間の正確性を向上させることを目的としている。

信頼区間カバレッジ確率(coverage probability)は、統計的推定において繰り返しサンプリングを行った場合に、真のパラメータがその信頼区間内に含まれる割合を指す。例えば、95%の信頼区間の場合、理論上は100回のサンプリングのうち95回は真のパラメータがその区間内に含まれることを期待する。カバレッジ確率が高いほど、信頼区間が真のパラメータを含む確率が高いことを示す。

Partial Ridge回帰(Partial Ridge Regression)は、リッジ回帰(Ridge Regression)を用いて回帰モデルの一部の係数を再推定する手法である。これは、通常のリッジ回帰の全体的な収縮効果を避け、特定の変数に対する信頼性の高い推定を提供することを目的としているようだ。特に、高次元スパースデータに対するLasso回帰の後処理として使用され、Lassoによって選択された変数に対してリッジ回帰を適用することで、バイアスを減少させ、より安定した推定を実現するというのがHDCIモデルの特徴である。

Rのレファレンス

https://cran.r-project.org/web/packages/HDCI/index.html

https://www.quantargo.com/help/r/latest/packages/HDCI/DESCRIPTION

データの読み込みと前処理

library(haven)
library(tidyr)
auto <- haven::read_dta("http://www.stata-press.com/data/r9/auto.dta")
auto <- auto %>% drop_na()  # 全カラムに対してNAがない行を抽出
auto$foreign <- as.integer(auto$foreign)

data <- auto[, c("price", "mpg", "rep78", "headroom", "trunk", "weight", "length", "turn", "displacement", "gear_ratio", "foreign")]
X <- as.matrix(data[, -1])
y <- data$price

1. 変数選択: Lasso回帰を用いて重要な変数を選択する。

library(glmnet)
library(HDCI)
set.seed(0)
lasso_model <- Lasso(x = X, y = y, fix.lambda = FALSE)
selected_variables <- which(lasso_model$beta != 0)

Lassoで選択された変数

print(colnames(X)[selected_variables])

"rep78" "headroom" "weight" "displacement" "foreign"

Lassoでは5つの変数が選ばれている。

2. 偏り補正: Partial Ridge回帰を使用して選択された変数の係数を再推定し、偏りを補正する。

Lasso回帰は、ペナルティ項を追加して多くの特徴量(変数)の係数をゼロに近づけ、重要な変数を選択する。しかし、このプロセスで選択されなかった変数も、完全に無視されるわけではない。Partial Ridge回帰では、選択された変数に対してはペナルティを加えずに再度回帰を行い、選択されなかった変数に対してL2正則化(リッジ回帰)を適用する。この手法により、選択されなかった変数も含めて推定が行われるため、結果としてLassoで選択されなかった変数の係数も非ゼロとなることがある。

partial_ridge_model <- PartRidge(x = as.matrix(data[, -1]), y = data$price, varset = selected_variables, lambda2 = 1/nrow(data))

Partial Ridge回帰の係数

print(partial_ridge_model$beta) 
         mpg        rep78     headroom        trunk       weight       length         turn displacement   gear_ratio      foreign 
    2.400132   221.874261  -672.899304     1.815335     2.397789    -6.081407   -32.300492    15.196493  -110.433805  3476.013031 

Lassoでは選ばれていなかった変数にも係数がついている。HDCIでPartial Ridge回帰による補正とはこのようなことである。

Partial Ridge回帰の切片

print(partial_ridge_model$beta0)
[1] -1242.533

3. 信頼区間の構築: 推定された係数に基づいた信頼区間

ci_result <- bootLOPR(x = as.matrix(data[, -1]), y = data$price, lambda2 = 1/nrow(data), B = 500, type.boot = "residual")

信頼区間の上下のみを表示する。

variables <- c("mpg", "rep78", "headroom", "trunk", "weight", "length", "turn", "displacement", "gear_ratio", "foreign")
colnames(ci_result$interval.LPR) <- variables
print(ci_result$interval.LPR)

結果。

           mpg     rep78     headroom     trunk    weight      length        turn displacement gear_ratio foreign
[1,] -34.42574 -403.9487 -1346.343118 -170.2001 0.8143867 -66.7195611 -236.258021     1.656846  -479.0976 1953.75
[2,] 132.72553  580.0130     5.245458  169.9695 4.4926832  -0.6497598    8.228028    27.383460  1721.1198 5140.85

95%信頼区間の真ん中にPartial Ridge回帰の推定量がくるわけではない点には多少注意が必要かもしれない。