상자그림¶

import matplotlib.pyplot as plt
import pandas as pd

# 그림에 한글을 쓰므로 한글 글꼴을 지정한다. 맥이면 'Apple SD Gothic Neo',
# 윈도우면 'Malgun Gothic', 리눅스면 'NanumGothic' 정도가 무난하다.
plt.rcParams['font.family'] = 'Apple SD Gothic Neo'
plt.rcParams['axes.unicode_minus'] = False

# 자료를 인터넷에서 내려받으므로 실행에 연결이 필요하다.
url = "https://raw.githubusercontent.com/datasciencedojo/datasets/f0ccab6a7ceafdff780052166fb6fab3311398eb/titanic.csv"
df = pd.read_csv(url, index_col='PassengerId')

fig, ax = plt.subplots(figsize=(5, 3))

# kind='box'로 pandas Series에서 곧바로 상자그림을 그린다.
# vert=False는 눕혀 그린다는 뜻. 눕히면 축 이름이 길어도 읽기 편하다.
# 결측값(Age가 비어 있는 177명)은 pandas가 알아서 제외한다.
df['Age'].plot(kind='box', ax=ax, vert=False)

ax.set_title("타이타닉 승객 나이의 가로 상자그림")
ax.set_xlabel("나이 (세)")
# 상자와 수염 자체에 눈이 가도록 불필요한 테두리를 지운다
ax.spines[["top", "left", "right"]].set_visible(False)
plt.show()

# 상자를 이루는 세 숫자(Q1, 중앙값, Q3)와 자료의 최솟값·최댓값을 확인한다.
# 주의: min 과 max 는 '수염의 끝'이 아니다. 수염은 울타리 안의 가장 극단적인
# 자료점까지만 뻗고, 그 바깥은 이상치 점으로 따로 찍힌다.
print(df['Age'].describe()[['min', '25%', '50%', '75%', 'max']].round(2))

min     0.42
25%    20.12
50%    28.00
75%    38.00
max    80.00
Name: Age, dtype: float64

import numpy as np
from scipy.stats import norm

z = norm.ppf(0.75)
fence = z + 1.5 * (2 * z)
p = 2 * norm.cdf(-fence)
print(f"z_0.75 = {z:.6f},  정규자료의 IQR = 2z = {2 * z:.6f} sigma")
print(f"울타리 = {z:.5f} + {1.5 * 2 * z:.5f} = {fence:.6f} sigma")
print(f"울타리 밖 확률 p = 2 Phi(-{fence:.5f}) = {p:.6f}  (1000 개 중 {p * 1000:.2f} 개)")


def flagged(x):
    q1, q3 = np.percentile(x, [25, 75])
    iqr = q3 - q1
    return (x < q1 - 1.5 * iqr) | (x > q3 + 1.5 * iqr)


rng = np.random.default_rng(3)
print(f"\n{'n':>8}{'반복':>7}{'모의 평균':>11}{'이론 p*n':>11}{'비':>7}")
for n in [100, 1000, 10000, 100000]:
    reps = max(20, 2_000_000 // n)
    m = np.mean([flagged(rng.standard_normal(n)).sum() for _ in range(reps)])
    print(f"{n:>8}{reps:>7}{m:>11.3f}{p * n:>11.3f}{m / (p * n):>7.3f}")

ns = np.arange(1000, 10001, 1000)
means = []
for n in ns:
    means.append(np.mean([flagged(rng.standard_normal(n)).sum() for _ in range(60)]))
means = np.array(means)
print(f"\nn 에 대한 회귀기울기 = {np.sum(ns * means) / np.sum(ns * ns):.6f}   "
      f"(이론 {p:.6f})")

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 3.6),
                               gridspec_kw={"width_ratios": [1, 1.2]})
x = rng.standard_normal(10000)
f = flagged(x)
ax1.boxplot(x, vert=False, widths=0.5, flierprops=dict(marker='o', ms=4,
            markerfacecolor="#D32F2F", markeredgecolor="#D32F2F"))
