14 Sep 2017 » statistics

2017-09-14-probability

수식이 깨질 경우 새로고침을 눌러주세요.


가설과 검정 ③

One-Tailed Test(단측 검정) vs Two-Tailed Test(양측 검정)

One-Tailed Test(단측 검정)

부호가 결정되어 있는 경우 $$H_0\;:\; \mu \leq 0\;:\;\; \mu가\; 양수이면서\; 희귀한\; 검정통계량\; 값이\; 나오는\; 경우에만\; 기각 $$

$$\alpha = 0.05\;인\; 경우 \;\; critical\; value는\; \text{CDF}(1-\alpha) = \text{CDF}(0.95) $$
$$H_0\;:\; \mu \geq 0\;:\;\; \mu가\; 음수이면서\; 희귀한\; 검정통계량\; 값이\; 나오는\; 경우에만\; 기각 $$
$$\alpha = 0.05\;인\; 경우 \;\; critical\; value는\; \text{CDF}(\alpha) = \text{CDF}(0.05) $$

Two-Tailed Test(양측 검정)

$$ H_0\;:\;\; \mu = 0\;:\;\; 부호와\; 상관없이\; 크기가\; 더\; 희귀한\; 검정통계량\; 값이\; 나오는\; 경우\; 기각 $$$$\alpha = 0.05\;인\; 경우 \;\; critical\; value는\; \text{CDF}(\alpha/2) = \text{CDF}(0.025)\; 와 \;\text{CDF}(1-\alpha/2) = \text{CDF}(0.975) $$
In [1]:
# -*- coding: utf-8 -*-

import seaborn as sns
import pandas as pd
import scipy as sp
import matplotlib as mpl
import matplotlib.pylab as plt
import numpy as np
%matplotlib inline
from matplotlib.patches import Polygon


x = np.linspace(-5, 5, 100)
y = sp.stats.norm.pdf(x)

def plot_tailed(f, subplot, alpha, title=None):
    ax = f.add_subplot(subplot)
    ax.plot(x, y, linewidth=2)
    plt.hold(True)
    
    # critical point
    cr1 = sp.stats.norm.ppf(alpha)
    
    ## np.sign(x) -> -1, 0 으로 출력(x < 0 → -1, 나머지 경우는 0)
    cr2 = cr1 + np.sign(cr1) * 10
    
    ix = np.linspace(cr1, cr2, 100)
    iy = sp.stats.norm.pdf(ix)
    
    verts = [(cr1, 0)] + list(zip(ix, iy)) + [(cr2, 0)]
    poly = Polygon(verts)
    ax.add_patch(poly)
    plt.xlim(-5, 5)
    if title is not None:
        plt.title(title)
In [2]:
alpha = 0.05
f = plt.figure(figsize=(20, 5))
plot_tailed(f, 131, alpha, "Negative One-Tailed")
plot_tailed(f, 133, 1 - alpha, "Positve One-Tailed")
plot_tailed(f, 132, alpha/2)
plot_tailed(f, 132, 1 - alpha/2, "Two-Tailed")
/Users/Leo/.pyenv/versions/anaconda3-4.0.0/envs/code_study/lib/python3.6/site-packages/ipykernel_launcher.py:19: MatplotlibDeprecationWarning: pyplot.hold is deprecated.
    Future behavior will be consistent with the long-time default:
    plot commands add elements without first clearing the
    Axes and/or Figure.
/Users/Leo/.pyenv/versions/anaconda3-4.0.0/envs/code_study/lib/python3.6/site-packages/matplotlib/__init__.py:917: UserWarning: axes.hold is deprecated. Please remove it from your matplotlibrc and/or style files.
  warnings.warn(self.msg_depr_set % key)
/Users/Leo/.pyenv/versions/anaconda3-4.0.0/envs/code_study/lib/python3.6/site-packages/matplotlib/rcsetup.py:152: UserWarning: axes.hold is deprecated, will be removed in 3.0
  warnings.warn("axes.hold is deprecated, will be removed in 3.0")

자주 쓰이는 검정 유형

