안녕하세요!

오늘은 “다층 연구 설계의 최적화(Sample Size and Power Analysis in Multilevel Designs)”에 대해 살펴보겠습니다. “학교 현장의 데이터”를 예시로 들어 직관적인 설명과 수리적 엄밀함을 모두 갖춘 형태로 재구성해 드리겠습니다.

분석 도구로는 jamovi의 기반이 되는 R 코드를 함께 제시하여 ‘최적 표본 크기 산출’부터 ‘모의 데이터 생성’, 그리고 ‘분석’까지 완벽하게 구현해 드리겠습니다. (참고: jamovi는 데이터 분석에는 강력하지만, 연구 설계 단계의 복잡한 ‘최적 표본 산출’ 기능은 제한적이므로, 이 부분은 R의 수리적 계산 기능을 활용하고 분석은 jamovi의 혼합모형 논리로 설명하겠습니다.)

1. 서론: 연구자의 영원한 딜레마, “돈이냐, 정확성이냐?”

교육 연구를 진행할 때 우리는 항상 두 가지 제약 조건 사이에서 고민합니다.

1. 통계적 검증력(Power): 효과가 있다면 있다고 말할 수 있는 힘 (높을수록 좋음).
2. 예산(Budget): 연구비와 시간은 한정되어 있음.

이 장에서는 다층 모형(Multilevel Modeling) 상황, 즉 학생이 학교에 소속되어 있는 구조에서 어떻게 하면 가장 적은 비용으로 가장 높은 정확성을 얻을 수 있는지, 그 최적 설계(Optimal Design) 방법을 알려드리겠습니다.

우리의 가상 시나리오를 소개합니다.

[시나리오: 프로젝트 ‘독서왕’]

교육심리학자 김 교수는 새로운 독서 프로그램이 초등학생의 ‘문해력’을 높이는지 검증하려 합니다.

• 총 예산: 1,000만 원 (가상의 화폐 단위)
• 비용 구조:
  • 학교 하나를 섭외하는 비용(행정 절차, 학교 보상 등): 200만 원 (c)
  • 학생 한 명을 검사하는 비용(검사지, 간식 등): 10만 원 (s)
• 연구 질문: 몇 개의 학교를 섭외하고, 학교당 몇 명의 학생을 뽑아야 내 돈 1,000만 원 안에서 가장 정확한 결과를 얻을까요?

2. 군집 무작위 배정(Cluster Randomized Trial)의 설계

2.1 왜 다층 설계인가?

가장 쉬운 방법은 전국의 학생 명부에서 무작위로 학생을 뽑는 것입니다. 하지만 현실적으로 불가능합니다.

• 실행 가능성: 학교 단위로 프로그램을 돌려야 합니다.
• 오염(Contamination): 한 반에서 철수는 실험집단, 영희는 통제집단이면 서로 이야기하며 효과가 섞여버립니다.

그래서 우리는 학교(Cluster)를 통째로 실험군 혹은 대조군으로 배정하는 군집 무작위 배정(Cluster Randomized Trial)을 사용합니다.

2.2 급내상관계수(ICC)와 설계 효과

문제는 같은 학교 아이들끼리는 서로 비슷하다는 점입니다(학교 분위기, 선생님의 영향 등). 이를 급내상관계수(ICC, $\rho$)라고 합니다.

• ICC가 높다 = 학교 간 차이가 크다 = 같은 학교 아이들은 매우 비슷하다.
• ICC가 높으면, 학생을 100명 더 뽑는 것보다 학교를 1개 더 섭외하는 게 훨씬 중요해집니다.

2.3 최적 표본 크기 공식 (The Magic Formula)

주어진 예산($B$) 하에서 처치 효과($\beta_1$)의 분산($Var(\hat{\beta}_1)$)을 최소화하는 최적의 학교 수($K$)와 학교당 학생 수($n$)는 다음과 같습니다.

$n_{optimal} = \sqrt{\frac{c(1-\rho)}{s\rho}}$

$K_{optimal} = \frac{B}{c + s \cdot n}$

• $c$: 학교당 비용 (200)
• $s$: 학생당 비용 (10)
• $\rho$: ICC (가정된 값, 보통 0.05~0.10)

이 공식을 보면, ICC($\rho$)가 커질수록 학교당 학생 수($n$)는 줄여야 합니다. 왜냐하면 같은 학교에서 많이 뽑아봤자 정보가 중복되기 때문입니다.

3. R과 jamovi를 활용한 최적 설계 및 데이터 생성

이제 김 교수의 ‘독서왕’ 프로젝트를 위해 R을 사용하여 최적 표본을 계산하고, 이를 분석할 수 있는 모의 데이터를 생성해 보겠습니다.

3.1 최적 표본 크기 계산 (R Code)

김 교수의 상황: 예산 10,000, $c=200, s=10$, 그리고 ICC($\rho$)는 선행연구를 통해 0.05로 가정합니다.

R

# [R Code] 최적 표본 크기 산출
# 파라미터 설정
B <- 10000   # 총 예산
c <- 200     # 학교(Cluster)당 비용
s <- 10      # 학생(Person)당 비용
rho <- 0.05  # 급내상관계수 (ICC)

# 1. 최적의 학생 수 (n) 계산 [cite: 120]
n_opt <- sqrt((c * (1 - rho)) / (s * rho))

# 2. 최적의 학교 수 (K) 계산 [cite: 119]
K_opt <- B / (c + s * n_opt)

# 결과 출력
cat("최적의 학교당 학생 수 (n):", round(n_opt, 2), "명\n")
cat("최적의 학교 수 (K):", round(K_opt, 2), "개교\n")

