안녕하세요!

오늘은 비선형 혼합 효과 모델(Nonlinear Mixed Effects Models)에 대해 살펴보겠습니다. “학교 현장의 데이터”를 예시로 들어 직관적인 설명과 수리적 엄밀함을 모두 갖춘 형태로 재구성해 드리겠습니다.

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

1. 왜 ‘비선형(Nonlinear)’ 모델이 필요할까요?

우리가 흔히 쓰는 ‘선형 모델(Linear Models)’은 계산이 쉽고 해석이 명확하다는 장점이 있습니다. 하지만 학교 현장에서 일어나는 실제 현상들은 직선으로 설명되지 않는 경우가 훨씬 많습니다.

• 선형 모델(Empirical Models): 관찰된 데이터를 잘 설명하기 위해 직선을 긋는 ‘경험적 모델’입니다. 데이터 범위 안에서는 잘 맞을지 몰라도, 범위를 벗어나면 예측력이 떨어질 수 있습니다.
• 비선형 모델(Mechanistic/Scientific Models): 데이터가 왜 그렇게 만들어졌는지, 즉 ‘생성 메커니즘’에 집중하는 ‘과학적 모델’입니다. 성장 곡선이나 학습 속도처럼 자연스러운 변화 과정을 반영하므로, 데이터가 없는 구간에 대해서도 더 정확한 예측이 가능합니다.

비선형 모델의 장점:

1. 실제 관계를 더 잘 근사합니다.
2. 관찰 데이터 범위를 벗어난 값에 대해서도 더 신뢰할 수 있는 예측을 제공합니다.
3. 파라미터(Parameter) 자체가 자연적인 물리적/교육적 의미를 갖습니다 (예: 학생의 최대 학습 잠재력, 학습 속도 등).

2. 학교 현장 스토리: “우리 아이들의 단어 암기량은 어떻게 변할까?”

초등학교 1학년 학생들이 매일 아침 10분씩 영어 단어를 외운다고 가정해 봅시다. 시간이 지날수록 외운 단어 수는 늘어나겠지만, 처음에는 천천히 늘다가 어느 순간 급격히 늘고, 결국 인간의 한계 때문에 정체기에 접어들 것입니다. 이러한 ‘S자 곡선(Logistic Curve)’은 직선인 선형 모델로는 설명할 수 없습니다.

[모의 데이터 생성 및 분석을 위한 R 코드]

먼저 학교 현장과 유사한 데이터를 만들어 보겠습니다.

R

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

# 1. 모의 데이터 생성 (학생 50명, 10회 측정)
set.seed(123)
n_students <- 50
n_days <- 10
data <- data.frame(
  student = rep(1:n_students, each = n_days),
  day = rep(1:n_days, n_students)
)

# 비선형 파라미터 설정 (Asym: 최대 암기량, xmid: 가속 지점, scal: 증가 속도)
# 학생마다 개인차(Random Effects) 부여
asym_pop <- 100  # 평균 최대 암기량 100개
xmid_pop <- 5    # 평균 5일째에 급격히 상승
scal_pop <- 1    # 증가 기울기 계수

# 학생별 랜덤 효과 생성
random_effects <- rnorm(n_students, 0, 5) 

# 로지스틱 성장 곡선 적용
data$word_count <- with(data, {
  asym_i <- asym_pop + random_effects[student]
  asym_i / (1 + exp((xmid_pop - day) / scal_pop)) + rnorm(nrow(data), 0, 2)
})

# 데이터 확인
head(data)

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

3. jamovi에서의 분석 방법

jamovi에서 비선형 혼합 모델을 직접 수행하려면 GAMLj 모듈이나 Rj Editor를 활용하는 것이 좋습니다. 기본 메뉴에 비선형 혼합 모델이 없을 경우, 다음과 같이 Rj Editor를 통해 분석합니다.

jamovi 분석 단계:

1. Library 설치: Rj - Editor에서 nlme 패키지를 호출합니다.
2. 모델 설정: 학생별로 다른 성장 곡선을 갖도록 ‘혼합 효과(Mixed Effects)’를 설정합니다.
3. 수식 입력:

# nlme 패키지를 이용한 비선형 혼합 모델 분석
library(nlme)
model <- nlme(word_count ~ SSlogis(day, Asym, xmid, scal),
              data = data,
              fixed = Asym + xmid + scal ~ 1,
              random = Asym ~ 1 | student,
              start = c(Asym = 100, xmid = 5, scal = 1))
summary(model)

# Asym: 학생이 도달할 수 있는 '최종 암기량' (상한선)
# xmid: 학습 효율이 최대가 되는 '시점'
# scal: 암기량이 증가하는 '가속도'

4. 결과 해석 및 시각화

분석 결과, 우리는 학생 전체의 평균적인 학습 곡선뿐만 아니라, 특정 학생이 다른 학생에 비해 얼마나 더 빨리 배우는지(xmid), 혹은 잠재력이 얼마나 큰지(Asym)를 개별적으로 파악할 수 있습니다.

R

# 시각화 코드
library(tidyverse)
ggplot(data, aes(x = day, y = word_count, group = student)) +
  geom_line(alpha = 0.3, color = "blue") + # 개별 학생 곡선
  stat_summary(aes(group = 1), fun = mean, geom = "line", size = 1.5, color = "red") + # 전체 평균
  labs(title = "학생별 영어 단어 암기량 성장 곡선",
       x = "학습 일수(Day)",
       y = "암기한 단어 수") +
  theme_minimal()

5. 주의할 점 (Limitations)

비선형 모델은 강력하지만 몇 가지 주의사항이 있습니다:

• 복잡성: 데이터 생성 메커니즘을 정확히 모르면 모델을 세우기 어렵습니다.
• 계산의 어려움: 선형 모델과 달리 반복적인 계산 과정이 필요하며, 때로는 수렴(Convergence)되지 않아 결과가 나오지 않을 수도 있습니다.
• 샘플 사이즈: 파라미터를 안정적으로 추정하기 위해 충분한 양의 데이터가 필요합니다.

6. 참고문헌 (APA Style)

Wu, L., & Liu, W. (2010). Nonlinear models. In J. J. Hox & J. K. Roberts (Eds.), The SAGE Handbook of Multilevel Modeling (pp. 249-266). SAGE Publications. +1

안녕하세요!

오늘은 “다층모형과 인과추론(Multilevel Models and Causal Inference)”에 대해 살펴보겠습니다. “학교 현장의 데이터”를 예시로 들어 직관적인 설명과 수리적 엄밀함을 모두 갖춘 형태로 재구성해 드리겠습니다.

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

1. 들어가며: 왜 다층모형인가?

전통적인 회귀분석은 모든 학생이 서로 독립적이라고 가정합니다. 하지만 교육 현장은 그렇지 않습니다. 같은 학교, 같은 반 학생들은 급훈, 담임 선생님, 학교 분위기 등을 공유합니다. 이를 위계적 구조(Hierarchical Structure) 또는 집단 의존성(Group Dependencies)이라고 합니다.

인과추론(Causal Inference)의 관점에서, 데이터가 이러한 계층 구조를 가질 때 다층모형(Multilevel Model)을 사용하는 것은 단순한 통계적 선호가 아니라, 편향(Bias)을 줄이고 정확한 표준오차를 추정하기 위한 필수 전략입니다.

2. 인과추론의 기초 개념과 교육적 예시

본격적인 분석에 앞서, 인과추론의 핵심 개념을 학교 상황에 빗대어 정의해 봅시다.

2.1 잠재적 결과 (Potential Outcomes)

어떤 학생 철수($i$)가 있습니다.

• $y_i(1)$: 철수가 ‘방과후 보충수업($z=1$)’을 들었을 때의 성적
• $y_i(0)$: 철수가 ‘방과후 보충수업($z=0$)’을 듣지 않았을 때의 성적

인과 효과(Causal Effect)는 이 둘의 차이 $y_i(1) – y_i(0)$입니다. 하지만 현실에서 우리는 철수가 수업을 듣거나, 듣지 않거나 둘 중 하나의 결과만 볼 수 있습니다. 이를 “인과추론의 근본적인 문제(Missing Data Problem)”라고 합니다.

2.2 SUTVA (Stable Unit Treatment Value Assumption)