ax1.set_title(f"완전한 정규자료 10000개 — 그래도 {f.sum()}개가 찍힌다", fontsize=11)
ax1.set_xlabel("$z$")
ax1.set_yticks([])
ax1.spines[["top", "left", "right"]].set_visible(False)

ax2.plot(ns, means, "o", color="#1565C0", ms=6, label="모의 (60회 평균)")
ax2.plot(ns, p * ns, "-", color="#D32F2F", lw=1.8,
         label=f"이론 $0.006977\\,n$")
ax2.set_xlabel("표본 크기 $n$")
ax2.set_ylabel("이상치로 찍힌 개수")
ax2.set_title("정규자료에서도 개수는 $n$ 에 비례해 늘어난다", fontsize=11)
ax2.legend(fontsize=9)
ax2.spines[["top", "right"]].set_visible(False)
plt.tight_layout()
plt.show()

z_0.75 = 0.674490,  정규자료의 IQR = 2z = 1.348980 sigma
울타리 = 0.67449 + 2.02347 = 2.697959 sigma
울타리 밖 확률 p = 2 Phi(-2.69796) = 0.006977  (1000 개 중 6.98 개)

       n     반복      모의 평균     이론 p*n      비
     100  20000      1.002      0.698  1.436
    1000   2000      7.186      6.977  1.030
   10000    200     70.760     69.766  1.014
  100000     20    705.200    697.660  1.011

n 에 대한 회귀기울기 = 0.007055   (이론 0.006977)

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

# 그림에 한글이 들어가므로 한글 글꼴을 지정한다(지정하지 않으면 네모로 깨진다).
plt.rcParams['font.family'] = 'Apple SD Gothic Neo'
plt.rcParams['axes.unicode_minus'] = False

np.random.seed(0)

# 오른쪽으로 치우친 자료를 만드는 요령:
# 중심이 다른 정규분포 셋을 개수를 줄여 가며 겹친다.
#   0 부근에 1000개  (본체)
#   2 부근에  200개  (오른쪽 어깨)
#   4 부근에  100개  (오른쪽 꼬리)
# 왼쪽에는 대응하는 덩어리가 없으므로 분포가 오른쪽으로 길어진다.
main_data = stats.norm().rvs(1_000)
right_1 = stats.norm(loc=2).rvs(200)
right_2 = stats.norm(loc=4).rvs(100)
combined = np.concatenate((main_data, right_1, right_2))

# 같은 자료를 위아래로 나란히 놓아 두 그림을 대응시킨다
fig, (ax_hist, ax_box) = plt.subplots(2, 1, figsize=(12, 6))

# 위: 히스토그램. 분포의 모양이 그대로 보인다.
ax_hist.hist(combined, density=True, bins=30)
ax_hist.set_title('오른쪽으로 치우친 자료의 히스토그램')

# 아래: 같은 자료의 상자그림. 다섯 숫자로 압축된 모습이다.
ax_box.boxplot(combined, vert=False)
ax_box.set_title('같은 자료의 상자그림')

plt.tight_layout()
plt.show()

# 치우침이 숫자로도 드러나는지 확인한다.
# 오른쪽으로 치우치면 평균 > 중앙값 이고, Q3-Q2 가 Q2-Q1 보다 크다.
q1, q2, q3 = np.percentile(combined, [25, 50, 75])
print(f"평균 {combined.mean():.3f}  중앙값 {q2:.3f}")
print(f"Q2-Q1 = {q2-q1:.3f}   Q3-Q2 = {q3-q2:.3f}  (오른쪽이 길다)")

평균 0.595  중앙값 0.314
Q2-Q1 = 0.794   Q3-Q2 = 1.098  (오른쪽이 길다)

from scipy.optimize import brentq

