Lasso / Elastic Net / glmmLasso — 종단 데이터의 변수 선택

고차원 종단 데이터에서 중요한 변수를 자동으로 골라내기

종단 데이터 분석에서 변수 선택이 왜 중요한지, L1/L2 정규화의 수학적 원리, glmmLasso로 혼합 모델과 Lasso를 결합하는 방법, AI Agent 개인화 실험의 고차원 피처 선택 실무 예시를 Python(sklearn, statsmodels)과 R(glmnet, glmmLasso)로 구현한다.

Statistics
Machine Learning
Regularization
저자

Kwangmin Kim

공개

2026년 03월 08일

1 Lasso / Elastic Net / glmmLasso

시리즈 위치: 12 ML Overview에서 네 번째 세부 주제. 종단 데이터에서 고차원 피처 중 “진짜 중요한 변수”를 자동으로 골라내는 정규화 기법을 다룬다.


1.1 1. 왜 정규화인가: 고차원 종단 데이터의 저주

1.1.1 1.1 종단 데이터가 고차원이 되는 이유

전통적인 횡단면(cross-sectional) 분석에서는 변수 수 \(p\)가 적당하고, 표본 크기 \(n\)이 충분해서 OLS나 MLE로 안정적인 추정이 가능하다. 그러나 종단 데이터에서는 상황이 급격히 달라진다.

피처 폭발의 원인:

원인 예시 피처 증가량
시간 변동 공변량 (time-varying covariates) 세션별 응답 속도, 클릭 수, 체류 시간 \(p \times T\)
교호작용 (interactions) 피처 \(\times\) 시간, 피처 \(\times\) 그룹 \(\binom{p}{2}\)
시차 변수 (lag terms) \(x_{t-1}, x_{t-2}, \ldots\) \(p \times L\)
다항 시간 항 (polynomial time) \(t, t^2, t^3\) \(3 \times p\)
파생 피처 (derived features) 이동 평균, 변화율, 누적합 \(p \times k\)

AI Agent 개인화 실험 예시를 생각해 보자. 사용자 200명이 각각 10개 세션에 참여하고, 각 세션에서 50개 피처를 수집한다. 여기에 교호작용과 시차 변수를 추가하면:

  • 원본 피처: 50개
  • 1차 교호작용: \(\binom{50}{2} = 1{,}225\)
  • 시차 1~2: \(50 \times 2 = 100\)
  • 시간 다항항: \(50 \times 2 = 100\)

총 후보 피처가 1,475개에 달하지만, 유효 표본 크기(독립 군집 수)는 200명에 불과하다. \(p \gg n_{\text{cluster}}\)인 상황에서 모든 변수를 모델에 넣으면:

  1. 과적합(overfitting): 학습 데이터에서는 좋지만 새 사용자에게 일반화되지 않음
  2. 다중공선성(multicollinearity): 시차 변수끼리, 파생 피처끼리 높은 상관 → 추정치 불안정
  3. 해석 불가: 1,475개 계수를 보고 “어떤 피처가 중요한가?”에 답할 수 없음
  4. 추정 불가: 순수 OLS는 \(p > n\)이면 유일해가 존재하지 않음

1.1.2 1.2 전통적 변수 선택의 한계

Stepwise selection (전진/후진/양방향)은 가장 널리 쓰이는 변수 선택 방법이었지만:

  • 탐욕적(greedy): 전역 최적이 아닌 지역 최적에 빠짐
  • 불안정: 데이터를 약간만 바꿔도 선택 결과가 크게 변함
  • p-value 남용: 다중 비교 문제를 무시
  • 고차원에서 불가: \(p > n\)이면 아예 시작할 수 없음

정규화(regularization)는 이 모든 문제를 목적 함수에 페널티 항을 추가하는 단일 프레임워크로 해결한다.

1.1.3 1.3 NYC Airbnb 예시에서의 고차원

NYC Airbnb 가격 분석에서도 비슷한 상황이 발생한다. 숙소별로 시간에 따라 가격, 리뷰 수, 점유율 등이 변하고, 여기에 동네(neighbourhood) 특성, 숙소 유형 더미, 시즌 효과, 교호작용까지 추가하면 후보 피처가 수백 개에 이른다. 동네(group) 수가 200여 개이고, 시간 포인트가 월별 12개라면 유효 군집 수 대비 피처 수가 과다해진다.


1.2 2. L1 / L2 정규화 수학

1.2.1 2.1 Ridge Regression (L2 정규화)

일반 선형 회귀의 목적 함수에 계수의 제곱합 페널티를 추가한다:

\[ \hat{\beta}^{\text{Ridge}} = \arg\min_{\beta} \left\{ \sum_{i=1}^{n} (y_i - X_i\beta)^2 + \lambda \sum_{j=1}^{p} \beta_j^2 \right\} \]

행렬 형태로 쓰면:

\[ \hat{\beta}^{\text{Ridge}} = \arg\min_{\beta} \left\{ \|y - X\beta\|_2^2 + \lambda \|\beta\|_2^2 \right\} \]

닫힌 해(closed-form solution)가 존재한다:

\[ \hat{\beta}^{\text{Ridge}} = (X^\top X + \lambda I_p)^{-1} X^\top y \]

핵심 특성:

  • \(\lambda \to 0\): OLS와 동일
  • \(\lambda \to \infty\): 모든 계수가 0에 수렴 (하지만 정확히 0이 되지는 않음)
  • 다중공선성이 있어도 \((X^\top X + \lambda I_p)\)는 항상 역행렬이 존재
  • 변수 선택은 하지 않음 — 모든 변수가 모델에 남아 있음

편향-분산 트레이드오프: Ridge는 OLS 대비 편향(bias)이 증가하지만 분산(variance)이 감소한다. 총 MSE = Bias\(^2\) + Variance이므로, 적절한 \(\lambda\)에서 MSE가 최소가 된다.

1.2.2 2.2 Lasso Regression (L1 정규화)

Tibshirani(1996)가 제안한 Lasso는 계수의 절대값 합 페널티를 사용한다:

\[ \hat{\beta}^{\text{Lasso}} = \arg\min_{\beta} \left\{ \frac{1}{2n}\|y - X\beta\|_2^2 + \lambda \sum_{j=1}^{p} |\beta_j| \right\} \]

또는 등가 제약 형태(Lagrangian duality):

\[ \hat{\beta}^{\text{Lasso}} = \arg\min_{\beta} \|y - X\beta\|_2^2 \quad \text{s.t.} \quad \sum_{j=1}^{p} |\beta_j| \leq t \]

핵심 특성:

  • 일부 계수를 정확히 0으로 만듦 → 자동 변수 선택
  • 닫힌 해 없음 → coordinate descent, LARS 등의 알고리즘 필요
  • \(p > n\)에서도 작동 (단, 최대 \(n\)개 변수만 선택 가능)

