🎨 直感で掴む — 「真ん中」 vs 「上限・下限」

例え話：通勤時間の予測

あなたは「家から駅まで何分かかるか」を予測したいとします。 普通の回帰なら「平均 12 分」と返してきます。 でも、 本当に知りたいのは：

• 「電車に遅れずに着くには何分前に出ればよいか？」 → 90% 分位点（最悪寄りの値）。 例：17 分
• 「だいたい何分かかる？」 → 中央値（50% 分位点）。 例：11 分
• 「最短はどれくらい？」 → 10% 分位点。 例：8 分

これらはすべて「同じデータから抽出できる別々の予測値」です。 分位点回帰は 3 本まとめて出してくれる 強力なツール。

不等分散（heteroskedastic）データに特に強い

「所得が高いほど消費のばらつきも大きい」「人口の多い県ほど県内総生産のばらつきも大きい」 のような 不等分散 データでは、 通常の OLS が引いた 1 本の直線は 「真ん中だけ」 しか見せてくれません。

OLS との対比

📐 数式 — Koenker & Bassett (1978) の定式化

分位点の定義

pinball loss（チェック関数）

分位点回帰の最適化問題

OLS との対比：損失関数の違い

🔬 数式を「言葉」で読み解く

🧮 SSDSE-B 都道府県データで計算してみる

SSDSE-B 2026 の 2023 年データを用い、 「県内総生産 $y$ ~ 人口総数 $x$」 の分位点回帰を τ = 0.10, 0.50, 0.90 の 3 本で推定します。

STEP 0: 何が知りたいか

普通の回帰なら 1 本の直線で「人口が 1 万人増えると県内総生産は何億円増えるか（平均）」を見ます。 しかし都道府県では 大都市と地方で生産性が大きく違う ので、 「平均だけ」では情報を失います。

STEP 1: 数値の準備（標本値の代表例）

説明変数 $x$ = 人口総数（万人）、 被説明変数 $y$ = 県内総生産（10 億円）。 5 県を抜粋：

STEP 2: 3 本の分位点回帰直線（47 県全部で推定した典型結果）

STEP 3: 解釈

• 傾きが τ で増加（4.2 → 4.9 → 6.3） → 「高生産地域ほど、 人口 1 万人あたりの生産が大きい」。 規模の経済が上位県で強い ことを示唆。
• 3 本の 幅（90% − 10%）が右ほど広がる → 不等分散の証拠。 大都市ほど予測の不確実性が大きい。
• τ=0.50 と τ=0.10/0.90 を比較すると、 OLS（平均直線）はだいたい τ=0.55 付近 に位置する（東京の極端値に引っ張られて平均が中央値より高い）。

🐍 Python 実装 — statsmodels と LightGBM の 2 通り

1. statsmodels (古典的・線形分位点回帰)

import pandas as pd
import numpy as np
import statsmodels.formula.api as smf
import matplotlib.pyplot as plt

# SSDSE-B 2026 を読み込み、2023 年に絞る
df = pd.read_csv('data/raw/SSDSE-B-2026.csv', encoding='utf-8-sig', header=1)
df = df[df['年度'] == 2023].copy()
df = df.rename(columns={'A1101': '人口', 'A1102': '県内総生産'})

# 単位を万人・10億円に揃える（数値の桁を整える）
df['人口万人'] = df['人口'] / 1e4
df['生産10億'] = df['県内総生産'] / 1e3

# 3 つの分位点で回帰
taus = [0.10, 0.50, 0.90]
models = {}
for tau in taus:
    mod = smf.quantreg('生産10億 ~ 人口万人', data=df)
    res = mod.fit(q=tau)
    models[tau] = res
    print(f'τ={tau}: intercept={res.params[0]:.2f}, slope={res.params[1]:.3f}')

2. 結果の可視化（扇形に広がる予測区間）

fig, ax = plt.subplots(figsize=(9, 6))

# 散布図
ax.scatter(df['人口万人'], df['生産10億'],
           s=70, alpha=0.6, edgecolor='white', linewidth=1.2,
           color='#1565C0', label='47 都道府県 (2023)')

