Application of minimum relative entropy theory for streamflow forecasting

Abstract

This paper develops and applies the minimum relative entropy (MRE) theory with spectral power as a random variable for streamflow forecasting. The MRE theory consists of (1) hypothesizing a prior probability distribution for the random variable, (2) determining the spectral power distribution, (3) extending the autocorrelation function, and (4) doing forecasting. The MRE theory was verified using streamflow data from the Mississippi River watershed. The exponential distribution was chosen as a prior probability in applying the MRE theory by evaluating the historical data of the Mississippi River. If no prior information is given, the MRE theory is equivalent to the Burg entropy (BE) theory. The spectral density obtained by the MRE theory led to higher resolution than did the BE theory. The MRE theory did not miss the largest peak at 1/12th frequency, which is the main periodicity of streamflow of the Mississippi River, but the BE theory sometimes did. The MRE theory was found to be capable of forecasting monthly streamflow with a lead time from 12 to 48 months. The coefficient of determination (r ²) between observed and forecasted stream flows was 0.912 for Upper Mississippi River and was 0.855 for Lower Mississippi River. Both MRE and BE theories were generally more reliable and had longer forecasting lead times than the autoregressive (AR) method. The forecasting lead time for MRE and BE could be as long as 48–60 months, while it was less than 48 months for the AR method. However, BE was comparable to MRE only when observations fitted the AR process well. The MRE theory provided more reliable forecasts than did the BE theory, and the advantage of using MRE is more significant for downstream flows with irregular flow patterns or where the periodicity information is limited. The reliability of monthly streamflow forecasting was the highest for MRE, followed by BE followed by AR.

Appendices

Appendix 1: Posterior spectral density without prior information given

If no prior information is given, one can consider the uniform distribution as $p(\mathop x\limits^{ \to } ) = C$and let $A = p(\mathop x\limits^{ \to } )\exp [ - \lambda_{0} ]$ in Eq. (14), then it becomes

$$q(\mathop x\limits^{ \to } ) = A\exp \left[ { - \sum\limits_{r = - N}^{N} {\lambda_{r} } \sum\limits_{k = - n}^{n} {x_{k} c_{rk} } } \right]$$

(35)

To simplify Eq. (35), let $u_{k} = \sum\limits_{r = - N}^{N} {\lambda_{r} c_{rk} }$. Then, Eq. (35) becomes

$$q(\mathop x\limits^{ \to } ) = A\exp \left[ { - \sum\limits_{k = - n}^{n} {x_{k} \sum\limits_{r = - m}^{m} {\lambda_{r} } c_{rk} } } \right] = A\exp \left[ { - \sum\limits_{k = - n}^{n} {u_{k} x_{k} } } \right]$$

(36)

Integrating Eq. (36) over 0 to infinite yields

$$1 = \int\limits_{0}^{\infty } {q(\mathop x\limits^{ \to } )d} \mathop x\limits^{ \to } = \iiint {A\exp \left[ { - \sum\limits_{k = - n}^{n} {u_{k} x_{k} } } \right]}dx_{k} = A\prod\limits_{k = 1}^{n} {\int {\exp ( - u_{k} x_{k} )} dx_{k} } = A\prod\limits_{k = 1}^{n} {\frac{1}{{u_{k} }}}$$

(37)

for $\int {\exp ( - u_{k} x_{k} )} dx_{k} = - \frac{1}{{u_{k} }}\exp ( - u_{k} x)\left| {_{0}^{\infty } } \right. = \frac{1}{{u_{k} }}$. Rearranging Eq. (37) yields

$$A = \prod\limits_{k} {u_{k} }$$

(38)

Thus, inputting Eq. (38) into Eq. (35), one obtains

$$q(\mathop x\limits^{ \to } ) = \prod\limits_{k = - n}^{n} {u_{k} \exp \left[ { - \sum\limits_{k = - n}^{n} {u_{k} x_{k} } } \right]} = \prod\limits_{k = - n}^{n} {u_{k} \exp [ - u_{k} x_{k} ]}$$

(39)

It is shown that Eq. (39) is a multivariate exponential distribution, where the probability of the spectral power at each frequency is $q(x_{k} ) = u_{k} \exp [ - u_{k} x_{k} ]$. Thus, the expected spectral power at each frequency can be obtained as