[분석 결과 해석]

• 계산 결과, $n \approx 19.49$, $K \approx 25.32$가 나옵니다.
• 현실적으로 반올림하여 학교당 19명, 총 26개 학교(실험 13, 통제 13)를 섭외하는 것이 최적입니다.

3.2 비용 효율성 시각화 (Design Efficiency Plot)

아래 그래프는 학교당 학생 수($n$)를 변화시킬 때, 추정의 오차(분산)가 어떻게 변하는지 보여줍니다. 우리는 분산이 가장 낮은 지점을 찾아야 합니다.

R

# [R Code] 시각화 생성
library(ggplot2)

n_seq <- seq(5, 50, by = 1) # 학생 수를 5명에서 50명까지 변화시킴
K_seq <- B / (c + s * n_seq) # 예산 제약에 따른 학교 수

# 분산 계산 공식 (단순화된 형태) [cite: 128]
# g_rho 부분과 예산 부분을 결합
design_var <- function(n, K, rho) {
  # Standard Error calculation based on Eq 11.3 & 11.6 logic relation
  # Here we look at relative variance proportional to the function
  design_effect <- 1 + (n - 1) * rho
  total_N <- n * K
  return(design_effect / total_N) 
}

var_values <- design_var(n_seq, K_seq, rho)

df_plot <- data.frame(n = n_seq, Variance = var_values)

ggplot(df_plot, aes(x = n, y = Variance)) +
  geom_line(color = "blue", linewidth = 1.2) +
  geom_vline(xintercept = 19, linetype = "dashed", color = "red") +
  annotate("text", x = 19, y = max(var_values), label = "Optimal n = 19", vjust = -1) +
  labs(title = "최적 표본 설계: 학교당 학생 수 vs. 추정 오차",
       subtitle = "예산 고정 시, 학생을 19명 뽑을 때 오차가 최소화됨",
       x = "학교당 학생 수 (n)", y = "치료 효과의 분산 (낮을수록 좋음)") +
  theme_minimal()

이 그래프를 통해 김 교수는 무작정 학생을 많이 뽑는다고 좋은 게 아니라, 학교 수와 학생 수의 황금 비율을 맞춰야 함을 알 수 있습니다.

4. 모의 데이터 생성 및 분석 (Linear Mixed Model)

이제 최적 설계($K=26, n=19$)에 따라 데이터를 수집했다고 가정하고, 이를 jamovi(또는 R의 lmer)에서 분석하는 방법을 보여드리겠습니다.

4.1 데이터 생성 (Backstory 포함)

• 상황: 26개 학교를 선정, 절반은 ‘독서왕 프로그램(Treatment)’, 절반은 ‘기존 수업(Control)’ 진행.
• 효과: 프로그램은 문해력 점수를 평균 5점 정도 높여줄 것으로 기대.
• 변산: 학교 간 차이(ICC) 반영.

R

# [R Code] 모의 데이터 생성
set.seed(123) # 재현성을 위해 시드 설정

K <- 26       # 학교 수
n <- 19       # 학교당 학생 수
N <- K * n    # 총 학생 수

# 학교 ID 및 치료 집단 배정 (0: 통제, 1: 처치)
school_id <- rep(1:K, each = n)
treatment <- rep(c(rep(0, K/2), rep(1, K/2)), each = n)

# 랜덤 효과 생성 (학교 간 차이)
u0j <- rep(rnorm(K, mean = 0, sd = sqrt(5)), each = n) # 학교 분산 = 5

# 오차항 생성 (학생 간 차이)
# ICC = 0.05 이려면, 학교분산/(학교분산+오차분산) = 0.05
# 5 / (5 + 95) = 0.05 -> 오차 분산은 95, SD는 약 9.75
eij <- rnorm(N, mean = 0, sd = sqrt(95))

# 고정 효과 (진짜 치료 효과 = 5점)
beta0 <- 50  # 평균 점수
beta1 <- 5   # 치료 효과

# 종속 변수 (문해력 점수) 생성 
# y_ij = beta0 + beta1*x_j + u_0j + e_ij
y <- beta0 + beta1 * treatment + u0j + eij

# 데이터 프레임 생성
data_sim <- data.frame(
  SchoolID = factor(school_id),
  StudentID = 1:N,
  Treatment = factor(treatment, labels = c("Control", "Program")),
  Score = y
)

# 데이터 확인
head(data_sim)

# CSV로 저장 (jamovi에서 불러오기 위함)
write.csv(data_sim, "chap11.csv", row.names = FALSE)

4.2 jamovi 식 분석 절차 (R lme4 문법 활용)

jamovi의 Linear Mixed Models (GAMLj 모듈) 메뉴에서 분석하는 것과 동일한 R 코드입니다.

데이터 구조는 다음과 같습니다.

• Level 1 (Person): 학생별 문해력 점수 (Score)
• Level 2 (Cluster): 학교 (SchoolID)
• Predictor: 처치 여부 (Treatment, Level 2 변수)

수식 모형은 다음과 같습니다15:

$y_{ij} = \beta_{0} + \beta_{1}Treatment_{j} + u_{0j} + \epsilon_{ij}$

R

# [R Code] 다층 모형 분석
library(lme4)
library(lmerTest)

# 모형 적합: 절편에 대해서만 학교별 무선 효과(Random Intercept) 허용
model <- lmer(Score ~ Treatment + (1 | SchoolID), data = data_sim)

# 결과 요약
summary(model)

[결과 해석 방법]