# 3 本の分位点直線
x_grid = np.linspace(df['人口万人'].min(), df['人口万人'].max(), 100)
colors = {0.10: '#D32F2F', 0.50: '#388E3C', 0.90: '#D32F2F'}
styles = {0.10: '--', 0.50: '-', 0.90: '--'}
for tau, res in models.items():
    a, b = res.params
    ax.plot(x_grid, a + b * x_grid,
            color=colors[tau], linestyle=styles[tau], linewidth=2.2,
            label=f'τ = {tau:.2f}')

# 比較のため OLS も
import statsmodels.formula.api as smf
ols = smf.ols('生産10億 ~ 人口万人', data=df).fit()
ax.plot(x_grid, ols.params[0] + ols.params[1] * x_grid,
        color='#1976D2', linewidth=2.2, linestyle=':', label='OLS (平均)')

ax.set_xlabel('人口総数（万人）')
ax.set_ylabel('県内総生産（10 億円）')
ax.set_title('分位点回帰：県内総生産 ~ 人口（10/50/90 パーセンタイル）')
ax.legend(loc='upper left')
plt.tight_layout()
plt.show()

3. LightGBM での非線形分位点回帰

大きなデータや非線形関係には、 LightGBM の objective='quantile' が高速で便利です。

from lightgbm import LGBMRegressor

# 同じデータで 3 本の分位点モデルを学習
X = df[['人口万人']].values
y = df['生産10億'].values

quantile_models = {}
for tau in [0.10, 0.50, 0.90]:
    m = LGBMRegressor(objective='quantile', alpha=tau,
                       n_estimators=300, learning_rate=0.05,
                       num_leaves=15, random_state=42, verbose=-1)
    m.fit(X, y)
    quantile_models[tau] = m

# 各 τ で予測（同じ x_grid 上）
X_grid = x_grid.reshape(-1, 1)
preds = {tau: m.predict(X_grid) for tau, m in quantile_models.items()}

# 予測区間として fill_between
fig, ax = plt.subplots(figsize=(9, 6))
ax.fill_between(x_grid, preds[0.10], preds[0.90],
                alpha=0.18, color='#D32F2F', label='80% 予測区間')
ax.plot(x_grid, preds[0.50], color='#388E3C', linewidth=2.3, label='中央値')
ax.scatter(df['人口万人'], df['生産10億'], s=60, alpha=0.6, color='#1565C0')
ax.set_xlabel('人口（万人）'); ax.set_ylabel('県内総生産（10億円）')
ax.legend(); plt.tight_layout(); plt.show()

4. pinball loss を直接計算して評価

from sklearn.metrics import mean_pinball_loss

# τ=0.9 モデルの平均 pinball loss
y_pred_90 = quantile_models[0.90].predict(X)
loss_90 = mean_pinball_loss(y, y_pred_90, alpha=0.90)
print(f'τ=0.90 の pinball loss: {loss_90:.2f}')

# カバレッジ（予測区間に入った観測値の割合）
y_lo = quantile_models[0.10].predict(X)
y_hi = quantile_models[0.90].predict(X)
coverage = ((y >= y_lo) & (y <= y_hi)).mean()
print(f'80%予測区間のカバレッジ: {coverage*100:.1f}%')  # 理論値 80%

5. 多変量への拡張（複数説明変数）

# 人口 + 高齢化率 + 出生数 で 3 つの分位点を推定
df = df.rename(columns={'A6101': '高齢化率', 'A4101': '出生数'})

for tau in [0.10, 0.50, 0.90]:
    res = smf.quantreg('生産10億 ~ 人口万人 + 高齢化率 + 出生数',
                       data=df).fit(q=tau)
    print(f'τ={tau}:')
    print(res.params.round(4))
    print()

⚠️ 5 つの「分位点回帰 落とし穴」

🌐 関連手法・派生

🔗 関連用語（前提・並列・発展）

🟦 前提となる用語

🟨 並列・比較対象

🟩 発展・応用

📚 関連グループ教材

• 回帰分析の総合教材 — OLS、 ロジスティック、 分位点の比較
• 線形回帰入門 — 平均を予測する古典的回帰
• 損失関数 — pinball loss を含む損失関数群の整理
• 外れ値処理 — 分位点回帰を含む頑健統計の章

