시간에 대한 신호를 주파수 영역으로 변환하여 분석하기 위해 통상 푸리에 변환을 많이 사용하여 왔다(Cohen, 1995; Hlawatsch and Auger, 2013). 푸리에 변환은 사인 함수 및 코사인 함수를 기저 함수로 사용하기 때문에, 국소적인 변화보다는 전체적인 주기성을 분석하는 경우 유용하게 사용된다. 그러나 시계열 자료가 비정상적인 경우 주어진 자료 기간 동안 동일한 주기특성을 가지는 것을 기대하기는 어렵다. 사실 이러한 문제는 자연에서 관측되는 대부분의 시계열이 가지고 있는 문제이기도 하다. 이런 문제를 극복하고자 시간-주파수 분석이 태동되었다(Wigner, 1932; Smith et al., 1998). 푸리에 변환에 기반한 시간-주파수 분석은 매 시점 사인 또는 코사인 함수의 상대적인 비중이 어떻게 달라지는지를 보여준다.

푸리에 변환에 기반한 시간-주파수 분석의 유용성도 크지만 한계도 있다. 가장 기본적인 것이 이 기저 함수들이 무한 반복된다는 점이다. 즉, 그 영향이 시계열 전체에 미치는 것이다. 이런 단점을 극복하고자 제안된 것이 웨이블릿 분석이다(Gabor, 1946; Misiti et al., 2007). 작은 파동이라는 의미의 웨이블릿은 유한 길이를 갖는다. 따라서 웨이블릿 변환은 유한한 길이의 기저 함수(basis function)를 이용하여 기존의 시계열 자료를 시간-주파수 영역으로 변환한다. 따라서 푸리에 변환과 달리 웨이블릿 변환에서는 국부적인 변화를 더욱 생생하게 나타낼 수 있다. 따라서, 신호의 불연속적인 지점을 식별하는 데에 유용하다.

웨이블릿 변환의 기저 함수를 모 웨이블릿(mother wavelet)이라 하며, 모 웨이블릿의 척도 변환(scaling) 및 전이(transition)를 통해 시간 도메인에 있는 신호가 주파수 영역으로 변환된다. Eq. (1)은 연속 웨이블릿 변환(continuous wavelet transform, CWT)의 기본적인 공식이다.

여기서, W_x(a,b)는 a와b를 매개변수로 가지는 연속 웨이블릿 변환을 의미한다. x(t)는 주어진 시계열을 의미한다. ψ는 모 웨이블릿 함수를 나타내며, 위첨자*는 복소 공액(complex conjugate)를 의미한다. a는 척도 매개변수(scale parameter)로 모 웨이블릿의 크기를 조정하며, b는 전이 매개변수(transition parameter)로 모 웨이블릿을 시간축 상 위치를 이동시킨다. 위 식으로 산정되는W_x(a,b) 은 통상 복소수의 값을 가지기 때문에 크기 비교를 위해| W_x(a,b)|²로 계산되는 wavelet power (WP)가 주로 사용된다.

웨이블릿 분석을 위해서는 모 웨이블릿이 분석의 목적에 맞게 선택되어야 한다. 다양한 종류의 모 웨이블릿이 존재하며, 동일한 자료를 분석하더라도 모 웨이블릿의 선택에 따라 그 결과가 확연하게 달라진다(Ngui et al., 2013). 단변량 모 웨이블릿 분석 시, Morlet 모 웨이블릿이 타 모 웨이블릿에 비해 시간과 주파수 영역에서의 조화가 좋다고 알려져 있다(Aguiar-Conraria et al., 2008; Rua, 2010; Lee et al., 2019). 뿐만 아니라, Lee et al. (2019)에서는 이변량 웨이블릿 분석에서도 Morlet 모 웨이블릿이 타 모 웨이블릿에 비해 많은 장점을 가진다는 것을 언급한 바 있다. 이에 본 연구에서도 Morlet 모 웨이블릿을 이용하여 분석을 진행하였다. Morlet 모 웨이블릿은 Eq. (2)로 표현된다.