1.2.3 2.3 기하학적 직관: 왜 L1은 희소(sparse)한가?

이것이 정규화 이해의 핵심이다. 2차원(\(\beta_1, \beta_2\))으로 시각화하면:

  • OLS 등고선: 타원형 (\((y - X\beta)^\top(y - X\beta) = c\))
  • L2 제약 영역: 원 (\(\beta_1^2 + \beta_2^2 \leq t\))
  • L1 제약 영역: 마름모 (\(|\beta_1| + |\beta_2| \leq t\))

등고선이 제약 영역에 처음 접하는 점이 해이다.

  • L2 (원): 등고선이 원의 매끄러운 곡면에 접함 → 축 위에 떨어질 확률이 0 (measure zero)
  • L1 (마름모): 등고선이 마름모의 꼭짓점에 접할 확률이 높음 → 꼭짓점은 한 좌표가 0 → 희소

\(p\)차원에서 마름모(cross-polytope)의 꼭짓점, 모서리, 면은 대부분 좌표축에 정렬되어 있으므로, 고차원일수록 희소 해를 얻을 확률이 더 높아진다.

1.2.4 2.4 KKT 조건 직관

Lasso의 최적화 조건(Karush-Kuhn-Tucker)을 통해 왜 계수가 0이 되는지 정확히 이해할 수 있다.

\(\beta_j\)에 대한 subgradient 조건:

\[ -\frac{1}{n} X_j^\top (y - X\hat{\beta}) + \lambda \cdot s_j = 0 \]

여기서 \(s_j \in \partial |\hat{\beta}_j|\) (subgradient):

\[ s_j = \begin{cases} \text{sign}(\hat{\beta}_j) & \text{if } \hat{\beta}_j \neq 0 \\ \in [-1, 1] & \text{if } \hat{\beta}_j = 0 \end{cases} \]

\(\hat{\beta}_j = 0\)이 되는 조건:

\[ \left| \frac{1}{n} X_j^\top (y - X\hat{\beta}_{-j}) \right| \leq \lambda \]

즉, 변수 \(j\)와 현재 잔차의 상관(내적)이 \(\lambda\)보다 작으면 해당 변수는 탈락한다. \(\lambda\)를 키울수록 탈락 기준이 엄격해지므로 더 많은 변수가 0이 된다.

Soft-thresholding 연산자: Lasso의 coordinate descent 업데이트는 다음과 같다:

\[ \hat{\beta}_j \leftarrow S\left(\frac{1}{n}X_j^\top r_{-j},\; \lambda\right) = \text{sign}\left(\frac{1}{n}X_j^\top r_{-j}\right) \max\left(\left|\frac{1}{n}X_j^\top r_{-j}\right| - \lambda,\; 0\right) \]

여기서 \(r_{-j} = y - X_{-j}\hat{\beta}_{-j}\)는 변수 \(j\)를 제외한 잔차이다.

1.2.5 2.5 Elastic Net

Zou & Hastie(2005)가 제안한 Elastic Net은 L1과 L2를 혼합한다:

\[ \hat{\beta}^{\text{EN}} = \arg\min_{\beta} \left\{ \frac{1}{2n}\|y - X\beta\|_2^2 + \lambda \left[ \alpha \|\beta\|_1 + \frac{1-\alpha}{2} \|\beta\|_2^2 \right] \right\} \]

  • \(\alpha = 1\): 순수 Lasso
  • \(\alpha = 0\): 순수 Ridge
  • \(0 < \alpha < 1\): 변수 선택 + 안정성의 절충

왜 Elastic Net이 필요한가?

Lasso의 세 가지 약점:

  1. Grouping effect 부재: 상관이 높은 변수 그룹에서 임의로 하나만 선택하고 나머지를 버림. 예: response_timeresponse_time_lag1의 상관이 0.8이면, Lasso는 둘 중 하나만 살림.
  2. \(p > n\) 한계: 이론적으로 최대 \(n\)개의 변수만 선택 가능. \(p = 1{,}475\)이고 \(n_{\text{effective}} = 200\)이면 최대 200개.
  3. 불안정성: Bootstrap으로 반복하면, 상관된 변수 중 어떤 것이 선택될지 매번 바뀜.

Elastic Net의 L2 항은 이 문제를 해결한다:

\[ |\hat{\beta}_i^{\text{EN}} - \hat{\beta}_j^{\text{EN}}| \leq \frac{1}{\lambda(1-\alpha)} \|X_i - X_j\|_2 \cdot \sqrt{2(1-\rho_{ij})} \]

상관 \(\rho_{ij}\)가 1에 가까운 변수 쌍의 계수는 비슷한 크기가 되므로, 상관된 변수들이 함께 선택되거나 함께 제거된다 (grouping effect).


1.3 3. 종단 데이터에서의 특수 고려

1.3.1 3.1 왜 나이브 Lasso가 위험한가

종단 데이터 \((y_{ij}, X_{ij})\)를 단순히 풀어서(unstack) 독립 관측치로 취급하고 Lasso를 적용하면:

\[ \hat{\beta}^{\text{naive}} = \arg\min_{\beta} \sum_{i=1}^{N} \sum_{j=1}^{T_i} (y_{ij} - X_{ij}\beta)^2 + \lambda \|\beta\|_1 \]

문제 1 — 유효 표본 크기 과대추정:

동일 사용자의 반복 측정치는 독립이 아니다. ICC = 0.5일 때, 200명 \(\times\) 10세션 = 2,000개 관측치의 유효 표본 크기는:

\[ n_{\text{eff}} = \frac{N \cdot T}{1 + (T-1) \cdot \text{ICC}} = \frac{2000}{1 + 9 \times 0.5} = 364 \]

그런데 나이브 Lasso는 \(n = 2{,}000\)으로 취급하므로 \(\lambda\)가 너무 작아지고, 노이즈 변수까지 선택해 버린다.

문제 2 — 개체 이질성(heterogeneity) 혼동:

사용자마다 기저 만족도가 다른데(랜덤 절편), 이를 모델링하지 않으면 개체 효과와 상관된 공변량이 과대 추정된다. 예: “자주 사용하는 사용자”가 “만족도가 높은” 경향이 있으면, session_count 변수의 효과가 개체 효과와 혼동되어 부풀려진다.

문제 3 — 표준 오차 과소 추정:

실제로는 200명인데, 2,000개 독립 관측치로 취급하므로 표준 오차가 \(\sqrt{10}\)배 정도 과소추정된다.

1.3.2 3.2 종단 데이터에 맞는 정규화 전략

