import random
from pathlib import Path
import matplotlib.pylab as plt
import numpy as np
import pandas as pd
import statsmodels.api as sm
import statsmodels.formula.api as smf
from scipy import stats
%matplotlib inline파이썬 statsmodels로 통계분석
Python
Statistics
Stat
통계적 추론이라는 것은 제한된 실험 데이터에서 얻은 결과를 모집단에도 적용하려는 것입니다. 이번 포스트에서는 통계적 추론에 사용되는 검정법을 배워봅니다.
여기에 사용된 코드 출처
1 0. 통계적 추론과 유의성 검정
try:
import common
DATA = common.dataDirectory()
except ImportError:
DATA = Path().resolve() / "data"WEB_PAGE_DATA_CSV = DATA / "web_page_data.csv"
FOUR_SESSIONS_CSV = DATA / "four_sessions.csv"
CLICK_RATE_CSV = DATA / "click_rates.csv"
IMANISHI_CSV = DATA / "imanishi_data.csv"2 1. A/B 검정
A/B 검정은 두 가지 중에 어떤 것이 더 좋은지를 증명하기 위해 실험군을 나누어 진행하는 실험입니다. 일반적으로 대조군(Negative control)과 처리군으로 나눕니다.
예시로 웹페이지에서 고객들이 얼마나 머무는지 측정(session_time)하고 A/B 검정을 해보겠습니다.
session_times = pd.read_csv(WEB_PAGE_DATA_CSV)
session_times.Time = 100 * session_times.Time먼저 상자 그림 시각화를 통해 비교해보자
ax = session_times.boxplot(by="Page", column="Time", figsize=(4, 4))
ax.set_xlabel("")
ax.set_ylabel("Time (in seconds)")
plt.suptitle("")
plt.tight_layout()
plt.show()
그림으로 봐서는 페이지B의 방문객들이 좀 더 오래 머무는 것 같다. 각각의 평균값을 구해 수치로 비교해보자.
mean_a = session_times[session_times.Page == "Page A"].Time.mean()
mean_b = session_times[session_times.Page == "Page B"].Time.mean()
print(mean_b - mean_a)35.66666666666667대략적으로 36초정도 차이가 난다. 문제는 이 차이가 우연에 의한 것인지 아니면 통계적으로 중요한 것인지 판단하는 것이다.
2.1 1.1. 순열 검정을 통한 A/B 검정
순열검정(permutaion test): 두 개 이상의 표본을 함께 결합하여 관측값들을 무작위로 재표본으로 추출하는 과정을 말한다.
파이썬에서 순열 검정을 구현하기 위해 아래와 같이 perm_fun 함수를 정의합니다.
# Permutation test example with stickiness
def perm_fun(x, nA, nB):
n = nA + nB
idx_B = set(random.sample(range(n), nB))
idx_A = set(range(n)) - idx_B
return x.loc[idx_B].mean() - x.loc[idx_A].mean()
nA = session_times[session_times.Page == "Page A"].shape[0]
nB = session_times[session_times.Page == "Page B"].shape[0]
print(perm_fun(session_times.Time, nA, nB))24.23809523809524단 한번의 계산을 통해서 24초라는 차이가 발생하였습니다. 계산을 반복해서 페이지A와 페이지B 사이의 시간 차이에 대한 도수 분포표를 그려봅시다.
random.seed(1)
perm_diffs = [perm_fun(session_times.Time, nA, nB) for _ in range(1000)]
fig, ax = plt.subplots(figsize=(5, 5))
ax.hist(perm_diffs, bins=11, rwidth=0.9)
ax.axvline(x=mean_b - mean_a, color="black", lw=2)
ax.text(50, 190, "Observed\ndifference", bbox={"facecolor": "white"})
ax.set_xlabel("Session time differences (in seconds)")
ax.set_ylabel("Frequency")
plt.tight_layout()
plt.show()
위 그림에서 수직선은 관측된 차이입니다. 이것을 통해 순열 검정에서 가끔 실제 관찰된 차이를 넘어가는 것을 알 수 있습니다. 그렇다면 어느정도의 확률로 그런 일이 벌어질까요?
print(np.mean(perm_diffs > mean_b - mean_a))0.121답은 12.1% 입니다. 이것을 통해 페이지A와 페이지B의 차이인 36초가 통계적으로 봤을때는 차이가 없어도 약 12% 확률로 발생할 수 있다는 결론을 얻었습니다.
2.2 1.1.2. t-Test를 사용한 A/B 검정
t-테스트(t-test) 또는 t-검정 또는 스튜던트 t-테스트(Student’s t-test)는 검정통계량이 귀무가설 하에서 t-분포를 따르는 통계적 가설 검정법이다. t-테스트는 일반적으로 검정통계량이 정규 분포를 따르며 분포와 관련된 스케일링 변숫값들이 알려진 경우에 사용한다. 이때 모집단의 분산과 같은 스케일링 항을 알 수 없으나, 이를 데이터를 기반으로 한 추정값으로 대체하면 검정통계량은 t-분포를 따른다. 예를 들어 t-테스트를 사용하여 두 데이터 세트(집단)의 평균이 서로 유의하게 다른지 여부를 판별할 수 있다. -wiki
파이썬에서는 ttest_ind 함수를 사용하면 손쉽게 구해볼 수 있습니다.
res = stats.ttest_ind(
session_times[session_times.Page == "Page A"].Time,
session_times[session_times.Page == "Page B"].Time,
equal_var=False,
)
print(f"p-value for single sided test: {res.pvalue / 2:.4f}")p-value for single sided test: 0.1408결과 값이 앞서 구한 순열 검정의 확률인 12%과 유사한 수치인 14%임을 확인 할 수 있습니다. 컴퓨터가 보급되기 전에 순열 검정은 실용적이지 않았고 그래서 통계학자들에게 t-Test가 널리 사용되었습니다.
tstat, pvalue, df = sm.stats.ttest_ind(
session_times[session_times.Page == "Page A"].Time,
session_times[session_times.Page == "Page B"].Time,
usevar="unequal",
alternative="smaller",
)
print(f"p-value: {pvalue:.4f}")p-value: 0.14083 2. 통계적 유의성과 P-Values
통계적 유의성이란 실험 결과가 우연하게 일어난 것인지 아닌지 판단하는 것입니다. 유의 확률(p-value)는 통계적 유의성을 판단하기 위해 가정한 귀무가설에서 극단적인 결과를 얻을 확률입니다. 따라서 p-value가 낮을수록 우연한 일이 아니라고 생각 할 수 있습니다.
다음의 예시 데이터를 가지고 예를 들어보도록 하겠습니다. 동일한 내용이지만 구성이 다른 상품판매 페이지A, 페이지B가 있고 방문자가 구매한 횟수와 그렇지 않은 데이터를 모았다고 생각합시다.
| 결과 | 페이지A | 페이지B |
|---|---|---|
| 구매 | 200 | 182 |
| 구매하지 않음 | 23539 | 22406 |
페이지A가 페이지B보다 약 5% 정도 우수한 결과를 보여주었는데 이것이 통계적으로 유의성을 갖는지 순열 검정을 통해 살펴봅시다.
random.seed(1)
obs_pct_diff = 100 * (200 / 23739 - 182 / 22588)
print(f"Observed difference: {obs_pct_diff:.4f}%")
conversion = [0] * 45945
conversion.extend([1] * 382)
conversion = pd.Series(conversion)
perm_diffs = [100 * perm_fun(conversion, 23739, 22588) for _ in range(1000)]
fig, ax = plt.subplots(figsize=(5, 5))
ax.hist(perm_diffs, bins=11, rwidth=0.9)
ax.axvline(x=obs_pct_diff, color="black", lw=2)
ax.text(0.06, 200, "Observed\ndifference", bbox={"facecolor": "white"})
ax.set_xlabel("Conversion rate (percent)")
ax.set_ylabel("Frequency")
plt.tight_layout()
plt.show()Observed difference: 0.0368%
약 1000번의 순열 검정을 통해 얻은 히스토그램을 살펴보니 관측된 값은 충분히 우연에 의한 것으로 보입니다. 그러나 좀 더 명확하게 하기 위해 p-value를 계산해봅시다.
3.1 2.1. P-Value 구하기
p-value는 순열 검정에서 얻은 결과 중에 관찰된 차이와 같거나 더 큰 차이를 보이는 경우의 비율이라고 할 수 있기에 다음과 같이 추정할 수 있습니다.
print(np.mean([diff > obs_pct_diff for diff in perm_diffs]))0.332예상한 것처럼 30%의 확률로 우연에 의해서 나타날 수 있는 차이였습니다. 따라서 페이지A와 페이지B의 차이는 통계적으로 유의미하지 않다고 말할 수 있겠습니다.
4 3. 일원 분산분석(ANOVA )
A/B 검정이 아닌 여러 그룹간의 통계적 유의미한 차이를 검정하는 절차를 분산분석이라 합니다. 예시로 페이지 1,2,3,4에 방문자들이 머문 시간에 대한 데이터를 살펴보겠습니다. 먼저 데이터셋을 불러옵니다.
print(pd.read_csv(FOUR_SESSIONS_CSV).head()) Page Time
0 Page 1 164
1 Page 2 178
2 Page 3 175
3 Page 4 155
4 Page 1 172네 그룹의 상자그림을 시각화해서 어떤 차이가 있는지 살펴봅시다.
four_sessions = pd.read_csv(FOUR_SESSIONS_CSV)
ax = four_sessions.boxplot(by="Page", column="Time", figsize=(4, 4))
ax.set_xlabel("Page")
ax.set_ylabel("Time (in seconds)")
plt.suptitle("")
plt.title("")
plt.tight_layout()
plt.show()
4.1 3.1. 순열검정을 통한 one-way ANOVA
파이썬에서는 다음 코드를 사용해 순열 검정을 통해 ANOVA 분석을 진행할 수 있습니다.
observed_variance = four_sessions.groupby("Page").mean().var()[0]
print("Observed means:", four_sessions.groupby("Page").mean().values.ravel())
print("Variance:", observed_variance)
# Permutation test example with stickiness
def perm_test(df):
df = df.copy()
df["Time"] = np.random.permutation(df["Time"].values)
return df.groupby("Page").mean().var()[0]
print(perm_test(four_sessions))Observed means: [172.8 182.6 175.6 164.6]
Variance: 55.426666666666655
57.02666666666678random.seed(1)
perm_variance = [perm_test(four_sessions) for _ in range(3000)]
print("Pr(Prob)", np.mean([var > observed_variance for var in perm_variance]))
fig, ax = plt.subplots(figsize=(5, 5))
ax.hist(perm_variance, bins=11, rwidth=0.9)
ax.axvline(x=observed_variance, color="black", lw=2)
ax.text(60, 200, "Observed\nvariance", bbox={"facecolor": "white"})
ax.set_xlabel("Variance")
ax.set_ylabel("Frequency")
plt.tight_layout()
plt.show()Pr(Prob) 0.07766666666666666
Pr(Prob)의 값은 p-value이며 결과는 0.07입니다. 통상적인 임계 p-value 값인 0.05이상임으로 네 페이지간의 차이가 우연히 발생할 수 있다고 결론 내릴 수 있습니다.
4.2 3.2. F-통계량을 통한 one-way ANOVA
F-통계량은 잔차 오류에 인한 분산과 그룹 평균의 분산에 대한 비율을 기초로 합니다. 비율이 높으면 통계적으로 유의미 하다고 할 수 있고 이를 토대로 p-value를 계산할 수 있습니다.
model = smf.ols("Time ~ Page", data=four_sessions).fit()
aov_table = sm.stats.anova_lm(model)
print(aov_table) df sum_sq mean_sq F PR(>F)
Page 3.0 831.4 277.133333 2.739825 0.077586
Residual 16.0 1618.4 101.150000 NaN NaNdf는 자유도, sum_sq는 제곱합, mean_sq는 평균제곱, F는 F-통계량을 나타냅니다.
res = stats.f_oneway(
four_sessions[four_sessions.Page == "Page 1"].Time,
four_sessions[four_sessions.Page == "Page 2"].Time,
four_sessions[four_sessions.Page == "Page 3"].Time,
four_sessions[four_sessions.Page == "Page 4"].Time,
)
print(f"F-Statistic: {res.statistic / 2:.4f}")
print(f"p-value: {res.pvalue / 2:.4f}")F-Statistic: 1.3699
p-value: 0.0388F-통계량을 사용한 방법은 p-value 값이 더 적게나와 임계값인 0.05 이하입니다. 그러나 ANOVA 분석의 p-value가 낮게 나왔다고 해서 모든 그룹에서 통계적으로 차이가 있다고 할 수는 없습니다. 추가적인 Ad hoc 분석을 진행해 어떤 그룹에서 차이가 있는지 확인해보아야 합니다.
5 4. 마무리하며
현재 실험 데이터를 분석할때는 통계적 추론과 유의성 검증은 필수적인 절차입니다. 그렇기 때문에 반대로 실험에 대조군을 반드시 포함하며 데이터에 인위적인 조작이 없도록 하는 것이 더욱 중요해 졌습니다.