[확률] 랜덤 숫자를 생성하는 방법 (1) Uniform distribution

우리가 MATLAB에서 랜덤함수를 쓸 때 별 생각없이 쓴다. 코드 한 줄이면 원하는 사이즈의 랜덤한 숫자를 원하는 분포(ex. 가우시안 분포)에서 추출할 수 있다.

 

그런데 실제로 이 값들이 어떻게 생성되는지는 잘 모른다.

 

임의의 확률분포가 있을 때 그 확률분포를 따르는 난수를 생성하기 위해서는 기본적으로 uniform distribution이 필요하다.

그래서 처음에는 uniform distribution을 만드는 방법을 배우고 확장해나가는 흐름이다.

 

 

Uniform distribution

 

interval $(0,m)$에서 uniformly distributed integer $z_{i}$를 구하기 위한 알고리즘은 다음과 같다.

 

general한 알고리즘으로는

$z_n = f(z_{n-1},\ldots,z_{n-r}) \mod m$의 형태이다.

 

여기서 $z_{n}$은 n으로부터 n-r까지의 값에 영향을 받을 수 있다. 또한 mod 함수를 사용해서 $m$으로 나눠서 사용한다.

그리고 $f$는 알고리즘마다 다르다. 대표적인 알고리즘을 소개하면 Lehmer's algorithm이 있다.

 

Lehmer's algorithm

 

$z_{n}=az_{n-1} \mod m \quad z_{0}=1, n\geq 1$

사실 상 다음과 동일하다.

$z_{n}=a^{n}\mod m$

 

다음과 같은 알고리즘으로 구하며, 보통 사용할 때는 prime number m은 $10^9$정도 수준으로 사용하고 a값을 찾는다.

Lehmer 씨는 prime number m으로 $2^31-1$을 제안하였다. 또한 위의 integer $a$의 경우에는 primitive root of m을 쓰는데 이는 다음을 만족한다.

 

An integer a is called the primitive root of m if the smallest n such that $a^n = 1 \mod m \text{ is } n = m-1$

 

 

 

왜 이런 조건을 만족하는 값을 쓰는 것인지는 정확하게 나오지는 않았다.

 

실제로 위의 알고리즘이 잘 동작하는지 보기 위해 책에 나온대로

$a=7^5, m = 2^31-1$인 경우에 계산을 해서 분포를 구해보았다.

 

참고로 위의 $z_{i}$는 $(0,m)$ 사이의 정수값이기 때문에 m으로 나눠줌으로써 MATLAB에 있는 rand 함수와 비교할 수 있다.

 

10000개의 random variable을 추출했을 때 다음과 같이 균일하게 분포가 이루어짐을 알 수 있다.

아래는 위 그림을 그리기 위한 소스코드이다.

 

a = 7^5 ;
m = 2^31 - 1 ;

N = 10000 ;
rng(1) ;
test = rand(N,1) ;

figure ;
histogram(test) ;

compare = zeros(N,1) ;
compare(1,1) = 1 ;

for i = 2:N
    compare(i) = mod( compare(i-1) * a, m ) ;
end

compare = compare / m ;

figure ;
histogram(test) ; hold on;
histogram(compare) ;
legend('MATLAB', 'Algorithm') ;
title('Uniform distribution') ;

 

 

이와 같이 uniform distribution을 추출할 때 다음에는 다른 분포에 대해 구하는 방법을 알아볼 생각이다.