전략 개체 내 상관 변수 선택 구현 난이도 소프트웨어
나이브 Lasso (무시) X O 쉬움 sklearn, glmnet
군집 인식 CV + Lasso 부분적 O 쉬움 sklearn + GroupKFold
2단계: Lasso → LMM 부분적 O 중간 sklearn → statsmodels
GEE + L1 O (working corr) O 중간 penalized GEE
glmmLasso O (랜덤 효과) O 중간 glmmLasso (R)
Group Lasso + 군집 부분적 O (그룹 단위) 높음 gglasso, grpreg
Bayesian Lasso + RE O O 높음 brms, rstanarm

1.3.3 3.3 시간 변동 계수와 그룹 구조

종단 데이터에서는 피처가 자연스러운 그룹 구조를 가진다:

  • 시차 그룹: \(x_1\)의 현재값, 시차 1, 시차 2가 한 그룹
  • 카테고리 그룹: 더미 변수 세트 (직업군 5개 → 4개 더미)
  • 교호작용 그룹: 주효과 + 교호작용이 한 그룹 (hierarchy 원칙)

일반 Lasso는 개별 변수 단위로 선택하므로, “시차 2는 남기고 현재값은 제거”하는 비논리적 결과가 나올 수 있다. Group Lasso는 그룹 단위로 선택/제거한다:

\[ \hat{\beta}^{\text{GL}} = \arg\min_{\beta} \left\{ \|y - X\beta\|_2^2 + \lambda \sum_{g=1}^{G} \sqrt{p_g} \|\beta_g\|_2 \right\} \]

여기서 \(\beta_g\)는 그룹 \(g\)에 속한 계수 벡터, \(p_g\)는 그룹 크기, \(\sqrt{p_g}\)는 그룹 크기에 따른 보정 가중치이다.

Group Lasso의 페널티는 각 그룹의 L2 노름(유클리드 길이)에 작용하므로, 그룹 전체의 노름이 충분히 작으면 그룹 전체가 0이 된다. 그룹 내부에서는 L2 수축이 적용되어 개별 변수가 0이 되지는 않는다.


1.4 4. glmmLasso: 혼합 모델 + Lasso

1.4.1 4.1 모델 정의

Groll & Tutz(2014)가 제안한 glmmLasso는 GLMM에 L1 페널티를 결합한다. 연속형 결과 변수의 경우:

\[ y_{ij} = X_{ij}\beta + Z_{ij} b_i + \varepsilon_{ij} \]

여기서:

  • \(y_{ij}\): 사용자 \(i\)의 세션 \(j\) 만족도 점수
  • \(X_{ij}\): 고정 효과 설계 행렬 (고차원, \(p\)개 열)
  • \(\beta\): 고정 효과 계수 벡터 → L1 페널티 대상
  • \(Z_{ij}\): 랜덤 효과 설계 행렬 (보통 절편 + 시간)
  • \(b_i \sim N(0, D)\): 사용자별 랜덤 효과 → 페널티 없음
  • \(\varepsilon_{ij} \sim N(0, \sigma^2)\): 잔차

페널티 로그 우도(penalized log-likelihood):

\[ \ell_p(\beta, D, \sigma^2) = \ell(\beta, D, \sigma^2) - \lambda \sum_{j=1}^{p} |\beta_j| \]

여기서 로그 우도:

\[ \ell(\beta, D, \sigma^2) = \sum_{i=1}^{N} \log \int \prod_{j=1}^{T_i} f(y_{ij} | b_i; \beta, \sigma^2) \cdot f(b_i; D) \, db_i \]

핵심 설계 원리: 랜덤 효과 \(b_i\)는 페널티를 주지 않는다. 세 가지 이유:

  1. 자연 수축(natural shrinkage): 랜덤 효과는 이미 분포 가정 \(b_i \sim N(0, D)\)으로 0 방향으로 수축됨. 이것이 mixed model의 BLUP(Best Linear Unbiased Predictor)이 shrinkage estimator인 이유.
  2. 군집 구조 보존: 랜덤 효과를 0으로 만들면 개체 내 상관 구조가 깨지고, 나이브 Lasso와 다를 바 없어짐.
  3. 목적 분리: 변수 선택은 고정 효과(집단 수준 관계)에서 “어떤 공변량이 중요한가?”를 묻는 것이고, 랜덤 효과는 “개체마다 얼마나 다른가?”를 묻는 것. 두 질문을 분리해야 한다.

1.4.2 4.2 추정 알고리즘

glmmLasso는 다음과 같은 반복 알고리즘을 사용한다:

Step 1 — 고정 효과 업데이트: 현재 \(D\), \(\sigma^2\) 고정 → \(\beta\)에 대해 gradient descent + L1 페널티 (proximal gradient descent)

\[ \beta^{(k+1)} = \text{prox}_{\lambda \eta \|\cdot\|_1}\left(\beta^{(k)} - \eta \nabla_\beta \left[-\ell(\beta^{(k)})\right]\right) \]

여기서 proximal operator는 원소별로 적용되는 soft-thresholding:

\[ [\text{prox}_{\lambda \eta \|\cdot\|_1}(u)]_j = \text{sign}(u_j) \max(|u_j| - \lambda \eta, 0) \]

이것이 바로 soft-thresholding이다. 그래디언트로 한 걸음 이동한 뒤, \(\lambda \eta\)보다 작은 계수는 0으로 잘라낸다.

Step 2 — 분산 성분 업데이트: 현재 \(\beta\) 고정 → \(D\), \(\sigma^2\)를 EM 알고리즘 또는 Newton-Raphson으로 업데이트

Step 3: 수렴할 때까지 Step 1-2 반복. 수렴 기준: \(\|\beta^{(k+1)} - \beta^{(k)}\| < \epsilon\)

1.4.3 4.3 \(\lambda\) 선택: BIC 기반

교차 검증(CV)은 종단 데이터에서 주의가 필요하다. 단순 랜덤 분할은 같은 사용자의 관측치가 학습/검증에 모두 포함되어 정보 누수(leakage)가 발생한다. 군집 단위 CV(GroupKFold)가 올바르지만, 계산 비용이 높다.

glmmLasso는 BIC\(\lambda\)를 선택하는 것을 권장한다:

\[ \text{BIC}(\lambda) = -2\ell(\hat{\beta}_\lambda, \hat{D}_\lambda, \hat{\sigma}^2_\lambda) + \text{df}(\lambda) \cdot \log(N_{\text{total}}) \]

여기서:

  • \(\ell\): 최대화된 로그 우도
  • \(\text{df}(\lambda)\): 0이 아닌 고정 효과 계수의 수 + 분산 성분 수
  • \(N_{\text{total}} = \sum_{i=1}^N T_i\): 총 관측치 수

