ラテン超方格サンプリング（LHS）

理論と物理

LHSの基本原理

先生、ラテン超方格サンプリングって普通のランダムサンプリングより何がいいんですか？ 名前だけ聞くとめちゃくちゃ難しそうなんですけど…

名前は仰々しいけど、やっていることはシンプルだよ。各入力変数の値域を均等に $N$ 分割して、それぞれの区間から 1点ずつ 選ぶ。ランダムサンプリングだと偶然ある領域に点が偏る——いわゆるクランピング——が起きるけど、LHSはそれを構造的に防ぐんだ。

へぇ、均等に区切るだけなのに効果があるんですか？

それが大きいんだ。同じサンプル数でモンテカルロ法と比べると、分散が30〜50%も小さくなることが多い。つまり少ないシミュレーション回数で精度の良い統計情報が得られる。CAEでは1回の解析に数時間かかることもあるから、サンプル数を減らせるのは計算コスト的に非常に大きいよ。

なるほど！ 例えばどういう場面で使うんですか？

典型的なのは3つだね。

• 不確かさ伝播（UQ）: 材料物性のバラツキが応力分布にどう影響するか
• 感度分析: どの入力パラメータが出力に一番効くかを特定する
• DOE（実験計画法）の初期点配置: サロゲートモデル（Kriging等）を作るときの学習データ点を決める

例えば自動車のクラッシュ解析で、板厚・降伏応力・摩擦係数の3変数にばらつきがあるとき、LHSで50点のサンプルを配置すれば、頭部傷害値（HIC）の分布をモンテカルロの100点分に匹敵する精度で推定できたりする。

層化サンプリングの数理

「均等に分割する」っていうのを数式で書くとどうなるんですか？

$d$ 次元の入力変数 $\mathbf{X} = (X_1, X_2, \ldots, X_d)$ に対して、$N$ 個のサンプルを生成することを考えよう。各変数 $X_j$ の累積分布関数（CDF）を $F_j$ とする。LHSの核心は、各変数の周辺分布を等確率で $N$ 層に分割する ことだ。

$j$ 番目の変数の $k$ 番目の層（$k = 1, 2, \ldots, N$）は以下の区間で定義される：

各層から1点ずつサンプリングするので、任意のサンプル点が第 $k$ 層に入る確率は：

具体的には、層内の位置を一様乱数 $U_{j,k} \sim \text{Uniform}(0,1)$ で決め、$k$ 番目のサンプル値を以下で生成する：

ここで $\pi_j$ は $\{1, 2, \ldots, N\}$ のランダム置換（permutation）であり、変数ごとに独立に生成する。この置換によって各変数の層の組み合わせがランダムになる。

つまり「均等な棚」を作ってから各棚の中でランダムに点を選ぶ、みたいなイメージですか？

そのイメージで完璧だ。2次元の場合を考えると分かりやすい。$N \times N$ のチェス盤を思い浮かべて、各行・各列に ちょうど1つ だけ点を置く——これがラテン方格の条件だ。「ラテン方格」という名前は、数学者レオンハルト・オイラーが研究した「ラテン方陣（Latin Square）」から来ている。各行・各列に同じ記号が1回だけ現れる配列のことだよ。

分散低減効果の証明

LHSだとなぜ分散が下がるのか、数学的に説明できますか？

出力 $Y = g(\mathbf{X})$ の期待値推定量 $\bar{Y}_N = \frac{1}{N}\sum_{i=1}^{N} g(\mathbf{x}^{(i)})$ の分散を比較しよう。McKay, Beckman & Conover (1979) の古典的結果によると：

第1項 $\sigma^2 / N$ は単純モンテカルロの分散そのもの。第2項がLHSによる分散低減量であり、$g$ が各入力変数に対して単調性を持つ（= 加法的成分が大きい）ほどこの項が大きくなる。つまり：

等号は $g$ が全入力変数に対して完全に非加法的（純粋な交互作用のみ）な場合にしか成り立たない。実際のCAE問題では加法的成分がほぼ必ず存在するので、LHSは常にモンテカルロより有利だ。

なるほど、入出力の関係が「素直」なほどLHSの恩恵が大きいってことですね。

その通り。構造解析の応力応答なんかは変数に対してだいたい単調だから、分散低減効果がものすごく効く。一方、カオス的な系（乱流の瞬時値とか）だと恩恵は薄くなることもある。ただそれでもLHSがモンテカルロより悪くなることは理論上ないから、デフォルトでLHSを使っておけば間違いない。

