상관의 이해¶

import matplotlib.pyplot as plt
import numpy as np
from scipy import stats

np.random.seed(0)

def generate_samples(mu_1, mu_2, sigma_1, sigma_2, rho, n):
    """이변량 정규분포에서 표본을 뽑는다.

    상관계수 rho 는 공분산행렬의 비대각 원소로 들어간다.
    공분산 = rho * sigma_1 * sigma_2 이므로, 표준편차를 1 로 두면
    공분산이 곧 상관계수가 된다.

    돌려주는 것은 모양 (n, 2) 인 배열이다.
    """
    covariance_matrix = [
        [sigma_1**2, rho * sigma_1 * sigma_2],
        [rho * sigma_1 * sigma_2, sigma_2**2]
    ]
    return stats.multivariate_normal([mu_1, mu_2], covariance_matrix).rvs(size=n)

def plot_correlations():
    """rho 를 0 에서 0.95 까지 키우며 점구름이 어떻게 좁아지는지 본다.

    0.4 와 0.6 의 그림 차이가 생각보다 작다는 점을 눈여겨볼 만하다.
    상관계수는 눈에 보이는 것보다 느리게 움직인다.
    """
    fig, axes = plt.subplots(1, 6, figsize=(15, 4))
    correlation_coefficients = (0.00, 0.40, 0.60, 0.80, 0.90, 0.95)

    for ax, rho in zip(axes, correlation_coefficients):
        xy = generate_samples(0, 0, 1, 1, rho, 100)
        ax.plot(xy[:, 0], xy[:, 1], 'ok')
        ax.set_title(f'rho = {rho}')
        ax.axis('off')
        ax.axis('equal')
        for loc in ('left', 'right', 'top', 'bottom'):
            ax.spines[loc].set_visible(False)
    plt.show()

if __name__ == "__main__":
    plot_correlations()

import numpy as np
from scipy import stats

def generate_samples(rho, n=100):
    return stats.multivariate_normal([0, 0], [[1, rho], [rho, 1]]).rvs(size=n)

rhos = (0.00, 0.40, 0.60, 0.80, 0.90, 0.95)

# (1)(2) 축비는 고윳값의 비의 제곱근이다
print("rho    축비   (이론)")
for rho in rhos:
    print(f"{rho:5.2f}  {np.sqrt((1 + abs(rho)) / (1 - abs(rho))):6.3f}")
for a, b in [(0.00, 0.05), (0.40, 0.45), (0.90, 0.95)]:
    f = lambda r: np.sqrt((1 + r) / (1 - r))
    print(f"  {a:.2f} -> {b:.2f}:  {f(a):.3f} -> {f(b):.3f}"
          f"   폭 {f(b) - f(a):+.3f}")

# (3) 표본 r 과 표본 축비
np.random.seed(0)
print("\n rho    표본 r      차      SE   차/SE   표본축비  이론축비")
for rho in rhos:
    xy = generate_samples(rho)
    r = np.corrcoef(xy.T)[0, 1]
    ev = np.linalg.eigvalsh(np.cov(xy.T))
    se = (1 - rho**2) / np.sqrt(len(xy))
    print(f"{rho:5.2f}  {r:+.4f}  {r - rho:+.4f}  {se:.4f}  {(r - rho) / se:+6.2f}"
          f"    {np.sqrt(ev[1] / ev[0]):6.3f}   {np.sqrt((1 + abs(rho)) / (1 - abs(rho))):6.3f}")

rho    축비   (이론)
 0.00   1.000
 0.40   1.528
 0.60   2.000
 0.80   3.000
 0.90   4.359
 0.95   6.245
  0.00 -> 0.05:  1.000 -> 1.051   폭 +0.051
  0.40 -> 0.45:  1.528 -> 1.624   폭 +0.096
  0.90 -> 0.95:  4.359 -> 6.245   폭 +1.886

 rho    표본 r      차      SE   차/SE   표본축비  이론축비
 0.00  -0.0240  -0.0240  0.1000   -0.24     1.025    1.000
 0.40  +0.4377  +0.0377  0.0840   +0.45     1.611    1.528
 0.60  +0.6116  +0.0116  0.0640   +0.18     2.045    2.000
 0.80  +0.7572  -0.0428  0.0360   -1.19     2.690    3.000
 0.90  +0.8858  -0.0142  0.0190   -0.75     4.067    4.359
 0.95  +0.9591  +0.0091  0.0098   +0.93     6.923    6.245

