Comparative assessment and uncertainty analysis of ensemble-based hydrologic data assimilation using airGRdatassim

Garim Lee; Songhee Lee; Bomi Kim; Dong Kook Woo; Seong Jin Noh

doi:10.3741/JKWRA.2022.55.10.761

Preview

Research Article

Journal of Korea Water Resources Association. 31 October 2022. 761-774
https://doi.org/10.3741/JKWRA.2022.55.10.761

Comparative assessment and uncertainty analysis of ensemble-based hydrologic data assimilation using airGRdatassim

airGRdatassim을 이용한 앙상블 기반 수문자료동화 기법의 비교 및 불확실성 평가

Garim Lee^a

Songhee Lee^b

Bomi Kim^c

Dong Kook Woo^d

Seong Jin Noh^e^*

이 가림^a

이 송희^b

김 보미^c

우 동국^d

노 성진^e^*

^aPh.D. Course Student, Department of Civil Engineering, Kumoh National Institute of Technology, Gumi, Korea

^bMaster’s Course Student, Department of Civil Engineering, Kumoh National Institute of Technology, Gumi, Korea

^cPh.D. Course Student, Department of Civil Engineering, Kumoh National Institute of Technology, Gumi, Korea

^dAssistant Professor, Department of Civil Engineering, Keimyung University, Daegu, Korea

^eAssistant Professor, Department of Civil Engineering, Kumoh National Institute of Technology, Gumi, Korea

^a금오공과대학교 토목공학과 석박사 통합과정

^b금오공과대학교 토목공학과 석사과정

^c금오공과대학교 토목공학과 석박사 통합과정

^d계명대학교 토목공학전공 조교수

^e금오공과대학교 토목공학과 조교수

^{*Corresponding Author}

ABSTRACT

Accurate hydrologic prediction is essential to analyze the effects of drought, flood, and climate change on flow rates, water quality, and ecosystems. Disentangling the uncertainty of the hydrological model is one of the important issues in hydrology and water resources research. Hydrologic data assimilation (DA), a technique that updates the status or parameters of a hydrological model to produce the most likely estimates of the initial conditions of the model, is one of the ways to minimize uncertainty in hydrological simulations and improve predictive accuracy. In this study, the two ensemble-based sequential DA techniques, ensemble Kalman filter, and particle filter are comparatively analyzed for the daily discharge simulation at the Yongdam catchment using airGRdatassim. The results showed that the values of Kling-Gupta efficiency (KGE) were improved from 0.799 in the open loop simulation to 0.826 in the ensemble Kalman filter and to 0.933 in the particle filter. In addition, we analyzed the effects of hyper-parameters related to the data assimilation methods such as precipitation and potential evaporation forcing error parameters and selection of perturbed and updated states. For the case of forcing error conditions, the particle filter was superior to the ensemble in terms of the KGE index. The size of the optimal forcing noise was relatively smaller in the particle filter compared to the ensemble Kalman filter. In addition, with more state variables included in the updating step, performance of data assimilation improved, implicating that adequate selection of updating states can be considered as a hyper-parameter. The simulation experiments in this study implied that DA hyper-parameters needed to be carefully optimized to exploit the potential of DA methods.

Keywords

Data assimilation

Ensemble Kalman filter

Particle filter

airGRdatassim

Uncertainty analysis

가뭄과 홍수의 예측, 기후변화가 유역 유출량, 더 나아가 수질 및 생태계에 미치는 영향의 정확한 분석을 위해서는 수문 모의 과정의 불확실성을 정량화하고 최소화하기 위한 노력이 필요하다. 수문자료동화는 수문모형의 상태량이나 매개변수를 갱신(update)하여 모의 초기 조건의 가장 가능성 있는 추정치를 생성하는 기법으로, 실시간 관측 정보를 이용하여 예측 정확도를 향상시킬 수 있는 방법이다. 본 연구에서는 airGRdatassim 모형을 이용하여 앙상블 기반 순차 자료동화 기법인 앙상블 칼만 필터와 파티클 필터로 용담댐 유역에 대해 일 유출을 모의하고, 자료동화 기법별 특성을 비교 및 분석하였다. 모의 결과, Kling-Gupta efficiency (KGE) 지표가 자료동화 적용 전 0.799에서 앙상블 칼만 필터와 파티클 필터 적용시 각각 0.826, 0.933으로 향상되었다. 또한 기상 강제력 노이즈의 범위, 갱신 대상 상태량 설정, 앙상블 수 등 수문자료동화의 설정과 관련된 하이퍼-매개변수(hyper-parameter)의 불확실성이 모의 예측 성능에 미치는 영향을 분석하였다. 강수 및 잠재 증발산 강제력의 오차 범위에 대한 민감도 분석 결과, 모든 모의 범위에서 파티클 필터가 앙상블 칼만 필터보다 예측 성능이 우수하였다. 파티클 필터는 기상 강제력 오차 크기가 작을수록 모의 성능이 향상되었으며, 앙상블 칼만 필터는 상대적으로 오차가 큰 경우 최적 성능이 확인되었다. 한편, 자료동화시 갱신되는 상태량의 종류를 줄일수록 자료동화에 의한 모의 성능은 감소하였다. 본 연구의 모의 실험 결과는 앙상블 자료동화를 이용하여 일 유출 모의 정확도 향상이 가능하지만, 최적 성능을 발휘하기 위해서는 수문자료동화 기법별 하이퍼-매개변수의 적정한 조정이 필요함을 함의한다.

키워드

자료동화

앙상블 칼만 필터

파티클 필터

airGRdatassim

불확실성 분석

MAIN

1. 서 론
2. 연구방법 및 자료
2.1 순차 자료동화 방법
2.2 수문 모의 및 자료동화 소프트웨어
2.3 적용 지역 및 입력자료
2.4 모의실험 및 평가 기준 설정
3. 결과 및 분석
3.1 수문자료동화 모의 결과 비교
3.2자료동화 하이퍼-매개변수 영향 분석
4. 결 론

1. 서 론

수문모형은 물순환 과정에 대한 정확한 이해와 지속가능한 수자원의 관리를 위해 필수적인 소프트웨어 도구이다. 가뭄과 홍수의 예측, 기후변화가 유량, 수질 및 생태에 미치는 영향 평가 등을 위해서는 수문모형을 통한 수자원 요소의 정확한 해석이 필요하다. 하지만, 수문 모의의 정확성은 1) 수문모형의 구조 및 모형 매개변수(model parameters)의 추정, 2) 예측 시작 시점에서의 수문모형 상태량(state variables) 초기조건, 3) 기상 강제력(weather forcings) 조건, 4) 인간 활동의 물순환 개입(예: 저수지 운영, 관개 등) 등 다양한 요소로부터 영향을 받는 것으로 알려져 있다(Noh et al., 2014). Blöschl et al. (2019)은 21세기 수문학에서의 주요 문제 중 하나로 수문 모의의 다양한 불확실성 원인을 규명하고 완화해야 할 필요성을 언급한 바 있다. 자료동화(Data Assimilation, DA)는 수문모형의 상태량이나 매개변수를 갱신(update)하여 수문모형 초기조건의 가장 가능성 있는 추정치를 생성하는 기법으로(Boucher et al., 2020; Piazzi et al., 2021), 수문 모의의 불확실성을 최소화하고 예측 정확도를 향상시킬 수 있는 것으로 알려져 있다. 예를 들어, 기상 예보나 홍수 예보 등을 위해서는 현재 상태의 정보를 이용하여 미래 시점의 예측 값 추정이 필요하다. 이 때 예측이 시간 단위로 진행됨에 따라 오차도 누적되어 특정 시간 이후의 예측 정확도가 감소하게 되는데, 이를 방지하기 위해 최신의 관측 정보를 이용하여 현재 시간에서의 모형 상태 정보를 갱신하는 기법을 자료동화라 한다.