空間充填基準（Space-Filling Criteria）

LHSでサンプルを取ったけど、なんか偏ってる気がする…ってこともあるんですか？

あるよ。基本LHSは周辺分布の層化は保証するけど、多次元空間での均一性までは保証しない。例えば2変数で「右上と左下に点が集中して、真ん中がスカスカ」みたいなことが起きうる。そこで空間充填基準を使って配置を最適化するんだ。

代表的な空間充填基準は以下の3つ：

1. Maximin距離基準

全サンプル点ペア間の最小距離を最大化する：

ここで $\| \cdot \|_p$ は $L_p$ ノルム（通常は $p=2$ のユークリッド距離）。直感的には「点同士をできるだけ離す」配置を求める。

2. Maximin距離の $\phi_p$ 基準

Morris & Mitchell (1995) の連続的緩和：

$p$ を大きくすると maximin 基準に漸近する。$p = 15 \sim 50$ が実用的。微分可能なので勾配ベースの最適化が使える利点がある。

3. ディスクレパンシー基準

点集合の一様分布からの乖離度を測る：

このスターディスクレパンシー $D_N^{*}$ が小さいほど一様分布に近い配置になる。Koksma-Hlawka不等式により、数値積分の誤差上界と直結する。

maximin距離とディスクレパンシー、実務ではどっちを使えばいいんですか？

サロゲートモデル（Kriging, RBF等）を作るならmaximin距離が鉄板だ。予測精度が予測点と最近傍学習点の距離に依存するから。一方、数値積分（期待値の推定）が目的ならディスクレパンシーのほうが理論的に正しい。迷ったらmaximin距離を選べば大きく外さないよ。

相関制御（Correlation Control）

LHSで生成した変数同士に意図しない相関が出ることってあるんですか？

鋭い質問だね。各変数に独立にランダム置換を適用するから、偶然に変数間に疑似相関（spurious correlation）が生じることがある。$N$ が小さいと結構大きな相関が出ることもあるんだ。これを防ぐのがIman & Conover (1982) の相関制御法だ。

アルゴリズムの核心は以下の通り：

1. 目標相関行列 $\mathbf{C}_{\text{target}}$ を設定する（独立なら単位行列 $\mathbf{I}$）
2. LHSのランク行列 $\mathbf{R}$ に対して、現在のランク相関行列 $\mathbf{C}_R$ を計算する
3. Cholesky分解 $\mathbf{C}_{\text{target}} = \mathbf{L}\mathbf{L}^T$ を行う
4. 補正行列 $\mathbf{Q} = \mathbf{L} \cdot \mathbf{C}_R^{-1/2}$ を用いてランクを再配列する

補正後のランク相関行列は目標相関行列にほぼ一致する：

入力変数同士に相関がある場合（例えばヤング率とポアソン比が相関するとか）にも使えるんですね！

その通り。材料物性の相関を反映させたサンプリングができるのは実務的にすごく重要だ。Nataf変換やCopulaと組み合わせて使うことも多いよ。

アルゴリズムと実装

基本LHS生成アルゴリズム

LHSを実際にプログラムで生成するには、具体的にどんな手順を踏めばいいんですか？

シンプルだよ。$N$ サンプル、$d$ 変数の場合：

1. 各変数 $j = 1, \ldots, d$ に対して $\{1, 2, \ldots, N\}$ のランダム置換 $\pi_j$ を生成する
2. 各サンプル $i$ の各変数 $j$ について一様乱数 $u_{ij} \sim U(0,1)$ を生成する
3. 単位超立方体上の座標を計算する：$s_{ij} = \frac{\pi_j(i) - 1 + u_{ij}}{N}$
4. 逆CDF変換で目的の分布に写す：$x_{ij} = F_j^{-1}(s_{ij})$

ステップ3がLHSの核で、$\pi_j(i) - 1$ で層を決め、$u_{ij}$ で層内の位置をランダムに選んでいる。

思ったよりシンプルですね！ 計算量はどのくらいですか？

基本LHSの生成は $O(Nd)$ で、モンテカルロと同じだ。最適化LHS（maximin距離基準）は最適化ループが加わるので $O(N^2 d \cdot T)$（$T$ は反復回数）になるけど、サンプル生成のコストはCAE解析1回分のコストに比べれば無視できるレベルだよ。

最適化LHS

maximin距離を最大化するにはどうやるんですか？ 全探索は無理ですよね？

