[논문 리뷰] An efficient algorithm for sampling from $\sin^k(x)$ for generating random correlation matrices

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

Cellular Automata and Applications참고 문헌 1인용 수 4

한 줄 요약

이 논문은 (0, π)에서 $\sin^k(x)$ 비례하는 분포에서 랜덤 변수를 생성하기 위한 매우 효율적인 거부 샘플링 알고리즘을 제시한다. 이는 베타 분포를 엔벨로프로 사용한다. 이 방법은 모든 $k \geq 1$ 에서 최대 평균 수용률이 $\pi/(2\sqrt{2}) \approx 1.11$ 번의 반복을 필요로 하여, 초구면 매개변수화를 통해 고차원 랜덤 상관계수 행렬을 빠르게 생성할 수 있다.

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)$ 비례하는 확률 밀도 함수에서 계산적으로 효율적인 샘플링 방법을 개발하는 것.
• 역변환 샘플링의 비효율성 — 특히 큰 $k$ 에서는 $O(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)$ 를 보장한다.
• $f_X(x)$ 와 $g_Y(x)$ 의 정규화 상수는 감마 함수와 베타 함수를 사용하여 표현되어 정확한 계산이 가능하다.
• 알고리즘은 Pourahmadi와 Wang의 초구면 매개변수화를 사용하여 $p \times p$ 상관계수 행렬을 생성하는 데 적용되어 검증되었다.
• 공개 사용을 위해 R(패키지 randcorr를 통해)과 MATLAB(File ID: 68810) 모두에서 알고리즘을 구현하였다.

실험 결과

연구 질문

• RQ1 모든 $k \geq 1$ 에 대해 $(0, \pi)$에서 $\sin^k(x)$ 에서 효율적인 거부 샘플링 방법을 설계할 수 있는가?
• RQ2 어떤 엔벨로프 분포가 수용 비율에 대한 날카운 경계를 제공하면서도 신속한 샘플링을 가능하게 하는가?
• RQ3 모든 $k \geq 1$ 에서 이론적으로 도달할 수 있는 최대 평균 반복 수는 얼마인가?
• RQ4 고차원 상관계수 행렬에서 제안된 방법이 역변환 샘플링보다 성능 면에서 어떻게 비교되는가?
• RQ5 큰 $p \times p$ 상관계수 행렬을 $p$ 가 최대 1000까지 효율적으로 생성할 수 있는가?

주요 결과

• 최대 평균 반복 수는 $\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$ 에서 수치적 적분의 정밀도 문제를 피한다.