최근 들어 하천유량이나 수위, 토양수분 등 실시간으로 가용한 수문관측 정보가 늘어남에 따라 수문자료동화(hydrologic data assimilation) 기법의 활용 또한 늘어나고 있다. Lee and Bae (2011)의 연구에서는 안동댐 상류 유역의 실시간 하천 유량 예측시스템에 자료동화기법을 연계할 경우, 자료동화를 적용하지 않은 결정론적 모의(deterministic simulation)에 비해 정확한 홍수량 예측이 가능하다고 보고하였다. Yoo et al. (2012)는 충주댐 유역을 대상으로 확장 칼만 필터(extended Kalman filter)를 적용하여 결/오측이나 포함된 실시간 유량이나 댐 유입량 자료의 품질향상에 효과적으로 작동함을 확인하였다. Lee et al. (2013)는 앙상블 칼만 필터 기법(Ensemble Kalman filter, EnKF), 연속형 강우-유출 모형 SURF모형, AutoRO을 결합한 실시간 댐 수문량 예측모형(DHVPM)을 개발하고 그 적용성을 평가하였다. 그 결과 실시간 예측유입량의 정확도를 향상시키기 위해서는 관측 유입량의 실시간 자료동화기법을 연계한 유입량 예측이 바람직하다고 평가하였다. Choi and Kim (2021)은 월단위 유량에 대한 합성 관측자료(synthetic observations)와 파티클 필터(Particle Filter, PF)를 이용하여 수문모형의 시간변동 매개변수를 추정하고 월 유량(monthly discharge) 예측 성능을 평가하였다. 그 결과, 월 유량의 예측 성능은 파티클 필터 방법과 Shuffled Complex Evolution Metropolis (SCEM)방법이 서로 비슷하였고, 우기를 제외한 시기의 월 유출 예측 성능은 자료 동화 방법을 이용한 주기적으로 변하는 매개변수에 의한 모의가 더 우수하였다. 또한 합성 관측 자료의 불확실성은 파티클 필터의 예측 성능에 중요한 역할을 한다고 보고하였다. Leach and Coulibaly (2019)는 캐나다 Ontario주의 Don강과 Black Creek강 유역에 대해 자료동화를 적용하여 하천유량 예측 신뢰도가 최대 78%까지 향상됨을 보여주었다. Piazzi et al. (2021)는 프랑스 내 232개 유역을 대상으로 집중형 강우-유출 모형 GR5J에 대해 앙상블 칼만 필터와 파티클 필터를 적용하고, 상태량 및 매개변수 동시 갱신(dual state-parameter updating) 기법의 최적 구성 조건에 대해 연구하였다.

자료동화 기법은 새로운 관측정보가 가용할 경우 모형의 상태량(혹은 매개변수)을 갱신(update)하는 방식에 따라 1) 순차 자료동화(sequential data assimilation)와 2) 변분법(variational assimilation)의 크게 두가지 종류로 나뉜다. 한편, 순차 자료동화는 사후 분포(posterior distribution)를 정규 분포(normal distribution) 형식으로 가정하는 칼만 필터와 사후 분포의 형식에 대한 제약이 없는 것이 특징인 파티클 필터의 두가지 종류로 나눌 수 있으며, 예측 모형(자료동화 일반 이론에서의 forward operator)의 비선형성을 고려하기 위해 모형 앙상블(ensemble) 혹은 파티클(particle)을 사용한다. 이 두가지 순차 자료동화 방법은 통계분포에 대한 가정이나 상태량 갱신 방식에서 분명한 차이가 있으므로, 수문 예측 정확도 향상을 위해서는 다양한 불확실성 조건에 대한 평가를 통해 최적 자료동화 기법 구성하는 것이 중요하다. 앙상블 칼만 필터와 파티클 필터에 대해 집중형 수문모형을 대상으로 비교한 국외 연구 사례는 Piazzi et al. (2021)의 연구가, 공간분포형 수문모형을 대상으로는 Noh et al. (2012)의 연구가 있다. 국내 사례로는 파티클 필터를 이용하여 토양침식모형의 매개변수를 추정한 Kim et al. (2015)의 연구, 집중형 강우-유출 모형의 매개변수 추정에 적용한 Lee et al. (2017)의 연구를 들 수 있다. 다만, 파티클 필터와 앙상블 칼만필터 등 순차 자료동화 기법의 특성을 비교한 국내 연구 사례는 드물다.

대부분의 앙상블 방식 순차 자료동화 기법에서는 모의와 관측의 불확실성에 대한 통계 모형을 상정하고, 이 통계 모형의 분포형 내에서 난수를 생성하여 여러가지 다른 초기 및 매개변수 조건을 생성하고, 이를 비선형 수문모형에 적용하여 불확실성을 추정하는 것이 일반적이다. 이 과정에서 모형 앙상블의 생성 방식과 앙상블의 불확실성 범위가 자료동화 예측 성능에 큰 영향을 주는 것으로 알려져 있으나(Noh et al., 2014; Shen et al., 2022), 수문모형의 불확실성 요소가 자료동화에 미치는 영향에 대한 연구는 비교적 많지 않은 편이다.

단일한 입력조건에 대해 단일한 예측 결과만 도출 가능한 결정론적 방식으로 구성된 대부분의 수문모형을 앙상블 방식으로 확장하고, 순차적으로 모형의 초기값 및 매개변수를 갱신할 수 있도록 모형의 작동방식을 변환해야 하는 점은 자료동화 기술의 적용을 어렵게 하는 요소 중 하나이다(Liu et al., 2012). 이러한 난점을 해소하기 위해 기존 모형의 프로그램 코드를 변경하지 않고 자료동화 기법을 적용할 수 있도록 지원하는 범용 자료동화 소프트웨어 프레임워크가 개발되었다. 대표적인 사례로 Data Assimilation Research Testbed (DART), Parallel DA Framework (PDAF), open source initiative for DA (OpenDA) 등이 있으나, 기상이나 해양 분야에 비해 수자원 분야의 적용 사례는 상대적으로 적은 편이다. 반면, 최근 개발된 airGRdatassim (Piazzi et al., 2021)은 GR5J 등 GR 계열 일 유출 해석 수문모형에 대해 앙상블 자료동화를 적용할 수 있도록 지원하는 오픈소스 R package로, 비교적 간단히 수문자료동화 기법을 적용할 수 있는 장점이 있다.

본 연구에서는 앙상블 방식의 순차 자료동화의 대표적 기법인 앙상블 칼만 필터와 파티클 필터를 유역단위 일 유출 모의에 적용하여, 자료동화 기법에 따른 차이점과 특징을 분석한다. 또한, 앙상블 생성 기법의 통계모형의 가정, 섭동(perturbation) 및 갱신(update) 대상 상태량 선정 등 자료동화 방법과 관련된 하이퍼-매개변수(hyper-parameter)의 영향을 검토하여, 순차 자료동화 기법의 최적화 방법을 제시한다. 구체적으로는 연구대상지역인 용담댐 유역에 대해 앙상블 칼만 필터와 파티클 필터를 적용하여 비교하고, 강수 및 잠재 증발산의 불확실성 범위가 일 유출 해석 정확성에 미치는 영향을 분석한다.

2. 연구방법 및 자료

2.1 순차 자료동화 방법

칼만 필터나 파티클 필터와 같은 순차 자료동화 방법에서는 시간에 따른 상태량의 변화, 관측과 모형과의 관계를 다음과 같은 상태-공간방정식(state-space equations)으로 일반화하여 표현한다.

(1a)

x_{t} = F (x_{t - 1}, I_{t}, v_{t})

(1b)

y_{t} = H (x_{t}, w_{t})

여기서, F와 H는 각각 상태방정식과 관측방정식이며, $x_{t}, I_{t}, v_{t}, y_{t}, w_{t}$ 는 각각 시간 t에서의 상태변수, 입력자료, 상태변수의 노이즈(noise), 관측변수, 관측변수의 노이즈이다. 수문자료동화에서 상태방정식 F는 강수나 기상자료와 같은 입력자료 I_t에 따른 유출량, 토양수분, 증발 등과 같은 수문 상태량 x_t의 예측 과정을, 관측방정식 H는 여러가지 세부 수문요소로부터 관측가능한 수문요소를 추정하는 과정을 각각 의미한다. 즉, 상태-공간방정식은 모두 수문모형 해석 과정을 통해 구현되며, 상태량(이나 매개변수를) 갱신하는 방법에 따라 칼만 필터나 파티클 필터 계열로 구분할 수 있다. 자료동화에서 모형의 구조, 입력자료, 매개변수의 불확실성은 상태변수 노이즈 v_t로 표현되며, 노이즈를 기상자료와 같은 입력자료에 적용하거나(Piazzi et al., 2021), 상태량에 노이즈를 직접 더하거나 곱하는 방법(Noh et al., 2011b)이 주로 사용된다.