이산 분포 파라미터 테스트

  • Bernoilli 분포의 파라미터 값이 theha인가 아닌가
    • Binomial Test
  • Categorical 분포의 파라미터 값이 theta vector인가 아닌가
    • Chi-Square Test

정규 분포 파라미터 테스트

  • 정규 분포의 mean 값이 0인가 아닌가
    • One-sample z-test (분산값을 아는 경우)
    • One-sample t-test (분산값을 모르는 경우)
  • 두 독립 정규 분포의 mean 값이 같은가 다른가
    • Independent-Two-sample t-test
  • 두 Paired 정규 분포의 mean 값이 같은가 다른가
    • Paired-Two-sample t-test
  • 두 정규 분포의 variance 값이 같은가 다른가
    • Equal-Variance test

분포 확인 테스트

  • 분포가 정규 분포인가 아닌가
    • Normality test
  • 두 분포가 같은 분포인가 아닌가 ← 기준이 되는 분포를 하나를 두고(대조군)
    • Kolmogorov Smirnov test

기초적인 검정을 실습해보자.

이항 검정 (Binomial test)

이항 검정은 이항 분포를 이용하여 Bernoulli 분포 모수 theta에 대한 가설을 조사하는 검정 방법이다.

$$ ex \;1)\; 데이터\; 갯수\; N=10,\; 실제\; 모수\; \theta = 0.5인\; 경우,\; 이항\; 검정을\; 해보자. $$
In [3]:
N = 10
theta = 0.5
np.random.seed(0)
x = sp.stats.bernoulli(theta).rvs(N)
n = np.count_nonzero(x)
n
Out[3]:
7
In [4]:
# N번 시행 중 n번 성공할 확률
sp.stats.binom_test(n, N)
Out[4]:
0.34374999999999989

유의 확률(p-value) = 34%로 높으므로 다음 귀무 가설을 기각할 수 없다.

$$ H_0 : \theta = 0.5 $$

$$ ex \;2)\; 데이터\; 갯수\; N=100,\; 실제\; 모수\; \theta = 0.5인\; 경우,\; 이항\; 검정을\; 해보자. $$
In [5]:
N = 100
theta = 0.5
np.random.seed(0)
x = sp.stats.bernoulli(theta).rvs(N)
n = np.count_nonzero(x)
n
Out[5]:
49
In [6]:
sp.stats.binom_test(n, N)
Out[6]:
0.92041076261282062

유의 확률(p-value) = 92%로 높으므로 다음 귀무 가설을 기각할 수 없다.

$$ H_0 : \theta = 0.5 $$

$$ ex \;3)\; 데이터\; 갯수\; N=100,\; 실제\; 모수\; \theta = 0.35인\; 경우,\; 이항\; 검정을\; 해보자. $$
In [7]:
N = 100
theta = 0.35
np.random.seed(0)
x = sp.stats.bernoulli(theta).rvs(N)
n = np.count_nonzero(x)
n
Out[7]:
31
In [8]:
sp.stats.binom_test(n, N)
Out[8]:
0.00018314322488235352

유의 확률(p-value)이 0.018%로 낮으므로 다음 귀무 가설을 기각한다.

$$ H_0 : 기각 → \theta \neq 0.5 $$

카이 제곱 검정 (Chi-square test)

카이 제곱 검정은 goodness of fit 검정이라고도 부른다.

Multinoulli 시도시 성공한 횟수로 파라미터 theta 벡터가 정해진 값을 만족하는지 검사한다. $$ 카테고리\; 분포의\; 모수\; \theta=(\theta_1, \ldots, \theta_K)에\; 대한\; 가설을\; 조사하는\; 검정\; 방법이다. $$

SciPy stats 서브패키지의 chisquare 명령을 사용한다.

  • 디폴트 귀무 가설은 $\theta = \left(\frac{1}{K}, \ldots, \frac{1}{K} \right)$이다.

  • scipy.stats.chisquare

