[R] 자기상관계수 (Autocorrelation Coefficients), 자기상관그림(Autocorrelation Plot)
R 분석과 프로그래밍/R 통계분석 2020. 3. 22. 23:58자기상관계수(Autocorrelation coefficients), 자기상관그림(Autocorrelation Plot)은 (a) 데이터가 시간에 의존하는 것 없이 무작위성 (randomness, no time dependence in the data) 을 띠는지 확인하는데 사용됩니다. 그리고 (b) 시계열분석에서 Box-Jenkins ARIMA (Autoregressive Integrated Moving Average) 모형 식별 단계(model identification stage)에서 AR 차수를 구하는데도 사용됩니다.
이번 포스팅에서는 R을 사용해서
(1) 자기상관계수(Autocorrelation coefficients)를 계산하고,
(2) 자기상관계수의 95% 신뢰구간(95% confidence interval)을 계산하고,
(3) 자기상관그림(Autocorrelation Plot)을 그리는 방법
을 소개하겠습니다.
본론에 들어가기 전에 먼저, 상관계수(correlation coefficients)와 자기상관계수(autocorrelation coefficients)를 비교해보겠습니다.
상관계수와 자기상관계수의 유사한 점은 두 연속형 변수 간의 관계를 분석하다는 점입니다. 기본 개념은 공분산을 표준편차로 나누어서 표준화해주며, (자기)상관계수가 -1 ~ 1 사이의 값을 가지게 됩니다.
상관계수와 자기상관계수가 서로 다른 점은, 상관계수는 특정 동일 시점을 횡단면으로 해서 Y와 다른 X1, X2, ... 변수들 간의 관계를 분석합니다. 반면에 자기상관계수는 동일한 변수(Yt, Yt-1, Yt-2, ...)의 서로 다른 시간 차이 (time lag) 를 두고 관계를 분석하는 것입니다.
기존에 cross-sectional 관점의 Y와 X 변수들 간의 상관관계 분석에 많이 익숙해져 있는 경우에, 시계열 분석을 공부할 때 보면 자기 자신의 시간 차이에 따른 자기상관관계 관점이 처음엔 헷갈리고 잘 이해가 안가기도 합니다. 관점이 바뀐것일 뿐 어려운 개념은 아니므로 정확하게 이해하고 가면 좋겠습니다.
(1) 자기상관계수(Autocorrelation coefficients) 구하고, 자기상관그림 그리기 |
자기상관계수는 아래의 공식과 같이 로 계산합니다.
간단한 예로서, 1 ~ 50까지의 시간 t에 대해서 싸인곡선 형태의 주기적인 파동을 띠는 값에 정규확률분포 N(0, 1) 에서 추출한 난수를 더하여 생성한 Y 데이터셋에 대해서 R로 위의 공식을 이용해서 각 시차(time lag)별로 자기상관계수를 구해보겠습니다.
# sine wave with noise set.seed(123) noise <- 0.5*rnorm(50, mean=0, sd=1) t <- seq(1, 50, 1) y <- sin(t) + noise plot(y, type="b") |
자기상관계수를 구하고 자기상관그림을 그리는 가장 간단한 방법은 R stats 패키지에 들어있는 acf() 함수를 사용하는 것입니다.
z <- acf(y, type=c("correlation"), plot=TRUE) # autocorrelation coefficients per time lag k round(z$acf, 4) [,1] [1,] 1.0000 [2,] 0.3753 [3,] -0.3244 [4,] -0.5600 [5,] -0.4216 [6,] 0.2267 [7,] 0.5729 [8,] 0.3551 [9,] -0.0952 [10,] -0.5586 [11,] -0.4582 [12,] 0.0586 [13,] 0.3776 [14,] 0.4330 [15,] 0.0961 [16,] -0.4588 [17,] -0.4759 # autocorrelation plot |
좀 복잡하기는 합니다만, 위의 자기상관계수를 구하는 공식을 이용해서 R로 직접 사용자 정의 함수를 만들고, for loop 반복문을 이용해서 시차(time lag) k 별로 자기상관계수를 구할 수도 있습니다. 자기공분산을 분산으로 나누어서 표준화해주었습니다. 시차 0(time lag 0) 일 경우에는 자기자신과의 상관관계이므로 자기상관계수는 '1'이 되며, 표준화를 해주었기 때문에 자기상관계수는 -1 <= autocorr(Y, Yt-k) <= 1 사이의 값을 가집니다.
# User Defined Function of Autocorrelation acf_func <- function(y, lag_k){ # y: input vector # lag_k : Lag order of k
N = length(y) # total number of observations y_bar = mean(y)
# Variance var = sum((y - y_bar)^2) / N
# Autocovariance auto_cov = sum((y[1:(N-lag_k)] - y_bar) * (y[(1+lag_k):(N)] - y_bar)) / N
# Autocorrelation coefficient = Autocovariance / Variance r_k = auto_cov / var
return(r_k) } # Compute Autocorrelation per lag (from 0 to 9 in this case) acf <- data.frame() for (k in 0:(length(y)-1)){ acf_k <- round(acf_func(y, lag_k = k), 4) acf[k+1, 'lag'] = k acf[k+1, 'ACF'] = acf_k } > print(acf) lag ACF 1 0 1.0000 2 1 0.3753 3 2 -0.3244 4 3 -0.5600 5 4 -0.4216 6 5 0.2267 7 6 0.5729 8 7 0.3551 9 8 -0.0952 10 9 -0.5586 11 10 -0.4582 12 11 0.0586 13 12 0.3776 14 13 0.4330 15 14 0.0961 16 15 -0.4588 17 16 -0.4759 18 17 -0.0594 19 18 0.2982 20 19 0.3968 21 20 0.1227 22 21 -0.2544 23 22 -0.3902 24 23 -0.1981 25 24 0.1891 26 25 0.3701 27 26 0.1360 28 27 -0.1339 29 28 -0.2167 30 29 -0.1641 31 30 0.0842 32 31 0.2594 33 32 0.1837 34 33 0.0139 35 34 -0.1695 36 35 -0.1782 37 36 -0.0363 38 37 0.1157 39 38 0.1754 40 39 0.0201 41 40 -0.1428 42 41 -0.1014 43 42 0.0000 44 43 0.0750 45 44 0.0641 46 45 -0.0031 47 46 -0.0354 48 47 -0.0405 49 48 -0.0177 50 49 -0.0054
|
이번 예의 경우 시간의 흐름에 따라 싸인 파동 (sine wave) 형태를 보이는 데이터이므로, 당연하게도 시차를 두고 주기적으로 큰 양(+)의 자기상관, 큰 음(-)의 자기 상관을 보이고 있습니다.
(2) 자기상관계수의 95% 신뢰구간(95% confidence interval) 계산 |
자기상관계수의 95% 신뢰구간은 아래의 공식으로 간단하게 구할 수 있습니다.
자기상관계수의 95% 신뢰구간 =
이때 N은 샘플 개수, z 는 표준정규분포의 누적분포함수, 는 유의수준 5% 입니다.
이미 익숙하게 알고 있다시피 유의수준 =0.05 에서의 z = 1.96 입니다. 이번 예에서의 관측치 개수 N = 50 입니다. 따라서 자기상관계수의 95% 신뢰구간은 이 됩니다.
> # z quantile of 95% confidence level > qnorm(0.975, mean=0, sd=1, lower.tail=TRUE) [1] 1.959964 > qnorm(0.025, mean=0, sd=1, lower.tail=TRUE) [1] -1.959964 > > # sample size > N <- length(y) > print(N) [1] 50 > > # 95% confidence interval of autocorrelation > qnorm(0.975, mean=0, sd=1, lower.tail=TRUE)/ sqrt(N) [1] 0.2771808 > qnorm(0.025, mean=0, sd=1, lower.tail=TRUE)/ sqrt(N) [1] -0.2771808 |
위의 acf() 함수로 그린 자기상관그림(autocorrelation plot)에서 ACF 축 +0.277 과 -0.277 위치에 수평으로 그어진 점선이 바로 95% 신뢰구간의 상, 하한 선입니다.
위의 자기상관그림에서 보면 95% 신뢰구간의 상, 하한 점선을 주기적으로 벗어나므로 이 Y 데이터셋은 무작위적이지 않으며 (not random), 시계열적인 자기상관(autocorrelation)이 존재한다고 판단할 수 있습니다.
만약, 무작위적인 데이터셋에 대해서 자기상관계수를 계산하고 자기상관그림을 그린다면, 시차(time lag)에 따른 자기상관계수가 모두 95% 신뢰구간 상, 하한 점선의 안쪽에 위치할 것입니다.
많은 도움이 되었기를 바랍니다.
이번 포스팅이 도움이 되었다면 아래의 '공감~'를 꾹 눌러주세요. :-)
'R 분석과 프로그래밍 > R 통계분석' 카테고리의 다른 글
[R] 빈도 데이터 분석을 위한 포아송 회귀모형과 대안 회귀모형에 대한 소개 (9) | 2019.10.23 |
---|---|
[R] 로지스틱 회귀분석을 통한 유방암 예측(분류) (4/4): 로지스틱 회귀모형 적합 및 모델평가, 해석 (19) | 2018.12.18 |
[R] 로지스틱 회귀분석을 통한 유방암 예측(분류) (3/4): 1차 변수 선택 및 목표변수와 설명변수 간 관계 분석 (8) | 2018.12.17 |
[R] 로지스틱 회귀분석을 통한 유방암 예측(분류) (2/4): 탐색적 데이터 분석 및 전처리 (32) | 2018.12.16 |
[R] 로지스틱 회귀분석을 통한 유방암 예측(분류) (1/4): WDBC Data 소개 (4) | 2018.12.16 |