Fig. 1은 순차 자료동화 방법 중 앙상블 칼만 필터(Fig. 1(a))와 파티클 필터(Fig. 1(b))의 개념을 비교한 것이다(노란색선: 관측에 대한 우도, 검은색선: 모의 상태량의 사전 정규분포, 파란색선: 모의 상태량의 사후 분포). Fig. 1(a)의 앙상블 칼만 필터는 자료동화 전, 후의 상태량의 분포, 즉 사전 분포(prior distribution)와 사후 분포(posterior distribution)를 정규 분포(normal distribution)로 근사(approximation)시키지만, Fig. 1 (b)의 파티클 필터의 사전, 사후 분포는 정규 분포가 아닌 것이 두 방법의 가장 큰 차이점이다. 이에 따라 관측에 대한 우도를 정규분포로 동일하게 가정하더라도 사후 분포의 형태가 두 자료동화 기법에서 달라지게 된다. Fig. 1(b)는 파티클 필터 중Sequential Importance Resampling (SIR) 방법을 도시화한 것으로, 관측과 비례한 가중치를 할당하여 리샘플링 과정을 통해 앙상블을 생성한다. 각 자료동화 방법에서 상태량의 사후 분포를 추정하는 이론적 방법에 대해서는 다음 절에서 설명한다.

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F1.jpg

Fig. 1.

Schematic diagram of data assimilation methods: (a) the ensemble Kalman filter and (b) the particle filter

2.1.1 앙상블 칼만 필터

앙상블 칼만 필터에서는 새로운 관측 정보(observation information)가 가용해질 때마다 관측과 모의와의 차이에 칼만 이득(Kalman gain)을 곱하여, 최적의 상태량을 추정한다(Evensen, 1994).

시간 t의 i번째 앙상블( $X_{t, i}^{a}$ )의 분석 상태(analysis state) 추정식은 다음 식과 같다.

(2)

x_{t, i}^{a} = x_{t, i}^{b} + K_{t} [y_{t, i} - H [x_{t, i}^{b}]]

$x_{t, i}^{b}$ 는 i번째 앙상블의 백그라운도 상태량(갱신전 상태량, 변수 길이 nstate), y_t,i는 i번째 앙상블의 관측정보(변수 길이 nobs), H는 Eq. (1b)의 관측방정식이다(Evensen, 1994). 칼만 이득(Kalman gain, K_t)은 최적의 업데이트를 위한 가중치로 모의와 관측정보로부터 산정된 오차 공분산 행렬의 조합이다.

(3)

K_{t} = P_{t} H^{T} {(H P_{t} H^{T} + R_{t})}^{- 1}

R_t는 관측 값의 오차 공분산 행렬(크기 nobs × nobs)이고, P_t는 앙상블로부터 생성되며 각 동화단계에서 동적으로 업데이트되는 모의 값의 오차 공분산 행렬로 자세한 내용은 Evensen (2003) 등을 참고할 수 있다.

앙상블 칼만 필터에서는 관측 및 모의 상태량의 사전 및 사후 분포를, 평균과 분산의 두가지 모멘트만으로 표현가능한 정규 분포(normal distribution)로 근사(approximation)하기 때문에 비교적 적은 수의 앙상블로도 다차원 문제에 대한 추정이 가능하여 여러 분야의 예측 응용 문제에 적용되고 있다.

상태량 뿐 아니라 매개변수 또한 칼만 필터를 통해 갱신 하는 방법론이 개발되어 활용되고 있으나(Hendricks Franssen and Kinzelbach, 2008), 일 유출 모의에 대해 자료동화를 적용한 기존 연구사례(Piazzi et al., 2021)에서 상태량-매개변수 동시 갱신의 효과가 그다지 크지 않아 본 연구에서는 상태량 갱신 만을 다룬다.

2.1.2 파티클 필터

파티클 필터(Particle Filter, PF)는 순차 몬테카를로(Sequential Monte Carlo, SMC) 방법이라고도 불리며, 비선형(non-linear) 비가우시안(non-Gaussian) 모형에 대한 자료동화에 적용된다. 특히, 기존의 칼만 필터 계열 자료동화 방법과 달리 정규 분포나 선형화 가정이 요구되지 않는 장점을 갖는다. 베이즈 이론(Baye’s theory), 마르코프 특성(Marov property)과 순차 중요도 샘플링(Sequential Importance Sampling, SIS) 기법(Arulampalam et al., 2002)에 의해 유도된 파티클 필터에서의 시간 t의 사후 확률은 다음과 같이 근사될 수 있다(Doucet et al., 2000).

(4a)

p (x_{t} | y_{t}) = \sum_{i = 1}^{n} w_{t, i} δ (x_{t} - x_{t, i})

(4b)

w_{t, i} \propto \frac{p (x_{t, i} | y_{t})}{q (x_{t, i} | y_{t})}

여기서, w_t,i는 가중치, $δ$ 는 Dirac delta함수, q(x_t,i | y_t)는 주요 확률 밀도이다. 다만, Ristic et al. (2004)이 지적한 바와 같이 SIS방법은 파티클의 퇴화현상(degeneracy; 하나의 파티클이 모든 가중치를 갖는 현상)이 발생할 수 있어, SIS에 리샘플링(resampling) 과정을 추가한 SIR (Sequential Importance Resampling)방법이 파티클 필터의 일반적인 방법으로 적용된다. 리샘플링은 세부 과정에 따라 multinomial resampling, stratified resampling, systematic resampling, residual resampling 등이 있으며, 리샘플링의 주기 또한 유효 파티클 수(effective particle number)를 기준으로 수행하거나, 매 자료동화 단계마다 수행하는 방법이 있다. 파티클 필터의 세부 기법과 수문 모의 적용 사례는 Noh (2013)의 연구를 참고할 수 있다. 본 연구에서는 multinomial resampling을 매 자료동화 단계마다 적용하였다. 또한, 각 파티클의 가중치 $w_{t, i}$ 는 다음의 식과 같이 관측에 대한 우도(likelihood)로부터 산정하였다.

(5)

w_{t, i} = p (y_{t} | x_{t, i}) = \frac{1}{σ \sqrt{2 π}} e^{- \frac{{(x_{t, i} - y_{t})}^{2}}{2 σ^{2}}}

여기서, $σ$ 는 관측자료의 표준편차로, Eq. (1b) 관측방정식에서 관측변수의 오차 w_t의 특성과 관련된다. 본 연구에서는 관측유량의 10%를 관측분포의 표준편차로 적용하였다.

2.1.3 기상 강제력 앙상블 생성 기법

대부분의 앙상블 자료동화 기법에서는 상태량, 매개변수, 기상 강제력 등을 확률적으로 섭동하여 입력자료와 모형의 불확실성을 추정하는 방식을 사용한다. airGRdatassim에서는 Clark et al. (2008)의 연구에서 제안된 기법으로 기상 강제력(강수 및 잠재 증발산) 앙상블을 생성하며 그 과정을 요약하면 다음과 같다. 개별 앙상블별 강수 강제력 p'는 다음 식과 같이 강수 강제력 p에 승수(multiplier) $φ_{p}$ 를 곱하여 생성한다.

(6a)

p' = p φ_{P}

(6b)

φ_{P} = (1 - ε_{p}) + 2 u_{p} ε_{p}

여기서, u_p는 [0,1] 범위의 균일 분포(uniform distribution)이고, 따라서 $φ_{p}$ 는 [1‒ $ε_{p}$ , 1+ $ε_{p}$ ] 범위의 균일 분포로부터 생성된 난수이다. Eq. (6)의 변수들은 매 모의시간 마다 변하는 값이지만, 식의 가독성과 설명의 편이를 위해 시간변수 t는 표시하지 않았다.