\(\lambda\) 그리드 \(\{\lambda_1 > \lambda_2 > \cdots > \lambda_K\}\)에 대해 각각 모델을 적합하고, BIC가 최소인 \(\lambda\)를 선택한다. 큰 \(\lambda\)에서 시작하여(모든 계수 0) 점점 줄여가면 계산 효율이 좋다(warm start).

1.4.4 4.4 GLMM으로의 확장

glmmLasso는 이진(binomial), 카운트(Poisson) 결과에도 적용된다:

\[ g(\mu_{ij}) = X_{ij}\beta + Z_{ij}b_i \]

\[ \ell_p(\beta, D) = \sum_{i=1}^{N} \log \int \prod_{j=1}^{T_i} f(y_{ij} | b_i; \beta) \cdot f(b_i; D) \, db_i - \lambda \|\beta\|_1 \]

  • 이진 결과: \(g = \text{logit}\), \(f = \text{Bernoulli}\)
  • 카운트 결과: \(g = \log\), \(f = \text{Poisson}\)

적분은 Laplace approximation 또는 adaptive Gauss-Hermite quadrature로 근사한다.


1.5 5. Python 실무 예시

1.5.1 5.1 시나리오 설정

AI Agent 개인화 실험: 200명의 사용자가 각 10세션 참여. 세션별 만족도(1-10)를 예측. 50개 원본 피처 + 교호작용/시차까지 약 150개 후보 변수. 진짜 영향을 주는 변수는 5개뿐.

1.5.2 5.2 시뮬레이션 데이터 생성

import numpy as np
import pandas as pd

np.random.seed(42)

n_users = 200
n_sessions = 10
n_features = 50

# --- 사용자별 랜덤 효과 ---
user_intercept = np.random.normal(0, 1.5, n_users)  # 기저 만족도 차이

# --- 세션 레벨 피처 생성 ---
rows = []
for i in range(n_users):
    for j in range(n_sessions):
        x = np.random.randn(n_features)
        # 진짜 효과: 변수 0, 3, 7, 12, 25만 영향
        y = (5.0 + user_intercept[i]
             + 0.8 * x[0]       # 응답 속도
             - 0.5 * x[3]       # 에러 횟수
             + 0.3 * x[7]       # 개인화 수준
             + 0.4 * x[12]      # 세션 길이
             - 0.2 * x[25]      # 광고 노출 수
             + 0.1 * j           # 시간 추세 (학습 효과)
             + np.random.normal(0, 0.8))
        row = {'user_id': i, 'session': j, 'satisfaction': y}
        for k in range(n_features):
            row[f'x{k}'] = x[k]
        rows.append(row)

df = pd.DataFrame(rows)
feature_cols = [f'x{k}' for k in range(n_features)]

print(f"데이터 크기: {df.shape}")
print(f"사용자 수: {df['user_id'].nunique()}")
print(f"피처 수: {len(feature_cols)}")
print(f"진짜 영향 변수: x0(+0.8), x3(-0.5), x7(+0.3), x12(+0.4), x25(-0.2)")

1.5.3 5.3 교호작용 및 시차 변수 추가

# --- 시차 1 변수 추가 (이전 세션 피처) ---
df_sorted = df.sort_values(['user_id', 'session'])
for col in feature_cols[:10]:  # 처음 10개만 시차
    df_sorted[f'{col}_lag1'] = df_sorted.groupby('user_id')[col].shift(1)

# 첫 세션은 시차 없음 → 제거
df_clean = df_sorted.dropna().reset_index(drop=True)

# 주요 교호작용 추가
df_clean['x0_x3'] = df_clean['x0'] * df_clean['x3']
df_clean['x0_x7'] = df_clean['x0'] * df_clean['x7']
df_clean['x3_x12'] = df_clean['x3'] * df_clean['x12']
df_clean['x7_x25'] = df_clean['x7'] * df_clean['x25']

# 전체 후보 피처 목록
lag_cols = [c for c in df_clean.columns if '_lag1' in c]
inter_cols = ['x0_x3', 'x0_x7', 'x3_x12', 'x7_x25']
all_features = feature_cols + lag_cols + inter_cols + ['session']

print(f"전체 후보 피처: {len(all_features)}개")
print(f"  원본: {len(feature_cols)}, 시차: {len(lag_cols)}, "
      f"교호작용: {len(inter_cols)}, 시간: 1")

1.5.4 5.4 나이브 Lasso (군집 무시)

from sklearn.linear_model import LassoCV, ElasticNetCV
from sklearn.preprocessing import StandardScaler

X = df_clean[all_features].values
y = df_clean['satisfaction'].values

# 표준화 (정규화 전 필수)
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# --- Lasso with 5-fold CV (군집 무시) ---
lasso_cv = LassoCV(
    alphas=np.logspace(-4, 1, 100),
    cv=5,              # 주의: 같은 사용자가 train/test에 섞임
    max_iter=10000,
    random_state=42
)
lasso_cv.fit(X_scaled, y)

print(f"최적 alpha: {lasso_cv.alpha_:.4f}")
print(f"0이 아닌 계수: {np.sum(lasso_cv.coef_ != 0)}개 / {len(all_features)}개")

# 선택된 변수 확인
selected = [(name, coef) for name, coef in zip(all_features, lasso_cv.coef_) if coef != 0]
selected_sorted = sorted(selected, key=lambda x: abs(x[1]), reverse=True)

print("\n--- 선택된 변수 (|계수| 내림차순) ---")
for name, coef in selected_sorted[:15]:
    marker = " *** TRUE" if name in ['x0','x3','x7','x12','x25','session'] else ""
    print(f"  {name:20s}  {coef:+.4f}{marker}")

1.5.5 5.5 Elastic Net

# --- Elastic Net with CV ---
enet_cv = ElasticNetCV(
    l1_ratio=[0.1, 0.3, 0.5, 0.7, 0.9, 0.95, 0.99],
    alphas=np.logspace(-4, 1, 50),
    cv=5,
    max_iter=10000,
    random_state=42
)
enet_cv.fit(X_scaled, y)

print(f"최적 alpha: {enet_cv.alpha_:.4f}")
print(f"최적 l1_ratio: {enet_cv.l1_ratio_:.2f}")
print(f"0이 아닌 계수: {np.sum(enet_cv.coef_ != 0)}개")

# Lasso vs Elastic Net 비교
n_lasso = np.sum(lasso_cv.coef_ != 0)
n_enet = np.sum(enet_cv.coef_ != 0)
print(f"\nLasso 선택: {n_lasso}개")
print(f"Elastic Net 선택: {n_enet}개")
print(f"→ Elastic Net이 {'더 많은' if n_enet > n_lasso else '더 적은'} 변수를 선택")

1.5.6 5.6 Coefficient Path 시각화

import matplotlib.pyplot as plt
from sklearn.linear_model import lasso_path

