선형성 확인¶

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

# 이 페이지의 모든 진단은 아래 모형 하나를 놓고 수행한다.
# 집단마다 표준편차를 1.0, 1.3, 1.6으로 다르게 주어 진단이 무엇을 잡아내는지
# (그리고 무엇을 못 잡아내는지) 볼 수 있게 했다.
rng = np.random.default_rng(42)
n = 20
data = pd.DataFrame({
    "group": np.repeat(["A", "B", "C"], n),
    "response": np.concatenate([
        rng.normal(10.0, 1.0, n),
        rng.normal(10.8, 1.3, n),
        rng.normal(12.0, 1.6, n),
    ]),
})
# 선형성은 연속형 공변량이 있을 때 비로소 문제가 된다.
# 그래서 이 페이지에서는 공변량을 하나 넣은 공분산분석 모형을 쓴다.
data["covariate"] = rng.uniform(0, 10, 3 * n)
data["response"] = data["response"] + 0.4 * data["covariate"]

model = ols("response ~ C(group) + covariate", data=data).fit()

print(data.groupby("group").response.agg(["count", "mean", "std"]).round(3))
print(f"\nF = {model.fvalue:.4f}, p = {model.f_pvalue:.4f}")
print(f"covariate 계수 = {model.params['covariate']:.4f}")

# (1) 의 네 등식을 설계행렬에 직접 내적해 확인한다.
e = model.resid
x = data["covariate"].values
print(f"\n설계행렬 열이름: {model.model.exog_names}")
print("X^T e =", np.array2string(model.model.exog.T @ e.values, precision=3))
print(f"  sum e_i            = {e.sum():+.3e}")
print(f"  sum_(i in B) e_i   = {e[data['group'] == 'B'].sum():+.3e}")
print(f"  sum_(i in C) e_i   = {e[data['group'] == 'C'].sum():+.3e}")
print(f"  sum x_i e_i        = {(x * e.values).sum():+.3e}")

# 그래서 잔차 그림의 직선 기울기는 항등적으로 0 이다.
print(f"\n잔차를 공변량에 회귀한 기울기 = {stats.linregress(x, e.values).slope:+.3e}")
print(f"잔차를 적합값에 회귀한 기울기 = "
      f"{stats.linregress(model.fittedvalues.values, e.values).slope:+.3e}")

b = model.params["covariate"]
se = model.bse["covariate"]
print(f"\ncovariate 계수 {b:.4f},  표준오차 {se:.4f},  참값 0.4")
print(f"  z = (추정 - 참)/표준오차 = {(b - 0.4) / se:+.3f}")
print(f"  95% 신뢰구간 = [{model.conf_int().loc['covariate', 0]:.4f}, "
      f"{model.conf_int().loc['covariate', 1]:.4f}]")
print(f"\n적합값의 범위 = [{model.fittedvalues.min():.4f}, {model.fittedvalues.max():.4f}]"
      f"   (서로 다른 값 {model.fittedvalues.nunique()}개)")

       count    mean    std
group                      
A         20  11.820  1.205
B         20  12.679  1.243
C         20  14.351  1.513

F = 35.4445, p = 0.0000
covariate 계수 = 0.3275

설계행렬 열이름: ['Intercept', 'C(group)[T.B]', 'C(group)[T.C]', 'covariate']
X^T e = [-2.292e-13 -2.114e-13  4.974e-14 -5.962e-13]
  sum e_i            = -2.292e-13
  sum_(i in B) e_i   = -2.114e-13
  sum_(i in C) e_i   = +4.974e-14
  sum x_i e_i        = -5.969e-13

잔차를 공변량에 회귀한 기울기 = +1.214e-15
잔차를 적합값에 회귀한 기울기 = +3.007e-15

covariate 계수 0.3275,  표준오차 0.0506,  참값 0.4
  z = (추정 - 참)/표준오차 = -1.432
  95% 신뢰구간 = [0.2260, 0.4289]

적합값의 범위 = [10.4040, 15.7564]   (서로 다른 값 60개)

import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns
from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm

