Welch의 이원배치 분산분석¶

import pingouin as pg
import pandas as pd

# 예시 자료
data = {
    "Temperature": ["High", "High", "High", "Low", "Low", "Low", "Medium", "Medium", "Medium"],
    "Fertilizer": ["A", "B", "C", "A", "B", "C", "A", "B", "C"],
    "Growth": [12, 15, 14, 10, 13, 11, 14, 16, 15],
}
df = pd.DataFrame(data)

# welch_anova는 between에 요인을 **하나만** 받는다. 그래서 따로 두 번 돌린다.
# 이렇게 하면 각 주효과는 다른 요인을 무시한 채 계산되며, 교호작용은 볼 수 없다.
print(pg.welch_anova(dv="Growth", between="Temperature", data=df))
print(pg.welch_anova(dv="Growth", between="Fertilizer", data=df))

        Source  ddof1     ddof2         F     p_unc       np2
0  Temperature      2  3.819209  5.173804  0.081821  0.645833
       Source  ddof1     ddof2         F     p_unc       np2
0  Fertilizer      2  3.915497  1.466347  0.334719  0.333333

import numpy as np
from scipy import stats


def welch_f_by_hand(groups):
    k = len(groups)
    n = np.array([len(g) for g in groups], float)
    m = np.array([g.mean() for g in groups])
    s2 = np.array([g.var(ddof=1) for g in groups])
    w = n / s2
    W = w.sum()
    Lam = (((1 - w / W) ** 2) / (n - 1)).sum()
    num = (w * (m - (w * m).sum() / W) ** 2).sum() / (k - 1)
    return num / (1 + 2 * (k - 2) / (k ** 2 - 1) * Lam), (k ** 2 - 1) / (3 * Lam)


for fac in ("Temperature", "Fertilizer"):
    gs = [v.values.astype(float) for _, v in df.groupby(fac)['Growth']]
    F, d2 = welch_f_by_hand(gs)
    print(f"{fac:>12}: 손계산 F = {F:.6f}, ddof2 = {d2:.6f}, "
          f"p = {stats.f(2, d2).sf(F):.6f}")

# 반복 없는 3x3 가법모형
M = df.pivot(index='Temperature', columns='Fertilizer', values='Growth').astype(float)
print("\n자료를 3x3 표로 보면")
print(M)
row, col = M.mean(axis=1), M.mean(axis=0)
grand = M.values.mean()
SSA = 3 * ((row - grand) ** 2).sum()
SSB = 3 * ((col - grand) ** 2).sum()
SST = ((M.values - grand) ** 2).sum()
SSE = SST - SSA - SSB
MSE = SSE / 4
print(f"\nSS(온도) = {SSA:.4f} (df 2),  SS(비료) = {SSB:.4f} (df 2),  "
      f"SS(오차) = {SSE:.4f} (df 4)")
print(f"SST = {SST:.4f},  MSE = {MSE:.6f}")
for name, SS in (("온도", SSA), ("비료", SSB)):
    F = (SS / 2) / MSE
    print(f"{name}: F = {F:.4f},  p = {stats.f(2, 4).sf(F):.6f}")

 Temperature: 손계산 F = 5.173804, ddof2 = 3.819209, p = 0.081821
  Fertilizer: 손계산 F = 1.466347, ddof2 = 3.915497, p = 0.334719

자료를 3x3 표로 보면
Fertilizer      A     B     C
Temperature
High         12.0  15.0  14.0
Low          10.0  13.0  11.0
Medium       14.0  16.0  15.0

SS(온도) = 20.6667 (df 2),  SS(비료) = 10.6667 (df 2),  SS(오차) = 0.6667 (df 4)
SST = 32.0000,  MSE = 0.166667
온도: F = 62.0000,  p = 0.000977
비료: F = 32.0000,  p = 0.003460

# Temperature에 대한 Games-Howell 사후검정
post_hoc_temp = pg.pairwise_gameshowell(dv="Growth", between="Temperature", data=df)
print(post_hoc_temp.round(4).to_string(index=False))

# Fertilizer에 대한 Games-Howell 사후검정
post_hoc_fert = pg.pairwise_gameshowell(dv="Growth", between="Fertilizer", data=df)
print(post_hoc_fert.round(4).to_string(index=False))

   A      B  mean_A  mean_B    diff     se       T     df   pval  hedges