$$ ex 1)\;\;데이터\; 개수\; N=10,\; 실제\; 모수\; \theta\;=\;(0.25,\; 0.25,\; 0.25,\; 0.25)인\; 경우\; 카이\; 제곱\; 검정을\; 해보자. $$
In [9]:
N = 10
K = 4
theta = np.ones(K)/K
np.random.seed(0)

# k개의 카테고리에서 p의 확률로, 이를 10번 수행해라.
x = np.random.choice(K, N, p=theta)
n = np.bincount(x, minlength=K)
n
Out[9]:
array([0, 3, 5, 2])
In [10]:
# 카테고리로 나눠진 어레이를 인풋으로 넣어준다.
sp.stats.chisquare(n)
Out[10]:
Power_divergenceResult(statistic=5.1999999999999993, pvalue=0.157724450396663)

유의 확률(p-value)이 15.7%로 높으므로 귀무 가설을 기각할 수 없다.

$$ 따라서\; H_0 : \theta \;=\; (0.25,\; 0.25,\; 0.25,\; 0.25)이다. $$

$$ ex 2)\;\;데이터\; 개수\; N=100,\; 실제\; 모수\; \theta\;=\;(0.35,\; 0.20,\; 0.20,\; 0.15)인\; 경우\; 카이\; 제곱\; 검정을\; 해보자. $$
In [11]:
N = 100
K = 4
theta = np.array([0.35, 0.30, 0.20, 0.15])
np.random.seed(0)

# k개의 카테고리에서 p의 확률로, 이를 10번 수행해라.
x = np.random.choice(K, N, p=theta)
n = np.bincount(x, minlength=K)
n
Out[11]:
array([37, 32, 20, 11])
In [12]:
sp.stats.chisquare(n)
Out[12]:
Power_divergenceResult(statistic=16.559999999999999, pvalue=0.00087034719789121269)

유의 확률(p-value)이 0.087%이므로 귀무 가설을 기각할 수 있다.

$$ 따라서\;H_0 : 기각 → \;\theta \neq (0.25,\; 0.25,\; 0.25,\; 0.25)이다. $$

단일 표본 z-검정 (One-sample z-test)

단일 표본 z-검정은 모분산의 값을 정확히 알고 있는 정규 분포의 표본에 대해 기대값을 조사하는 검정 방법이다.

  • 단, 실제로 이런 경우는 없어 사용할 일이 없다.

단일 표본 z-검정의 경우에는 SciPy에 별도의 함수가 준비되어 있지 않으므로 norm 명령의 cdf 메서드를 사용하여 직접 구현해야 한다.

$$ ex 1)\;\;데이터\; 개수\; N=10,\; 실제\; 모수\; \mu = 0인\; 경우\; 단일\; 표본\; z검정을\; 해보자. $$
In [13]:
N = 10
mu = 0
np.random.seed(0)
x = sp.stats.norm(mu).rvs(N)
x
Out[13]:
array([ 1.76405235,  0.40015721,  0.97873798,  2.2408932 ,  1.86755799,
       -0.97727788,  0.95008842, -0.15135721, -0.10321885,  0.4105985 ])
In [14]:
def z_test_1_sampling(x, sigma_2=1, mu=0):
    z = (x.mean() - mu)/ np.sqrt(sigma_2 / len(x))
    p_value = 2 * sp.stats.norm().sf(np.abs(z))
    return (z, p_value)
In [15]:
z_test_1_sampling(x)
Out[15]:
(2.3338341854824276, 0.019604406021683538)

유의 확률(p-value) = 0.0196 이므로 유의 수준 = 0.05 이라면 귀무 가설을 기각할 수 있다.

$$H_0 : 기각 → \mu \neq 0$$

이 경우는 검정 결과가 오류인 예라고 볼 수 있다.

  • 검정 결과가 오류로 나온 이유는 데이터 수가 10개로 부족하기 때문이다.
  • 제 1종 오류(Type 1 Error) : 귀무 가설이 참임에도 기각하는 경우