全探索は $(N!)^d$ の組み合わせがあるから当然無理だ。実用的なのは以下の方法：

• ESE法（Enhanced Stochastic Evolutionary）: 2行をランダムに入れ替えて改善する。Jin et al. (2005) の手法で実装が簡単
• シミュレーテッドアニーリング: Morris & Mitchell (1995) が提案。局所解からの脱出が可能
• 遺伝的アルゴリズム: 大規模問題に有効だが計算コストが高い
• 列ペアワイズ交換: Viana et al. (2009) の高速アルゴリズム。2変数ずつ最適化

実務ではESE法かシミュレーテッドアニーリングが定番だ。pyDOEやSciPyにも実装されているよ。

Pythonによる実装例

Pythonだとどう書くんですか？

SciPyの qmc.LatinHypercube を使うのが今の標準だ。最適化LHS（maximin基準）もワンライナーで生成できる。

import numpy as np
from scipy.stats.qmc import LatinHypercube
from scipy.stats import norm, uniform

# --- 基本LHS（5変数, 50サンプル）---
d, N = 5, 50
sampler = LatinHypercube(d=d, seed=42, optimization="random-cd")
unit_samples = sampler.random(n=N)  # [0,1]^d 上のLHS

# --- 逆CDF変換で目的の分布に写す ---
# 変数1: 正規分布 N(210e3, 5e3)  ← ヤング率 [MPa]
# 変数2: 正規分布 N(0.3, 0.02)   ← ポアソン比
# 変数3: 一様分布 U(1.5, 2.5)    ← 板厚 [mm]
x1 = norm.ppf(unit_samples[:, 0], loc=210e3, scale=5e3)
x2 = norm.ppf(unit_samples[:, 1], loc=0.3, scale=0.02)
x3 = uniform.ppf(unit_samples[:, 2], loc=1.5, scale=1.0)
# ...以下同様

optimization="random-cd" って何ですか？

Centered Discrepancy（中心ディスクレパンシー）を最小化する最適化オプションだ。他に "Lloyd" もある。デフォルト（None）だと最適化なしの基本LHSになる。実務では必ず最適化オプションを入れること。素のLHSだと空間充填性が悪いことがあるからね。

収束性とサンプル数の決定

LHSの収束速度ってモンテカルロより速いんですか？

オーダーとしてはどちらも $O(1/\sqrt{N})$ なんだけど、定数項が違う。LHSのほうが定数が小さいから、同じ誤差水準に達するのに必要なサンプル数が少ない。具体的な目安を表にまとめよう。

実践ガイド

CAE解析でのLHS適用フロー

実際のCAEプロジェクトでLHSを使うときの全体の流れを教えてください。

典型的なフローは5ステップだ：

1. 入力変数の選定: どのパラメータを変動させるかを決める（材料物性、形状寸法、荷重条件など）。感度が低い変数はスクリーニング（Morris法等）で事前に除外すると効率が上がる
2. 確率分布の設定: 各変数の分布型（正規、対数正規、一様、Weibull等）とパラメータを決める。実測データがあれば分布フィッティング、なければ専門家判断
3. LHSサンプルの生成: 最適化LHS（maximin基準推奨）で $N$ 点を生成。相関がある変数ペアにはIman-Conover法を適用
4. CAE解析の実行: 各サンプル点に対してソルバーを実行。パラメトリックスタディの自動化スクリプトを組む
5. 後処理: 出力の統計量（平均、標準偏差、パーセンタイル）を計算。必要に応じてサロゲートモデルを構築し追加解析

例えば自動車の部品設計だとどんな変数を入れるんですか？

よくあるのは：

• 板厚（公差 ±0.1mm → 正規分布 $N(t_0, 0.033)$ で3σ範囲）
• 降伏応力（ロットばらつき → 対数正規分布）
• スポット溶接の位置（±2mm → 一様分布）
• 摩擦係数（0.1〜0.3 → 一様分布）

変数数が4〜8個なら $N = 50 \sim 100$ でかなり良い精度が出る。衝突解析1ケースに4時間かかるとしても、50ケースなら200時間。クラスタで10並列すれば20時間で終わる。モンテカルロだと同精度に200〜300ケース必要になるから、コスト差は歴然だ。

サンプル数の決め方

「$10d$ 以上」とは言いますけど、もっと実践的な判断基準はありますか？

以下のアプローチがおすすめだ：