import matplotlib.pyplot as plt
import numpy as np
from scipy import stats

np.random.seed(0)

def generate_samples(mu_1, mu_2, sigma_1, sigma_2, rho, n):
    covariance_matrix = [
        [sigma_1**2, rho * sigma_1 * sigma_2],
        [rho * sigma_1 * sigma_2, sigma_2**2]
    ]
    return stats.multivariate_normal([mu_1, mu_2], covariance_matrix).rvs(size=n)

def plot_negative_correlations():
    """같은 일을 음의 상관에서 되풀이한다. 모양은 같고 기울기만 뒤집힌다."""
    fig, axes = plt.subplots(1, 6, figsize=(15, 4))
    correlation_coefficients = (0.00, -0.40, -0.60, -0.80, -0.90, -0.95)

    for ax, rho in zip(axes, correlation_coefficients):
        xy = generate_samples(0, 0, 1, 1, rho, 100)
        ax.plot(xy[:, 0], xy[:, 1], 'ok')
        ax.set_title(f'rho = {rho}')
        ax.axis('off')
        ax.axis('equal')
        for loc in ('left', 'right', 'top', 'bottom'):
            ax.spines[loc].set_visible(False)
    plt.show()

if __name__ == "__main__":
    plot_negative_correlations()

import numpy as np
from scipy import stats

def generate_samples(rho, n=100):
    return stats.multivariate_normal([0, 0], [[1, rho], [rho, 1]]).rvs(size=n)

pos = (0.00, 0.40, 0.60, 0.80, 0.90, 0.95)
neg = tuple(-r for r in pos)

np.random.seed(0)
A = [generate_samples(r) for r in pos]
np.random.seed(0)
B = [generate_samples(r) for r in neg]

print("|rho|    양수 판 r             음수 판 r             합        90도 회전인가")
for rho, a, b in zip(pos, A, B):
    ra, rb = np.corrcoef(a.T)[0, 1], np.corrcoef(b.T)[0, 1]
    rot = np.allclose(b[:, 0], a[:, 1]) and np.allclose(b[:, 1], -a[:, 0])
    print(f"{rho:5.2f}  {ra:+.16f}  {rb:+.16f}  {ra + rb:+.1e}   {rot}")

|rho|    양수 판 r             음수 판 r             합        90도 회전인가
 0.00  -0.0240052014940784  -0.0240052014940784  -4.8e-02   False
 0.40  +0.4376998281800900  -0.4376998281800900  +0.0e+00   True
 0.60  +0.6116463353944882  -0.6116463353944882  +0.0e+00   True
 0.80  +0.7571693721126749  -0.7571693721126749  +0.0e+00   True
 0.90  +0.8858200908150257  -0.8858200908150257  +0.0e+00   True
 0.95  +0.9590898844142728  -0.9590898844142730  -2.2e-16   True

import matplotlib.pyplot as plt
import numpy as np
from scipy import stats

np.random.seed(0)

def generate_samples(mu_1, mu_2, sigma_1, sigma_2, rho, n):
    covariance_matrix = [
        [sigma_1**2, rho * sigma_1 * sigma_2],
        [rho * sigma_1 * sigma_2, sigma_2**2]
    ]
    return stats.multivariate_normal([mu_1, mu_2], covariance_matrix).rvs(size=n)

def plot_all_correlations():
    """-0.95 부터 0.95 까지 열한 칸에 늘어놓아 한눈에 견준다."""
    fig, axes = plt.subplots(1, 11, figsize=(20, 2))
    rhos = (-0.95, -0.90, -0.80, -0.60, -0.40,
             0.00,  0.40,  0.60,  0.80,  0.90, 0.95)

    for ax, rho in zip(axes, rhos):
        xy = generate_samples(0, 0, 1, 1, rho, 100)
        ax.plot(xy[:, 0], xy[:, 1], 'ok')
        ax.set_title(f'rho = {rho}')
        ax.axis('off')
        ax.axis('equal')
        for loc in ('left', 'right', 'top', 'bottom'):
            ax.spines[loc].set_visible(False)
    plt.show()

