11.3 선형 회귀 모형의 분산 분석
관찰값 \(y_i\) 들의 전체 평균을 \(\bar{y}\)라 하고, (선형이든 비선형이든) 모형에 의한 y값 (즉, 잔차/오차 부분이 없는 값) 을 \(\hat{y}\) 라 했을 때 다음의 식이 성립한다.
\(y_i - \bar{y} = y_i - \hat{y}_i + \hat{y}_i - \bar{y}\)
신기하게도 다음도 성립한다.
\(\sum (y_i - \bar{y})^2 = \sum (y_i - \hat{y}_i)^2 + \sum (\hat{y}_i - \bar{y})^2\)
위의 것은
총변동 = (회귀모형이 설명하지 못하는 변동) + (회귀 모형이 설명하는 변동)
이라고 할 수 있으며 영어로는
SST (sum of squares total) = SSE (sum of squares residual/error) + SSR(sum of squares regression)
Error(오차)는 주로 확률변수를 가리키고, 잔차(residual)은 주로 값(value)을 의미하지만, 경우에 따라서 위 처럼 혼용하고 있다.
설명력(\(R^2\))을 다음과 같이 정의하여 모형이 얼마나 잘 적합되었는지를 나타내는 지표로 쓸 수 있다.
\(R^2 = \frac{SSR}{SST} = 1 - \frac{SSE}{SST}\)
제곱합(sum of square)을 자신의 자유도로 나누어 준 것을 평균제곱합(mean square)이라고 하는데, 이것의 비가 F 분포를 따른다.
\(\frac{SSR/df_{SSR}}{SSE/df_{SSE}} = MSR/MSE \sim F(df_{SSR}, df_{SSE})\)
따라서, SSR이 SSE에 비해 통계적으로 유의하게 큰 지(비율이 1보다 큰지)를 통계적으로 (단측) 검정을 할 수 있다.
일반적으로 전체평균(\(\bar{Y}\))으로부터의 deviation을 제곱합한 총변동을 corrected sum of square (corrected SST)라 하고, 그렇지 않은 \(\sum y_i^2\)을 uncorrected sum of square (uncorrected SST)라 한다. 분산 분석에서는 전자를 주로 사용하기 때문에 (전체 평균 또는 절편(intercept)가 자유도를 하나 소모하여) 분산 분석에서 자유도가 n - 1이 된다.
Corrected sum of squared total을 사용하는 분산분석표는 다음과 같다.
Source: Source of Variability
SS : Sum of Square
DF : Degree of Freedom
MS : Mean Square
n : count of observations, \(y_1, y_2, \cdots, y_n\)
p : count of parameters, \(\beta_0, \beta_1, \cdots, \beta_p\)
\(R^2\)는 설명력을 나타내므로 모형 간에 더 우수한 모형을 찾고자 할 때 쓰기도 한다. (모형 간의 비교 방법으로 다른 방법들이 있어서, 현재는 이 방법은 추천하지 않는다.) 그런데, \(\beta\) parameter의 개수를 늘리면 \(R^2\)는 항상 더 커지므로, 관찰값 개수와 parameter 개수로 보정해 준 것을 더 많이 쓰는데, 그것을 adjusted \(R^2\)라고 부른다.
Adjusted \(R^2\)를 구하는 공식은 다음과 같으며, \(R^2\)는 항상 0이상 1이하의 값을 갖지만, adjusted \(R^2\)는 음수가 나오기도 한다.
\(R_{adjusted}^2 = 1 - (1 - R^2)\frac{(DF_{corrected SST})}{DF_{residual}}\)
앞의 예제 (\(y \sim x_1 + x_2\) model)의 \(R^2\)와 분산 분석 결과는 다음과 같다.
SST = sum((y - mean(y))^2) # corrected sum of square total
SSE = sum((y - yhat)^2) # sum of square error
SSR = SST - SSE # sum of square regression
Rsq = SSR/SST ; Rsq # R-square[1] 0.9745628n = length(y) # number of observation
p = length(b) # number of parameters
R2adj = 1 - (1 - Rsq)*(n - 1)/(n - p) ; R2adj # Adjusted R-square[1] 0.9236884SS = c(SSR, SSE, SST)
DF = c(p - 1, n - p, n - 1)
MS = c(SS[1:2]/DF[1:2], NA)
Fval = MS[1]/MS[2]
alpha = 0.05
Fcrit = qf(1 - alpha, p - 1, n - p)
pval = 1 - pf(Fval, p - 1, n - p)
r3 = cbind(SS, DF, MS, c(Fval, NA, NA), c(Fcrit, NA, NA), c(pval, NA, NA))
rownames(r3) = c("Regression", "Error", "Total")
colnames(r3) = c("Sum Sq", "Df", "Mean Sq", "F value", "F crit", "Pr(>F)")
class(r3) = "anova" # to see in a better way
r3 # compare with anova(r2) and summary(r2) Sum Sq Df Mean Sq F value F crit Pr(>F)
Regression 613 2 306.5 19.156 199.5 0.1595
Error 16 1 16.0
Total 629 3