1. Fixed Effects: TreatmentProgram의 Estimate가 5에 가까운지, p-value가 0.05보다 작은지 확인합니다. (우리가 5로 설정했으므로 유의하게 나올 것입니다.)
2. Random Effects: SchoolID의 Variance가 설계 시 가정한 분산과 유사한지 확인합니다.

5. 심화: 더 복잡한 상황들

5.1 다기관 임상시험 (Multisite Trials)

만약 김 교수가 학교 전체를 배정하는 게 아니라, 각 학교 안에서 철수는 실험반, 영희는 통제반으로 나눌 수 있다면 어떨까요?

이를 Multisite Trial이라고 합니다.

• 장점: 학교 효과($u_{0j}$)가 치료 효과 추정에서 사라지므로 훨씬 강력한 검증력(Power)을 가집니다.
• 설계: 이 경우 $Var(\hat{\beta}_1)$은 학교 분산에 영향을 받지 않으므로, 더 적은 예산으로도 유의한 결과를 얻을 수 있습니다. 하지만 ‘오염’ 문제가 없어야만 가능합니다.

5.2 반복 측정 (Longitudinal Design)

만약 김 교수가 프로그램을 1년 동안 진행하면서 학생들의 변화를 보고 싶다면 몇 번 측정해야 할까요?

• 선형 변화(직선 성장)를 가정할 때: 시작(Baseline)과 끝(End), 딱 2번 측정하거나, 중간 지점 하나를 추가하는 것이 비용 대비 가장 효율적입니다.
• 이차 함수(곡선 성장)를 가정할 때: 최소 3번의 측정이 필요하며, 등간격으로 측정하는 것이 좋습니다.

5.3 현실적인 문제와 해결책

1. ICC를 모를 때: 보통 문헌 연구를 통해 보수적으로(약간 높게) 잡습니다. ICC를 실제보다 절반 정도로 낮게 잘못 예측했더라도, 최적 설계 대비 효율성 손실은 약 10% 내외로 크지 않다는 연구가 있습니다.
2. 학교 크기가 다를 때: 모든 학교가 19명일 수는 없습니다. 학교 간 크기 편차(CV)가 0.5 정도라면, 계산된 학교 수($K$)보다 약 11% 정도 더 많은 학교를 섭외하여 이를 보정해야 합니다.

6. 결론 및 제언

오늘 살펴본 내용을 요약하면 다음과 같습니다.

1. 학교 현장 연구에서는 학생 수만큼이나 학교 수(Cluster Number)가 중요하다.
2. 비용 함수와 ICC를 고려하면, 무조건 많은 표본보다 최적의 비율을 찾는 것이 경제적이다.
3. 학교를 통째로 배정(CRT)하는 것보다 학교 내 무선 배정(Multisite)이 통계적으로는 더 유리하나, 오염 가능성을 고려해야 한다.

참고문헌 (APA Style)

• Maas, C. J. M., & Hox, J. J. (2005). Sufficient sample sizes for multilevel modeling. Methodology, 1(3), 86-92.
• Moerbeek, M., & Teerenstra, S. (2011). Optimal design in multilevel experiments. In J. J. Hox & J. K. Roberts (Eds.), Handbook of Advanced Multilevel Analysis (pp. 257-281). New York: Routledge.
• Raudenbush, S. W. (1997). Statistical analysis and optimal design for cluster randomized trials. Psychological Methods, 2(2), 173-185.
• Snijders, T. A. B., & Bosker, R. J. (1993). Standard errors and sample sizes for two-level research. Journal of Educational Statistics, 18(3), 237-259.
• Van Breukelen, G. J. P., & Moerbeek, M. (2013). Design considerations in multilevel studies. In M. A. Scott, J. S. Simonoff, & B. D. Marx (Eds.), The SAGE Handbook of Multilevel Modeling (pp. 183-200). SAGE Publications.

안녕하세요!

오늘은 개체 내 오차 구조(Within-Individual Error Structures)의 복잡성과 모델링에 대해 살펴보겠습니다. “학교 현장의 데이터”를 예시로 들어 직관적인 설명과 수리적 엄밀함을 모두 갖춘 형태로 재구성해 드리겠습니다.

분석 도구로는 jamovi의 사용법을 설명하되, jamovi의 기반이 되는 R 코드를 함께 제시하여 모의 데이터 생성부터 분석, 시각화까지 완벽하게 구현해 드리겠습니다.

우리가 학교에서 학생들의 성장을 추적할 때(종단 연구), 단순히 “평균 점수가 올랐나?”만 보는 것은 반쪽짜리 분석입니다. 학생 개개인이 어떻게 변화하는지, 그 변화의 폭(분산)은 일정한지, 어제의 성적이 오늘의 성적에 얼마나 영향을 미치는지(상관)를 파악해야 합니다. 제공해주신 텍스트는 이러한 공분산 구조(Covariance Structure)를 어떻게 모델링할 것인가에 대한 깊이 있는 통찰을 제공합니다.

1. 교육 현장의 예시: “읽기 유창성 성장 프로젝트”

이해를 돕기 위해 가상의 시나리오를 설정하겠습니다.

상황: A 초등학교에서 3학년 학생 50명을 대상으로 ‘읽기 유창성(1분당 읽은 단어 수)’을 1년 동안 4회(3월, 6월, 9월, 12월) 측정했습니다.

핵심 질문:

1. 모든 학생의 읽기 실력 격차(분산)는 3월이나 12월이나 똑같을까요? (등분산성)
2. 3월 성적이 좋은 학생은 6월에도 좋을까요? 12월까지 그 영향이 갈까요? (계열 상관)

2. 기본 모델과 개념 (The General Model)

