拉丁超立方抽样(Latin Hypercube Sampling, LHS)MATLAB实现

"超拉丁立方抽样"通常指拉丁超立方抽样(LHS) 或其高级变种(如正交拉丁超立方抽样)。LHS是一种高效的实验设计方法,通过分层抽样确保样本在空间中的均匀覆盖,适用于高维参数空间的抽样问题。以下是完整的MATLAB实现,包括基础LHS、优化变种(最大化最小距离、最小化相关性)、正交LHS及可视化分析。

一、基础拉丁超立方抽样(LHS)

1.1 核心原理

将每个维度的取值范围划分为n个等概率区间,在每个区间内随机抽取1个样本,再对各维度样本随机排列,确保样本点在空间中均匀分布。

1.2 基础实现代码

function samples = lhs_basic(n, d, bounds)
    % 基础拉丁超立方抽样
    % 输入:
    %   n: 样本数量
    %   d: 维度
    %   bounds: 变量边界 (d×2矩阵,每行[min, max])
    % 输出:
    %   samples: n×d样本矩阵
    
    % 处理边界(默认[0,1])
    if nargin < 3, bounds = zeros(d,2) + [0,1]; end
    if size(bounds,1)==1, bounds = repmat(bounds, d, 1); end
    
    samples = zeros(n, d);
    for i = 1:d
        % 1. 将区间[0,1]分为n等份
        intervals = linspace(0, 1, n+1);
        % 2. 每个区间随机取1点
        points = intervals(1:end-1) + diff(intervals)/2 + (rand(1,n)-0.5)*diff(intervals);
        % 3. 随机排列
        samples(:,i) = points(randperm(n))';
        % 4. 缩放到实际边界
        samples(:,i) = bounds(i,1) + samples(:,i)*(bounds(i,2)-bounds(i,1));
    end
end

二、优化拉丁超立方抽样(OLHS)

2.1 最大化最小距离(Maximin LHS)

通过优化样本点布局,最大化任意两点间的最小欧氏距离,提升空间填充性。

function samples = lhs_maximin(n, d, bounds, max_iter)
    % 最大化最小距离的LHS
    if nargin < 4, max_iter = 1000; end
    samples = lhs_basic(n, d, bounds);  % 初始样本
    [~, d_min] = min(pdist(samples));  % 初始最小距离
    
    for iter = 1:max_iter
        % 随机交换两个样本的某个维度值
        i1 = randi(n); i2 = randi(n); j = randi(d);
        temp = samples(i1,j);
        samples(i1,j) = samples(i2,j);
        samples(i2,j) = temp;
        
        % 计算新最小距离
        new_d_min = min(pdist(samples));
        if new_d_min > d_min
            d_min = new_d_min;  % 接受改进
        else
            % 恢复原样本
            samples(i2,j) = samples(i1,j);
            samples(i1,j) = temp;
        end
    end
end

2.2 最小化相关性(MinCorr LHS)

通过优化样本点布局,最小化各维度间的 Pearson 相关性,提升样本独立性。

function samples = lhs_mincorr(n, d, bounds, max_iter)
    % 最小化相关性的LHS
    if nargin < 4, max_iter = 1000; end
    samples = lhs_basic(n, d, bounds);
    corr_init = mean(abs(corr(samples)-eye(d)), 'all');  % 初始相关性
    
    for iter = 1:max_iter
        i1 = randi(n); i2 = randi(n); j = randi(d);
        temp = samples(i1,j);
        samples(i1,j) = samples(i2,j);
        samples(i2,j) = temp;
        
        new_corr = mean(abs(corr(samples)-eye(d)), 'all');
        if new_corr < corr_init
            corr_init = new_corr;
        else
            samples(i2,j) = samples(i1,j);
            samples(i1,j) = temp;
        end
    end
end

2.3 正交拉丁超立方抽样(OLHS)

在LHS基础上引入正交性约束,确保样本点在各维度上的投影均匀且相互独立,适用于高维空间的高效抽样。