🧪 SSDSE-B 実データで「OLS と分位点回帰の差」を実感する

「分位点回帰は外れ値に強い」「不等分散データで威力を発揮する」と説明されても、 実感が湧きにくいですよね。 ここでは SSDSE-B の 47 都道府県データで 5 つの比較実験 を行い、 OLS との具体的な差を数値で確認します。

実験 1: 東京を含めるか除くかで OLS の係数はどう変わる？

SSDSE-B 2023 年の「県内総生産 ~ 人口」回帰で：

解釈：OLS の傾きは大都市圏の「梃の効いた点」に大きく引っ張られている。 中央値回帰（τ=0.5）はその影響を受けにくく、 「典型的な都道府県での関係」 をより正しく表す。

実験 2: 説明変数を増やしたとき τ ごとに係数はどう動く？

$y$ = 県内総生産、 $x_1$ = 人口、 $x_2$ = 高齢化率、 $x_3$ = 大学進学率 で重回帰：

解釈：

• 人口係数は τ で単調増加（規模の経済が上位ほど強い）
• 高齢化率係数は τ で減少（負の効果が上位では弱い）
• 大学進学率は τ=0.9 で τ=0.1 の 3 倍 の効果 → 高生産県では教育の生産性貢献がさらに大きい

OLS 1 本だけでは、 これらの「分布のなかでの違い」は 絶対に 見えません。

実験 3: 残差プロットの読み方

OLS で残差プロット（x 軸：予測値、 y 軸：残差）を描いたとき、 残差の広がりが 「漏斗状」 なら不等分散の証拠。 これは分位点回帰の出番です。

# OLS 残差プロット
import matplotlib.pyplot as plt
ols_resid = ols.resid
fig, ax = plt.subplots(figsize=(8, 5))
ax.scatter(ols.fittedvalues, ols_resid, alpha=0.6, s=60)
ax.axhline(0, color='red', linestyle='--')
ax.set_xlabel('予測値'); ax.set_ylabel('残差')
ax.set_title('OLS 残差プロット — 漏斗状なら不等分散')
plt.show()

実験 4: Breusch-Pagan 検定で不等分散を統計的に確認

from statsmodels.stats.diagnostic import het_breuschpagan

bp = het_breuschpagan(ols.resid, ols.model.exog)
print(f'BP statistic: {bp[0]:.2f}, p-value: {bp[1]:.4f}')
# p < 0.05 なら不等分散の証拠 → 分位点回帰の出番

実験 5: ブートストラップで分位点係数の信頼区間

import numpy as np

n_boot = 1000
boot_slopes = []
for _ in range(n_boot):
    idx = np.random.choice(len(df), len(df), replace=True)
    sample = df.iloc[idx]
    res = smf.quantreg('生産10億 ~ 人口万人', data=sample).fit(q=0.9)
    boot_slopes.append(res.params[1])

ci_lo, ci_hi = np.percentile(boot_slopes, [2.5, 97.5])
print(f'τ=0.9 傾きの 95% 信頼区間: [{ci_lo:.2f}, {ci_hi:.2f}]')

📊 応用シナリオ集 — 分位点回帰が真価を発揮する場面

シナリオ 1: 公共政策の格差分析

「教育予算が増えたら学力は上がるか？」を考えるとき、 OLS の平均効果だけ見ても 下位層の子どもに効いたのか、 上位層がさらに伸びたのか が分からない。 分位点回帰で τ=0.1（学力下位 10%）と τ=0.9（学力上位 10%）の係数を比較すると、 教育政策の 平等化効果か格差拡大効果か を識別できる。

SSDSE-B では 47 都道府県の「教育費」と「学力テスト得点」の関係を τ で分解し、 政策の 「底上げ効果」 vs 「トップ層伸長効果」を可視化できる。

シナリオ 2: 金融リスク管理 (VaR の推定)

Value-at-Risk (VaR) は「99% の確率でこれ以上は損しない」金額。 これは τ=0.01 の分位点 そのもの。 リターン分布の左裾を 説明変数（マクロ経済指標）の関数として 推定できる。