우리는 데이터를 다음의 선형 회귀 모델로 표현할 수 있습니다.

$y = X\beta + Zu + \epsilon$

• $y$: 학생들의 읽기 점수 벡터
• $X\beta$ (고정 효과): 전체 학생들의 평균적인 성장 곡선 (예: 시간이 지날수록 점수가 오른다).
• $Zu$ (임의 효과): 학생 개인별 특성 (예: 어떤 학생은 시작부터 잘하고, 어떤 학생은 성장 속도가 빠르다).
• $\epsilon$ (오차 항): 설명되지 않는 나머지 변동. 오늘의 핵심 주제는 바로 이 $\epsilon$의 구조인 $\Sigma$를 파헤치는 것입니다.

3. 데이터 시각화와 진단 (Diagnostics)

모델링을 하기 전에 눈으로 확인해야 합니다. 텍스트에서는 프로파일 도표(Profile Plot), OSM, PRISM을 추천합니다.

3.1 프로파일 도표 (Profile Plot)

학생 개개인의 성장 궤적을 그린 그래프입니다.

• 해석: 선들이 서로 꼬이지 않고 나란히 간다면? $\rightarrow$ 학생 내 상관이 높음(잘하던 애가 계속 잘함).
• 분산: 시간이 갈수록 선들의 폭이 넓어진다면? $\rightarrow$ 이분산성(Heterogeneity) 존재.

3.2 OSM (Ordinary Scatterplot Matrix)

모든 시점 간의 산점도 행렬입니다.

• 특징: 대각선에서 멀어질수록(시간 격차가 클수록) 상관관계가 어떻게 변하는지 보여줍니다. 읽기 데이터에서는 보통 시간이 지날수록(대각선에서 멀어질수록) 상관이 낮아지는 패턴을 보입니다.

3.3 PRISM (Partial Regression-on-Intervenors Scatterplot Matrix)

이것은 저자(Zimmerman)가 강조하는 고급 도구입니다.

• 개념: $t_1$과 $t_3$의 관계를 볼 때, 중간 시점인 $t_2$의 영향을 제거하고(통제하고) 봅니다.
• 해석: 만약 $t_2$를 통제했을 때 $t_1$과 $t_3$의 관계가 사라진다면(무작위 산포), 이는 “조건부 독립”을 의미하며, Antedependence 모델이나 AR(1) 모델이 적합하다는 강력한 신호입니다.

4. 공분산 구조 모델링: 유연성 vs 간명성

어떤 모델을 선택해야 할까요? 너무 단순하면 현실을 반영 못 하고(과소적합), 너무 복잡하면 추정이 어렵습니다(과대적합). 이 둘 사이의 균형(Balance)이 핵심입니다.

4.1 가장 단순한 모델들 (Parsimonious)

1. 복합 대칭 (Compound Symmetry, CS):
   • 가정: 3월의 분산 = 12월의 분산. 3월-6월 상관 = 3월-12월 상관. 모든 것이 일정함.
   • 교육적 현실: 글쎄요, 보통 고학년이 될수록 격차(분산)가 커지고, 먼 시점 간의 상관은 떨어지기 때문에 현실적이지 않을 수 있습니다.
   • 구조: $\sigma^2 \begin{pmatrix} 1 & \rho & \rho & \rho \\ \rho & 1 & \rho & \rho \\ \rho & \rho & 1 & \rho \\ \rho & \rho & \rho & 1 \end{pmatrix}$
2. 1차 자기회귀 (AR(1)):
   • 가정: 바로 전 시험 성적이 다음 시험에 가장 큰 영향을 줌. 시간이 멀어질수록 상관이 지수적으로 감소 ($\rho, \rho^2, \rho^3…$).
   • 교육적 현실: 성장 데이터에 매우 적합합니다. 단, 분산이 일정하다고 가정하는(Homogeneous) 한계가 있습니다.

4.2 중간 단계의 모델들 (Heterogeneous & Banded)

3. 이분산 자기회귀 (Heterogeneous AR, ARH(1)):
   • AR(1)과 같지만, 시점별로 분산(시험 점수의 퍼짐 정도)이 다를 수 있게 허용합니다. 학년이 올라갈수록 실력 격차가 벌어지는 현상을 설명하기 좋습니다.
4. 토플리츠 (Toeplitz, TOEP):
   • AR(1)처럼 시간이 지날수록 상관이 변하지만, 그 패턴을 $\rho^k$처럼 강제로 수학적 공식에 넣지 않고 자유롭게 둡니다(Band 구조). 더 유연하지만 파라미터가 많아집니다.

4.3 가장 유연한 모델들 (Flexible)

5. 무구조 (Unstructured, UN):
   • 가정: 아무런 제약이 없습니다. 모든 시점의 분산과 상관을 각각 다 추정합니다.
   • 장단점: 가장 정확하지만, 측정 시점이 많아지면 추정해야 할 파라미터가 기하급수적으로 늘어나($n(n+1)/2$) 모델이 터질 수 있습니다.
6. Antedependence (AD) 및 Structured AD (SAD) 모델:
   • 텍스트의 핵심 제안: 1차 자기회귀(AR1)는 너무 딱딱하고, 무구조(UN)는 너무 무겁습니다.
   • AD 모델: 시점 간의 “조건부 독립”을 가정하여 파라미터를 줄이면서도, 비정상성(Non-stationarity: 분산과 상관이 변하는 성질)을 허용하는 똑똑한 모델입니다.
   • SAD 모델은 이를 더 구조화하여 파라미터 수를 $O(1)$ 수준으로 줄여 간명성을 확보합니다.

5. 분석 실습: R과 jamovi 활용

