💡 この記事でわかること

・ ベイズ統計における解析の流れ（事前分布 → 尤度 → 事後分布）とMCMCアルゴリズムの仕組み
・ Rの brms パッケージを使ったベイズ回帰モデルの実装方法（コード・出力・解釈セット）
・ SASの PROC MCMC を使ったベイズ解析の実装方法
・ 製薬業界でベイズ統計が活用される具体的な場面（アダプティブデザイン・PK/PD・バスケット試験など）
・ 規制当局（FDA・ICH）のガイダンスにおけるベイズ手法の位置づけ

はじめに

ベイズ統計は、従来の頻度論的アプローチとは根本的に異なる統計的推論の枠組みです。事前情報をデータで更新し、パラメータの確率的な表現（事後分布）を得ることができるため、情報の蓄積と柔軟な意思決定が求められる製薬業界との相性が非常に高いと言えます。

当ブログでは、これまでにベイズ統計の基礎理論として、ベイズ統計入門：頻度論との違いとベイズの定理や事前分布・事後分布・尤度の関係と共役事前分布を解説してきました。本記事では「実際にどう解析するのか」「製薬業界のどの場面で使われているのか」に焦点を当て、RとSASのコードを交えて体系的に解説します。

対象読者は、ベイズ統計の基本概念は理解しているが実装経験が少ない製薬企業の生物統計担当者や、統計を学ぶ社会人・学生です。

ベイズ統計解析の基本的な流れ

ベイズ統計の解析は、次の4つのステップで進みます。

ステップ1：事前分布の設定
解析前の知識・仮定をもとに、パラメータ θ に対する事前分布 p(θ) を設定します。過去の試験結果や文献情報を反映させることができます。

ステップ2：尤度の定式化
データ x が生成されるモデル p(x|θ) を設定します。通常の統計モデリングと同様の作業です。

ステップ3：事後分布の計算
ベイズの定理に従い、事後分布を計算します。

p(θ|x) ∝ p(x|θ) · p(θ)

ここで p(θ|x) は事後分布（データ観測後のパラメータの確率）、p(x|θ) は尤度、p(θ) は事前分布を表します。

ステップ4：事後分布からの推論
事後分布から、事後平均・事後中央値・信用区間（Credible Interval）・仮説の確率などを計算します。

💡 ポイント
実際の問題では事後分布が解析的に求められないケースが大半です。その場合、MCMC（マルコフ連鎖モンテカルロ法）などの数値的サンプリング手法を用いて事後分布からサンプルを生成し、推論を行います。

MCMCアルゴリズムの概要

MCMCは事後分布から直接サンプリングできない場合に、事後分布に収束するマルコフ連鎖を構成してサンプルを生成する手法です。主要なアルゴリズムを以下の表で整理します。

アルゴリズム概要主な実装
Metropolis-Hastings法提案分布から候補サンプルを生成し、採択率で受け入れを判断する最も基本的な手法汎用（SAS PROC MCMCでも利用可）
ギブスサンプリング他のパラメータを固定した条件付き分布から各パラメータを順番にサンプリングSAS PROC MCMC（デフォルト）
HMC / NUTS勾配情報を利用したハミルトニアンモンテカルロ法。高次元パラメータで高効率Stan（R: brms / rstanarm のバックエンド）
📝 解釈・補足
MCMCでは最初のサンプルはマルコフ連鎖が定常分布に収束する前の「バーンイン（warm-up）期間」のものであり、推論には使用しません。Stanでは warmup 引数、SASでは nbi（Number of Burn-In iterations）で指定します。

Rによるベイズ解析の実装（brms）

ここでは、製薬研究でよく見られる「用量と薬効の関係」をベイズ線形回帰でモデル化する例を示します。brms パッケージはStanをバックエンドに使用し、比較的シンプルな構文でベイズ回帰モデルを実装できます。

パッケージの準備とデータ生成

# インストール（初回のみ）
# install.packages("brms")

library(brms)
library(bayesplot)

set.seed(2024)

# 用量-反応データの生成（n=60）
n <- 60
dose <- runif(n, min = 0, max = 100)
response <- 3 + 0.6 * dose + rnorm(n, mean = 0, sd = 12)
df_pharma <- data.frame(dose = dose, response = response)

事前分布の設定とモデルの実行

事前分布は以下のように設定します。切片・傾きには広い正規分布（弱情報事前分布）、残差標準偏差には半コーシー分布を使用します。

priors_custom <- c(
  prior(normal(0, 100), class = "Intercept"),
  prior(normal(0, 10),  class = "b"),
  prior(cauchy(0, 5),   class = "sigma")
)

fit_bayes <- brm(
  formula  = response ~ dose,
  data     = df_pharma,
  family   = gaussian(),
  prior    = priors_custom,
  chains   = 4,
  iter     = 4000,
  warmup   = 2000,
  seed     = 2024,
  cores    = 4
)