if __name__ == "__main__":
    plot_all_correlations()

import numpy as np
from scipy import stats

rhos = (-0.95, -0.90, -0.80, -0.60, -0.40,
         0.00,  0.40,  0.60,  0.80,  0.90, 0.95)

np.random.seed(0)
print("  rho     표본 r    rho^2")
for rho in rhos:
    xy = stats.multivariate_normal([0, 0], [[1, rho], [rho, 1]]).rvs(size=100)
    print(f"{rho:+6.2f}   {np.corrcoef(xy.T)[0, 1]:+.3f}   {rho**2:.2f}")

# 이웃한 두 판의 순서가 뒤집힐 확률 (피셔 z)
se = 1 / np.sqrt(100 - 3)
print("\n이웃 판의 순서가 뒤집힐 확률")
for a, b in [(0.00, 0.40), (0.40, 0.60), (0.60, 0.80), (0.80, 0.90), (0.90, 0.95)]:
    d = (np.arctanh(b) - np.arctanh(a)) / (se * np.sqrt(2))
    print(f"  P(r_{a:.2f} > r_{b:.2f}) = {stats.norm.sf(d):.4f}")

  rho     표본 r    rho^2
 -0.95   -0.951   0.90
 -0.90   -0.909   0.81
 -0.80   -0.807   0.64
 -0.60   -0.526   0.36
 -0.40   -0.338   0.16
 +0.00   +0.107   0.00
 +0.40   +0.340   0.16
 +0.60   +0.653   0.36
 +0.80   +0.805   0.64
 +0.90   +0.879   0.81
 +0.95   +0.956   0.90

이웃 판의 순서가 뒤집힐 확률
  P(r_0.00 > r_0.40) = 0.0016
  P(r_0.40 > r_0.60) = 0.0303
  P(r_0.60 > r_0.80) = 0.0024
  P(r_0.80 > r_0.90) = 0.0046
  P(r_0.90 > r_0.95) = 0.0061

import matplotlib.pyplot as plt
import pandas as pd

def plot_height_weight_scatter():
    # openintro의 bdims 자료: 성인 507명의 신체 치수.
    # hgt(cm), wgt(kg), sex(1 = 남성, 0 = 여성) 열을 쓴다.
    url = ("https://raw.githubusercontent.com/vincentarelbundock/Rdatasets/"
           "1dcc2bf5f955cc1224a3e1307256e1fe86b68dae/csv/openintro/bdims.csv")
    data = pd.read_csv(url).rename(columns={"hgt": "Height", "wgt": "Weight"})
    data["Gender"] = data["sex"].map({1: "Male", 0: "Female"})
    filtered = data[data.Gender == "Male"][:300]

    fig, ax = plt.subplots(figsize=(10, 10))
    ax.plot(filtered.Height, filtered.Weight, '.k')
    ax.set_xlabel('Height', fontsize=15)
    ax.set_ylabel('Weight', fontsize=15)
    ax.spines['top'].set_visible(False)
    ax.spines['right'].set_visible(False)
    plt.show()

if __name__ == "__main__":
    plot_height_weight_scatter()

import numpy as np
import pandas as pd

url = ("https://raw.githubusercontent.com/vincentarelbundock/Rdatasets/"
       "1dcc2bf5f955cc1224a3e1307256e1fe86b68dae/csv/openintro/bdims.csv")
data = pd.read_csv(url).rename(columns={"hgt": "Height", "wgt": "Weight"})
data["Gender"] = data["sex"].map({1: "Male", 0: "Female"})

# (2) [:300] 은 몇 개를 남기는가
male = data[data.Gender == "Male"]
print(f"남성 전체 {len(male)}명   [:300] 뒤 {len(male[:300])}명"
      f"   잘린 수 {len(male) - len(male[:300])}")

