潜在クラス分析（LCA）｜観測指標から“見えないクラス”を推定（混合モデル）

2026年10月8日にscikit-learnのGaussianMixture、EM反復、共分散型、AIC/BICによるモデル選択本文を再確認した。潜在クラスは観測された回答から確率的な所属を推定するモデルで、クラス数、初期値、局所解、共分散、局所独立の仮定に依存する。所属確率やラベルを顧客の固定的な人格とみなさず、適合度、安定性、再現、施策の不利益を別々に確認する。

アンケートをいくつかの顧客タイプに分けたいとき、分類の数値が出たことと、施策に使える集団が見つかったことは別です。潜在クラス分析（LCA）は回答パターンを確率で表す方法であり、個人の本質や診断名を確定する道具ではありません。

30秒で要点

• 基本のLCAは、クラス内では各回答が独立になると仮定した、カテゴリ分布の混合モデルです。
• クラス別の回答確率、混合比、個人の事後所属確率を分けて読みます。最大の確率だけを残すと不確かさが消えます。
• 複数初期値で推定を確かめ、同じデータと尤度の下でBIC/AICを比較します。少数例の最小値を「真のタイプ数」とはしません。

最初に分けたい用語

表が収まらない場合は、横にスクロールできます

poLCAの論文抄録は、カテゴリ回答を条件付き独立の分布の混合として扱う基本モデルを説明しています。公式マニュアルの実装は名義カテゴリとして扱い、順序を自動利用しません。「満足・普通・不満」の順番をモデルに反映したい場合は、別の順序制約を持つ仕様かどうかを確認します。

同じ回答でも、確率を1人1ラベルに潰さない

二つのクラスA・Bを仮定した、説明専用の値を置きます。下の値は推定結果ではなく、著者が計算のために指定した架空のパラメータです。

表が収まらない場合は、横にスクロールできます

回答が(1, 1, 1)なら、混合比込みの確率はA側が0.6×0.8×0.7×0.6＝0.2016、B側が0.4×0.2×0.3×0.4＝0.0096です。和0.2112で割ると、Aへの事後確率は21/22、約95.45%になります。これは「q1の1がクラスAを引き起こす」という因果の計算ではありません。

回答(0, 1, 0)なら、A側は0.6×0.2×0.7×0.4＝0.0336、B側は0.4×0.8×0.3×0.6＝0.0576です。合計0.0912で割り、Aは7/19で約36.84%、Bは約63.16%。最大値だけならBですが、Aにも36.84%を残します。自動で別の価格や重要な取扱いを割り当てる前に、この不確かさと利用目的を検討します。

元の10人を残し、二値回答のモデルへ直す

元稿の教材にあったq1・q2・q3の全10行を使います。ここでも実データではなく架空の回答です。旧コードは0/1を連続値へ変換してGaussianMixtureを当てていましたが、同APIの公式説明が扱うのはガウス分布の密度です。カテゴリの確率質量を使うLCAとは異なるため、例の中心を二値のベルヌーイ混合へ変更します。

表が収まらない場合は、横にスクロールできます

1件の回答確率は、クラスごとに「混合比×各回答の確率の積」を求め、それを全クラスで足します。EMのE段階ではこの値から各人の所属確率を求め、M段階ではその重み付きの回答平均とクラス割合へ更新します。NumPy 2.1.3で動く小さな教材です。欠測、標準誤差、重み付き標本、一般の多値回答には対応していません。

import numpy as np

# 元の架空10人・3回答をそのまま保持。実顧客の観測値ではない。
X = np.array(list(zip(
    [0,1,1,0,1,0,1,1,0,0],
    [1,1,0,0,1,0,1,0,0,1],
    [0,1,0,0,1,0,1,1,0,1])), dtype=float)
n, j = X.shape
eps = 1e-6  # log(0)回避の境界制約。無制約の大域最適解の証明ではない。

def expectation(pi, theta):
    log_joint = (X @ np.log(theta).T
                 + (1-X) @ np.log1p(-theta).T + np.log(pi))
    shift = log_joint.max(axis=1, keepdims=True)
    weight = np.exp(log_joint-shift)
    z = weight / weight.sum(axis=1, keepdims=True)
    ll = np.sum(shift[:, 0] + np.log(weight.sum(axis=1)))
    return z, float(ll)

def fit(k, seed):
    rng = np.random.default_rng(seed)
    pi, theta = np.full(k, 1/k), rng.uniform(.15, .85, (k, j))
    previous = -np.inf
    for iteration in range(1, 20001):
        z, ll = expectation(pi, theta)
        assert ll >= previous-1e-8
        if ll-previous < 1e-9:
            return ll, pi, theta, z, iteration
        previous = ll
        nk = z.sum(axis=0)
        pi = nk/n
        theta = np.clip((z.T @ X)/nk[:, None], eps, 1-eps)
    raise RuntimeError('収束条件に未到達')