# 그림에서 읽으려는 것을 먼저 수로 적어 둔다.
print(f"{'기울기':>8s} {'표준오차':>9s} {'상관 r':>8s} {'x 범위':>16s}  집단")
for name, d in data.groupby("group"):
    sl = stats.linregress(d["covariate"], d["response"])
    print(f"{sl.slope:8.4f} {sl.stderr:9.4f} {sl.rvalue:8.4f}"
          f" [{d['covariate'].min():5.2f},{d['covariate'].max():6.2f}]  {name}")

# 평행성(기울기 동질성)은 눈이 아니라 교호작용 항으로 판정한다.
m_int = ols("response ~ C(group) * covariate", data=data).fit()
print("\n교호작용 검정 (기울기가 집단마다 같은가)")
print(anova_lm(m_int, typ=2).round(4).to_string())

# 휘어짐도 항을 넣어 재 본다.
m_quad = ols("response ~ C(group) + covariate + I(covariate ** 2)", data=data).fit()
print(f"\n이차항 계수 {m_quad.params['I(covariate ** 2)']:+.5f}, "
      f"p = {m_quad.pvalues['I(covariate ** 2)']:.4f}")

# 집단마다 색을 달리해 흩뿌린다. 각 색의 점들이 직선을 이루는지,
# 그리고 그 직선들의 기울기가 서로 비슷한지를 본다. 기울기가 다르면
# 공변량과 집단 사이에 교호작용이 있는 것이라 공분산분석의 전제가 깨진다.
sns.scatterplot(data=data, x='covariate', y='response', hue='group', alpha=0.6)
plt.title("Response vs. Covariate by Group")
plt.show()

     기울기      표준오차     상관 r             x 범위  집단
  0.3567    0.0866   0.6966 [ 0.31,  8.53]  A
  0.2677    0.0760   0.6391 [ 0.23,  9.62]  B
  0.3816    0.1039   0.6544 [ 0.96,  9.69]  C

교호작용 검정 (기울기가 집단마다 같은가)
                     sum_sq    df        F  PR(>F)
C(group)            51.5285   2.0  24.6561  0.0000
covariate           42.9301   1.0  41.0835  0.0000
C(group):covariate   1.0613   2.0   0.5078  0.6046
Residual            56.4271  54.0      NaN     NaN

이차항 계수 -0.00556, p = 0.7870

import numpy as np
import pandas as pd
from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm

# 평행성 검정이 기울기 차이를 얼마나 잡아내는가. n = 20 에서 재 본다.
rng_sim = np.random.default_rng(2026)
B = 2000
grp = np.repeat(["A", "B", "C"], 20)
print(f"반복 {B}회, 집단당 n = 20, 명목수준 0.05, 몬테카를로 오차 "
      f"{np.sqrt(0.05 * 0.95 / B):.4f}")
print(f"\n{'참 기울기':>22s} {'기각률':>8s}")
for slopes in [(0.4, 0.4, 0.4), (0.4, 0.5, 0.6), (0.2, 0.4, 0.6), (0.0, 0.4, 0.8)]:
    s = dict(zip("ABC", slopes))
    hit = 0
    for _ in range(B):
        xv = rng_sim.uniform(0, 10, 60)
        yv = np.array([10.0 + s[g] * xx for g, xx in zip(grp, xv)]) \
            + rng_sim.normal(0, 1, 60)
        d = pd.DataFrame({"group": grp, "covariate": xv, "response": yv})
        m = ols("response ~ C(group) * covariate", data=d).fit()
        hit += anova_lm(m, typ=2).loc["C(group):covariate", "PR(>F)"] < 0.05
    print(f"{str(slopes):>22s} {hit / B:8.4f}")

반복 2000회, 집단당 n = 20, 명목수준 0.05, 몬테카를로 오차 0.0049

                 참 기울기      기각률
       (0.4, 0.4, 0.4)   0.0550
       (0.4, 0.5, 0.6)   0.3220
       (0.2, 0.4, 0.6)   0.8650
       (0.0, 0.4, 0.8)   1.0000

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from statsmodels.formula.api import ols

# 그림에서 읽으려는 것을 먼저 수로 적어 둔다.
def read_curve(v, e, label):
    """사분위 구간별 잔차 평균과 이차항 검정으로 휘어짐을 잰다."""
    q = np.quantile(v, [0, .25, .5, .75, 1.0])
    b = np.clip(np.digitize(v, q[1:-1]), 0, 3)
    means = " ".join(f"{e[b == j].mean():+.3f}" for j in range(4))
    m = ols("e ~ v + I(v ** 2)", data=pd.DataFrame({"e": e, "v": v})).fit()
    print(f"{label}  구간 평균 {means}   "
          f"이차항 {m.params['I(v ** 2)']:+.5f} (p = {m.pvalues['I(v ** 2)']:.4f})")

