library(tidyverse)
library(here)
library(knitr)
library(CCA)
#아래 3 문장은 한글을 포함한 ggplot 그림이 포함된 HTML, PDF로 만드는 경우 사용
library(showtext)
#font_add_google("Nanum Pen Script", "gl")
font_add_google(name = "Noto Sans KR", family = "noto")
showtext_auto()7 정준상관분석
7.1 상관계수
두 개의 확률변수의 상관계수(Correlation Cefficient)는 두 변수의 선형 관계의 정도를 나타내는 측도이다. 두 개의 확률변수 \(X_1\)과 \(X_2\) 의 상관계수 \(\rho\) 는 다음과 같이 정의된다.
\[ \rho =\rho(X_1, X_2) = cor(X_1, X_2) = \frac{ Cov(X_1,X_2)} { \sqrt{ Var(X_1) Var(X_2)}} \tag{7.1}\]
상관계수 \(\rho\)는 -1과 1 사이의 값을 가지며 상관계수가 0 이면 두 확률변수의 선형관계는 존재하지 않는다. 상관계수가 1 에 가까울수록 두 변수는 양의 선형관계가 강해지며 반대로 -1 에 가까울수록 두 변수는 음의 선형관계가 강해진다. 두 변수의 상관계수가 0이라고 해서 관계(relation)가 없다고 단정할 수 없다. 왜냐하면 상관계수는 두 변수의 선형관계(linear relationship)만을 나타내는 측도이고 비선형 관계 등 다른 특별한 관계를 반영하지는 못한다.
두 개의 변수에 대한 \(n\)개의 독립표본 \((X_{11},X_{12}),(X_{21},X_{22}), \dots, (X_{n1},X_{n2})\)이 주어지면 표본으로 부터 모집단의 상관계수를 추정할 수 있는 표본 상관계수 \(\hat \rho\) 을 다음과 같이 계산할 수 있다.
\[ \hat \rho = \frac{ \sum_{i=1}^n (X_{i1} -\bar X_1)(X_{i2}-\bar X_2)} { \sqrt{ \sum_{i=1}^n (X_{i1} - \bar X_1)^2 (X_{i2} - \bar X_2)^2}} \tag{7.2}\]
여기서 \(\bar X_1\)과 \(\bar X_2\)는 각각 첫 번째 변수의 표본 \(X_{11},X_{21},\dots X_{n1}\)과 두 번째 변수의 표본 \(X_{12},X_{22},\dots X_{n2}\)의 평균이다.
보기 7.1 (스위스의 47개 주 자료) R 패키지에 내장된 swiss 자료는 1988년 스위스의 47개 주에 대한 출산율과 사회경제변수를 모아놓은 자료이다(\(n=47\)). 6개의 변수에 대한 설명과 자료의 일부는 다음과 같다.
Fertility: \(I_g\), common standardized fertility measureAgriculture: % of males involved in agriculture as occupationExamination: % draftees receiving highest mark on army examinationEducation: % education beyond primary school for draftees.Catholic: % catholic (as opposed to protestant).Infant.Mortality: live births who live less than 1 year.
head(swiss) Fertility Agriculture Examination Education Catholic
Courtelary 80.2 17.0 15 12 9.96
Delemont 83.1 45.1 6 9 84.84
Franches-Mnt 92.5 39.7 5 5 93.40
Moutier 85.8 36.5 12 7 33.77
Neuveville 76.9 43.5 17 15 5.16
Porrentruy 76.1 35.3 9 7 90.57
Infant.Mortality
Courtelary 22.2
Delemont 22.2
Franches-Mnt 20.2
Moutier 20.3
Neuveville 20.6
Porrentruy 26.6자료에서 두개의 변수에 대한 상관계수는 cor 함수로 다음과 같이 계산할 수 있고 더 나아가 데이터프레임 swiss 에 있는 모든 변수들에 대한 상관계수행렬도 동시에 계산할 수 있다.
cor(swiss) Fertility Agriculture Examination Education Catholic
Fertility 1.0000000 0.35307918 -0.6458827 -0.66378886 0.4636847
Agriculture 0.3530792 1.00000000 -0.6865422 -0.63952252 0.4010951
Examination -0.6458827 -0.68654221 1.0000000 0.69841530 -0.5727418
Education -0.6637889 -0.63952252 0.6984153 1.00000000 -0.1538589
Catholic 0.4636847 0.40109505 -0.5727418 -0.15385892 1.0000000
Infant.Mortality 0.4165560 -0.06085861 -0.1140216 -0.09932185 0.1754959
Infant.Mortality
Fertility 0.41655603
Agriculture -0.06085861
Examination -0.11402160
Education -0.09932185
Catholic 0.17549591
Infant.Mortality 1.000000007.2 다중상관계수
상관계수는 두 개의 확률변수에 대한 선형관계를 나타내는 측도이다. 두 개의 변수에 대한 측도인 상관계수를 를 두 개의 확률벡터에 대한 관계를 나타내는 측도로 확장할 수 있다.
이렇게 두 개의 확률벡터에 대한 상관관계를 나타내는 측도를 정준상관계수(Canonical Correlation Coeffcient)라고 한다. 이 절에서는 두 개의 확률벡터에 대한 관계를 측정하는 정준상관계수를 정의하기 전에 하나의 확률변수와 여러 개의 변수를 포함하는 확률벡터의 관계를 나타내는 다중상관계수(Multiple Correlation Coefficient)을 먼저 정의하고자 한다.
확률벡터 \(\pmb X\)를 \(p\)개의 확률변수로 이루어졌다고 하고 하나의 확률변수 \(X_1\) (편의상 첫 번째 확률변수를 선택하였다)와 나머지 \(p-1\)개의 변수로 구성된 확률벡터를 \(\pmb X_*\)라고 하자.
\[ \pmb X = \begin{bmatrix} X_1 \\ X_2 \\ X_3 \\ \vdots \\ X_p \end{bmatrix} = \begin{bmatrix} X_1 \\ \pmb X_* \end{bmatrix}, \quad Cov(\pmb X) = \begin{bmatrix} \sigma_{11} & \pmb \sigma_{12}^t \\ \pmb \sigma_{12} & \pmb \Sigma_{22} \end{bmatrix} \] 위의 식에서 \(\sigma_{11}\)은 \(X_1\)의 분산, \(\pmb \sigma_{12}\)은 \(X_1\)과 \(\pmb X_*\)의 (\(p-1\))-차원 공분산 벡터, \(\pmb \Sigma_{22}\)는 \(\pmb X_*\)의 \((p-1)\times (p-1)\)-차원 공분산행렬이다.
하나의 확률변수와 여러 개의 확률변수로 구성된 확률벡터간의 선형관계는 \(X_1\)과 \(\pmb X_*\)에 포함된 각각의 확률변수에 대하여 각각 \(p-1\)개의 상관계수를 구하여 따로 따로 파악할 수 있다.
\[ \rho_{12}=cor(X_1,X_2), ~ \rho_{13}=cor(X_1,X_3), \dots, ~ \rho_{1p}=cor(X_1,X_p) \]
하지만 이러한 관계는 각각 두 개의 변수들에 대한 관계로서 하나의 확률변수와 확률벡터의 관계를 하나의 측도로서 종합적으로 반영하는 것은 아니다. 이렇게 확률변수와 확률벡터의 관계를 하나의 측도로 나타내기 위하여 주성분분석이나 회귀분석과 유사하게 확률벡터 \(\pmb X_*\)에 포함된 확률변수들의 선형조합을 생각하고 확률변수 \(X_1\)과 선형조합으로 만들어진 새로운 확률변수 \(\pmb a^t \pmb X_*\)의 상관관계를 생각해 보자.
\[ cor(X_1, a_2 X_2 + a_3 X_3 +\dots + a_p X_p) = cor(X_1,\pmb a^t \pmb X_*) \]
위의 식에서 선형조합으로 만들어진 새로운 확률변수와의 상관계수는 계수벡터 \(\pmb a\)의 계수값에 따라 달라진다. 이렇게 무수히 많은 값이 가능한 경우 최대의 상관계수를 가지는 계수를 고려하는 것이 다중상관계수의 정의이다.
\[ \rho(X_1, \pmb X_*) = \max_{\pmb a} cor(X_1, \pmb a^t \pmb X_*) \tag{7.3}\]
최대의 상관계수를 가지는 계수를 구하기 위하여 먼저 임의의 벡터 \(\pmb a\)에 대하여 \(X_1\)과 \(\pmb a^t \pmb X_*\)의 상관계수를 유도해보자. 먼저
\[ \begin{aligned} Cov(X_1,\pmb a^t \pmb X_*) & = E [(X_1 -E(X_1)][\pmb a^t \pmb X_*-E(\pmb a^t \pmb X_*) ] \\ &= \pmb a^t E [(X_1 -E(X_1)] [ \pmb X_*-E(\pmb X_*) ] \\ &= \pmb a^t \pmb \sigma_{12} \\ Var(\pmb a^t \pmb X_*) &= \pmb a^t Var( \pmb X_*) \pmb a \\ &= \pmb a^t \pmb \Sigma_{22} \pmb a \end{aligned} \]
따라서
\[ cor(X_1,\pmb a^t \pmb X_*) = \frac{ Cov(X_1,\pmb a^t \pmb X_*)} { \sqrt{ Var(X_1) Var(\pmb a^t \pmb X_*)}} = \frac{\pmb a^t \pmb \sigma_{12} } { \sqrt{\sigma_{11} (\pmb a^t \pmb \Sigma_{22} \pmb a ) }} \] 여기서 한 가지 유의해야할 점은 계수벡터 \(\pmb a\)의 계수들의 부호만을 바꾸면 다중상관계수의 부호가 바뀐다는 점이다.
\[ cor(X_1,-\pmb a^t \pmb X_*) = \frac{-\pmb a^t \pmb \sigma_{12} } { \sqrt{\sigma_{11} [(-\pmb a^t) \pmb \Sigma_{22} (-\pmb a )] }} = - cor(X_1,\pmb a^t \pmb X_*) \] 따라서 다중상관계수를 구하는 경우에는 그 부호에 관계없이 절대값이 가장 큰 경우를 고려해야 한다. 이러한 점을 고려하여 다중상관계수를 다음과 같이 자신의 제곱값을 통하여 유도할 수 있다. 이제 다중상관계수의 계수벡터를 구하기 위한 다중상관계수의 제곱에 대한 최대값을 구해보자.
\[ \begin{aligned} cor^2(X_1,\pmb a^t \pmb X_*) &= \frac{ (\pmb a^t \pmb \sigma_{12})^2 } { \sigma_{11} (\pmb a^t \pmb \Sigma_{22} \pmb a ) } \\ &= \frac{(\pmb a^t \pmb \Sigma_{22}^{1/2} \pmb \Sigma_{22}^{-1/2} \pmb \sigma_{12})^2 } { \sigma_{11} (\pmb a^t \pmb \Sigma_{22} \pmb a ) } \\ & \le \frac{ (\pmb a^t \pmb \Sigma_{22}^{1/2}\pmb \Sigma_{22}^{1/2} \pmb a)(\pmb \sigma_{12}^t \pmb \Sigma_{22}^{-1/2} \pmb \Sigma_{22}^{-1/2} \pmb \sigma_{12}) } { \sigma_{11} (\pmb a^t \pmb \Sigma_{22} \pmb a ) } \\ &= \frac{ (\pmb a^t \pmb \Sigma_{22} \pmb a)(\pmb \sigma_{12}^t \pmb \Sigma_{22}^{-1} \pmb \sigma_{12}) } { \sigma_{11} (\pmb a^t \pmb \Sigma_{22} \pmb a ) } \\ &= \frac{ \pmb \sigma_{12}^t \pmb \Sigma_{22}^{-1} \pmb \sigma_{12} } { \sigma_{11} } \\ \end{aligned} \]
위의 유도식에서 공분산 행렬 \(\pmb \Sigma\) 의 제곱근 행렬 \({\pmb \Sigma}^{1/2}\) 는 양정치 행렬의 성질 식 C.17 을 이용하였다.
또한 위의 유도식에서 사용한 부등식은 코쉬-쉬바르쯔(Cauchy-Schwarz) 부등식을 적용한 결과이다 (식 5.18 참조)
\[ (\pmb \alpha^t \pmb \beta)^2 \le (\pmb \alpha^t \pmb \alpha)(\pmb \beta^t \pmb \beta) \] 또한 위의 유도식에서 부등식의 등식이 성립하는 경우는 임의의 상수 \(\lambda\)에 대하여 다음과 같은 관계가 성립하는 경우이다.
\[ \Sigma_{22}^{1/2} \pmb a = \lambda \pmb \Sigma_{22}^{-1/2} \pmb \sigma_{12} \quad \rightarrow \quad \pmb a = \lambda \pmb \Sigma_{22}^{-1} \pmb \sigma_{12} \tag{7.4}\]
이제 \(X_1\)과 \(\pmb X_*\) 의 다중상관계수는 다음과 같이 구해진다.
\[ \rho(X_1, \pmb a^t \pmb X_*) = \frac{ (\pmb \sigma_{12}^t \pmb \Sigma_{22}^{-1} \pmb \sigma_{12})^{1/2} } { (\sigma_{11})^{1/2} } \tag{7.5}\]
여기서 \(\pmb a = \lambda \pmb \Sigma_{22}^{-1} \pmb \sigma_{12}\) , \(0 \le \rho(X_1, \pmb a^t \pmb X_*) \le 1\)
다중상관계수의 계산은 다음과 같이 상관계수들로 구할 수 있다.
\[ \rho(X_1, \pmb a^t \pmb X_*) = (\pmb \rho_{12}^t \pmb R_{22}^{-1} \pmb \rho_{12})^{1/2} \] 여기서 \(\pmb \rho_{12}\)는 \(X_1\)과 \(\pmb X_*\)의 상관계수 벡터이며 \(\pmb R_{22}\)는 \(\pmb X_*\)의 상관계수행렬이다. 이러한 계산식의 유도는 아래와 같은 공분산 행렬과 상관계수 행렬의 관계(2차원 확률벡터의 예)로부터 유도할 수 있다.
\[ \pmb R = \begin{bmatrix} 1 & \rho_{12} \\ \rho_{12} & 1 \\ \end{bmatrix} = \begin{bmatrix} 1/\sqrt{\sigma_{11}} & 0 \\ 0 & 1/\sqrt{\sigma_{22}} \\ \end{bmatrix} \begin{bmatrix} \sigma_{11} & \sigma_{12} \\ \sigma_{12} & \sigma_{22} \\ \end{bmatrix} \begin{bmatrix} 1/\sqrt{\sigma_{11}} & 0 \\ 0 & 1/\sqrt{\sigma_{22}} \\ \end{bmatrix} \]
표본으로서 자료가 주어진 경우에는 식 7.5 에서 각각의 모수에 대하여 그 추정량을 구하면 표본 다중상관계수를 구할 수 있다.
\[ \hat \rho(X_1, \pmb a^t \pmb X_*) = \frac{ ( \hat {\pmb \sigma}_{12}^t \hat {\pmb \Sigma}_{22}^{-1} \hat {\pmb \sigma}_{12})^{1/2} } { (\hat \sigma_{11})^{1/2} } = ( \hat {\pmb \rho}_{12}^t \hat {\pmb R}_{22}^{-1} \hat {\pmb \rho}_{12} )^{1/2} \]
swiss 자료에서 출산율(\(X_1\))과 나머지 5개의 사회경제변수(\(\pmb X_*\))의 표본 다중상관계수는 다음과 같이 구할 수 있다.
n<-dim(swiss)[1] ; p <- dim(swiss)[2]
dim(swiss)[1] 47 6R <- cor(swiss)
R12 <- matrix(R[2:p,1],p-1,1)
R12 [,1]
[1,] 0.3530792
[2,] -0.6458827
[3,] -0.6637889
[4,] 0.4636847
[5,] 0.4165560R22 <- matrix(R[2:p,2:p],p-1,p-1)
R22 [,1] [,2] [,3] [,4] [,5]
[1,] 1.00000000 -0.6865422 -0.63952252 0.4010951 -0.06085861
[2,] -0.68654221 1.0000000 0.69841530 -0.5727418 -0.11402160
[3,] -0.63952252 0.6984153 1.00000000 -0.1538589 -0.09932185
[4,] 0.40109505 -0.5727418 -0.15385892 1.0000000 0.17549591
[5,] -0.06085861 -0.1140216 -0.09932185 0.1754959 1.00000000mulcor <- sqrt(t(R12) %*% solve(R22) %*% R12)
mulcor [,1]
[1,] 0.8406753참고로 표본 다중상관계수의 제곱값은 는 \(X_1\)을 종속변수, \(\pmb X_*\)을 독립변수로 선형회귀직선(linear regression)을 적합하였을 경우 결정계수(coefficient of determination) \(R^2\) 와 같다.
\[
\hat \rho^2 (X_1, \pmb a^t \pmb X_*) = R^2
\] 아래 R 프로그램의 결과에서 확인할 수 있다.
res <- lm(Fertility~ Agriculture + Examination + Education + Catholic + Infant.Mortality,data=swiss)
summary(res)
Call:
lm(formula = Fertility ~ Agriculture + Examination + Education +
Catholic + Infant.Mortality, data = swiss)
Residuals:
Min 1Q Median 3Q Max
-15.2743 -5.2617 0.5032 4.1198 15.3213
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 66.91518 10.70604 6.250 1.91e-07 ***
Agriculture -0.17211 0.07030 -2.448 0.01873 *
Examination -0.25801 0.25388 -1.016 0.31546
Education -0.87094 0.18303 -4.758 2.43e-05 ***
Catholic 0.10412 0.03526 2.953 0.00519 **
Infant.Mortality 1.07705 0.38172 2.822 0.00734 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 7.165 on 41 degrees of freedom
Multiple R-squared: 0.7067, Adjusted R-squared: 0.671
F-statistic: 19.76 on 5 and 41 DF, p-value: 5.594e-10mulcor^2 [,1]
[1,] 0.7067357.3 정준상관계수
7.3.1 정준상관계수의 정의
이제 두 개의 확률벡터에 대한 상관관계를 나타내는 측도인 정준상관계수(Canonical Correlation Coeffcient)에 대하여 알아보자. 확률벡터 \(\pmb X\)가 두 개의 확률벡터 \(\pmb X_1\)과 \(\pmb X_2\)로 나누어져 있다고 가정하자. 두 확률 벡터의 차원은 각각 \(p\)와 \(q\) 라고 하자 (편의상 \(p \le q\)라고 가정한다.)
\[ \pmb X = \begin{bmatrix} \pmb X_1 \\ \pmb X_2 \\ \end{bmatrix} \quad Cov(\pmb X) = \begin{bmatrix} \pmb \Sigma_{11} & \pmb \Sigma_{12} \\ \pmb \Sigma_{12}^t & \pmb \Sigma_{22} \end{bmatrix} \]
앞 절에서와 유사한 방법으로 각 확률벡터의 선형조합으로 만들어진 두 개의 새로운 확률변수를 이용하여 정준상관계수를 다음과 같이 정의한다.
\[ \rho(\pmb X_1, \pmb X_2) = \max_{\pmb a, \pmb b}~~ cor(\pmb a^t \pmb X_1, \pmb b^t \pmb X_2) \tag{7.6}\]
앞절에서 공분산행렬을 구할 때 사용한 같은 방법을 적용하면
\[ \begin{aligned} Cov(\pmb a^t \pmb X_1,\pmb b^t \pmb X_2) &= E [(\pmb a^t \pmb X_1 -E(\pmb a^t \pmb X_1)][\pmb b^t \pmb X_2-E(\pmb b^t \pmb X_2) ] \\ &= \pmb a^t E [(\pmb X_1 -E(\pmb X_1)] [ \pmb X_2-E(\pmb X_2) ]^t \pmb b \\ &= \pmb a^t \pmb \Sigma_{12} \pmb b \\ \end{aligned} \]
따라서
\[ cor(\pmb a^t \pmb X_1,\pmb b^t \pmb X_2) = \frac{\pmb a^t \pmb \Sigma_{12} \pmb b } { \sqrt{(\pmb a^t \pmb \Sigma_{11} \pmb a ) (\pmb b^t \pmb \Sigma_{22} \pmb b ) }} = (\pmb a^t \pmb \Sigma_{11} \pmb a )^{-1/2} (\pmb a^t \pmb \Sigma_{12} \pmb b)(\pmb b^t \pmb \Sigma_{22} \pmb b )^{-1/2} \]
이제 새로운 두 변수의 상관계수를 최대로 하는 벡터 \(\pmb a\)와 \(\pmb b\)를 찾으면 정준상관계수가 구해진다. 여기서 정준상관계수행렬(canonical correlation matrix) \(\pmb C\)을 다음과 같이 정의하자.
\[ \pmb C = \pmb \Sigma_{11}^{-1/2} \pmb \Sigma_{12} \pmb \Sigma_{22}^{-1/2} \tag{7.7}\]
또한 정준상관계수행렬 \(\pmb C\) 의 Sigular Value Decomposition(SVD)를 다음과 같이 고려하자 (SVD 에 대한 자세한 내용은 섹션 C.4 를 참조)
\[ \pmb C = \pmb U \pmb S \pmb V^t \tag{7.8}\]
여기서 \(\pmb U\)와 \(\pmb V\) 는 각각 차원이 \(p \times p\), \(q \times q\)인 정규직교행렬(orthonormal matrix)이고 \(\pmb S\)는 \(p \times q\) 행렬로 대각원소는 \(\pmb C \pmb C^t\)의 고유값 \(\lambda_i\)의 제곱근을 대각원소로 하고 비대각원소는 0이다.
\[ \pmb U^t \pmb U = \pmb I_p ,\quad \pmb V^t \pmb V = \pmb I_q \] 만약 \(\pmb u_1, \pmb u_2, \dots ,\pmb u_p\)와 \(\pmb v_1, \pmb v_2, \dots ,\pmb v_q\)를 각각 \(\pmb U\)와 \(\pmb V\)의 정규직교벡터라고 놓으면 다음과 같이 표시할 수 있다.
\[ \pmb C = \pmb U \pmb S \pmb V^t = [ \pmb u_1~ \pmb u_2~ \dots ~\pmb u_p] \begin{bmatrix} \sqrt{\lambda_1} & 0 & \dots & 0 & 0 & \dots & 0\\ 0 & \sqrt{\lambda_2} & \dots & 0 & 0 & \dots & 0\\ 0 & 0 & \dots & 0 & 0 & \dots & 0\\ 0 & 0 & \dots & \sqrt{\lambda_p} & 0 & \dots & 0\\ \end{bmatrix} \begin{bmatrix} \pmb v_1^t \\ \pmb v_2^t \\ \vdots \\ \pmb v_q^t \\ \end{bmatrix} \]
여기서 참고로 \(\pmb u_1, \pmb u_2, \dots ,\pmb u_p\)는 \(\pmb C \pmb C^t\)의 고유벡터이고 \(\pmb v_1, \pmb v_2, \dots ,\pmb v_q\)는 \(\pmb C^t \pmb C\)의 고유벡터이다. 또한 두 행렬 \(\pmb C \pmb C^t\)와 \(\pmb C^t \pmb C\)는 0이 아닌 고유값이 \(\lambda_1 \ge \lambda_2 \ge \dots \ge \lambda_p\)로 같다.
\[ \pmb C \pmb C^t = \pmb \Sigma_{11}^{-1/2} \pmb \Sigma_{12} \pmb \Sigma_{22}^{-1} \pmb \Sigma_{12}^t \pmb \Sigma_{11}^{-1/2} \quad \pmb C^t \pmb C = \pmb \Sigma_{22}^{-1/2} \pmb \Sigma_{12}^t \pmb \Sigma_{11}^{-1} \pmb \Sigma_{12} \pmb \Sigma_{22}^{-1/2} \] SVD 분해에 대한 자세한 내용은 섹션 C.4 를 참조하자.
7.3.2 정준상관계수의 유도
다시 정준상관계수의 정의에 따른 벡터 \(\pmb a\)와 \(\pmb b\)를 찾는 문제로 돌아가서 두 벡터 \(\pmb X_1\)과 \(\pmb X_2\)을 표준화한 후에 선형조합을 고려한다.
\[ Z = \pmb a^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1), \quad W = \pmb b^t \pmb \Sigma_{22}^{-1/2}(\pmb X_2-\pmb \mu_2) \]
여기서 상관계수는 상수를 더하거나 빼도 변함이 없으므로 평균벡터를 빼는것은 영향이 없다. 또한 두 벡터 \(\pmb a\)와 \(\pmb b\) 대신에 \(\pmb \Sigma_{11}^{-1/2} \pmb a\)와 \(\pmb \Sigma_{22}^{-1/2} \pmb b\)를 고려해도 어차피 두 개의 임의의 벡터이므로 정준상관계수의 정의에는 영향을 미치지 않는다.
\[ \begin{aligned} \rho(\pmb X_1, \pmb X_2) &= \max_{\pmb a, \pmb b}~~ cor[\pmb a^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1), \pmb b^t \pmb \Sigma_{22}^{-1/2}(\pmb X_2-\pmb \mu_2) ]\\ &= \max_{\pmb a_*, \pmb b_*}~~ cor(\pmb a_*^t \pmb X_1, \pmb b_*^t \pmb X_2) \end{aligned} \]
두 확률변수 \(Z\)와 \(W\)의 공분산과 분산은 다음과 같이 주어지고
\[ \begin{aligned} Cov(Z,W) &= Cov(\pmb a^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1), \pmb b^t \pmb \Sigma_{22}^{-1/2}(\pmb X_2-\pmb \mu_2) ) \\ &= \pmb a^t \pmb \Sigma_{11}^{-1/2} \pmb \Sigma_{12} \pmb \Sigma_{22}^{-1/2} \pmb b \\ &= \pmb a^t \pmb C \pmb b \\ Var(Z) &= \pmb a^t \pmb a \\ Var(W) &= \pmb b^t \pmb b \end{aligned} \]
따라서
\[ cor(Z,W) = \frac{ \pmb a^t \pmb C \pmb b } { \sqrt{(\pmb a^t \pmb a)(\pmb b^t \pmb b)}} = \frac{\pmb a^t} {\sqrt{(\pmb a^t \pmb a)}} \pmb C \frac{\pmb b} {\sqrt{(\pmb b^t \pmb b)}} \]
위의 식에서 확률변수 \(Z\)와 \(W\)의 상관계수는 두 벡터 \(\pmb a/\sqrt{\pmb a^t \pmb a}\)와 \(\pmb b/\sqrt{\pmb b^t \pmb b}\)의 길이가 1이므로 두 벡터 \(\pmb a\)와 \(\pmb b\)의 길이를 1로 가정해도 무방하다 (\(\pmb a^t \pmb a=\pmb b^t \pmb b=1\)). 따라서 다음과 같이 표시할 수 있다.
\[ cor(Z,W) = \pmb a^t \pmb C \pmb b, \quad \pmb a^t \pmb a=\pmb b^t \pmb b=1 \tag{7.9}\]
여기서 벡터공간의 알려진 사실을 이용한다. 두 벡터 \(\pmb a\)와 \(\pmb b\)는 각각 \(p\)와 \(q\)-차원 벡터이므로 다음과 같은 두 벡터 $$와 \(\pmb \beta\)가 존재하고 다음과 같이 표시할 수 있다.
\[ \begin{aligned} \pmb a &= \pmb U \pmb \alpha = \alpha_1 \pmb u_1 + \alpha_2 \pmb u_2 + \dots + \alpha_p \pmb u_p \\ \pmb b &= \pmb V \pmb \beta = \beta_1 \pmb v_1 + \beta_2 \pmb v_2 + \dots + \beta_q \pmb v_q \end{aligned} \tag{7.10}\]
여기서 \(\pmb \alpha^t \pmb \alpha = \pmb \beta^t \pmb \beta =1\).
이제 식 7.9 에 식 7.10 를 대입하고 식 식 7.8 의 SVD를 이용하면
\[ \begin{aligned} cor(Z,W) &= \pmb a^t \pmb C \pmb b \\ &= \pmb \alpha^t \pmb U^t \pmb U \pmb S \pmb V^t \pmb V \pmb \beta \\ &= \pmb \alpha^t \pmb S \pmb \beta \\ &= \sum_{i=1}^ p \alpha_i \beta_i \sqrt{\lambda_i} \\ &\le (\max_i \sqrt{\lambda_i} ) \sum_{i=1}^ p \alpha_i \beta_i \\ & = \sqrt{\lambda_1} \sum_{i=1}^ p \alpha_i \beta_i \end{aligned} \]
위의 식에서 \(cor(Z,W)\)가 상한값(upper bound)와 같아지려면 \(\pmb \alpha\)와 \(\pmb \beta\)의 계수는 다음과 같은 조건일 때이다.
\[ \alpha_1=1,~ \alpha_2=\dots=\alpha_p=0, \quad \beta_1=1,~ \beta_2=\dots=\beta_p=0 \]
따라서 표준화된 두 벡터의 정준상관계수는 \(\pmb C \pmb C^t\)의 최대고유값의 제곱근이 되고 선형조합을 만드는 두 벡터 \(\pmb a\) 와 \(\pmb b\)는 각각 \(\pmb C \pmb C^t\) 와 \(\pmb C^t \pmb C\)의 첫번쨰 고유벡터가 된다.
\[ \begin{aligned} \rho(\pmb X_1,\pmb X_2) &= Cov(Z,W) \\ &= \max_{\pmb a, \pmb b} ~~cor(\pmb a^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1) , \pmb b^t \pmb \Sigma_{22}^{-1/2}(\pmb X_2-\pmb \mu_2)) \\ &= cor(\pmb u_1^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1) , \pmb v_1^t \pmb \Sigma_{22}^{-1/2}(\pmb X_2-\pmb \mu_2)) \\ &= \sqrt{\lambda_1} \end{aligned} \]
여기서 \(\pmb a =\pmb u_1, ~ \pmb b= \pmb v_1\)이다.
정준상관계수를 정의하고 유도하는 과정을 보면 주성분분석과 유사하게 서로 공분산이 0인 확률변수들를 만들고 대응하는 두 변수들간의 상관계수를 큰 순서대로 만들 수 있다. 예를 들어 \(Z_2= \pmb u_2^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1)\)는 \(Z_1=\pmb u_1^t \pmb \Sigma_{11}^{-1/2} (\pmb X_1 - \pmb \mu_1)\)과 공분산이 0이고 또한 \(W_2=\pmb v_2^t \pmb \Sigma_{22}^{-1/2} (\pmb X_2 - \pmb \mu_2)\)는 \(W_1=\pmb v_1^t \pmb \Sigma_{22}^{-1/2} (\pmb X_2 - \pmb \mu_2)\)과 공분산이 0 이다 (정규분포인경우 독립이다). 더 나아가 \(Z_2\)와 \(W_2\)의 상관계수는 각각 \(Z_1\)과 \(W_1\)의 공분산이 0인 변수들중에서 최대의 상관계수 \(\sqrt{\lambda_2}\) 를 가지게 된다. 유사한 방법으로 \(Z_i\)와 \(W_i\)들을 정의할 수 있다.
7.3.3 표본 정준상관계수
이제 모집단이 아닌 표본이 주어졌을 경우 표본 정준상관계수를 구하는 방법은 다음과 같다. 이번에도 \(p \le q\)를 가정하자.
확률벡터 \(\pmb X\)의 표본 상관계수행렬을 두 벡터 \(\pmb X_1\), \(\pmb X_2\)의 차원에 맞게 다음과 같이 표시하고 구한다.
\[ \pmb R = \begin{bmatrix} \pmb R_{11} & \pmb R_{12} \\ \pmb R_{12}^t & \pmb R_{22} \end{bmatrix} \]
다음과 같은 두 행렬을 구한다.
\[ \pmb E_1 = \pmb R_{11}^{-1} \pmb R_{12} \pmb R_{22}^{-1} \pmb R_{12}^t \]
\[ \pmb E_2 = \pmb R_{22}^{-1} \pmb R_{12}^t \pmb R_{11}^{-1} \pmb R_{12} \]
두 벡터 \(\pmb X_1\), \(\pmb X_2\)의 표본 정준상관계수는 \(\pmb E_1\) 또는 \(\pmb E_2\)의 최대 고유치의 제곱근이다. 참고로 \(\pmb E_1\) 와 \(\pmb E_2\)은 이 아닌 고유값이 같다.
swiss 자료에서 출산율(Fertility)과 영아사망율(Infant.Mortality)을 하나의 벡터로 나머지 4개의 사회경제변수를 다른 벡터로 하여 표본 정준상관계수는 다음과 같이 구할 수 있다.
n <-dim(swiss)[1]
p <- 2
q <- 4
swiss0 <- swiss[,c(1,6,2,3,4,5)] # 순서를 바꾸는 작업
R <- cor(swiss0)
R11 <- matrix(R[1:p,1:p],p,p)
R22 <- matrix(R[(p+1):(p+q),(p+1):(p+q)],q,q)
R12 <- matrix(R[1:p,(p+1):(p+q)],p,q)
R11 [,1] [,2]
[1,] 1.000000 0.416556
[2,] 0.416556 1.000000R22 [,1] [,2] [,3] [,4]
[1,] 1.0000000 -0.6865422 -0.6395225 0.4010951
[2,] -0.6865422 1.0000000 0.6984153 -0.5727418
[3,] -0.6395225 0.6984153 1.0000000 -0.1538589
[4,] 0.4010951 -0.5727418 -0.1538589 1.0000000R12 [,1] [,2] [,3] [,4]
[1,] 0.35307918 -0.6458827 -0.66378886 0.4636847
[2,] -0.06085861 -0.1140216 -0.09932185 0.1754959E1 <- solve(R11) %*% R12 %*% solve(R22) %*% t(R12)
E2 <- solve(R22) %*% t(R12) %*% solve(R11) %*% R12
E1.eigen <- eigen(E1)
E2.eigen <- eigen(E2)
rho <- sqrt(E1.eigen$value[1])
rho # sample CCA[1] 0.8142291sqrt(E2.eigen$value[1])[1] 0.8142291+0iR 패키지 CCA의 함수 cc 를 이용하면 표본 정준상관계수와 선형변환을 위한 벡터를 구할 수 있다.
library(CCA)
X1 <- swiss[,c(1,6)]
X2 <- swiss[,c(2,3,4,5)]
res1 <- cc(X1,X2)
res1$cor # sample CCA[1] 0.8142291 0.2222637res1$xcoef # vectors for linear transformation for X1 [,1] [,2]
Fertility -0.08456464 -0.02455185
Infant.Mortality 0.05534823 0.37357101res1$ycoef # vectors for linear transformation for X2 [,1] [,2]
Agriculture 0.01985292 -0.05136137
Examination 0.02690124 0.02476760
Education 0.09414900 -0.03527373
Catholic -0.01164052 0.01793981