# Lasso path 계산
alphas_path, coefs_path, _ = lasso_path(
    X_scaled, y,
    alphas=np.logspace(-4, 1, 200),
    max_iter=10000
)

# --- Coefficient Path Plot ---
fig, ax = plt.subplots(figsize=(12, 6))

# 진짜 변수를 강조 (색상 + 굵은 선)
true_vars = {
    'x0': all_features.index('x0'),
    'x3': all_features.index('x3'),
    'x7': all_features.index('x7'),
    'x12': all_features.index('x12'),
    'x25': all_features.index('x25'),
    'session': all_features.index('session')
}

# 나머지 변수는 회색으로
for idx in range(coefs_path.shape[0]):
    if idx not in true_vars.values():
        ax.plot(np.log10(alphas_path), coefs_path[idx],
                color='grey', alpha=0.15, linewidth=0.5)

# 진짜 변수는 색상으로
colors = ['#e41a1c', '#377eb8', '#4daf4a', '#984ea3', '#ff7f00', '#a65628']
for (name, idx), color in zip(true_vars.items(), colors):
    ax.plot(np.log10(alphas_path), coefs_path[idx],
            linewidth=2.5, label=f'{name} (TRUE)', color=color)

ax.axhline(y=0, color='black', linestyle='--', linewidth=0.5)
ax.axvline(x=np.log10(lasso_cv.alpha_), color='red', linestyle=':',
           linewidth=1.5, label=f'CV optimal (alpha={lasso_cv.alpha_:.4f})')
ax.set_xlabel('log10(alpha)', fontsize=12)
ax.set_ylabel('Coefficient value', fontsize=12)
ax.set_title('Lasso Coefficient Path — AI Agent Personalization (50+ features)',
             fontsize=14)
ax.legend(loc='upper right', fontsize=9, framealpha=0.9)
plt.tight_layout()
plt.savefig('lasso_path_agent.png', dpi=150, bbox_inches='tight')
plt.show()

1.5.7 5.7 GroupKFold로 군집 인식 교차 검증

나이브 5-fold CV 대신 사용자 단위로 fold를 나누어 정보 누수를 방지한다.

from sklearn.model_selection import GroupKFold
from sklearn.linear_model import Lasso, ElasticNet
from sklearn.metrics import mean_squared_error

groups = df_clean['user_id'].values
gkf = GroupKFold(n_splits=5)

# --- 군집 인식 CV로 최적 alpha 탐색 ---
alphas = np.logspace(-4, 1, 50)
cv_scores_naive = []   # 나이브 CV
cv_scores_group = []   # 군집 인식 CV

for alpha in alphas:
    # 군집 인식 CV
    fold_scores = []
    for train_idx, test_idx in gkf.split(X_scaled, y, groups):
        model = Lasso(alpha=alpha, max_iter=10000)
        model.fit(X_scaled[train_idx], y[train_idx])
        y_pred = model.predict(X_scaled[test_idx])
        fold_scores.append(mean_squared_error(y[test_idx], y_pred))
    cv_scores_group.append(np.mean(fold_scores))

best_alpha_group = alphas[np.argmin(cv_scores_group)]
print(f"군집 인식 CV 최적 alpha: {best_alpha_group:.4f}")
print(f"나이브 CV 최적 alpha:   {lasso_cv.alpha_:.4f}")
print(f"→ 군집 인식 alpha가 {'더 큼' if best_alpha_group > lasso_cv.alpha_ else '더 작음'} "
      f"(더 강한 정규화)")

# 군집 인식 alpha로 재적합
lasso_group = Lasso(alpha=best_alpha_group, max_iter=10000)
lasso_group.fit(X_scaled, y)
n_selected_group = np.sum(lasso_group.coef_ != 0)
print(f"\n군집 인식 선택 변수: {n_selected_group}개")
print(f"나이브 선택 변수:   {np.sum(lasso_cv.coef_ != 0)}개")

1.5.8 5.8 2단계 접근: Lasso → statsmodels MixedLM

Lasso로 변수를 선택하고, 그 결과를 Mixed Model에 넣는 실용적 접근법.

import statsmodels.formula.api as smf

# Step 1: Lasso로 변수 선택
selected_vars = [name for name, coef in zip(all_features, lasso_group.coef_)
                 if coef != 0]
print(f"Step 1 — Lasso 선택 변수 ({len(selected_vars)}개): {selected_vars}")

# Step 2: 선택된 변수로 MixedLM 적합
formula = f"satisfaction ~ {' + '.join(selected_vars)}"
mixed_model = smf.mixedlm(
    formula,
    data=df_clean,
    groups=df_clean['user_id'],
    re_formula='~session'   # 랜덤 절편 + 랜덤 기울기(시간)
)
result = mixed_model.fit(reml=True)
print("\n--- Mixed Model with Lasso-selected variables ---")
print(result.summary())

# 랜덤 효과 분산 확인
print(f"\n랜덤 절편 SD: {np.sqrt(result.cov_re.iloc[0, 0]):.3f}")
print(f"잔차 SD: {np.sqrt(result.scale):.3f}")

# --- 비교: 전체 변수 MixedLM ---
try:
    formula_full = f"satisfaction ~ {' + '.join(all_features)}"
    mixed_full = smf.mixedlm(
        formula_full,
        data=df_clean,
        groups=df_clean['user_id'],
        re_formula='~1'   # 랜덤 절편만 (수렴 위해 단순화)
    )
    result_full = mixed_full.fit(reml=True)
    print(f"\n전체 변수 모델 AIC: {result_full.aic:.1f}")
    print(f"선택 변수 모델 AIC: {result.aic:.1f}")
    print(f"→ AIC 차이: {result_full.aic - result.aic:.1f}")
except Exception as e:
    print(f"전체 변수 모델 실패: {e}")
    print("→ 이것이 변수 선택이 필요한 이유")

1.6 6. R 실무 예시

1.6.1 6.1 데이터 생성

library(tidyverse)
set.seed(42)

n_users <- 200
n_sessions <- 10
n_features <- 50

# --- 사용자별 랜덤 효과 ---
user_intercept <- rnorm(n_users, 0, 1.5)

# --- 시뮬레이션 데이터 ---
df <- expand_grid(user_id = 1:n_users, session = 0:(n_sessions - 1))

# 50개 피처 생성
X_mat <- matrix(rnorm(nrow(df) * n_features), ncol = n_features)
colnames(X_mat) <- paste0("x", 0:(n_features - 1))
df <- bind_cols(df, as_tibble(X_mat))

# 진짜 효과 부여
df <- df %>%
  mutate(
    satisfaction = 5.0 + user_intercept[user_id] +
      0.8 * x0 - 0.5 * x3 + 0.3 * x7 + 0.4 * x12 - 0.2 * x25 +
      0.1 * session +
      rnorm(n(), 0, 0.8)
  )

