이번 포스팅에서는 모집단에서 유한한 갯수의 표본을 추출해서 평균값을 추정하는 상황에 대해 다뤄보겠습니다. 유한한 크기의 표준편차를 가지는 모집단이 있을 때, 표본평균의 값이 모집단의 평균에 비해 얼마나 차이가 날 수 있는지 살펴볼텐데요. 이를 결정하는 중요한 이론인 중심 극한 정리 (central limit theorem, 줄여서 CLT)에 대해서 알아봅시다.
먼저 평균, 분산 및 표준편차의 수학적인 정의를 알아둘 필요가 있습니다. 그리고 특정한 확률 분포를 따르는 모집단으로부터 유한한 갯수의 표본을 추출한 뒤 평균값을 구하면, 그게 모집단의 실제 평균과는 일반적으로 다르다는 사실도 인지할 필요가 있죠. 여기에 대해서는 다음 포스팅에 더 자세하게 소개되어 있습니다.
확률밀도함수 (probability density function)를 통해 표본평균의 값이 특정 범위 내에 있을 확률을 특정하고, 확률밀도함수의 형태를 알아내기 위해서 특성함수 (characteristic function)라는 것을 사용할텐데요. 이들의 수학적 정의와 평균, 표준편차와의 관계에 대한 더 자세한 내용은 다음 포스팅에 소개되어 있습니다.
뒤에서 더 자세히 알아보겠지만, 표본의 갯수가 매우 많은 경우에 표본평균의 값은 정규 분포 (normal distribution)라고 하는 확률 분포를 따르게 됩니다. 그리고 확률밀도함수는 가우스 (Gaussian) 함수의 형태로 주어지게 되죠. 가우스 함수가 생소하게 느껴지신다면 시작하기에 앞서서 다음 포스팅을 읽어보시면 큰 도움이 되리라 생각합니다.
CLT의 개념과 증명
중심 극한 정리가 시사하는 바를 간단히 말하자면, 표본의 갯수가 많아질수록 표본평균의 확률 분포가 정규 분포에 가까워진다는 것입니다. 이말인즉슨 모집단의 확률 분포가 정규 분포를 따르지 않더라도, 표준편차의 크기가 유한하고 표본의 갯수가 많다면 표본평균이 특정한 값을 가질 확률은 정규 분포를 따른다는 것입니다.
정규 분포를 따르는 변수의 확률밀도함수는 가우스 함수로 주어지고, 평균과 표준편차를 알면 그 형태를 완벽하게 특정할 수 있다는 특징이 있습니다. 표본평균의 확률밀도함수의 평균은 모집단의 평균과 동일합니다. 그리고 표본평균의 표준편차는 모집단의 표준편차에서 표본의 갯수의 제곱근을 나눈 것과 같죠. 여기서는 흔히 통용되는 방식에 따라 모집단의 평균과 표준편차를 각각 그리스 문자 뮤 (mu)와 시그마 (sigma)로 표기하겠습니다.
중심 극한 정리의 취지를 알았으니, 이를 어떻게 도출해낼 수 있는지 살펴봅시다. 표본평균의 특성함수를 통해서 증명할 수 있는데요. 특성함수는 확률밀도함수의 푸리에 변환으로 주어집니다. 표본의 갯수가 매우 많은 경우에 특성함수가 어떤 형태를 띠는지 살펴보고, 푸리에 역변환을 통해 확률밀도함수를 계산하는 것이 핵심입니다.
여기서는 연속적인 값을 가질 수 있는 변수에 대한 증명을 다루지만, 불연속적인 값을 가지는 변수의 확률 분포에 대해서도 같은 논리를 적용할 수 있습니다. 확률밀도함수가 불연속적인 위치에서 피크를 가지고, 피크의 높이가 확률에 비례하는 형태를 생각해보면 되겠죠. 그렇게 되면 정적분을 통해 서로 다른 불연속적인 값을 가질 확률을 도출할 수 있게 되겠습니다.
모집단이 가지는 변수의 값이 특정한 확률밀도함수를 따르고 표본을 추출하는 것이 서로 독립적이라면, 모든 표본들에 대한 확률밀도함수는 각 표본에 대한 확률밀도함수들의 곱으로 주어지게 됩니다. 이러한 확률 분포에 입각해서 표본평균의 특성함수를 계산할 수 있는데요. 여기서 표본들이 가질 수 있는 모든 값들을 고려하기 위한 다차원 정적분이 들어갑니다.
한 가지 더 주목할 점은 표본평균의 특성함수가 모집단의 특성함수들의 곱으로 주어질 수 있다는 것입니다. 정확히 말해서 특성함수에 인자로 들어가는 값을 표본의 갯수로 나누고, 특성함수를 표본의 갯수만큼 거듭제곱한 것과 동일합니다. 여기서 표본의 갯수가 매우 많은 경우의 극한값을 취하게 되면, 특성함수가 가우스 함수와 삼각함수의 곱으로 나타나는 복소수 함수의 형태를 띠게 된다는 결론에 도달하게 됩니다.
특성함수를 알았으니 확률밀도함수 역시 계산할 수 있습니다.
표본평균의 확률밀도함수가 가우스 함수의 형태를 띠고 있는것을 볼 수 있고, 이로부터 표본평균이 얼마나 큰 편차를 가지고 분포하는지에 대해서도 알 수 있습니다.
예시
주사위 던지기
중심 극한 정리를 확인할 수 있는 가장 간단한 예시로서 주사위를 여러번 던지는 상황을 생각해볼 수 있습니다. 주사위를 던지는 것은 서로 독립적이죠. 다시 말해서 첫번째로 주사위를 던져서 특정한 갯수의 눈이 나올 확률은 두번째로 던졌을 때와 관련이 없다는 것입니다.
서로 다른 갯수의 눈이 나올 확률은 6분의 1로 모두 동일하기 때문에, 모집단의 평균과 표준편차는 쉽게 계산할 수 있습니다. 다만 주사위를 던지는 것은 무한히 가능하므로 모집단의 크기가 무한하다고 볼 수 있는데요. 여기서 알아볼 것은 유한한 횟수만큼 주사위를 던져서 나온 눈의 갯수의 평균값을 구했을 때, 그 표본평균이 어떻게 분포하는가 하는 것입니다.
표본평균의 확률분포를 알아보기 위해서는 유한한 갯수의 표본을 추출하는 이벤트를 상정하고, 이 이벤트를 여러번 발생시켜야 합니다. 예를 들어서 주사위를 10번 던졌을 때의 표본평균의 분포를 알고 싶다면, 주사위를 10번 던지는 이벤트를 여러번 발생시켜야겠죠. 이벤트를 10만번 발생시킨다면, 주사위를 던지는 총 횟수는 100만이 되어야 하는 것입니다.
그런데 제가 직접 주사위를 100만번씩이나 던지기는 현실적으로 어렵기 때문에 시뮬레이션을 통해서 이를 알아볼까 합니다. 컴퓨터를 이용해서 사고실험을 해 보는 것인데요. 난수생성 함수를 포함한 C++ 프로그램을 통해 계산을 해 봅시다.
C++ 소스 코드입니다.
DiceRoller.h [다운로드]
#ifndef _DICEROLLER_H_
#define _DICEROLLER_H_
#include<stdlib.h>
#include<math.h>
class DiceRoller {
private :
unsigned long int n_sample_;
unsigned int *list_number_;
double mean_sample_;
double sigma_sample_;
double mean_pdf_;
double sigma_pdf_;
double sigma_mean_sample_;
unsigned int (*ptr_random_dice_)();
bool initialized_;
public :
DiceRoller() {
n_sample_ = 2;
initialized_ = false;
}
~DiceRoller() {
free_array();
}
void init(unsigned long int n_sample_in,
unsigned int (*ptr_in_random_dice)());
void next();
double get_mean_sample() {
return mean_sample_;
}
double get_sigma_sample() {
return sigma_sample_;
}
double get_sigma_mean_sample() {
return sigma_mean_sample_;
}
void free_array() {
if (initialized_) {
delete [] list_number_;
}
}
};
#endif
DiceRoller.cpp [다운로드]
#include"DiceRoller.h"
void DiceRoller::init(unsigned long int n_sample_in,
unsigned int (*ptr_in_random_dice)()) {
free_array();
n_sample_ = n_sample_in;
ptr_random_dice_ = ptr_in_random_dice;
list_number_ =
new unsigned int [n_sample_];
mean_pdf_ = 3.5;
sigma_pdf_ = sqrt(35. / 12.);
sigma_mean_sample_ =
sigma_pdf_ / sqrt((double)n_sample_);
initialized_ = true;
}
void DiceRoller::next() {
unsigned long int i;
mean_sample_ = 0.;
for (i = 0; i < n_sample_; i++) {
unsigned int num_now =
1 + (*ptr_random_dice_)() % 6;
list_number_[i] = num_now;
mean_sample_ += (double)num_now;
}
mean_sample_ =
mean_sample_ / (double)n_sample_;
sigma_sample_ = 0.;
for (i = 0; i < n_sample_; i++) {
double diff =
(double)list_number_[i] -
mean_sample_;
sigma_sample_ += fabs(diff * diff);
}
sigma_sample_ =
sigma_sample_ / ((double)n_sample_ - 1.);
sigma_sample_ = sqrt(sigma_sample_);
}
주사위를 던지는 시뮬레이션을 하기 위한 DiceRoller 클래스를 도입했습니다. 표본의 갯수와 난수생성을 위한 함수의 포인터를 인자로 주고 init 함수를 호출해서 초기화를 합니다. 그러면 next 함수를 호출함으로써 표본을 추출하고 표본평균과 표준편차 등을 구할 수 있게 되죠. 한 가지 짚고 넘어갈 점이 있다면 표본평균의 확률분포를 알아내기 위해서 next 함수를 여러번 호출하여 여러개의 이벤트를 발생시켜야 한다는 것입니다.
test1_rng_dice.cpp [다운로드]
#include<stdio.h>
#include<string.h>
#include<chrono>
#include<random>
#include"DiceRoller.h"
unsigned int seed_mt = 0;
std::mt19937 rng_mt(seed_mt);
unsigned int rng_uni() {
unsigned int num =
static_cast<unsigned int>(rng_mt());
return num;
}
// 표본의 갯수
unsigned long int n_sample = 64;
// 이벤트를 발생시키는 횟수
unsigned long int n_events = 1000000;
int nbin_x = 25;
double xmin = 1.;
double xmax = 6.;
double delta_x = (xmax - xmin) / (double)nbin_x;
int main(int argc, char *argv[]) {
seed_mt =
std::chrono::system_clock::now().time_since_epoch().count();
rng_mt.seed(seed_mt);
// 시뮬레이터의 선언과 초기화
DiceRoller roll;
roll.init(n_sample, &rng_uni);
double *pdf_x = new double [nbin_x];
for (int ix = 0; ix < nbin_x; ix++) {
pdf_x[ix] = 0.;
}
// 이벤트를 발생시키는 시뮬레이션 반복문
unsigned long int iev;
for (iev = 0; iev < n_events; iev++) {
roll.next();
double x_now = roll.get_mean_sample();
int ix = (int)floor((x_now - xmin) / delta_x);
if (ix < 0 || ix > nbin_x) {
continue;
}
pdf_x[ix] +=
1. / (delta_x * (double)n_events);
}
char filename[200];
char buffer[20];
strcpy(filename, "pdf_mean_dice_");
sprintf(buffer, "%lu", n_sample);
strcat(filename, buffer);
strcat(filename, "x");
sprintf(buffer, "%lu", n_events);
strcat(filename, buffer);
strcat(filename, ".dat");
// 표본평균의 확률분포 출력
FILE *fout;
fout = fopen(filename, "w");
fprintf(fout, "# n_sample = %lu\n", n_sample);
fprintf(fout, "# n_events = %lu\n", n_events);
fprintf(fout,
"# expected mean = 3.5 / sigma = %e\n",
roll.get_sigma_mean_sample());
for (int ix = 0; ix < nbin_x; ix++) {
double x1 = xmin + delta_x * (double)ix;
double x2 = x1 + delta_x;
fprintf(fout,
" %e %e %e\n", x1, x2, pdf_x[ix]);
}
fclose(fout);
delete [] pdf_x;
}
프로그램의 main 함수가 있는 소스 코드입니다. 여기서 표본의 갯수을 정하고, 이 이벤트를 몇 번 발생시킬지도 전역 변수 n_sample_ 및 n_events_를 통해 결정하게 됩니다. 앞서 언급한대로 DiceRoller 클래스의 init 함수와 next 함수를 호출하여 시뮬레이션을 하게 됩니다.
GNU 컴파일러를 사용하는 경우, 다음과 같이 코드를 빌드해서 실행파일을 얻을 수 있습니다.
- 컴파일
g++ DiceRoller.cpp -c
g++ test1_rng_dice.cpp -c - 링크
g++ test1_rng_dice.o DiceRoller.o -lm -o [실행파일 이름]
프로그램을 실행시켜 보면 시뮬레이션을 거쳐서, 주사위를 던졌을 때의 표본평균이 어떤 확률 분포를 따르는지를 볼 수 있게 됩니다.
여기서는 16개와 64개의 표본에 대해서 시뮬레이션을 하고, 이를 정규분포와 비교했는데요. 표본의 갯수가 16개인 경우에는 시뮬레이션을 통해 도출한 확률밀도함수와 정규 분포 사이에 차이가 좀 있습니다. 반면에 표본의 갯수를 64개로 늘리게 되면, 시뮬레이션을 통해 파악한 표본평균의 확률 분포가 가우스 함수에 매우 근접해 있는것을 볼 수 있죠. 중심 극한 정리가 시사하는 바가 잘 드러나는 부분입니다.
앞에서 C++ 프로그램을 하나 소개했는데, 소스 코드와 헤더 파일들을 이용해서 프로그램을 만드는 과정이 궁금하시다면 다음 포스팅이 큰 도움이 되리라 생각합니다.