잠재 증발산 강제력 앙상블 E' 또한 강수 강제력과 동일한 방법으로, [1‒ $ε_{e}$ , 1+ $ε_{e}$ ] 범위의 균일 분포로부터 추출된 난수인 승수 $φ_{e}$ 를 잠재 증발산 입력자료 E'에 곱하여 생성한다. 즉, 기상 강제력 앙상블의 다양성은 강수 및 잠재 증발산 강제력 승수의 범위에 대한 변수인 $ε_{p}$ , $ε_{e}$ 에 의해 조절될 수 있다. 이러한 앙상블 생성 방식은 기존의 수문자료동화 연구에서 가장 활발히 적용되는 방법론 중 하나로, 하이퍼-매개변수 $ε_{p}$ , $ε_{e}$ 는 기상 자료의 불확실성을 표현할 뿐 아니라 수문모형 앙상블의 확률분포, 궁극적으로 자료동화의 성능에 영향을 미치는 요소이다. 본 연구에서는 기상 강제력 앙상블의 범위 관련 하이퍼-매개변수가 자료동화 성능에 변화를 분석하였으며, 모의 실험의 구성과 결과는 각각 2.4절과 3절에서 설명한다. 한편, 기상 강제력 시계열의 불확실성은 매시간 독립적이지 않고 시간 축에서 서로 연관될 수 있다. Clark et al. (2008)의 방법론에서는 다음의 식과 같은 1차 자기 회귀 변형 모형을 통해 앙상블에 적용되는 난수의 시공간적 상관관계를 고려한다.

(7a)

s_{t} = ρ s_{t - 1} + \sqrt{1 - ρ^{2}} W_{t - 1}

(7b)

u_{t} = \frac{1}{2} e r f c (\frac{s_{t}}{\sqrt{2}})

여기서 s_t는 정규분포에서 발생된 현재 시간 t의 난수로, 시간 상관관계에 대한 가중치 $ρ$ 를 통해 이전 시간 t‒1의 난수 s_t_‒1과 연계된다. W_t_‒1은 공간 상관도에 관한 계수이다. s_t가 결정되면, Eq. (7b)의 error function에 의해 균일 분포로 변환되어, Eq. (6b)의 u_p 혹은 u_e로 적용된다.

시간 상관관계 가중치 $ρ$ 는 다음과 같이 정의될 수 있다.

(8)

ρ = 1 - \frac{△ t}{τ}

여기서 $τ$ 는 시간 역상관 길이, ∆t는 모의 시간 간격이다(Clark et al., 2008). 본 연구에서는 집중형 수문모형을 이용하여 일단위 기상 강제력 자료를 입력자료로 일유량을 모의하였기 때문에, 공간 상관 관계는 고려하지 않았으며, 시간 상관 가중치는 기존 연구(Piazzi et al., 2021)에서 제시한 값(시간 역상관 길이: 강수 1(day), 잠재 증발산 2(day))을 적용하였다.

2.2 수문 모의 및 자료동화 소프트웨어

2.2.1 airGRdatassim