cat("데이터 크기:", nrow(df), "x", ncol(df), "\n")
cat("진짜 영향 변수: x0(+0.8), x3(-0.5), x7(+0.3), x12(+0.4), x25(-0.2)\n")

1.6.2 6.2 glmnet: 기본 Lasso / Elastic Net

library(glmnet)

feature_names <- paste0("x", 0:(n_features - 1))
X <- as.matrix(df[, feature_names])
y <- df$satisfaction

# --- 군집 단위 fold 생성 ---
# 같은 사용자가 같은 fold에 들어가야 함
user_fold <- rep(1:5, length.out = n_users)  # 사용자를 5개 fold에 배정
foldid <- user_fold[df$user_id]               # 관측치별 fold ID

# --- Lasso (alpha=1) with 군집 인식 CV ---
cv_lasso <- cv.glmnet(X, y, alpha = 1, foldid = foldid)
plot(cv_lasso, main = "Lasso CV — 군집 인식 fold")

cat("최적 lambda (min):", cv_lasso$lambda.min, "\n")
cat("최적 lambda (1se):", cv_lasso$lambda.1se, "\n")

# 선택된 변수 (lambda.1se = 더 보수적, 간결한 모델)
coef_lasso <- coef(cv_lasso, s = "lambda.1se")
selected <- rownames(coef_lasso)[which(coef_lasso != 0)]
selected <- selected[selected != "(Intercept)"]
cat("\n선택된 변수 (", length(selected), "개):\n")
for (v in selected) {
  cat("  ", v, ":", round(coef_lasso[v, ], 4), "\n")
}

# --- Elastic Net (alpha=0.5) ---
cv_enet <- cv.glmnet(X, y, alpha = 0.5, foldid = foldid)
coef_enet <- coef(cv_enet, s = "lambda.1se")
selected_enet <- rownames(coef_enet)[which(coef_enet != 0)]
selected_enet <- selected_enet[selected_enet != "(Intercept)"]
cat("\nElastic Net 선택 변수 (", length(selected_enet), "개)\n")

1.6.3 6.3 Coefficient Path 시각화 (R)

library(ggplot2)

# glmnet path
fit_lasso <- glmnet(X, y, alpha = 1)

# ggplot2로 coefficient path
coef_mat <- as.matrix(fit_lasso$beta)
df_path <- as.data.frame(t(coef_mat)) %>%
  mutate(lambda = fit_lasso$lambda) %>%
  pivot_longer(-lambda, names_to = "variable", values_to = "coefficient")

true_vars <- c("x0", "x3", "x7", "x12", "x25")
df_path <- df_path %>%
  mutate(is_true = variable %in% true_vars)

ggplot(df_path, aes(x = log10(lambda), y = coefficient, group = variable)) +
  geom_line(data = filter(df_path, !is_true),
            color = "grey80", alpha = 0.5, linewidth = 0.3) +
  geom_line(data = filter(df_path, is_true),
            aes(color = variable), linewidth = 1.2) +
  geom_vline(xintercept = log10(cv_lasso$lambda.1se),
             linetype = "dashed", color = "red") +
  annotate("text", x = log10(cv_lasso$lambda.1se) + 0.1, y = max(coef_mat) * 0.8,
           label = "lambda.1se", color = "red", hjust = 0) +
  scale_color_brewer(palette = "Set1") +
  labs(title = "Lasso Coefficient Path — AI Agent 개인화 실험",
       subtitle = "색상 = 진짜 영향 변수, 회색 = 노이즈 변수",
       x = "log10(lambda)", y = "Coefficient") +
  theme_minimal(base_size = 12)

1.6.4 6.4 glmmLasso: 혼합 모델 + L1 페널티

library(glmmLasso)

# glmmLasso는 formula 기반 → 피처 수를 관리 가능한 수준으로
# 방법 1: 전체 50개를 넣고 glmmLasso가 선택하게 함
# 방법 2: glmnet에서 사전 선별(30개) 후 glmmLasso에 넣음 (더 빠름)

# 여기서는 방법 2를 사용 (실무에서 권장)
coef_screen <- coef(cv_lasso, s = "lambda.min")  # 관대한 lambda로 사전 선별
nonzero_idx <- which(abs(coef_screen[-1, ]) > 0)
screened_vars <- feature_names[nonzero_idx]
cat("사전 선별된 변수:", length(screened_vars), "개\n")

# 선별된 변수가 30개 이하가 되도록 조정
if (length(screened_vars) > 30) {
  top_idx <- order(abs(coef_screen[-1, ]), decreasing = TRUE)[1:30]
  screened_vars <- feature_names[top_idx]
}

# glmmLasso formula 구성
fixed_formula <- as.formula(
  paste("satisfaction ~", paste(screened_vars, collapse = " + "))
)
cat("Formula:", deparse(fixed_formula), "\n")

# --- Lambda 그리드 탐색 (BIC 기준) ---
lambda_grid <- seq(5, 500, by = 5)
bic_results <- data.frame(lambda = lambda_grid, bic = NA, n_nonzero = NA)

for (i in seq_along(lambda_grid)) {
  tryCatch({
    fit <- glmmLasso(
      fix = fixed_formula,
      rnd = list(user_id = ~ 1),  # 랜덤 절편
      data = df,
      lambda = lambda_grid[i],
      family = gaussian(link = "identity"),
      control = list(
        start = NULL,
        center = TRUE,
        standardize = TRUE
      )
    )
    bic_results$bic[i] <- fit$bic
    bic_results$n_nonzero[i] <- sum(fit$coefficients[-1] != 0)
  }, error = function(e) {
    bic_results$bic[i] <<- NA
    bic_results$n_nonzero[i] <<- NA
  })
}

# 최적 lambda
best_row <- bic_results[which.min(bic_results$bic), ]
cat("\n최적 lambda (BIC):", best_row$lambda, "\n")
cat("최적 모델의 0이 아닌 계수:", best_row$n_nonzero, "개\n")

# BIC 곡선 시각화
bic_clean <- na.omit(bic_results)
ggplot(bic_clean, aes(x = lambda, y = bic)) +
  geom_line(linewidth = 0.8) +
  geom_point(data = filter(bic_clean, lambda == best_row$lambda),
             color = "red", size = 3) +
  geom_text(data = filter(bic_clean, lambda == best_row$lambda),
            aes(label = paste0("lambda=", lambda)),
            vjust = -1.5, color = "red") +
  labs(title = "glmmLasso: BIC vs Lambda",
       subtitle = "AI Agent 개인화 실험 (랜덤 절편 모델)",
       x = "Lambda (페널티 강도)", y = "BIC") +
  theme_minimal(base_size = 12)

1.6.5 6.5 최적 모델 적합 및 결과 해석

