「100年に一度の洪水」「観測史上最高気温」「想定される最大地震動」——私たちはニュースで、平均値ではなく極端な値にしばしば直面します。ダムの堤防は平均的な雨量ではなく、何十年に一度の豪雨に耐えるよう設計しなければなりません。保険会社は平均的な損害額ではなく、ごく稀に起こる巨額損害に備えて保険料を決めます。こうした場面で本当に知りたいのは「平均がどうか」ではなく「最大値がどこまで大きくなりうるか」です。
平均値の振る舞いは中心極限定理が見事に説明してくれます。たくさんの確率変数を足し合わせれば正規分布に近づく——これは統計学の土台です。では、足し算ではなく最大値を取ったらどうなるのでしょうか？ 驚くことに、最大値にも中心極限定理に対応する深い法則があります。適切に標準化した最大値の分布は、母集団がどんな形であっても、たった3種類の分布のいずれかに収束するのです。その代表格が本記事の主役、ガンベル分布（Gumbel distribution） です。
この理論は実務で広く使われています。河川工学では年最大流量をガンベル分布に当てはめて「T年に一度の流量」を推定し、堤防の高さを決めます。気象・地震・風工学でも極端現象のリスク評価に必須です。さらに意外な応用先として、深層学習があります。カテゴリカル分布から離散サンプリングしつつ勾配を流したい——この難題を解く Gumbel-Max トリック と Gumbel-Softmax は、まさにガンベル分布の性質を巧みに使った技法です。
本記事の内容
- 最大値の分布がどう振る舞うかの直感と、極値理論（Fisher-Tippett-Gnedenkoの定理）
- ガンベル分布の密度・平均・分散・分布関数の導出
- 一般化極値分布（GEV）の3型 — ガンベル・フレシェ・ワイブル
- ブロック最大値法と「再現期間（100年に一度）」の数学
- 機械学習のGumbel-MaxトリックとGumbel-Softmaxの導出
- Pythonでの実装、洪水データへのフィット、離散サンプリングの実演
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
最大値の分布を直感でつかむ
まず、最大値を取るとは何が起きるのかを実験で見てみましょう。標準正規分布から $n$ 個の値をサンプルし、その最大値 $M_n = \max(X_1, \dots, X_n)$ を記録します。これを何十万回も繰り返すと、$M_n$ 自体の分布が描けます。
この図から3つのことが読み取れます。第一に、$n$ が増えるほど最大値の分布は右へ移動します。たくさん引けばより大きな外れ値に当たりやすいので当然です。第二に、分布の形が左右非対称になります。最大値は下に裾を引かず、右側（大きい方）に長い裾を持ちます。第三に、$n$ が大きくなっても分布の幅はさほど広がらず、形がある一定の非対称な山に落ち着いていくように見えます。
この「落ち着いていく形」こそがガンベル分布です。平均値が正規分布に収束するように、最大値も（標準化すれば）ある決まった分布に収束する——この事実が極値統計の出発点になります。
なぜ最大値に注目するのでしょうか。それは、リスク評価のほとんどが「最悪のケース」に関わるからです。1年間の日降水量の平均を知っても堤防は設計できません。知りたいのは「過去最大の日降水量」やそれを超える確率です。次に、この最大値の分布を数学的に特徴づける理論——極値理論——を見ていきます。
極値理論とFisher-Tippett-Gnedenkoの定理
中心極限定理の主張を思い出しましょう。独立同分布な確率変数 $X_1, \dots, X_n$ の和を適切に標準化すると、$n \to \infty$ で正規分布に収束します。極値理論は、これの「最大値版」です。
最大値 $M_n = \max(X_1, \dots, X_n)$ をそのまま見ても、$n\to\infty$ で母集団の上限（無限大のこともある）に発散してしまうので意味がありません。そこで、和のときに $\frac{1}{\sqrt n}$ で縮めて平均を引いたように、最大値も適切な定数 $a_n>0,\ b_n$ で標準化します。
$$ \frac{M_n – b_n}{a_n} $$
このとき次が成り立ちます。これが極値理論の中心定理です。
Fisher-Tippett-Gnedenkoの定理: ある正の定数列 $a_n$ と実数列 $b_n$ が存在して、標準化した最大値 $(M_n – b_n)/a_n$ が $n\to\infty$ で退化しない（一点に潰れない）分布 $G$ に収束するとき、$G$ は次の3型のいずれかに限られる。 1. ガンベル型（Type I） 2. フレシェ型（Type II） 3. ワイブル型（Type III、上限つき）
これは驚くべき主張です。母集団の分布が正規だろうと指数だろうとパレートだろうと、最大値の極限分布はたった3種類しかありえない。中心極限定理で「足し算の極限は正規分布だけ」となるのと同じ普遍性が、最大値にも成り立つのです。
なぜ3型に限られるのか、その核心は「最大値の安定性」にあります。$2n$ 個の最大値は、$n$ 個ずつの最大値2つの、そのまた最大値に等しい。
$$ \max(X_1,\dots,X_{2n}) = \max\bigl(\max(X_1,\dots,X_n),\ \max(X_{n+1},\dots,X_{2n})\bigr) $$
この「最大値を取っても同じ族に留まる」という自己再現性（max-stability）を満たす分布を分類すると、ちょうど3型が出てきます。導出の詳細は本記事の範囲を超えますが、和に対する安定分布（正規分布など）の議論と完全に対応します。
この3型は1つの式にまとめられます。それが次に見る一般化極値分布です。まずは3型の代表であるガンベル分布の性質を、丁寧に調べていきましょう。
ガンベル分布の定義と性質
ガンベル分布は「指数分布のように軽い裾を持つ母集団（正規・指数・対数正規など）」の最大値が収束する先です。日常的には、年最大降水量や年最高気温のモデルとして最もよく使われます。
定義式に入る前に、形のイメージを持っておきましょう。ガンベル分布は右に長い裾を引く非対称な山です。最大値という「上に飛び出した値」を表すので、上側（大きい方）に裾を伸ばすのは自然です。
左のパネルは位置パラメータ $\mu$ を動かしたもので、形を変えずに横へスライドします。右は尺度パラメータ $\beta$ を動かしたもので、$\beta$ が大きいほど山が低く広くなります。どちらの場合も「右に長い裾」という非対称性は保たれます。正規分布の左右対称な釣鐘型とは明確に違うことが見て取れます。
累積分布関数と密度関数
ガンベル分布（最大値型、$\mathrm{Gumbel}_{\max}$）の累積分布関数（CDF） は次で定義されます。位置パラメータ $\mu \in \mathbb{R}$、尺度パラメータ $\beta>0$ です。
$$ \begin{equation} F(x) = \exp\!\left(-\exp\!\left(-\frac{x-\mu}{\beta}\right)\right) \end{equation} $$
二重指数の形が特徴的で、ガンベル分布が「二重指数分布」とも呼ばれる所以です。この式が「最大値の分布」として自然に出てくることは、後ほど正規分布の最大値から導いて確認します。密度関数（PDF）は CDF を微分して得ます。標準化変数 $z = (x-\mu)/\beta$ と置くと、合成関数の微分で
$$ \begin{align} f(x) = \frac{dF}{dx} &= \frac{d}{dx}\exp\!\left(-e^{-z}\right) \\ &= \exp\!\left(-e^{-z}\right) \cdot \left(-\frac{d}{dx}e^{-z}\right) \end{align} $$
ここで $\dfrac{d}{dx}e^{-z} = e^{-z}\cdot\dfrac{d(-z)}{dx} = e^{-z}\cdot\left(-\dfrac{1}{\beta}\right)$ を代入すると、
$$ \begin{equation} f(x) = \frac{1}{\beta}\exp\!\left(-\frac{x-\mu}{\beta} – \exp\!\left(-\frac{x-\mu}{\beta}\right)\right) \end{equation} $$
が得られます。$x \to +\infty$ では指数項 $e^{-z}\to 0$ となり、密度はおよそ $\frac{1}{\beta}e^{-z}$ のように指数的に減衰します。一方 $x\to-\infty$ では二重指数 $e^{-e^{-z}}$ が猛烈に速くゼロになります。この非対称な減衰が「右に長い裾」を作っているわけです。
平均・分散・最頻値
ガンベル分布の主要な統計量は次のとおりです。導出にはガンマ関数の微分が絡むので結果のみ示します。
$$ \begin{align} \mathbb{E}[X] &= \mu + \gamma\beta, \quad \gamma \approx 0.5772 \ (\text{オイラー・マスケローニ定数}) \\ \mathrm{Var}[X] &= \frac{\pi^2}{6}\beta^2 \\ \text{最頻値（モード）} &= \mu \end{align} $$
平均がモード $\mu$ より $\gamma\beta$ だけ右にずれている点に注目してください。これは分布が右に裾を引いているため、平均が裾に引っ張られることの現れです。左右対称な正規分布では平均=モードですが、ガンベルではそうなりません。分散が尺度 $\beta$ の2乗に比例し、$\pi^2/6 \approx 1.645$ という係数が付くのも特徴です。
ここまでで、3型の代表であるガンベル分布の素性がわかりました。次に、3型を1本の式で統一する一般化極値分布を導入し、ガンベルがその特別な場合であることを見ます。
一般化極値分布（GEV）の3型
Fisher-Tippettの3型は、別々に覚える必要はありません。形状パラメータ $\xi$ を1つ導入すると、3型はすべて次の1本の式（一般化極値分布、GEV）に統合されます。
$$ \begin{equation} F(x) = \exp\!\left\{-\left[1 + \xi\left(\frac{x-\mu}{\sigma}\right)\right]^{-1/\xi}\right\}, \quad 1 + \xi\frac{x-\mu}{\sigma} > 0 \end{equation} $$
ここで $\mu$ は位置、$\sigma>0$ は尺度、$\xi$ は形状（裾の重さ）を決めるパラメータです。$\xi$ の符号によって3型に分かれます。
- $\xi > 0$：フレシェ型 — べき乗則のような重い裾。下限はあるが上限なし。パレート分布やコーシー分布、$t$分布など重い裾を持つ母集団の最大値はこれに収束します。
- $\xi < 0$：ワイブル型 — 有限の上限を持つ。一様分布のように上が打ち切られた母集団の最大値はこれに収束します。
- $\xi = 0$：ガンベル型 — その中間。指数的な軽い裾。$\xi\to 0$ の極限を取ると、先ほどのガンベルCDFに一致します。
$\xi\to0$ の極限がガンベルになることを確認しましょう。$\xi\to0$ のとき、内側のべき乗を対数で評価します。$y = 1 + \xi z$（ただし $z=(x-\mu)/\sigma$）として、
$$ \begin{align} \left(1 + \xi z\right)^{-1/\xi} &= \exp\!\left(-\frac{1}{\xi}\ln(1+\xi z)\right) \end{align} $$
ここで $\xi\to0$ のとき $\ln(1+\xi z) \approx \xi z – \frac{(\xi z)^2}{2} + \cdots$ とテイラー展開できるので、
$$ \begin{align} -\frac{1}{\xi}\ln(1+\xi z) &\approx -\frac{1}{\xi}\left(\xi z – \cdots\right) \xrightarrow{\xi\to0} -z \end{align} $$
したがって $(1+\xi z)^{-1/\xi} \to e^{-z}$ となり、GEVのCDFは $\exp(-e^{-z})$、つまりガンベル分布のCDFに一致します。これで「ガンベルはGEVの $\xi=0$ という特別な場合」だと確認できました。
3型の密度を並べて見比べてみましょう。
この図から、$\xi$ の符号が裾の性質を支配することがわかります。赤のフレシェ型（$\xi>0$）は右に最も重い裾を引き、極端に大きな値が出やすいことを示します。緑のワイブル型（$\xi<0$）は右側がスパッと打ち切られ、ある上限を超える値は決して出ません。青のガンベル型（$\xi=0$）はちょうど中間で、裾は指数的に減ります。
裾の重さの違いは片対数プロットにするとさらに鮮明です。
縦軸を対数にした生存関数 $1-F(x)$ を見ると、フレシェ型（赤）はほぼ直線的にしか減らず、極端値の確率が桁違いに大きいことがわかります。ガンベル型（青）は下に凸の曲線で、フレシェより速く減衰します。ワイブル型（緑）はある $x$ で確率が完全にゼロになり、上限の存在がはっきり見えます。
実務では、対象がどの型に従うかを $\xi$ の推定を通じて判断します。河川流量や風速は弱い正の $\xi$（軽いフレシェ寄り）になることが多く、$\xi$ を固定せず推定するのが安全です。では、観測データからこのGEV（やガンベル）をどう当てはめるのか——その標準手法を次に見ます。
ブロック最大値法と再現期間
実データからガンベル/GEVを推定する最も基本的な手法がブロック最大値法（Block Maxima） です。アイデアは極めてシンプルです。長い時系列を等しい長さのブロック（例: 1年ごと）に区切り、各ブロックの最大値だけを取り出します。極値理論より、その最大値たちはGEVに従うので、それらにGEVを当てはめます。
左のパネルでは、灰色の日次観測値を30日ごとのブロックに区切り、各ブロックで最も大きい値（赤点）を1つずつ拾っています。右のパネルは、こうして集めた最大値のヒストグラムで、ガンベル分布（青線）がよく当てはまっていることがわかります。「全データではなく各ブロックの最大値だけを使う」——これが極値理論を正しく適用するための鍵です。中央付近の平凡なデータを混ぜると、それは最大値の分布ではなくなってしまうからです。
再現期間：100年に一度の意味
防災で頻出する「100年に一度の洪水」という言い回しを、数式で正確に定義しましょう。年最大値が分布関数 $F$ に従うとき、ある水準 $z$ を1年間に超える確率は $1-F(z)$ です。「平均して $T$ 年に一度だけ $z$ を超える」とは、この超過確率が $1/T$ であることを意味します。
$$ 1 – F(z_T) = \frac{1}{T} $$
この $z_T$ をリターンレベル（return level）、$T$ を再現期間（return period） と呼びます。ガンベル分布の場合、CDF を代入して $z_T$ について解けます。$F(z_T) = 1 – 1/T$ より、
$$ \begin{align} \exp\!\left(-\exp\!\left(-\frac{z_T-\mu}{\beta}\right)\right) &= 1 – \frac{1}{T} \end{align} $$
両辺の対数を2回取ります。まず1回目の対数で
$$ -\exp\!\left(-\frac{z_T-\mu}{\beta}\right) = \ln\!\left(1 – \frac{1}{T}\right) $$
符号を整理して、もう一度対数を取ると（左辺の指数が外れます）、
$$ -\frac{z_T-\mu}{\beta} = \ln\!\left(-\ln\!\left(1-\frac{1}{T}\right)\right) $$
最後に $z_T$ について解くと、ガンベル分布のリターンレベルが求まります。
$$ \begin{equation} z_T = \mu – \beta\,\ln\!\left(-\ln\!\left(1 – \frac{1}{T}\right)\right) \end{equation} $$
この式を $T$ の関数として描いたのが再現期間プロットです。
横軸（再現期間）を対数軸にすると、リターンレベルはほぼ直線的に増えます。これは $T$ が大きいとき $\ln(1-1/T)\approx -1/T$ より $z_T \approx \mu + \beta\ln T$ となるためです。図の例（$\mu=100,\ \beta=20$）では、10年に一度が約145、100年に一度が約192、500年に一度が約227です。注目すべきは、再現期間を10倍にしてもリターンレベルは比例的には増えず、$\beta\ln 10 \approx 46$ ずつ足し算的に増える点です。「1000年に一度」の備えが「100年に一度」のたかだか1.2倍程度で済むのは、この対数の効き方ゆえです。
ここまでで、ガンベル/GEVを使ってリスクを定量化する流れがつかめました。ここで話を大きく変えます。実はこのガンベル分布、災害統計とはまったく無縁に見える深層学習でも主役を演じます。次の章でその意外なつながりを見ていきましょう。
Gumbel-Max トリック
深層学習でカテゴリカル分布からサンプリングしたい場面を考えます。たとえば言語モデルが次の単語を選ぶとき、強化学習でエージェントが行動を選ぶとき、あるいは離散的な潜在変数を持つVAEを学習するとき。確率 $p_1, \dots, p_K$（$\sum_k p_k=1$）からカテゴリ $k$ を1つサンプルしたいわけです。
通常これは累積確率と一様乱数で実装できますが、ニューラルネットの学習では困った問題が起きます。サンプリング（argmax や離散選択）は微分不可能なので、誤差逆伝播で勾配を流せないのです。この壁を鮮やかに越えるのが Gumbel-Max トリックです。
トリックの主張はこうです。カテゴリ $k$ のロジット（対数確率に相当）を $\ell_k = \log p_k$ とします。各カテゴリに独立な標準ガンベル雑音 $g_k$ を足し、その最大値のインデックスを取ると、それはちょうど確率 $p_k$ でカテゴリ $k$ を選んだことになります。
$$ \begin{equation} k^\ast = \operatorname*{argmax}_{k}\ \bigl(\ell_k + g_k\bigr), \quad g_k \overset{\text{iid}}{\sim} \mathrm{Gumbel}(0,1) \end{equation} $$
これが本当に確率 $p_k$ を再現することを実験で確かめましょう。
左のパネルは1回のサンプリングの内訳です。各カテゴリのロジット（青）にガンベル雑音（赤）を足し、その和（黒点）が最大のカテゴリ（星印）を選びます。右のパネルは、これを4万回繰り返した選択頻度（オレンジ）と、ロジットから計算したsoftmax確率（青）の比較です。両者がぴたりと一致しています。つまり Gumbel-Max トリックは、softmax 分布からの厳密なサンプリングになっているのです。
なぜ成り立つのか — 導出
この一致は偶然ではありません。導出してみましょう。$\mathrm{Gumbel}(0,1)$ のCDFは $F(g)=e^{-e^{-g}}$、PDFは $f(g)=e^{-g}e^{-e^{-g}}$ です。$Y_k = \ell_k + g_k$ と置くと、$Y_k$ は位置 $\ell_k$ のガンベル分布で、CDFは $\Pr(Y_k \le y)=\exp(-e^{-(y-\ell_k)})$ です。カテゴリ $k$ が選ばれる確率は「$Y_k$ が他のすべての $Y_j$ より大きい確率」なので、$y=Y_k$ で条件づけて積分します。
$$ \begin{align} \Pr(k^\ast = k) &= \int_{-\infty}^{\infty} f_{Y_k}(y) \prod_{j\ne k}\Pr(Y_j < y)\, dy \end{align} $$
各因子を代入します。$f_{Y_k}(y) = e^{-(y-\ell_k)}\exp(-e^{-(y-\ell_k)})$、$\Pr(Y_j $$
\begin{align}
\prod_{j} \exp\!\left(-e^{-(y-\ell_j)}\right)
&= \exp\!\left(-\sum_{j} e^{-y}e^{\ell_j}\right)
= \exp\!\left(-e^{-y}\sum_j e^{\ell_j}\right)
\end{align}
$$ ここで $S = \sum_j e^{\ell_j} = \sum_j p_j = 1$（ロジットが対数確率なら）ですが、一般性のため $S$ のまま進めます。被積分関数は $$
\Pr(k^\ast=k) = \int_{-\infty}^{\infty} e^{-(y-\ell_k)} \exp\!\left(-e^{-y} S\right) dy
$$ となります。$u = e^{-y}$ と置換すると $du = -e^{-y}dy$、すなわち $e^{-y}dy = -du$ で、積分範囲は $y:-\infty\to\infty$ が $u:\infty\to0$ に変わります。$e^{\ell_k}$ は定数として外に出して、 $$
\begin{align}
\Pr(k^\ast=k)
&= e^{\ell_k}\int_{0}^{\infty} e^{-uS}\, du
= e^{\ell_k}\cdot\frac{1}{S}
= \frac{e^{\ell_k}}{\sum_j e^{\ell_j}}
\end{align}
$$ きれいに softmax が出てきました。すなわち $\Pr(k^\ast=k) = \mathrm{softmax}(\boldsymbol{\ell})_k = p_k$ です。ガンベル分布の二重指数という特殊な形が、積分を $u$ の単純な指数積分に潰してくれたのがポイントです。これがガンベル分布が深層学習で重宝される理由です。 ただし argmax 自体はまだ微分できません。雑音の足し算で「サンプリング元の乱数性」を入力側に追い出せた（再パラメータ化）のは大きな前進ですが、勾配を流すには argmax を滑らかにする必要があります。それを実現するのが次の Gumbel-Softmax です。 Gumbel-Max の argmax を、微分可能な softmax で置き換えたのが Gumbel-Softmax（別名 Concrete 分布）です。温度パラメータ $\tau>0$ を導入し、one-hot ベクトルの代わりに次の「ソフトな one-hot」を出力します。 $$
\begin{equation}
y_k = \frac{\exp\bigl((\ell_k + g_k)/\tau\bigr)}{\sum_{j}\exp\bigl((\ell_j + g_j)/\tau\bigr)}, \quad g_k \sim \mathrm{Gumbel}(0,1)
\end{equation}
$$ この $\boldsymbol{y}=(y_1,\dots,y_K)$ は要素が非負で和が1のベクトル、つまり確率単体の上の点です。温度 $\tau$ が出力の「鋭さ」を制御します。 この図で温度を $\tau=2.0\to0.5\to0.1$ と下げると、出力が次第に1つのカテゴリに集中し、one-hot ベクトル（離散サンプル）に近づくのが見て取れます。$\tau\to0$ の極限では softmax は argmax に一致し、Gumbel-Max と同じ厳密なサンプルになります。逆に $\tau$ が大きいと出力は一様に近づき、すべてのカテゴリに薄く広がります。 実用上の使い方は「温度を徐々に下げていくアニーリング」です。学習初期は $\tau$ を大きめにして勾配を安定に流し、訓練が進むにつれ $\tau$ を下げて離散サンプルに近づけます。これにより、離散的なサンプリングを含むモデルを誤差逆伝播で end-to-end に学習できます。離散潜在変数VAE、微分可能なアーキテクチャ探索、ハードアテンションなどで広く使われています。 理論からPythonに移り、ここまでの内容を実際に動かして確かめましょう。 まず、ガンベル分布の密度・CDF・統計量を scipy で確認します。 出力は次のようになります。平均が 平均がモード（=0）より右にずれて約0.577になっているのが、右に裾を引く非対称性の現れです。分散の係数 $\pi^2/6\approx1.645$ も理論どおりで、scipy の実装が定義式と完全に整合していることが確認できました。 次に、極値理論の核心——正規母集団の最大値が標準化するとガンベルに収束する——を数値実験で確かめます。正規分布の標準化定数は理論的に $a_n = 1/\sqrt{2\ln n}$、$b_n = \sqrt{2\ln n} – \frac{\ln\ln n + \ln 4\pi}{2\sqrt{2\ln n}}$ と知られています。 出力例は次のとおりです。 標準化した最大値の平均がガンベルの理論平均 $\gamma\approx0.577$ にほぼ一致し、KS統計量も小さく、分布が標準ガンベルとよく合っていることがわかります。これを図にすると一致が一目瞭然です。 灰色のヒストグラム（標準化した最大値）が、赤い標準ガンベル曲線にぴたりと重なっています。$n=1000$ という有限のサンプル数でも収束がよく見えます（正規分布の収束は理論的には遅いことが知られていますが、標準化定数を正しく取れば実用上は十分です）。これで「正規分布の最大値はガンベルに収束する」という極値理論の主張を、自分の手で確認できました。 実務的な応用として、年最大流量データにGEVを当てはめ、再現期間ごとのリターンレベルを推定してみます。ここでは真のGEV（$\xi=0.15$）から50年分の合成データを生成し、それを推定します。 出力例は次のとおりです。 推定された $\xi=0.121$ は真値0.15に近く、わずか50年分のデータからでも形状パラメータがそれなりに推定できています。「100年に一度」の流量が約195と推定され、防災計画の根拠になります。フィットの様子を図で見てみましょう。 オレンジのGEVフィットがヒストグラムによく沿っています。青破線は $\xi=0$ に固定したガンベルフィットで、こちらも悪くはありませんが、右裾でわずかにGEVと差が出ます。実務では $\xi$ を固定せず推定し、$\xi$ がほぼ0なら裾は指数的、正なら重い裾と判断します。$\xi$ を不用意に0と決め打ちすると、極端なリターンレベルを過小評価する危険がある点に注意が必要です。 最後に、機械学習側の応用を実装します。まず Gumbel-Max トリックが本当に softmax からのサンプリングになっているかを確認します。 出力例は次のとおりです。 経験頻度がsoftmax理論値に小数点以下3桁まで一致しており、先ほど導出した「Gumbel-Max は softmax からの厳密サンプリング」が数値的に裏付けられました。続いて、微分可能な Gumbel-Softmax を実装し、温度の効果を見ます。 出力例は次のとおりです。 温度 $\tau$ を下げるほど出力が1つの成分に集中し、$\tau=0.1$ ではほぼ完全な one-hot になっています。高温では確率が複数カテゴリに分散します。重要なのは、この出力が 本記事では、ガンベル分布と極値統計を、災害リスク評価から深層学習まで横断して解説しました。 ガンベル分布は「最大値」という単純な操作から、防災工学・保険・気候科学・深層学習をつなぐ普遍的な道具です。次のステップとして、重い裾そのものを扱う分布や、離散サンプリングの理論を深めると理解が広がります。 次に読むと理解が深まる記事も挙げておきます。Gumbel-Softmax（Concrete分布）
Pythonでの実装
ガンベル分布の基本量
scipy.stats.gumbel_r が最大値型のガンベル分布です（gumbel_l は最小値型）。import numpy as np
from scipy import stats
mu, beta = 0.0, 1.0
dist = stats.gumbel_r(loc=mu, scale=beta)
# 理論値と数値積分の比較
gamma = 0.5772156649 # オイラー・マスケローニ定数
print("平均 理論:", mu + gamma * beta, " scipy:", dist.mean())
print("分散 理論:", (np.pi**2 / 6) * beta**2, " scipy:", dist.var())
print("最頻値（モード）= mu =", mu, " （密度最大の x）")
# CDF が二重指数になっているか
x = 1.5
print("CDF(1.5):", dist.cdf(x), " 手計算:", np.exp(-np.exp(-(x - mu) / beta)))
0.5772（=$\gamma\beta$）、分散が 1.6449（=$\pi^2/6$）、CDF が手計算の二重指数と一致します。平均 理論: 0.5772156649 scipy: 0.5772156649015329
分散 理論: 1.6449340668482264 scipy: 1.6449340668482264
最頻値（モード）= mu = 0.0 （密度最大の x）
最適値: CDF(1.5): 0.8003194996137927 手計算: 0.8003194996137927
正規分布の最大値がガンベルに収束する
import numpy as np
from scipy import stats
rng = np.random.default_rng(17)
n = 1000 # 1ブロックのサンプル数
trials = 300000 # 最大値を取る試行回数
maxima = rng.standard_normal((trials, n)).max(axis=1)
# 正規分布の極値正規化定数
ln_n = np.log(n)
a_n = 1.0 / np.sqrt(2 * ln_n)
b_n = np.sqrt(2 * ln_n) - (np.log(ln_n) + np.log(4 * np.pi)) / (2 * np.sqrt(2 * ln_n))
z = (maxima - b_n) / a_n # 標準化
# 標準ガンベルとの一致度（コルモゴロフ-スミルノフ検定）
ks = stats.kstest(z, stats.gumbel_r(loc=0, scale=1).cdf)
print(f"標準化した最大値の平均: {z.mean():.3f} (ガンベル理論 0.577)")
print(f"KS統計量: {ks.statistic:.4f}")
標準化した最大値の平均: 0.575 (ガンベル理論 0.577)
KS統計量: 0.0193
洪水データへのGEVフィット（応用）
import numpy as np
from scipy import stats
rng = np.random.default_rng(7)
# 真のGEV（やや重い裾）から年最大流量を50年分生成
true_c, true_loc, true_scale = -0.15, 90.0, 22.0 # scipyのc = -xi
data = stats.genextreme(true_c, loc=true_loc, scale=true_scale).rvs(50, random_state=rng)
# GEVを最尤推定で当てはめ
c_hat, loc_hat, scale_hat = stats.genextreme.fit(data)
xi_hat = -c_hat
print(f"推定: xi={xi_hat:.3f}, mu={loc_hat:.1f}, sigma={scale_hat:.1f}")
# 再現期間ごとのリターンレベル（GEVの逆関数 = ppf）
fitted = stats.genextreme(c_hat, loc=loc_hat, scale=scale_hat)
for T in [10, 50, 100, 500]:
zT = fitted.ppf(1 - 1.0 / T) # 超過確率 1/T となる水準
print(f"再現期間 {T:>4}年 のリターンレベル: {zT:.1f}")
推定: xi=0.121, mu=89.6, sigma=21.4
再現期間 10年 のリターンレベル: 144.3
再現期間 50年 のリターンレベル: 178.6
再現期間 100年 のリターンレベル: 194.8
再現期間 500年 のリターンレベル: 237.2
Gumbel-Max トリックとGumbel-Softmax
import numpy as np
rng = np.random.default_rng(11)
logits = np.array([1.0, 2.5, 0.5, 2.0])
probs = np.exp(logits) / np.exp(logits).sum() # 目標のsoftmax確率
# Gumbel-Maxトリック: ロジット + ガンベル雑音の argmax
N = 200000
g = rng.gumbel(size=(N, len(logits))) # 標準ガンベル雑音
choices = np.argmax(logits + g, axis=1)
freq = np.bincount(choices, minlength=len(logits)) / N
print("softmax 理論:", np.round(probs, 4))
print("経験頻度 :", np.round(freq, 4))
softmax 理論: [0.1131 0.5069 0.0686 0.3114]
経験頻度 : [0.1138 0.5066 0.0683 0.3113]
import numpy as np
def gumbel_softmax(logits, tau, rng):
"""温度tauのGumbel-Softmaxサンプルを返す（微分可能な緩和）"""
g = rng.gumbel(size=logits.shape)
z = (logits + g) / tau
z -= z.max() # 数値安定化
y = np.exp(z)
return y / y.sum()
rng = np.random.default_rng(13)
logits = np.array([1.0, 2.5, 0.5, 2.0])
for tau in [2.0, 0.5, 0.1]:
# 同じ乱数列を使うため毎回シードを揃える
sub = np.random.default_rng(13)
y = gumbel_softmax(logits, tau, sub)
print(f"tau={tau}: {np.round(y, 3)} (最大成分 {y.max():.3f})")
tau=2.0: [0.21 0.327 0.157 0.306] (最大成分 0.327)
tau=0.5: [0.041 0.62 0.012 0.327] (最大成分 0.620)
tau=0.1: [0. 1. 0. 0. ] (最大成分 1.000)
logits について微分可能だという点です。argmax なら勾配がゼロですが、softmax 緩和のおかげで勾配が流れ、ニューラルネットの学習に組み込めます。これが離散変数を含むモデルを誤差逆伝播で学習する仕組みの核心です。まとめ