High    Low 13.6667 11.3333  2.3333 1.2472  1.8708 4.0000 0.2604  1.2220
High Medium 13.6667 15.0000 -1.3333 1.0541 -1.2649 3.4483 0.4915 -0.8262
 Low Medium 11.3333 15.0000 -3.6667 1.0541 -3.4785 3.4483 0.0660 -2.2722
A B  mean_A  mean_B    diff     se       T     df   pval  hedges
A B 12.0000 14.6667 -2.6667 1.4530 -1.8353 3.7409 0.2767 -1.1988
A C 12.0000 13.3333 -1.3333 1.6667 -0.8000 3.9936 0.7230 -0.5226
B C 14.6667 13.3333  1.3333 1.4907  0.8944 3.6697 0.6740  0.5842

from itertools import combinations

g = df.groupby('Temperature')['Growth'].agg(['mean', 'var', 'count'])
k = 3
print(f"{'pair':<14}{'diff':>9}{'se':>8}{'T':>9}{'df':>9}{'pval':>8}{'문턱 T*':>10}")
for a, b in combinations(['High', 'Low', 'Medium'], 2):
    ma, va, na = g.loc[a]
    mb, vb, nb = g.loc[b]
    se = np.sqrt(va / na + vb / nb)
    T = (ma - mb) / se
    nu = (va / na + vb / nb) ** 2 / ((va / na) ** 2 / (na - 1)
                                     + (vb / nb) ** 2 / (nb - 1))
    pval = stats.studentized_range.sf(np.sqrt(2) * abs(T), k, nu)
    thr = stats.studentized_range.ppf(0.95, k, nu) / np.sqrt(2)
    print(f"{a + '-' + b:<14}{ma - mb:>9.4f}{se:>8.4f}{T:>9.4f}{nu:>9.4f}"
          f"{pval:>8.4f}{thr:>10.4f}")

# 수준당 관측값을 늘리면 문턱이 얼마나 내려가는가 (분산은 그대로 둔다).
va, vb = 2.3333333, 1.0
print(f"\n{'수준당 n':>9}{'nu':>9}{'문턱 T*':>10}{'필요한 평균차':>14}")
for n in (3, 4, 5, 8, 20):
    se = np.sqrt(va / n + vb / n)
    nu = (va / n + vb / n) ** 2 / ((va / n) ** 2 / (n - 1) + (vb / n) ** 2 / (n - 1))
    thr = stats.studentized_range.ppf(0.95, k, nu) / np.sqrt(2)
    print(f"{n:>9}{nu:>9.3f}{thr:>10.4f}{thr * se:>14.4f}")

pair               diff      se        T       df    pval     문턱 T*
High-Low         2.3333  1.2472   1.8708   4.0000  0.2604    3.5640
High-Medium     -1.3333  1.0541  -1.2649   3.4483  0.4915    3.8425
Low-Medium      -3.6667  1.0541  -3.4785   3.4483  0.0660    3.8425

    수준당 n       nu     문턱 T*       필요한 평균차
        3    3.448    3.8425        4.0504
        4    5.172    3.2155        2.9353
        5    6.897    2.9558        2.4134
        8   12.069    2.6658        1.7208
       20   32.759    2.4546        1.0021

import numpy as np
import pandas as pd
import pingouin as pg
import statsmodels.api as sm
from statsmodels.formula.api import ols

# 칸마다 반복이 있고 분산도 칸마다 다른 2×3 자료
rng = np.random.default_rng(202)
a, b = 2, 3
mu = np.array([[10.0, 13.0, 16.0],
               [11.0, 14.0, 24.0]])       # A1·B2 에서 교호작용
sd = np.array([[1.0, 2.0, 6.0],
               [1.5, 3.0, 7.0]])          # 칸마다 다른 분산
ns = np.array([[12, 10, 8],
               [9, 12, 11]])              # 불균형

rows = []
for i in range(a):
    for j in range(b):
        rows.append(pd.DataFrame({"y": rng.normal(mu[i, j], sd[i, j], ns[i, j]),
                                  "A": f"A{i}", "B": f"B{j}"}))
df = pd.concat(rows, ignore_index=True)

print(df.groupby(["A", "B"]).y.agg(["size", "mean", "var"]).round(3).to_string())

print("\npingouin 주변 웰치 (다른 요인 무시)")
for f in ["A", "B"]:
    r = pg.welch_anova(dv="y", between=f, data=df)
    print(f"  {f}: F = {r.F[0]:8.4f}, df2 = {r.ddof2[0]:7.3f}, "
          f"p = {r.p_unc[0]:.4f}")