w = np.array([1000., 200., 100.]); w /= w.sum()
mu = np.array([0., 2., 4.])
M = (w * mu).sum()
V = (w * (mu ** 2 + 1.0)).sum() - M ** 2
dlt = mu - M
m3 = (w * (dlt ** 3 + 3 * dlt)).sum()
print(f"혼합분포의 참값:  평균 {M:.6f},  분산 {V:.6f},  표준편차 {np.sqrt(V):.6f}")
print(f"  3차 중심적률 {m3:.6f},  적률 왜도 = m3/sd^3 = {m3 / V ** 1.5:.6f}")

cdf = lambda x: (w * stats.norm.cdf(x, loc=mu)).sum()
Q = [brentq(lambda x: cdf(x) - q, -6, 12) for q in [0.25, 0.5, 0.75]]
bowley = ((Q[2] - Q[1]) - (Q[1] - Q[0])) / (Q[2] - Q[0])
print(f"  참 사분위수 Q1={Q[0]:.4f}  Q2={Q[1]:.4f}  Q3={Q[2]:.4f}")
print(f"  참 Q2-Q1={Q[1]-Q[0]:.4f}  Q3-Q2={Q[2]-Q[1]:.4f}  보울리 왜도 {bowley:.4f}")
print(f"  참 피어슨 3(평균-중앙값)/sd = {3 * (M - Q[1]) / np.sqrt(V):.4f}")

q1, q2, q3 = np.percentile(combined, [25, 50, 75])
bw_s = ((q3 - q2) - (q2 - q1)) / (q3 - q1)
print(f"\n표본 (n = {len(combined)}):")
print(f"  평균 {combined.mean():.6f}  중앙값 {q2:.6f}  표준편차 {combined.std(ddof=1):.6f}")
print(f"  표본 Q1={q1:.4f}  Q2={q2:.4f}  Q3={q3:.4f}")
print(f"  적률 왜도 {stats.skew(combined):.6f}   "
      f"(몬테카를로 오차 ~ sqrt(6/n) = {np.sqrt(6/len(combined)):.4f})")
print(f"  보울리 왜도 {bw_s:.4f}   피어슨 "
      f"{3 * (combined.mean() - q2) / combined.std(ddof=1):.4f}")
print(f"\n상자가 덮는 가로 폭 = IQR / 자료 범위 = "
      f"{(q3-q1)/(combined.max()-combined.min()):.4f}")

fig, ax = plt.subplots(figsize=(11, 3.4))
ax.hist(combined, bins=40, density=True, color="#CFD8DC", zorder=0)
ax.axvspan(q1, q3, color="#DCEBFB", zorder=1, label=f"상자 (IQR = {q3-q1:.2f})")
ax.axvline(q2, color="#1565C0", lw=2, zorder=3, label=f"중앙값 {q2:.2f}")
ax.axvline(combined.mean(), color="#D32F2F", lw=2, ls="--", zorder=3,
           label=f"평균 {combined.mean():.2f}")
ax.axvline(q3 + 1.5 * (q3 - q1), color="#E65100", lw=1.6, ls=":", zorder=3,
           label=f"위 울타리 {q3 + 1.5*(q3-q1):.2f}")
ax.set_xlabel("값")
ax.set_ylabel("밀도")
ax.set_title(f"적률 왜도 {stats.skew(combined):.2f},  피어슨 "
             f"{3*(combined.mean()-q2)/combined.std(ddof=1):.2f},  보울리 {bw_s:.2f}"
             "  — 같은 자료, 다른 수", fontsize=11)
ax.legend(fontsize=9)
ax.spines[["top", "right"]].set_visible(False)
plt.tight_layout()
plt.show()

혼합분포의 참값:  평균 0.615385,  분산 2.467456,  표준편차 1.570814
  3차 중심적률 3.211652,  적률 왜도 = m3/sd^3 = 0.828618
  참 사분위수 Q1=-0.4577  Q2=0.3583  Q3=1.4020
  참 Q2-Q1=0.8159  Q3-Q2=1.0438  보울리 왜도 0.1225
  참 피어슨 3(평균-중앙값)/sd = 0.4911