$$ ex 2)\;\;데이터\; 개수\; N=100,\; 실제\; 모수\; \mu = 0인\; 경우\; 단일\; 표본\; z검정을\; 해보자. $$
In [16]:
N = 100
mu = 0
np.random.seed(0)
x = sp.stats.norm(mu).rvs(N)
In [17]:
z_test_1_sampling(x)
Out[17]:
(0.59808015534484993, 0.54978645086241684)

유의 확률(p-value) = 0.5497 이므로 귀무 가설을 기각할 수 없다.

$$ H_0 : \mu = 0 $$

단일 표본 t-검정 (One-sample t-test)

정규 분포의 표본에 대해 기대값을 조사하는 검정 방법이다.

scipy의 stats 서브 패키지의 ttest_1samp을 사용한다.

  • ttest_1samp 명령의 경우에는 디폴트 모수가 없으므로 popmean 인수를 사용하여 직접 지정해야 한다.
  • scipy.stats.ttest_1samp
$$ ex 1)\;\;데이터\; 개수\; N=10,\; 실제\; 모수\; \mu = 0인\; 경우\; 단일\; 표본\; t-test를\; 해보자. $$
In [18]:
N = 10
mu = 0
np.random.seed(0)
x = sp.stats.norm(mu).rvs(N)
sns.distplot(x, kde=False, fit=sp.stats.norm)
plt.show()
In [19]:
sp.stats.ttest_1samp(x, popmean=0)
Out[19]:
Ttest_1sampResult(statistic=2.2894396723896699, pvalue=0.047818464908570578)

p-value = 0.0468 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 있다.

$$H_0 : 기각 → \mu \neq 0 $$

마찬가지로 제 1종 오류(귀무가설이 참일 때 기각할 오류, 샘플수가 적기 때문이다.)

$$ ex 2)\;\;데이터\; 개수\; N=100,\; 실제\; 모수\; \mu = 0인\; 경우\; 단일\; 표본\; t-test를\; 해보자. $$
In [20]:
N = 100
mu = 0
np.random.seed(0)
x = sp.stats.norm(mu).rvs(N)
sns.distplot(x, kde=False, fit=sp.stats.norm)
plt.show()
In [21]:
sp.stats.ttest_1samp(x, popmean=0)
Out[21]:
Ttest_1sampResult(statistic=0.59042834028516977, pvalue=0.55624891586946745)

p-value = 0.5562 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 없다.

$$H_0 : \mu = 0 $$

독립 표본 t-검정 (Independent-two-sample t-test)

Independent-two-sample t-test는 간단하게 two sample t-test라고도 한다.

  • 두 개의 독립적인 정규 분포에서 나온 두 개의 데이터 셋을 사용하여 두 정규 분포의 기대값이 동일한지를 검사한다.
  • scipy stats 서브패키지의 ttest_ind를 사용한다.
  • 독립 표본 t-검정은 두 정규 분포의 분산값이 같은 경우와 같지 않은 경우에 사용하는 검정 통계량이 다르다.
    • equal_var 인수를 사용하여 이를 지정해 주어야 한다.
  • scipy.stats.ttest_ind
$$ ex)\;\;데이터\; 개수\; N_1 = N_2 = 10,\; 두\; 정규\; 분포의\; 기대값이\; \mu_1 = 0,\; \mu_2 = 1로\; 다르고,\\ \\ 분산은\; \sigma_1 = \sigma_2 = 1으로\; 같을 때,\; two-sample-t-test를\; 해보자. $$
In [22]:
N_1 = 10; mu_1 = 0; sigma_1 = 1
N_2 = 10; mu_2 = 0.5; sigma_2 = 1
np.random.seed(0)
x1 = sp.stats.norm(mu_1, sigma_1).rvs(N_1)
x2 = sp.stats.norm(mu_2, sigma_2).rvs(N_2)

plt.figure(figsize=(12,6))
sns.distplot(x1, kde=False, fit=sp.stats.norm)
sns.distplot(x2, kde=False, fit=sp.stats.norm)
plt.show()
In [23]:
sp.stats.ttest_ind(x1, x2, equal_var=True)
Out[23]:
Ttest_indResult(statistic=-0.41399685269886549, pvalue=0.68376768941164268)