# (1) 남성의 상관, 설명된 변동, 기울기
r = male.Height.corr(male.Weight)
b = r * male.Weight.std() / male.Height.std()
print(f"\nr = {r:.4f}   r^2 = {r**2:.4f}   기울기 = {b:.4f} kg/cm")
print(f"키 {male.Height.min()}~{male.Height.max()} cm   "
      f"몸무게 {male.Weight.min()}~{male.Weight.max()} kg")
print(f"sd  키 {male.Height.std():.2f}   몸무게 {male.Weight.std():.2f}")

# (3) 집단별과 전체
print()
for label, sub in [("남성", male),
                   ("여성", data[data.Gender == "Female"]),
                   ("전체", data)]:
    rr = sub.Height.corr(sub.Weight)
    print(f"{label}: n = {len(sub):3d}   r = {rr:.4f}   r^2 = {rr**2:.4f}"
          f"   키평균 {sub.Height.mean():6.2f}   몸무게평균 {sub.Weight.mean():5.2f}")

# 성별 평균을 뺀 잔차의 상관 = 성별을 통제한 편상관
res_h = data.Height - data.groupby("Gender").Height.transform("mean")
res_w = data.Weight - data.groupby("Gender").Weight.transform("mean")
c = data[["Height", "Weight"]].assign(S=data.sex).corr()
rxy, rxs, rws = c.loc["Height", "Weight"], c.loc["Height", "S"], c.loc["Weight", "S"]
print(f"\n잔차의 상관   {res_h.corr(res_w):.4f}")
print(f"편상관 공식   "
      f"{(rxy - rxs * rws) / np.sqrt((1 - rxs**2) * (1 - rws**2)):.4f}")

남성 전체 247명   [:300] 뒤 247명   잘린 수 0

r = 0.5347   r^2 = 0.2859   기울기 = 0.7826 kg/cm
키 157.2~198.1 cm   몸무게 53.9~116.4 kg
sd  키 7.18   몸무게 10.51

남성: n = 247   r = 0.5347   r^2 = 0.2859   키평균 177.75   몸무게평균 78.14
여성: n = 260   r = 0.4311   r^2 = 0.1858   키평균 164.87   몸무게평균 60.60
전체: n = 507   r = 0.7173   r^2 = 0.5145   키평균 171.14   몸무게평균 69.15

잔차의 상관   0.4863
편상관 공식   0.4863

import numpy as np
from scipy import stats

x1 = np.array([10, 8, 13, 9, 11, 14, 6, 4, 12, 7, 5], float)
y1 = np.array([8.04, 6.95, 7.58, 8.81, 8.33, 9.96, 7.24, 4.26, 10.84, 4.82, 5.68])
y2 = np.array([9.14, 8.14, 8.74, 8.77, 9.26, 8.10, 6.13, 3.10, 9.13, 7.26, 4.74])
y3 = np.array([7.46, 6.77, 12.74, 7.11, 7.81, 8.84, 6.08, 5.39, 8.15, 6.42, 5.73])
x4 = np.array([8, 8, 8, 8, 8, 8, 8, 19, 8, 8, 8], float)
y4 = np.array([6.58, 5.76, 7.71, 8.84, 8.47, 7.04, 5.25, 12.50, 5.56, 7.91, 6.89])

for i, (x, y) in enumerate([(x1, y1), (x1, y2), (x1, y3), (x4, y4)], 1):
    b = np.polyfit(x, y, 1)
    print(f"  {i}: 평균 x={x.mean():.2f}, 평균 y={y.mean():.2f}, "
          f"r={np.corrcoef(x, y)[0, 1]:.4f}, 기울기={b[0]:.3f}, "
          f"절편={b[1]:.3f}, 스피어만={stats.spearmanr(x, y).statistic:+.4f}")

  1: 평균 x=9.00, 평균 y=7.50, r=0.8164, 기울기=0.500, 절편=3.000, 스피어만=+0.8182
  2: 평균 x=9.00, 평균 y=7.50, r=0.8162, 기울기=0.500, 절편=3.001, 스피어만=+0.6909
  3: 평균 x=9.00, 평균 y=7.50, r=0.8163, 기울기=0.500, 절편=3.002, 스피어만=+0.9909
  4: 평균 x=9.00, 평균 y=7.50, r=0.8165, 기울기=0.500, 절편=3.002, 스피어만=+0.5000