표본 (n = 1300):
  평균 0.595182  중앙값 0.313942  표준편차 1.587044
  표본 Q1=-0.4800  Q2=0.3139  Q3=1.4115
  적률 왜도 0.847545   (몬테카를로 오차 ~ sqrt(6/n) = 0.0679)
  보울리 왜도 0.1605   피어슨 0.5316

상자가 덮는 가로 폭 = IQR / 자료 범위 = 0.2019

import numpy as np
import matplotlib.pyplot as plt

# 그림에 한글이 들어가므로 한글 글꼴을 지정한다(지정하지 않으면 네모로 깨진다).
plt.rcParams['font.family'] = 'Apple SD Gothic Neo'
plt.rcParams['axes.unicode_minus'] = False

# 표본 크기를 10^4, 5*10^4, 10^5 로 늘려 가며 얻은 몬테카를로 추정 오차를 흉내 낸 자료다.
# 표본이 커질수록 퍼짐이 줄어드는 모습을 만들기 위해 비슷한 값들의 묶음에
# 0.5, 0.25를 곱했다. 극단값은 일부러 서로 다르게 두어, 표본이 커져도
# 이상치는 남을 수 있다는 점을 보인다.
data_a = np.array([1, 2, 0, 0, 0, 1, 3, 1, 2, 1, 2, 4, 5, -1, -2, 0, 8])
data_b = np.array([1, 2, 0, 0, 0, 1, 3, 1, 2, 1, 2, 4, 5, -1, -2, 0, -8]) * 0.5
data_c = np.array([1, 2, 0, 0, 0, 1, 3, 1, 2, 1, 2, 4, 5, -1, -2, 0, 10, -7]) * 0.25

fig, ax = plt.subplots()

# 리스트를 넘기면 상자를 나란히 그린다. 이것이 상자그림의 가장 큰 쓸모다.
ax.boxplot([data_a, data_b, data_c])

# 각 상자 아래에 이름을 붙인다.
# boxplot에 직접 주는 인자는 matplotlib 버전에 따라 이름이 다르다
# (3.9 미만은 labels=, 3.9 이상은 tick_labels=). 아래처럼 축에 직접 주면
# 버전에 상관없이 동작한다.
ax.set_xticklabels(["$10^4$", "$5 \\cdot 10^4$", "$10^5$"])

# 비교 기준선. 이론값 1을 가로선으로 깔아 두면
# 각 상자가 그 선을 얼마나 감싸는지 눈으로 볼 수 있다.
ax.plot([0, 1, 2, 3, 4], [1, 1, 1, 1, 1],
        label="이론값 = 1", linestyle="--", color="r", alpha=0.7)

ax.legend()
ax.set_ylim(-10.0, 10.0)      # 세 상자를 같은 눈금에 두어야 비교가 성립한다
ax.set_xlabel('표본 크기')
ax.set_ylabel('몬테카를로 추정값')
plt.show()

# 퍼짐이 실제로 줄어드는지 IQR로 확인한다
for name, d in [("10^4", data_a), ("5*10^4", data_b), ("10^5", data_c)]:
    q1, q3 = np.percentile(d, [25, 75])
    print(f"{name:>7}: 중앙값 {np.median(d):6.2f}  IQR {q3-q1:.2f}")

   10^4: 중앙값   1.00  IQR 2.00
 5*10^4: 중앙값   0.50  IQR 1.00
   10^5: 중앙값   0.25  IQR 0.50

from scipy.stats import norm

z = norm.ppf(0.75)
sizes = [10_000, 50_000, 100_000]
names = ["$10^4$", "$5\\cdot 10^4$", "$10^5$"]
sigma = 1.0                      # Exp(1) 의 표준편차
print(f"이론: IQR = 2 z_0.75 sigma / sqrt(n) = {2*z:.6f} / sqrt(n)")
print(f"{'n':>8}{'이론 IQR':>12}{'모의 IQR':>12}{'비':>8}{'첫 IQR 대비(이론)':>18}"
      f"{'(모의)':>10}")