이제 가상의 “읽기 유창성 데이터”를 생성하고, 위에서 배운 모델들을 비교 분석해 보겠습니다.

(참고: jamovi의 GAMLj 모듈이나 Linear Mixed Models는 CS, AR1, UN, Toeplitz 등을 지원하지만, PRISM 도표나 고급 AD 모델은 R이 필요합니다. 따라서 R 코드를 중심으로 설명하고 jamovi 활용법을 덧붙이겠습니다.)

5.1 R을 이용한 모의 데이터 생성 및 시각화

R

# 필요한 라이브러리 로드
library(MASS)
library(nlme)
library(lattice)
library(ggplot2)

# --- 1. 가상의 읽기 유창성 데이터 생성 (Backstory) ---
# 50명의 학생, 4회 측정 (t=1, 2, 3, 4)
# 가정: 시간이 갈수록 분산이 커지고(이분산), 상관은 멀어질수록 낮아짐(AR 구조)
set.seed(2024)
N <- 50
Time <- 4
Mean_Vector <- c(40, 55, 65, 72) # 시간별 평균 점수 상승
Sigma <- matrix(0, Time, Time)

# 공분산 행렬 생성 (ARH(1) 구조와 유사하게 설정)
# 분산은 100, 140, 180, 220으로 증가
vars <- c(100, 140, 180, 220)
cor_rho <- 0.7 # 1-lag 상관계수
for(i in 1:Time){
  for(j in 1:Time){
    Sigma[i,j] <- sqrt(vars[i]) * sqrt(vars[j]) * (cor_rho^abs(i-j))
  }
}

# 데이터 생성
data_matrix <- mvrnorm(n = N, mu = Mean_Vector, Sigma = Sigma)
colnames(data_matrix) <- c("T1", "T2", "T3", "T4")
df_wide <- data.frame(ID = 1:N, data_matrix)

# Long format으로 변환 (분석용)
df_long <- reshape(df_wide, direction = "long", varying = list(2:5), 
                   v.names = "Score", timevar = "Time", times = 1:4)
df_long <- df_long[order(df_long$ID, df_long$Time), ]

# --- 2. 진단 시각화 ---

# (1) 프로파일 도표 (Profile Plot) [cite: 65]
ggplot(df_long, aes(x = Time, y = Score, group = ID)) +
  geom_line(alpha = 0.5) +
  geom_smooth(aes(group = 1), method = "loess", color = "red", size = 1.5) +
  labs(title = "학생 읽기 유창성 프로파일 도표", subtitle = "시간에 따른 개별 성장 곡선") +
  theme_minimal()

# (2) OSM (Ordinary Scatterplot Matrix) 
splom(df_wide[,2:5], 
      main = "OSM: 시점 간 산점도 행렬",
      pscales = 0, varname.cex = 0.8)

# (3) PRISM (간략화된 개념적 구현) [cite: 127]
# T1과 T3의 관계에서 T2를 통제한 편상관을 시각적으로 확인
res_T1_T2 <- residuals(lm(T1 ~ T2, data = df_wide))
res_T3_T2 <- residuals(lm(T3 ~ T2, data = df_wide))
plot(res_T1_T2, res_T3_T2, 
     main = "PRISM 예시 (T1 vs T3 | T2)", 
     xlab = "T1 | T2 residuals", ylab = "T3 | T2 residuals")
abline(lm(res_T3_T2 ~ res_T1_T2), col="blue")
# 설명: 이 점들이 무작위로 퍼져있으면(기울기가 0에 가까우면) 조건부 독립 -> AD 모델 적합

# --- 3. 모델 적합 및 비교 (gls 함수 사용) ---

# Model 1: 복합 대칭 (Compound Symmetry) - 너무 단순함 [cite: 232]
fit_cs <- gls(Score ~ factor(Time), data = df_long,
              correlation = corCompSymm(form = ~ 1 | ID))

# Model 2: AR(1) - 등분산 가정 [cite: 250]
fit_ar1 <- gls(Score ~ factor(Time), data = df_long,
               correlation = corAR1(form = ~ 1 | ID))

# Model 3: Unstructured (무구조) - 가장 유연함 [cite: 498]
fit_un <- gls(Score ~ factor(Time), data = df_long,
              correlation = corSymm(form = ~ 1 | ID),
              weights = varIdent(form = ~ 1 | Time)) # 이분산 허용

# Model 4: ARH(1) (이분산 AR1) - 현실적 대안 
fit_arh1 <- gls(Score ~ factor(Time), data = df_long,
                correlation = corAR1(form = ~ 1 | ID),
                weights = varIdent(form = ~ 1 | Time))

# --- 4. 모델 선택 (AIC/BIC 비교) [cite: 47, 55] ---
anova(fit_cs, fit_ar1, fit_un, fit_arh1)

5.2 jamovi에서의 분석 절차

jamovi에서는 [Linear Mixed Models] 메뉴를 사용합니다.

1. 데이터 불러오기: 위에서 생성한 데이터를 csv로 저장하여 jamovi에서 엽니다.
2. 분석 설정:
   • Dependent Variable: Score
   • Covariates: Time (또는 Factors로 설정)
   • Cluster: ID
3. Random Effects & Covariance Structure:
   • jamovi 옵션 중 Residual Covariance (혹은 Random Effect의 상관 구조) 설정 탭을 찾습니다.
   • 여기서 Unstructured, AR(1), Compound Symmetry 등을 선택할 수 있습니다.
4. 모델 비교: 각 모델별로 AIC/BIC 값을 확인하여 가장 작은 값을 가진 모델을 선택합니다. (보통 교육 성장 데이터는 ARH(1)이나 Unstructured가 잘 나옵니다).