print("\n표준 이원배치 (제II형)")
print(sm.stats.anova_lm(ols("y ~ C(A)*C(B)", data=df).fit(),
                        typ=2).round(4).to_string())

       size    mean     var
A  B
A0 B0    12  10.109   1.272
   B1    10  13.713   1.814
   B2     8  16.181  21.283
A1 B0     9  11.282   6.237
   B1    12  13.790   5.857
   B2    11  24.174  71.906

pingouin 주변 웰치 (다른 요인 무시)
  A: F =   6.1709, df2 =  44.378, p = 0.0168
  B: F =  24.6111, df2 =  33.987, p = 0.0000

표준 이원배치 (제II형)
              sum_sq    df        F  PR(>F)
C(A)        124.1164   1.0   6.8635  0.0113
C(B)        987.4697   2.0  27.3031  0.0000
C(A):C(B)   178.9091   2.0   4.9468  0.0105
Residual   1012.6747  56.0      NaN     NaN

import numpy as np
import pandas as pd
from scipy import stats
import pingouin as pg

def welch_james(m, v, n, Cm):
    """요한센(1980)의 웰치-제임스 검정. Cm 은 칸 평균에 대한 대비 행렬."""
    m, v, n = np.asarray(m, float), np.asarray(v, float), np.asarray(n, float)
    V = np.diag(v / n)
    Q = Cm.T @ np.linalg.inv(Cm @ V @ Cm.T) @ Cm
    T = m @ Q @ m
    q = np.linalg.matrix_rank(Cm)
    VQ = V @ Q
    A = 0.0
    for i in range(len(n)):
        E = np.zeros_like(V)
        E[i, i] = 1
        M = VQ @ E
        A += 0.5 * (np.trace(M)**2 + np.trace(M @ M)) / (n[i] - 1)
    c = q + 2 * A - 6 * A / (q + 2)
    f = q * (q + 2) / (3 * A)
    return T / c, q, f, stats.f.sf(T / c, q, f)

def factorial_contrasts(a, b):
    """a×b 요인설계의 주효과·교호작용 대비 행렬."""
    Ja, Jb = np.ones(a) / a, np.ones(b) / b
    Ca, Cb = np.zeros((a - 1, a)), np.zeros((b - 1, b))
    for i in range(a - 1):
        Ca[i, i], Ca[i, i + 1] = 1, -1
    for j in range(b - 1):
        Cb[j, j], Cb[j, j + 1] = 1, -1
    return np.kron(Ca, Jb), np.kron(Ja, Cb), np.kron(Ca, Cb)

# --- 검증 1: 일원배치에서 pingouin 과 일치하는가 ---
rng = np.random.default_rng(42)
ns1, mus, sds = [12, 18, 9], [0, 0.4, 0.9], [1.0, 2.0, 3.0]
gs = [rng.normal(u, s, n) for n, u, s in zip(ns1, mus, sds)]
C1 = np.array([[1.0, -1, 0], [0, 1, -1]])
F, q, f, p = welch_james([g.mean() for g in gs],
                         [g.var(ddof=1) for g in gs], ns1, C1)
print(f"웰치-제임스: F = {F:.6f}, df = ({q}, {f:.6f}), p = {p:.6f}")
d1 = pd.DataFrame({"y": np.concatenate(gs),
                   "g": np.concatenate([[f"G{i}"] * n
                                        for i, n in enumerate(ns1)])})
r = pg.welch_anova(dv="y", between="g", data=d1)
print(f"pingouin  : F = {r.F[0]:.6f}, df = ({r.ddof1[0]}, "
      f"{r.ddof2[0]:.6f}), p = {r.p_unc[0]:.6f}")

# --- 적용 2: 연습문제 4 의 2×3 자료 ---
rng = np.random.default_rng(202)
a, b = 2, 3
mu = np.array([[10.0, 13.0, 16.0], [11.0, 14.0, 24.0]])
sd = np.array([[1.0, 2.0, 6.0], [1.5, 3.0, 7.0]])
ns = np.array([[12, 10, 8], [9, 12, 11]])
rows = []
for i in range(a):
    for j in range(b):
        rows.append(pd.DataFrame({"y": rng.normal(mu[i, j], sd[i, j], ns[i, j]),
                                  "A": f"A{i}", "B": f"B{j}"}))
