HTML

This study provides an assessment of the soil moisture predictions based on SiB2 and EnKF utilizing the hourly insitu soil moisture observations from 11 August to 14 September 2007 (days 223–257) at Point A and from 1 December 2008 to 14 February 2009 (days 335–410) at Point B. Considering the short temporal coverage of the experiment and the observations, the meteorological data (from 22 June to 10 August 2007 and from 14 September to 30 November 2008) were used to spinup the model for 10 times at Points A and B, respectively, before undertaking the assimilation. The daily observed precipitation of the corresponding assimilation and prediction period is shown in Fig. 2, and averagely interpolated into each model running time. Once the observations are obtained, the observations can be assimilated into the model by EnKF to improve the soil moisture simulations, so the assimilation procedure can be conducted according to Fig. 3. The simulated, assimilated, and predicted soil moistures were assessed at the soil depths of 5, 50, and 100 cm at Point A, and 5, 30, and 100 cm at Point B.
Figure 2. Precipitation at Points A and B in Meilin watershed, eastern China: (a) 11 August to 14 September 2007 (days 223–257) at Point A and (b) 1 December 2008 to 14 February 2009 (days 335–410) at Point B.
The time interval ∆t of the soil moisture simulation is set to be 1 h, the assimilation frequencies are set to be every one hour (hourly), every three hours, and once a day (daily). During the prediction period, the model runs without assimilation, which is shown in Fig. 3. Thus, the assimilated soil moisture, with different assimilation frequencies, at the end of assimilation period is set to be the initial state value for the soil moisture prediction period. Additionally, the precipitation in the prediction period uses the measured precipitation, the randomly generated precipitation [i.e., measured precipitation value multiplied by a random value generated from the uniform distribution U (0, 2)], and zero precipitation. There are sparse in situ soil moisture observations, and remote soil moisture is limited to the soil surface, so setting the observed soil moisture to be the initial values is only used to test whether the model is stable at the beginning of prediction period. Based on the initial value and precipitation data, six cases for Point A are considered as listed in Table 1, and only four cases for Point B are considered as listed in Table 2.
Case Initial value for prediction period Precipitation in prediction period Assimilation frequency (h) Soil moisture result 0 Measured value Measured value Fig. 4 1 Assimilated_1h Measured value 1 Figs. 5, 8 2 Assimilated_3h Measured value 3 Fig. 9 3 Assimilated_1d Measured value 24 Fig. 10 4 Assimilated_1h Measured × U (0, 2) 1 Fig. 12 5 Assimilated_1h 0 1 Fig. 13 Note: Assimilated_1h, Assimilated_3h, and Assimilated_1d represent cases with an assimilation frequency of 1 h (hourly), 3 h, and 1 day (daily), respectively. Case 0 is used to test the model stability in the prediction period. Table 1. Case studies designed for soil moisture predictions at Point A in Meilin watershed

