R lm 함수 뜯어보기

Published

March 24, 2024

다음과 같은 회귀식을 생각한다.

$y = X β + ϵ$

$X$ 는 n $\times$ p의 행렬이고, $y$ 는 n $\times$ 1의 행렬이다.

우리는 보통 R의 lm 함수를 사용할 것이다. lm 함수를 직접 구현해보면서, 선형회귀의 핵심을 알아보자.

OLS 추정에 따르면, 오차를 최소화 하는 $\hat{β}$ 는 아래와 같다.

$\hat{β} = (X^{T} X)^{- 1} X^{T} y$ $\hat{β}$ 가 구해지면 예측값 $X \hat{β}$ 를 구할 수 있고, 그러면 예측잔차 $\hat{ϵ} = y - X \hat{β}$ 도 구해진다.

또한 회귀분석에서는 $\hat{β}$ 의 효과가 0과 다르다는 것을 검정하기 위해서, 회귀계수의 분산 $s d (\hat{β})$ 을 알 필요가 있다. 이는 다음과 같다.

$s d (\hat{β}) = \sqrt{σ^{2} d i a g (X^{T} X)^{- 1}}$ $d i a g$ 는 대각행렬의 원소를 의미한다.

다만, $σ^{2}$ 가 알려지지 않았을 때는 $s^{2} = \frac{Σ_{i = 1}^{n} {\hat{ϵ}}^{2}}{(n - p)}$ 로 대체 가능하다. 즉, 회귀계수의 분산식은 $s d (\hat{β}) = \sqrt{s^{2} d i a g (X^{T} X)^{- 1}}, s^{2} = \frac{Σ_{i = 1}^{n} {\hat{ϵ}}^{2}}{(n - p)}$ 로 대체된다.

예제

$X$ 와 $y$ 가 아래와 같이 주어졌다고 하자.

n <- 7
p <- 3

X <- matrix(c(rep(1,n), 
              -floor(n/2):floor(n/2),
              (0:(n-1))^2),
            nrow=n, ncol=p); X

     [,1] [,2] [,3]
[1,]    1   -3    0
[2,]    1   -2    1
[3,]    1   -1    4
[4,]    1    0    9
[5,]    1    1   16
[6,]    1    2   25
[7,]    1    3   36

y <- c(1,2,3,5,7,9,11); y

[1]  1  2  3  5  7  9 11

$(X^{T} X)^{- 1} X^{T} y$

beta.hat <- solve(t(X) %*% X) %*% t(X) %*% y; beta.hat

          [,1]
[1,] 3.8809524
[2,] 1.0000000
[3,] 0.1190476

$X \hat{β}$

X %*% beta.hat

           [,1]
[1,]  0.8809524
[2,]  2.0000000
[3,]  3.3571429
[4,]  4.9523810
[5,]  6.7857143
[6,]  8.8571429
[7,] 11.1666667

$\hat{ϵ} = y - X \hat{β}$

residuals <- y - X %*% beta.hat; residuals

              [,1]
[1,]  1.190476e-01
[2,] -1.776357e-15
[3,] -3.571429e-01
[4,]  4.761905e-02
[5,]  2.142857e-01
[6,]  1.428571e-01
[7,] -1.666667e-01

$s d (\hat{β}) = \sqrt{s^{2} (X^{T} X)^{- 1}}$

s.squared <- sum((y - X %*% beta.hat)^2)/(n-p)
sqrt(s.squared * diag(solve(t(X) %*% X)))

[1] 0.35813355 0.16624095 0.02661986

lm 함수의 summary와 비교하여 결과가 동일한지 확인한다.

summary(lm(y ~ X[,(2:3)]))


Call:
lm(formula = y ~ X[, (2:3)])

Residuals:
         1          2          3          4          5          6          7 
 1.190e-01  2.137e-15 -3.571e-01  4.762e-02  2.143e-01  1.429e-01 -1.667e-01 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  3.88095    0.35813  10.837 0.000411 ***
X[, (2:3)]1  1.00000    0.16624   6.015 0.003846 ** 
X[, (2:3)]2  0.11905    0.02662   4.472 0.011056 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.244 on 4 degrees of freedom
Multiple R-squared:  0.9972,    Adjusted R-squared:  0.9957 
F-statistic: 701.2 on 2 and 4 DF,  p-value: 8.089e-06

회귀계수(Coefficients)의 추정치와 표준오차, 그리고 잔차(Residuals)가 동일함을 확인할 수 있다.