出力結果

> summary(fit_bayes)

 Family: gaussian
  Links: mu = identity; sigma = identity
Formula: response ~ dose
   Data: df_pharma (Number of observations: 60)
  Draws: 4 chains, each with iter = 4000; warmup = 2000

Population-Level Effects:
          Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept     3.47      3.21    -2.83     9.81 1.00     7842     6015
dose          0.59      0.05     0.49     0.70 1.00     7621     5988

Family Specific Parameters:
      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma    11.84      1.13     9.79    14.26 1.00     7294     5722
📝 解釈
dose の事後平均は 0.59（95%信用区間：0.49〜0.70）で、用量が1単位増えると薬効が約0.59増加することが示されています。
Rhat はすべて 1.00 に近く、4本のMCMCチェーンが収束していることを確認できます。Bulk_ESS・Tail_ESS（有効サンプルサイズ）もいずれも5,700以上と十分で、推定の精度が高いと判断できます。

収束診断とビジュアライゼーション

# トレースプロット（収束確認）
mcmc_trace(fit_bayes, pars = c("b_Intercept", "b_dose", "sigma"))

# 事後分布の可視化
mcmc_areas(fit_bayes, pars = c("b_dose"), prob = 0.95)

# 予測分布の可視化
plot(conditional_effects(fit_bayes))
⚠️ 注意
Rhat > 1.01 の場合、MCMCチェーンが収束していない可能性があります。iter を増やす・adapt_delta を上げる・モデルのスケーリング見直し等を試みてください。収束診断なしの結果をそのまま使用することは避けましょう。

SASによるベイズ解析の実装（PROC MCMC）

SASでは PROC MCMC プロシジャを使って同様のベイズ線形回帰を実装します。

/* 用量-反応データの生成 */
data pharma_data;
  call streaminit(2024);
  do i = 1 to 60;
    dose     = rand('uniform') * 100;
    response = 3 + 0.6 * dose + rand('normal') * 12;
    output;
  end;
  drop i;
run;

/* ベイズ線形回帰（PROC MCMC） */
proc mcmc data     = pharma_data
          nmc      = 8000
          nbi      = 2000
          thin     = 1
          seed     = 2024
          plots    = all
          statistics = all;

  parms beta0  0;
  parms beta1  0;
  parms sigma2 100;

  prior beta0  ~ normal(0, var=10000);
  prior beta1  ~ normal(0, var=100);
  prior sigma2 ~ igamma(shape=0.001, scale=0.001);

  mu = beta0 + beta1 * dose;
  model response ~ normal(mu, var=sigma2);

run;

出力例

                   Posterior Summaries and Intervals

Parameter    N       Mean     Std Dev      5%       25%      50%      75%      95%
beta0       8000     3.421     3.184     -1.89     1.289    3.418    5.536    8.718
beta1       8000     0.597     0.054      0.51     0.561    0.597    0.633    0.686
sigma2      8000   143.275    27.614     103.4   123.831  140.158  159.211  192.451
📝 解釈
beta1 の事後中央値は 0.597（95%信用区間：0.51〜0.69）で、Rの結果と整合しています。
sigma2 は残差分散で、事後中央値は約140.2。標準偏差に換算すると √140.2 ≈ 11.8 となり、真の値（12）に近い推定が得られています。
📝 補足：PROC MCMCとPROC BGLIMMの使い分け
SAS 9.4 以降では PROC BGLIMM も利用可能で、一般化線形混合モデル（GLMM）のベイズ推定をよりシンプルな構文で実行できます。定型モデルには PROC BGLIMM を、カスタムモデルや非線形モデルには PROC MCMC を使い分けると効率的です。

製薬業界でのベイズ統計の活用場面

製薬業界では、以下の領域でベイズ統計が積極的に活用されています。

アダプティブデザイン

ベイズ流アダプティブデザインでは、試験途中のデータを逐次的に更新しながら試験設計を変更できます。代表的な適用例は次のとおりです。

  • サンプルサイズ再推定：中間解析で事後確率を評価し、必要に応じてサンプルサイズを調整します
  • 早期成功・早期中止の判断：「治療効果あり」の事後確率が95%を超えれば早期成功、「無効」の事後確率が高ければ早期中止（フュータリティ）と判断します
  • 用量選択の最適化：複数の用量群で並行試験を行い、中間解析の結果をもとに有望な用量に被験者配分を増やします（Response-Adaptive Randomization）

ベイズ統計を用いたサンプルサイズ設定では、頻度論とベイズの方法を比較しながら詳しく解説していますので、ぜひ参照ください。

希少疾患・単群臨床試験

希少疾患や小児適応では十分な症例数を確保することが困難です。ベイズ統計は事前情報を積極的に活用することで、少ない症例数でも信頼性の高い推論を可能にします。過去の試験データや文献情報を情報的事前分布（Informative Prior）として組み込む「Borrowing（情報借用）」の考え方は、希少疾患治験で特に有効です。二値エンドポイントの単群臨床試験におけるベイズ流デザイン入門でも具体的な手法を解説しています。