In this study, the SiB2 model (Sellers et al., 1996) is adopted as the model operator, which can simulate the water–energy transfer process in land–atmosphere interaction. To simulate the soil moisture, three soil layers in the soil column are considered, which are the soil surface, root zone, and deep soil layers. The node for soil surface, root zone, and deep soil layers are set to be 5, 50 or 30, and 100 cm, respectively, as the soil moisture observation collected at the corresponding depths. The governing equation in each soil layer is as follows:
$$ {\rm Soil \; surface \; layer: }\;\; \frac{{\partial {W_1}}}{{\partial t}} = \frac{1}{{{\theta _{\rm{s}}}{D_1}}}\left[ {P  {Q_{1, 2}}  \frac{1}{\rho }{E_{\rm{g}}}} \right]; \;$$ (1) $${\rm Root \; zone \; layer:}\;\;\frac{{\partial {W_2}}}{{\partial t}} = \frac{1}{{{\theta _{\rm{s}}}{D_2}}}\left[ {{Q_{1, 2}}  {Q_{2, 3}}  \frac{1}{\rho }{E_{\rm{c}}}} \right];$$ (2) $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!{\rm Deep \; soil \; layer:}\;\;\frac{{\partial {W_3}}}{{\partial t}} = \frac{1}{{{\theta _{\rm{s}}}{D_3}}}\left[ {{Q_{2, 3}}  {Q_3}} \right]. $$ (3) Here, W_{1}, W_{2}, and W_{3} are the soil moisture wetness (state variables) at each soil layer (m^{3} m^{–3}); θ_{s} is the saturated soil water content (m^{3} m^{–3}); D_{i} (i = 1, 2, 3) is the ith soil layer thickness (m); Q_{i,i+1} (i = 1, 2) is the water flow between layers i and i + 1 (m s^{–1}); Q_{3} is the gravitational drainage from the recharge soil moisture store (m s^{–1}); P is the precipitation infiltration into the soil surface layer (m s^{–1}); E_{g} is the evaporation rate from the soil surface layer (m s^{–1}); E_{c} is the canopy transpiration rate (m s^{–1}); and ρ is the water density.
In this land surface model, the input variables include the precipitation, air temperature, relative humidity, wind speed, shortwave radiation, and longwave radiation. The longwave radiation is calculated by using a standard method, i.e., Brunt’s equation (Monteith, 1973), in SiB2.

