ラテン超方格サンプリング（LHS）の実装（Python）

研究でラテン超方格サンプリング（LHS）を実装する機会があったのでまとめてみました。
間違いやご提案があればご連絡ください。

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

実験計画法の一つでできるだけ探索空間から均一にサンプリングする手法の一つ
単純なランダムサンプリングよりも効率的にサンプリングできるらしい。

具体的には、サンプリング点数で探索空間を分割する。
例えば2次元空間から10点をサンプリングするなら、10×10の格子を生成する。
そこから、行と列が一致しないように格子を10個選び、その中でランダムサンプリングする。
以下は、2次元空間から5点サンプリングするときの例。

実装

ここでは二次元空間からのサンプリングを例として実装を行う。

import numpy as np

n_dim = 2  # number of dimension
n_sample = 10  # number of sample

lhs_idx = np.zeros([n_sample, n_dim])
for i in range(n_dim):
    lhs_idx[:, i] = np.random.permutation(grid_array)

# array for storing values
lhs_result = np.zeros([n_sample, n_dim])
for i in range(n_sample):
    for j in range(n_dim):
        lhs_result[i, j] = (np.random.rand() + lhs_idx[i, j])/n_sample

これはPythonのpyDOEライブラリを使うと、

import pyDOE2 as doe

n_dim = 2  # number of dimension
n_sample = 10  # number of sample
lhs = doe.lhs(n_dim, samples=n_sample )

となる。

追加でサンプリングしたいとき

最初にn点サンプリングした後にm点サンプリングに変更したい場合について考える（n < m）。
この場合は、m×mでグリッドを生成し、n点サンプリングで追加されているグリッド以外から行と列が一致しないようにグリッドを選んでランダムサンプリングする。

下に、先ほど示した5点から追加で5点サンプリングする場合を示している。
この場合は下の青い領域から行と列が一致しないように格子を決めればいい。

# 'init_array' is initial lhs array.
old_lhs_idx = np.floor(init_array / (1/n_sample)).astype('int32')
old_n_sample = init_array.shape[0]

new_n_sample = n_sample-old_n_sample
new_lhs_idx = np.zeros([new_n_sample, n_dim])

for i in range(n_dim):
    # remove the old index and permutation
    new_lhs_idx[:, i] = np.random.permutation(np.delete(grid_array, old_lhs_idx[:, i], 0))

new_lhs_result = np.zeros([new_n_sample, n_dim])
for i in range(new_n_sample):
    for j in range(n_dim):
        new_lhs_result[i, j] = (np.random.rand() + new_lhs_idx[i, j])/n_sample

lhs_idx = np.concatenate([old_lhs_idx, new_lhs_idx], axis=0)
lhs_result = np.concatenate([init_array, new_lhs_result], axis=0)

まとめたコードはGithubにあげときます。
こちらから