$$\overline{{{\text{x}}_{k} }} = \int {x_{k} q(} x_{k} )dx_{k} = \int {x_{k} } u_{k} \exp [ - u_{k} x_{k} ]dx_{k} = \frac{1}{{u_{k} }} = \frac{1}{{\sum\limits_{r = - N}^{N} {\lambda_{r} c_{rk} } }}$$

(40)

which is equivalent to the Burg entropy.

Appendix 2: Posterior spectral density with exponential prior distribution

If prior information points to an exponential distribution, with the expected spectrum power S _k for each frequency, then the prior becomes

$$p(\mathop x\limits^{ \to } ) = \prod\limits_{k} {\frac{1}{{S_{k} }}\exp \left( { - \frac{{x_{k} }}{{S_{k} }}} \right)}$$

(41)

Inputting Eq. (41) into the posterior distribution Eq. (14), the result becomes

$$\begin{aligned} q(\mathop x\limits^{ \to } ) & = p(\mathop x\limits^{ \to } )\exp \left[ { - \lambda_{0} - \sum\limits_{r = - N}^{N} {\lambda_{r} } \sum\limits_{k = - n}^{n} {x_{k} c_{rk} } } \right] \\ & = \exp ( - \lambda_{0} )\prod\limits_{k} {\frac{1}{{S_{k} }}\exp \left( { - \frac{{x_{k} }}{{S_{k} }}} \right)} \exp \left[ { - \sum\limits_{r = - N}^{N} {\lambda_{r} c_{rk} x_{k} } } \right] \\ & = \exp ( - \lambda_{0} )\prod\limits_{k} {\frac{1}{{S_{k} }}\exp \left( { - \left( {u_{k} + \frac{1}{{S_{k} }}} \right)x_{k} } \right)} \\ \end{aligned}$$

(42)

where $u_{k} = \sum\limits_{r = - N}^{N} {\lambda_{r} c_{rk} }$, the same as previously defined. It can be seen from Eq. (42) that the posterior distribution of spectral power is in the form that is transformable to the exponential distribution of $\eta \exp ( - \eta x)$, where $u_{k} + \frac{1}{{S_{k} }}$ can be considered as the exponential parameter of each distribution. Thus, λ₀ needs to satisfy that $\exp ( - \lambda_{0} )\prod\limits_{k} {\frac{1}{{S_{k} }}} = \prod\limits_{k} {u_{k} + \frac{1}{{S_{k} }}}$ so that Eq. (42) forms the multivariate exponential distribution as

$$q(\mathop x\limits^{ \to } ) = \prod\limits_{k = 1}^{n} {\left( {u_{k} + \frac{1}{{S_{k} }}} \right)\exp } \left[ { - \left( {u_{k} + \frac{1}{{S_{k} }}} \right)x_{k} } \right]$$

(43)

and the posterior distribution of each frequency is

$$q(x_{k} ) = \left( {u_{k} + \frac{1}{{S_{k} }}} \right)\exp \left[ { - \left( {u_{k} + \frac{1}{{S_{k} }}} \right)x_{k} } \right]$$

(44)

Thus, the expected posterior power spectral power at each frequency becomes

$$\overline{{x_{k} }} = \int {x_{k} q(} x_{k} )dx_{k} = \int {x_{k} } \left( {u_{k} + \frac{1}{{S_{k} }}} \right)\exp \left[ { - \left( {u_{k} + \frac{1}{{S_{k} }}} \right)x_{k} } \right]dx_{k} = \frac{1}{{\frac{1}{{S_{k} }} + \sum\limits_{r = - N}^{N} {\lambda_{r} c_{rk} } }}$$

(45)

Appendix 3: Linear extension of autocovariance

Recall that in time series, the AR process or the linear prediction process has the spectral density in the form of

$$T(f) = \frac{{\sigma^{2} }}{{\left| {1 + \sum\limits_{n = - N}^{N} {a_{n} z^{n} } } \right|^{2} }} = \frac{{\sigma^{2} }}{{\sum\limits_{n = 0}^{N} {a_{s} z^{n} } \sum\limits_{n = 0}^{N} {a_{n}^{*} z^{ - n} } }}$$

(46)