6. 결론 및 제언 (Conclusions)

이 장의 저자들(Núñez-Antón & Zimmerman)은 다음과 같이 결론짓습니다.

1. 평균과 공분산을 함께 보라: 평균적인 성장만 보지 말고, 오차의 구조(변동성, 상관성)를 함께 모델링해야 정확한 추론이 가능합니다.
2. 균형을 찾아라: 무조건 복잡한 모델(Unstructured)이 좋은 것이 아닙니다. 데이터의 특성(비정상성, 계열 상관 등)을 반영하되, 가능한 한 파라미터 수가 적은 모델(예: Structured Antedependence 등)을 찾는 것이 최선입니다.
3. 시각화를 믿어라: 통계 수치(AIC/BIC)도 중요하지만, 프로파일 도표와 PRISM 등을 통해 데이터의 실제 패턴을 눈으로 확인하는 것이 강력한 진단 도구가 됩니다.

[User를 위한 다음 단계]

위의 R 코드를 복사하여 RStudio에서 실행해 보시겠습니까? 생성된 그래프(프로파일 도표, PRISM)를 통해 방금 설명한 “조건부 독립” 개념이 시각적으로 어떻게 나타나는지 직접 확인해 보시면 이해가 훨씬 빠를 것입니다.

참고문헌 (References)

• Akaike, H. (1974). A new look at the statistical model identification. IEEE Transactions on Automatic Control, 19, 716-723.
• Diggle, P. J., Heagerty, P. J., Liang, K.-Y., & Zeger, S. L. (2002). Analysis of longitudinal data (2nd ed.). Oxford University Press.
• Gabriel, K. R. (1962). Ante-dependence analysis of an ordered set of variables. Annals of Mathematical Statistics, 33, 201-212.
• Goldstein, H. (2003). Multilevel statistical model (3rd ed.). Arnold.
• Kenward, M. G. (1987). A method for comparing profiles of repeated measurements. Applied Statistics, 36, 296-308.
• Laird, N. M., & Ware, J. H. (1982). Random-effects models for longitudinal data. Biometrics, 38, 963-974.
• Núñez-Antón, V., & Zimmerman, D. L. (2001). Parametric modeling of growth curve data: An overview (with discussion). Test, 10, 1-73.
• Zimmerman, D. L. (2000). Viewing the correlation structure in longitudinal data through a PRISM. The American Statistician, 54, 310-318.
• Zimmerman, D. L., & Núñez-Antón, V. (2010). Antedependence models for longitudinal data. Chapman & Hall/CRC.

안녕하세요!

오늘은 “고정 효과(Fixed Effects)와 무선 효과(Random Effects)의 선택, 그리고 하이브리드 모델”에 대해 살펴보겠습니다. “학교 현장의 데이터”를 예시로 들어 직관적인 설명과 수리적 엄밀함을 모두 갖춘 형태로 재구성해 드리겠습니다.

분석 도구로는 jamovi의 사용법을 설명하되, jamovi의 기반이 되는 R 코드를 함께 제시하여 모의 데이터 생성부터 분석, 시각화까지 완벽하게 구현해 드리겠습니다.

1. 들어가며: 우리는 왜 고민하는가?

연구자 여러분, 우리가 다층모형(Multilevel Modeling)을 사용할 때 가장 흔하게 마주치는 질문이 있습니다.

“교사(또는 학교) 효과를 고정 효과(Fixed Effect)로 볼 것인가, 무선 효과(Random Effect)로 볼 것인가?”

이 선택은 단순히 통계적 취향의 문제가 아닙니다. 이 선택에 따라 모수의 해석이 달라지고, 추정의 효율성(Efficiency)과 편향(Bias) 사이의 중대한 트레이드오프(Trade-off)가 발생하기 때문입니다.

• 고정 효과(FE): 각 그룹(학교)을 고유한 특성을 가진 개별적 존재로 보고, 그 자체를 변수로 투입합니다.
• 무선 효과(RE): 각 그룹을 모집단에서 추출된 하나의 표본으로 보고, 그룹 효과가 정규분포를 따른다고 가정합니다.

이 글에서는 이 두 모델의 차이를 명확히 하고, 하이브리드 모델(Hybrid Model)이라는 아주 매력적인 대안까지 함께 알아보겠습니다.

2. 시나리오: “방과 후 자습 시간은 성적을 올리는가?”

이해를 돕기 위해 가상의 교육 데이터를 만들어 봅시다.

• 연구 문제: 학생의 ‘방과 후 자습 시간($x$)’이 ‘수학 성취도($y$)’에 미치는 영향.
• 데이터 구조: 학생($i$)들은 여러 학교($j$)에 소속되어 있음 (2수준 데이터).
• 숨겨진 함정(Confounder): 사실 ‘학교의 면학 분위기($u$)’가 좋은 학교일수록, 학생들의 자습 시간도 길고 성적도 높습니다. 만약 이 학교 효과를 제대로 통제하지 않으면, 자습 시간의 순수한 효과를 과대평가할 수 있습니다.

먼저, R을 사용하여 이 시나리오에 맞는 모의 데이터를 생성해보겠습니다. (이 코드를 R Studio나 jamovi의 Rj Editor에서 실행하면 데이터를 얻을 수 있습니다.)

R

# [R Code] 데이터 생성
set.seed(123)
n_schools <- 50    # 학교 수
n_students <- 20   # 학교당 학생 수

# 학교 효과 (면학 분위기): 학교마다 다름
school_effect <- rnorm(n_schools, mean = 0, sd = 10)

