「査読者から『多変量モデルでのステップワイズ法による自動変数選択は過学習やP-hackingを招くため不可』と指摘されました。そこで機械学習で流行りのLasso(ラッソ)を試したのですが、Rの glmnet パッケージで解析したら……あれ? P値も標準誤差(SE)も95%信頼区間も一切出力されません! オッズ比だけ出ても論文にどう書けばいいのでしょうか?」
近年、論文査読におけるステップワイズ法(Stepwise selection)への風当たりは急速に厳しくなっている。その代替案として、データサイエンティストや機械学習の分野で標準的に使われる「Lasso正則化ロジスティック回帰」を導入する研究者が増えている。
しかし、従来の多変量解析と同じ感覚で glmnet を実行すると、画面に表示されるのは「非ゼロの係数(オッズ比)」だけであり、医学論文の標準フォーマットである「$P$ 値」や「95%信頼区間」が影も形もないことにパニックに陥る。
Lassoの本質は、未知のデータに対する「予測精度(Prediction)」を高め「過学習(Overfitting)」を防ぐための正則化手法であり、医学統計が得意とする「因果推論(Causality)・主曝露因子のオッズ比算出」とは目的が根本的に異なる。
本記事では、なぜLassoでは $P$ 値が出ないのかという理論的背景、選択後推論(Post-selection inference)のバイアス、目的別の正しい使い分け、$\lambda$(ラムダ)の選び方(lambda.min vs lambda.1se)、主曝露変数を守る実践的なR実装、そして論文で使える英文テンプレートまでを徹底解説する。
1. 「ステップワイズはダメと言われたのでLassoを使いました」臨床医が直面するP値消失パニック
臨床研究や論文執筆の現場で、候補となる背景因子や臨床検査値が20〜30項目以上に及ぶ際、どの変数をモデルに入れるべきかという「変数選択」は避けて通れない課題である。
かつて多用されたステップワイズ法は、検定を何重にも繰り返すことで第一種の誤り率(偽陽性)が跳ね上がり、サンプルサイズに対して過剰に適合した過学習モデルを作り出すとして、現代の主要な医学誌の査読(Peer Review)ではほぼ禁止されつつある。
そこで「ステップワイズの代わりにLassoを使おう」とRの glmnet パッケージを実行すると、以下のような壁にぶつかる:
summary()を実行しても $P$ 値が出力されない。- 各変数の標準誤差(SE)や95%信頼区間(95% CI)が出ない。
- 変数が自動で削られて「0」になる仕組みは分かるが、残った因子の $P$ 値をどう論文に書けばいいか分からない。
Lassoは $P$ 値を出すためのツールではなく、「変数の係数をあえて縮小(シュリンク)させてモデルを簡素化する手法」である。まずは過学習と正則化の基本原理を直感的に理解しよう。
【コラム1】直感で理解する!「過学習」と「正則化」(過去問丸暗記 vs スパルタ指導者)
統計学や機械学習で登場する「過学習」と「正則化」は、受験勉強に例えると一発で理解できる。
- 過学習(Overfitting)=「過去問の答えを丸暗記した受験生」手元にある過去問(既存データ)を隅々まで解き明かし、問題のちょっとした印刷のカスや偶然の傾向まで100%完璧に覚えた状態。しかし、出題傾向が少し変わった本番の試験(未知の臨床データ)を受けると、応用が一切効かずにボロボロの点数を取る。これが過学習である。
- 正則化(Regularization)=「無駄な丸暗記を許さないスパルタ指導者」過去問を解かせるときに、「あまりにも複雑な解き方や丸暗記をした場合は、テストの点数を減点する(ペナルティを与える)」という制約を課す手法。あえてモデルを自由にフィットさせず足枷(ペナルティ)をはめることで、細かいノイズの丸暗記を防ぎ、本番の試験(未知のデータ)でも通用する「本質的な解法(シンプルで汎用的なモデル)」を学ばせる。
【コラム2】L1ペナルティの正体(減点主義の仕分け人)
正則化の中でも、Lassoが採用しているのが L1ペナルティ(各変数の回帰係数$\beta$の絶対値の合計に対するペナルティ)である。
通常のロジスティック回帰(
glm)は、データへのあてはまりの良さ(尤度)だけをひたすら追求する「甘やかしモデル」であり、変数を入れれば入れるほど過学習を起こす。これに対しLassoは、モデルの複雑さに対して「重税(L1ペナルティ)」を課す。モデルの中に大して役に立たない微小な効果の変数を抱えていると、その変数は「重税」に耐えきれなくなり、係数をちょうど「ゼロ($0$)」に潰されてモデルから強制退場させられる。
これが、Lassoによって特定の変数が自動的に消去され、重要な変数だけが手元に残る「自動変数選択(スパース性)」の仕組みである。
2. 臨床医が絶対に知っておくべきLassoの「4つの罠」
Lassoを臨床研究に適用する際、その特性を正しく理解していないと査読で重大なバイアスを指摘されることになる。絶対に知っておくべき「4つの罠」を整理する。
罠1:P値が出ない罠(理論的制約)
Lassoで $P$ 値が出ないのは glmnet パッケージの不備ではなく、統計理論上の必然である。
通常の $P$ 値や標準誤差は「無バイアス推定(真の値に平均して一致する推定値)」を前提として計算される。しかしLassoは、ペナルティを課すことで推定値にあえて「バイアス」を注入し、その代わりに関数のバリアンス(ばらつき)を激減させる手法である。意図的にバイアスを含ませているため、従来の「無バイアスを前提とした標準誤差や $P$ 値」は数学的に直接定義されない。
罠2:「Lassoでスクリーニング $\rightarrow$ 普通の glm()」の罠(選択後推論)
実務や一部の論文で頻繁に見かけるのが、「Lassoで残った非ゼロの変数だけを抽出し、その変数セットを使って改めて普通の glm() を回して $P$ 値やオッズ比を算出する」という2段階アプローチ(Two-stage approach)である。
一見スマートに見えるが、これは 選択後推論(Post-selection inference)のバイアス を生む。「変数選択に使ったデータ」と「$P$ 値の計算に使ったデータ」が同一であるため、得られる $P$ 値は不当に小さく(過剰に有意に)なり、オッズ比は過大評価される。この限界(Limitation)を理解した上で感度分析として提示する姿勢が必要となる。
罠3:目的の混同の罠(「予測 Prediction」 vs 「因果推論 Causality」)
データサイエンス(機械学習)と医学統計では、モデル構築の目的が根本的に異なる。
- データサイエンティスト(機械学習): 未知データに対する「予測精度(Out-of-sample prediction / C-index向上)」が目的。
- 臨床医(医学統計): 「特定のリスク因子(主曝露)と疾病の因果関係や多変量調整オッズ比の解釈」が目的。
目的が「因果推論」であるにもかかわらず無条件にLassoを回すと、最も評価したい主曝露変数(Main exposure)や重要な交絡要因(Confounder)が、Lassoのペナルティによって $0$ に潰されて消滅する悲劇が起きる。
罠4:標準化(Standardization)未考慮の罠
Lassoは係数 $\beta$ の絶対値の大きさにペナルティを課すため、変数の単位(例:年齢[歳]と性別[0/1]、血圧[mmHg]など)のスケールが揃っていないと、単位の大きい変数が不当に強くペナルティを受けて消されてしまう。解析前の標準化(Standardization)または glmnet 内部での自動標準化(standardize = TRUE)が必須である。
※補足(関連記事):
リッジ(Ridge)、ラッソ(Lasso)、エラスティックネット(Elastic Net)の違いや、Rでの基本的な実行比較については Ridge・Lasso・Elastic NetのRでの実行方法ガイド を参照されたい。
3. 目的で180度変わる!Lassoの正しい使い分けと実践戦略
Lassoを使用する際は、研究の目的が「予測モデルの構築」なのか「特定因子のオッズ比推定(因果推論)」なのかを明確に区別しなければならない。
【研究の目的はどちらか?】
│
┌──────────────────────┴──────────────────────┐
▼ ▼
【A: 予後予測モデル構築】 【B: 主曝露因子のオッズ比推定】
(Prognostic Risk Score) (Causal Inference)
│ │
・P値に固執しない ・無条件のLassoはNG
・C-index / Calibration で評価 ・主曝露変数は penalty.factor = 0 で保護
・変数選択は lambda.1se を採用 ・Double Selection Lasso を検討
シナリオ A:目的が「予後予測モデル構築(Prognostic Model / Risk Score)」の場合
- Lassoは最強の正攻法: 過学習を防ぎ、スコア化に適したコンパクトなモデルを作ることができる。
- 評価軸: $P$ 値の提示に執着せず、C-index(AUC) や Calibration plot(校正プロット) でモデルの予測精度を評価・提示する。
- 閾値の選択: 予測誤差が最小となる
lambda.minよりも、1標準誤差ルールを適用したlambda.1seを選ぶのが医学論文での世界標準アプローチである。
シナリオ B:目的が「因果推論(特定因子のオッズ比推定)」の場合
- 主暴露変数の保護: 評価したい主曝露因子(例:新治療の有無)がLassoで消されないよう、特定変数のペナルティ係数を $0$ に設定する(
penalty.factor = 0)。 - Double Selection Lasso の検討: 曝露変数とアウトカムの両方に関連する共変量をLassoでそれぞれ選出し、最終的な多変量モデルに投入する統計的手法(Double Selection)を検討する。
【コラム3】$\lambda$(ラムダ)の選び方:
lambda.minvslambda.1seLassoを実行する際、ペナルティの強さを決定するハイパーパラメータが$\lambda$(ラムダ)である。クロスバリデーション(
cv.glmnet)を実行すると、主に2つの$\lambda$が算出される。
lambda.min(ギリギリまで欲張る派):クロスバリデーションで予測誤差(Deviance)が最も小さくなったポイント。予測精度を極限まで追求するため、変数は多めに残りやすい。lambda.1se(安全第一・シンプル派):最小予測誤差から「1標準誤差(1 Standard Error)ぶんペナルティを強めた(変数をさらに絞り込んだ)」ポイント。予測精度をほとんど落とすことなく、より簡素で過学習しにくい堅牢なモデルが得られる。医学論文ではモデルの簡素性と再現性を重視するためlambda.1seの採用が推奨されることが多い。
4. 【実践Rコード】glmnet パッケージによるLasso正則化ロジスティック回帰(Step by Step)
ここからは、Rの glmnet パッケージを用い、主曝露変数を保護しながら lambda.1se でモデルを構築し、オッズ比と予測精度(C-index)を評価する手順を Step by Step で解説する。
Step 1: 疑似臨床データ(高次元・多重共線性を含む)の作成
患者300人分、アウトカム $Y$(発症: 1/0)、主曝露変数 $X_1$(exposure)、および共変量20項目のデータを生成する。
library(glmnet)
library(pROC)
library(tidyverse)
# 乱数シードの設定
set.seed(123)
n <- 300
p <- 20 # 共変量の数
# 説明変数の行列を作成
X_cov <- matrix(rnorm(n * p), nrow = n, ncol = p)
colnames(X_cov) <- paste0("X", 1:p)
# 主曝露変数 (exposure: 0/1)
exposure <- rbinom(n, size = 1, prob = 0.5)
# 真のロジット式(exposure, X1, X2 のみが真の影響を持つ)
logit_p <- -1.0 + 0.8 * exposure + 0.6 * X_cov[, 1] - 0.7 * X_cov[, 2]
prob <- 1 / (1 + exp(-logit_p))
y <- rbinom(n, size = 1, prob = prob)
# データフレームの統合
df_clinical <- as_tibble(X_cov) %>%
mutate(y = y, exposure = exposure)
# データの確認
head(df_clinical)
Step 2: 行列形式への整形と penalty.factor の設定
glmnet 関数は formula 表記(y ~ x)ではなく、説明変数の行列(matrix)と目的変数のベクトル(vector)を要求する。
また、主曝露変数 exposure がLassoによって $0$ に潰されないよう、penalty.factor を設定する(保護したい変数は $0$、ペナルティを課す変数は $1$ とする)。
【査読コメント】model.matrix() で行列を作成する際、model.matrix が出力する列の並び順と penalty.factor ベクトルの順序が寸分違わず一致していることを必ず確認してください。位置がズレていると、保護するつもりのない変数のペナルティが0になってしまいます。
# 説明変数の行列 (x) と目的変数ベクトル (y) の準備
# model.matrix でダミー変数化された行列を作成(切片列は除く)
x_mat <- model.matrix(y ~ . - 1, data = df_clinical)
y_vec <- df_clinical$y
# penalty.factor の設定(exposure の位置を 0 にしてペナルティ除外)
p_fac <- rep(1, ncol(x_mat))
p_fac[colnames(x_mat) == "exposure"] <- 0
# 設定の確認
p_fac
Step 3: クロスバリデーション(cv.glmnet)による最適 $\lambda$ の探索
10-fold Cross-Validation を実行し、最適な $\lambda$(lambda.min および lambda.1se)を探索する。
cv.glmnet は分割の仕方(10-fold)によって選ばれる lambda や残る変数がわずかに変動する。論文用の解析の再現性を保証するために、cv.glmnet を呼び出す直前に set.seed() を明記しておくことが必須である。
# クロスバリデーションの実行 (ロジスティック回帰のため family = "binomial")
set.seed(123) # 再現性のため必須
cv_fit <- cv.glmnet(
x = x_mat,
y = y_vec,
family = "binomial",
alpha = 1, # alpha = 1 は Lasso (alpha = 0 は Ridge)
penalty.factor = p_fac, # 主曝露変数の保護
standardize = TRUE # 内部での標準化を実行
)
# CVプロットの描画
plot(cv_fit)
プロット上の左側の点線が lambda.min、右側の点線が lambda.1se を示す。