# --- 최적 lambda로 최종 모델 ---
final_fit <- glmmLasso(
  fix = fixed_formula,
  rnd = list(user_id = ~ 1),
  data = df,
  lambda = best_row$lambda,
  family = gaussian(link = "identity"),
  control = list(center = TRUE, standardize = TRUE)
)

# --- 고정 효과 분석 ---
coef_fixed <- final_fit$coefficients
nonzero_mask <- abs(coef_fixed) > 1e-8
nonzero_coefs <- coef_fixed[nonzero_mask]
nonzero_coefs <- sort(nonzero_coefs[names(nonzero_coefs) != "(Intercept)"],
                      decreasing = TRUE)

cat("\n===== glmmLasso 최종 결과 =====\n")
cat("\n--- 선택된 변수 (계수 크기 내림차순) ---\n")
for (v in names(nonzero_coefs)) {
  true_marker <- ifelse(v %in% c("x0","x3","x7","x12","x25"), " *** TRUE", "")
  cat(sprintf("  %-15s %+.4f%s\n", v, nonzero_coefs[v], true_marker))
}

cat("\n--- 제거된 변수 (계수 = 0) ---\n")
zero_vars <- names(coef_fixed[!nonzero_mask])
cat("  총", length(zero_vars), "개:", paste(head(zero_vars, 10), collapse=", "), "...\n")

# 랜덤 효과 분산
cat("\n--- 랜덤 효과 ---\n")
cat("  랜덤 절편 SD:", round(sqrt(final_fit$StdDev), 3), "\n")

# --- lme4로 검증 ---
library(lme4)
survived_vars <- names(coef_fixed[nonzero_mask])
survived_vars <- survived_vars[survived_vars != "(Intercept)"]

lmm_formula <- as.formula(
  paste("satisfaction ~", paste(survived_vars, collapse = " + "),
        "+ (1 | user_id)")
)
lmm_fit <- lmer(lmm_formula, data = df, REML = TRUE)

cat("\n===== lme4 검증 (glmmLasso 선택 변수) =====\n")
print(summary(lmm_fit))

cat("\n--- 모델 비교 ---\n")
cat("glmmLasso BIC:", final_fit$bic, "\n")
cat("lme4 BIC:     ", BIC(lmm_fit), "\n")

1.7 7. Elastic Net vs Group Lasso vs glmmLasso 비교

1.7.1 7.1 종합 비교표

기준 Lasso (L1) Elastic Net (L1+L2) Group Lasso glmmLasso
변수 선택 단위 개별 변수 개별 변수 그룹 단위 개별 변수
상관 변수 처리 하나만 선택 그룹 유지 (grouping effect) 그룹 유지 개별 선택
개체 내 상관 처리 X X X (별도 처리 필요) O (랜덤 효과)
\(p > n\) 지원 O (최대 \(n\)개) O (제한 없음) O O (느릴 수 있음)
계산 속도 매우 빠름 빠름 중간 느림
GLMM 확장 X X X O (binomial, poisson)
Lambda 선택 CV (GroupKFold 권장) CV CV BIC (권장)
이론적 정당성 Oracle property (irrepresentable cond.) Oracle property 개선 그룹 Oracle 일관성 증명 있음
R 패키지 glmnet glmnet gglasso, grpreg glmmLasso
Python 패키지 sklearn sklearn group-lasso 없음 (2단계 사용)

1.7.2 7.2 상세 사용 시나리오

Lasso를 쓸 때:

  • 변수 간 상관이 낮고, 진짜 변수가 소수일 때
  • 빠른 탐색적 분석이 목적일 때
  • \(p\)가 중간 크기(\(< 1{,}000\))이고 군집 구조가 약할 때

Elastic Net을 쓸 때:

  • 상관된 변수 그룹이 있고, 그룹 전체가 중요할 때 (예: 유전체 데이터)
  • \(p \gg n\)이고 Lasso의 \(n\)개 제한이 문제될 때
  • Lasso 결과가 불안정할 때 (Bootstrap에서 선택 변수가 계속 바뀜)

Group Lasso를 쓸 때:

  • 자연적 그룹 구조: 카테고리 더미 세트, 시차 변수 세트, 다항 기저
  • “주효과가 없으면 교호작용도 없다” (hierarchy) 원칙을 강제하고 싶을 때
  • 그룹 내부 변수를 개별 선택할 필요가 없을 때

glmmLasso를 쓸 때:

  • 종단/반복 측정 데이터에서 변수 선택과 랜덤 효과를 동시에 필요로 할 때
  • 개체 수 대비 피처 수가 많지만, \(p\)가 수백 이하일 때
  • 최종 보고용 모델에서 해석 가능한 변수 집합이 필요할 때
  • GLMM으로 확장해야 할 때 (이진, 카운트 결과)

2단계 접근법 (Lasso/Elastic Net → LMM/GLMM):

  • 피처가 매우 많아서 (\(p > 500\)) glmmLasso가 비현실적으로 느릴 때
  • Python 환경에서 작업할 때 (glmmLasso는 R 전용)
  • 빠르게 후보 변수를 줄이고, 정밀 모델은 lme4나 statsmodels로 적합할 때
  • GroupKFold로 군집 인식 CV를 수행하면 합리적 결과를 얻을 수 있음

1.8 8. 실무 팁

1.8.1 8.1 표준화가 필수인 이유

Lasso의 페널티 \(\lambda |\beta_j|\)계수의 크기에 작용한다. 변수의 스케일이 다르면:

  • 단위가 큰 변수 (예: 연봉 4,000만 원) → 계수가 작아서 페널티를 적게 받음
  • 단위가 작은 변수 (예: 클릭률 0.03) → 계수가 커서 페널티를 많이 받음

결과적으로 스케일에 따라 선택/제거가 결정되는 부당한 결과가 나온다.

# Python: 반드시 표준화 후 적합
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)  # 평균=0, 분산=1

# 주의: 표준화된 계수를 원래 스케일로 변환하려면
# coef_original = coef_scaled / scaler.scale_
# R: glmnet은 내부적으로 표준화 (standardize=TRUE가 기본값)
# glmmLasso는 control에서 standardize = TRUE 설정
# 직접 표준화하려면:
X_scaled <- scale(X)  # center=TRUE, scale=TRUE

1.8.2 8.2 Lambda 그리드 전략

# --- 나쁜 예: 임의의 값 ---
alphas = [0.001, 0.01, 0.1, 1, 10]  # 5개뿐, 최적 lambda를 놓칠 가능성 높음

# --- 좋은 예: 로그 등간격 ---
alphas = np.logspace(-4, 1, 100)  # 0.0001 ~ 10, 100개