airGRdatassim은 프랑스 INRAE-Antony 연구소의 Catchment Hydrology 연구팀에서 개발한 R패키지로, airGR 수문학적 모델링 패키지를 기반으로 GR 계열 일 유출 수문모형(GR4J, GR5J, GR6J)을 이용한 앙상블 기반 자료동화 기능을 제공한다(https://gitlab.irstea.fr/HYCAR-Hydro/airgrdatassim). 모형 상태량, 기상강제력, 매개변수 등 다양한 조건에 대한 불확실성을 고려하여 앙상블을 생성하고 자료동화 기법을 적용할 수 있다. 자료동화 방법으로 앙상블 칼만 필터와 파티클 필터, 두가지 방법을 선택하여 사용할 수 있으며, 파티클 필터 기법은 SIR방법이 적용된다. airGRdatassim은 RunModel_DA (DA_EnKF, DA_PF 모듈 포함), CreateInputsPert, StatePert (StateEnKF, StatePF 모듈 포함) 크게 세가지 모듈로 구성된다. StatePert 모듈을 통해 모형의 섭동(perturbation) 대상을 지정하고, StateEnKF모듈은 앙상블 칼만 필터, StatePF모듈은 파티클 필터의 갱신(update) 대상을 지정할 수 있다. CreateInputsPert 모듈은 1차 자기 회귀 변형 모델에 의해 임의의 섭동을 생성하며, 기상강제력(강수 및 잠재 증발), 예측선행시간, 앙상블 크기 등을 설정하여 교란시킨다. 강수 및 잠재 증발산 오차 매개변수의 기본값은 0.65로 설정되어 있으며 시간 역상관 길이의 기본값은 강수 1일(day), 잠재 증발 2일(day)로 설정되어 있다. 본 연구에서는 GR 계열 일 유출 해석 수문모형 중 GR5J를 사용하여 앙상블 칼만 필터와 파티클 필터 자료동화를 수행한다.

2.2.2 GR5J 수문모형

GR5J 모형은 유역단위 물순환 과정을 토양수분 및 하도의 저류, 추적 과정을 개념화한 집중형 수문모형이다(Lavenne et al., 2016; Moine, 2008). GR5J 모형은 두 개의 저류지 모듈과 한 개의 단위유량도를 적용하며(Fig. 2), X1부터 X5까지 5개의 매개변수가 있다(Table 1). X1(mm)과 X3(mm)는 각 저류지에 저장될 수 있는 최대저류량을 의미한다. 입력자료인 강수와 잠재 증발산 중 일부는 첫번째 토양수분 저류지에서 유효강수량으로 환산된다. 두 번째 저류지에서 단위유량도를 통해 흘러온 유효 강수량 중 일부는 하천으로 유입되고 나머지는 지표면을 통해 유출된다. 하천으로부터 유출되는 유량과 지표면을 통해 유출되는 유량은 유역 간 교환계수인 X2 (mm/d)와 유역 교환 임계값인 X5에 의해 유역 간의 지하수 흐름이 정량화 된다(Adeyeri et al., 2020). 최종적으로 하도와 지하수 유출을 합산하여 유역의 유출량이 산정되는 모형 구조를 가지고 있다. 본 연구에서는 airGRdatassim을 개발하여 자료동화를 적용한 기존 연구(Piazzi et al., 2021)의 연장선상에서 국내 적용성 및 불확실성을 분석하기 위해서 GR계열 모형 중 GR5J를 선정하였다.

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F2.jpg

Fig. 2.

Structure of GR5J model

Table 1.

Descriptions of GR5J model parameters (Adeyeri et al., 2020)

Model parameters	Description
X1	Production store capacity (mm)
X2	Inter-catchment exchange coefficient (mm/day)
X3	Routing store capacity (mm)
X4	Unit hydrograph time constant (day)
X5	Inter-catchment exchange threshold (-)

2.3 적용 지역 및 입력자료

본 연구는 금강 유역의 상류에 위치한 용담댐을 대상 유역으로 선정하였다. 용담댐은 전라북도 진안군에 위치하고 있으며, 댐 상류 유역은 무주군, 진안군, 장수군을 포함한 충청남도, 전라북도, 경상남도에 걸쳐 있고 유역 면적은 930 km²로 금강 유역 면적의 약 9.5%를 차지하고 있다.

GR5J 모형의 입력자료는 기상자료개방포털(http://data.kma.go.kr)과 국가수자원관리종합 정보시스템(WAMIS)에서 제공하는 일 단위 평균 기온, 강우, 댐의 유입량 자료를 사용하였다. 모형의 입력자료 구축을 위해 유역 반경 5 km 이내에 있는 12 개의 강우관측소와 4 개의 기상관측소 자료를 활용하였다(Fig. 3(a)). 2001 년부터 2014 년까지의 일 단위 평균 강우량 및 기온 자료 구축에는 Thiessen 가중평균법을 사용하였다(Figs. 3(b) and 3(c)). 용담댐 유역 강우의 연 평균과 표준편차는 각각 1409.6 mm, 303.2 mm이고 연 평균 기온은 11.3℃로 산정되었다.

모형 입력자료 중 잠재증발산은 Oudin et al. (2005)의 연구에서 제시한 기온 기반 방정식을 사용하였으며 그 식은 다음과 같다.

(9)

P E \{\begin{cases} 0.408 R_{e} (T + 5) & if (T + 5) > 0 \\ 0 & o t h e r w i s e \end{cases}

여기서, R_e는 태양복사에너지(MJ m^-2d^-1)이고, T는 2 m 높이에서의 평균 기온(℃)이다(Oudin et al., 2005; 2010). 본 연구에서는 Eq. (9)의 평균 기온은 Thiessen 가중평균에 의해 산정된 면적 평균값을 적용하였다.

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F3.jpg

Fig. 3.

Study area and observation network: (a) Yongdam-dam catchment (b) Thiessen network of rainfall stations, and (c) Thiessen network of weather stations. Associated Thiessen polygons are represented in black lines

2.4 모의실험 및 평가 기준 설정

2.4.1 모의 성능 평가 기준

수문 모의의 성능평가 지표로는 Kling-Gupta Efficiency (KGE)를 사용하였으며, 식은 다음과 같다.

(10)

K G E = 1 - \sqrt{{(γ - 1)}^{2} + {(α - 1)}^{2} + {(β - 1)}^{2}}

여기서 $γ$ 은 관측 값과 모의 값 사이의 선형 상관계수, $α$ 는 관측 값의 표준편차에 대한 모의 값의 표준편차로 관측과 모의 값의 사이 변동성에 대한 오차의 척도를 나타내며, $β$ 는 관측 평균에 대한 모의 평균 비이다. KGE가 1에 가까울수록 관측과 모의 사이의 일치도가 높다는 것을 의미한다(Gupta et al., 2009; Knoben et al., 2019). Nash-Sutcliffe Efficiency (NSE) 등의 지표가 전체 자료기간 중 일부 홍수 기간의 모의 결과에 지배적으로 영향을 받는데 비해, KGE는 $α, β$ 등에 대한 고려를 통해 비강우 기간에 대한 모의 성능 비교가 가능하기 때문에 평가 지표로 적용하였다.

2.4.2 모의 실험 조건

모의 실험은 용담댐 유역에 대해 구축된 총 14년의 일 단위 시계열 자료(2001~2014년)를 세가지 기간으로 구분하여, 모형 매개변수 보정은 2001~2004년(4년), 자료동화 하이퍼-매개변수 보정은 2005~2008년(4년), 자료동화 검정은 2009~2014년(6년)의 기간에 대해 수행하였다(Table 2). 모의초기 조건이 open loop 결과에 미치는 영향을 최소화하기 위해 자료동화 검정 기간 모의시 1년간(2008년)을 warm-up 기간으로 설정하여 추가 모의하고, 모의 성능은 검정 기간에 대해서만 평가하였다.

목적함수를 KGE로 적용하여 최적 보정된 GR5J의 모형 매개변수는 각각 X1=88, X2=1.61, X3=37, X4=0.9, X5=0.41였으며(Table 3), 이 때 KGE 지표는 0.94로 매개변수 보정기간(2001~2004년)의 수문사상을 재현하기에 적절하였다. 다만, 매개변수 보정 기간 외에 대해서는 open loop모의의 정확도가 KGE 지표 기준 0.8 이하로 저하되었으며, 이에 자료동화 기법을 적용하여 모의 정확도 변동의 정도를 평가하였다. 앙상블 칼만 필터와 파티클 필터의 성능 비교 시, 모형 상태량 중 토양 수분 저류량(production store), 하도 저류량(routing store), 단위유량도(unit-hydro graph)를 모두 섭동 및 갱신하는 조건을 적용하였다(GR5J의 수문 요소에 대해서는 Fig. 2 참조). 이는 Piazzi et al. (2021)의 상태량 및 매개변수 동시 갱신 기법의 최적 구성 조건 연구결과, 가장 성능이 높았던 구성 조건을 적용한 것이다. 모의 앙상블의 수는 두 자료동화 기법에 대해 각각 100 개를 적용하였다. 기존 연구(Noh et al., 2011a)에서 파티클 필터의 경우 앙상블 수에 따라 모의 성능이 달라지는 현상이 보고되어, 본 연구에서도 앙상블 수에 따른 모의 정확도 변동을 비교하였으나, 큰 변화가 없어 결과 분석에는 포함하지 않았다.

자료동화 방법과 관련된 하이퍼-매개변수로는 강수 및 잠재 증발산 강제력 앙상블 생성, 섭동 및 갱신 대상 상태량 선정 등의 영향을 모의 실험을 통해 검토하였다. 시간 상관 가중치는 Piazzi et al. (2021)에서 제시한 값(시간 역상관 길이: 강수 1(day), 잠재 증발산 2(day))을 적용하였다. 또한, 섭동 및 갱신 대상 상태량의 선택에 따른 모형 정확도의 변화를 모의하여 비교하였다. 하이퍼-매개변수의 민감도 모의 결과는 3.2절에서 설명한다.

Table 2.

Modeling periods for parameter calibration, DA calibration and validation

Parameter calibration	Data assimilation calibration	Data assimilation validation
2001~2004	2005~2008	2009~2014

Table 3.

Optimized parameter values of GR5J

Study area	Objective function	Optimized parameter values
Study area	Objective function	X1	X2	X3	X4	X5
Yongdam	KGE	88	1.61	37	0.9	0.41

3. 결과 및 분석

3.1 수문자료동화 모의 결과 비교

자료동화 검정기간(2009~2014년)에 대해 앙상블 칼만 필터, 파티클 필터를 적용하여 모의된 일 유량 시계열을 관측 자료와 연도별로 비교하였다(Fig. 4). 자료동화 결과는 자료동화 보정 기간(2005~2008년)에 대해 민감도 평가 후, 각각의 자료동화 기법에 대한 최적 하이퍼-매개변수를 선정하여, 검정 기간에 대해 모의하는 방식으로 계산되었다. Fig. 4에서 앙상블 칼만 필터(EnKF), 파티클 필터(PF)는 각각 초록색, 파란색 실선으로, Open Loop (OL), 관측 유량(Obs)은 각각 빨간색 실선, 검은색 점선으로 표시하였다. 자료동화 유무에 따른 일 유량 모의 결과는 비홍수기에 현격한 차이를 보였기 때문에 y축을 대수축으로 도시하였다. 모든 연도별 일 유량 시계열 그래프에서 파티클 필터에 의한 모의 시계열이 관측과 가장 잘 일치하는 것을 시각적으로 확인할 수 있으며, 그 다음으로 앙상블 칼만 필터, open loop의 순으로 관측과의 차이가 늘어났다. 특히, warm-up 기간(2008년)을 추가로 모의하였음에도 불구하고, 2009 년 3~6 월 기간의 open loop 모의 결과(Fig. 4(a))는 자료동화 결과에 비해 관측 값과 차이가 크게 나타나는 것을 확인할 수 있었다. 이는 매개변수 보정에도 불구하고, 기상수문 조건에 따라 수문 모의의 정확도가 급격히 감소하는 사례이며, 자료동화를 적용하면 관측자료 정보를 실시간으로 반영하여 수문모형 해석의 한계를 보완할 수 있음을 보여준다.

Fig. 4에서 파티클 필터에 의한 모의 시계열은 연중 1~3 월, 10 월~12 월에서 관측 값과 가장 가까운 결과를 보이며, 반면 앙상블 칼만 필터는 open loop 보다 일부 향상되었으나 파티클 필터에 비해 정확도 향상이 제한적이었다. 이는 Table 4에서 KGE기준으로 모의 정확도를 정량 평가한 결과에서도 확인된다. 자료동화 보정 기간에 대한 KGE값은 open loop, 앙상블 칼만 필터, 파티클 필터에서 각각 0.736, 0.798, 0.885 이며, 자료동화 검정 기간에 대한 KGE값은 각각 0.799, 0.826, 0.933으로 분석되었다. 즉, 자료동화를 적용한 KGE값이 자료동화 적용 전보다 높게, 특히 자료동화 기법 중에서는 파티클 필터를 적용했을 때 KGE 값이 가장 높은 것으로 분석되었다.

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F4.jpg

Fig. 4.

Daily time-series of observed and simulated discharge from 2009 to 2014

Table 4.

Comparison of performance measure (KGE) of the open loop, the ensemble Kalman filter, and the particle filter with optimized hyper-parameters for DA calibration (2005~2008) and validation (2009~2014) periods

Simulations (hyper-parameters)	Data assimilation calibration (2005~2008)	Data assimilation validation (2009~2014)
Open loop	0.736	0.799
EnKF ( $ε_{p}$ = 0.4, $ε_{e}$ = 0.3)	0.798	0.826
PF ( $ε_{p}$ = 0.1, $ε_{e}$ = 0.2)	0.885	0.933

Fig. 5에서는 관측 및 모의 일 유량의 1:1그래프를 통해 앙상블 칼만 필터, 파티클 필터open loop의 성능을 비교하였다. Fig. 5(a)는 자료동화 보정기간, Fig. 5(b)는 자료동화 검정기간에 대한 관측 모의 1:1 그래프 이다. 초록색은 앙상블 칼만 필터, 파란색은 파티클 필터, 빨간색은 open loop를 나타낸다. Fig. 5(a)에서 open loop는 약 50~90 mm/day 모의 유량 사이에서, Fig. 5(b)에서 open loop는 약 40~200 mm/day 모의 유량 사이에서 관측 모의 일대일 선과 가장 떨어져 있는 것을 확인할 수 있다. 반면 앙상블 칼만 필터와 파티클 필터는 open loop의 불안정한 모의 성능을 자료동화를 통해 보완하여 관측과 더 가깝게 모의된 것을 볼 수 있다.

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F5.jpg

Fig. 5.

Observed and simulated discharges via the open loop (red dot), the ensemble Kalman filter (green dot), and the particle filter (blue dot)

3.2자료동화 하이퍼-매개변수 영향 분석

자료동화 하이퍼-매개변수 중 강수 및 잠재 증발산의 불확실성 범위가 모의 성능에 미치는 영향을 분석하였다. 구체적으로는, 기상 강제력 중 강수( $ε_{p}$ ) 및 잠재 증발산( $ε_{e}$ )의 오차 범위 관련 매개변수를 0.1~1 사이에서, 0.1 단위로 변경하며 모의 정확도가 어떻게 변화되는지 분석하였다(Fig. 6). Fig. 6(a)는 앙상블 칼만 필터에 대하여 강수 및 잠재 증발 오차 매개변수 별 KGE값을 분석한 결과를, Fig. 6(b)는 파티클 필터에 대하여 분석한 결과를 나타낸다. 앙상블 칼만 필터 모의 결과 $ε_{p}$ = 0.4, $ε_{e}$ = 0.3 일 때 KGE 값이 0.798로 가장 높았고, 파티클 필터로 모의 결과 $ε_{p}$ = 0.1, $ε_{e}$ = 0.2 일 때 KGE 값이 0.885 로 가장 높았다. 앙상블 칼만 필터와 파티클 필터 모두 기상 강제력 중에서도 강우 오차의 매개변수가 상대적으로 작을수록 KGE지표가 향상되는 것으로 분석되었다. 반면 자료동화한 결과 KGE값이 가장 낮을 때는 앙상블 칼만 필터가 $ε_{p}$ = 1.0, $ε_{e}$ = 0.8 일 때 KGE 값이 0.744 로 가장 낮은 것을, 파티클 필터가 $ε_{p}$ = 1.0, $ε_{e}$ = 0.4 일 때 KGE값이 0.801 로 가장 낮은 것을 확인 할 수 있다.

두 자료동화기법을 비교하면, 앙상블 칼만 필터와 파티클 필터 모의 결과의 KGE지표 범위는 각각 0.744~0.798, 0.801~0.885로 파티클 필터의 모의 성능이 모든 기상 강제력 오차 매개변수 모의 결과에서 앙상블 칼만 필터보다 우수하였다. 앙상블 칼만 필터는 최적 강수 및 잠재 증발 오차 매개변수의 범위가 전반적으로 넓었지만, 파티클 필터는 상대적으로 좁았다. 즉, 앙상블 칼만 필터는 파티클 필터 기법에 비해 상대적으로 강수 및 잠재 증발 오차 매개변수의 영향을 적게 받으며, 파티클 필터 기법은 강수 오차 매개변수가 작을수록 모의 정확도가 향상되는 경향이 있는 것으로 분석되었다.

본 연구에서 비교한 두가지 순차 자료동화 기법은 앙상블에 기반한 방법이므로, 모의된 앙상블의 확률 분포에 대한 분석이 필요하다. 하지만, KGE 지표는 모의 평균값만을 관측과 비교하기 때문에 모의 앙상블 통계적 특성을 분석하기 어렵다. 이에 앞서 Fig. 6에서 분석한 기상 강제력 오차 매개변수 민감도 모의 결과 중 일부 사례에 대해 랭크히스토그램(rank histogram)을 이용하여 앙상블 확률 분포의 특징을 분석하였다(Fig. 7). Fig. 7의 랭크히스토그램은 자료동화 보정기간(2005~2008년)에 대해 앙상블 칼만 필터는 0.3, 파티클 필터는 0.2로 잠재 증발 오차 매개변수를 고정하고, 강수 오차 매개변수를 0.1, 0.3, 0.5, 0.7, 0.9 간격으로 확률적 섭동을 다르게 설정하여 모의한 앙상블을 비교한 것이다. 이는 민감도 분석(Fig. 6)에서 모의 성능이 가장 우수했을 때의 잠재 증발산 오차 매개변수를 기준으로, 매개변수 변화에 따라 모의 정확도의 변동성이 큰 강수 오차 매개변수에 대한 변화를 분석한 것이다. 랭크히스토그램은 관측 값이 앙상블 확률 분포 범위에서 얼마나 떨어져 있는지에 따라 순서를 매긴다. 예를 들어, 100개의 모형으로 구성된 앙상블의 경우, 모의 유량의 확률 분포 범위에서 관측 값의 상대적 위치에 따라 1~101위의 순서를 매기는데, 1위는 앙상블의 최소값보다 관측값이 작은 것이고, 반대로 101위는 해당 시간 앙상블의 최대값보다 관측이 더 큰 경우이다. 이와 같은 방식으로 모든 관측 시계열에 대해 순위를 계산하여 각 순위별 도수(counts)의 분포를 도시한 것이 랭크히스토그램이다. 앙상블 기반 모의는 랭크히스토그램을 통해 모의된 확률 분포 범위의 적절성 여부를 평가할 수 있다. 가장 이상적인 모의 앙상블은 균등 분포에 가깝고, 앙상블 범위가 과대 추정되면 정규 분포와 비슷하게 가운데로 치우치게 되며, 앙상블 범위가 과소 추정되면 양 끝 부분의 도수가 늘어난다.

Fig. 7에서 왼쪽 패널의 칼만 필터의 앙상블은 전반적으로 과대 추정된 분포로, 오른쪽 패널의 파티클 필터는 양 끝 쪽의 도수가 대부분 높아 과소 추정된 분포인 것으로 분석되었다. 강수 강제력 오차 매개변수를 확률적으로 섭동의 크기가 증가하면, 앙상블 칼만 필터는 랭크히스토그램 분포의 가운데 부분 도수가 증가하며 과대 추정 양상이 확대되지만, 파티클 필터는 앙상블 확률 분포 범위에 큰 영향을 받지 않았다. 파티클 필터에서는 사전 분포에서 모의 앙상블의 범위가 넓게 계산되더라도 리샘플링 과정에서 우도가 큰 앙상블 위주로 사후 분포가 생성되었기 때문으로 추정된다.

한편, 랭크히스토그램으로부터 평가된 앙상블 분포의 적정성과 앙상블 평균의 모의 정확도와의 상관관계를 찾기는 어려웠다. 일례로, 파티클 필터 기법은 $ε_{p}$ = 0.1, $ε_{e}$ = 0.2 일 때 KGE값이 0.885 로 가장 높았지만(Fig. 6(b)), 이때 랭크히스토그램의 앙상블 확률 분포 범위(Fig. 7(b))는 과소추정에 가까우며, 파티클 필터의 다른 기상 강제력 하이퍼-매개변수를 적용한 결과와 분포 양상이 크게 다르지 않았다. 앙상블 칼만 필터에서도 $ε_{p}$ = 0.1, $ε_{e}$ = 0.3 일 때 KGE값이 0.79 로 높게 산정 되었지만(Fig. 6(a)), 이때의 앙상블 확률 분포 범위(Fig. 7(a))는 과대추정 영역이었다.

Table 5는 섭동 및 갱신 대상 상태량을 다르게 선택한 경우의 모의 성능을 분석하였다. 자료동화 보정기간(2005~2008년)에 대해 상태량, 즉 토양 수분 저류량(production store, P), 하도 저류량(routing store, R), 단위유량도(unit-hydro graph, UH)의 종류를 다르게 섭동 및 갱신한 KGE 결과 값이다(Table 5). 앙상블 칼만 필터는 상태량 별 종류를 변경하면서 모두 섭동 및 갱신하였다. 다만, 파티클 필터는 리샘플링 과정에서 상태량을 모두 갱신하기 때문에 P, R, UH 모두를 갱신하였고, 섭동은 상태량 별 종류를 변경하면서 모의 하였다. 상태량 별 모의 결과, 앙상블 칼만 필터는 P, R, UH를 모두 섭동 및 갱신 하였을 때 KGE값이 가장 높은 결과인 0.798, 파티클 필터 또한 P, R, UH를 모두 섭동 및 갱신 하였을 때 KGE값이 가장 높은 결과인 0.885로 분석되었다. 이로부터 일 유량에 대한 자료동화 기법 적용시에는 토양 수분 저류량, 하도 저류량, 단위 유량도 등 대부분의 수문과정에 대한 섭동 및 갱신을 모두 고려하는 것이 적절한 것으로 판단된다.

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F6.jpg

Fig. 6.

Comparison of performance measure (KGE) according to varying weather forcing uncertainty parameters

https://cdn.apub.kr/journalsite/sites/kwra/2022-055-10/N0200551005/images/kwra_55_10_05_F7.jpg

Fig. 7.

Rank histogram plots of (a) the ensemble Kalman filter and (b) the particle filter with varying weather forcing uncertainty parameters (selective cases shown in Fig. 6)

Table 5.

Comparison of performance measure (KGE) with varying perturbed (and updated) state variables (P: production store, R: routing store, UH: unit hydrograph)

DA methods	Perturbed (and updated) state variables
DA methods	P, R, UH	P, R	P
EnKF ( $ε_{p}$ = 0.4, $ε_{e}$ = 0.3)	0.798	0.791	0.747
PF ( $ε_{p}$ = 0.1, $ε_{e}$ = 0.2)	0.885	0.879	0.862

4. 결 론

본 연구에서는 용담댐 유역에 대해 앙상블 칼만 필터, 파티클 필터 기법을 일 유출 모의에 적용하여 두 기법을 차이를 비교하고, 자료동화 방법과 관련된 하이퍼-매개변수의 불확실성이 순차 자료동화 기법의 모의 성능에 미치는 영향을 분석하였다.

본 연구의 결과를 요약하면 다음과 같다.

1) 자료동화 적용 전 수문 모의 결과, GR5J는 KGE값이 0.799이며, 자료동화 적용 후 앙상블 칼만 필터는 KGE값이 0.826, 파티클 필터는 KGE 값이 0.933으로 자료동화를 적용 후 모의성능이 개선되었다. 또한 두 자료동화 기법 중 파티클 필터에 의한 모의가 정확성이 높은 것으로 분석되었다.

2) 기상 강제력 중 강수 및 잠재 증발산 오차 매개변수의 불확실성 영향을 검토한 결과, 앙상블 칼만 필터는 $ε_{p}$ = 0.4, $ε_{e}$ = 0.3 일 때, 파티클 필터는 $ε_{p}$ = 0.1, $ε_{e}$ = 0.2 일 때 KGE 값이 0.798, 0.885 로 가장 높았다. 이는 자료동화 기법에 따라 최적 성능을 얻을 수 있는 불확실성 조건이 다른 것을 의미하며, 파티클 필터를 이용하여 일 유출을 모의할 경우, 과도한 기상 강제력 오차를 설정하지 않는 것이 권고된다.

3) 랭크히스토그램을 이용하여 앙상블 확률 분포의 범위를 비교한 결과, 전반적으로 앙상블 칼만 필터는 과대 추정 분포였고, 파티클 필터는 과소 추정 분포였다. 또한 강수 오차 매개변수 변화에 따른 앙상블의 확률 분포의 변동을 비교한 결과, 앙상블 분포의 적정성과 앙상블 평균의 모의 정확도와의 상관관계는 찾을 수 없었다.