p-value = 0.6837 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 없다.

$$H_0 : \mu_1 = \mu_2 $$

이 경우는 검정 결과가 오류다. 제 2종 오류(Type 2 Error) 귀무 가설이 거짓임에도 참으로 나온 경우

  • 데이터 수를 더 늘려야 한다.

샘플 수를 더 늘려서 다시 해보면

In [24]:
N_1 = 50; mu_1 = 0; sigma_1 = 1
N_2 = 100; mu_2 = 0.5; sigma_2 = 1
np.random.seed(0)
x1 = sp.stats.norm(mu_1, sigma_1).rvs(N_1)
x2 = sp.stats.norm(mu_2, sigma_2).rvs(N_2)

plt.figure(figsize=(12,6))
sns.distplot(x1, kde=False, fit=sp.stats.norm)
sns.distplot(x2, kde=False, fit=sp.stats.norm)
plt.show()
In [25]:
sp.stats.ttest_ind(x1, x2, equal_var=True)
Out[25]:
Ttest_indResult(statistic=-2.6826951236616963, pvalue=0.0081339709157226582)

p-value = 0.0081 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 있다.

$$H_0 : 기각 → \mu_1 \neq \mu_2 $$

대응 표본 t-검정 (Paired-two-sample t-test)

대응 표본 t-검정은 독립 표본 t-검정을 두 집단의 샘플이 1대1 대응하는 경우에 대해 수정한 것이다.

즉, 독립 표본 t-검정과 마찬가지로 두 정규 분포의 기대값이 같은 지 확인하기 위한 검정

  • 예를 들어 어떤 반의 학생들이 특강을 수강하기 전과 수강한 이후에 각각 시험을 본 시험 점수의 경우에는 같은 학생의 두 점수는 대응할 수 있다.
  • 이 대응 정보를 알고 있다면 보통의 독립 표본 t-검정에서 발생할 수 있는 샘플 간의 차이의 영향을 없앨 수 있기 때문에 특강 수강의 영향을 보다 정확하게 추정할 수 있다.
  • scipy.stats.ttest_rel
$$ ex)\;\;데이터\; 개수\; N = 5,\; 두\; 정규\; 분포의\; 기대값이\; \mu_1 = 0,\; \mu_2 = 0.5로\; 달라진 경우,\; Paired-two-sample \;t-test를\; 해보자. $$
In [26]:
N = 5
mu_1 = 0
mu_2 = 0.5

np.random.seed(1)
x1 = sp.stats.norm(mu_1).rvs(N)
# 표준편차 = 0.1 추가로 설정
x2 = x1 + sp.stats.norm(mu_2, 0.1).rvs(N)

plt.figure(figsize=(12,6))
sns.distplot(x1, kde=False, fit=sp.stats.norm)
sns.distplot(x2, kde=False, fit=sp.stats.norm)
plt.show()
In [27]:
sp.stats.ttest_rel(x1, x2)
Out[27]:
Ttest_relResult(statistic=-7.1723380661732756, pvalue=0.0020008849290622677)

p-value = 0.002 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 있다.

$$H_0 : 기각 → \mu_1 \neq \mu_2 $$

5 개의 데이터만으로도 두 평균이 다르다는 것을 알아낼 수 있다.


카이 제곱 분산 검정 (Chi-Square Test for the Variance)

지금까지는 정규 분포의 기대값을 비교하는 검정을 살펴보았다. 이제는 정규 분포의 분산에 대해 살펴보자.

카이 제곱 분산 검정(Chi-Square Test for the Variance)은 정규 분포의 샘플 분산값은 정규화하면 카이 제곱 분포를 따른다는 점을 이용하는 검정 방법이다.

scipy는 카이 제곱 분산 검정에 대한 명령이 없으므로 chi2 클래스를 사용하여 직접 구현해야 한다.