• $y$ = ポートフォリオの日次リターン
• $x$ = VIX、 金利スプレッド、 為替変動など
• $\tau = 0.01, 0.05$ の分位点回帰で「悪条件下での最悪リターン」を予測

シナリオ 3: 気象・気候予測の不確実性

「明日の最高気温は 23 ℃」と単一値で予測するより、 「80% の確率で 19〜26 ℃ に入る」と区間で予測する方が意思決定に有用。 分位点回帰で τ=0.1, 0.9 を同時推定すれば、 即座に 80% 予測区間が手に入る。

シナリオ 4: 医療における体重曲線

WHO の 成長曲線（growth chart） は本質的に分位点回帰の応用。 「年齢 $x$ の子どもの体重」の τ=3%, 50%, 97% を引いて、 「あなたの子は同年齢の何 % 位にいる」と判断する。

シナリオ 5: 不平等指標との接続

所得分布の τ=0.1 と τ=0.9 の差は 「10/90 比」（不平等指標）。 これを 地域・年齢・教育水準 の関数として推定すれば、 不平等が「誰のあいだで」「いつ」拡大したかを分解できる。 ピケティ的な分析。

📋 比較チートシート — どの回帰を選ぶ？

分位点回帰を「使うべきサイン」 — 5 つの兆候

1. 残差プロットで 「漏斗状」（広がりが x で変化）の不等分散が見える
2. OLS の残差が 非対称な分布（歪んでいる）
3. 「平均」よりも「下位／上位の動き」 が政策・実務的に重要
4. 「予測値だけでなく 予測区間」を出す必要がある
5. 少数の外れ値が OLS 結果を大きく動かしている

🎓 理論的補足 — 推定量の性質と歴史

歴史的経緯

• 1755 年：Roger Boscovich が天文観測データに 絶対誤差最小化（中央値回帰の祖）を提案。 同年代に Laplace も独立に研究。
• 1789 年：Laplace が「絶対誤差和最小化」の解法を整備。 これが現代の τ=0.5 分位点回帰の原型。
• 1805 年：Legendre が 最小二乗法（OLS）を発表。 計算容易性から L1 系を圧倒し、 19〜20 世紀の主流に。
• 1978 年：Roger Koenker と Gilbert Bassett が論文 “Regression Quantiles”（Econometrica）を発表。 任意 τ の分位点を線形計画問題として推定する現代的フレームワークを確立。
• 1990s：計算機の発達で実用化。 経済学・公衆衛生・気象で急速に普及。
• 2005〜2015：分位点ランダムフォレスト（Meinshausen 2006）、 ニューラル分位点回帰、 LightGBM/XGBoost の objective 統合などで 機械学習版 が確立。

分位点回帰推定量の漸近的性質

適切な正則条件下で、 分位点回帰推定量 $\hat{\boldsymbol\beta}(\tau)$ は次の性質を持つ：

• 一致性：$\hat{\boldsymbol\beta}(\tau) \overset{p}{\to} \boldsymbol\beta(\tau)$
• 漸近正規性：$\sqrt{n}\bigl(\hat{\boldsymbol\beta}(\tau) - \boldsymbol\beta(\tau)\bigr) \overset{d}{\to} \mathcal{N}(0, \Sigma(\tau))$
• 分散共分散行列 $\Sigma(\tau)$ は 誤差の条件付き密度 $f_{u|x}(0)$ に依存。 これが推定誤差の計算を難しくする。

OLS との効率比較

誤差が正規分布のとき、 中央値回帰の漸近相対効率（OLS 比）は：

線形計画問題としての定式化

pinball loss の最小化は、 補助変数 $u_i^+, u_i^- \ge 0$ を導入することで 線形計画問題（LP） として書ける：

これによりシンプレックス法・内点法など 確立された LP ソルバー がそのまま使え、 大規模問題でも実用的に解ける。

🔢 statsmodels.QuantReg / LightGBM 主要パラメータ表

statsmodels.regression.quantile_regression.QuantReg.fit() の主要引数

LightGBM の objective='quantile' 関連パラメータ

📚 ミニケーススタディ — SSDSE-B 死亡率の分位点分析

もう 1 つ実データで使い方を確認しましょう。 SSDSE-B 2026 の 「死亡率 ~ 高齢化率」 の関係を、 通常の回帰と分位点回帰で比較します。