4) 상태량 별 대상을 다른 조건으로 섭동 및 갱신했을 때의 예측 성능을 분석하였다. 분석 결과, 두가지 자료동화기법 모두에서 저류량(S), 하도 저류량(R), 단위유량도(UH)를 모두 섭동 및 갱신하였을 때 KGE값이 가장 높았다.

본 연구의 모의 실험 결과는 앙상블 자료동화를 적용하여 일 유출 모의의 정확도 향상이 가능하며, 최적 성능을 발휘하기 위해서는 수문자료동화 기법 별 하이퍼-매개변수의 적정한 조정이 필요함을 함의한다. 오픈 소스 소프트웨어인 airGRdatassim은 집중형 모형에 의한 결정론적 모의(deterministic predictions)를 앙상블 기반 확률 모의로 확장하고 순차 자료동화 기술에 대한 활용성을 향상시킬 수 있는 유용한 도구로 평가된다.

자료동화 적용 시 홍수기 보다 갈수 기간에 대한 수문 모의 정확도 향상 효과가 뛰어났는데, 이는 수문 반응시간이 상대적으로 짧은 대상 유역에서 일 유량 모의 시 상태량 갱신으로 인한 홍수기 예측 정확도 향상이 제한적인 것인 이유로 추정되나, 후속 연구가 필요한 부분이다. 또한, 본 고의 일 유량 앙상블 예측 기술을 확장하여 시간 단위 이하 유량 예측에 대한 자료동화 적용 및 불확실성 분석에 대한 연구가 가능할 것으로 판단된다.