이 가정은 “철수가 보충수업을 받았는지 여부가, 옆 짝꿍 영희의 성적에 영향을 주지 않아야 한다(상호간섭 없음)”는 것입니다.

• 문제점: 학교에서는 이 가정이 자주 깨집니다. 철수가 보충수업에서 배운 내용을 영희에게 알려줄 수 있기 때문입니다. 이를 해결하기 위해 집단(학교/학급) 단위 무선화가 권장되기도 합니다.

3. 연구 설계에 따른 다층모형 적용

3.1 무선화 실험 (Randomized Experiments)

가장 이상적인 상황입니다. 처치(Treatment)가 무작위로 배정되면, 평균적으로 두 집단은 성향이 비슷해집니다($E[y_i(1)] = E[y_i|z_i=1]$).

A. 개인 단위 무선배정 (학생별 제비뽑기)

학생들에게 무작위로 새로운 ‘독서 프로그램’을 배정했습니다. 하지만 학생들은 학교($j$)라는 집단에 속해 있습니다. 학교마다 평균 독서 능력이 다를 수 있으므로, 이를 반영한 다층모형(Random Intercept Model)이 필요합니다.

$y_{ij} = \mu + \alpha_j + \tau z_{ij} + \epsilon_{ij}$

• $\alpha_j$: $j$번째 학교의 고유한 특성(학교 효과, 랜덤 절편)
• $\tau$: 독서 프로그램의 효과 (우리가 알고 싶은 값)

이 모형을 쓰면 학교 간 차이($\alpha_j$)를 통제하고 순수한 프로그램 효과($\tau$)를 더 정밀하게 추정할 수 있습니다.

B. 집단 단위 무선배정 (학교별 제비뽑기)

교육 정책 연구에서는 흔히 “A학교는 실험군, B학교는 대조군”으로 배정합니다. 이를 군집 무선화(Cluster Randomized Experiments)라고 합니다.

• 이유: ‘학교 폭력 예방 캠페인’처럼 학교 전체 분위기를 바꾸는 처치는 학생 개인별로 쪼개서 적용할 수 없기 때문입니다.
• 분석: 처치 변수($z_j$)가 학생 수준($i$)이 아닌 학교 수준($j$)에 들어갑니다.

3.2 관찰 연구 (Observational Studies)

현실적으로 무선 배정이 불가능할 때(예: 사립학교 진학 효과), 우리는 무시가능성(Ignorability) 가정을 도입합니다. 즉, “부모의 소득, 지능 등 공변량($x$)이 같다면, 사립학교와 공립학교 학생은 비교 가능하다”고 가정하는 것입니다.

• 성향점수(Propensity Score) 활용: 다층 구조에서는 성향점수를 추정할 때도 다층모형을 사용하는 것이 좋습니다.

4. [실습] jamovi & R을 활용한 다층 인과 분석

이제 가상의 시나리오를 통해 실제 데이터를 생성하고 분석해 보겠습니다.

4.1 시나리오: “아침 독서 마라톤” 효과 분석

연구 배경: 경기도 교육청은 초등학생의 어휘력 향상을 위해 매일 아침 20분간 책을 읽는 ‘아침 독서 마라톤’ 프로그램을 개발했습니다.

연구 설계:

• 총 20개 학교, 학교당 30명의 학생(총 600명).
• 군집 무선화(Cluster RCT): 학교 단위로 제비뽑기를 하여 10개 학교는 ‘프로그램 시행(Treatment)’, 10개 학교는 ‘기존 자습(Control)’을 하도록 했습니다.
• 데이터 구조:
  • Level 1: 학생 (사후 어휘력 점수 score)
  • Level 2: 학교 (school_id)
  • 처치: program (1=시행, 0=미시행)

4.2 R을 이용한 모의 데이터 생성

jamovi는 R 기반이므로, 아래 코드로 데이터를 생성하여 CSV로 저장한 뒤 jamovi에서 불러오면 됩니다.

R

# 필수 라이브러리 로드
library(lme4)
library(tidyverse)

set.seed(2026) # 재현성을 위한 시드 설정