rng = np.random.default_rng(11)
reps = 400
est = []
for n in sizes:
    # 한 번에 큰 배열을 잡지 않도록 반복마다 따로 뽑는다. 참값은 E[X] = 1.
    est.append(np.array([rng.exponential(1.0, n).mean() for _ in range(reps)]))
first_t = 2 * z * sigma / np.sqrt(sizes[0])
first_s = np.subtract(*np.percentile(est[0], [75, 25]))
for n, e in zip(sizes, est):
    t = 2 * z * sigma / np.sqrt(n)
    s = np.subtract(*np.percentile(e, [75, 25]))
    print(f"{n:>8}{t:>12.6f}{s:>12.6f}{s/t:>8.3f}{t/first_t:>18.4f}{s/first_s:>10.4f}")

big = np.array([rng.exponential(1.0, sizes[0]).mean() for _ in range(10_000)])
print(f"\nn = 10^4 에서 반복을 10000 회로 늘리면 비 = "
      f"{np.subtract(*np.percentile(big, [75, 25])) / first_t:.4f}")
print(f"이론 비:  sqrt(1/5) = {np.sqrt(1/5):.4f},  sqrt(1/10) = {np.sqrt(1/10):.4f}")
print(f"손으로 만든 자료의 비: 1.00 -> 0.50 -> 0.25")
print(f"  0.50 에 맞는 n 배수 = 1/0.50^2 = {1/0.5**2:.0f},  "
      f"0.25 에 맞는 n 배수 = 1/0.25^2 = {1/0.25**2:.0f}")

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 3.6))
ax1.boxplot(est, widths=0.5)
ax1.axhline(1.0, color="#D32F2F", ls="--", lw=1.4, label="참값 = 1")
ax1.set_xticklabels(names)
ax1.set_xlabel("표본 크기 $n$")
ax1.set_ylabel("몬테카를로 추정값")
ax1.set_title("진짜 몬테카를로 추정값 (400 회씩)", fontsize=11)
ax1.legend(fontsize=9)

iqr_t = [2 * z * sigma / np.sqrt(n) for n in sizes]
iqr_s = [np.subtract(*np.percentile(e, [75, 25])) for e in est]
ax2.plot(sizes, np.array(iqr_t) / iqr_t[0], "-", color="#D32F2F", lw=1.8,
         label="이론 $\\sqrt{n_1/n}$")
ax2.plot(sizes, np.array(iqr_s) / iqr_s[0], "o", color="#1565C0", ms=7,
         label="모의")
ax2.plot(sizes, [1.0, 0.5, 0.25], "s", color="#E65100", ms=7,
         label="쪽의 예시 자료")
ax2.set_xticks(sizes)
ax2.set_xticklabels(names)
ax2.set_xlabel("표본 크기 $n$")
ax2.set_ylabel("IQR (첫 값 대비)")
ax2.set_title("예시 자료는 이론보다 빨리 줄어든다", fontsize=11)
ax2.legend(fontsize=9)
for ax in (ax1, ax2):
    ax.spines[["top", "right"]].set_visible(False)
plt.tight_layout()
plt.show()

이론: IQR = 2 z_0.75 sigma / sqrt(n) = 1.348980 / sqrt(n)
       n      이론 IQR      모의 IQR       비      첫 IQR 대비(이론)      (모의)
   10000    0.013490    0.013231   0.981            1.0000    1.0000
   50000    0.006033    0.005954   0.987            0.4472    0.4500
  100000    0.004266    0.004197   0.984            0.3162    0.3172