import numpy as np

print("r_xy=a, r_yz=b 일 때 r_xz 의 가능 범위")
print(f"{'r_xy':>6s} {'r_yz':>6s} {'하한':>8s} {'상한':>8s} {'음수 가능?':>10s}")
for a, b in [(0.5, 0.5), (0.7, 0.7), (0.8, 0.8), (0.9, 0.9),
             (0.71, 0.71), (0.6, 0.9)]:
    s = np.sqrt((1 - a**2) * (1 - b**2))
    lo, hi = a * b - s, min(a * b + s, 1.0)
    print(f"{a:6.2f} {b:6.2f} {lo:8.4f} {hi:8.4f} "
          f"{'예' if lo < 0 else '아니오':>10s}")

r_xy=a, r_yz=b 일 때 r_xz 의 가능 범위
  r_xy   r_yz       하한       상한     음수 가능?
  0.50   0.50  -0.5000   1.0000          예
  0.70   0.70  -0.0200   1.0000          예
  0.80   0.80   0.2800   1.0000        아니오
  0.90   0.90   0.6200   1.0000        아니오
  0.71   0.71   0.0082   1.0000        아니오
  0.60   0.90   0.1913   0.8887        아니오

rng = np.random.default_rng(24001)
print("\nr_xy=0.7, r_yz=0.7 을 고정하고 r_xz 를 바꿔 본다")
for target in [-0.2, 0.0, 0.5, 0.95]:
    R = np.array([[1, 0.7, target], [0.7, 1, 0.7], [target, 0.7, 1]])
    ev = np.linalg.eigvalsh(R)
    if ev.min() > 0:
        d = rng.multivariate_normal([0, 0, 0], R, size=200_000)
        print(f"  목표 r_xz={target:+.2f}: 가능. 실현값 "
              f"r_xy={np.corrcoef(d[:, 0], d[:, 1])[0, 1]:+.4f} "
              f"r_yz={np.corrcoef(d[:, 1], d[:, 2])[0, 1]:+.4f} "
              f"r_xz={np.corrcoef(d[:, 0], d[:, 2])[0, 1]:+.4f}")
    else:
        print(f"  목표 r_xz={target:+.2f}: 불가능 "
              f"(최소 고윳값 {ev.min():.4f} < 0)")

r_xy=0.7, r_yz=0.7 을 고정하고 r_xz 를 바꿔 본다
  목표 r_xz=-0.20: 불가능 (최소 고윳값 -0.0950 < 0)
  목표 r_xz=+0.00: 가능. 실현값 r_xy=+0.6996 r_yz=+0.6995 r_xz=-0.0013
  목표 r_xz=+0.50: 가능. 실현값 r_xy=+0.7013 r_yz=+0.7008 r_xz=+0.4999
  목표 r_xz=+0.95: 가능. 실현값 r_xy=+0.7001 r_yz=+0.7001 r_xz=+0.9505

Y = X 성분 + 잔여 성분
Z 는 Y 의 "잔여 성분"과만 연결될 수 있다
  → X 와 Z 는 무관해진다

import warnings
warnings.filterwarnings("ignore")

import numpy as np