여기서, η는 척도 및 전이 매개변수에 의해 변환된 무차원 시간 매개변수(non-dimensional time parameter)를 나타낸다. 아울러 w₀ 은 파수(wave number)를 의미한다. 이는 무차원 주파수(nondimensional frequency)라고 불리기도 한다. 본 연구에서는 해당 값으로 6을 적용하였다. 이러한 경우 Morlet 모 웨이블릿의 형태는 6번의 진동(oscillation)을 가지게 된다. 아울러, w₀ 으로 6을 사용하는 경우가 웨이블릿 변환 시 시간-주파수 영역의 지역화(localization) 측면에서 최적의 균형을 준다고 알려져 있다(Grinsted et al., 2004; Aguiar-Conraria et al., 2008; Rua, 2010).

모 웨이블릿이 결정되면, 분석 대상이 되는 주파수가 설정되어야 한다. 통상 이를 위해 주파수의 역수인 스케일이 활용된다. 이 스케일은 Eq. (3)으로 결정되며, 가장 큰 스케일의 index를 의미하는J는 Eq. (4)를 통해 얻어진다.

여기서, s_j는j 에서의 스케일 크기이며, j 는 스케일의 index를 의미한다. 이는 0부터 J 까지 존재한다. 따라서 분석 대상 스케일은 총 J+1 개가 된다. s₀ 는 분석 대상 스케일 중 가장 작은 값이며, δ_s 는 스케일의 조밀함을 결정하는 값이다. δ_s 는 상수이며 이를 작게 설정할수록 웨이블릿 스펙트럼을 보다 정밀하게 표현할 수 있다. δ_s 는 일반적으로 0.4875가 사용된다(Misiti et al., 2007). 아울러, δ_t 는 웨이블릿 분석을 위한 시간 단위를 의미한다. 그리고fix[] 는 숫자의 정수 부분만을 반환하는 함수이다. 아울러, 작은 스케일(s₀)은 일반적으로 2δ_t 로 설정된다.

웨이블릿 변환은 시간 영역의 자료를 시간-주파수 영역으로 변환하기 때문에, 그 해석이 푸리에 변환처럼 단순하지 않다. 특히, 웨이블릿 스펙트럼이 시간-주파수 영역에서 보여주는 주파수가 원자료의 주파수를 의미하는 것은 아니라는 점에 주의할 필요가 있다. 주파수와 역수 관계를 가지는 스케일 또한 마찬가지이다. 이러한 현상이 나타나는 이유는 일정한 주파수를 가지는 사인 또는 코사인 함수를 기저함수로 사용하는 푸리에 변환과는 달리, 웨이블릿 분석에서는 기저함수로 사용되는 모 웨이블릿이 가지는 중심 주파수(center frequency)가 모두 다르기 때문이다. 다행히 Morlet 모 웨이블릿(m=6)을 사용하는 경우에는 웨이블릿의 분석 결과에 나타난 주파수가 원 자료의 주파수와 거의 유사한 값(λ_s=1.03s)을 가지게 된다. 따라서, 본 연구에서는 이후의 설명에서 스케일과 주기를 동일시하였다. 추가로, 웨이블릿 분석 결과에서는 통상 가장자리 효과(edge effect)가 무시될 수 없는 영역을 COI (Cone of Influence)로 나타내 준다. 이 영역에서는 불확실성이 상대적으로 크기 때문에 결과의 해석에 주의할 필요가 있다. Morlet 모 웨이블릿을 사용하는 경우, 스케일 s 에서의 COI는 √2 s 로 계산된다(Aguiar-Conraria et al., 2008).