models = {}
for k in (1, 2):
    trials = [fit(k, seed) for seed in range(20)]
    best = max(trials, key=lambda x: x[0])
    models[k] = best
    ll, pi, theta, z, iteration = best
    parameters = (k-1) + k*j
    aic = -2*ll + 2*parameters
    bic = -2*ll + np.log(n)*parameters
    print(k, 'll/AIC/BIC:', np.round([ll, aic, bic], 6))
    print('start range:', np.round([min(t[0] for t in trials), ll], 6))
    print('mixing:', np.round(pi, 6), 'posterior sums:', z.sum(axis=1))
    assert np.allclose(z.sum(axis=1), 1)
# 3クラス11・4クラス15パラメータに対し、3二値回答の自由なセル確率は7。
# この10件で1〜4クラスの順位を「真のタイプ数」とは解釈しない。

この境界制約と停止条件での結果は次のとおりです。20初期値は6桁表示で同じ最良尤度へ達しましたが、全パラメータ空間の大域最適性を証明したわけではありません。

表が収まらない場合は、横にスクロールできます

Rの情報量規準の説明に沿えば、AICは「-2×対数尤度＋2×パラメータ数」、BICは後半を「log(人数)×パラメータ数」にしたものです。この例では1クラスの値が小さくなります。ただし10人しかおらず、2クラスの推定では一部の回答確率が設定した境界へ寄ります。通常の大標本近似を安心して使える状況と見なさず、数値比較の教材にとどめます。

3二値回答には8通りのセルがあり、自由なセル確率は合計1という制約を除いて7個です。Kクラスの名目パラメータ数は(K-1)＋3Kで、3クラスなら11、4クラスなら15。旧稿の1〜4クラス比較をそのまま実施して順位を採用する根拠にはなりません。2クラスでも名目数が7と一致するだけで、識別可能性や安定推定が保証されるわけではありません。連続密度のGMMと離散確率のLCAのBICを同じ順位表へ入れることも避けます。

実務ではクラス数より、仮定と再現性を点検する

表が収まらない場合は、横にスクロールできます

旧稿の「好み・頻度・不満を6〜10設問で尋ねる」は設計例として残せますが、必要項目数の万能な基準にはなりません。尺度、対象者、欠測、回答負担を決めてから項目を選びます。poLCAならカテゴリを1始まりの連番で渡し、欠測の扱いと複数初期値を明示します。mclustの公式概要もガウス混合を中心とした別のモデルであり、二値LCAと交換可能な名前として挙げません。本稿のPythonは専用パッケージによる実調査の推定を代行するものではありません。

ラベルの入替えを「クラスが不安定」と誤認しない

同じ4人を二つに分けた結果が、Aでは[0, 0, 1, 1]、Bでは[1, 1, 0, 0]だったとします。番号だけが逆で、分割は同じです。クラス0同士の集合をそのまま比べたJaccard係数（共通部分÷和集合）は0ですが、対応するクラスを合わせると1になります。

from sklearn.metrics import adjusted_rand_score

# 同一4人の分割。Bはクラス番号だけを逆にしたもの。
a = [0, 0, 1, 1]
b = [1, 1, 0, 0]
left = {i for i, label in enumerate(a) if label == 0}
wrong = {i for i, label in enumerate(b) if label == 0}
matched = {i for i, label in enumerate(b) if label == 1}
def jaccard(x, y):
    return len(x & y) / len(x | y)
print('unmatched/matched Jaccard:', jaccard(left, wrong), jaccard(left, matched))
print('ARI:', adjusted_rand_score(a, b))
assert jaccard(left, wrong) == 0
assert jaccard(left, matched) == 1
assert adjusted_rand_score(a, b) == 1

scikit-learn 1.5.2の調整ランド指数（ARI）は、このラベル置換に影響されず1を返します。別期間の「別の人」をそのまま同じ行番号で比較してはいけません。同じ対象への再割当なら対応を保ち、異なる対象群ならクラス別回答確率や構成比を比較するなど、問いを分けます。

三つの練習を、採用判断に結び付ける

1. BICは3クラスが最小、2クラスの方が説明しやすい。 まず同じデータ・尤度で収束と複数初期値を確認します。差、識別可能性、不確かさ、別標本での再現、運用目的を合わせて判断し、説明しやすさだけで統計上の問題を隠しません。
2. 1クラスが5%未満になった。 5%だけを削除基準にせず、実人数・境界推定・異常回答・局所独立の不適合を調べます。希少だが重要な回答を消す可能性もあります。この10人教材で5%という推定割合は、確定した半人の所属を意味しません。
3. 期間を変えると割当が違う。 同一人かを確認し、ラベル対応と所属確率をそろえてから比較します。そのうえで標本構成の変化、測定の変化、推定の局所解、実際の行動変化を切り分けます。分類の安定性と施策効果は別の問いです。