10.1 단순선형모형의 표현
좁은 의미의 선형모형(Linear Model)은 연속형 결과변수(반응변수, 종속변수, y 변수)를 연속형 또는 범주형 설명변수(독립변수, x 변수)의 선형식으로 설명하고자 하는 모형이다.
설명변수(x 변수)가 연속형인 경우의 선형분석방법이 선형회귀분석(regression analysis)이고, 설명변수가 범주형인 경우의 분석방법이 분산분석(analysis of variance, ANOVA)이다. 설명변수에 연속형과 범주형이 섞여 있는 경우 공분산분석(analysis of covariance, ANCOVA)을 한다. 실제로 분산분석이나 공분산분석을 하기 위해서는 (내부적으로) 회귀적합(regression fitting)을 먼저 한다.
일반적으로는 t-검정, 카이제곱검정, 선형회귀분석, 분산분석, 공분산분석, 이분형 로지스틱 회귀분석의 모형을 통칭해서 선형모형(linear model)이라 한다.
혹자는 분포를 지수족 분포함수로 확대한 일반화 선형 모형(generalized linear model)과 무작위 효과(random effect, 변량효과, 임의효과)를 함께 고려한 혼합효과 모형(mixed effects model)까지 포함해서 선형모형이라고 생각한다.
연구자가 상정한 어떠한 모형이 다음과 같은 수식 관계를 가지고 있다고 하자.
\[y = \beta_0 + \beta_1 \cdot x_1 + \beta_2 \cdot x_2 + \cdots + \beta_k \cdot x_k\]
이 때 \(\beta_0\)를 절편(intercept)이라고 한다. \(\beta_1, \cdots, \beta_k\)는 그에 상응하는 \(x_i\)들의 계수(coefficient)이다.
위의 모형수식에는 확률변수를 포함시키지 않았고, 이렇게 확률변수가 포함되지 않은 모형을 구조모형(structural model)이라고 한다.
만약 위의 수식을 행렬로 나타내면 다음과 같다. \[y = X\beta\]
이때 X를 설계행렬(design matrix)라고 하며, 이제 y와 \(\beta\)는 이제 scalar가 아닌 열벡터(column vector)를 의미한다.
절편이 있는 모형의 경우 절편(\(\beta_0\))과 곱해질 값이 필요해서 설계행렬의 첫 열(column)은 모두 1이 된다.
위의 구조모형 수식만으로 어떻게 독립변수(x 변수)로 결과변수(y 변수)를 계산하는지 알 수 있다. 하지만, 실제 자료에는 잡음(noise)이 섞여 있고, 이를 통계에서는 확률변수 \(\epsilon\)을 이용해서 표현한다.
따라서 위의 구조모형이 예측하는 값을 F라고 할 때, 가장 단순한 통계모형은 additive homoscedastic error model로서 다음과 같이 쓸 수 있다.
\[Y = F + \epsilon\]
결과변수를 대문자 Y로 표기한 것은, 결과변수가 확률변수가 되었기 때문이다. 즉, \(\epsilon\)이 확률변수이며, 여기에다 일반적인 사칙연산을 한 결과는 새로운 확률변수이기 때문이다.
이를 행렬이 포함된 수식으로 표현하자면 다음과 같다.
\[Y = X\beta + \epsilon\]
이제 \(\epsilon\)도 Y와 차원(dimension)이 같은 column vector가 된다.
이제 k가 3까지 있는 경우를 R의 formula로 나타내면 다음과 같다.
y ~ x1 + x2 + x3
위 식에서 \(\beta\)들은 따로 표현하지 않는다. 절편을 따로 표시하지 않아도 절편이 포함된 모형을 나타낸다. 굳이 절편이 없는 모형을 표현하고 싶다면 -1이나 +0을 다음과 같이 추가하여야 한다.
y ~ x1 + x2 + x3 - 1
y ~ x1 + x2 + x3 + 0
선형모형에서는 절편이 있는 모형을 주로 사용하기 때문에 따로 명시하지 않는 한, 절편이 있는 모형을 의미한다.
\(x_i\)가 연속형 변수일 때는 연속형 변수 하나는 설계행렬에서 하나의 열을 차지하며, 값은 그대로 dataset과 같은 값을 가진다.
R에서는 설계 행렬(design matix)을 모형 행렬(model matrix)라고 부르며, model.matrix라는 함수로 볼 수 있다.
= data.frame(Age = c(25, 51, 75, 34, 21, 48))
d1 model.matrix(~ Age, d1)
(Intercept) Age
1 1 25
2 1 51
3 1 75
4 1 34
5 1 21
6 1 48
attr(,"assign")
[1] 0 1
하지만, \(x_i\)가 범주형 변수일 때는 그 범주형 변수가 갖는 (원래는) 수준(level)의 수만큼 열(column)을 갖게 된다.
= data.frame(Treatment = c("A", "A", "B", "B", "C", "C"))
d2 = model.matrix(~ Treatment, d2) ; x1 x1
(Intercept) TreatmentB TreatmentC
1 1 0 0
2 1 0 0
3 1 1 0
4 1 1 0
5 1 0 1
6 1 0 1
attr(,"assign")
[1] 0 1 1
attr(,"contrasts")
attr(,"contrasts")$Treatment
[1] "contr.treatment"
범주형 자료의 설계행렬은 각 수준(level)별로 하나의 열(column)을 가지며, 그 값은 0 (FALSE의 의미) 또는 1 (TRUE의 의미)만 갖는다. 그런데, 위의 결과를 보면 Treatment A를 받은 경우를 표시하는 컬럼이 나오지 않는다. 그러면, 왜 R의 내장 함수인 model.matrix는 열(column)이 하나 빠진 설계행렬을 return할까?
이는 절편(intercept)이 포함된 모형에서는 모든 수준을 나타내려면 \(\beta\) coefficient를 구하는 중간 과정에 필요한 \(X'X\)의 역행렬이 존재하지 않기 때문이다. 역행렬이 존재하지 않는 행렬을 선행대수(linear algebra)에서는 ’비정칙 행렬(singular matrix)’라고 하며, 비정칙(singular)임을 확인하는 방법으로 대표적인 것은 행렬식(determinant)가 0인지 확인하는 것이다.
만약 모든 수준(level)이 포함된 설계행렬을 model.matrix로 보고 싶다면 절편을 제외한 모형을 써서 볼 수 있다.
= model.matrix(~ Treatment - 1, d2) ; x0 x0
TreatmentA TreatmentB TreatmentC
1 1 0 0
2 1 0 0
3 0 1 0
4 0 1 0
5 0 0 1
6 0 0 1
attr(,"assign")
[1] 1 1 1
attr(,"contrasts")
attr(,"contrasts")$Treatment
[1] "contr.treatment"
이러한 경우들에서는 \(X'X\) 행렬이 역행렬을 가짐을 알 수 있다.
solve(crossprod(x1))
(Intercept) TreatmentB TreatmentC
(Intercept) 0.5 -0.5 -0.5
TreatmentB -0.5 1.0 0.5
TreatmentC -0.5 0.5 1.0
solve(crossprod(x0))
TreatmentA TreatmentB TreatmentC
TreatmentA 0.5 0.0 0.0
TreatmentB 0.0 0.5 0.0
TreatmentC 0.0 0.0 0.5
만약 절편과 함께 모든 수준이 포함된 설계행렬을 보고 싶다면 sasLM package의 ModelMatrix 함수를 이용하여 다음과 같이 볼 수 있다.
require(sasLM)
= ModelMatrix(~ Treatment, d2) ; x2$X x2
(Intercept) TreatmentA TreatmentB TreatmentC
1 1 1 0 0
2 1 1 0 0
3 1 0 1 0
4 1 0 1 0
5 1 0 0 1
6 1 0 0 1
그런데, 이 경우에는 역행렬이 존재하지 않음을 다음과 같이 확인할 수 있다.
solve(crossprod(x2$X))
Error in solve.default(crossprod(x2$X)) :
system is computationally singular: reciprocal condition number = 1.38778e-17
참고로 ModelMatrix 함수는 설계행렬 뿐 아니라, 다른 정보도 return하기 때문에 설계행렬만 보려면 뒤에 ’$X’를 붙여주어야 한다.