function samples = lhs_orthogonal(n, d, bounds)
    % 正交拉丁超立方抽样(简化实现)
    % 基于Plackett-Burman设计的改进LHS
    if d < 2, error('维度至少为2'); end
    if n < d, error('样本数至少为维度'); end
    
    % 1. 生成基础LHS
    base_samples = lhs_basic(n, d, [0,1]);
    
    % 2. 正交化处理(Gram-Schmidt)
    ortho_samples = zeros(n, d);
    ortho_samples(:,1) = base_samples(:,1);
    for i = 2:d
        proj = ortho_samples(:,1:i-1) * (ortho_samples(:,1:i-1)' * base_samples(:,i)) / ...
               diag(ortho_samples(:,1:i-1)' * ortho_samples(:,1:i-1));
        ortho_samples(:,i) = base_samples(:,i) - proj;
    end
    
    % 3. 归一化并缩放到边界
    samples = ortho_samples ./ vecnorm(ortho_samples, 2, 2);  % 单位化
    for i = 1:d
        samples(:,i) = bounds(i,1) + (samples(:,i)+1)/2*(bounds(i,2)-bounds(i,1));
    end
end

三、可视化与评估

3.1 样本空间分布可视化

function plot_lhs(samples, bounds, title_str)
    % 可视化2D/3D LHS样本
    [n, d] = size(samples);
    figure;
    if d == 2
        scatter(samples(:,1), samples(:,2), 50, 'filled');
        xlabel('维度1'); ylabel('维度2');
        title(title_str);
        grid on; axis equal;
        xlim(bounds(1,:)); ylim(bounds(2,:));
    elseif d == 3
        scatter3(samples(:,1), samples(:,2), samples(:,3), 50, 'filled');
        xlabel('维度1'); ylabel('维度2'); zlabel('维度3');
        title(title_str);
        grid on; axis equal;
    else
        % 高维:显示前2维+相关矩阵
        subplot(1,2,1);
        scatter(samples(:,1), samples(:,2), 50, 'filled');
        xlabel('维度1'); ylabel('维度2'); title('前2维分布');
        subplot(1,2,2);
        imagesc(corr(samples)); colorbar; title('维度相关性');
    end
end

3.2 抽样质量评估

function evaluate_lhs(samples)
    % 评估LHS样本质量
    [n, d] = size(samples);
    metrics = struct();
    
    % 1. 最小距离(空间填充性)
    metrics.min_dist = min(pdist(samples));
    % 2. 平均距离
    metrics.mean_dist = mean(pdist(samples));
    % 3. 相关性(越小越好)
    metrics.corr = mean(abs(corr(samples)-eye(d)), 'all');
    % 4. 边际分布均匀性(K-S检验)
    metrics.ks_p = mean(cellfun(@(x) kstest(x), num2cell(samples,1)));
    
    disp('LHS质量评估:');
    disp(['最小距离: ', num2str(metrics.min_dist)]);
    disp(['平均距离: ', num2str(metrics.mean_dist)]);
    disp(['平均相关性: ', num2str(metrics.corr)]);
    disp(['K-S检验p值: ', num2str(metrics.ks_p)]);
end

参考代码 超拉丁立方抽样matlab-latin-sampling www.youwenfan.com/contentcnt/160563.html

四、应用示例

4.1 基础LHS使用

% 参数设置
n = 50;          % 样本数
d = 2;           % 维度
bounds = [0,10;   % 维度1范围
          0,5];   % 维度2范围

% 生成LHS样本
samples = lhs_basic(n, d, bounds);

% 可视化与评估
plot_lhs(samples, bounds, '基础LHS样本分布');
evaluate_lhs(samples);

4.2 优化LHS对比

% 生成不同LHS样本
samples_basic = lhs_basic(n, d, bounds);
samples_maximin = lhs_maximin(n, d, bounds, 1000);
samples_mincorr = lhs_mincorr(n, d, bounds, 1000);
samples_olhs = lhs_orthogonal(n, d, bounds);

% 对比可视化
figure;
subplot(2,2,1); plot_lhs(samples_basic, bounds, '基础LHS');
subplot(2,2,2); plot_lhs(samples_maximin, bounds, 'Maximin LHS');
subplot(2,2,3); plot_lhs(samples_mincorr, bounds, 'MinCorr LHS');
subplot(2,2,4); plot_lhs(samples_olhs, bounds, '正交LHS');

4.3 高维LHS应用(5维参数空间)

% 5维参数,100样本
n = 100; d = 5;
bounds = [0,1; 0,1; 0,1; 0,1; 0,1];  % 5维[0,1]边界
samples = lhs_orthogonal(n, d, bounds);

% 评估高维LHS
evaluate_lhs(samples);

五、关键优势与应用场景

  1. 高效性:用少量样本覆盖高维空间,优于随机抽样;
  2. 均匀性:通过分层抽样避免样本聚集;
  3. 灵活性:可结合优化准则(如Maximin、MinCorr)提升性能。

应用场景

  • 实验设计(DOE)、不确定性量化(UQ)、机器学习超参数调优、蒙特卡洛模拟等。
posted @ 2026-04-07 11:34  kang_ms  阅读(124)  评论(0)    收藏  举报