read_curve(model.fittedvalues.values, model.resid.values, "적합값에 대해     ")
read_curve(data["covariate"].values, model.resid.values, "공변량에 대해     ")

# 비교: 참 관계가 이차인 자료를 선형으로 적합하면 두 그림이 어떻게 달라지는가
dq = data.copy()
dq["response"] = (10.0 + 0.4 * dq["covariate"]
                  + 0.10 * (dq["covariate"] - 5) ** 2
                  + np.random.default_rng(2026).normal(0, 1, 60))
mq = ols("response ~ C(group) + covariate", data=dq).fit()
print()
read_curve(mq.fittedvalues.values, mq.resid.values, "[이차] 적합값에 대해")
read_curve(dq["covariate"].values, mq.resid.values, "[이차] 공변량에 대해")
print(f"\n[이차] 적합값과 공변량의 상관 = "
      f"{np.corrcoef(mq.fittedvalues.values, dq['covariate'].values)[0, 1]:.4f}")

# 잔차 대 적합값 그림은 진단의 출발점이다. 점들이 0 선 둘레에 폭을 일정하게
# 유지하며 흩어져 있으면 좋다. 깔때기 모양이면 등분산이 깨진 것이고, 굽은
# 모양이면 모형이 놓친 구조가 남아 있다는 뜻이다.
plt.scatter(model.fittedvalues, model.resid, alpha=0.6)
plt.axhline(y=0, color='r', linestyle='--')
plt.xlabel("Fitted Values")
plt.ylabel("Residuals")
plt.title("Residuals vs. Fitted Values")
plt.show()

적합값에 대해       구간 평균 +0.015 +0.051 -0.090 +0.023   이차항 +0.01221 (p = 0.8597)
공변량에 대해       구간 평균 -0.058 +0.038 +0.247 -0.227   이차항 -0.00532 (p = 0.7877)

[이차] 적합값에 대해  구간 평균 +0.274 -0.224 -0.068 +0.017   이차항 +0.30729 (p = 0.0245)
[이차] 공변량에 대해  구간 평균 +0.365 -0.494 -0.305 +0.434   이차항 +0.06945 (p = 0.0004)

[이차] 적합값과 공변량의 상관 = 0.9257

import warnings
warnings.filterwarnings("ignore")

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

rng = np.random.default_rng(7001)
B = 4_000
print("y = (x-1.5)² + ε,  집단 효과는 실제로 0 (명목 0.05)")
print(f"{'공변량 분포':>26s} {'선형 적합 오류율':>14s} {'이차항 포함':>12s}")
for lab, same in [("모든 집단 x ~ U(0,3) (동일)", True),
                  ("집단마다 x 의 중심이 다름", False)]:
    a = b = 0
    for _ in range(B):
        rows = []
        for g in range(3):
            x = rng.uniform(0, 3, 30) if same else rng.uniform(g * 0.6,
                                                               g * 0.6 + 2, 30)
            y = (x - 1.5)**2 + rng.normal(0, 1.0, 30)
            rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
        df = pd.concat(rows, ignore_index=True)
        t1 = sm.stats.anova_lm(ols("y ~ x + C(g)", data=df).fit(), typ=2)
        t2 = sm.stats.anova_lm(ols("y ~ x + I(x**2) + C(g)", data=df).fit(),
                               typ=2)
        a += t1.loc["C(g)", "PR(>F)"] < 0.05
        b += t2.loc["C(g)", "PR(>F)"] < 0.05
    print(f"{lab:>26s} {a / B:14.4f} {b / B:12.4f}")

y = (x-1.5)² + ε,  집단 효과는 실제로 0 (명목 0.05)
                    공변량 분포      선형 적합 오류율       이차항 포함
     모든 집단 x ~ U(0,3) (동일)         0.0530       0.0508
           집단마다 x 의 중심이 다름         0.1650       0.0498