The EnKF (Evensen, 1994; Heemink et al., 2001; Vrugt et al., 2006; Weerts and Serafy, 2006; Xie and Zhang, 2010; Dumedah and Coulibaly, 2013; Liu H. R. et al., 2016; Fu et al., 2018b) has achieved a great success and is frequently applied to improve the simulation accuracy in the fields of atmosphere, hydrology, and ocean science (Whitaker and Hamill, 2002). The EnKF includes a state function M and an observation function H. M is the aforementioned SiB2 model for simulating soil moisture at different soil layers, and H is used to link the soil moisture simulations to observations in this study as follows:
$${ X_k} = { M}\left({{{ X}_{k  1}}, {{ V}_{k  1}}} \right), $$ (4) $${{ Y}_k} = { H}\left({{{ X}_k}, {{ U}_k}} \right), \quad\quad$$ (5) where X_{k} is the simulated value of the state variable at time k, and X_{k} = [W_{1}, W_{2}, W_{3}]
$ _k^{\rm {T}} $ is the soil moisture content at each soil layer at time k in this study; V_{k–1} and U_{k–1} are the model and observation errors that conform to Gaussian distributions; M() is the model operator (i.e., SiB2 model).(1) Forecasting: The initial ensembles
$ { X}_{0,i} $ (i=1, …, N) are generated by adding random perturbations to${ X} $ _{0} at the initial time and are used to forecast the state ensemble based on the model operator:$${{ X}_{k, i}} = { M}\left({{{ X}_{k  1, i}}, {{ V}_{k  1, i}}} \right)\;\;\;\;{{ V}_{k, i}} \sim N\left({0, { Q}} \right), $$ (6) where X_{k,i} is the forecasted soil moisture for each member at time k; N is the ensemble size, and is set to 80 in this study, as Moradkhani et al. (2005) and Vrugt et al. (2006) pointed out that the assimilated results tend to be stable when the ensemble size is larger than 50.
(2) Updating: For updating the state variable (soil moisture), the observed soil moisture is assimilated into the model by using the following equation:
$${ X}_{k, i}^{\rm up} = {{ X}_{k, i}} + {{ K}_k}\left[ {{{ Y}_{k, i}}  {{H}}\left({ {{\bar{ X}}_k} } \right)} \right].$$ (7) The Kalman gain K_{k} and generated measurements ensemble Y_{k,i} at time k (Weerts and Serafy, 2006) are obtained based on the following equations:
$$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!{{ K}_k} = {{ P}_k}{{ H}^{\rm T}}{\left({{ H}{{ P}_k}{{ H}^{\rm T}} + { R}} \right)^{  1}}, $$ (8) $$\!{{ Y}_{k, i}} = { H}\left({{{ X}_{k, i}}, {{ U}_{k, i}}} \right)\;\;\;\;\;{{ U}_{k, i}} \sim N\left({0, { R}} \right), $$ (9) $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!{{H}}\left({{{\bar{ X}}_k} } \right) = \frac{1}{N}\mathop \sum \nolimits_{i = 1}^N { H} \left({{{ X}_{k, i}}} \right), \;$$ (10) where R is the measurement noise variance, which is set to be the unit matrix multiplied by 0.005, and Q is set to equal to R.
To obtain the Kalman gain, the error covariance matrix P_{k} can be calculated as follow:
$${{ P}_k} = \frac{1}{{N  1}}{{ E}_k}{ E}_k^{\rm T}.$$ (11) The error expectation in each ensemble member is obtained as follows:
$${{ E}_k} = \left[ {{{{ X}}_{k, 1}}  {{{\bar{ X}}_k}}, \cdots, {{{ X}}_{k, N}}  {{{{\bar{ X}}}_k}} } \right], $$ (12) $$ {{{\bar { X}}_k}} = \frac{1}{N}\mathop \sum \nolimits_{i = 1}^N {{ X}_{k, i}}. \;\quad\quad\quad\quad\quad\quad$$ (13) Finally, the updated state variable value at time k can be obtained as follow:
$${ X}_k^{\rm up} = \frac{1}{N}\mathop \sum \nolimits_{i = 1}^N { X}_{k, i}^{\rm up}.$$ (14) When the measurements are a nonlinear combination of state variables, PH^{T} and HPH^{T} can be directly calculated by using the following equations (Houtekamer and Mitchell, 2001; Huang et al., 2008):
$$\!\;\;{{P}}{{{ H}}^{\rm T}} = \frac{1}{{N  1}}\mathop \sum \nolimits_{i = 1}^N \left[ {{{ X}_{k, i}}  { {\bar X}}} \right]{\left[ {{ H} \left({{{ X}_{k, i}}} \right)  { H}\left({ {\bar X}}\right)} \right]^{\rm T}}, $$ (15) $${{HP}}{{{H}}^{\rm T}} = \frac{1}{{N  1}}\mathop \sum \nolimits_{i = 1}^N \left[ {{ H}\left({{{ X}_{k, i}}} \right)  { H} \left({ {\bar X} }\right)} \right]{\left[ {{ H}\left({{{ X}_{k, i}}} \right)  { H}\left({ {\bar X}} \right)} \right]^{\rm T}}.$$ (16) 
To analyze how the initial state value and precipitation affect the soil moisture prediction, the RootMeanSquare Error (RMSE), Absolute Error (AE), and daily average bias (Bias) are used and calculated as follows:
$${\rm{RMSE}} = \sqrt {\frac{1}{n}\mathop \sum \nolimits_{k = 1}^N {{\left({{{ X}_{{\rm{pred}}, k}}  {{ X}_{{\rm{obs}}, k}}} \right)}^2}}, $$ (17) $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!{\rm{AE}} = \left {{{ X}_{{\rm{pred}}, k}}  {{ X}_{{\rm{obs}}, k}}} \right, \quad\quad\quad\quad$$ (18) $$\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!{\rm{Bias}} = {{{\bar { X}}_{{\rm{ass}}, d}}}  {{{\bar { X}}_{{\rm{sim}}, d}}}, \quad\quad\quad\quad$$ (19) where X_{pred,k} is the simulated/assimilated/predicted soil moisture at time k; X_{obs,k} is the observed soil moisture at time k; n is the total number of time steps;
${{{\bar { X}}_{{\rm{ass}}, d}}} $ and$ {{{\bar { X}}_{{\rm{sim}}, d}}} $ are the daily averaged soil moisture following the assimilated and simulated values at the end of assimilation period, respectively.