# 데이터 프레임 생성
data <- data.frame()

for(j in 1:n_schools) {
  # 자습 시간(X): 학교 효과와 상관이 있도록 설정 (중요! Assumption 위반 상황 연출)
  # 면학 분위기가 좋은 학교 학생들이 공부를 더 많이 함
  study_hours <- rnorm(n_students, mean = 5 + 0.1 * school_effect[j], sd = 1)
  
  # 수학 성취도(Y): 기본점수 + 자습효과(2점) + 학교효과 + 오차
  math_score <- 50 + 2 * study_hours + school_effect[j] + rnorm(n_students, mean = 0, sd = 5)
  
  temp <- data.frame(school_id = factor(j), study_hours, math_score)
  data <- rbind(data, temp)
}

# CSV로 저장 (jamovi에서 불러오기 위함)
write.csv(data, "chap05.csv", row.names = FALSE)

3. 두 모델의 핵심 개념 비교

3.1. 무작위 효과 모델 (Random Effects Model: RE)

무작위 효과 모델은 우리가 흔히 쓰는 다층모형의 기본 형태입니다.

수식은 다음과 같습니다.

$y_{ij} = \beta_{0} + \beta_{1}x_{1ij} + u_{0j} + \epsilon_{ij}$

여기서 $u_{0j}$는 학교별 효과(Intercept)인데, RE 모델은 이 $u_{0j}$가 평균이 0이고 분산이 $\sigma_{u0}^2$인 정규분포를 따른다고 가정합니다.

• 특징:

1. 부분 풀링(Partial Pooling): 전체 학교의 정보를 빌려와서 추정하므로, 데이터가 적은 학교의 추정치도 안정적입니다(수축 효과, Shrinkage).
2. 효율성: 추정해야 할 파라미터가 적어(분산만 추정하면 됨) 통계적으로 효율적입니다.
3. 치명적 가정: 학교 효과($u_{0j}$)와 설명변수($x_{1ij}$, 자습시간)가 서로 독립이어야(상관이 없어야) 합니다.
   • 문제점: 위 시나리오처럼 ‘면학 분위기가 좋은 학교($u$↑)’일수록 ‘자습 시간($x$↑)’이 길다면, 이 가정이 깨지고 결과는 편향(Bias)됩니다.

3.2. 고정 효과 모델 (Fixed Effects Model: FE)

고정 효과 모델은 각 학교를 고유한 더미 변수(Dummy Variable)로 취급하거나, 학교 평균을 빼버리는 방식(De-meaning)을 사용합니다.

$y_{ij} – \bar{y}_{j} = \beta_{1}(x_{1ij} – \bar{x}_{1j}) + (\epsilon_{ij} – \bar{\epsilon}_{j})$

• 직관적 이해: 학교 간의 차이는 아예 보지 않겠다는 뜻입니다. 오로지 “같은 학교 내에서(Within-school)” 자습 시간이 늘어날 때 성적이 오르는지만 봅니다.
• 특징:
  1. 편향 제거: 학교 효과($u_{0j}$)가 자습 시간($x$)과 상관이 있든 없든 상관없습니다. 학교 고유의 특성을 완벽히 통제합니다.
  2. 비효율성: 학교 수만큼의 더미 변수를 만드는 셈이므로 자유도(df)를 많이 잡아먹습니다.
  3. 한계: 시간(또는 그룹)에 따라 변하지 않는 변수(예: 학교의 설립 유형, 지역)의 효과는 추정할 수 없습니다. 다 삭제되기 때문입니다.

4. 분석 실습: jamovi & R

이제 jamovi를 이용해 두 모델을 분석하고 비교해 보겠습니다.

4.1. 고정 효과(FE) 분석 (in jamovi)

jamovi에는 ‘Panel Fixed Effects’ 전용 버튼은 없지만, 일반 선형 회귀(Linear Regression)에서 학교 ID를 더미 변수로 넣거나, GAMLj 모듈을 사용하여 구현할 수 있습니다. 가장 교과서적인 ‘Within Estimator’ 방식은 변수를 중심화(Centering)하여 분석하는 것입니다.

[jamovi 절차]

1. 데이터 열기: 위에서 만든 chap05.csv를 엽니다.
2. 변수 계산 (Compute):
   • mean_study_hours = VMEAN(study_hours, group_by=school_id) (학교별 평균 자습시간)
   • within_study_hours = study_hours – mean_study_hours (학교 평균 중심화)
3. 분석 (Linear Regression):
   • 종속변수: math_score
   • 공변량: within_study_hours (이때 학교 간 차이는 이미 제거되었습니다.)
   • (참고: 엄밀한 FE 추정치는 더미변수를 넣어야 하지만, 계수값 $\beta_1$은 중심화된 변수 회귀와 동일한 원리를 갖습니다.)

4.2. 무작위 효과(RE) 분석 (in jamovi)

[jamovi 절차]

1. 메뉴: Linear Models > Mixed Model 선택.
2. 설정:
   • Dependent Variable: math_score
   • Covariates: study_hours
   • Cluster: school_id
   • Random Effects: Intercept 체크 (학교별 절편을 무작위로 가정).

4.3. 결과 비교 및 하우스만 검정 (Hausman Test)

전통적으로 이 두 모델 중 무엇을 쓸지 결정할 때 하우스만 검정을 사용합니다.

• 귀무가설($H_0$): RE와 FE의 추정치 차이가 없다 (즉, RE가 효율적이고 편향도 없으니 RE를 써라).
• 대립가설($H_1$): 차이가 있다 (즉, RE는 편향되었으니 FE를 써라).