Acknowledgements

이 연구는 금오공과대학교 학술연구비로 지원되었음(과제 번호: 2019104170).

Conflicts of Interest

The authors declare no conflict of interest.

References

Adeyeri, O.E., Laux, P., Arnault, J., Lawin, A.E., and Kunstmann, H. (2020). "Conceptual hydrological model calibration using multi-objective optimization techniques over the transboundary Komadugu-Yobe basin, Lake Chad Area, West Africa." Journal of Hydrology: Regional Studies, Vol. 27, 100655. 10.1016/j.ejrh.2019.100655

Arulampalam, M.S., Maskell, S., Gordon, N., and Clapp, T. (2002). "A tutorial on particle filters for online nonlinear/non-Gaussian Bayesian tracking." IEEE Transactions on Signal Processing, Vol. 50, No. 2, pp. 174-188. 10.1109/78.978374

Blöschl, G., Bierkens, M.F.P., Chambel, A., Cudennec, C., Destouni, G., Fiori, A., Kirchner, J.W., McDonnell, J.J., Savenije, H.H.G., Sivapalan, M. et al. (2019). "Twenty-three unsolved problems in hydrology (UPH) - a community perspective." Hydrological Sciences Journal, Vol. 64, No. 10, pp. 1141-1158. 10.1080/02626667.2019.1620507

Boucher, M.-A., Quilty, J., and Adamowski, J. (2020). "Data assimilation for streamflow forecasting using extreme learning machines and multilayer perceptrons." Water Resources Research, Vol. 56, No. 6, e2019WR026226. 10.1029/2019WR026226

Choi, J.-H., and Kim, S.-D. (2021). "Estimating time-varying parameters for monthly water balance model using particle filter: assimilation of stream flow data." Journal of Korea Water Resources Association, Vol. 54, No. 6, pp. 365-379.

Clark, M.P., Rupp, D.E., Woods, R.A., Zheng, X., Ibbitt, R.P., Slater, A.G., Schmidt, J., and Uddstrom, M.J. (2008). "Hydrological data assimilation with the ensemble Kalman filter: Use of streamflow observations to update states in a distributed hydrological model." Advances in Water Resources, Vol. 31, No. 10, pp. 1309-1324. 10.1016/j.advwatres.2008.06.005

Doucet, A., Godsill, S., and Andrieu, C. (2000). "On sequential Monte Carlo sampling methods for Bayesian filtering." Statistics and Computing, Vol. 10, No. 3, pp. 197-208. 10.1023/A:1008935410038

Evensen, G. (1994). "Sequential data assimilation with a nonlinear quasi-geostrophic model using Monte Carlo methods to forecast error statistics." Journal of Geophysical Research: Oceans, Vol. 99, No. C5, pp. 10143-10162. 10.1029/94JC00572

Evensen, G. (2003). "The Ensemble Kalman Filter: theoretical formulation and practical implementation." Ocean Dynamics, Vol. 53, No. 4, pp. 343-367. 10.1007/s10236-003-0036-9

Gupta, H.V., Kling, H., Yilmaz, K.K., and Martinez, G.F. (2009). "Decomposition of the mean squared error and NSE performance criteria: Implications for improving hydrological modelling." Journal of Hydrology, Vol. 377, No. 1-2, pp. 80-91. 10.1016/j.jhydrol.2009.08.003

Hendricks Franssen, H.J., and Kinzelbach, W. (2008). "Real-time groundwater flow modeling with the Ensemble Kalman Filter: Joint estimation of states and parameters and the filter inbreeding problem." Water Resources Research, Vol. 44, No. 9, W09408. doi: 10.1029/2007WR006505. 10.1029/2007WR006505