PK/PD（薬物動態・薬力学）解析

非線形混合効果モデル（NLME）によるPopPK（集団薬物動態）解析では、ベイズ推定が広く使われています。

  • 個体ベイズ推定（EBE: Empirical Bayes Estimates）：集団パラメータを事前分布として活用し、個体パラメータを推定します
  • 小サンプルPK解析：限られたサンプリングポイントしかない小児・高齢者データへの適用
  • NONMEMやMonolixなどのPopPKソフトウェアはベイズ推定を標準機能として搭載しています

バスケット試験・プラットフォーム試験

がん領域では複数の適応症を同時に評価するバスケット試験が増加しています。階層ベイズモデルを用いて適応症間で情報を共有しながら効率的に評価する手法（EXNEXモデル等）が普及しています。詳細はがん臨床試験におけるベイズ流バスケットデザインの理論と実装をご参照ください。

規制当局のガイダンスと対応

近年、規制当局もベイズ統計手法の活用を明確に認めるガイダンスを相次いで発出しています。

機関ガイダンス名主な内容
FDAAdaptive Designs for Clinical Trials of Drugs and Biologics（2019年）ベイズ流アダプティブデザインの適用基準と提出資料の考え方を整理
FDABayesian Statistics for Medical Device Clinical Trials（2010年）医療機器の臨床試験にベイズ統計を適用する際の考え方・事前分布の設定指針
ICHICH E9(R1) 統計的原則（2019年）推定目標（Estimand）の概念の中で、頻度論・ベイズを含む推定手法の選択が整理された

これらのガイダンスからわかるように、ベイズ統計を使う場合は事前分布の選択根拠・感度分析の実施・シミュレーションによる操作特性の確認を規制当局に示すことが求められます。

実務でのポイント

🔑 まとめ・実務ポイント

① 収束診断は必須
Rhat（<1.01）、有効サンプルサイズ（ESS）、トレースプロットを必ず確認してください。収束していない結果をそのまま使用するのは統計的根拠に欠けます。

② 事前分布の感度分析を実施する
特に規制当局への提出を想定した解析では、複数の事前分布（弱情報事前分布・情報的事前分布）でモデルを回し、結果の安定性を確認する感度分析が求められます。

③ シミュレーションで操作特性を確認する
アダプティブデザインに適用する場合、事前に大規模シミュレーション（10,000〜100,000回）で第I種過誤・検出力・サンプルサイズの操作特性を確認し、プロトコルに記載することが重要です。

④ バーンインサンプルは除外する
Rの warmup、SASの nbi で指定したバーンイン期間のサンプルは推論に使用しないよう注意してください。

📚 この記事をより深く理解するための参考書籍

ベイズ統計をさらに深く学びたい方に、おすすめの書籍をご紹介します。

『RとStanではじめるベイズ統計モデリングによるデータ分析入門』馬場真哉（講談社）
本記事で扱ったbrmsパッケージ・Stanの実装をより体系的に学びたい方に最適な一冊です。MCMCアルゴリズムの仕組みから実践的なモデリングまで、豊富なRコードとともに丁寧に解説されています。ベイズ解析を業務に取り入れたい実務家に強くおすすめします。

Amazon

楽天市場

関連記事

まとめ

本記事では、ベイズ統計の解析手法と製薬業界での活用について、RとSASの実装例を交えて解説しました。

ベイズ統計の解析はMCMCによる数値サンプリングを中心に進み、Rでは brms（Stanバックエンド）、SASでは PROC MCMC が主要なツールです。収束診断（Rhat・ESS・トレースプロット）と事前分布の感度分析は、いずれの実装においても欠かせないプロセスです。

製薬業界では、アダプティブデザイン・希少疾患試験・PK/PD解析・バスケット試験と多岐にわたる場面でベイズ統計が活用されており、FDAやICHのガイダンスでもその重要性が認められています。特に情報の乏しい状況や柔軟な意思決定が求められる試験設計において、ベイズ統計は頻度論にはない大きな強みを発揮します。

ぜひ本記事のコードを参考に、ご自身の解析にベイズ統計を取り入れていただければと思います。次回は、ベイズ流の信用区間と頻度論の信頼区間について、シミュレーションを交えて掘り下げる予定です。

MCMC R SAS ベイズ統計 生物統計 臨床試験 製薬業界
ABOUT ME
tomokichi
外資系製薬会社で生物統計家として働ている1児のパパ。生物統計家とは何か、どのようなスキルが必要か、何を行っているのかを共有していきたいと思っております！生物統計に関する最新情報を皆様にお届けすべく、日々奮闘中です。趣味は筋トレ、温泉巡り、家族と散歩。
BLOG：https://www.biostat4869.com