rng = np.random.default_rng(24002)
n = 200_000
for r in [0.3, 0.6, 0.9]:
    z = rng.standard_normal((n, 2))
    x = z[:, 0]
    y = r * z[:, 0] + np.sqrt(1 - r**2) * z[:, 1]
    print(f"\n  r = {r}")
    for lo, hi, lab in [(2, 99, "상위 2.3% (z>2)"), (1, 2, "z ∈ (1,2)"),
                        (-99, -2, "하위 2.3% (z<-2)")]:
        m = (x > lo) & (x < hi)
        print(f"    {lab:>18s}: 평균 x={x[m].mean():+.3f} → "
              f"평균 y={y[m].mean():+.3f} (예측 r·x̄={r * x[m].mean():+.3f})")

  r = 0.3
         상위 2.3% (z>2): 평균 x=+2.376 → 평균 y=+0.720 (예측 r·x̄=+0.713)
             z ∈ (1,2): 평균 x=+1.382 → 평균 y=+0.416 (예측 r·x̄=+0.414)
        하위 2.3% (z<-2): 평균 x=-2.376 → 평균 y=-0.719 (예측 r·x̄=-0.713)

  r = 0.6
         상위 2.3% (z>2): 평균 x=+2.369 → 평균 y=+1.438 (예측 r·x̄=+1.422)
             z ∈ (1,2): 평균 x=+1.385 → 평균 y=+0.830 (예측 r·x̄=+0.831)
        하위 2.3% (z<-2): 평균 x=-2.378 → 평균 y=-1.454 (예측 r·x̄=-1.427)

  r = 0.9
         상위 2.3% (z>2): 평균 x=+2.381 → 평균 y=+2.143 (예측 r·x̄=+2.143)
             z ∈ (1,2): 평균 x=+1.385 → 평균 y=+1.246 (예측 r·x̄=+1.246)
        하위 2.3% (z<-2): 평균 x=-2.367 → 평균 y=-2.132 (예측 r·x̄=-2.130)

from scipy import stats

rng = np.random.default_rng(24002)
n, r = 100_000, 0.7
z = rng.standard_normal((n, 2))
pre = z[:, 0] * 10 + 100
post = (r * z[:, 0] + np.sqrt(1 - r**2) * z[:, 1]) * 10 + 100

print("아무 처치도 하지 않았다 (사전·사후 상관 0.7, 평균 100, SD 10)")
for lab, m in [("하위 20% 선발", pre < np.quantile(pre, 0.2)),
               ("상위 20% 선발", pre > np.quantile(pre, 0.8)),
               ("전체", np.ones(n, bool))]:
    t = stats.ttest_rel(post[m], pre[m])
    print(f"  {lab:>14s}: 사전 {pre[m].mean():6.2f} → 사후 {post[m].mean():6.2f}"
          f"  변화 {post[m].mean() - pre[m].mean():+6.2f}  p={t.pvalue:.2e}")

아무 처치도 하지 않았다 (사전·사후 상관 0.7, 평균 100, SD 10)
       하위 20% 선발: 사전  86.06 → 사후  90.16  변화  +4.10  p=0.00e+00
       상위 20% 선발: 사전 114.11 → 사후 109.90  변화  -4.22  p=0.00e+00
              전체: 사전 100.07 → 사후 100.02  변화  -0.06  p=1.93e-02

print("r 을 두 집단의 성공률 차이로 읽으면")
print(f"{'r':>6s} {'r²':>8s} {'처치군':>8s} {'대조군':>8s} {'차이':>7s}")
for r in [0.05, 0.1, 0.2, 0.3, 0.5]:
    print(f"{r:6.2f} {r**2:8.4f} {0.5 + r / 2:8.2f} {0.5 - r / 2:8.2f} {r:7.2f}")

r 을 두 집단의 성공률 차이로 읽으면
     r       r²      처치군      대조군      차이
  0.05   0.0025     0.53     0.47    0.05
  0.10   0.0100     0.55     0.45    0.10
  0.20   0.0400     0.60     0.40    0.20
  0.30   0.0900     0.65     0.35    0.30
  0.50   0.2500     0.75     0.25    0.50

print("\n실제 연구의 효과크기")
for lab, r in [("아스피린과 심근경색 예방", 0.034), ("흡연과 폐암", 0.29),
               ("수면제와 단기 수면 개선", 0.30), ("심리치료 효과", 0.32)]:
    print(f"  {lab:>22s}: r={r:.3f}, r²={r**2 * 100:5.2f}%, "
          f"BESD {0.5 - r / 2:.3f} → {0.5 + r / 2:.3f}")

실제 연구의 효과크기
           아스피린과 심근경색 예방: r=0.034, r²= 0.12%, BESD 0.483 → 0.517
                  흡연과 폐암: r=0.290, r²= 8.41%, BESD 0.355 → 0.645
           수면제와 단기 수면 개선: r=0.300, r²= 9.00%, BESD 0.350 → 0.650
                 심리치료 효과: r=0.320, r²=10.24%, BESD 0.340 → 0.660