웨이블릿 변환은 주어진 함수(모 웨이블릿)에 대한 대역통과필터(bandpass filter)의 개념이므로, 역회선적분(deconvolution) 또는 역필터(inverse filter) 개념을 이용하면 원 시계열 자료를 복원(reconstruction)할 수 있다. 이를 웨이블릿 역변환(inverse wavelet transform)이라고 하는데, 연속형 웨이블릿 변환의 경우 위 과정이 매우 복잡하여 통상 모 웨이블릿과 완전히 다른 함수를 도입하여 역변환을 수행한다(Georgiou and Cohen, 2001). 가장 쉬운 방법은 Dirac delta 함수를 도입하는 것이다. 원자료가 실수로 이루어진 경우, Eqs. (5) and (6)을 활용하여 웨이블릿 스펙트럼의 실수 부분을 통해 원자료 시계열에 근접한 값으로 복원할 수 있다(Farge, 1992).

여기서, x^’(t)는 복원된 시계열을 의미한다. δ_s, δ_s 그리고s_j는 앞서 정의된 것과 같다. Ψ₀(o)는 원점에서의 모 웨이블릿 함수 값을 의미하며, 이는 스케일의 영향을 제거하기 위해 사용된다. W_x(s_j,t) 는 스케일 s_j , 시간 t에서의 연속 웨이블릿 변환 값을 의미하며, R[] 은 웨이블릿 스펙트럼의 실수 부분을 반환하는 역할을 한다. 또한 C_δ 는 모 웨이블릿 함수에 따라 달라지는 상수이며, Eq. (6)과 같이 Dirac delta 함수를 이용하여 계산된다. 이때 W_δ 는 Dirac delta 함수에 대한 연속 웨이블릿 변환을 의미한다. 이때, Morlet 모 웨이블릿(ω₀=6)의 경우, C_δ 값은 0.776으로 계산된다.

본 연구에서는 웨이블릿 스펙트럼의 확장을 위해 block resampling 방법을 사용하였다. Block resampling은 Bullmore et al. (2001)에 의해 제안된 ‘random (free) permutation’ 방법을 개선하기 위해 제안된 방법이다. Random permutation은 단순하게 동일한 스케일 내의 스펙트럼을 임의로 섞는 방법을 의미한다. 이러한 경우 결과적으로 스펙트럼이 가지는 선형적인 관계는 모두 제거된다. 이와 같은 한계점을 보완하기 위해 Breakspear et al. (2003)은 일련의 블록 크기로 스펙트럼을 나눈 뒤 이를 재추출하는 방식인 block resampling을 제안하였다.

Block resampling 방법의 적용을 위해서는 스펙트럼 재추출을 위한 시간 단위를 스케일별로 설정할 필요가 있다. 이러한 시간 단위를 블록 길이(block length)라고 한다(Breakspear et al., 2003). Breakspear et al. (2003)은 동일한 크기의 블록을 모든 스케일에 적용하였다. 본 연구에서는 이를 개선한 기법인 (1) 각 스케일에 따라 다른 블록 길이를 고려하는 방법(block resampling with scale-dependent block length)과 Breakspear et al. (2003)과 동일한 기법인 (2) 모든 스케일에 동일한 블록 길이를 설정하는 방법(block resampling with constant block length)을 모두 고려하였다. 편의상 각 방법을 Method (1), Method (2)라고 칭하고자 한다. 추가로, 모든 스케일에 대해 같은 블록이 선정되게 하는 방법(block resampling with same block length at the same time)을 고려하였다. 이를 Method (3)이라 부르고자 한다. Method (3)은 동일한 시점의 스케일 간 상관성을 보존하는 경우를 고려하기 위한 방법으로 Method (2)와는 구분된다. Fig. 2는 각 방법에 대한 개념도를 보여준다.

