Introduction

ANOVA는 Analysis of Variance의 약자로서, 2개 이상의 그룹의 평균을 비교하고자 할때 Anova를 사용합니다.
T-test와 유사한데 차이는 T-test는 2개의 그룹일때만 사용하고, 2개 이상일때는 ANOVA를 사용합니다.

  • Tests
    • Continuous VS Category : ANOVA 사용
    • Category VS Category : Chi-sqaure test 사용
  • Effect Size (얼마나 관련성이 있는지..)
    • Continous VS Category: intraclass correlation 사용
    • Category VS Category: Cramer’s V 사용
  • One-way ANOVA
    • 모든 그룹의 평균값은 동일한가? (예. 연령별 (20대, 30대, 40대, 50대) 에서 커피에 대한 선호도 차이 유무를 검정
    • t-test와 유사하지만, 2개 이상의 그룹을 비교할때 사용할 수 있다. (t-test는 2개만 가능)

Anova

Grand Mean

\[\mu_{GM} = \frac{\sum x}{N} = \frac{\sum n \bar{x}}{\sum n}\]

ANOVA에서는 2종류의 평균값이 사용이 됩니다.
첫번째는 각 그룹들의 표본평균 \(\mu_1, \mu_2, \mu_3, ...\) 그리고 grand mean \(\mu\) 를 사용합니다.
grand mean은 각각의 표본평균들에 대한 평균값이라고 보면 됩니다.

Null Hypothesis

  • \(H_0\) : \(\mu_1 = \mu_2 = .. = \mu_k\) 즉 귀무가설은 모든 그룹의 평균은 정확하게 일치한다 있다.
  • \(H_A\) : \(\mu_l \ne \mu_m\) : 최소한 한개의 그룹은 \(\mu_i\) 는 나머지 그룹과 다은 평균값을 갖고 있다.

\(\mu\) 는 해당 그룹의 평균이고, k는 그룹의 갯수입니다.
\(\mu_l\) 그리고 \(\mu_m\) 은 전체 샘플 중에 아무거나 특정 표본평균을 뜻 합니다.

Between Group Variability

\[SS_{\text{between}} = \sum^{k}_{i=1} n_i ( \bar{x}_i - \bar{X}_{GM} )^2\]
  • 서로 다른 두개의 그룹을 SS 로 표현을 합니다. SS는 Sum of squares의 약자.
  • 만약 각 그룹간의 표본평균이 서로 가깝다면, Grand Mean 과의 차이도 작을 것 입니다.
  • \(n_i\) 는 각 그룹의 갯수를 사용해서 각각의 제곱편차(squared deviation)에 weight값으로 사용을 합니다.

쉽게 공식을 풀어 쓰면 다음과 같습니다.

\[SS_{\text{between}} = n_1 (\bar{x}_1 - \bar{x}_{GM})^2 + n_2 (\bar{x}_2 - \bar{x}_{GM})^2 + n_3 (\bar{x}_3 - \bar{x}_{GM})^2 + ... + n_k (\bar{x}_k - \bar{x}_{GM})^2\]

아래의 두 표본의 분포를 보면, 서로 겹치는 부분이 있습니다. 각각의 표본평균과 grand mean의 차이는 크지 않을 것 입니다.

def plot_normal_dist(means, sigmas):
    for i, (sample_mean, sigma) in enumerate(zip(means, sigmas)):
        normal_data = np.random.normal(sample_mean, scale=sigma, size=100)
        sns.distplot(normal_data)

means = [4, 3.5]
sigmas = [0.2, 0.5]
plot_normal_dist(means, sigmas)

반면에 아래의 두 표본분포를 서로 크게 떨어져 있습니다. 이 경우 각각의 표본 평균과 grand mean은 서로 크게 차이가 날 것입니다.

plot_normal_dist([2, 5], [0.2, 0.3])

이렇게 각 표본분포간의 차이를 Between-group variability 라고 합니다.

between-group variability의 경우 자유도는 표본 평균의 갯수 k에 -1 을 한 것입니다.

\[MS_{\text{between}} = \frac{\sum^{k}_{i=1} n_i(\bar{x}_i - \bar(x)_{GM})}{k-1}\]
  • k : 그룹의 갯수

Within-Group Variability

\[\begin{align} SS_{\text{within}} &= \sum (n_{1i}-\bar{x}_1) + \sum(n_{2i}-\bar{x}_2) + ... + \sum(n_{ki}-\bar{x}_k) \\ &= \sum (x_{ij} - \bar{x}_j)^2 \end{align}\]

\(n_{1i}\) 라는 뜻은 첫번째 샘플에서의 i번째의 값을 의미 합니다.
\(n_{2i}\) 라는 뜻은 두번째 샘플에서의 i번째의 값을 의미 합니다.
\(n_{ji}\) 라는 뜻은 j번째 샘플에서의 i번째 값을 의미 합니다.

예를 들어서 아래 세개의 표본에 대한 분포를 그렸습니다.
각 그룹들이 옆으로 더 퍼질수록 (variability), 각 그룹들의 분포는 서로 겹치게 되고, 하나의 큰 population이 되게 됩니다.

plot_normal_dist([1, 20, 40], [7, 7, 7])

이번에는 위와 동일한 평균값을 갖고 있지만, 적은 변동성 (variability)을 갖고 있는 분포를 보도록 하겠습니다.
위와 동일하게 샘플들간의 표본평균은 동일하지만, 서로 완전히 다른 population에 종속이 되어 있습니다.

plot_normal_dist([1, 20, 40], [2, 2, 2])

위와 같이 각 샘플안에서의 variation을 Within-group Variation이라고 합니다.
즉 Between-group variability는 각 샘플들 사이의 표본평균의 차를 구했다면..
Within-group variation에서는 각각의 샘플들 마다 존재하는 편차를 구하게 됩니다.

Between-group variability에서 했듯이, 제곱편차의 합을 degrees of freedom으로 나누어 줄수 있습니다.
이렇게 나누어 주면 좀더 덜 편향된 (less biased) estimator를 구할 수 있습니다.

Within-group variability의 degrees of freedom은 전체 샘플 크기들의 합 (sum of the sample sizes) N 에다 샘플듯의 갯수 k 를 빼주게 됩니다.

\[\text{MS}_{\text{within}} = \frac{\sum (x_{ji} - \bar{x}_i)^2}{N-k}\]
  • MS: mean sum of square
  • N : 전체 데이터 갯수
  • k : 샘플들의 갯수

Degree of Freedom

각 공식별 degree of freedoms은 다음과 같이 정의 할 수 있습니다.

  • \(DOF(SS_{between})\) = \(k -1\)
  • \(DOF(SS_{within})\) = \(N - k\)
  • \(DOF(total)\) = \(N - 1\)

F Test

독립변수와 종속변수 사이의 관계성을 테스트하기 위해서는 F distribution이 사용이 됩니다.

\[F = \frac{MS_{between}}{MS_{within}}\]
  • \(MS_{\text{between}}\) : mean sum of squares due to treatment
  • \(MS_{\text{within}}\) : mean sum of squares due to error

표본평균들간의 차이를 나타낸 것을 F-Ratio 라고 합니다.
F-Ratio의 값이 1 또는 그 이하라면, 표본평균들은 서로 유사하다고 볼수 있으며,
표본평균들이 모두 모평균과 동일하다는 귀무가설이 맞게 됩니다..

만약 F-Ratio의 값이 매우 높다면, 귀무가설은 기각이 됩니다.

앞서 언급했듯이 between-group variability(그룹 간 변동성) 증가할수록, 표본평균들간의 거리는 서로 더 멀리 떨어져 있다고 볼 수 있습니다.
즉 멀어지면 멀어질수록 다른 모집단에 속해 있을 가능성이 큽니다.

Source: Dr. Asim’s Anatomy Cafe

z 그리고 t-distributions과는 다르게, F-distribution은 negative value가 없습니다.
이유는 between group이나 within group 모두 각 편차를 제곱(squaring each deviation)하기 때문입니다.

따라서 아래 그래프에서 보듯이 임계 영역(critical region)은 오른쪽 꼬리(right tail)에 한군데 존재합니다.

만약 F-statistic값이 해당 임계 영역안에 있다면 서로간의 표본평균들은 매우 다르기 때문에 귀무가설을 기각할 수 있습니다.

Source: Statistics How To

임계 영역을 나누는 critical value를 찾아야 하는데, 이때 F-table을 사용하게 됩니다.

One Way ANOVA

예제를 갖고서 설명을 하겠습니다.

그룹 X, Y, Z그룹으로 나뉘어 지며, 공부 방법에 따른 성적이라고 생각하면 됩니다.

data = pd.DataFrame([[7, 9, 5, 8, 6, 8, 6, 10, 7, 4], 
                     [4, 3, 6, 2, 7, 5, 5, 4, 1, 3],
                     [6, 1, 3, 5, 3, 4, 6, 5, 7, 3]]).T

data.columns = ('X', 'Y', 'Z')

print('[표본평균]')
print(data.mean(axis=0))
print()

print('[Grand Mean]')
print(data.mean(axis=0).mean())
[표본평균]
X    7.0
Y    4.0
Z    4.3
dtype: float64

[Grand Mean]
5.1000000000000005
  • 각 그룹의 표본평균: \(\mu_1 = 7\), \(\mu_2 = 4\), \(\mu_3 = 4.3\)
  • 전체 평균 (Grand Mean): \(\mu = 5.1\)

Boxplot 을 통한 분석

간단하게는 boxplot을 통해서 대략의 인싸이트를 얻을 수 있습니다.

아래에서도 언급하겠지만, boxplot만 갖고서 결론을 내리기에는 부족합니다.

sns.boxplot(data=data)
sns.swarmplot(data=data, color='.2')

Boxplot을 통해서 분산분석을 하는 것은 꽤 유용합니다.
다만 boxplot만 보고서 판단은 할 수가 없습니다.
\(SS_{between}\) 에서 보면, n * square of deviation 을 하게 되는데.. 여기서 n값이 weight값으로 사용이 됩니다.
즉 boxplot자체에는 각 샘플마다 얼마나 중요도가 있는지 알 수가 없으며,
또한 각 표본평균들간에 실제로 얼마나 차이가 나고 의미가 있는지, 분산이 큰것과는 어떤 관계가 있는지 등등 유의성을 알 수가 없습니다.
이 경우 일원분산분석이나 이원분산분석등을 사용하게 됩니다.

With Scipy

Scipy를 이용해서 두개 이상의 그룹이 동일한 모평균을 갖고 있다는 귀무가설을 테스트 할 수 있습니다.
이때 P-value가 유효하기 위해서는 다음의 가정이 충족되어야 합니다.

  • 표본들은 모두 독립적이다.
  • 각각의 표본들은 정규분포를 따른다
  • 모든 그룹들의 모표준편차 (population standard deviations of the groups) 는 모두 동일하다

만약 위의 가정들이 맞지 않는다면, Kruskal Wallis H-test 같은 다른 방법을 사용해야 합니다.

from scipy import stats

F_statistic, p_value = stats.f_oneway(data['X'], data['Y'], data['Z'])
print('F={0:.4f}, p={1:.5f}'.format(F_statistic, p_value))
F=8.1809, p=0.00167

위에서 P-value 값이 0.05보다 작기 때문에 통계적으로 유의미함을 나타냅니다.
또한 F-statistic값이 8.18이라는 뜻은 between-group variance가 within-group variance보다 약 8.18배 더 크다는 뜻입니다.
결론적으로 표본평균은 모두 동일하다는 귀무가설을 기각할수 있습니다.

Scipy에서 제공하는 라이브러리의 문제점은 F-critical value 를 따로 계산을 해주어야 합니다.
이때 scipy.stats.f 안의 ppf를 사용해서 critical value를 계산할 수 있습니다.

from scipy.stats import f, norm

N = data.shape[0] # 10
k = data.shape[1] # 3

dfn = k - 1 # degree of freedom numerator -> k - 1 
dfd = N - k # degree of freedom denominator -> N - K 
critical_value = f.ppf(q=0.95, dfn=dfn, dfd=dfd)

print('Degree of freedom Numerator  : {0} - 1 = {1}'.format(k, dfn) )
print('Degree of freedom Denominator: {0} - {1} = {2}'.format(N, k, dfd))
print('Critical Value (rejection)   :', critical_value)
Degree of freedom Numerator  : 3 - 1 = 2
Degree of freedom Denominator: 10 - 3 = 7
Critical Value (rejection)   : 4.73741412777588
def plot_f_distribution(dfn, dfd, f_value=None):
    critical_value = f.ppf(q=.95, dfn=dfn, dfd=dfd)
    mean, var, skew, kurt = f.stats(dfn, dfd, moments='mvsk')
    x = np.linspace(f.ppf(0.01, dfn, dfd), f.ppf(0.99, dfn, dfd), 100)
    
    # Plotting
    plt.plot(x, f.pdf(x, dfn, dfd), alpha=0.6, label='X~ F({0}, {1})'.format(dfn, dfd))
    plt.vlines(critical_value, 0.0, 1.0, linestyles='dashdot', 
               label='Crit. Value: {:.2f}'.format(critical_value), color='red')
    
    if f_value is not None:
        plt.vlines(f_value, 0., 1., label='F-value')
    plt.legend()
    plt.grid()
    

plot_f_distribution(dfn, dfd, 8.1809)

with Statsmodel

import statsmodels.api as sm
from statsmodels.formula.api import ols

moore = sm.datasets.get_rdataset("Moore", "carData", cache=True).data
moore.rename(columns={'partner.status': 'partner_status'}, inplace=True) # make name pythonic
display(moore.head())

moore_lm = ols('conformity ~ C(fcategory, Sum) * C(partner_status, Sum)', data=moore).fit()
display(moore_lm.summary())

sm.stats.anova_lm(moore_lm, typ=2)

MSB = 11.6147 / 2.
MSE = 817.763961/39.

print('F-statistic for fcategory:', MSB/MSE)
print('F-statistic for partner_status:', (212.213778/1.)/MSE)
print('F-statistic for fcategory + partner_status:', (175.488928/2.)/MSE)
F-statistic for fcategory: 0.276958463323624
F-statistic for partner_status: 10.120692200570085
F-statistic for fcategory + partner_status: 4.1846232644140695