rng = np.random.default_rng(7005)
B = 3_000
print("집단별 추정 효과의 평균 (참값 모두 0), x 중심이 집단마다 다름")
print(f"{'적합':>18s} {'G1-G0':>9s} {'G2-G0':>9s}")
for lab, form in [("선형 (y ~ x + g)", "y ~ x + C(g)"),
                  ("이차항 포함", "y ~ x + I(x**2) + C(g)")]:
    a, b = [], []
    for _ in range(B):
        rows = []
        for g in range(3):
            x = rng.uniform(g * 0.6, g * 0.6 + 2, 30)
            y = (x - 1.5)**2 + rng.normal(0, 1.0, 30)
            rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
        df = pd.concat(rows, ignore_index=True)
        f = ols(form, data=df).fit()
        a.append(f.params["C(g)[T.G1]"])
        b.append(f.params["C(g)[T.G2]"])
    print(f"{lab:>18s} {np.mean(a):9.4f} {np.mean(b):9.4f}")

집단별 추정 효과의 평균 (참값 모두 0), x 중심이 집단마다 다름
                적합     G1-G0     G2-G0
    선형 (y ~ x + g)   -0.3599    0.0039
            이차항 포함    0.0059    0.0068

import warnings
warnings.filterwarnings("ignore")

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

rng = np.random.default_rng(7009)
rows = []
for g in range(3):
    x = rng.uniform(0, 3, 40)
    y = 2.0 + 1.5 * (x - 1.5)**2 + 0.8 * g + rng.normal(0, 0.8, 40)
    rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
df = pd.concat(rows, ignore_index=True)

fit = ols("y ~ x + C(g)", data=df).fit()
b = fit.params["x"]
df["cr"] = fit.resid + b * df["x"]          # 성분 + 잔차
print(f"선형 적합의 x 계수 = {b:.4f}")

print("\nx 를 6구간으로 나눠 성분+잔차의 평균을 본다")
df["bin"] = pd.cut(df.x, 6)
t = df.groupby("bin", observed=True).agg(
    x평균=("x", "mean"), 성분잔차평균=("cr", "mean"),
    잔차평균=("y", lambda s: fit.resid[s.index].mean()), n=("x", "size"))
print(t.round(4).to_string())

q = np.polyfit(df.x, df.cr, 2)
print(f"\n성분+잔차에 2차 곡선을 맞추면: "
      f"{q[0]:.4f}·x² + {q[1]:.4f}·x + {q[2]:.4f}")
print("  (참 곡률 계수는 1.5)")

f2 = ols("y ~ x + I(x**2) + C(g)", data=df).fit()
print(f"\n이차항 모형의 x² 계수 = {f2.params['I(x ** 2)']:.4f} "
      f"(p = {f2.pvalues['I(x ** 2)']:.6f})")
print(f"집단 효과 추정: 선형 {fit.params['C(g)[T.G2]']:.4f} / "
      f"이차항 포함 {f2.params['C(g)[T.G2]']:.4f}  (참값 1.6)")

선형 적합의 x 계수 = 0.2489

x 를 6구간으로 나눠 성분+잔차의 평균을 본다
                   x평균  성분잔차평균    잔차평균   n
bin                                       
(0.104, 0.588]  0.3747  1.0736  0.9804  21
(0.588, 1.069]  0.8376 -0.0565 -0.2650  16
(1.069, 1.55]   1.3165 -0.2544 -0.5821  22
(1.55, 2.031]   1.7797 -0.4981 -0.9410  22
(2.031, 2.513]  2.2216  0.0673 -0.4857  15
(2.513, 2.994]  2.7283  1.6975  1.0185  24

성분+잔차에 2차 곡선을 맞추면: 1.3570·x² + -3.9902·x + 2.3871
  (참 곡률 계수는 1.5)

이차항 모형의 x² 계수 = 1.3602 (p = 0.000000)
집단 효과 추정: 선형 1.6839 / 이차항 포함 1.5862  (참값 1.6)

x 평균   0.37   0.84   1.32   1.78   2.22   2.73
성분잔차  1.07  -0.06  -0.25  -0.50   0.07   1.70
          ↓                    ↓                ↑
        높음 ────── 낮아짐 ──── 최저 ──── 다시 높아짐

import warnings
warnings.filterwarnings("ignore")

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