Method (1)과 같이 스케일별로 다른 크기의 블록 길이를 설정하는 경우, 각 스케일별 블록 길이를 결정해야 한다. 이는 시간-주파수 영역에서의 자기상관성과 교차상관성을 어느 정도 유지시킬 것이냐의 문제와도 직결된다. 이러한 고려 없이 해당 스케일의 전⋅후 스케일 평균값을 사용한 사례가 있기도 하나(Erkyihun et al., 2016), 이는 본 연구의 목적과는 부합하지 않는다. Breakspear et al. (2003)도 블록 길이는 원자료 스펙트럼의 상관성(correlation)에 대한 특성 길이(characteristic length), 즉, 상관거리(correlation length) 이상으로 설정하는 것을 권장한 바 있다. 시간축에서 분석이라면 상관거리는 상관시간(correlation time)으로 부를 수도 있다. 본 연구에서도 불록 길이의 결정에 상관거리의 개념을 고려하였다. 상관거리는 통상 자기상관함수(auto-correlation function)의 분석을 통해 추정된다. 간단히 상관거리는 자기상관함수가 포괄하는 면적으로 정의된다. 자기상관이 양(+)과 음(-)을 반복하는 경우에는 0으로 수렴하는 함수로 근사화하여 상관거리를 추정한다. 만일 AR (1) 모형의 특성을 보인다면 그 면적은 상관계수e^-1 에 해당하는 길이와 같다. 이런 이유로 상관거리를 e-folding distance라고 부르기도 한다(Rodríguez-Iturbe and Mejía, 1974).

본 연구에서는 앞서 설명된 방법들을 적용해보고 평가하기 위해 비정상 시계열을 모의 생성하여 이용하였다. 모의자료를 생성을 위해 비선형 시계열 모형을 활용하였다. 비선형 모형은 autoregressive conditional Herteroskedasticity (ARCH) 모형(Engle, 1982), exponetial autoregressive (EXPAR) 모형(Tong, 1990), threshold autoregressive (TAR) 모형(Tong, 1990), smooth transition autoregressive (STAR) 모형(Terasvirta and Anderson, 1992) 등 다양한 형태가 제시된 바 있다. 본 연구에서는 이러한 모형들 중 특정한 추세를 가지지 않으며 비교적 간단한 식을 통해 생성 가능한 NMA (Nonlinear moving average) 모형과 NAR (Nonlinear Autoregressive) 모형(Zhang et al., 2001)을 활용하였다. 이들 두 모형은 각각 Eqs. (7) and (8)로 표현된다.

여기서, x_t 는 시간 t에 대해 생성되는 값이며, ∈ 는 평균 0, 표준편차 1인 백색잡음이다.

본 연구에서는 이 모형들로 모의된 비정상 시계열을 각각 NMA, NAR로 명명하였다. 모의 시계열의 시간 단위는 월(month)로 설정하였으며, 자료의 길이는 672개월(56년)로 가정하였다. 최종적으로 모의된 시계열은 Fig. 3과 같다. 생성된 시계열은 시각적으로 어떠한 주기성이나 경향성이 확인되지 않는다. NMA와 NAR의 평균은 각각 -0.222, 0.243이며, 표준편차는 1.186, 0.996이다.

웨이블릿 변환을 위한 주요 변수들의 설정 과정은 다음과 같다. 먼저, 본 연구에서는 Eqs. (3) and (4)를 참고하여 분석 대상 스케일을 설정하였다. 웨이블릿 분석을 위한 시간 단위(δ_t)는 자료의 시간 단위인 1개월(1/12년)로 설정하였다. 따라서 최소 스케일 s₀ 는 2δ_t 인 2개월(1/6년)이 된다. Δ_s 는 Misiti et al. (2007)를 참고하여 일반적으로 사용되는 값인 0.4875를 적용하였다. 이러한 경우, J는 17로 계산되며, 이를 통해 산정되는 최대 스케일은 625개월(약 52년)이 된다. 하지만 이때 √2 s 로 계산되는 COI가 가장 큰 스케일에서 884개월(약 74년)이 된다. 따라서 COI가 전체 자료의 길이를 초과하게 되어, 해당 스케일의 결과가 큰 의미를 가지지 못하게 된다. 이에 본 연구에서는 최대 스케일 index를 16까지만 고려하였다. 이러한 경우 스케일은 최소(s₀) 2개월(1/12년)부터 최대(s₁₆) 446개월(약 37년)까지 총 17개로 설정된다.