# 1. 파라미터 설정
n_schools <- 20       # 학교 수
n_students <- 30      # 학교당 학생 수
n_total <- n_schools * n_students

# 2. 학교 수준 효과 (Level 2)
# 학교마다 평균 어휘력이 다름 (표준편차 5)
school_intercept <- rnorm(n_schools, mean = 0, sd = 5)

# 처치 배정 (학교 단위 무선화)
# 1~10번 학교: 통제군(0), 11~20번 학교: 실험군(1)
school_treatment <- c(rep(0, 10), rep(1, 10))

# 학교 데이터 프레임
school_data <- data.frame(
  school_id = 1:n_schools,
  school_eff = school_intercept,
  program = school_treatment
)

# 3. 학생 수준 데이터 생성 (Level 1)
data <- data.frame(
  student_id = 1:n_total,
  school_id = rep(1:n_schools, each = n_students)
)

# 학교 정보 병합
data <- left_join(data, school_data, by = "school_id")

# 4. 결과 변수 생성 (어휘력 점수)
# 기본 점수 70점 + 프로그램 효과 8점 + 학교 효과 + 개인 오차(sd=8)
# y_ij = 70 + 8 * z_j + u_j + e_ij
data <- data %>%
  mutate(
    error = rnorm(n_total, mean = 0, sd = 8),
    score = 70 + 8 * program + school_eff + error
  )

# 팩터 변환
data$school_id <- as.factor(data$school_id)
data$program <- factor(data$program, levels = c(0, 1), labels = c("Control", "Treatment"))

# 데이터 확인
head(data)

4.3 jamovi 분석 절차

이 데이터는 학교 간 차이(School Effect)가 존재하고, 처치가 학교 단위로 부여되었으므로 다층모형(Linear Mixed Model)을 사용해야 정확합니다.

Step 1: 데이터 탐색 및 시각화

분석 전에 데이터의 구조를 눈으로 확인해야 합니다.

1. jamovi 메뉴: Exploration > Descriptives
2. Variables에 score를 넣고, Split by에 program을 넣습니다.
3. Box Plot: 학교별 차이를 보기 위해 Box plot을 체크하고, X축에 program을 둡니다. (※ jamovi 기본 기능으로는 학교별 boxplot을 한 번에 그리기 어려우므로 R 모듈인 seolmatrix나 scatr 모듈을 설치하여 시각화하면 좋습니다.)

[R 시각화 코드]

R

# 학교별 점수 분포 시각화
ggplot(data, aes(x = school_id, y = score, fill = program)) +
  geom_boxplot() +
  theme_minimal() +
  labs(title = "학교별 어휘력 점수 분포", y = "어휘력 점수", x = "학교 ID")

해석: 상자 그림을 보면 같은 처치 집단 내에서도 학교마다 점수의 높낮이가 다름을 알 수 있습니다. 이것이 바로 $\alpha_j$ (학교 효과)입니다.

Step 2: 다층모형 분석 (Linear Mixed Models)

1. 모듈 선택: 상단 메뉴에서 Linear Models > Mixed Model을 클릭합니다. (보이지 않으면 jamovi Library에서 GAMLj 모듈을 설치하는 것을 강력 추천합니다. 여기서는 기본 Mixed Model 기준으로 설명합니다.)
2. 변수 설정:
   • Dependent Variable (종속변수): score
   • Covariates (공변량) 또는 Factors: program (처치 변수)
   • Cluster (군집 변수): school_id
3. Random Effects (랜덤 효과) 설정:
   • 왼쪽의 program을 오른쪽으로 옮기지 않고, Intercept만 Random Coefficients에 둡니다. (기본적으로 (Intercept | school_id)로 설정됨)
   • 이는 학교마다 평균 점수(절편)가 다름을 허용하는 것입니다.
4. Fixed Effects (고정 효과) 설정:
   • program을 Model Terms에 넣습니다. 이것이 우리가 알고 싶은 ‘독서 마라톤 효과’입니다.

Step 3: 결과 해석

jamovi의 결과표(Estimates)는 다음과 유사하게 나옵니다.