import warnings
warnings.filterwarnings("ignore")

import numpy as np
from scipy import stats

rng = np.random.default_rng(24003)
n = 300
target = 0.60
print(f"목표 r = {target} 를 갖는 서로 다른 자료 다섯 가지")
print(f"{'구조':>22s} {'r':>8s} {'스피어만':>9s} {'최대 쿡 D':>10s} {'y 왜도':>8s}")

def rep(lab, x, y):
    b = np.polyfit(x, y, 1)
    e = y - np.polyval(b, x)
    h = 1 / len(x) + (x - x.mean())**2 / ((x - x.mean())**2).sum()
    s2 = (e**2).sum() / (len(x) - 2)
    cook = e**2 * h / (2 * s2 * (1 - h)**2)
    print(f"{lab:>22s} {np.corrcoef(x, y)[0, 1]:8.4f} "
          f"{stats.spearmanr(x, y).statistic:9.4f} {cook.max():10.4f} "
          f"{stats.skew(y):8.3f}")

z = rng.standard_normal((n, 2))                      # (1) 이변량 정규
x = z[:, 0]
y = target * z[:, 0] + np.sqrt(1 - target**2) * z[:, 1]
rep("이변량 정규", x, y)

x2 = rng.standard_normal(n)                          # (2) 분산이 커지는 구조
y2 = 0.78 * x2 + rng.standard_normal(n) * (0.4 + 0.5 * np.abs(x2))
y2 = (y2 - y2.mean()) / y2.std()
rep("이분산 (부채꼴)", x2, y2)

x3 = rng.standard_normal(n)                          # (3) 두 덩어리
g = rng.random(n) < 0.5
y3 = np.where(g, x3 + 1.2, x3 - 1.2) * 0.90 + rng.standard_normal(n) * 0.55
rep("두 군집", x3, y3)

x4 = rng.standard_normal(n - 1)                      # (4) 이상점 하나가 주도
y4 = rng.standard_normal(n - 1) * 1.0 + 0.44 * x4
x4 = np.append(x4, 12.0)
y4 = np.append(y4, 12.0)
rep("이상점 하나", x4, y4)

x5 = rng.uniform(-2, 2, n)                           # (5) 곡선 + 잡음
y5 = 0.47 * x5 + 0.30 * x5**2 + rng.standard_normal(n) * 0.7
rep("곡선 (2차항 포함)", x5, y5)

목표 r = 0.6 를 갖는 서로 다른 자료 다섯 가지
                    구조        r      스피어만     최대 쿡 D     y 왜도
                이변량 정규   0.6338    0.6291     0.0413   -0.083
             이분산 (부채꼴)   0.6082    0.6462     0.2856   -0.275
                  두 군집   0.6071    0.5871     0.0332   -0.106
                이상점 하나   0.6139    0.4235     6.2910    2.611
           곡선 (2차항 포함)   0.6144    0.6174     0.0306    0.134

1. 산점도               ─ 모양, 군집, 이상점
2. 잔차 대 적합값 그림   ─ 곡률, 이분산
3. 피어슨 vs 스피어만    ─ 이상점, 단조 비선형
4. 쿡 D / 지렛값         ─ 영향점

□ 산점도를 그렸는가
□ 표본 크기 n 을 밝혔는가
□ 신뢰구간을 함께 썼는가 (피셔 z)
□ 어느 상관계수인지 명시했는가 (피어슨/스피어만/켄들)
□ 결측 처리 방식을 밝혔는가
□ 이상점 유무와 그 영향을 점검했는가
□ 인과를 암시하는 표현을 피했는가

학습 시간과 시험 점수 (n = 120)

  피어슨 r = 0.48,  95% CI [0.33, 0.61],  p < 0.001
  스피어만 r_s = 0.46 (두 값이 가까워 이상점의 영향은 작다)

  산점도에서 선형 관계가 확인되었고 잔차의 등분산도 만족했다.
  쿡의 거리가 0.1 을 넘는 관측은 없었다.

  r² = 0.23 이지만, 관찰연구이므로 "학습 시간이 점수를
  높인다"고 해석할 수 없다. 능력·동기 등 공통 원인이 있을 수 있다.