n = 10^4 에서 반복을 10000 회로 늘리면 비 = 1.0097
이론 비:  sqrt(1/5) = 0.4472,  sqrt(1/10) = 0.3162
손으로 만든 자료의 비: 1.00 -> 0.50 -> 0.25
  0.50 에 맞는 n 배수 = 1/0.50^2 = 4,  0.25 에 맞는 n 배수 = 1/0.25^2 = 16

import numpy as np
import matplotlib.pyplot as plt

# 그림에 한글이 들어가므로 한글 글꼴을 지정한다(지정하지 않으면 네모로 깨진다).
plt.rcParams['font.family'] = 'Apple SD Gothic Neo'
plt.rcParams['axes.unicode_minus'] = False

rng = np.random.default_rng(0)
n = 2000
unimodal = rng.normal(0, 1, n)
bimodal = np.concatenate([rng.normal(-1.6, 0.35, n // 2),
                          rng.normal(1.6, 0.35, n // 2)])

# 중앙값 0, IQR 1 로 맞춘다
norm = lambda x: (x - np.median(x)) / np.subtract(*np.percentile(x, [75, 25]))
unimodal, bimodal = norm(unimodal), norm(bimodal)

print("다섯 수치 요약")
for label, x in [("단봉", unimodal), ("이봉", bimodal)]:
    print(f"  {label}: {np.round(np.percentile(x, [0, 25, 50, 75, 100]), 3)}")

fig, axes = plt.subplots(1, 3, figsize=(13, 4))
axes[0].boxplot([unimodal, bimodal])
axes[0].set_xticklabels(["단봉", "이봉"])
axes[0].set_title("상자그림 — 상자가 거의 같다", fontsize=10)

axes[1].violinplot([unimodal, bimodal], showmedians=True)
axes[1].set_xticks([1, 2]); axes[1].set_xticklabels(["단봉", "이봉"])
axes[1].set_title("바이올린 — 구조가 드러난다", fontsize=10)

# 점이 2000개면 너무 빽빽하므로 5개마다 하나씩 400개만 뽑아 찍는다.
# 앞에서 400개를 잘라 쓰면 이봉 자료의 '앞쪽 봉우리'만 뽑히므로 안 된다.
for i, (label, x) in enumerate([("단봉", unimodal), ("이봉", bimodal)]):
    sub = x[::5]
    axes[2].scatter(np.full(len(sub), i + 1) + rng.normal(0, 0.06, len(sub)), sub,
                    s=4, alpha=0.3)
axes[2].set_xticks([1, 2]); axes[2].set_xticklabels(["단봉", "이봉"])
axes[2].set_title("점 흩뿌리기 — 원자료를 그대로", fontsize=10)
fig.tight_layout()
plt.show()

다섯 수치 요약
  단봉: [-2.843 -0.491 -0.     0.509  2.279]
  이봉: [-0.823 -0.512  0.     0.488  0.834]

import numpy as np

rng = np.random.default_rng(0)
n = 200_000

print(f"{'분포':>10}{'표시 비율':>12}{'위쪽':>10}{'아래쪽':>10}")
for label, x in [("정규", rng.normal(0, 1, n)),
                 ("지수", rng.exponential(1, n)),
                 ("로그정규", rng.lognormal(0, 1, n))]:
    q1, q3 = np.percentile(x, [25, 75])
    iqr = q3 - q1
    lo, hi = q1 - 1.5 * iqr, q3 + 1.5 * iqr
    print(f"{label:>10}{np.mean((x < lo) | (x > hi)):>12.4f}"
          f"{np.mean(x > hi):>10.4f}{np.mean(x < lo):>10.4f}")

        분포       표시 비율        위쪽       아래쪽
        정규      0.0071    0.0035    0.0036
        지수      0.0481    0.0481    0.0000
      로그정규      0.0774    0.0774    0.0000

import numpy as np

print(f"{'n':>5s}{'노치 반폭/상자 반높이':>22s}{'모양':>10s}")
for n in [4, 6, 8, 9, 10, 12, 20, 50]:
    ratio = 1.57 / np.sqrt(n) / 0.5
    print(f"{n:>5d}{ratio:>22.4f}{'모래시계' if ratio > 1 else '정상':>10s}")

print(f"\n경계: n = {3.14 ** 2:.4f}  ->  n <= 9 이면 모래시계")

    n          노치 반폭/상자 반높이        모양
    4                1.5700      모래시계
    6                1.2819      모래시계
    8                1.1102      모래시계
    9                1.0467      모래시계
   10                0.9930        정상
   12                0.9064        정상
   20                0.7021        정상
   50                0.4441        정상

경계: n = 9.8596  ->  n <= 9 이면 모래시계

import numpy as np
import matplotlib.pyplot as plt

# 그림에 한글이 들어가므로 한글 글꼴을 지정한다(지정하지 않으면 네모로 깨진다).
plt.rcParams['font.family'] = 'Apple SD Gothic Neo'
plt.rcParams['axes.unicode_minus'] = False

rng = np.random.default_rng(7)
fig, ax = plt.subplots(1, 2, figsize=(10, 4.4))

for k, n in enumerate([6, 60]):
    data = [rng.normal(0, 1, n) for _ in range(3)]
    bp = ax[k].boxplot(data, notch=True, patch_artist=True, widths=0.5)
    for b in bp["boxes"]:
        b.set_facecolor("lightsteelblue"); b.set_edgecolor("black")
    for med in bp["medians"]:
        med.set_color("crimson"); med.set_linewidth(2)
    ratio = 1.57 / np.sqrt(n) / 0.5
    ax[k].set_title(f"n = {n}   (노치 반폭 / 상자 반높이 = {ratio:.2f})")
    ax[k].set_xlabel("집단"); ax[k].grid(alpha=0.25, axis="y")
ax[0].set_ylabel("값")
fig.suptitle("표본이 작으면 노치가 상자를 넘어 접힌다 — 모래시계 모양", y=1.00)
fig.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# 그림에 한글이 들어가므로 한글 글꼴을 지정한다(지정하지 않으면 네모로 깨진다).
plt.rcParams['font.family'] = 'Apple SD Gothic Neo'
plt.rcParams['axes.unicode_minus'] = False

rng = np.random.default_rng(5)
groups = {"A (n=8)": rng.normal(10, 2, 8),
          "B (n=40)": rng.normal(10, 2, 40),
          "C (n=500)": rng.normal(10, 2, 500)}

print("세 집단은 모두 같은 모집단 N(10, 4) 에서 나왔다")
for name, x in groups.items():
    q1, q3 = np.percentile(x, [25, 75])
    print(f"  {name:<10} 중앙값 {np.median(x):>6.3f}   IQR {q3 - q1:>6.3f}   "
          f"중앙값의 SE {1.253 * x.std(ddof=1) / np.sqrt(len(x)):>6.3f}")

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
data = list(groups.values())
axes[0].boxplot(data)
axes[0].set_xticklabels(list(groups))
axes[0].set_title("보통 상자그림 — 표본 크기를 알 수 없다", fontsize=10)

axes[1].boxplot(data, widths=[0.25 * np.sqrt(len(d) / 500) + 0.12 for d in data],
                notch=True, bootstrap=2000)
for i, d in enumerate(data):
    axes[1].scatter(np.full(len(d), i + 1) + rng.normal(0, 0.04, len(d)), d,
                    s=6, alpha=0.25, color="black", zorder=3)
axes[1].set_xticklabels(list(groups))
axes[1].set_title("가변 너비 + 노치 + 점", fontsize=10)
fig.tight_layout()
plt.show()

세 집단은 모두 같은 모집단 N(10, 4) 에서 나왔다
  A (n=8)    중앙값  9.199   IQR  1.953   중앙값의 SE  0.697
  B (n=40)   중앙값  9.439   IQR  2.661   중앙값의 SE  0.364
  C (n=500)  중앙값  9.978   IQR  2.719   중앙값의 SE  0.109