• ブートストラップ法: 現在の $N$ 点からリサンプリングして統計量のばらつきを評価。目標精度に達していれば打ち切り
• 逐次追加: まず $N_0 = 10d$ で開始し、サロゲートモデルの交差検証誤差が閾値以下になるまで点を追加（SLHS — 後述）
• 計算予算ベース: 「1ケースあたり計算時間 × サンプル数 ÷ 並列数 ≤ 許容日数」で $N$ を逆算する。現実的にはこれが最も多い

よくある失敗と対策

LHSを使って失敗するパターンってありますか？

めちゃくちゃあるよ。現場で多いのは：

分布のテールからサンプルされると解析が壊れるのは盲点でした…！

LHSは層化により確実にテール領域からもサンプルを引くから、純粋ランダムよりこの問題が起きやすい。正規分布なら区間 $[\mu - 3\sigma, \mu + 3\sigma]$ でトランケートするのが定石だ。対数正規分布なら0.1%点〜99.9%点で切る。

ソフトウェア比較

商用・OSSツール比較

LHSを使えるツールにはどんなものがありますか？

商用・OSSの両方に豊富な選択肢があるよ。

無料で始めるならどれがおすすめですか？

PythonならSciPyの qmc.LatinHypercube が今のベストチョイスだ。NumPy/SciPyだけで完結するし、最適化オプションもある。本格的なUQフレームワークが欲しければDAKOTA（Sandia国立研究所が開発）が業界標準。GUIはないけど、コマンドラインで大規模なパラメトリックスタディを回すのに最適だよ。

先端技術

適応的LHS（Adaptive LHS）

最近の研究でLHSはどう進化しているんですか？

大きなトレンドは適応的サンプリングだ。従来のLHSは事前にサンプル点を全部決めてから解析を回す「ワンショット設計」だったけど、適応的LHSでは解析結果を見ながら次のサンプル点を追加する。具体的には、サロゲートモデルの予測不確かさが大きい領域に重点的に点を足すんだ。

• EGO（Efficient Global Optimization）との統合: Krigingの予測分散が大きい点にLHS点を追加
• EIVAR（Expected Improvement of Variance）: 分散低減量を最大化する位置に追加
• Multi-fidelity LHS: 低精度モデル（粗メッシュ）と高精度モデル（密メッシュ）で異なるサンプル数を使い分ける

逐次LHS（SLHS: Sequential Latin Hypercube Sampling）

さっき出た「追加サンプルでLHS性が崩れる問題」、SLHSってどう解決するんですか？

SLHSのアイデアは、$N_1$ 点の初期LHSを $N_1 + N_2$ 点のLHSの「部分集合」になるようにあらかじめ設計しておくことだ。Sheikholeslami & Razavi (2017) のPLHS（Progressive LHS）が代表的で、初期配置の段階で将来の拡張を見据えた層分割をしておく。こうすれば追加しても全体のLHS性が保たれる。

高次元問題への拡張

変数が数十〜数百になるとLHSでも限界がありますか？

高次元（$d > 20$）だと空間充填の最適化が困難になるし、$10d$ の法則に従うとサンプル数が膨大になる。対策としては：

• 次元削減との併用: 活性部分空間（Active Subspace）やPCAで次元を落としてからLHSを適用する
• スクリーニングとの2段階設計: Morris法やグループSobol法で重要変数を絞り込んでからLHSを実行する
• スパースグリッドとのハイブリッド: 低次元方向にスパースグリッド、高次元にLHSを組み合わせる
• 準モンテカルロ法（QMC）との比較検討: Sobol列やHalton列は高次元でもディスクレパンシーが低い。$d > 50$ ではQMCがLHSより有利な場合もある

QMCとLHS、どっちが「正解」なんですか？

正解はケースバイケースだ。$d \leq 20$ で入力の確率分布が明確ならLHSがほぼ最良。$d$ が大きくて確率分布の情報が乏しい場合はQMCが効率的なことがある。最近のDAKOTAやUQLabはどちらもサポートしているから、両方試して比べるのが実務的だよ。

トラブルシューティング

よくあるトラブルと対策

LHSを使っているときに困ったら見返せるチェックリストみたいなの、ありませんか？

まとめておこう。

これは保存版ですね。最後に、LHSを一言で要約するとどうなりますか？

「入力空間を均等にカバーする最も実用的なサンプリング戦略」だ。モンテカルロより常に効率的で、実装も簡単。CAEの不確かさ定量化・DOE・サロゲートモデル構築のどれにおいても第一選択肢になる。理論で迷ったらこの記事に戻ってきてくれ。