jamovi의 기본 메뉴에는 하우스만 검정이 없으므로, R 코드를 통해 수행하거나 뒤에 소개할 하이브리드 모델로 대체하여 판단합니다.

R

# [R Code] FE vs RE 및 하우스만 검정
library(plm)

# 1. 고정 효과 모델 (Fixed Effects)
fe_model <- plm(math_score ~ study_hours, data=data, index=c("school_id"), model="within")

# 2. 무작위 효과 모델 (Random Effects)
re_model <- plm(math_score ~ study_hours, data=data, index=c("school_id"), model="random")

# 3. 하우스만 검정
phtest(fe_model, re_model)


	Hausman Test

data:  math_score ~ study_hours
chisq = 192.49, df = 1, p-value < 2.2e-16
alternative hypothesis: one model is inconsistent

해석: 만약 p-value가 0.05보다 작다면, RE 모델의 가정(학교효과와 $X$가 독립)이 기각된 것입니다. 즉, 편향이 발생했으므로 FE를 써야 한다는 신호입니다.

5. 최선의 대안: 하이브리드 모델 (Hybrid Model)

많은 학자들은 FE와 RE 중 하나만 고르는 이분법 대신, 두 모델의 장점을 합친 하이브리드 모델을 추천합니다. 이 모델은 ‘Mundlak 모델’ 또는 ‘Group-Mean Centering’ 방법으로도 불립니다.

5.1. 하이브리드 모델의 원리

설명변수 $x_{ij}$를 두 부분으로 쪼개서 모델에 넣습니다.

$y_{ij} = \beta_{0} + \beta_{W}(x_{ij} – \bar{x}_{j}) + \beta_{B}(\bar{x}_{j}) + u_{0j} + \epsilon_{ij}$

• $\beta_{W}$ (Within effect): 학교 내 효과. 학생이 자기 학교 평균보다 더 공부했을 때 성적이 얼마나 오르는가? (= 고정 효과 추정치와 동일).
• $\beta_{B}$ (Between effect): 학교 간 효과. 공부를 많이 시키는 학교가 성적이 더 높은가?
• 장점:
  1. 변수 간 상관으로 인한 편향 문제 해결 (FE의 장점).
  2. 학교 수준 변수나 학교 간 차이도 추정 가능 (RE의 장점).
  3. $\beta_{W}$와 $\beta_{B}$가 같은지 검정하여(Wald test), 맥락 효과(Contextual Effect)가 있는지 볼 수 있음.

5.2. jamovi에서 하이브리드 모델 구현하기

이것이 오늘 강의의 핵심 꿀팁입니다. 별도의 코딩 없이 jamovi 메뉴만으로 가능합니다.

[jamovi 실습 절차]

1. 변수 생성:
   • 앞서 만든 mean_study_hours (학교 평균, Between 성분)
   • 앞서 만든 within_study_hours (개인 편차, Within 성분)
2. Mixed Model 분석:
   • Linear Models > Mixed Model
   • Dependent Variable: math_score
   • Covariates: within_study_hours 그리고 mean_study_hours 두 개를 모두 넣습니다.
   • Cluster: school_id
   • Random Effects: Intercept

[결과 해석]

• within_study_hours의 계수: 이것이 바로 순수한 개인 노력의 효과입니다. 학교 분위기(교란변수)가 통제된 FE 추정치와 같습니다.
• mean_study_hours의 계수: 학교 간의 차이 효과입니다. 만약 Within 계수와 Between 계수가 크게 다르다면, 단순히 개인 노력이 아니라 학교 분위기가 성적에 영향을 미치고 있음을 시사합니다.

6. 결론 및 제언

우리가 살펴본 내용을 요약하면 다음과 같습니다.

1. 무작위 효과(RE)는 효율적이지만, 그룹 효과와 설명변수가 관련이 있을 경우 편향될 위험이 있습니다 (교육 데이터에서는 흔한 일입니다).
2. 고정 효과(FE)는 편향을 제거해주지만, 그룹 수준의 변수(예: 사립/공립 여부)를 분석할 수 없고 비효율적일 수 있습니다²¹.
3. 하이브리드 모델은 설명변수를 ‘그룹 내 편차(Within)’와 ‘그룹 평균(Between)’으로 분해하여 모델에 투입함으로써, 편향 제거와 정보 활용이라는 두 마리 토끼를 모두 잡을 수 있는 강력한 방법입니다.

WaurimaL의 조언:

학교 데이터를 분석할 때, 무조건 RE만 돌리지 마세요. 설명변수를 학교 평균 중심으로(Group-mean centering) 변환하여 투입하는 하이브리드 접근법을 사용한다면, 훨씬 더 풍부하고 정확한 교육적 시사점을 얻을 수 있습니다.

📚 참고문헌 (APA Style)

• Townsend, Z., Buckley, J., Harada, M., & Scott, M. A. (2013). The choice between fixed and random effects. In The SAGE Handbook of Multilevel Modeling (pp. 73-88). SAGE Publications.
• Allison, P. D. (2009). Fixed effects regression models. SAGE.
• Bafumi, J., & Gelman, A. (2006). Fitting multilevel models when predictors and group effects correlate. Paper presented at the annual meeting of the Midwest Political Science Association, Chicago, IL.
• Hausman, J. A. (1978). Specification tests in econometrics. Econometrica, 46, 1251–1271.
• Raudenbush, S. W., & Bryk, A. S. (2002). Hierarchical linear models: Applications and data analysis methods (2nd ed.). SAGE Publications.
• Wooldridge, J. (2010). Econometric analysis of cross section and panel data. MIT Press.