• Fixed Effects: program의 Estimate가 약 10.162입니다. 즉, 독서 마라톤을 한 학교 학생들이 하지 않은 학교보다 평균적으로 약 10.162점 더 높은 어휘력을 보입니다. $p < .05$이므로 통계적으로 유의합니다.
• Random Components (Variance):
  • $\sigma^2_{school}$ (School Intercept): 학교 간 분산. 이 값이 0보다 크다면 학교 효과가 존재한다는 뜻입니다.
  • ICC (Intraclass Correlation Coefficient): 전체 분산 중 학교가 설명하는 비율입니다.

5. 심화: 불응(Noncompliance)과 도구변수(IV)

실험을 했는데, 독서 프로그램을 하라고 배정받은 학교의 일부 학생이 땡땡이를 쳤다면(Noncompliance) 어떻게 될까요? 이때는 “배정된 상태($z$)”를 도구변수(Instrument)로 사용하여, 실제 “참여한 상태($d$)”의 효과를 추정해야 합니다.

jamovi/R 구현 (2단계 최소자승법 개념)

1. 1단계: 실제 참여 여부($d$)를 배정 여부($z$)로 예측합니다.
   $d_{ij} = \gamma_j + \delta z_{ij} + \nu_{ij}$
2. 2단계: 1단계에서 예측된 참여값($\hat{d}$)을 사용하여 점수($y$)를 예측합니다.
   $y_{ij} = \alpha_j + \tau \hat{d}_{ij} + \epsilon_{ij}$

이 분석은 jamovi의 sem (구조방정식) 모듈이나 R의 AER 패키지(ivreg)를 통해 수행할 수 있습니다. 중요한 건 배정($z$)은 오직 참여($d$)를 통해서만 결과($y$)에 영향을 미쳐야 한다(배제 제한)는 가정입니다.

6. 결론

다층모형을 활용한 인과추론은 교육 현장과 같이 “집단 속에 개인이 속한 데이터”를 분석할 때 가장 강력한 도구입니다.

1. 설계: 가능하다면 학교 단위 무선화(Cluster RCT)가 상호간섭(SUTVA 위배) 문제를 피하는 데 유리합니다.
2. 분석: 단순히 평균을 비교하는 t-test 대신, 학교의 무선 절편(Random Intercept)을 포함한 혼합 모형을 사용해야 표준오차의 과소추정을 막을 수 있습니다.
3. 해석: 결과는 “개인 수준의 효과”인지 “학교 수준의 효과”인지 명확히 구분하여 해석해야 합니다.

이 장의 내용이 여러분의 연구에 튼튼한 방법론적 기초가 되기를 바랍니다.

참고문헌 (APA Style)

• Almond, D., Chay, K., & Lee, D. (2005). The costs of low birth weight. The Quarterly Journal of Economics, 120(3), 1031-1083.
• Angrist, J. D., Imbens, G. W., & Rubin, D. B. (1996). Identification of causal effects using instrumental variables. Journal of the American Statistical Association, 91(434), 444-472.
• Cornfield, J. (1978). Randomization by group: A formal analysis. American Journal of Epidemiology, 108(2), 100-102.
• Gelman, A., & Hill, J. (2007). Data analysis using regression and multilevel/hierarchical models. Cambridge University Press.
• Hill, J. (2013). Multilevel models and causal inference. In The SAGE Handbook of Multilevel Modeling (Chapter 12, pp. 201-219).
• Hong, G., & Raudenbush, S. W. (2006). Evaluating kindergarten retention policy: A case study of causal inference for multilevel observational data. Journal of the American Statistical Association, 101(475), 901-910.
• Kim, J., & Seltzer, M. (2007). Causal inference in multilevel settings in which selection processes vary across schools (Tech. Rep.). CRESST, UCLA.
• Rubin, D. B. (1978). Bayesian inference for causal effects: The role of randomization. The Annals of Statistics, 6(1), 34-58.
• Rubin, D. B. (1990). Formal modes of statistical inference for causal effects. Journal of Statistical Planning and Inference, 25(3), 279-292.
• Slavin, R. E., Madden, N. A., Dolan, L. J., & Wasik, B. A. (1996). Every child, every school: Success for all. Corwin Press.

안녕하세요!

오늘은 “다층 연구 설계의 최적화(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.