背景

「高齢化率が上がると死亡率も上がる」は当然の関係ですが、 「どの程度の傾きで」「下位／上位の県では別の傾向か」を分位点回帰で詳細に見ます。

OLS の結果（参考）

$y$ = 死亡率（人口千あたり）、 $x$ = 高齢化率（%）。 47 県で OLS：

分位点回帰の結果

解釈

• 傾きが τ で 緩やかに増加（0.532 → 0.561、 約 5% 増）。 高齢化に対する死亡率の感応度は 高死亡率群でやや強い が、 全体としては均一（等分散に近い）。
• OLS と τ=0.5 の傾きはほぼ一致（0.547 vs 0.546）。 これは 外れ値・歪みが小さく、 線形関係が極めて頑健 なことを示す。
• これは 分位点回帰を使ってもあまり情報が増えない ケース。 逆にこういう状況なら OLS を安心して使ってよいと判定できる。

逆に分位点回帰が威力を発揮するケース

SSDSE-B の 「県内総生産 ~ 人口」 では、 τ=0.10 と τ=0.90 で傾きが 4.18 vs 6.31 と 50% も違う。 これは 明確な不等分散 の証拠で、 OLS の 1 本の傾きでは情報を 1/2 失う。 こういうデータには分位点回帰が必須です。

使い分けのまとめ

❓ よくある質問 (FAQ)

Q1. τ=0.5 と OLS、 結局どっちを使うべき？

A1. データに 明らかな外れ値や歪み があるなら τ=0.5（中央値回帰）の方が頑健。 誤差が正規・等分散で外れ値もない理想的状況なら OLS の方が 1.57 倍効率がよい（漸近相対効率）。 実務では 両方推定して比較 するのが正攻法。 大きく違うなら外れ値や不等分散を疑う合図。

Q2. 何個 τ を推定すべき？

A2. 探索目的なら τ ∈ {0.05, 0.10, 0.25, 0.50, 0.75, 0.90, 0.95} の 7 点が定番。 予測区間が目的なら 90% 区間 → τ=0.05, 0.95、 80% 区間 → τ=0.1, 0.9。 多すぎると分位点交差問題が増える。

Q3. p 値や信頼区間はどう計算？

A3. statsmodels の summary() はカーネル密度を使った漸近近似で SE を返すが、 小標本では信頼性が低い。 ブートストラップ法（1000 リサンプル）で信頼区間を求めるのが推奨。 上記コード例参照。

Q4. 多変量・カテゴリ変数も扱える？

A4. 完全に扱える。 statsmodels の formula API なら y ~ x1 + x2 + C(category) のように OLS と同じ書き方でカテゴリダミー化も自動。 LightGBM 版は数値・カテゴリ・欠損値もそのまま処理可能。

Q5. 時系列データでも使える？

A5. 使える。 ただし誤差の独立性が崩れるので、 信頼区間にブロックブートストラップを使うこと。 ボラティリティ予測やリスク管理で CAViaR（Conditional Autoregressive VaR、 Engle & Manganelli 2004）という拡張モデルが有名。

Q6. ニューラルネットでもできる？

A6. できる。 損失関数を pinball loss に変えるだけ。 PyTorch なら：

import torch

def pinball_loss(y_pred, y_true, tau):
    err = y_true - y_pred
    return torch.max(tau * err, (tau - 1) * err).mean()

# 学習ループの中で：
#   loss = pinball_loss(model(x), y, tau=0.9)
#   loss.backward(); optimizer.step()

Q7. 解釈は OLS と同じでいい？

A7. ほぼ同じ。 「$x_1$ が 1 単位増えると、 $y$ の τ 分位点 が $\beta_1(\tau)$ だけ増える」と読む。 OLS の「平均が増える」を「τ 分位点が増える」に置き換えるだけ。

Q8. 計算コストは？

A8. 線形分位点回帰は LP として解けるので、 n=10,000、 p=50 程度なら 1 秒以内。 OLS の正規方程式解法に比べれば遅いが、 LightGBM 版なら GBDT と同等。

🗺 概念マップ — 分位点回帰の位置づけ