def reset_p(fit, df, power=3):
    """램지 RESET 검정: 적합값의 거듭제곱을 넣어 본다."""
    d2 = df.copy()
    for p in range(2, power + 1):
        d2[f"yh{p}"] = fit.fittedvalues**p
    big = ols("y ~ x + C(g) + "
              + " + ".join(f"yh{p}" for p in range(2, power + 1)),
              data=d2).fit()
    idx = list(big.params.index)
    R = np.zeros((power - 1, len(idx)))
    for r, p in enumerate(range(2, power + 1)):
        R[r, idx.index(f"yh{p}")] = 1
    return float(big.f_test(R).pvalue)

rng = np.random.default_rng(7003)
B = 3_000
print("k=3, 집단당 n=30, 오차 SD 1.0, 명목 0.05")
print(f"{'참 관계':>18s} {'RESET 탐지율':>12s} {'이차항 t 검정':>13s}")
for lab, f in [("선형 (H0 참)", lambda x: 2 * x),
               ("약한 곡률", lambda x: 2 * x + 0.3 * (x - 1.5)**2),
               ("뚜렷한 곡률", lambda x: 2 * x + 1.0 * (x - 1.5)**2),
               ("로그형", lambda x: 4 * np.log(x + 0.5))]:
    a = b = 0
    for _ in range(B):
        rows = []
        for g in range(3):
            x = rng.uniform(0, 3, 30)
            y = f(x) + rng.normal(0, 1.0, 30)
            rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
        df = pd.concat(rows, ignore_index=True)
        fit = ols("y ~ x + C(g)", data=df).fit()
        a += reset_p(fit, df) < 0.05
        f2 = ols("y ~ x + I(x**2) + C(g)", data=df).fit()
        b += f2.pvalues["I(x ** 2)"] < 0.05
    print(f"{lab:>18s} {a / B:12.4f} {b / B:13.4f}")

k=3, 집단당 n=30, 오차 SD 1.0, 명목 0.05
              참 관계    RESET 탐지율      이차항 t 검정
         선형 (H0 참)       0.0410        0.0450
             약한 곡률       0.3310        0.4590
            뚜렷한 곡률       0.9997        1.0000
               로그형       0.9710        0.9833

import warnings
warnings.filterwarnings("ignore")

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

forms = [("선형", "y ~ x + C(g)"),
         ("이차항", "y ~ x + I(x**2) + C(g)"),
         ("자연 스플라인(df=4)", "y ~ cr(x, df=4) + C(g)"),
         ("x 를 4분위로 층화", "y ~ C(xq) + C(g)")]

rng = np.random.default_rng(7011)
B = 2_000
print("참 관계는 로그형 y = 4·log(x), 집단 효과 = (0, 0.4, 0.8), 집단당 n=40")
print(f"{'방법':>22s} {'검정력':>8s} {'잔차 SD':>9s}")
res = {k: [0, []] for k, _ in forms}
for _ in range(B):
    rows = []
    for g in range(3):
        x = rng.uniform(0.2, 5, 40)
        y = 4 * np.log(x) + 0.4 * g + rng.normal(0, 1.0, 40)
        rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
    df = pd.concat(rows, ignore_index=True)
    df["xq"] = pd.qcut(df.x, 4, labels=False).astype(str)
    for k, form in forms:
        f = ols(form, data=df).fit()
        t = sm.stats.anova_lm(f, typ=2)
        res[k][0] += t.loc["C(g)", "PR(>F)"] < 0.05
        res[k][1].append(np.sqrt(f.mse_resid))
for k, (c, s) in res.items():
    print(f"{k:>22s} {c / B:8.4f} {np.mean(s):9.4f}")

참 관계는 로그형 y = 4·log(x), 집단 효과 = (0, 0.4, 0.8), 집단당 n=40
                    방법      검정력     잔차 SD
                    선형   0.5915    1.4243
                   이차항   0.8215    1.0975
         자연 스플라인(df=4)   0.8530    1.0451
           x 를 4분위로 층화   0.5250    1.5074

import warnings
warnings.filterwarnings("ignore")

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

rng = np.random.default_rng(7008)
B = 3_000
print("참 모형: y = (1 + g·Δ)·x + ε,  절편은 세 집단 모두 0 (집단 효과 없음)")
print(f"{'Δ':>5s} {'공통기울기 C(g) 오류율':>20s} {'교호작용 탐지율':>14s} "
      f"{'G2-G0 @x=0.5':>13s} {'@x=1.5':>9s} {'@x=2.5':>9s}")