In [28]:
def chi2var_test(x, sigma2=1):
    v = x.var(ddof=1)
    t = (len(x) - 1)*v/sigma2
    p_value = sp.stats.chi2(df=len(x)-1).sf(np.abs(t))
    return (t, p_value)
In [29]:
N = 10
mu_0 = 0
sigma_0 = 1.1

np.random.seed(0)
x = sp.stats.norm(mu_0, sigma_0).rvs(N)

plt.figure(figsize=(12,6))
sns.distplot(x, kde=False, fit=sp.stats.norm)
plt.show()
x.std()
Out[29]:
1.0637871321863899
In [30]:
chi2var_test(x)
Out[30]:
(11.316430626053437, 0.25464123584764542)

p-value = 0.2546 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 없다.

$$H_0 : 채택 → \sigma_0 = 1.1 $$

등분산 검정 (Equal-variance test)

등분산 검정은 두 정규 분포로부터 생성된 두 개의 데이터 집합으로부터 두 정규 분포의 분산 모수가 같은지 확인하기 위한 검정이다.

  • 가장 기본적인 방법은 F분포를 사용하는 것이지만 실무에서는 이보다 더 성능이 좋은 bartlett, fligner, levene 방법을 주로 사용한다.
  • scipy의 stats 서브패키지는 이를 위한 bartlett, fligner, levene 명령을 제공한다.
  • scipy.stats.bartlett
  • scipy.stats.fligner
  • scipy.stats.levene
In [31]:
N1 = 100
N2 = 100
sigma_1 = 1
sigma_2 = 1.2

np.random.seed(0)
x1 = sp.stats.norm(0, sigma_1).rvs(N1)
x2 = sp.stats.norm(0, sigma_2).rvs(N2)

plt.figure(figsize=(12,6))
sns.distplot(x1, kde=False, fit=sp.stats.norm)
sns.distplot(x2, kde=False, fit=sp.stats.norm)
plt.show()
x1.std(), x2.std()
Out[31]:
(1.0078822447165796, 1.2416003969261071)
In [32]:
sp.stats.bartlett(x1, x2)
Out[32]:
BartlettResult(statistic=4.2534738372322662, pvalue=0.039170128783651344)
In [33]:
sp.stats.fligner(x1, x2)
Out[33]:
FlignerResult(statistic=7.2248419904094572, pvalue=0.0071901501067483673)
In [34]:
sp.stats.levene(x1, x2)
Out[34]:
LeveneResult(statistic=7.6807089476794372, pvalue=0.0061135154970207925)
bartlett fligner levene
p-value 0.03917 0.0071 0.0061

위에서 보듯이 세 가지의 등분산 검정 결과 모두 귀무 가설을 기각할 수 있다.

$$H_0 : 기각 → \sigma_1 \neq \sigma_2 $$

정규성 검정

회귀 분석 등에서는 확률 분포가 가우시안 정규 분포를 따르는지 아닌지를 확인하는 것이 중요하다.

이러한 검정을 정규성 검정(Normality Test)이라고 한다.

statsmodels에서 제공하는 정규성 검정 명령어

scipy 에서 제공하는 정규성 검정 명령어

이 중에서 Kolmogorov-Smirnov 검정은 정규 분포에 국한되지 않고 두 샘플이 같은 분포를 따르는지 확인할 수 있는 방법이다.

In [35]:
np.random.seed(0)
N1 = 50
N2 = 100

x1 = sp.stats.norm(0, 1).rvs(N1)
x2 = sp.stats.norm(0.5, 1.5).rvs(N2)

plt.figure(figsize=(12,6))
sns.distplot(x1)
sns.distplot(x2)
plt.show()
In [36]:
sp.stats.ks_2samp(x1, x2)
Out[36]:
Ks_2sampResult(statistic=0.23000000000000004, pvalue=0.049516112814422863)

p-value = 0.04951 이므로 유의 수준 = 0.05일 때 귀무가설을 기각할 수 있다.

$$H_0 : 기각 → X_1,\; X_2가\; 따르는\; 분포는\; 서로\; 다른 \;분포라고\; 볼\; 수\; 있다. $$


Related Posts