Taking the inverse Fourier transformation of Eq. (46), one obtains the autocovariance function like in Eq. (7) as

$$\begin{aligned}R_{r} = &\int {T(f)} \exp (2\pi ir\Delta tf) = \frac{{\sigma^{2} }}{2\pi i}\oint {\frac{{z^{ - r - 1} dz}}{{\sum\limits_{n = 0}^{N} {a_{n} z^{n} } \sum\limits_{n = 0}^{N} {a_{n}^{*} z^{ - n} } }}},\\ &- N < n < N \end{aligned}$$

(47)

Multiplying Eq. (47) by a ^* _r and summing up for n from 0 to N,

$$\sum\limits_{r = 0}^{N} {a_{r}^{*} R_{N - r} = \frac{{\sigma^{2} }}{2\pi i}\oint {\frac{{z^{r - 1} dz\sum\limits_{n = 0}^{N} {a_{n}^{*} z^{ - n} } }}{{\sum\limits_{n = 0}^{N} {a_{n} z^{n} } \sum\limits_{n = 0}^{N} {a_{n}^{*} z^{ - n} } }}} } = \frac{{\sigma^{2} }}{2\pi i}\oint {\frac{{z^{r - 1} dz}}{{\sum\limits_{n = 0}^{N} {a_{n} z^{n} } }}}$$

(480)

It is noted from Cauchy’s integral formula that

$$\oint {n\frac{{z^{r - 1} dz}}{{\sum\limits_{s = 0}^{N} {a_{n} z^{n} } }}} = \left\{ {\begin{array}{*{20}c} {2\pi i} \\ 0 \\ \end{array} } \right.\begin{array}{*{20}c} , & {r = 0} \\ , & {r \ge 1} \\ \end{array}$$

(49)

Thus, integrating of Eq. (49) using Cauchy’s formula, one obtains

$$\sum\limits_{r = 0}^{N} {a_{r}^{*} R_{N - r} = \frac{{\sigma^{2} }}{2\pi i}\oint {\frac{{z^{r - 1} dz}}{{\sum\limits_{n = 0}^{N} {a_{n} z^{n} } }}} = \frac{{\sigma^{2} }}{2\pi i}2\pi i = \sigma^{2} \begin{array}{*{20}c} , & {r = 0} \\ \end{array} }$$

(50a)

$$\sum\limits_{r = 0}^{N} {a_{r}^{*} R_{N - r} = \frac{{\sigma^{2} }}{2\pi i}\oint {\frac{{z^{r - 1} dz}}{{\sum\limits_{n = 0}^{N} {a_{n} z^{n} } }}} = \frac{{\sigma^{2} }}{2\pi i}0 = 0\begin{array}{*{20}c} , & {r \ge 1} \\ \end{array} }$$

(50b)

Taking the conjugate of Eq. (50a, b) yields

$$\sum\limits_{n = 0}^{N} {R_{r - n} a_{n} = \sigma^{2} \begin{array}{*{20}c} , & {r = 0} \\ \end{array} }$$

(51a)

$$\sum\limits_{n = 0}^{N} {R_{r - n} a_{n} = 0\begin{array}{*{20}c} , & {r \ge 1} \\ \end{array} }$$

(51b)

which can be written in the matrix form as:

$$\left[ {\begin{array}{*{20}c} {R_{0} } & {R_{ - 1} } & \cdots & {R_{ - N} } \\ {R_{1} } & {R_{0} } & {} & {R_{1 - N} } \\ \vdots & \vdots & {} & \vdots \\ {R_{N} } & {R_{N - 1} } & \cdots & {R_{0} } \\ \end{array} } \right]\left\{ {\begin{array}{*{20}c} 1 \\ {a_{1} } \\ \vdots \\ {a_{N} } \\ \end{array} } \right\} = \left\{ {\begin{array}{*{20}c} {\sigma^{2} } \\ 0 \\ \vdots \\ 0 \\ \end{array} } \right\}$$

(52)

Equation (52) is the Yule-Walker equation (Yule 1927), with coefficient a _n and gain σ ², and according to the linear extension, the N + 1th autocovariance can be obtained by

$$R_{N + 1} = a_{1} R_{N} + a_{2} R_{N - 1} + \cdots + a_{m} R_{N + m - 1}$$

(53)