# --- 더 좋은 예: lambda_max에서 시작 ---
# lambda_max: 모든 계수가 0이 되는 최소 lambda
lambda_max = np.max(np.abs(X_scaled.T @ y)) / (2 * len(y))
alphas = np.logspace(np.log10(lambda_max), np.log10(lambda_max * 1e-4), 100)
# → lambda_max에서 시작하여 0.01%까지 100단계로 줄여감
# R: glmnet은 자동으로 적절한 lambda 그리드를 생성
# 수동 지정도 가능:
lambda_grid <- 10^seq(log10(max_lambda), log10(min_lambda), length.out = 100)

# glmmLasso는 수동 지정 필수:
lambda_grid <- seq(5, 500, by = 5)  # 넓은 범위에서 시작
# BIC 최소 근처에서 더 세밀하게 탐색
lambda_fine <- seq(best_lambda - 20, best_lambda + 20, by = 1)

1.8.3 8.3 수렴 문제 대처

glmmLasso에서 흔한 문제와 해결법:

문제 원인 해결법
수렴 안 됨 (max iterations) lambda가 너무 작음, 또는 피처 간 높은 상관 lambda 키우기; 상관 > 0.9인 변수 사전 제거
모든 계수가 0 lambda가 너무 큼 lambda 줄이기; lambda_max 이하에서 시작
느린 수렴 (100+ iterations) 피처 수 과다 2단계: glmnet 선별 → glmmLasso
메모리 부족 군집 수 \(\times\) 피처 수가 큼 랜덤 효과 구조 단순화 (random slope 제거)
랜덤 효과 분산 = 0 (boundary) 군집 수 부족 또는 ICC가 매우 낮음 랜덤 효과 필요성 재검토; LRT로 확인
BIC가 단조 감소 lambda 범위가 너무 좁음 lambda 범위를 넓혀서 U자형이 보일 때까지

1.8.4 8.4 결과 보고 방법

학술 논문이나 사내 보고서에서 정규화 결과를 보고할 때 포함해야 할 내용:

1. 방법론 기술:

“고차원 공변량에서 핵심 예측 변수를 선별하기 위해 glmmLasso를 적용하였다. 사용자별 랜덤 절편을 포함하는 선형 혼합 모델에 고정 효과의 L1 페널티를 결합하였으며, BIC 기준으로 \(\lambda = 85\)를 선택하였다.”

2. 선택된 변수 테이블:

변수 glmmLasso 계수 lme4 계수 (SE) p-value 해석
x0 (응답 속도) 0.78 0.81 (0.04) < 0.001 빠를수록 만족도 높음
x3 (에러 횟수) -0.48 -0.52 (0.05) < 0.001 에러 많을수록 불만족
session 0.09 0.10 (0.02) < 0.001 학습 효과

3. Coefficient Path 그래프: Lambda에 따라 변수가 진입/퇴출하는 과정을 시각화

4. 민감도 분석: lambda.min vs lambda.1se 결과 비교, 또는 Elastic Net 대비

5. Post-selection inference 주의 명시:

“주의: glmmLasso로 선택된 변수에 대한 추론(p-value, 신뢰구간)은 post-selection inference 문제가 있으므로, 선택된 변수로 별도의 LMM(lme4)을 적합하여 추론을 수행하였다.”

1.8.5 8.5 Post-Selection Inference 주의

Lasso로 변수를 선택한 뒤, 같은 데이터로 LMM을 적합하면 p-value와 신뢰구간이 낙관적으로 편향된다. 이는 “선택 과정에서 이미 데이터를 사용”했기 때문이다 (double dipping).

대안:

  • Sample splitting: 데이터를 반으로 나눠서 선택(전반)/추론(후반) 분리
  • Selective inference: Lee et al.(2016)의 조건부 추론. selectiveInference R 패키지
  • Stability selection: Meinshausen & Buhlmann(2010). 여러 subsample에서 반복 선택하여 안정적인 변수만 보고. 거짓 양성 수의 상한을 이론적으로 보장
# --- Stability Selection 예시 ---
library(stabs)

stab_fit <- stabsel(
  x = X, y = y,
  fitfun = glmnet.lasso,
  cutoff = 0.75,       # 선택 확률 75% 이상인 변수만
  PFER = 1             # 기대 거짓 양성 수(Per-Family Error Rate) <= 1
)

cat("안정적으로 선택된 변수:\n")
print(names(stab_fit$selected))
plot(stab_fit, main = "Stability Selection — 변수별 선택 확률")

# 선택 확률이 높은 변수일수록 "진짜"일 가능성이 높음
# 보통 cutoff=0.6~0.9, PFER=1~2로 설정

1.9 9. 핵심 요약

1.9.1 9.1 한 줄 정리

기법 핵심 아이디어 페널티 수식 변수 선택
Ridge (L2) 계수 수축, 분산 감소 \(\lambda \|\beta\|_2^2\) X
Lasso (L1) 희소 해, 자동 선택 \(\lambda \|\beta\|_1\) O
Elastic Net L1+L2 혼합, grouping \(\lambda[\alpha\|\beta\|_1 + (1-\alpha)\|\beta\|_2^2/2]\) O
Group Lasso 그룹 단위 선택 \(\lambda \sum_g \sqrt{p_g}\|\beta_g\|_2\) O (그룹)
glmmLasso LMM + L1, RE 보존 \(-\ell(y|\beta, b) + \lambda\|\beta\|_1\) O

1.9.2 9.2 핵심 체크리스트

  • L1 페널티는 마름모 제약의 꼭짓점에서 해를 찾아 희소성(sparsity)을 유도한다
  • KKT 조건: 변수와 잔차의 상관이 \(\lambda\)보다 작으면 해당 변수의 계수는 0
  • 종단 데이터에 나이브 Lasso를 쓰면 유효 표본 크기를 과대추정하고, 개체 효과와 혼동된 변수를 잘못 선택할 수 있다
  • glmmLasso는 고정 효과에만 L1 페널티를 주고, 랜덤 효과는 분포 가정 \(N(0, D)\)으로 자연 수축시킨다
  • Lambda 선택: 종단 데이터에서는 BIC가 CV보다 안전하다 (정보 누수 방지). CV를 쓸 때는 반드시 GroupKFold
  • 표준화 없이 정규화를 쓰면 스케일이 큰 변수가 부당하게 살아남는다
  • 변수 선택 후 추론(p-value)을 할 때는 post-selection inference 문제를 인식하고, sample splitting 또는 stability selection을 고려한다
  • 피처가 매우 많으면 2단계 접근 (glmnet/sklearn 선별 → glmmLasso/lme4 정밀 모델)이 현실적이다
  • Elastic Net은 상관 변수의 grouping effect를 제공하여 Lasso의 불안정성을 보완한다

1.9.3 9.3 관련 파일

Subscribe

Enjoy this blog? Get notified of new posts by email: