능형회귀

choyunjeong·2025년 1월 8일

회귀분석 제 3판 - 박성현

목록 보기
6/8

17.2 능형회귀 (Ridge regression)

중선형회귀모형에서 행렬 XX와 yy가 표준화되었다고 하고 XX의 열의 수를 pp라고 한다면 X′XX'X와 X′yX'y의 각각 상관계수행렬이 될 것이다.

X′X=[x11⋯x1nx21⋯x2n]⋅[x11x21⋮⋮x1nx2n]=[∑x1i2∑x1ix2i∑x2ix1i∑x2i2]=[1r12r121]X′y=[x11⋯x1nx21⋯x2n]⋅[y1⋮yn]=[∑x1iyi∑x2iyi]=[r1yr2y]\begin{aligned} X'X &= \begin{bmatrix} x_{11} & \cdots & x_{1n} \\[10pt] x_{21} & \cdots & x_{2n} \end{bmatrix} \cdot \begin{bmatrix} x_{11} & x_{21} \\[10pt] \vdots & \vdots \\[10pt] x_{1n} & x_{2n} \end{bmatrix} =\begin{bmatrix} \sum x_{1i}^2 & \sum x_{1i}x_{2i} \\[10pt] \sum x_{2i}x_{1i} & \sum x_{2i}^2 \end{bmatrix} =\begin{bmatrix} 1 & r_{12} \\[10pt] r_{12} & 1 \end{bmatrix} \\[30pt] X'y &= \begin{bmatrix} x_{11} & \cdots & x_{1n} \\[10pt] x_{21} & \cdots & x_{2n} \end{bmatrix} \cdot \begin{bmatrix} y_{1} \\[10pt] \vdots \\[10pt] y_{n} \end{bmatrix} =\begin{bmatrix} \sum x_{1i}y_i \\[10pt] \sum x_{2i}y_i \end{bmatrix} =\begin{bmatrix} r_{1y} \\[10pt] r_{2y} \end{bmatrix} \end{aligned}

이 때 최소제곱추정량 b=(X′X)−1X′yb=(X'X)^{-1}X'y와 β\beta간의 거리의 제곱을

L=(b−β)′(b−β)L=(b-\beta)'(b-\beta)

으로 정의하고 λi\lambda_i를 행렬 X′XX'X의 고유값이라 한다면 L2L^2의 기대값은 (총)편균제곱오차 (MSE(b))(MSE(b))라고도 부르며 다음과 같다.

E(L2)=MSE(b)=E[(b−β)′(b−β)]=tr{E[(b−β)′(b−β)]}(∵tr(a)=a; scalar)=E{tr[(b−β)(b−β)′]}(p×p ∵ outer product)=∑i=1pVar(bi)=σ2∑i=1pλi−1((Var(bi)=σ2(X′X)−1=σ2QΛ−1Q′=σ2Λ−1)\begin{aligned} E(L^2) &= MSE(b) \\[10pt] &=E[(b-\beta)'(b-\beta)] \\[10pt] &= \text{tr}\{E[(b-\beta)'(b-\beta)]\}\quad(\because tr(a)=a;\ scalar) \\[10pt] &= E\{\text{tr}[(b-\beta)(b-\beta)']\}\quad(p\times p\ \because \text{ outer product})\\[10pt] &= \sum_{i=1}^{p}\text{Var}(b_i) \\[15pt] &= \sigma^2\sum_{i=1}^{p}\lambda_i^{-1} \quad ((\text{Var}(b_i)=\sigma^2(X'X)^{-1}=\sigma^2Q\Lambda^{-1}Q'=\sigma^2\Lambda^{-1}) \\[10pt] \end{aligned}

이 때 Λ=diag(λ1,…,λp))\Lambda=\text{diag}(\lambda_1,\ldots,\lambda_p)) QQ 는 고유벡터이다. 따라서 독립변수들간에 완전에 가까운 다중공선성이 존재하면 λi\lambda_i중에서 거의 0에 가까운 값이 있으며 bb는 β\beta로부터 대단히 멀리 떨어져 있을 가능성이 크며 bb는 β\beta의 좋은 추정량이라고 볼 수 없다.

\\[40pt]

17.2.0 능형 회귀 (Ridge regression)

위의 단점을 보완하기 위하여 다음과 같은 추정을 고안하였으며 이를 능형추정량 (ridge estimator)이라고 한다.

b(k)=(X′X+kI)−1X′y,k>0b(k)=(X'X+kI)^{-1}X'y,\quad k>0

여기서 kk는 양의 상수로서 보통 0<k<10<k<1범위 내에 있다. b(k)b(k)의 기댓값과 분산은 다음과 같다.

  • 기댓값
E[b(k)]=(X′X+kI)−1X′E(y)=(X′X+kI)−1X′Xβ≠β\begin{aligned} E[b(k)] &= (X'X+kI)^{-1}X'E(y) \\[10pt] &=(X'X+kI)^{-1}X'X\beta\\[10pt] &\neq \beta \end{aligned}

가 되어, 편향(biased) 추정량이다.

  • 분산
Var[b(k)]=(X′X+kI)−1X′X(X′X+kI)−1σ2\text{Var}[b(k)]=(X'X+kI)^{-1}X'X(X'X+kI)^{-1}\sigma^2 \\[15pt]

잔차제곱합 (SSE(k))(SSE(k))은 다음과 같으므로 bb의 잔차제곱합(SSE)(SSE) 보다 크다.

[y−Xb(k)]′[y−Xb(k)]=[y−X{b+(b(k)−b)}]′[y−X{b+(b(k)−b)}]=[(y−Xb)+{X(b(k)−b)}]′[(y−Xb)+{X(b(k)−b)}]=(y−Xb)′(y−Xb)+[b(k)−b]′X′X[b(k)−b](∵(y−Xb)′(X(b(k)−b))+(X(b(k)−b))′(y−Xb)=0)\begin{aligned} [y-Xb(k)]'[y-Xb(k)] &= [y-X\{b+(b(k)-b)\}]'[y-X\{b+(b(k)-b)\}] \\[15pt] &= [(y-Xb)+\{X(b(k)-b)\}]'[(y-Xb)+\{X(b(k)-b)\}] \\[15pt] &= (y-Xb)'(y-Xb)+[b(k)-b]'X'X[b(k)-b] \\[15pt] (\because(y&-Xb)'(X(b(k)-b))+(X(b(k)-b))'(y-Xb)=0) \end{aligned}

마지막 식이 0이 되는 이유는 y−Xby−Xb잔차 벡터로 회귀 계수 추정값의 변화에 직교(orthogonal) 한다는 성질을 가지기 때문입니다. 즉, 회귀 계수 b(k)b(k)와 b(k)−bb(k)−b는 잔차 벡터 y−Xby−Xb와 직교하므로 두 항의 내적 결과는 0이 됩니다.

능형추정량 b(k)b(k)와 β\beta간의 거리의 제곱을 L2(k)L^2(k)이라 정의하고, 그 기대값은 (총)평균제곱오차 (MSE(b(k)))(MSE(b(k)))라고도 부르며 다음과 같다.

E(L2(k))=MSE(b(k))=E[(b(k)−β)′(b(k)−β)]=E{[b(k)−E(b(k))]′[b(k)−E(b(k))]}+[E(b(k))−β]′[E(b(k))−β](∵MSE=Var(y^)+(E(y^−y))2)=σ2∑i=1pλi(λi+k)−2+k2β(X′X+kI)−2β\begin{aligned} E(L^2(k)) &= MSE(b(k)) \\[10pt] &= E[(b(k)-\beta)'(b(k)-\beta)] \\[15pt] &= E\{[b(k)-E(b(k))]'[b(k)-E(b(k))]\} \\[15pt] &\quad + [E(b(k))-\beta]'[E(b(k))-\beta]\quad(\because MSE =\text{Var}(\hat{y})+(E(\hat{y}-y))^2) \\[10pt] &= \sigma^2\sum_{i=1}^{p}\lambda_i(\lambda_i+k)^{-2}+k^2\beta(X'X+kI)^{-2}\beta \end{aligned}

이고 첫번째 항은 b(k)b(k)의 총분산이며 둘째 항은 편향의 제곱에 해당된다. 이제 E[L2(k)]E[L^2(k)]을 E(L2)E(L^2)과 비교하였을 때

E[L2(k)]<E(L2)E[L^2(k)] < E(L^2)

이 성립되는 k>0k>0의 값이 항상 존재함을 증명함으로써 능형추정량이 최소제곱추정량보다 거리의 제곱이라는 관점에서 볼 때 더 우수할 수 있음을 입증하였다. 일반적으로 kk가 대단히 작은 값에서 총분산은 급격히 감소하고, 편향의(bias) 제곱은 서서히 증가하기 때문에 E[L2(k)]E[L^2(k)]가 E(L2)E(L^2)보다 작게되는 k>0k>0의 값이 존재한다. 그러나 kk가 너무 커지면 편의의 제곱이 급격히 증가하여 E[L2(k)]>E(L2)E[L^2(k)]>E(L^2)이 되게 될것이다.

\\[20pt]

참고로 머신러닝분야에서 릿지 회귀를

SSE(L2)=(y−Xb)′(y−Xb)+{b(k)−b}′X′X{b(k)−b}=∑(y−y^)2+λ∣∣w∣∣2(∣∣w∣∣2=∑b2)\begin{aligned} SSE(L_2) &=(y-Xb)'(y-Xb) +\{b(k)-b\}'X'X\{b(k)-b\} \\[15pt] &= \sum (y-\hat{y})^2 + \lambda||w||^2\quad(||w||^2=\sum b^2 )\\[15pt] \end{aligned}

으로 표기한다. 라소 회귀는 제곱을 절대값으로 변환한 식이다.

\\[40pt]

17.2.1 능형추정량 성질

행렬 X′XX'X의 고유값들의 행렬을 Λ\Lambda라고 하고 고유벡터들을 열로 하는 행렬을 PP라고 하자. 즉,

Λ=[λ10…00λ2…0⋮⋮⋱⋮00…λn]P=[p1,p2,⋯ ,pp]\Lambda= \begin{bmatrix} \lambda_1 & 0 & \ldots & 0 \\[10pt] 0 & \lambda_2 & \ldots & 0 \\[10pt] \vdots & \vdots & \ddots & \vdots \\[10pt] 0 & 0 & \ldots & \lambda_n \end{bmatrix}\quad P=[p_1,p_2,\cdots,p_p]

이며

(X′X−λiI)pi=0,i=1,2,…,p(X'X-\lambda_iI)p_i=0,\quad i=1,2,\ldots,p

를 만족하고 P′P=PP′=IP'P=PP'=I가 되며 또한

X′X=PΛPX'X=P\Lambda P

성립한다. 선형모형을 Z=XP, α=P′βZ=XP,\ \alpha=P'\beta로 변환시키면

y=Xβ+ϵ=XPP′β+ϵ=Zα+ϵ,ϵ∼N(0,Iσ2)\begin{aligned} y &= X\beta+\epsilon \\[10pt] &= XPP'\beta+\epsilon \\[10pt] &=Z\alpha+\epsilon,\quad \epsilon\sim N(0,I\sigma^2) \end{aligned}

으로 표현된다. 이 때 α\alpha의 최소제곱추정량은

α^=(Z′Z)−1Z′y=Λ−1Z′y(∵Z′Z=P′X′XP=Λ)\begin{aligned} \hat{\alpha} &= (Z'Z)^{-1}Z'y \\[10pt] &= \Lambda^{-1}Z'y\quad (\because Z'Z=P'X'XP=\Lambda) \end{aligned}

이며, 적합된 모형 y^=Zα^\hat{y}=Z\hat{\alpha}의 잔차제곱합은

SSE(α^)=(y−Zα^)′(y−Zα^)=y′y−α^Z′y\begin{aligned} SSE(\hat{\alpha}) &= (y-Z\hat{\alpha})'(y-Z\hat{\alpha}) \\[10pt] &=y'y-\hat{\alpha}Z'y \end{aligned}

이다. α^\hat{\alpha}의 분산과 분산-공분산 행렬은 다음과 같다.

  • 분산
Var(α^)=Var(Λ−1Z′y)=Λ−1Z′yVar(y)ZΛ−1=σ2Λ−1\begin{aligned} \text{Var}(\hat{\alpha}) &= \text{Var}(\Lambda^{-1}Z'y) \\[10pt] &= \Lambda^{-1}Z'y\text{Var}(y)Z\Lambda^{-1} \\[10pt] &=\sigma^2\Lambda^{-1} \end{aligned}
  • 분산-공분산 행렬

앞에서 정의한 바와 같이 α\alpha의 능형추정량이

α^(k)=(Z′Z+kI)−1Z′y=(Λ+kI)−1Z′y=(Λ+kI)−1Λα^\begin{aligned} \hat{\alpha}(k) &= (Z'Z+kI)^{-1}Z'y \\[10pt] &= (\Lambda+kI)^{-1}Z'y \\[10pt] &= (\Lambda+kI)^{-1}\Lambda\hat{\alpha} \\[10pt] \end{aligned}

이므로, α^(k)\hat{\alpha}(k)의 분산-공분산행렬은

Var[α^(k)]=Var[(Λ+kI)−1Z′y]=(Λ+kI)−1Z′⋅Var[y]⋅Z(Λ+kI)−1=σ2Λ(Λ+kI)−2\begin{aligned} \text{Var}[\hat{\alpha}(k)] &= \text{Var}[(\Lambda+kI)^{-1}Z'y] \\[10pt] &= (\Lambda+kI)^{-1}Z'\cdot\text{Var}[y]\cdot Z(\Lambda+kI)^{-1} \\[10pt] &= \sigma^2\Lambda(\Lambda+kI)^{-2} \end{aligned}

로서 대각 행렬이 되며 α^(k)\hat{\alpha}(k)가 갖는 편향의 크기는

E[α^(k)−α]=[(Λ+kI)−1Λ−I]αE[\hat{\alpha}(k)-\alpha]=[(\Lambda+kI)^{-1}\Lambda-I]\alpha

이다. 따라서, 총평균제곱오차는

E[L2(k)]=E[(α^(k)−α)′(α^(k)−α)]=tr{Var[α^(k)]}+([(Λ+kI)−1Λ−I]α)2=σ2tr[Λ(Λ+kI)−2]+α′[(Λ+kI)−1Λ−I]′[(Λ+kI)−1Λ−I]α=σ2∑i=1pλi(λi+k)2+∑i=1pk2αi2(λi+k)2\begin{aligned} E[L^2(k)] &= E[(\hat{\alpha}(k)-\alpha)'(\hat{\alpha}(k)-\alpha)] \\[10pt] &= \text{tr}\{\text{Var}[\hat{\alpha}(k)]\}+([(\Lambda+kI)^{-1}\Lambda-I]\alpha)^2 \\[10pt] &= \sigma^2\text{tr}[\Lambda(\Lambda+kI)^{-2}] + \alpha'[(\Lambda+kI)^{-1}\Lambda-I]'[(\Lambda+kI)^{-1}\Lambda-I]\alpha \\[10pt] &= \sigma^2\sum_{i=1}^{p}\dfrac{\lambda_i}{(\lambda_i+k)^2} + \sum_{i=1}^{p}\dfrac{k^2\alpha_i^2}{(\lambda_i+k)^2} \end{aligned}

이 됨을 증명할 수 있다. 여기서 α′=(α1,α2,…,αn)\alpha'=(\alpha_1,\alpha_2,\ldots,\alpha_n)으로 나타내 준다.

E[L2(k)]E[L^2(k)]를 kk로 미분하면

dE[L2(k)]dk=−2σ2∑i=1pλi(λi+k)3+2k∑i=1pλiαi2(λi+k)3\dfrac{dE[L^2(k)]}{dk}=-2\sigma^2\sum_{i=1}^{p}\dfrac{\lambda_i}{(\lambda_i+k)^3} + 2k\sum_{i=1}^{p}\dfrac{\lambda_i\alpha_i^2}{(\lambda_i+k)^3}

이 되는데, k=0k=0에서 E[L2(0)]=E[L2]E[L^2(0)]=E[L^2]이며, L2(k)L^2(k)는 kk의 연속함수이고 k→0+k\rightarrow0+로 접근할 때에

lim⁡k→0+dE[L2(k)]dk=−2σ2∑i=1p1λi2<0\lim_{k\rightarrow0+}\dfrac{dE[L^2(k)]}{dk}=-2\sigma^2\sum_{i=1}^{p}\dfrac{1}{\lambda_i^2}<0

으로 음이므로 E[L2(k)]<E[L2]E[L^2(k)]<E[L^2]인 kk의 값이 존재한다.

E[L2(k)]=σ2∑i=1pλi(λi+k)2+∑i=1pk2αi2(λi+k)2\begin{aligned} E[L^2(k)] &= \sigma^2\sum_{i=1}^{p}\dfrac{\lambda_i}{(\lambda_i+k)^2} + \sum_{i=1}^{p}\dfrac{k^2\alpha_i^2}{(\lambda_i+k)^2} \end{aligned}

을 최소로 하는 kk의 값을 간단한 공식으로 구할 수는 없으나, 만약 X′X=IX'X=I라면 E[L2(k)]E[L^2(k)]가

k=pσ2α′αk=\dfrac{p\sigma^2}{\alpha'\alpha}

일 때에 최소가 된다. 만약 능형추정량의 정의를 보다 더 일반적으로

b=(X′X+K)−1X′y,K=[k10…00k2…0⋮⋮⋱⋮0k2…kp]b=(X'X+K)^{-1}X'y,\quad K=\begin{bmatrix} k_1 & 0 & \ldots 0 \\[10pt] 0 & k_2 & \ldots 0 \\[10pt] \vdots & \vdots & \ddots \vdots \\[10pt] 0 & k_2 & \ldots k_p \end{bmatrix}

으로 정의하여, pp개의 다른 kk 값을 허용한다면 α^(K)\hat{\alpha}(K)를

α^(K)=(Λ+K)−1Λα^\hat{\alpha}(K) = (\Lambda+K)^{-1}\Lambda\hat{\alpha}

으로 정의할 수 있다. 이 경우 총평균제곱오차는

E[L2(K)]=σ2∑i=1pλi(λi+ki)2+∑i=1pki2αi2(λi+ki)2\begin{aligned} E[L^2(K)] &= \sigma^2\sum_{i=1}^{p}\dfrac{\lambda_i}{(\lambda_i+k_i)^2} + \sum_{i=1}^{p}\dfrac{k_i^2\alpha_i^2}{(\lambda_i+k_i)^2} \end{aligned}

이 되며,

ki=σ2αi2k_i=\dfrac{\sigma^2}{\alpha_i^2}

에서 최소가 됨을 쉽게 보일 수 있다.

\\[40pt]

17.2.2 kk 추정

방법 1

능형추정값 b(k)=(X′X+kI)X′yb(k)=(X'X+kI)X'y을 kk를 0에서 1까지 증가하였을 때 각각의

b′(k)=(b1(k),b2(k),…,bp(k))b'(k)=(b_1(k),b_2(k),\ldots,b_p(k))

가 변화하는 것을 그림으로 도시하고 bi(k)b_i(k)들의 궤적을 능형 트레이스(trace)라 부르고, 이 그림으로부터 bi(k)b_i(k)값들이 급격히 변화하지 않고 안정되어 가는 kk의 값을 선택해 주어야 한다. 이렇게 선택된 kk의 값을 능형추정값 b(k)=(X′X+kI)X′yb(k)=(X'X+kI)X'y에 대입하면 다중공선성이 거의 제거된 상태가 되며, 능형회귀에서 가장 많이 사용되는 방법이다. 그러나 이 방법은 많은 계산이 필요하고 주관에 따라 kk의 값이 약간씩 다를 수 있다는 점이다.

\\[30pt]

방법 2

b(k)=(X′X+kI)X′yb(k)=(X'X+kI)X'y의 관계식으로부터 kk의 추정을 다음과 같이 한다.

k=pσ^2αi^′αi^k=\dfrac{p\hat{\sigma}^2}{\hat{\alpha_i}'\hat{\alpha_i}}

즉, σ2\sigma^2 대신 불편추정값 σ^2\hat{\sigma}^2을 구하여 대입시키고 α\alpha 대신에

α^=P′b=P′(X′X)−1X′y\hat{\alpha}=P'b=P'(X'X)^{-1}X'y

를 대입하여 kk의 값을 구하는 것이다. 이처럼 한 번에 kk의 값을 구하기도 하고, 반복적인 방법을 사용하여 다음과 같이 구하기도 한다.

k=pσ^2αi^′αi^k=\dfrac{p\hat{\sigma}^2}{\hat{\alpha_i}'\hat{\alpha_i}}

에 의하여 얻어진 kk를 k1k_1이라 하고 이 값으로

b(k1)=(X′X+k1I)X′yb(k_1)=(X'X+k_1I)X'y

을 구한 후에

α^(k1)=P′b(k1)\hat{\alpha}(k_1)=P'b(k_1)

을 구한다. 다음으로 이 α^(k1)\hat{\alpha}(k_1)을 사용하여 kk의 두번째 값

k2=pσ^2αi^(k1)′αi^(k1)k_2 =\dfrac{p\hat{\sigma}^2}{\hat{\alpha_i}(k_1)'\hat{\alpha_i}(k_1)}

을 구한 후

b(k2)=(X′X+k2I)X′yα^(k2)=P′b(k2)b(k_2)=(X'X+k_2I)X'y \\[10pt] \hat{\alpha}(k_2)=P'b(k_2)

을 얻고, 다음의 kk의 값

k3=pσ^2αi^(k2)′αi^(k2)k_3 =\dfrac{p\hat{\sigma}^2}{\hat{\alpha_i}(k_2)'\hat{\alpha_i}(k_2)}

을 구한다. 이와 같은 방법으로 반복적으로 하여 k1,k2,k3,…,kqk_1,k_2,k_3,\ldots,k_q를 구하여 가는 과정에서 만약 kq−1k_{q-1}의 값과 kqk_q의 값에 유의한 차이가 없으면 kq−1k_{q-1}의 값을 kk의 값으로 정하여준다

[참고문헌]

  • 회귀분석 제 3판 - 박성현
profile
choyunjeong
이전 포스트

편의 추정

다음 포스트

주성분 분석

0개의 댓글