[처음으로] [모든 글 보기] [같은 카테고리의 글 보기: 메모, 통계]
개인적으로 학부 졸업하고 나서 가장 아쉬운 것 중 하나가 학부 다닐 때 통계학을 제대로 공부하지 않은 것이다. 통계가 왜 중요한가? 어떤 것을 잘 하기 위해 필요한 첫 번째 단계는 그것을 얼마나 잘 하고 있는 지 재는 것이고, 이것을 제대로 재기 위해서는 통계가 필요하기 때문이다.
이 메모는 이 책을 읽으면서 가설 검증에 관련된 챕터 내용들을 정리한 것이다. 이 글은 남에게 보여주기 위해 쓴 것이 아니라 단지 개인적인 용도를 위해 정리한 것이고, 나는 통계학의 전문가와는 거리가 아주 멀기 때문에, 틀린 내용이 얼마든지 있을 수 있다. 틀린 점은 이메일이나 트위터로 알려주면 고맙겠다.
다양한 종류의 가설 검정 기법이 있으며 어떤 자료를 비교하는지, 샘플은 몇 개인지에 따라 각각 다른 기법을 써야 한다.
모든 관찰값이 쌍을 이루고 있는 경우이다. 이 경우 한 쌍은 테스트하고 싶은 변수 외의 모든 요소가 같아야 한다.
같은 가정에서 자란 일란성 쌍둥이인데, 한명은 정신분열이고 한명은 아닌 경우들을 찾는다. 이 때 이 둘의 왼쪽 해마의 크기를 MRI를 이용해 비교하자. 이 때 각 샘플은 정신분열이라는 점만 빼면 모든 점이 똑같다.
각 샘플마다, 정신분열이 아닌 사람의 해마 부피에서 정신분열인 사람의 해마 부피를 뺀다. 이 차이의 집합이 이 문제의 샘플이 된다. 따라서 단 1개의 샘플이 있게 된다.
해마 부피 차이의 기대값은 0일까?
scipy 패키지에 $t$ 분포의 구현이 있기 때문에 간단하게 구현할 수 있다.
from math import sqrt
from scipy.stats import t
def paired_t_test(diff):
n = len(diff)
mean = diff.mean()
std = diff.std()
se = std / sqrt(n)
df = n - 1
t_ratio = mean / se
print 't ratio (assuming mean is zero): %.8lf' % t_ratio
rv = t(df)
p_value = (1.0 - rv.cdf(abs(t_ratio))) * 2
print 'two-sided p value: %.8lf' % p_value
# t-stat = (mean - pop_mean) / se
# pop_mean = mean - t-stat * se
pop_mean = lambda t_stat: mean - t_stat * se
hi, lo = rv.ppf([0.025, 0.975])
print '95%% confidence interval: %.4lf ~ %.4lf' % (pop_mean(lo), pop_mean(hi))
가장 흔한 경우로, 두 개의 서로 다른 분포에서 얻은 두 개의 샘플이 있다고 하자. 이 두 분포의 평균이 다를까? 다르다면 얼마나 다를까?
혹독한 겨울 폭풍이 몰아친 이후 야생동물 보호 센터(?)로 보내진 59마리의 참새 중, 35마리는 살아남고 24마리는 죽었다. 이들의 덩치와 살아남을 가능성 사이에 관계가 있는지 알아보기 위해 이들의 날개뼈 길이를 쟀다.
Pooled sample standard deviation($s_p$)으로부터 차이의 standard error 구하기:
$$
SE(\bar{Y}_{2} - \bar{Y}_{1}) = s_{p} \sqrt{ \frac{1}{n_1} + \frac{1}{n_2}}
$$
def two_sample_t_test(sample1, sample2):
n1, n2 = len(sample1), len(sample2)
mean1, mean2 = sample1.mean(), sample2.mean()
var1, var2 = sample1.var(), sample2.var()
df = n1 + n2 - 2
pooled_sd = sqrt((var1 * (n1 - 1) + var2 * (n2 - 1)) / df)
print 'pooled SD = %.8lf' % pooled_sd
diff_se = pooled_sd * sqrt(1.0 / n1 + 1.0 / n2)
print 'Diff SE = %.8lf' % diff_se
t_ratio = (mean1 - mean2) / diff_se
print 't ratio (assuming diff is zero): %.8lf' % t_ratio
rv = t(df)
p_value = (1.0 - rv.cdf(abs(t_ratio))) * 2
print 'two-sided p value: %.8lf' % p_value
# t-stat = (diff - pop_mean) / diff_se
# pop_mean = diff - t-stat * se
pop_mean = lambda t_stat: (mean1 - mean2) - t_stat * diff_se
hi, lo = rv.ppf([0.025, 0.975])
print '95%% confidence interval: %.4lf ~ %.4lf' % (pop_mean(lo), pop_mean(hi))
$t$-test는 모집단이 정규분포를 띠며, 두 집단의 표준 편차가 같다고 가정한다. 샘플의 수가 크면 Central Limit Theorem에 의해 샘플의 평균값은 정규분포를 띠므로 이 가정이 깨져도 상관없다고 주장하고 싶지만, 모집단의 표준 편차를 샘플 표준 오차로 근사하는 과정이 이 가정을 사용하기 때문에 이 가정은 아주 중요하다.
다행히 모집단이 정규분포가 아니라도 $t$-test는 상당부분 성립함이 널리 알려져 있다. 안 되는 시나리오들을 정리해 보면
두 모집단의 표준 편차 | 샘플 크기 | 스큐 | 판정 |
---|---|---|---|
같을때 | 다를때 | 클때 | 망했다 |
같을때 | 아주 다를 때 | 망했다 | |
다를때 | 다를때 | 망했다 |
outlier의 존재와 serial/cluster effect의 존재 또한 $t$-test 의 가정을 망쳐놓을 수 있다.
ANOVA는 다음 글에서 따로 다루겠다.