Kim, Y.S., Lee, G.H., Lee, D.E., and Noh, S.J. (2015). "Parameter estimation and uncertainty assessment of a soil erosion model using data assimilation method." Journal of Korean Society of Hazard Mitigation, Vol. 15, No. 6, pp. 373-382. 10.9798/KOSHAM.2015.15.6.373

Knoben, W.J.M., Freer, J.E., and Woods, R.A. (2019). "Technical note: Inherent benchmark or not? Comparing Nash-Sutcliffe and Kling-Gupta efficiency scores." Hydrology and Earth System Sciences, Vol. 23, No. 10, pp. 4323-4331. 10.5194/hess-23-4323-2019

Lavenne, A., Thirel, G., Andréassian, V., Perrin, C., and Ramos, M.-H. (2016). "Spatial variability of the parameters of a semi-distributed hydrological model." Proceedings of the International Association of Hydrological Sciences, Vol. 373, pp. 87-94. 10.5194/piahs-373-87-2016

Leach, J.M., and Coulibaly, P. (2019). "An extension of data assimilation into the short-term hydrologic forecast for improved prediction reliability." Advances in Water Resources, Vol. 134, 103443. 10.1016/j.advwatres.2019.103443

Lee, B.J., and Bae, D.-H. (2011). "Development of real-time river flow forecasting model with data assimilation technique." Journal of Korea Water Resources Association, Vol. 44, No. 3, pp. 199-208. 10.3741/JKWRA.2011.44.3.199

Lee, B.J., Jung, I.-W., Jeong, H.-S., and Bae, D.-H. (2013). "Development of realtime dam's hydrologic variables prediction model using observed data assimilation and reservoir operation techniques." Journal of Korea Water Resources Association, Vol. 46, No. 7, pp. 755-765. 10.3741/JKWRA.2013.46.7.755

Lee, D.U., Kim, Y.S., Yu, W.S., and Lee, G.H. (2017). "Evaluation on applicability of on/off-line parameter calibration techniques in rainfall-runoff modeling." Journal of Korea Water Resources Association, Vol. 50, No. 4, pp. 241-252. 10.3741/JKWRA.2017.50.4.241

Liu, Y., Weerts, A.H., Clark, M., Hendricks Franssen, H.-J., Kumar, S., Moradkhani, H., Seo, D.-J., Schwanenberg, D., Smith, P., van Dijk, A.I.J.M., van Velzen, N., He, M., Lee, H., Noh, S.J., Rakovec, O., and Restrepo, P. (2012). "Advancing data assimilation in operational hydrologic forecasting: Progresses, challenges, and emerging opportunities." Hydrology and Earth System Sciences, Vol. 16, No. 10, pp. 3863-3887. 10.5194/hess-16-3863-2012

Le Moine, N. (2008). Le bassin versant de surface vu par le souterrain: Une voie d'amélioration des performances et du réalisme des modèles pluie-débit?. Ph. D. Dissertation, Université Pierre et Marie Curie Paris VI, Paris, France, pp. 149-152.

Noh, S.J. (2013). Sequential Monte Carlo methods for probabilistic forecasts and uncertainty assessment in hydrologic modeling, Ph. D. Dissertation, Kyoto University, Kyoto, Japan.

Noh, S.J., Rakovec, O., Weerts, A.H., and Tachikawa, Y. (2014). "On noise specification in data assimilation schemes for improved flood forecasting using distributed hydrological models." Journal of Hydrology, Vol. 519, pp. 2707-2721. 10.1016/j.jhydrol.2014.07.049

Noh, S.J., Tachikawa, Y., Shiiba, M., and Kim, S. (2011a). "Dual state-parameter updating scheme on a conceptual hydrologic model using sequential Monte Carlo Filters." Journal of Japan Society of Civil Engineers, Ser. B1 (Hydraulic Engineering), Vol. 67, No. 4, p. I_1-I_6. 10.2208/jscejhe.67.I_1

Noh, S.J., Tachikawa, Y., Shiiba, M., and Kim, S. (2011b). "Applying sequential Monte Carlo methods into a distributed hydrologic model: lagged particle filtering approach with regularization." Hydrology and Earth System Sciences, Vol. 15, No. 10, pp. 3237-3251. 10.5194/hess-15-3237-2011

Noh, S.J., Tachikawa, Y., Shiiba, M., and Kim, S. (2012). "Ensemble Kalman Filtering and particle filtering in a lag-time window for short-term streamflow forecasting with a distributed hydrologic model." Journal of Hydrologic Engineering, Vol. 18, No. 12, pp. 1684-1696. 10.1061/(ASCE)HE.1943-5584.0000751

Oudin, L., Hervieu, F., Michel, C., Perrin, C., Andréassian, V., Anctil, F., and Loumagne, C. (2005). "Which potential evapotranspiration input for a lumped rainfall - runoff model?: Part 2 - Towards a simple and efficient potential evapotranspiration model for rainfall - runoff modelling." Journal of Hydrology, Vol. 303, No. 1, pp. 290-306. 10.1016/j.jhydrol.2004.08.026

Oudin, L., Moulin, L., Bendjoudi, H., and Ribstein, P. (2010). "Estimating potential evapotranspiration without continuous daily data: possible errors and impact on water balance simulations." Hydrological Sciences Journal, Vol. 55, No. 2, pp. 209-222. 10.1080/02626660903546118

Piazzi, G., Thirel, G., Perrin, C., and Delaigue, O. (2021). "Sequential data assimilation for streamflow forecasting: Assessing the sensitivity to uncertainties and updated variables of a conceptual hydrological model at basin scale." Water Resources Research, Vol. 57, No. 4, e2020WR028390. doi: 10.1029/2020WR028390 10.1029/2020WR028390

Ristic, B., Arulampalam, S., and Gordon, N. (2004). Beyond the Kalman Filter: Particle filters for tracking applications. Artech House, Boston, MA, U.S. and London, UK. 10.1155/S1110865704405095

Shen, H., Seo, D.-J., Lee, H., Liu, Y., and Noh, S. (2022). "Improving flood forecasting using conditional bias-aware assimilation of streamflow observations and dynamic assessment of flow-dependent information content." Journal of Hydrology, Vol. 605, 127247. 10.1016/j.jhydrol.2021.127247

Yoo, C., Hwang, J.-H., and Kim, J. (2012). "Use of the extended Kalman Filter for the real-time quality improvement of runoff data: 1. Algorithm construction and application to one station." Journal of Korea Water Resources Association, Vol. 45, No. 7, pp. 697-711. 10.3741/JKWRA.2012.45.7.697

Journal of Korea Water Resources Association ISSN:2799-8746(Print) 2799-8754(Online) 한국수자원학회 논문집

Preview

Comparative assessment and uncertainty analysis of ensemble-based hydrologic data assimilation using airGRdatassim

ABSTRACT

MAIN

(1a)

(1b)

Fig. 1.

Schematic diagram of data assimilation methods: (a) the ensemble Kalman filter and (b) the particle filter

(2)

(3)

(4a)

(4b)

(5)

(6a)

(6b)

(7a)

(7b)

(8)

Fig. 2.

Structure of GR5J model

Table 1.

Descriptions of GR5J model parameters (Adeyeri et al., 2020)

(9)

Fig. 3.

Study area and observation network: (a) Yongdam-dam catchment (b) Thiessen network of rainfall stations, and (c) Thiessen network of weather stations. Associated Thiessen polygons are represented in black lines

(10)

Table 2.

Modeling periods for parameter calibration, DA calibration and validation

Table 3.

Optimized parameter values of GR5J

Fig. 4.

Daily time-series of observed and simulated discharge from 2009 to 2014

Table 4.

Comparison of performance measure (KGE) of the open loop, the ensemble Kalman filter, and the particle filter with optimized hyper-parameters for DA calibration (2005~2008) and validation (2009~2014) periods

Fig. 5.

Observed and simulated discharges via the open loop (red dot), the ensemble Kalman filter (green dot), and the particle filter (blue dot)

Fig. 6.

Comparison of performance measure (KGE) according to varying weather forcing uncertainty parameters

Fig. 7.

Rank histogram plots of (a) the ensemble Kalman filter and (b) the particle filter with varying weather forcing uncertainty parameters (selective cases shown in Fig. 6)

Table 5.

Comparison of performance measure (KGE) with varying perturbed (and updated) state variables (P: production store, R: routing store, UH: unit hydrograph)

Acknowledgements

Conflicts of Interest

References