Korean Perl Workshop 2012 발표자료 - Perl 101

Korean Perl Workshop 2012 발표자료

PDL, PDL::Stats 예제 t분포

t-분포는 정규분포와 유사한 확률분포이며, 모집단의 표준편차를 알지 못하며(표본의 표준편차로 모집단의 표준편차를 대용) 표본의 크기가 작을 때(n<30) 주로 이용한다. 자유도(n-1)가 낮을수록 봉이 낮아지고 꼬리가 높아지며 자유도가 높을수록 정규분포에 가까워진다.

#!/usr/bin/env perl
use 5.010;
use PDL;
use PDL::Stats;
use PDL::GSL::CDF;
use PDL::Graphics::PGPLOT::Window;
use Math::GSL::Randist qw/gsl_ran_tdist_pdf/;

# t 분포
my $win = pgwin();
my $x = zeroes(100)->xlinvals(-4, 4);
$win->env(-4,4,0,1,{Title=>'t distribution'});
$win->hold();
#CDF
$win->line($x, gsl_cdf_tdist_P($x, 29),{COLOR=>'RED'});
$win->line($x, gsl_cdf_tdist_P($x, 4),{COLOR=>'MAGENTA'});
#PDF
$win->line($x, [ map { gsl_ran_tdist_pdf($_, 29) } $x->list ],{COLOR=>'BLUE'});
$win->line($x, [ map { gsl_ran_tdist_pdf($_, 4) } $x->list ],{COLOR=>'CYAN'});

# t분포 자유도10, 우측꼬리확률 0.05일때 t확률변수 값은?
say gsl_cdf_tdist_Pinv(1-0.05, 10);
# t분포 자유도10, 양측꼬리확률 0.05일때 t확률변수 값은?
say gsl_cdf_tdist_Pinv(pdl(0.05/2, 1-0.05/2), 10);

1.81246112281168
[-2.2281389  2.2281389]

자유도29(표본수 30), 자유도4(표본수 5)에 대해서 PDF그래프를 보면 자유도가 낮을수록 펑퍼짐하며 높을수록 봉이 높아지며  정규분포에 가까워짐을 볼 수 있다. 분포가 펑퍼짐하다는 것은 표본수가 작을 수록 불확실성이 크기때문에 양쪽으로 치우칠 확률이 더 높다는 것을 의미한다.

PDL, PDL::Stats 예제 정규분포

정규분포(가우스분포)는 앞서 올린 이항분포의 n이 아주 커질 경우 근사적으로 따르게되는 확률분포이며(참고) 통계학에서 기본이 되는 연속확률변수 이다. 평균과 표준편차가 다른 많은 정규분포가 존재할 수 있으나 z점수를 통해 평균 0, 표준편차 1 인 표준정규분포로 변환해서 계산가능하다. 계산기 및 컴퓨터가 발달하지 않았던 옛날에는 주로 미리 계산된 표를 통해서 계산했지만 요즘에서는 굳이 그럴 필요없이 컴퓨터에 미리정의된 함수를 활용하면 간단하다.

다음은 표준정규분포의 확률밀도함수(PDF)와 누적분포함수(CDF)를 그리고 정규분포을 이용해서 몇가지 예에 대해 계산하는 코드이다.

#!/usr/bin/env perl
use 5.010;
use PDL;
use PDL::Stats;
use PDL::GSL::CDF;
use PDL::Graphics::PGPLOT::Window;
use Math::GSL::Randist qw/gsl_ran_gaussian_pdf/;

# 정규분포
my $win = pgwin();
my $x = zeroes(100)->xlinvals(-4, 4);
$win->env(-4,4,0,1,{Title=>'standard normal(gaussian) distribution'});
$win->hold();
#CDF
$win->line($x, gsl_cdf_gaussian_P($x,1), {COLOR=>'RED'});
#PDF
$win->line($x, [ map { gsl_ran_gaussian_pdf($_, 1) } $x->list ], {COLOR=>'BLUE'});
# PDL::Stats::Distr의 pdf_gaussian 함수 사용시
#$win->line($x, $x->pdf_gaussian(0, 1), {COLOR=>'BLUE'});

# 시험결과 학생들의 성적이 평균이 63점 분산이 100인 정규분포를 따른다고 하면
# 50점 미만의 학생은 몇 %인가 ?
# Z = (50-63)/sqrt(100) 이므로
say gsl_cdf_gaussian_P((50-63)/sqrt(100), 1)*100, " %";

# 상위 10%에게 A를 주면 A를 받기 위해서 몇점 이상이 되어야 하는가?
my $Z = gsl_cdf_gaussian_Pinv(0.9, 1); # 표준정규분포에서 하위(100-10)=90% 위치의 Z값
# Z = (X-63)/sqrt(100) 이므로 점수X는
say $Z*sqrt(100)+63

<결과>

9.68004845856103 %
75.815515655446

그래프에서 파란색이 PDF, 빨간색이 CDF의 그래프이다. 표준정규분포의 확률변수구간은 음의무한대에서 양의 무한대 까지이나 -4~4구간만 그렸으며 PDF아래의 면적이1 CDF는 PDF를 구간별로 적분한 값으로 0에서 1에 수렴한다.

PDL, PDL::Stats 예제 이항분포

앞서 올린 PDL, PDL::Stats 로 기본 통계량 계산하기 에 이어 이항분포에 대한 예를 들어보면 다음과 같다.

이항분포는  연속된 n번의 독립적 시행에서 각 시행이 확률 p를 가질 때의 확률분포를 말하는데 말로 하는 것 보다 다음 그림을 보면 그 의미가 명확해진다.

위 틀의 꼭대기에서 구슬을 떨어뜨리면 각 층에서 좌우로 갈 확률p는 반반(1/2)으로 각층을 통과 하며(독립시행) 구슬이 제각각의 칸으로 떨어지게 되는데 아래 구슬이 쌓인 모양이 바로 이항분포가 되는 것이다.

 그러면 n=30, p=0.5 일때 이항분포의 확률질량함수 PMF[Probability Mass Function](연속확률변수의 PDF에 해당), 누적분포함수 CDF[Cumulative Distribution Function] 의 예를 코딩해 보면 다음과 같다.

#!/usr/bin/env perl
use 5.010;
use PDL;
use PDL::Stats;
use PDL::GSL::CDF;
use Math::GSL::Randist qw/gsl_ran_binomial_pdf/;
use PDL::Graphics::PGPLOT::Window;

my $win1 = pgwin();
my $x = sequence(30);
$win1->env(0,30,0,1,{Title=>'binomial distribution'});
$win1->hold();

#CDF
# gsl_cdf_binomial_P는 PDL::Stats에 포함된 PDL::GSL::CDF 모듈의 함수
$win1->bin($x, gsl_cdf_binomial_P($x, 0.5, 30),{COLOR=>'RED'});

#PMF
# Math::GSL::Randist의 gsl_ran_binomial_pdf 함수를 사용
# PDL용 함수가 아니라서 piddle을 계산하지 못하므로 map으로 piddle $x내부의 값을
# map 함수를 사용하여 리스트를 순회하며 적용
$win1->bin($x, [ map { gsl_ran_binomial_pdf($_, 0.5, 30) } $x->list ], {COLOR=>'BLUE'});
# PDL::Stats::Distr의 pmf_binomial 함수를 사용할 시
# $win1->bin($x, $x->pmf_binomial(30, 0.5), {COLOR=>'BLUE'}); 

#동전을 30번 던졌을때 12번 앞면이 나올 확률
# PMF 그래프 가로축 12일때 세로축의 값
say gsl_ran_binomial_pdf(12, 0.5, 30) * 100, " %";
#say piddle(12)->pmf_binomial(30, 0.5) * 100, " %"; # pmf_binmoial 함수를 사용할 시
#동전을 30번 던졌을때 12번 이하로 앞면이 나올 확률
# CDF 그래프 가로축 12일때 세로축의 값
# 결론적으로 12까지의 PMF그래프 값을 누적한(적분) 것
say gsl_cdf_binomial_P(12, 0.5, 30) * 100, " %";

<결과>

8.05530929937956 %
18.0797304026779 %

위 그래프에서 파란색이 PMF이고 빨간색이 CDF이다. PMF그래프 아래의 면적을 구하면 1이 되고 CDF는 0에서 시작하여 1까지 증가한다는 사실을 보면 확률의 분포와 누적확률이 시각적으로 느껴질 것이라고 본다.

PDL, PDL::Stats 로 기본 통계량 계산하기

Perl에는 Matlab같은 PDL이란 것이 있다.

여기에 PDL::Stats 모듈을 같이 사용하면 오픈소스 통계패키지인 R같이 여러가지 통계계산도 할 수 있다.

다음은 어떤 데이터의 평균, 분산, 표준편차, 어떤 두 데이터간의 공분산과 상관계수를 구하는 예제 코드이다.

#!/usr/bin/env perl
use 5.010;
use PDL;
use PDL::Stats;
use PDL::Graphics::PGPLOT::Window;

my $data = pdl( [1,1,5,6,7,3,3,9,4,10,8,6,6,2,6,15] );
say "data: ",$data;
say "average: ",$data->average;
say "variance: ",$data->var;
say "standard deviation: ",$data->stdv;

my $x = pdl( [1,2,3,4,5,6] );
my $y = pdl( [1,1.5,2,2.25,2.5,3] );
my $y2 = pdl( [6,5,4,3,2,1] );

say "x:", $x;
say "y:", $y;
say "y2:", $y2;
say "covariance x,y: ",$x->cov($y);
say "correlation x,y: ",$x->corr($y);
                  
say "covariance x,y2: ",$x->cov($y2);
say "correlation: x,y2: ",$x->corr($y2);

my $win = pgwin();
$win->hold();
$win->env(0,10,0,10);
$win->points($x,$y, {SYMBOL=>'SQUARE', COLOR=>'BLUE'});
$win->points($x,$y2, {SYMBOL=>'PLUS', COLOR=>'CYAN'});

<결과>

data: [1 1 5 6 7 3 3 9 4 10 8 6 6 2 6 15]
average: 5.75
variance: 12.4375
standard deviation: 3.52668399491647
x:[1 2 3 4 5 6]
y:[1 1.5 2 2.25 2.5 3]
y2:[6 5 4 3 2 1]
covariance x,y: 1.10416666666667
correlation x,y: 0.991332701357794
covariance x,y2: -2.91666666666667
correlation: x,y2: -1

차후에는 여러가지 확률분포에 대해 테스트해보고자 한다.