ここまでのシリーズでは、制御入力やコマンドという「操作の情報」が観測できる前提でした。残差ベースは制御入力 $u$ を、状態遷移はアクチュエータ指令を、条件付きフローはモードラベル $c$ を使いました。
しかし現実には、いまどの運用モードにいるかのラベルが手に入らないことがよくあります。化学プラントのグレード切替、バッチ工程の段階、設備の負荷状態——センサ値の履歴はあっても、「この時刻はモードA」という正解ラベルは付いていません。
このとき登場するのが、サーベイで挙げた (い) マルチモード監視 の王道、有限混合ガウスモデル（FGMM）＋ベイズ推論です。Yu と Qin が AIChE Journal（2008）で確立したこの手法は、化学工学のプロセス監視で長く使われてきました。アイデアは明快です——データから運用モードを自動で推定し、モードごとに別々の正常範囲を当て、所属確率で重み付けて融合する。
上図の星印を見てください。これは3つのモードの谷間に落ちる異常です。どのモードの正常データからも遠いのに、全体をひとまとめにした単一の正常範囲では「真ん中だから正常」と判定されてしまいます。この見逃しを防ぐのが本手法です。
この手法が効く場面は明確です。
- モードラベルが無い多変量プロセス: 教師なしでモードを発見して監視したい
- モードが滑らかに切り替わる系: ソフトな所属確率で遷移を扱いたい
- モードごとに相関構造が違う系: 各モードに専用の共分散を当てたい
この記事の内容
- FGMMで運用モードを推定する — 各ガウス成分＝1モード
- 成分数 $K$ の自動選択（BIC / Figueiredo–Jain）
- 局所マハラノビス距離とカイ二乗による「正常確率」
- ベイズ推論による融合 — 所属事後確率で局所統計を重み付け（BIP）
- 遷移（過渡）モードの扱い
- scikit-learn + scipy 実装 — 無条件の単一ガウス（AUC0.08）に対しFGMM+ベイズ融合（0.96）が捉える様子を実測
前提知識
1. FGMMで運用モードを推定する
モードラベルが無いなら、データそのものからモードを見つけます。正常運転データが複数のクラスタからなると考え、それを有限個のガウス分布の混合で表現します。
$$ \begin{equation} p(x) = \sum_{k=1}^{K} \pi_k\, \mathcal{N}(x \mid \mu_k, \Sigma_k) \end{equation} $$
各ガウス成分 $\mathcal{N}(x\mid \mu_k, \Sigma_k)$ が1つの運用モードに対応し、混合比 $\pi_k$ がそのモードの出現頻度です。パラメータはEMアルゴリズムで推定します。
当てはめると、各モードのクラスタにぴったりガウスが乗ります。重要なのは、各モードが別々の平均と共分散（相関構造）を持てること。モードAは2変数が正相関、モードBは負相関、といった違いをそのまま表現できます。
ここで「成分数 $K$ をいくつにするか」という問題が残ります。
$K$ を増やせばデータへの当てはまりは良くなりますが、過剰適合します。そこでBIC（ベイズ情報量規準）のようなモデル選択規準を使い、当てはまりと複雑さのバランスが最良の $K$ を選びます。上図では BIC が $K=3$ で最小になり、正しくモード数を言い当てています。
ここで上のBIC法は「$K=1,2,3,\dots$ を別々に全部当てはめてから一番良いものを選ぶ」総当たりです。Yu & Qin の原論文はもっと洗練された Figueiredo–Jain（FJ）アルゴリズムを使います。これは「当てはめ」と「成分数選び」を一つのEMの中で同時にやってしまう点が新しいところです。
Figueiredo–Jain: 成分を「育てる」のでなく「間引く」
発想を逆転させます。少ない $K$ から増やすのではなく、過剰に多い成分数 $K_{\max}$ から出発し、要らない成分を学習の途中で消していくのです。判断基準は MML（最小記述長, Minimum Message Length）。モデルとデータを符号化するのに必要な「総ビット数」を最小化します。
$$ \begin{equation} \mathcal{L}(\theta, K) = \underbrace{\frac{N_p}{2}\sum_{k:\,\pi_k>0} \log\!\frac{n\,\pi_k}{12} + \frac{K_{\text{nz}}}{2}\log\frac{n}{12} + \frac{K_{\text{nz}}(N_p+1)}{2}}_{\text{モデルの符号化長(複雑さの罰則)}} \;\underbrace{-\,\log p(X\mid\theta)}_{\text{データの符号化長(あてはまり)}} \end{equation} $$
ここで $n$ はデータ数、$N_p$ は1成分あたりの自由パラメータ数、$K_{\text{nz}}$ は重み $\pi_k>0$ の生き残り成分数です。第1項に注目してください。混合比 $\pi_k$ が小さい成分は $\log(n\pi_k/12)$ が負に大きくなり、罰則を下げる方向に働きます——つまり重みの小さい成分は「消した方が記述長が短くなる」ように設計されています。
具体的な成分ごとのEM更新（component-wise EM）では、各成分の重みが
$$ \begin{equation} \hat{\pi}_k = \frac{\max\!\left(0,\ \sum_i r_{ik} – \tfrac{N_p}{2}\right)}{\sum_j \max\!\left(0,\ \sum_i r_{ij} – \tfrac{N_p}{2}\right)} \end{equation} $$
と更新されます（$r_{ik}$ は責任度＝所属事後確率）。分子の $-N_p/2$ がポイントで、責任度の総和がパラメータ数の半分に満たない成分は重みが $0$ になり、自動的に脱落します。EMを回すうちに余剰成分が次々と $\pi_k=0$ に落ち、最終的に必要なモード数だけが残ります。
左図は $K_{\max}=8$ から出発した様子です。青い実線が生き残った3つの成分（＝3モード）、赤い破線がプルーニングされた余剰成分。右図では、反復が進むにつれMMLコスト（青）が単調に下がり、それと連動して有効成分数 $K$（橙）が $8\to3$ へ階段状に減っていきます。BIC法のように $K$ を外から総当たりせず、一回の学習で「モードはいくつか」と「各モードの形」が同時に決まる——これがFJの強みです。なお scikit-learn の BayesianGaussianMixture（変分ベイズGMM, VBGMM）も同じ思想で、Dirichlet過程事前により不要成分の重みを自動でゼロ近傍に押し下げます。後半で触れる新しい手法はこのVBGMMを使います。
モードが推定できたら、新しいデータがどのモードに属するかを確率で測ります。
2. 所属事後確率 — どのモードにいるか
新しい観測 $x$ が来たら、各モードへの所属事後確率をベイズの定理で計算します。
$$ \begin{equation} P(\text{mode}_k \mid x) = \frac{\pi_k\, \mathcal{N}(x \mid \mu_k, \Sigma_k)}{\sum_{j} \pi_j\, \mathcal{N}(x \mid \mu_j, \Sigma_j)} \end{equation} $$
これは「この観測は、どのモードから来た可能性が高いか」を表します。
各モードの中心付近では所属確率がほぼ1（はっきりそのモード）、モードの境界付近では複数モードに確率が割れます（ソフト割当）。この「ソフトさ」が、後で遷移を扱うときに効いてきます。
ポイントは、いまどのモードにいるかを一つに決めつけないこと。所属確率という形で曖昧さを保ったまま、次の異常判定に渡します。
3. 局所統計とベイズ融合（BIP）
3.1 各モードでの正常確率
モード $k$ から見て、観測 $x$ がどれだけ典型的かをマハラノビス距離で測ります。
$$ \begin{equation} D_k(x) = (x – \mu_k)^\top \Sigma_k^{-1} (x – \mu_k) \end{equation} $$
正規分布のもとでは、この距離の二乗は自由度 $d$ のカイ二乗分布に従います。そこで「そのモードで、この距離以上に外れる確率」を計算します。
$$ \begin{equation} P^k(x) = P\!\left(\chi^2_d > D_k(x)\right) \end{equation} $$
$P^k$ が大きい（1に近い）なら「モード $k$ にとってありふれた値＝正常」、小さい（0に近い）なら「モード $k$ にとって極端＝異常」です。これがモードごとの局所的な正常度です。
3.2 ベイズ融合
最後に、所属事後確率で局所正常確率を重み付けて融合します。これが Yu & Qin のベイズ推論ベースの確率指標（BIP, Bayesian Inference Probability）です。
$$ \begin{equation} \mathrm{BIP}(x) = \sum_{k=1}^{K} P(\text{mode}_k \mid x)\, P^k(x) \end{equation} $$
この式の美しさは、「いまどのモードか」を陽に決めなくても、自動的に最も近いモードの基準で評価されることです。観測がモードAの中心にあれば $P(\text{mode}_A|x)\approx 1$ なのでモードAの正常確率がそのまま効き、境界付近なら複数モードの基準が滑らかに混ざります。BIP が小さければ「どのモードの基準でも珍しい」＝異常です。
この「滑らかに混ざる」性質が、マルチモード監視の最難所である遷移を緩和します。
4. 発展：VBGMM × CCA — 「相関構造の崩れ」まで見る
ここまでのFGMM+BIPは、各モードを1つのガウス（平均と共分散）で表し、マハラノビス距離で正常度を測りました。しかし産業プロセスでは、「操作変数 $x$ と品質変数 $y$ の対応関係」そのものが正常性の鍵になることが多い——温度を上げれば収率が上がる、という入出力の相関が崩れたら、各変数は範囲内でも異常です。この相関構造を陽にモデル化するのが、Jiang らの VBGMM-CCA（IEEE, 2019） に代表される新しい流れです。考え方は3段構えになります。
オフライン学習（上段）とオンライン監視（下段）の二段構成です。FGMM+BIPと骨格は同じ——モードを自動同定し、モードごとに局所モデルを当て、所属確率でベイズ融合する——ですが、局所モデルがマハラノビスから「局所CCAモデル」に置き換わるのが核心の工夫です。
4.1 モード同定：VBGMMで自動的に
第1段はモード同定。§1 のFJと同じ思想で、変分ベイズGMM（VBGMM）を使います。VBGMMは混合比 $\pi_k$ にDirichlet事前を置き、変分下界 $\mathcal{L}$（ELBO）を最大化する過程で不要成分の重みを自動的に $0$ 近くへ押し下げます。$K_{\max}$ から出発して有効モード数が自然に決まるので、ラベルも $K$ の事前指定も要りません。
4.2 局所モデル：モードごとにCCA
第2段が新しいところです。モード $k$ に属するデータだけを集め、プロセス変数ブロック $x$ と品質変数ブロック $y$ の間で正準相関分析（CCA）を行います。CCAは、$x$ の線形結合 $u=a^\top x$ と $y$ の線形結合 $v=b^\top y$ の相関 $\rho=\mathrm{corr}(u,v)$ を最大化する射影 $a,b$ を求めます。
左図のように、正常運転では正準変数が $u\approx v$（強く相関）となります。異常が起きて入出力の対応が崩れると、この 残差 $e=u-v$ が拡大します。各変数を単独で見る監視では捉えられない「関係の異常」を、残差1本に凝縮できるわけです。
4.3 監視統計量：T² と SPE
局所CCAモデルから、おなじみの2つの統計量を作ります（右図）。
$$ \begin{equation} T^2_k = z_k^\top \Lambda_k^{-1} z_k, \qquad \mathrm{SPE}_k = \lVert x – \hat{x}_k \rVert^2 \end{equation} $$
- $T^2$（Hotelling統計量）：モデルが捉えた正準部分空間の内側で、正常時のばらつき $\Lambda_k$ に照らしてどれだけ離れているか。右図の星印（横に飛び出した点）が $T^2$ 逸脱です。
- $\mathrm{SPE}$（二乗予測誤差 / Q統計量）：モデルでは説明できない直交残差の大きさ。右図のプラス印（部分空間の外へ立ち上がった点）が $\mathrm{SPE}$ 逸脱です。
両者は「想定内の方向での過大」と「想定外の方向への逸脱」という直交した異常の捉え方で、相補的に働きます。各統計量の管理限界は、正常時の分布（$T^2$ なら $F$ 分布や $\chi^2$、$\mathrm{SPE}$ なら重み付き $\chi^2$ 近似）から決めます。
4.4 ベイズ融合：BIPで1本に
第3段は §3 と同じBIPです。各局所CCAモデルの統計量を「その統計量を超える確率（正常確率）」に変換し、所属事後確率で加重平均して1つの監視指標にまとめます。モード数がいくつあっても、オペレータが見るチャートは1本——これがBIP融合の運用上の利点です。
上はBIPベースの異常スコアを時系列に並べた管理図です。前半（青）は正常運転で、途中にモード切替が含まれていても管理限界（橙の破線）を超えません——これがマルチモード監視の眼目で、切替を異常と誤らないことです。後半（赤）はあるモードから別モード方向への緩慢な逸脱（ステルス異常）で、各点が次々と限界を破って警報を上げています。「正常なモード切替は見逃さず通し、異常な逸脱だけを捉える」——マルチモード監視が達成したかった挙動が、1本のBIP管理図に表れています。
骨格（モード同定→局所モデル→ベイズ融合）を保ったまま、局所モデルを差し替えれば手法が拡張できる——このモジュール性こそ、FGMM+BIPの枠組みが長く使われてきた理由です。
5. 遷移（過渡）モードの扱い
起動・停止・グレード切替の最中、システムはどの定常モードにも当てはまりません。これがFGMMを含むマルチモード監視の弱点です。
ハードにモードを1つ選ぶ方式だと、遷移期は「どのモードにも合わない＝異常」と誤検知しがちです。一方、ベイズ融合のようにソフトな所属確率で複数モードの基準を混ぜると、遷移期にも「行き先のモードの基準」が一定の重みを持つため、誤報が緩和されます。
左のハード割当（argmax でモードを1つに決める）は、境界が断崖になっています。モードAからBへ移る途中の点は、どちらかのモードへ強制的に分類され、その瞬間に「選ばれなかった側のモード基準では大きく外れている」と判定されかねません。右のソフト割当では、境界付近で最大所属確率が $1$ から $1/K$ 近くまで滑らかに下がり（暗い帯）、複数モードの確率が割れています。BIPはこの割れた確率で各モードの正常確率を加重平均するので、遷移期でも「行き先モードの正常域に近づきつつある」評価が連続的に効き、誤報が出にくくなります。「いまどのモードか」を一つに決めつけないことが、そのまま遷移のロバスト性に直結しているわけです。
より積極的な対策として、遷移そのものを独立の「遷移モード」として別の成分でモデル化する、あるいは前記事までの制御入力ベースの予測と組み合わせる、といった拡張もあります。
理屈が揃ったので、実装で効果を確かめます。
6. 実装 — scikit-learn + scipy
3モードのデータ（モードごとに平均・相関が異なる）でFGMMを学習し、BIPで監視します。比較対象は、モードを無視した無条件の単一ガウス・マハラノビスです。
import numpy as np
from scipy.stats import chi2
from sklearn.mixture import GaussianMixture
rng = np.random.default_rng(0)
M, n = 3, 1200
means = np.array([[-3.0, -2.0], [0.5, 3.0], [3.2, -1.6]])
covs = [ (lambda A: A@A.T*0.4 + 0.2*np.eye(2))(rng.normal(0,1,(2,2))) for _ in range(M) ]
X = np.vstack([rng.multivariate_normal(means[k], covs[k], n) for k in range(M)])
perm = rng.permutation(len(X)); X = X[perm]
mu_all, sd_all = X.mean(0), X.std(0); Xs = (X - mu_all) / sd_all
n_tr = int(0.7*len(Xs)); Xtr, Xte_n = Xs[:n_tr], Xs[n_tr:]
異常は「自モードから別モード方向へ部分シフトした」点（＝モード間の谷や別モード寄りに紛れる、ステルスな逸脱）です。
na = len(Xte_n); Xan = []
for _ in range(na):
m = rng.integers(M); w = (m + rng.integers(1, M)) % M
tgt = means[m] + 0.5*(means[w] - means[m])
Xan.append((rng.multivariate_normal(tgt, covs[m], 1)[0] - mu_all) / sd_all)
Xan = np.vstack(Xan)
Xte = np.vstack([Xte_n, Xan]); yte = np.concatenate([np.zeros(len(Xte_n)), np.ones(len(Xan))])
成分数 $K$ をBICで自動選択してFGMMを学習し、BIPスコアを計算します。
Ks = range(1, 9)
gmms = [GaussianMixture(k, covariance_type="full", random_state=0).fit(Xtr) for k in Ks]
kbest = list(Ks)[int(np.argmin([g.bic(Xtr) for g in gmms]))]
gmm = gmms[list(Ks).index(kbest)]
d = Xtr.shape[1]
def bip_score(Xq):
post = gmm.predict_proba(Xq) # P(mode_k | x)
Dk = np.zeros((len(Xq), gmm.n_components))
for k in range(gmm.n_components):
diff = Xq - gmm.means_[k]; inv = np.linalg.inv(gmm.covariances_[k])
Dk[:, k] = np.einsum("ij,jl,il->i", diff, inv, diff) # マハラノビス距離^2
Pk = chi2.sf(Dk, df=d) # そのモードでの正常確率
return -(post * Pk).sum(1) # 異常スコア(大=異常)
s_bip = bip_score(Xte)
# 無条件の単一ガウス マハラノビス
mu0 = Xtr.mean(0); inv0 = np.linalg.inv(np.cov(Xtr.T))
s_glob = np.einsum("ij,jl,il->i", Xte - mu0, inv0, Xte - mu0)
def auc(score, y):
order = np.argsort(score); ranks = np.empty(len(score)); ranks[order] = np.arange(1, len(score)+1)
npos = y.sum(); nneg = len(y) - npos
return (ranks[y==1].sum() - npos*(npos+1)/2) / (npos*nneg)
print(f"自動選択 K = {kbest}")
print(f"無条件マハラノビス AUC={auc(s_glob, yte):.3f}")
print(f"FGMM+ベイズ融合(BIP) AUC={auc(s_bip, yte):.3f}")
実行結果（seed固定）：
自動選択 K = 3
無条件マハラノビス AUC=0.082
FGMM+ベイズ融合(BIP) AUC=0.958
左の無条件の単一ガウスは、3モードを覆うように楕円が大きく広がり、谷の異常を内側に取り込んでしまいます。右のFGMMは各モードを別々に囲み、谷の異常を正しく外に置きます。
数字の読み取りです。
- 無条件マハラノビス＝0.082：ランダム（0.5）を大きく下回っています。モードを無視すると、別モード方向へ紛れた異常は全体平均から見るとむしろ近く、正常な端のデータよりも「正常らしく」評価してしまう——順位づけが逆転しているのです。
- FGMM+ベイズ融合＝0.958：モードを推定して別々に評価することで、「このモードにしてはあり得ない」点をくっきり捉えます。
- K自動選択＝3：BICが正しくモード数を言い当て、ラベル無しでもモード構造を発見できています。
教訓は親サーベイの実験と同じ方向ですが、ここではモードラベルが一切無くても、データからモードを推定して同じ効果が得られる点が核心です。教師なしでモード構造を発見し、ベイズ融合で滑らかに監視する——これがマルチモード監視の王道です。
7. まとめ
モードラベルの無いマルチモードプロセスの異常検知を、FGMM＋ベイズ推論で見てきました。
- 核心：正常データをガウス混合で表し、各成分を運用モードとみなす。所属事後確率で局所統計を融合（BIP）。
- モード推定：EMで混合を学習、BIC / Figueiredo–Jain で成分数 $K$ を自動選択。
- 局所統計：各モードのマハラノビス距離 → カイ二乗で正常確率 $P^k$。
- ベイズ融合：$\mathrm{BIP} = \sum_k P(\text{mode}_k\mid x)\,P^k$。モードを陽に決めず最近モード基準で評価。
- 拡張：局所モデルをマハラノビスから局所CCA（残差 $e=u-v$ と $T^2$/$\mathrm{SPE}$）へ差し替えると、入出力の相関構造の崩れまで監視できる（VBGMM-CCA）。骨格（モード同定→局所モデル→BIP融合）は不変。
- 実証：無条件マハラノビス（0.082）が逆転する一方、FGMM+ベイズ融合（0.958）が捉える。K自動選択も成功。BIP管理図はモード切替を通し、緩慢な逸脱だけを警報。
FGMMは「モードを静的なクラスタとして空間的に分ける」アプローチです。次の記事では、モードを時間的に切り替わる潜在状態として追跡する、切替状態空間モデル／HMMに進みます。
主な参考文献
- J. Yu, S. J. Qin, “Multimode process monitoring with Bayesian inference-based finite Gaussian mixture models,” AIChE Journal, 54(7), 1811–1829, 2008.
- M. A. F. Figueiredo, A. K. Jain, “Unsupervised learning of finite mixture models,” IEEE Trans. PAMI, 24(3), 381–396, 2002.
- Z. Ge, Z. Song, “Multimode process monitoring based on Bayesian method,” Journal of Chemometrics, 23(12), 636–650, 2009.
- Q. Jiang, X. Yan, et al., “Multimode Process Monitoring Using Variational Bayesian Inference and Canonical Correlation Analysis,” IEEE Trans. Automation Science and Engineering, 16(4), 1814–1824, 2019.