df = pd.concat(rows, ignore_index=True)
cell = df.groupby(["A", "B"]).y.agg(["size", "mean", "var"])
m, v, n = cell["mean"].values, cell["var"].values, cell["size"].values

CA, CB, CAB = factorial_contrasts(a, b)
print("\n2×3 웰치-제임스")
for lab, Cm in [("A 주효과", CA), ("B 주효과", CB), ("A×B", CAB)]:
    F, q, f, p = welch_james(m, v, n, Cm)
    print(f"  {lab:8s} F = {F:8.4f}  df = ({q}, {f:7.3f})  p = {p:.4f}")

웰치-제임스: F = 2.062345, df = (2, 17.185272), p = 0.157479
pingouin  : F = 2.062345, df = (2, 17.185272), p = 0.157479

2×3 웰치-제임스
  A 주효과    F =   8.0112  df = (1,  21.182)  p = 0.0100
  B 주효과    F =  24.8803  df = (2,  26.461)  p = 0.0000
  A×B      F =   3.2044  df = (2,  26.461)  p = 0.0567

import warnings
warnings.filterwarnings("ignore")

import numpy as np
import pandas as pd
from scipy import stats
import statsmodels.api as sm
from statsmodels.formula.api import ols

def welch_james(m, v, n, Cm):
    m, v, n = np.asarray(m, float), np.asarray(v, float), np.asarray(n, float)
    V = np.diag(v / n)
    Q = Cm.T @ np.linalg.inv(Cm @ V @ Cm.T) @ Cm
    T = m @ Q @ m
    q = np.linalg.matrix_rank(Cm)
    VQ = V @ Q
    A = 0.0
    for i in range(len(n)):
        E = np.zeros_like(V); E[i, i] = 1
        M = VQ @ E
        A += 0.5 * (np.trace(M)**2 + np.trace(M @ M)) / (n[i] - 1)
    c = q + 2 * A - 6 * A / (q + 2)
    return stats.f.sf(T / c, q, q * (q + 2) / (3 * A))

def factorial_contrasts(a, b):
    Ja, Jb = np.ones(a) / a, np.ones(b) / b
    Ca, Cb = np.zeros((a - 1, a)), np.zeros((b - 1, b))
    for i in range(a - 1):
        Ca[i, i], Ca[i, i + 1] = 1, -1
    for j in range(b - 1):
        Cb[j, j], Cb[j, j + 1] = 1, -1
    return np.kron(Ca, Jb), np.kron(Ja, Cb), np.kron(Ca, Cb)

a, b = 2, 3
CA, CB, CAB = factorial_contrasts(a, b)
NS = np.array([[15, 12, 6], [14, 10, 5]])     # 작은 칸에 큰 분산 = 역페어링
SD = np.array([[1.0, 2.0, 5.0], [1.0, 2.0, 5.0]])

rng = np.random.default_rng(1212)
B = 1_500
cnt = {k: np.zeros(3) for k in ["표준 F", "웰치-제임스", "HC3"]}
for _ in range(B):
    rows = []
    for i in range(a):
        for j in range(b):
            rows.append(pd.DataFrame({"y": rng.normal(0, SD[i, j], NS[i, j]),
                                      "A": f"A{i}", "B": f"B{j}"}))
    df = pd.concat(rows, ignore_index=True)
    cell = df.groupby(["A", "B"]).y.agg(["size", "mean", "var"])
    m, v, n = cell["mean"].values, cell["var"].values, cell["size"].values

    t = sm.stats.anova_lm(ols("y ~ C(A)*C(B)", data=df).fit(), typ=2)
    cnt["표준 F"] += np.array([t.loc["C(A)", "PR(>F)"],
                               t.loc["C(B)", "PR(>F)"],
                               t.loc["C(A):C(B)", "PR(>F)"]]) < 0.05
    cnt["웰치-제임스"] += np.array([welch_james(m, v, n, CA),
                                    welch_james(m, v, n, CB),
                                    welch_james(m, v, n, CAB)]) < 0.05

    fit = ols("y ~ C(A, Sum)*C(B, Sum)", data=df).fit(cov_type="HC3")
    idx = list(fit.params.index)
    ps = []
    for pick in [lambda s: s.startswith("C(A, Sum)[") and ":" not in s,
                 lambda s: s.startswith("C(B, Sum)[") and ":" not in s,
                 lambda s: ":" in s]:
        nm = [x for x in idx if pick(x)]
        R = np.zeros((len(nm), len(idx)))
        for r, x in enumerate(nm):
            R[r, idx.index(x)] = 1
        ps.append(float(fit.f_test(R).pvalue))
    cnt["HC3"] += np.array(ps) < 0.05

