Distributional Regression

R研究集会2025
統計数理研究所共同利用研究集会 「データ解析環境Rの整備と利用」

Transcript

1. Confidential © Takehiko Yasuakwa All rights reserved. 0 Rで分布回帰をやってみる 〜Distributional Regression〜 R研究集会2025 統計数理研究所共同利用研究集会 「データ解析環境Rの整備と利用」 2 02 5年1 1月2 3日 安川武彦 株式会社JDSC * 本発表の内容は発表者が所属する組織の意見を代表するものではありません。
2. Confidential © Takehiko Yasuakwa All rights reserved. 1 目次 • はじめに • Distributional Regressionの分類 • GAMLSSの理論 • `gamlss package` & `gamlss2 package` • Rによる分析例 • まとめ
3. Confidential © Takehiko Yasuakwa All rights reserved. 2 平均値だけではなく、分布全体に対して回帰を行いたいときにどうすればよいかを調べてみた はじめに Stasinopoulos, D. M., Rigby, R. A., et al. (2017). Flexible Regression and Smoothing: Using GAMLSS in R. Chapman and Hall/CRC. FIGURE 2.13 Distributional Regression：共変量の関数として条件付き分布全体をモデル化することを目的としたア プローチの総称。 従来の回帰モデルが応答変数の条件付き平均のみをモデル化するのに対し、分布回帰は条件付き分布のすべ ての側面（平均、分散、歪度、尖度、分位点など）を共変量でモデル化する。 項目 通常の回帰 分布回帰 モデル化の対象 条件付き平均のみ 条件付き分布全体 定式化 （正規分布） 分散は定数 両方が共変量の関数 得られる情報 平均的な傾向のみ 平均、分散、歪度、尖度、分位点など 適用場面 平均への効果が主な関心 分散、裾、格差などへの関心
4. Confidential © Takehiko Yasuakwa All rights reserved. 3 Rage Against the Mean – A Review of Distributional Regression Approaches 分布回帰 • Kneib (2023) は、分布回帰のアプローチを以下のように分類し包括的にレビュー 分類 特徴 モデリング方法 Rパッケージ パラメトリック パラメトリックな分布のパラメータを 共変量でモデリング GAMLSS gamlss, bamlss セミパラメトリック 変換関数で条件付き分布を参照 分布にマッピング Conditional Transformation Models mlt, tram ノンパラメトリック 目的変数の条件付き密度を表現 密度回帰、ノンパラメトリックベ イズ logcondens 局所的 分位点やexpectileをモデル化 分位点回帰、期待値回帰 quantreg,expectreg, qgam Kneib, T., Silbersdorff, A., & Säfken, B. (2023). Rage Against the Mean – A Review of Distributional Regression Approaches. Econometrics and Statistics, 26, 99-123. Klein, N. (2024). Distributional Regression for Data Analysis. Annual Review of Statistics and Its Application, 11, 321-346.
5. Confidential © Takehiko Yasuakwa All rights reserved. 4 従来の平均中心の回帰分析では多くの重要な情報（分散、歪度、裾の形状など）が失わ れている 分布回帰 • 以下のような状況で分布回帰が必要であることが指摘されている（ Klein (2024) とKneib et al. (2023) ） Klein, N. (2024). Distributional Regression for Data Analysis. Annual Review of Statistics and Its Application, 11, 321-346. Kneib, T., Silbersdorff, A., & Säfken, B. (2023). Rage Against the Mean – A Review of Distributional Regression Approaches. Econometrics and Statistics, 26, 99-123. 平均だけでなく、極端な値への関心 1 分散・ばらつきへの関心 2 分布形状の複雑さ 3 予測区間の精度向上 4 • 児童栄養失調では極端な栄養不良のリスクが重要 • 市場：極端な変動が市場参加者の意識決定に直結 • 所得格差や教育効果のばあい、格差の拡大・縮小を捉える 必要がある • 多峰性、歪度 • 分散を評価して正確な予測区間がほしい 関心 例
6. Confidential © Takehiko Yasuakwa All rights reserved. 5 分布回帰に対するパラメトリックなアプローチがGAMLSS(Generalized Additive Models for Location, Scale and Shape)である GAMLSS GAMLSS：Rigby and Stasinopoulos (2005)で提案された、一般化線形モ デルに、応答変数の分布全体（Location, Scale, Shape)に対して回帰構造を 導入したモデル Rigby, R. A., & Stasinopoulos, D. M. (2005). Generalized additive models for location, scale and shape. Journal of the Royal Statistical Society: Series C (Applied Statistics), 54(3), 507-554.
7. Confidential © Takehiko Yasuakwa All rights reserved. 6 GAMLSSでは、分布のパラメータを共変量の関数としてモデル化することで、共変量に応じ た分布の変化を捉えた分析が可能となる GAMLSS • 各観測値は、パラメータベクトル （location, scale, skewness, kurtosis）を持つ確率分布に従うことを仮定
8. Confidential © Takehiko Yasuakwa All rights reserved. 7 GAMLSSは、平均だけでなくばらつきや歪みが説明変数によってどう変化するかを分析可能 な柔軟なモデリングフレームワーク GAMLSS 豊富な分布の選択肢 - 指数分部族に限らない多くの分布が利用可能（正規、対数正規、ガンマ、ベータ、t分布、BCPE、BCTなど） すべての分布パラメータをモデル化 - location, scale, shape のすべてを共変量の関数としてモデル化する 柔軟な予測子 - 線形効果だけでなく、非線形効果（スプライン）、ランダム効果、空間効果などを含めることが可能
9. Confidential © Takehiko Yasuakwa All rights reserved. 8 GAMLSSは、一般化線形モデルおよび一般化加法モデルの自然な拡張である GAMLSS 比較項目 GLM （一般化線形モデル） GAM （一般化加法モデル） GAMLSS （位置・スケール・形状モデル） モデル対象 位置（平均）のみ 位置（平均）のみ（非線形） 位置・スケール・形状（μ, σ, ν, τ） 説明変数の効果 線形 非線形（スプライン等） 非線形（各パラメータに個別指定可 能） 利用可能な分布 指数型分布（正規, ポアソン等） 指数型分布 指数型以外も可（t分布、Box-Cox型、 歪度付きなど） スケール・分散のモデリング 固定または関数形に限定 固定または簡単な構造のみ 自由に回帰構造を持たせることが可能 歪度・尖度のモデル化 不可 不可 可能（分布によるが多くが対応） 推定法 尤度最大化（MLE） MLE / 平滑化 MLE / ベイズ Rパッケージ glm(), stats::glm mgcv::gam, gam gamlss, bamlss, gamboostLSS 応用例 一般的な平均の回帰問題 非線形効果の分析（例：年齢効果） 成長曲線、家賃分布、リスク分布、歪度 付きデータなど
10. Confidential © Takehiko Yasuakwa All rights reserved. 9 利用可能な分布の例 GAMLSS • 以下の代表的な分布の他に、打ち切り、区間データ、混合分布なども用意されている。 https://www.gamlss.com/information/distributions/ 分布 名称 パラメータ数 用途 連 続 分 布 正規分布 NO 2 対称な連続データ 対数正規分布 LOGNO 2 正の歪んだデータ ガンマ分布 GA 2 正の連続データ ベータ分布 BE 2 (0,1) 区間のデータ t分布 TF 3 裾の重いデータ Box-Cox Power Exponential BCPE 4 柔軟な形状 Box-Cox t BCT 4 柔軟な形状と裾 離 散 分 布 ポアソン分布 PO 1 カウントデータ 負の二項分布 NBI 2 過分散のカウントデータ ゼロ過剰ポアソン ZIP 2 ゼロが多いカウントデータ ベータ二項分布 BB 2 二項データの過分散
11. Confidential © Takehiko Yasuakwa All rights reserved. 10 確率分布：BCCG (Box-Cox Cole and Green) 分布 GAMLSS • Box-Cox変換を施したデータが、正規分布に従うと仮定。3つのパラメータ mu (中央値の近似), sigma (変動係数の近 似), nu (歪度) を持つ。主に正の値をとる歪んだデータのモデリングに使用される
12. Confidential © Takehiko Yasuakwa All rights reserved. 11 確率分布：BCT (Box-Cox t)分布 GAMLSS • Box-Cox変換を施したデータが、t分布に従うと仮定。BCCGに加えて4つ目のパラメータ tau (尖度、裾の重さ) を持つ。これ により、歪みだけでなく、外れ値の多さ（裾の厚み）も調整できる、さらに柔軟。
13. Confidential © Takehiko Yasuakwa All rights reserved. 12 確率分布：ZAP (Zero-Altered Poisson) 分布 GAMLSS • 構造的なゼロ（Structural Zeros）と偶然のゼロ（Sampling Zeros）を区別してモデリング。ゼロかどうかの二値過程 （ベルヌーイ分布）と、ゼロではない場合のカウント過程（切断ポアソン分布）を組み合わせたモデル。
14. Confidential © Takehiko Yasuakwa All rights reserved. 13 パラメータの推定は対数尤度関数を最大化することで求め、これを計算するためにRigby- Stasinopoulos (RS) Algorithmが採用されている GAMLSS RSアルゴリズム：Backfittingの考え方にもとづき、すべて のパラメータを同時に最適化するのではなく、「一つのパラ メータベクトルに注目し、他は固定する」というサイクルを交 互に繰り返すアルゴリズム • 各パラメータを順番に更新 • 各ステップはFisher’s Scoringアルゴリズムに帰着 • 非線形項の推定では、ペナルティ項（平滑化パラメータ）を 導入しオーバーフィッティングを防いでいる • 指数型分布族以外では大域解があることが保証されていな いことに注意 尤度：
15. Confidential © Takehiko Yasuakwa All rights reserved. 14 推定は最尤法に根ざしており、漸近的な性質は保持。一方、指数型分布族以外に拡張した ためGLMで保証されている性質の一部は失われる GAMLSS 項目 GLM GAMLSS モデル化の対象 平均 (μ) のみ 分布全体 (μ, σ, ν, τ) 応答変数の分布 指数型分布族に限定 ほぼ制限なし 推定方法の名称 最尤法 (IWLSアルゴGAMLSSの推定と統計的性質についてまと めていきたいゴリズム) 罰則付き最尤法 (RS/CGアルゴリズム) 尤度関数の形状 正準リンクの場合、上に凸であることが保証される。 形状は保証されない。複数の極値を持つ可能性がある。 解の一意性 解は一意に定まることが多い。 局所解に陥る可能性があり、解の一意性は保証されない。 推定の安定性 非常に高い。 モデルによる。柔軟性が高い分、収束の確認が重要。 統計的性質の保証 指数型分布族の性質に強く依存する。 最尤法の一般的な漸近理論に依存する。 保証の前提条件 指数型分布族であること。 正則条件（識別可能性など）を満たすこと。 必要なサンプルサイズ 比較的小さくても安定しやすい。 複雑なモデル（特に形状パラメータを含む場合）ではより多くを要 求する。
16. Confidential © Takehiko Yasuakwa All rights reserved. 15 モデル診断の方法としてワームプロットと正規化分位点残差が推奨されている GAMLSS ワームプロット*（分布の適合度） 正規化分位点残差 QQプロットを拡張したもので残差のパターンを評価するもの。全て の点が水平な２本の点線で示される信頼区間の内側にランダムに 収まっていることが理想 • 水平：よくあてはまっている。傾いているときは分散のモデル化を疑う • U字・逆U字：分布の裾が軽い（平均のモデル化が不適切） • S字：歪度のモデル化を疑う • M字：尖度のモデル化を疑う モデルが正しければ以下で定義された正規化残差が標準正規分 布に従う * van Buuren and Fredriks M. (2001) Worm plot: simple diagnostic device for modelling growth reference curves. Statistics in Medicine, 20, 1259–1277. 例：胎児超音波検査 Stasinopoulos 他 (2024)のChapter 9の例から
17. Confidential © Takehiko Yasuakwa All rights reserved. 16 Rのgamlssパッケージは分布回帰を実装した包括的なフレームワークである gamlss package • 各パラメータに関するモデルを定義して推定する https://www.gamlss.com https://cran.r-project.org/web/packages/gamlss/index.html パラメータ 説明 デフォルト値 formula 平均(mu)のモデル式 必須 sigma.formula 分散(sigma)のモデル式 ~ 1（定数） nu.formula 歪度(nu)のモデル式（3 パラメータ分布以上） ~ 1 tau.formula 尖度(tau)のモデル式（4 パラメータ分布） ~ 1 data データフレーム 必須 family 分布族 NO（正規分布） method 推定アルゴリズム RS()（Rigby- Stasinopoulos） trace 収束過程の表示 TRUE control 制御パラメータ gamlss.control()
18. Confidential © Takehiko Yasuakwa All rights reserved. 17 加法（Additive）モデルを構成する上で、平滑化関数は極めて重要な役割を担う. GAMLSSでは多くの選択肢がある gamlss package 関数 名称 概要 主な用途 特徴・使い分け pb() P-spline (Penalized B- Spline) B-spline基底関数と差分罰則を組み合わせ た、柔軟で頑健な平滑化手法。 最も一般的。滑らかな非線形関係のモデリング。 デフォルトの選択肢。dfやlambdaで平滑度を 直感的に制御可能。過剰適合を効果的に抑 制 cs() Cubic Spline データを複数の区間に分割し、各区間を3次 多項式で近似し、それらを滑らかに接続する古 典的な手法。 pb()登場前の標準的な平滑化手法。比較 的単純な非線形関係。 pb()に比べて柔軟性が劣る場合がある。ノット の数と位置の指定が重要になるため、現在は pb()の使用が推奨 re() Random Effects カテゴリカル変数（個体、地域、グループなど） が持つばらつきを、ランダム効果としてモデル化 する。 パネルデータ、階層データ、反復測定データにお ける個体差やグループ差の吸収。 各カテゴリの水準ごとに独立したパラメータを推 定するのではなく、それらが共通の正規分布か ら生成されると仮定することで、パラメータ数を 節約し、頑健に推定 ga() GAM Additive Term mgcvパッケージのgam()関数で使われる平 滑化項（例: s(), te()）をGAMLSS内で利 用するためのインターフェース。 空間データのモデリング、複数変数の交互作用 の平滑化など、より高度な平滑化が必要な場 合。 mgcvパッケージの強力な平滑化機能（薄板 回帰スプライン、テンソル積スプラインなど）を活 用できる。GAMLSSの枠組みで、より複雑な 交互作用をノンパラメトリックに扱いたい場合 fp() Fractional Polynomials べき乗変換 (x^p) を組み合わせることで、多 様な曲線形状を表現するパラメトリックなアプ ローチ。 単調な増加・減少関係や、S字カーブなど、特 定の関数形で表現できる非線形関係。 ノンパラメトリックなスプラインとは異なり、陽な関 数形で関係を記述できるため、解釈が容易な 場合がある。しかし、表現できる形状は限定的 https://www.gamlss.com/additive-terms/
19. Confidential © Takehiko Yasuakwa All rights reserved. 18 モデルの要約、診断、予測は一般のRの作法で実行できる gamlss package summary(model object): - モデル全体の適合度指標（Global Deviance, AIC, BIC）。各パラメータ（Mu, Sigmaなど）のモデル式と、係数の推定値、 標準誤差、z値、p値を表示 plot(model object): - デフォルトで4つの診断プロットを生成。これにより、残差のパターンを視覚的に確認する wp(model object): - ワームプロット（Worm Plot）を生成。これはgamlssにおける最も重要な診断ツール。 term.plot(model object, what = "mu", pages = 1): - 指定したパラメータ（what引数）のモデルに含まれる各項（特に平滑化項）の効果をプロットします。 predict(model object, newdata = ...): - デフォルトでは、各パラメータの線形予測子の値を返す。type = "response" を指定すると、リンク関数を逆変換した後のパラメー タの値（例: μ や σ の予測値）を返す。what引数で、どのパラメータについて予測するかを指定（what = "mu", what = "sigma"など）。 centiles.pred(model object, xname = "fl", xvalues = ...): - 指定した説明変数 xname の値の変化に伴い、応答変数の予測パーセンタイル（デフォルトでは5, 25, 50, 75, 95パーセンタイ ル）がどのように変化するかを計算し、プロットする。
20. Confidential © Takehiko Yasuakwa All rights reserved. 19 機能を拡張するパッケージエコシステムが多く開発されている gamlss package https://www.gamlss.com/the-r-packages/ パッケージ名 主な機能と役割 gamlss GAMLSSモデルをフィッティングするための中心的なパッケージ。gamlss.distとgamlss.dataに依存している。 gamlss.dist gamlss.familyで利用可能な、100を超える連続・離散・混合分布の確率分布関数を提供する。 gamlss.data GAMLSSの書籍やドキュメントで使用されるサンプルデータセットを収録している。 gamlss.add モデルに追加できる、特殊な加法項（additive terms）を提供する。 gamlss.cens 打ち切り（Censored）データ（左側、右側、区間打ち切り）を扱うモデルのフィッティングを可能にする。 gamlss.mx 有限混合分布（Finite Mixture Distributions）をフィッティングするための機能を提供する。 gamlss.tr 切断分布（Truncated Distributions）を生成し、フィッティングするための機能を提供する。 gamlss.inf ゼロ過剰（Zero-inflated）およびゼロ調整（Zero-adjusted）連続分布のフィッティングをサポートする。 gamlss.spatial ガウス型マルコフ確率場（Gaussian Markov Random Fields, MRF）モデルをGAMLSS内でフィッティングする。 gamlss.nl 非線形モデル（Non-linear models）のフィッティングをサポートする。 gamlss.ggplots ggplot2を用いて、モデルの診断や解釈に役立つグラフィックスを作成する。 gamlss.foreach foreachパッケージを利用して、並列計算を可能にする。 gamlss.demo rpanelパッケージを使用して、分布や平滑化のデモンストレーションを行う。 gamlss.util GAMLSSに直接関連しないユーティリティ関数を含む。 gamboostLSS ブースティング（Boosting）法を用いてGAMLSSモデルをフィッティングする。 bamlss ベイズ的アプローチ（MCMC法など）を用いてGAMLSSモデルをフィッティングする。
21. Confidential © Takehiko Yasuakwa All rights reserved. 20 `gamlss`をオーバーホールした`gamlss2`が開発されている gamlss package • 推定速度と適応性の向上を目指して開発。ただし、 gamlss2はまだ開発中 - 速度と効率（大規模なデータへの対応） - `mgcv`との統合（平滑化項の拡張） - 構文の簡潔化 - 拡張性（他のモデリング手法との連携） https://gamlss-dev.github.io/gamlss2/
22. Confidential © Takehiko Yasuakwa All rights reserved. 21 ある地域の5つの小学校で、学年が上がるにつれて学力がどのように変化するかを分析する 分析例（1） • 分析シナリオ・データ：指導法の効果の測定 - 学校数: 5校 - 各学校の生徒数: 20名 - 学年: 1年生、2年生、3年生 - 応答変数: 学力スコア（0-100点） • 問い 1. 学年が上がると、平均学力は上昇するか？ 2. 学年が上がると、学力のばらつきは変化するか？ 3. 学校間の違いを考慮するには？ 変数名 説明 型 school_id 学校ID（1-5） カテゴリカル student_id 生徒ID（1-100） カテゴリカル grade 学年（1, 2, 3） 数値 score 学力スコア 連続値 • 平均的な効果は肯定的: 指導法は学力向上に 寄与する • 学校間の差: 学校の特性により効果が異なる可 能性がある • 学力格差の拡大: 下位層への追加支援が必要
23. Confidential © Takehiko Yasuakwa All rights reserved. 22 最も単純なモデルから開始 分析例（1） • 学校間の差を考慮していない • すべての生徒が同じ母集団から来ていると仮定 • 階層性のあるデータ構造を反映できていない
24. Confidential © Takehiko Yasuakwa All rights reserved. 23 学年が上がると分析も大きくなっている状況（分散不均一）をモデル化 分析例（1） • gradeの固定効果に加えて、school_idのランダム切片とランダム 係数を含む（学校間の違いを考慮） • 平均のモデルに加えて、分散のモデルを構築。分散のモデルは、 sigma.formula引数にて指定（分散不均一を考慮） * re()項の情報を取り出すときは：re_obj <- getSmo(final_model, "re(random = ~ grade | school_id)"); summary(re_obj) 学年が1つ上がると、平均学 力が約9.49点上昇 学年が1つ上がると、 log(σ)が約0.27増加
25. Confidential © Takehiko Yasuakwa All rights reserved. 24 モデル診断では、’gamlss_ggplots’ パッケージを使うのが便利 分析例（1） まだ裾の当てはまりが悪い
26. Confidential © Takehiko Yasuakwa All rights reserved. 25 モデル診断では、’gamlss_ggplots’ パッケージを使うのが便利 分析例（1） • 予測した平均と分散を利用して、予測分布のプロットで、学年が上がると分散も大きくなっているモデルの様子がわかる 学校ID=1の場合
27. Confidential © Takehiko Yasuakwa All rights reserved. 26 ニューヨークの大気データ：オゾン濃度は気温にどのような影響を与えるのか 分析例（2） • データ： - `airquality` • 問い 1. 気温は、日々のオゾン濃度の平均値だけでな く、そのばらつき（リスク）や分布の形状（歪 み）にどのように影響を与えるか？
28. Confidential © Takehiko Yasuakwa All rights reserved. 27 高温日に極端な高濃度オゾンが発生するリスクが高まるリスクがある 分析例（2） • 気温が上昇するにつれて、オゾン濃度の平均値は非線形に増加 • 気温が上昇するほど、オゾン濃度のばらつき（標準偏差）も大きくなる
29. Confidential © Takehiko Yasuakwa All rights reserved. 28 モデル診断からは裾のあてはまり、尖度・歪度に改善の余地がある 分析例（2） 裾のあてはまりに問題がある
30. Confidential © Takehiko Yasuakwa All rights reserved. 29 平均値のみの分析で「多くの情報を捨てているのではないか？」と感じたら、、、分布回帰をす ればいい まとめ 分布回帰：条件付き分布をモデリングするという発想 指数型分布族以外の分布の利用 平滑化関数の組み合わせにより現象をモデリング 柔軟であるが故に、分析者の力量が試される
31. Confidential © Takehiko Yasuakwa All rights reserved. 30 参考文献 • Stasinopoulos, M. D., Rigby, R. A., Heller, G. Z., Voudouris, V., De Bastiani, F . (2017). Flexible Regression and Smoothing: Using GAMLSS in R. Chapman & Hall/CRC - GAMLSSの入門書。 • Rigby, R. A., Stasinopoulos, D. M., Heller, G. Z., De Bastiani, F . (2019). Distributions for Modeling Location, Scale, and Shape: Using GAMLSS in R. Chapman & Hall/CRC - GAMLSSで利用可能な分布について網羅的に解説 • Stasinopoulos, M. D., Rigby, R. A., Heller, G. Z., Voudouris, V., De Bastiani, F . (2024). Generalized Additive Models for Location, Scale and Shape: A Distributional Regression Approach, with Applications. Cambridge University Press - GAMLSSの分布回帰としての位置づけをまとめた総説
32. Confidential © Takehiko Yasuakwa All rights reserved. 31 参考文献 • Rigby, R. A., & Stasinopoulos, D. M. (2005). Generalized additive models for location, scale and shape, (with discussion). Journal of the Royal Statistical Society: Series C (Applied Statistics), 54(3), 507-554. - GAMLSSを提案した論文 • Borghi, E., de Onis, M., Garza, C., Van den Broeck, J., Frongillo, E.A., Grummer-Strawn, L., Van Buuren, S., Pan, H., Molinari, L., Martorell, R., Onyango, A.W., Martines, J.C. and for the WHO Multicentre Growth Reference Study Group (2006), Construction of the World Health Organization child growth standards: selection of methods for attained growth curves. Statist. Med., 25: 247- 265. - WHO（世界保健機関）の子供の成長基準曲線の作成を解説した論文。GAMLSSが用いられている。 • Klein, N. (2024). Distributional Regression for Data Analysis. Annual Review of Statistics and Its Application, 11:321-46. • Kneib, T., Silbersdorff, A., & Säfken, B. (2023). Rage Against the Mean – A Review of Distributional Regression Approaches. Econometrics and Statistics, 26, 99-123.
33. Confidential © Takehiko Yasuakwa All rights reserved. 32 本ファイルの内容の一部、または全部を無断で転用・転載することを禁じます。 EOF