Fig. 4는 모의자료에 대한 단변량 웨이블릿 분석 결과를 보여준다. 먼저, NMA의 결과(Fig. 4(a))를 살펴보면, 대체로 짧은 주기(약 10년 이하)에서 강한 스펙트럼이 나타난다. 약 2년 이하 주기에서는 강한 스펙트럼이 강약을 반복하며, 약 6년 주기 성분은 중간 시점에서 시작되어 점점 주기가 짧아진다. 해당 스펙트럼은 약 30년 이상의 긴 시간에 걸쳐 있다. 아울러 GWP를 통해 NMA는 장주기보다는 단주기적 성향이 상대적으로 강하게 나타남을 알 수 있다. 한편 NAR의 결과(Fig. 4(b))를 보면, 보다 다양한 시간대 및 주기에서 강한 스펙트럼이 나타난다. 약 20년 주기에 나타나는 장주기 성분이 가장 크게 두드러진다. 이외에도 약 4년 주기 성분이 강하게 나타나는 것 또한 확인할 수 있다. 종합적으로 NAR은 NMA보다는 상대적으로 장주기적 경향이 강하다고 판단된다.

스펙트럼 재추출(resample)을 위한 블록 길이 결정을 위해, 스케일별로 시차에 따른 상관거리를 산정하였다. Fig. 5는 특정 스케일(j = 8)에 대한 상관거리 산정 예시를 보여준다. 참고로, 해당 스케일은 29.9개월(약 2.5년)에 해당한다. 자기상관함수의 면적을 통해 NMA, NAR에 대해 산정된 상관거리는 각각 35.9, 39.8개월이다. 이는 약 3.0, 3.3년에 해당한다. 예시로 산정된 스케일의 경우 두 시계열의 상관거리는 매우 유사하게 나타난다.

Table 1은 동일한 방식을 이용하여 각 스케일별로 산정된 상관거리를 보여준다. 전반적으로 스케일이 커짐에 따라 상관거리가 길어지는 것을 확인할 수 있다. 일정 스케일(약 10년) 이상에서 어느 정도 일정해진다. 이와 같이 산정된 상관거리가 Method (1)의 블록 길이로 고려되었다. 아울러 Method (2) 및 Method (3)의 스펙트럼 재추출(resample)을 위한 블록 길이는 스케일별 블록 길이 중 최댓값을 고려하였다. 이에 상응하는 NMA, NAR의 블록 길이는 각각 114.3개월, 114.0개월이다. 이는 연으로 환산하면 약 9.5년으로 두 시계열 모두 매우 유사한 규모의 상관거리를 가지는 것을 알 수 있다.

앞서 고려한 세 가지 방법을 적용하여 각 모의자료에 대한 웨이블릿 스펙트럼의 확장을 수행하였다. 확장된 웨이블릿 스펙트럼을 역변환하여 가상의 기간에 대한 시계열을 생성하였다. 자료 확장을 위한 기간은 원자료의 5배가 고려되었다. 즉, 원자료의 길이가 기간 T라고 하였을 때, 5T 의 길이를 가지는 확장 스펙트럼 및 확장 시계열이 생성되었다.

Fig. 6은 NMA에 대해 생성된 스펙트럼과 시계열을 보여준다. 방법에 따라 나타나는 스펙트럼의 차이를 살펴보면, Method (1)에서는 작은 스케일에 존재하는 스펙트럼이 단절되어 나타나는 경향이 관찰된다. 한편 Method (2)와 (3)은 매우 유사한 양상을 보여준다. 역변환된 자료를 보면 세 방법 모두 원자료와 유사한 형태로 확장된 것을 확인할 수 있다. 다만, Method (1)과 Method (2)는 Method (3)에 비해 변동 범위가 다소 작게 나타난다. 이는 스펙트럼 재추출 시 스케일별로 다른 시간대의 값들이 추출됨에 따라 스케일 간 상관성이 제거되어 나타나는 현상이라 판단된다.