for D in [0.0, 0.5, 1.0]:
    a = c = 0
    e = [[], [], []]
    for _ in range(B):
        rows = []
        for g in range(3):
            x = rng.uniform(0, 3, 30)
            y = (1 + g * D) * x + rng.normal(0, 1.5, 30)
            rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
        df = pd.concat(rows, ignore_index=True)
        t1 = sm.stats.anova_lm(ols("y ~ x + C(g)", data=df).fit(), typ=2)
        a += t1.loc["C(g)", "PR(>F)"] < 0.05
        f2 = ols("y ~ x * C(g)", data=df).fit()
        c += sm.stats.anova_lm(f2, typ=2).loc["x:C(g)", "PR(>F)"] < 0.05
        for i, xv in enumerate([0.5, 1.5, 2.5]):
            p0 = f2.predict(pd.DataFrame({"x": [xv], "g": ["G0"]}))[0]
            p2 = f2.predict(pd.DataFrame({"x": [xv], "g": ["G2"]}))[0]
            e[i].append(p2 - p0)
    print(f"{D:5.1f} {a / B:20.4f} {c / B:14.4f} "
          f"{np.mean(e[0]):13.4f} {np.mean(e[1]):9.4f} {np.mean(e[2]):9.4f}")

참 모형: y = (1 + g·Δ)·x + ε,  절편은 세 집단 모두 0 (집단 효과 없음)
    Δ       공통기울기 C(g) 오류율       교호작용 탐지율  G2-G0 @x=0.5    @x=1.5    @x=2.5
  0.0               0.0520         0.0490       -0.0069   -0.0063   -0.0058
  0.5               0.9173         0.4633        0.4997    1.4986    2.4975
  1.0               1.0000         0.9677        1.0047    3.0053    5.0058

import warnings
warnings.filterwarnings("ignore")

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

rng = np.random.default_rng(7007)
B = 3_000
print("참 관계는 이차 y=(x-1.5)², 집단 효과 0, 선형 공분산분석으로 적합")
print(f"{'겹침':>10s} {'조정평균 차이 추정':>16s} {'표준오차':>10s} {'오류율':>8s}")
for lab, (r0, r1) in [("완전 겹침", ((0, 3), (0, 3))),
                      ("절반 겹침", ((0, 2), (1, 3))),
                      ("접점만", ((0, 1.5), (1.5, 3))),
                      ("완전 분리", ((0, 1.2), (1.8, 3)))]:
    d, se, hit = [], [], 0
    for _ in range(B):
        rows = []
        for g, (lo, hi) in enumerate([r0, r1]):
            x = rng.uniform(lo, hi, 40)
            y = (x - 1.5)**2 + rng.normal(0, 1.0, 40)
            rows.append(pd.DataFrame({"y": y, "x": x, "g": f"G{g}"}))
        df = pd.concat(rows, ignore_index=True)
        f = ols("y ~ x + C(g)", data=df).fit()
        d.append(f.params["C(g)[T.G1]"])
        se.append(f.bse["C(g)[T.G1]"])
        hit += f.pvalues["C(g)[T.G1]"] < 0.05
    print(f"{lab:>10s} {np.mean(d):16.4f} {np.mean(se):10.4f} {hit / B:8.4f}")

참 관계는 이차 y=(x-1.5)², 집단 효과 0, 선형 공분산분석으로 적합
        겹침       조정평균 차이 추정       표준오차      오류율
     완전 겹침          -0.0112     0.2691   0.0527
     절반 겹침          -0.0066     0.3551   0.0253
       접점만           0.0110     0.5434   0.0347
     완전 분리          -0.0265     0.7445   0.0373

1. 집단별 x 의 분포를 그린다  ──→ 겹치는가?
     ↓  (겹치지 않으면 여기서 멈추고 범위를 제한)
2. y ~ x * C(g) 를 적합한다   ──→ 기울기가 같은가?
     ↓  (다르면 교호작용을 유지하고 특정 x 에서 보고)
3. 부분잔차 그림을 본다       ──→ 곡률이 있는가?
     ↓  (있으면 이차항이나 스플라인)
4. y ~ x + C(g) 로 조정평균을 보고