[論文レビュー] An efficient algorithm for sampling from $\sin^k(x)$ for generating random correlation matrices

Enes Makalic, Daniel F. Schmidt|arXiv (Cornell University)|Sep 14, 2018

Cellular Automata and Applications参考文献 1被引用数 4

ひとこと要約

この論文は、区間 $(0, \pi)$ における $\sin^k(x)$ の比例分布から乱数変量を生成する、非常に効率的な拒否採択サンプリングアルゴリズムを提示する。この手法では、包み込み密度としてベータ分布を用いる。すべての $k \geq 1$ に対して、平均して最大で $\pi/(2\sqrt{2}) \approx 1.11$ 回の反復でサンプルが得られ、高次元の乱数相関行列の生成が超球座標パrameterization を用いて高速に可能になる。

ABSTRACT

In this note, we develop a novel algorithm for generating random numbers from a distribution with a probability density function proportional to $\sin^k(x)$, $x \in (0,π)$ and $k \geq 1$. Our algorithm is highly efficient and is based on rejection sampling where the envelope distribution is an appropriately chosen beta distribution. An example application illustrating how the new algorithm can be used to generate random correlation matrices is discussed.

研究の動機と目的

• $(0, \pi)$ における $\sin^k(x)$ の比例する確率密度関数からサンプリングする計算的に効率的な手法の開発。この分布は、乱数相関行列の生成に現れる。
• 逆変換サンプリングの非効率性を克服すること。この手法は $O(k)$ の計算量を要し、$k$ が大きくなると実用的でなくなる。
• 大規模な $k$ に対して精度の問題を抱える数値積分の代替として、数値的に安定した手法を提供すること。
• 必要な角度 $\theta_{ij}$ を効率的にサンプリングすることで、高次元相関行列のスケーラブルな生成を可能とすること。
• R や MATLAB といった広く使われているプラットフォームでアルゴリズムを実装・ベンチマークし、実用的な利用を可能とすること。

提案手法

• アルゴリズムは、対称的な $\text{Beta}(k+1, k+1)$ 分布を包み込み密度として用いる拒否採択サンプリングを採用する。
• 受理基準は、比率 $f_X(x)/g_Y(x)$ に基づき、$\frac{1}{k}\log U \leq \log\left(\frac{\pi^2 \sin X}{4X(\pi - X)}\right)$ を満たす場合に受理する。
• 上界 $M_k = \frac{f_X(\pi/2)}{g_Y(\pi/2)}$ は解析的に導出され、包み込み条件 $f_X(x) \leq M_k g_Y(x)$ を保証する。
• 正規化定数はガンマ関数およびベータ関数を用いて表現され、正確な計算が可能になる。
• パウラハミディとワンの超球座標パrameterization を用いて、$p \times p$ 相関行列の生成にこの手法を適用することで検証された。
• R（パッケージ randcorr を介して）および MATLAB（File ID: 68810）の両方で実装され、一般公開されている。

実験結果

リサーチクエスチョン

• RQ1 すべての $k \geq 1$ に対して、$(0, \pi)$ における $\sin^k(x)$ から効率的にサンプリングできる拒否採択アプローチを設計できるか？
• RQ2 受理率の上限をタイトに保ちつつ、高速なサンプリングを可能にする包み込み分布は何か？
• RQ3 すべての $k \geq 1$ に対して、理論的な最大期待反復回数（1サンプルあたり）は何か？
• RQ4 高次元相関行列に対して、提案手法は逆変換サンプリングと比較して性能に優れているか？
• RQ5 $p$ が最大 1000 に達するまで、このアルゴリズムは大規模な $p \times p$ 相関行列の生成に効率的にスケーリングできるか？

主な発見

• 1サンプルあたりの最大期待反復回数は $\pi/(2\sqrt{2}) \approx 1.11$ であり、$k \to \infty$ の極限で達成される。
• $k=1$ の場合、期待反復回数は $\pi/3 \approx 1.047$ で、$k$ が増加するにつれて単調に増加する。
• このアルゴリズムは、一般のラップトップで $100 \times 100$ の乱数相関行列を約 0.02 秒で生成できる。
• $1000 \times 1000$ の相関行列は約 0.5 秒で生成され、行列サイズに比例して線形スケーリングしていることが示された。
• 提案手法は、逆変換サンプリングを上回る性能を示し、同等のハードウェアで $100 \times 100$ 行列には 0.9 秒、$600 \times 600$ 行列には 42 秒を要した。
• 拒否サンプラーは数値的に安定しており、大規模な $k$ に対して数値積分の精度問題を回避する。