Fig. 7은 NAR에 대해 생성된 스펙트럼과 시계열을 보여준다. 앞선 결과와 유사하게 나타나나, Method (2)와 Method (3)의 스펙트럼에서 또 다른 차이가 발견된다. 즉, Method (2)의 경우 스케일 간의 상관은 고려하지 않고 임의로 추출함에 따라 스펙트럼의 형태가 매우 다양하게 나타나는 반면, Method (3)의 경우 동시간대의 스펙트럼만을 추출함에 따라 확장 스펙트럼의 부분들이 원자료 스펙트럼과 매우 닮아 있다. 또한 역변환된 자료의 경우 NMA와 유사하게 Method (1), (2)의 경우 극단적인 값들이 다소 완화되어 나타난다.

Fig. 8은 방법별로 확장된 시계열 자료와 원자료의 box-plot을 비교한 것이다. 고려된 방법들 중 원자료와 가장 유사한 분포를 보여주는 것은 Method (3)이다. Method (1)과 (2)의 경우 전반적으로 자료가 중앙값 부근으로 다소 집중되어 나타난다.

추가적으로 본 연구에서는 웨이블릿 이변량 분석을 통한 검토를 수행하였다. 이변량 분석 결과를 검토한 이유는 모의되는 자료가 원자료들이 가지고 있는 스케일 간의 간섭 또는 소위 동조화(co-movement)라고 불리는 특성들을 그대로 유지하고 있는지를 확인하기 위함이다. 이와 관련된 여러 개의 지표가 존재하지만, 통상 cross wavelet power (XWP)와 wavelet coherence (WC)가 주로 사용된다. Cross wavelet transform은 두 시계열의 continuous wavelet transform 결과 중 하나를 복소 공액(complex conjugate) 형태로 변환하여 이를 서로 곱함으로써 계산된다. 크기 비교를 위해 cross wavelet transform에 절댓값을 취한 것을 cross wavelet power (XWP)라 한다. 한편, WC (wavelet coherence)는 XWP를 각 스펙트럼의 크기로 정규화한 지표라고 볼 수 있다. XWP는 범위를 제한하지 않기 때문에 어떤 특정 값이 너무 큰 경우에는 두 시계열이 동시에 발생하는 중요한 동조화를 놓칠 수 있을 가능성이 있는 반면, WC는 두 시계열의 스펙트럼 모두가 작은 경우에도 이를 포착할 수 있다(Grinsted et al., 2004). 즉, WC는 XWP의 보완적 성격을 가진다. WC는 두 스펙트럼이 가지는 시간-주파수 도메인에서의 국부적인 상관도와 유사한 개념으로 볼 수 있으므로, 이를 통해 두 시계열의 동시 발생 영역을 대략적으로 파악할 수 있다. WC는 0에서 1 사이의 값을 가지며, 0에 가까울수록 두 시계열 간의 선형적인 관계가 약하며, 1에 가까울수록 두 시계열 간의 선형적인 관계가 강하다는 것을 의미한다(Labat et al., 2005). 아울러 cross wavelet transform의 실수부와 복소부가 만들어내는 복소 편각(complex argument)은 두 시계열의 국부적인 상대적 위상(relative phase)을 나타내며, 이는 phase angle이라 한다. Phase angle은 통상 WC의 결과가 0.5 이상인 경우에 대해서 WC의 결과 위에 화살표로 표현된다(De Boer, 1985). 실제 자료들에 대한 적용 시 유의수준에 대한 고려를 위해 Monte Carlo 시뮬레이션 방법이 주로 활용되며, 유의수준 값으로 통상 0.05가 사용된다(Torrence and Compo, 1998).

Fig. 9는 원자료에 대한 XWP 및 WC를 보여준다. 이를 살펴보면 약 1-2년 주기 및 약 2-4년 주기에 존재하는 강한 XWP가 부각된다. 매우 짧은 주기에서도 일부 강한 XWP가 관찰된다. 한편, WC의 경우 XWP와는 다른 양상을 보여준다. 약 16년 이하에서 강한 WC가 존재하고, 4-8년 주기에 존재하는 성분은 매우 긴 시간에 걸쳐 있다. 아울러 phase는 초반부에서는 NAR이 선행하며, 중반부에서는 NMA가 선행하다 다시 후반부에 그 양상이 전환되는 것을 확인할 수 있다.

