拉丁超立方采样（LHS）
理论与物理
LHS的基本原理
老师，拉丁超立方采样比普通随机采样好在哪里？光听名字感觉超难的样子…
名字听起来唬人，但做的事情很简单。将每个输入变量的取值范围均匀分割成 $N$ 份，从每个区间里 各选1点。随机采样可能会偶然在某些区域聚集点——也就是所谓的聚类（clumping）——而LHS从结构上防止了这种情况。
诶，只是均匀分割就有效果吗？
效果很大。与相同样本数的蒙特卡洛法相比，方差通常能减小30%~50%。这意味着可以用更少的仿真次数获得精度更好的统计信息。在CAE中，一次分析有时需要数小时，因此能减少样本数在计算成本上意义重大。
原来如此！那具体在什么场景下使用呢？
典型的有三种。
- 不确定性传播（UQ）: 材料物性的波动对应力分布有何影响
- 敏感性分析: 确定哪个输入参数对输出影响最大
- DOE（实验设计法）的初始点配置: 构建代理模型（如Kriging）时确定学习数据点
例如在汽车碰撞分析中，当板厚、屈服应力、摩擦系数这三个变量存在波动时，用LHS配置50个样本点，就能以媲美蒙特卡洛100个样本点的精度来估计头部伤害值（HIC）的分布。
分层采样的数学原理
“均匀分割”用数学公式怎么写？
考虑对 $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$）由以下区间定义：
由于从每层采样一个点，任意样本点落入第 $k$ 层的概率为：
具体来说，用均匀随机数 $U_{j,k} \sim \text{Uniform}(0,1)$ 决定层内的位置，第 $k$ 个样本值由下式生成：
这里 $\pi_j$ 是 $\{1, 2, \ldots, N\}$ 的随机置换（permutation），每个变量独立生成。这个置换使得各变量的层组合变得随机。
也就是说，像是先做“均匀的架子”，然后在每个架子里随机选点，这样的感觉吗？
这个理解完全正确。考虑二维情况就容易明白了。想象一个 $N \times N$ 的棋盘，在每一行、每一列放置恰好一个点——这就是拉丁方格的条件。“拉丁方格”这个名字来源于数学家莱昂哈德·欧拉研究的“拉丁方阵（Latin Square）”，指的是每行每列中每个符号只出现一次的阵列。
方差缩减效果的证明
为什么LHS能降低方差，能从数学上解释吗？
我们来比较输出 $Y = g(\mathbf{X})$ 的期望值估计量 $\bar{Y}_N = \frac{1}{N}\sum_{i=1}^{N} g(\mathbf{x}^{(i)})$ 的方差。根据McKay, Beckman & Conover (1979) 的经典结果：
第一项 $\sigma^2 / N$ 是简单蒙特卡洛的方差本身。第二项是LHS带来的方差缩减量，当 $g$ 对每个输入变量具有单调性（即加性成分较大）时，此项会变大。也就是说：
等号仅在 $g$ 对所有输入变量完全非加性（纯交互作用）时才成立。实际的CAE问题中几乎总是存在加性成分，因此LHS总是比蒙特卡洛有利。
原来如此，输入输出关系越“顺滑”，LHS的收益就越大对吧。
没错。像结构分析的应力响应等，对变量大致是单调的，所以方差缩减效果非常显著。另一方面，对于混沌系统（如湍流的瞬时值等），收益可能会变薄。但即便如此，理论上LHS绝不会比蒙特卡洛更差，所以默认使用LHS是稳妥的。
空间填充准则（Space-Filling Criteria）
用LHS采样后，会不会感觉有点偏呢？
会的。基本LHS只保证边缘分布的分层，不保证多维空间中的均匀性。例如在二维情况下，可能出现“点集中在右上和左下，中间稀疏”的情况。因此需要使用空间填充准则来优化配置。
代表性的空间填充准则有以下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) 的相关性控制法。
算法的核心如下：
- 设定目标相关矩阵 $\mathbf{C}_{\text{target}}$（若独立则为单位矩阵 $\mathbf{I}$）
- 对LHS的秩矩阵 $\mathbf{R}$，计算当前的秩相关矩阵 $\mathbf{C}_R$
- 进行Cholesky分解 $\mathbf{C}_{\text{target}} = \mathbf{L}\mathbf{L}^T$
- 使用修正矩阵 $\mathbf{Q} = \mathbf{L} \cdot \mathbf{C}_R^{-1/2}$ 重新排列秩
修正后的秩相关矩阵几乎与目标相关矩阵一致：
当输入变量之间存在相关性时（例如杨氏模量和泊松比相关），也能用这个方法呢！
没错。能够进行反映材料物性相关性的采样，这在实践中非常重要。也常与Nataf变换或Copula结合使用。
算法与实现
基本LHS生成算法
实际用程序生成LHS，具体需要什么步骤？
很简单。对于 $N$ 个样本，$d$ 个变量的情况：
- 对每个变量 $j = 1, \ldots, d$，生成 $\{1, 2, \ldots, N\}$ 的随机置换 $\pi_j$
- 对每个样本 $i$ 的每个变量 $j$，生成均匀随机数 $u_{ij} \sim U(0,1)$
- 计算单位超立方体上的坐标：$s_{ij} = \frac{\pi_j(i) - 1 + u_{ij}}{N}$
- 通过逆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分析一次的成本相比，可以忽略不计。
优化LHS
要最大化maximin距离该怎么做？穷举搜索不可能吧？
穷举搜索有 $(N!)^d$ 种组合，当然不可能。实用的方法如下：
- ESE法（增强随机进化法）: 随机交换两行以改进。Jin et al. (2005) 的方法，实现简单
- 模拟退火: Morris & Mitchell (1995) 提出。能够跳出局部解
- 遗传算法: 对大规模问题有效，但计算成本高
- 列对交换: Viana et al. (2009) 的快速算法。每次优化两个变量
实践中，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的常数更小，因此达到相同误差水平所需的样本数更少。具体参考下表。
| 用途 | 维度 $d$ | 推荐样本数 $N$ | 备注 |
|---|---|---|---|
| 均值・标准差的估计 | 任意 | $N \geq 10d$ | 最低线 |
| 概率分布估计（直方图） | $\leq 10$ | $N \geq 50d$ | 确保尾部精度 |
| 代理模型构建（Kriging） | $\leq 20$ | $N = 10d \sim 50d$ | 使用空间填充准则优化 |
| Sobol敏感性指标估计 | 任意 | $N \geq 100(d+1)$ | 与Saltelli法并用 |
实践指南
CAE分析中的LHS应用流程
请告诉我实际CAE项目中使用LHS时的整体流程。
典型的流程有5个步骤：
- 输入变量的选定: 决定哪些参数需要变动（材料物性、形状尺寸、载荷条件等）。敏感性低的变量通过筛选（
この記事の評価ご回答ありがとうございます！参考に
なったもっと
詳しく誤りを
報告