print("2×3 완전 귀무, n=[[15,12,6],[14,10,5]], σ=[[1,2,5],[1,2,5]], B=1500")
print(f"{'방법':>12s} {'A 주효과':>9s} {'B 주효과':>9s} {'A×B':>9s}")
for k, x in cnt.items():
    print(f"{k:>12s} {x[0] / B:9.4f} {x[1] / B:9.4f} {x[2] / B:9.4f}")

2×3 완전 귀무, n=[[15,12,6],[14,10,5]], σ=[[1,2,5],[1,2,5]], B=1500
          방법     A 주효과     B 주효과       A×B
        표준 F    0.0687    0.2400    0.2113
      웰치-제임스    0.0467    0.0467    0.0447
         HC3    0.0487    0.0540    0.0493

import numpy as np
from statsmodels.stats.libqsturng import qsturng

rng = np.random.default_rng(9090)
M = 20_000
print("세 집단, 집단당 n, 등분산 σ=1 — 게임스-하웰이 한 쌍을 잡을 확률")
print(f"{'n':>4s} {'d=1':>8s} {'d=2':>8s} {'d=3':>8s} {'d=4':>8s}")
for n in [3, 5, 10, 20]:
    out = []
    for d in [1, 2, 3, 4]:
        hit = 0
        for _ in range(M // 4):
            x = rng.normal(0, 1, n)
            y = rng.normal(d, 1, n)
            s = x.var(ddof=1) / n + y.var(ddof=1) / n
            df = s**2 / ((x.var(ddof=1) / n)**2 / (n - 1)
                         + (y.var(ddof=1) / n)**2 / (n - 1))
            hit += abs(x.mean() - y.mean()) / np.sqrt(s / 2) > qsturng(0.95, 3, df)
        out.append(hit / (M // 4))
    print(f"{n:4d} " + " ".join(f"{o:8.4f}" for o in out))

print("\n최소 탐지 차이 (등분산 가정)")
for n in [3, 5, 10]:
    df = 2 * (n - 1)
    q = qsturng(0.95, 3, df)
    print(f"  n={n:2d}: q(0.95, 3, {df}) = {q:.4f} "
          f"→ 최소 탐지 차이 ≈ {q * np.sqrt(1 / n):.3f}σ")

세 집단, 집단당 n, 등분산 σ=1 — 게임스-하웰이 한 쌍을 잡을 확률
   n      d=1      d=2      d=3      d=4
   3   0.0562   0.2112   0.4624   0.7126
   5   0.1434   0.5858   0.9290   0.9972
  10   0.3922   0.9590   1.0000   1.0000
  20   0.7526   1.0000   1.0000   1.0000

최소 탐지 차이 (등분산 가정)
  n= 3: q(0.95, 3, 4) = 5.0332 → 최소 탐지 차이 ≈ 2.906σ
  n= 5: q(0.95, 3, 8) = 4.0375 → 최소 탐지 차이 ≈ 1.806σ
  n=10: q(0.95, 3, 18) = 3.6080 → 최소 탐지 차이 ≈ 1.141σ

import numpy as np
import pandas as pd

mu = np.array([[10.0, 13.0, 16.0],
               [11.0, 14.0, 24.0]])
ns = np.array([[12, 10, 8],
               [9, 12, 11]])

print("칸 평균")
print(pd.DataFrame(mu, index=["A0", "A1"],
                   columns=["B0", "B1", "B2"]).to_string())
print("\n칸 크기")
print(pd.DataFrame(ns, index=["A0", "A1"],
                   columns=["B0", "B1", "B2"]).to_string())

un = mu.mean(1)
wt = (mu * ns).sum(1) / ns.sum(1)
print(f"\n{'':4s} {'비가중 주변평균':>14s} {'가중 주변평균':>13s}")
for i, lab in enumerate(["A0", "A1"]):
    print(f"{lab:4s} {un[i]:14.4f} {wt[i]:13.4f}")
print(f"{'차이':4s} {un[1] - un[0]:14.4f} {wt[1] - wt[0]:13.4f}")

unb = mu.mean(0)
wtb = (mu * ns).sum(0) / ns.sum(0)
print(f"\nB 수준별  비가중 {np.round(unb, 3).tolist()}   "
      f"가중 {np.round(wtb, 3).tolist()}")

칸 평균
      B0    B1    B2
A0  10.0  13.0  16.0
A1  11.0  14.0  24.0

칸 크기
    B0  B1  B2
A0  12  10   8
A1   9  12  11

           비가중 주변평균       가중 주변평균
A0          13.0000       12.6000
A1          16.3333       16.5938
차이           3.3333        3.9938

B 수준별  비가중 [10.5, 13.5, 20.0]   가중 [10.429, 13.545, 20.632]

import numpy as np
import pandas as pd
from scipy import stats

def welch_james(m, v, n, Cm):
    m, v, n = np.asarray(m, float), np.asarray(v, float), np.asarray(n, float)
    V = np.diag(v / n)
    Q = Cm.T @ np.linalg.inv(Cm @ V @ Cm.T) @ Cm
    T = m @ Q @ m
    q = np.linalg.matrix_rank(Cm)
    VQ = V @ Q
    A = 0.0
    for i in range(len(n)):
        E = np.zeros_like(V); E[i, i] = 1
        M = VQ @ E
        A += 0.5 * (np.trace(M)**2 + np.trace(M @ M)) / (n[i] - 1)
    c = q + 2 * A - 6 * A / (q + 2)
    f = q * (q + 2) / (3 * A)
    return T / c, q, f, stats.f.sf(T / c, q, f)

rng = np.random.default_rng(202)
a, b = 2, 3
mu = np.array([[10.0, 13.0, 16.0], [11.0, 14.0, 24.0]])
sd = np.array([[1.0, 2.0, 6.0], [1.5, 3.0, 7.0]])
ns = np.array([[12, 10, 8], [9, 12, 11]])
rows = []
for i in range(a):
    for j in range(b):
        rows.append(pd.DataFrame({"y": rng.normal(mu[i, j], sd[i, j], ns[i, j]),
                                  "A": f"A{i}", "B": f"B{j}"}))
df = pd.concat(rows, ignore_index=True)
cell = df.groupby(["A", "B"]).y.agg(["size", "mean", "var"])
m, v, n = cell["mean"].values, cell["var"].values, cell["size"].values

print("B 수준별 A 의 단순효과")
print(f"{'B':>4s} {'A1-A0':>8s} {'F':>9s} {'df2':>8s} {'p':>9s} "
      f"{'본페로니 p':>10s}   대조(웰치 t 검정)")
for j in range(b):
    Cm = np.zeros((1, 6))
    Cm[0, j], Cm[0, 3 + j] = -1, 1
    F, q, f, p = welch_james(m, v, n, Cm)
    tt = stats.ttest_ind(df[(df.A == "A1") & (df.B == f"B{j}")].y,
                         df[(df.A == "A0") & (df.B == f"B{j}")].y,
                         equal_var=False)
    print(f"B{j:>3d} {m[3 + j] - m[j]:8.3f} {F:9.4f} {f:8.3f} {p:9.4f} "
          f"{min(1, 3 * p):10.4f}   t² = {tt.statistic**2:.4f}, p = {tt.pvalue:.4f}")

B 수준별 A 의 단순효과
   B    A1-A0         F      df2         p     본페로니 p   대조(웰치 t 검정)
B  0    1.173    1.7228   10.457    0.2174     0.6522   t² = 1.7228, p = 0.2174
B  1    0.077    0.0089   17.706    0.9258     1.0000   t² = 0.0089, p = 0.9258
B  2    7.993    6.9466   16.008    0.0180     0.0539   t² = 6.9466, p = 0.0180

칸별 n, 평균, 분산을 먼저 본다
    │
    ├─ 칸당 반복이 없다  ──→ 이원배치 웰치는 불가
    │                        (터키 비가법성 검정 또는 설계 재검토)
    │
    └─ 반복이 있다
          │
          ├─ 분산비 < 2, 균형 ──→ 표준 이원배치 F 로 충분
          │
          └─ 그 외 ──→ HC3 (기본)  또는  웰치-제임스 (칸이 작을 때)
                          │
                          └─ 유의한 교호작용 → 단순효과 + 보정

2×3 요인설계, 칸별 n = 8~12, 칸 분산 1.3~71.9 (57배 차이)
→ 등분산 가정 불가, 웰치-제임스(요한센 F 근사) 사용

  A 주효과   F(1, 21.2) =  8.01,  p = 0.010
  B 주효과   F(2, 26.5) = 24.88,  p < 0.001
  A×B       F(2, 26.5) =  3.20,  p = 0.057