Fig. 10은 확장된 자료에 대한 XWP 및 WC를 보여준다. 먼저 XWP를 살펴보면 Method (1)의 경우 전반적으로 약해진 XWP가 다소 퍼져서 나타나는 양상을 보인다. 이는 짧은 스케일에서 타방법에 비해 짧은 블록 길이가 고려됨에 따른 결과라 판단된다. Method (2)와 (3)의 경우 Method (1)에 비해서 XWP가 선명하게 나타난다. 확장된 자료의 길이(5T)를 고려한다면 강하게 나타나는 XWP의 빈도가 적절하게 표현된다고 생각된다. WC에 대한 결과를 살펴보면, Method (1)에서 짧은 주기 성분이 표현되지 않는다는 점이 관찰된다. 앞선 XWP의 결과와 마찬가지로 이는 블록 길이에 따른 영향으로 판단된다. 아울러 원자료 간의 결과에서는 확인하기 어려웠던 약 16년 이상의 장주기 WC 성분이 관찰된다. 이는 세 방법 모두에서 공통적으로 확인되는 결과이다. 원자료의 WC의 4-8년 주기에서 약 30년 이상에 걸쳐 나타났던 성분은 모의 자료 결과에서는 다소 단절되어 나타난다. Method (3)에서는 다른 방법에 비해 이러한 성분이 보다 길게 표현된다. 이는 동시간대의 스펙트럼이 추출됨에 따라 스펙트럼의 스케일 간 관계가 어느 정도는 보존되고 있음을 보여준다. 한편 phase angle의 경우 세 방법 모두 특별한 경향성을 보이고 있지는 않다. 이는 자료 확장 시 두 시계열에 대한 선행 관계를 특별히 고려하지는 않았기 때문이라 판단된다.

확장을 시도할 때 부여되는 난수에 따라 각 방법별로 모의 자료의 분포가 달라진다. 다시 말해, 이상의 결과는 특정 난수에 대한 하나의 사례에 불과하다. 이로 인해 관련 연구들(Erkyihun et al., 2016; Erkyihun et al., 2017)에서는 다수의 확장 자료를 생성하여 이들의 분포를 비교한 바 있다. 참고로 이들은 동일한 가상 기간에 대한 자료 확장을 위해 500개의 앙상블을 고려하였다. 이에 본 연구에서는 확장되는 자료의 길이 및 임의성을 충분히 고려하기 위해 총 5,000개의 확장 자료를 생성하였다.

Fig. 11은 모의자료에 대한 앙상블 결과를 확률분포함수(Probability Density Function, PDF) 및 박스플롯을 이용하여 나타낸 것이다. 여기서, PDF를 이용하여 확장된 자료의 분포를 원자료의 것과 비교하고자 하였으며, PDF 추정을 위해 kernel density estimation (Rosenblatt, 1956)을 이용하고, 전체 확장 자료들의 PDF를 세부 구간별 box-plot으로 나타내었다. 먼저, NMA의 결과를 살펴보면 Method (3)가 원자료를 가장 잘 모방하는 것을 확인할 수 있다. Method (1)과 (2)의 경우 앞서 Box-plot에서도 확인한 것과 같이 과도하게 크거나 작은 값들은 놓치게 된다. 아울러 이 두 방법은 확장되는 자료들의 분포가 정규분포에 가까운 모습을 보여준다. 한편, NAR의 경우 비슷한 결과를 보여주나, 전반적으로 NMA보다는 그 차이가 작게 나타난다. 이러한 결과를 통해, 원자료의 특성에 따라 방법별 자료 확장 성능이 달라질 수 있다는 점 또한 확인할 수 있다.