Step 4: lambda.1se での非ゼロ係数(オッズ比)の取り出し
lambda.1se のポイントで残った変数を抽出し、係数 $\beta$ およびオッズ比($\exp(\beta)$)を算出する。
# lambda.1se での非ゼロ係数を抽出
coef_1se <- coef(cv_fit, s = "lambda.1se")
# スパース行列からデータフレームへ変換
df_coef <- tibble(
var = rownames(coef_1se),
coef = as.vector(coef_1se)
) %>%
filter(coef != 0) %>% # 0でない(選択された)変数のみ抽出
mutate(odds_ratio = exp(coef))
print(df_coef)
出力結果
> print(df_coef)
# A tibble: 4 × 3
var coef odds_ratio
<chr> <dbl> <dbl>
1 (Intercept) -0.776 0.460
2 X1 0.0820 1.09
3 X2 -0.176 0.839
4 exposure 0.669 1.95
主曝露変数 exposure は penalty.factor = 0 により確実に保護されつつ、真に関連のある変数(X1, X2 )のみが正しく選択されたことが確認できる。
Step 5: C-index(AUC)によるモデル予測精度の評価
算出されたLassoモデルの予測確率を求め、pROC パッケージを用いて C-index(ROC-AUC)とその95%信頼区間を計算する。
# lambda.1se での予測確率の算出
pred_prob <- predict(cv_fit, newx = x_mat, s = "lambda.1se", type = "response")
# ROC曲線の作成と C-index (AUC) の算出
roc_obj <- roc(y_vec, as.vector(pred_prob))
auc_value <- auc(roc_obj)
auc_ci <- ci.auc(roc_obj)
cat(sprintf("C-index (AUC): %.3f (95%% CI: %.3f - %.3f)\n", auc_value, auc_ci[1], auc_ci[3]))
出力結果
> cat(sprintf("C-index (AUC): %.3f (95%% CI: %.3f - %.3f)\n", auc_value, auc_ci[1], auc_ci[3]))
C-index (AUC): 0.681 (95% CI: 0.621 - 0.741)
5. 論文でそのまま使える英文テンプレート(Methods & Results)
論文の Methods および Results セクションに記述する標準的な英文表現テンプレートを示す。
Methods セクション英文テンプレート
“To avoid overfitting and address potential multicollinearity among a large set of candidate covariates, variable selection was performed using Lasso (Least Absolute Shrinkage and Selection Operator) penalized logistic regression via the
glmnetpackage in R.To ensure that the primary exposure variable [e.g., treatment group] remained in the multivariable model, its penalty factor was set to zero (
penalty.factor = 0). Optimal penalization parameter ($\lambda$) was determined using 10-fold cross-validation, and the one-standard-error rule (lambda.1se) was selected to build the most parsimonious model. Model discrimination was evaluated using the C-index (area under the receiver operating characteristic curve [AUC]).”
Results セクション英文テンプレート
“Out of [20] candidate covariates, Lasso logistic regression with
lambda.1seselected [2] independent predictors ([X1] and [X2]) in addition to the primary exposure variable.After regularized adjustment, the primary exposure remained significantly associated with the outcome (regularized Odds Ratio [OR], [1.95]). The final parsimonious model demonstrated robust discriminative ability, achieving a C-index of [0.681] (95% CI, [0.621] to [0.741]).”
6. まとめ&一括コピペ用Rスクリプト
本記事のポイントおさらい
- LassoでP値が出ないのは当然の仕様: 正則化でバイアスを注入して過学習を防ぐ手法であるため、無バイアス前提の $P$ 値や95% CIは直接算出されない。
- 「Lasso $\rightarrow$ 通常の
glm()」の安易な実行に注意: 選択後推論(Post-selection bias)が生じるため、そのままの $P$ 値を真に受けてはならない。 - 主暴露変数は
penalty.factor = 0で保護する: 因果推論目的の解析では、消されては困る重要な変数のペナルティを解除して回す。 lambda.1seの採用が世界標準: 最小誤差のlambda.minよりも、1標準誤差ルールを適用したlambda.1seの方が簡素で堅牢な医学モデルになる。
一括実行用Rスクリプト
以下のコードをコピー&ペーストすることで、疑似データ作成から penalty.factor 設定、CV実行、lambda.1se でのオッズ比抽出、C-index評価までを一気通貫で実行できる。
# ==============================================================================
# Lasso正則化ロジスティック回帰(主曝露変数保護 + lambda.1se 選択)一括実行スクリプト
# パッケージ: glmnet, pROC, tidyverse
# ==============================================================================
if (!requireNamespace("glmnet", quietly = TRUE)) install.packages("glmnet")
if (!requireNamespace("pROC", quietly = TRUE)) install.packages("pROC")
if (!requireNamespace("tidyverse", quietly = TRUE)) install.packages("tidyverse")
library(glmnet)
library(pROC)
library(tidyverse)
# 1. 疑似データの作成
set.seed(123)
n <- 300
p <- 20
X_cov <- matrix(rnorm(n * p), nrow = n, ncol = p)
colnames(X_cov) <- paste0("X", 1:p)
exposure <- rbinom(n, size = 1, prob = 0.5)
logit_p <- -1.0 + 0.8 * exposure + 0.6 * X_cov[, 1] - 0.7 * X_cov[, 2]
prob <- 1 / (1 + exp(-logit_p))
y <- rbinom(n, size = 1, prob = prob)
df_clinical <- as_tibble(X_cov) %>%
mutate(y = y, exposure = exposure)
# 2. 行列整形と penalty.factor 設定(exposure を保護)
x_mat <- model.matrix(y ~ . - 1, data = df_clinical)
y_vec <- df_clinical$y
p_fac <- rep(1, ncol(x_mat))
p_fac[colnames(x_mat) == "exposure"] <- 0
# 3. cv.glmnet による 10-fold CV 実行
set.seed(123)
cv_fit <- cv.glmnet(
x = x_mat,
y = y_vec,
family = "binomial",
alpha = 1,
penalty.factor = p_fac,
standardize = TRUE
)
# 4. lambda.1se での非ゼロ係数・オッズ比の算出
coef_1se <- coef(cv_fit, s = "lambda.1se")
df_results <- tibble(
Variable = rownames(coef_1se),
Coefficient = as.vector(coef_1se)
) %>%
filter(Coefficient != 0) %>%
mutate(Odds_Ratio = exp(Coefficient))
cat("\n--- Lasso 選択変数およびオッズ比 (lambda.1se) ---\n")
print(df_results)
# 5. C-index (AUC) による予測精度の評価
pred_prob <- predict(cv_fit, newx = x_mat, s = "lambda.1se", type = "response")
roc_obj <- roc(y_vec, as.vector(pred_prob))
auc_ci <- ci.auc(roc_obj)
cat(sprintf("\nModel C-index (AUC): %.3f (95%% CI: %.3f - %.3f)\n", auc(roc_obj), auc_ci[1], auc_ci[3]))





コメント