Unscented weighted ensemble Kalman filter for soil moisture assimilation

Journal Pre-proofs Research papers Unscented weighted ensemble Kalman filter for soil moisture assimilation Xiaolei Fu, Zhongbo Yu, Yongjian Ding, Yu Qin, Lifeng Luo, Chuancheng Zhao, Haishen Lü, Xiaolei Jiang, Qin Ju, Chuanguo Yang PII: DOI: Reference:

S0022-1694(19)31087-X https://doi.org/10.1016/j.jhydrol.2019.124352 HYDROL 124352

To appear in:

Journal of Hydrology

Received Date: Accepted Date:

21 September 2019 11 November 2019

Please cite this article as: Fu, X., Yu, Z., Ding, Y., Qin, Y., Luo, L., Zhao, C., Lü, H., Jiang, X., Ju, Q., Yang, C., Unscented weighted ensemble Kalman filter for soil moisture assimilation, Journal of Hydrology (2019), doi: https:// doi.org/10.1016/j.jhydrol.2019.124352

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Elsevier B.V. All rights reserved.

1

Unscented weighted ensemble Kalman filter for soil moisture assimilation

2

Xiaolei Fu1,2,*, Zhongbo Yu3,*, Yongjian Ding2, Yu Qin2, Lifeng Luo4, Chuancheng Zhao2,

3

Haishen Lü3, Xiaolei Jiang3, Qin Ju3, Chuanguo Yang3

4

1College

5

2State

6

Resources, Chinese Academy of Sciences, Lanzhou 730000, China

7

3State

8

Nanjing 210098, China

9

4Department

of Civil Engineering, Fuzhou University, Fuzhou 350116, China

Key Laboratory of Cryospheric Science, Northwest Institute of Eco-Environment and

Key Laboratory of Hydrology-Water Resources and Hydraulic Engineering, Hohai University,

of Geography, Environment, and Spatial Sciences, Michigan State University, East

10

Lansing, MI, 48824, USA

11

Corresponding authors: Xiaolei Fu, [email protected]; Zhongbo Yu, [email protected]

12

Abstract

13

A new data assimilation technique, unscented weighted ensemble Kalman filter (UWEnKF) was

14

developed based on the scaled unscented transformation and ensemble Kalman filter (EnKF). In

15

UWEnKF, the individual members selected are unequally weighted and symmetric about the

16

expectation. To investigate the performance of UWEnKF, nine assimilation experiments with different

17

ensemble sizes (161, 1601, 16001) and different assimilation frequencies (every 6 h, every 12 h, every

18

24 h) were designed to assimilate soil surface (5 cm) moisture data observed at station HY in the upper

19

reaches of the Yellow River, in the northeastern of Tibetan plateau, China into the Richards equation.

20

The results showed that the performance of the filter was greatly affected by random noise, and the

21

filter was sensitive to ensemble size and assimilation frequency. Increasing the ensemble size reduced

22

the effects of random noise on filter performance in several independent assimilation runs (i.e., it

23

decreased the differences between the results of the several independent assimilation runs). Reducing

24

the assimilation frequency also reduced the effects of random noise on filter performance. UWEnKF

25

gave more accurate soil moisture model results than EnKF for all ensemble sizes and assimilation

26

frequencies at all soil depths. Additionally, EnKF may have different performances according to

27

different initial conditions, but not for UWEnKF. Precipitation and soil properties uncertainties had

28

some impact on filter performance. Thus, UWEnKF is a better choice than EnKF, while it is more

29

computationally demanding, for improving soil moisture predictions by assimilating data from many

30

sources, such as satellite-observed soil moisture data, at a low assimilation frequency (e.g., every 24

31

h).

32

Keywords: Soil moisture; Richards equation; ensemble Kalman filter (EnKF); unscented weighted

33

ensemble Kalman filter (UWEnKF)

34

35

1. Introduction

36

Soil moisture is an important state variable in land surface processes, meteorology, hydrology and

37

agriculture (Heathman et al., 2003). Soil moisture determines the proportions of rainfall directed into

38

surface runoff and soil infiltration (Yu et al., 1999; Silberstein et al., 1999; Western and Bloschl, 1999;

39

Koster et al., 2003), influences the distribution of sensible and latent heat (Koster et al., 2004), and

40

affects irrigation scheduling and crop production in agriculture (Hanson et al., 1998; Ma et al., 1998;

41

Lü et al., 2011; Fu et al., 2014). However, despite its importance, accurate long-term observations of

42

soil moisture over large areas and in deep soil layers are difficult to obtain (Yu, et al., 2001; Huang et

43

al., 2008).

44

In recent decades, land surface models and hydrological models have been used to provide

45

long-term soil moisture values for increased depths of soil at a regional scale (Yu et al., 2006). For

46

example, the Biosphere-Atmosphere Transfer Scheme (Dickinson et al., 1993) and the Simple

47

Biosphere Model, SiB2 (Sellers et al., 1996) were used to obtain soil moisture values at different soil

48

depths. The widely-used Common Land Model (Oleson et al., 2004) has been used to simulate

49

horizontal and vertical soil moisture and temperature behavior for a long period. The Variable

50

Infiltration Capacity model, VIC (Liang et al., 1994) and the TOPMODEL-based Land

51

Surface-Atmosphere Transfer Scheme (Famiglietti and Wood, 1994a, b; Fu et al., 2018a) also simulate

52

soil moisture behavior for different spatial scales. However, the results given by these models are

53

biased because of uncertainties in model parameters, model structure, and forcing data (Yu et al., 2001,

54

2002; Huang et al., 2008; Xu et al. 2017). Remote sensing can be used to obtain soil moisture data on a

55

large scale; for example, Parinussa et al. (2012) used a radiative-transfer-based model to derive surface

56

soil moisture data from the WindSat spaceborne polarimetric microwave radiometer. However, remote

57

sensing observation and data collection are limited to the soil surface (Jackson et al., 1999; Wang et

58

al., 2016; Balsamo et al., 2007).

59

Data assimilation is an effective technique to improve model predictions and simulations and has

60

been widely used to estimate soil moisture by combining multi-source datasets (e.g., in situ soil

61

moisture observations and remote sensing data) over the past few decades (Kalman, 1960; Han and Li,

62

2008; Yu et al., 2014a, 2014b; Fu et al., 2014; Dong et al., 2015, 2016). Sequential data assimilation

63

has often been used to improve the accuracy of soil moisture predictions and simulations since the

64

introduction of the Kalman filter, KF (Kalman, 1960). However, the application of KF is limited to

65

linear dynamic systems and the linear observation operator (Luo and Moroz, 2009). To overcome the

66

limitation of KF, the extended Kalman filter (EKF) was developed by restricting the Taylor series for

67

nonlinear functions to the second order (Evensen, 1992; Miller et al., 1994). However, EKF is not

68

effective for large-scale problems because the full forward covariance matrix greatly increases its

69

computational cost (Luo and Moroz, 2009). The ensemble Kalman filter (EnKF), introduced by

70

Evensen (1994), is a computationally efficient technique that can accommodate nonlinear dynamics

71

and which has been used in soil moisture assimilation experiments (Huang et al., 2008; Fu et al., 2014,

72

2019; Liu et al., 2016; Brandhorst et al., 2017). Huang et al. (2008) used EnKF to assimilate the

73

Tropical Rainfall Measuring Mission Microwave Imager brightness temperature data into the SiB2

74

model to improve soil moisture predictions in Tibetan Plateau. Lievens et al. (2016) used EnKF to

75

assimilate the Soil Moisture and Ocean Salinity brightness temperature data products into VIC model

76

to improve soil moisture predictions. Brandhorst et al. (2017) used EnKF to manage the uncertainty of

77

soil hydraulic parameters in predicting soil moisture.

78

In EnKF, all ensemble members are equally weighted, which reduces the effects of important

79

members and exaggerates the effects of unimportant members. Luo and Moroz (2009) demonstrated

80

that even randomly selected members can have the correct mean and covariance in EnKF, but still

81

introduce spurious modes in the transformed distribution, which affects the performance of EnKF (Han

82

and Li, 2008; De Lannoy et al., 2006). Sigma-points with different weights for different members,

83

which supersedes the scaled unscented transformation (Julier and Uhlmann, 2004; Han and Li, 2008;

84

Luo and Moroz, 2009; Fu et al., 2018b), can more accurately reflect the importance of each member,

85

and may lead to incorporate the scaled unscented transformation to benefit the EnKF performance.

86

Han and Li (2008) showed that filter performance is influenced by different equally-sized ensembles

87

of randomly selected members. Thus, we aim to investigate the impact of the different randomly and

88

equally-sized ensembles on filters performance, and lessen its influence in successive independently

89

simulation runs; then, combined the characteristics of the scaled unscented transformation with EnKF

90

to create a new highly effective data assimilation technique, the unscented weighted ensemble Kalman

91

filter (UWEnKF).

92

To achieve the objective, experiments were designed to assimilate observed soil surface moisture

93

data into the Richards equation with ensemble sizes 161, 1601 and 16001, and assimilation intervals Δ

94

𝑡𝑎𝑠𝑠 = 6 h, 12 h, and 24 h, using data observed at station HY in the upper reaches of the Yellow River,

95

in the northeastern of the Tibetan Plateau, China. The field site and available data are described in

96

section 2. Section 3 describes the model operator, filter structure, experimental design and evaluation

97

criteria. The results are given and discussed in section 4, and the summary and conclusions are given

98

in section 5.

99 100

2. Field site and data

101

The source of the Yellow River is in the northeast of the Tibetan Plateau. The area has been the

102

subject of much research because of its unique land surface processes which are due to the

103

geographical location, terrain, and altitude (Fig. 1). There are areas of permanent and seasonal

104

permafrost in this region, so the soil moisture observations and atmospheric forcing data (e.g.,

105

precipitation, wind speed, air temperature) from 2015-05-17 to 2015-09-21 recorded at experimental

106

site HY in Hongyuan county for 5 cm, 20 cm, 40 cm and 80 cm depths were selected to avoid the

107

effects of permafrost. Soil moisture was measured every hour using an auto-measurement system

108

(Decagon Inc., USA), which consisted of an EM50 and four 5TM sensors; precipitation was recorded

109

every hour using a rain gauge. Land cover at the site is mainly alpine meadow, and the average values

110

of soil properties at different depths are listed in Table 1. Average annual precipitation is 749.1 mm,

111

and average annual temperature is 1.4 °C.

112

[Insert: Fig. 1]

113

[Insert: Table 1]

114 115

3. Methodology

116

3.1 Richards equation

117

The Richards equation is an one-dimensional equation that describes soil water fluxes in the

118

unsaturated zone of a homogeneous and isotropic soil (Richards, 1931; Hoeben and Troch, 2000;

119

Oleson et al., 2004; Lü et al., 2011; Chirico et al., 2014; Medina et al, 2014a, b). It was used as the

120

model operator in the data assimilation scheme. The -based Richards equation can be written as:

121 122

∂𝜃 ∂𝑡

∂

[

∂𝜃

]

= ∂𝑧 𝐷(𝜃)∂𝑧 ― 𝐾(𝜃) ―𝑒

(1)

where θ is volumetric soil moisture content (m3/m3), t is time (s), z is soil depth (mm), D(θ) = K(θ)

123

𝑑𝜓 𝑑𝜃

is the unsaturated diffusivity, thus, equation (1) is equal to

∂𝜃 ∂𝑡

∂

[

∂𝜓

]

= ∂𝑧 𝐾(𝜃) ∂𝑧 ― 𝐾(𝜃) ―𝑒; 𝜓 is soil

124

water matric potential (mm), K is soil hydraulic conductivity (mm/s), and e is evapotranspiration loss

125

(mm/s). The Clapp-Hornberger relationships (Clapp and Hornberger, 1978) are used to specify the

126

dependence between K, 𝜓 and soil moisture content θ: 𝜃 ―𝐵

127

ψ(θ) = 𝜓𝑠(𝜃𝑠)

𝜃 2𝐵 + 3

128

K(θ) = 𝐾𝑠(𝜃𝑠)

(2) (3)

129

Where 𝜓𝑠, Ks and 𝜃𝑠 are the saturated soil matric potential (mm), saturated soil hydraulic

130

conductivity (mm/s) and saturated volumetric soil moisture content (m3/m3), respectively; B is the

131

distribution exponent of soil porosity.

132

In the numerical solution for the Richards equation (Hoeben and Troch, 2000; Lü et al., 2010;

133

Medina et al, 2014a, b), the upper boundary condition was determined by infiltration and precipitation,

134

and the lower boundary condition was defined as free drainage. The initial state values used were the

135

observed data for each layer to achieve the objective, and they were set 0.2 m3/m3 to analyze the

136

impact of initial value on filter performance. Soil hydraulic conductivity was calculated from the soil

137

properties in Table 1, and evapotranspiration was calculated by using the FAO-56 Penman-Monteith

138

equation (Allen, 2002).

139

3.2 Ensemble Kalman filter (EnKF)

140

The state function F, which maps the previous state 𝑋𝑡 ― 1 at time 𝑡 ― 1 to state 𝑋𝑡 at time t.

141

The observation function H, which indicates the deterministic relationship between the system states

142

and observations. Both of them are required for EnKF-based data assimilation (Evensen, 1994;

143

Moradkhani et al., 2005; Kumar et al., 2008; Xie and Zhang, 2010; Yu et al., 2014a,b; Fu et al., 2014;

144

Han et al., 2014). The Richards equation is used as the state function. The observation function H is

145

used to link the system state vector with the observations. The two functions are presented as:

146

𝑋𝑡 = 𝐹(𝑋𝑡 ― 1, 𝑣𝑡 ― 1)

(4)

147

𝑌𝑡 = 𝐻(𝑋𝑡, 𝑢𝑡)

(5)

148

where 𝑋𝑡 and 𝑌𝑡 are the state and measurement vectors of soil moisture at time t, and 𝑣𝑡 and 𝑢𝑡 are

149

state noise and measurement noise (or error) at time t.

150 151

There are two steps in EnKF: forecasting and updating the state. Forecasting

152

Initial ensemble members 𝑋0,𝑖(𝑖 = 1,⋯,𝑀) is obtained by adding the random perturbations to the

153

initial value 𝑋0 according to the Gaussian distribution with zero mean and error covariance matrix P

154

(Huang et al., 2008; Fu et al., 2014). The state variables can then be predicted by the state function: 𝑋𝑡,𝑖 = 𝐹(𝑋𝑢𝑝 𝑡 ― 1,𝑖, 𝑣𝑡 ― 1,𝑖)

155 156

𝑣𝑡 ― 1,𝑖~𝑁(0,𝑄)

where 𝑋𝑢𝑝 𝑡 ― 1,𝑖 is the updated state value at time t ―1, 𝑋𝑡,𝑖 is the predicted state value, Q is the

157

model error covariance matrix, and M is the EnKF ensemble size.

158

Updating the state

159 160 161 162 163 164 165 166

(6)

The state predictions can be updated with the observations by: 𝑋𝑢𝑝 𝑡,𝑖 = 𝑋𝑡,𝑖 + 𝐾𝑡[𝑌𝑡,𝑖 ―𝐻(𝑋𝑡,𝑖)]

(7)

where 𝐾𝑡 is the Kalman gain which is calculated by: 𝐾𝑡 = 𝑃𝑡𝐻𝑇(𝐻𝑃𝑡𝐻𝑇 + 𝑅)

―1

(8)

where R is the observation error covariance matrix. The error covariance matrix 𝑃𝑡 is calculated by: 1

𝑃𝑡 = 𝑀 ― 1𝐸𝑡𝐸𝑇𝑡

(9)

The error matrix 𝐸𝑡 for each ensemble member and the mean value of the forecast state are given

167

by:

168

𝐸𝑡 = [𝑋𝑡,1 ― 𝑋𝑡 ,⋯,𝑋𝑡,𝑀 ― 𝑋𝑡]

169

𝑋𝑡 = 𝑀∑𝑖 = 1𝑋𝑡,𝑖

1

𝑀

(11)

1

𝑀

(12)

The updated value at time t is:

170

𝑢𝑝 𝑋𝑢𝑝 𝑡 = 𝑀∑𝑖 = 1𝑋𝑡,𝑖

171

If the measurements are a nonlinear combination of state variables, the terms 𝑃𝑡𝐻𝑇 and 𝐻𝑃𝑡𝐻𝑇

172 173

(10)

can be approximated by (Houtekamer and Mitchell, 2001): 1

𝑀

𝑇

174

𝑃𝑡𝐻𝑇 = 𝑀 ― 1∑𝑖 = 1[𝑋𝑡,𝑖 ― 𝑋𝑡][𝐻(𝑋𝑡,𝑖) ― 𝐻(𝑋𝑡)]

175

𝐻𝑃𝑡𝐻𝑇 = 𝑀 ― 1∑𝑖 = 1[𝐻(𝑋𝑡,𝑖) ― 𝐻(𝑋𝑡)][𝐻(𝑋𝑡,𝑖) ― 𝐻(𝑋𝑡)]

176

1

𝑀

𝑇

(14)

where: 1

𝑀

177

𝐻(𝑋𝑡) = 𝑀∑𝑖 = 1𝐻(𝑋𝑡,𝑖)

178

3.3 Unscented weighted ensemble Kalman filter (UWEnKF)

179

3.3.1 Scaled unscented transformation

(15)

Given a n-dimensional random state variable X with mean 𝑋 and covariance 𝑃𝑥, and the

180 181

(13)

transformation of X is estimated using:

182

(16)

Y = G(X)

183

To calculate the mean and covariance of Y using scaled unscented transform, the 2n + 1 points

184

𝑋𝑖 with corresponding weight (Julier and Uhlmann, 2004, Van der Merwe, 2004) are drawn according

185

to: 𝜆

𝜔(𝑚) 0 = 𝑛+𝜆

𝑋0 = 𝑋 186

𝑋𝑖 = 𝑋 + ( (𝑛 + 𝜆)𝑃𝑥)𝑖 𝑋𝑖 = 𝑋 ― ( (𝑛 + 𝜆)𝑃𝑥)𝑖

187

𝑖 = 1,⋯,𝑛 𝑖 = 𝑛 + 1,⋯,2𝑛

𝑖=0

𝜆

2 𝜔(𝑐) 0 = 𝑛 + 𝜆 + (1 ― 𝛼 + 𝛽) 𝑖 = 0

𝜔(𝑚) = 𝜔(𝑐) 𝑖 𝑖 =

1 2(𝑛 + 𝜆)

(17)

𝑖 = 1,⋯,2𝑛

where: λ = 𝛼2(𝑛 + 𝑙) ―𝑛 is an adjustable scaling parameter, and those points are called the

188

sigma-points. ( (𝑛 + 𝜆)𝑃𝑥)𝑖 is the ith row or column of the matrix square root of (𝑛 + 𝜆)𝑃𝑥. In the

189

sigma-points, α(0 ≤ α ≤ 1) controls the distribution spread of the sigma-points. 𝑙(𝑙 ≥ 0) is a

190

parameter introduced to guarantee the positive semi-definiteness of the covariance matrix. 𝛽(𝛽 ≥ 0)

191

is used to include the higher order moments of the distribution (Julier and Uhlmann, 2004; Van der

192

Merwe, 2004).

193

After that, the transformation of each sigma-point can be propagated as:

194

𝑌𝑖 = 𝐺(𝑋𝑖)

195

The mean 𝑌, covariance 𝑃𝑦 and cross-covariance 𝑃𝑥𝑦 can be calculated:

𝑖 = 0,⋯,2𝑛

2𝑛

196

𝑌 ≈ ∑𝑖 = 0𝜔(𝑚) 𝑖 𝑌𝑖

197

𝑇 𝑃𝑦 ≈ ∑𝑖 = 0𝜔(𝑐) 𝑖 (𝑌𝑖 ― 𝑌)(𝑌𝑖 ― 𝑌)

198

𝑇 𝑃𝑥𝑦 ≈ ∑𝑖 = 0𝜔(𝑐) 𝑖 (𝑋𝑖 ― 𝑋)(𝑌𝑖 ― 𝑌)

199

(18)

(19)

2𝑛

(20)

2𝑛

(21)

3.3.2 The framework of UWEnKF

200

The ensemble of sigma-points in UWEnKF can be created as follows.

201

Errors 𝑞𝑗𝑗(𝑗𝑗 = 1,⋯,𝑁) are randomly generated using a Gaussian distribution with zero mean and

202

203

error covariance matrix 𝑃𝑥. The sigma-points ensemble with number 2𝑛𝑁 + 1 is defined by: 𝑋0 = 𝑋 𝑋𝑖 = 𝑋 + 𝑞𝑗𝑗( (𝑛 + 𝜆)𝑎𝑏𝑠(𝑃𝑥))𝑗 𝑋𝑖 = 𝑋 ― 𝑞𝑗𝑗( (𝑛 + 𝜆)𝑎𝑏𝑠(𝑃𝑥))𝑗

𝑖 = 1,⋯,𝑛𝑁 𝑖 = 𝑛𝑁 + 1,⋯,2𝑛𝑁

𝑗𝑗 = 1,⋯,𝑁 𝑗 = 1,⋯,𝑛 𝑗𝑗 = 1,⋯,𝑁 𝑗 = 1,⋯,𝑛

(22)

204

where ( (𝑛 + 𝜆)𝑎𝑏𝑠(𝑃𝑥))𝑗 is the jth row or column of the matrix square root of (𝑛 + 𝜆)𝑎𝑏𝑠(𝑃𝑥), and

205

𝑎𝑏𝑠(𝑃𝑥) = (|𝑃𝑜,𝑘|)𝑛 × 𝑛, 𝑜 = 1,⋯,𝑛

206

study, there are the same number of ensemble members (i.e., the number of sigma-points) in UWEnKF

207

and EnKF, therefore, 𝑀 = 2𝑛𝑁 + 1.

208

Forecasting

and 𝑘 = 1,⋯,𝑛; and N is the number of errors generated. In this

209

The state prediction for each ensemble member is given by the state function: 𝑋𝑡|𝑡 ― 1, 𝑖 = 𝐹(𝑋𝑡 ― 1,𝑖, 𝑣𝑡 ― 1,𝑖)

210 211 212

215

The error covariance matrix of state X at time t is calculated by: 𝑃𝑋𝑡|𝑡 ― 1 = 𝐴(𝐼. × 𝑊𝑇)𝐴𝑇

𝐴 = [𝑋𝑡|𝑡 ― 1, 0 ― 𝑋𝑡|𝑡 ― 1,⋯,𝑋𝑡|𝑡 ― 1,

222 223 224

― 𝑋𝑡|𝑡 ― 1]

2𝑛𝑁

𝑊𝑡,𝑖 =

218

221

2𝑛𝑁

𝑋𝑡|𝑡 ― 1 = ∑𝑖 = 0𝑊𝑡,𝑖𝑋𝑡|𝑡 ― 1,

217

220

(24)

where I is the identity matrix, 𝑊 = [𝑊𝑡,0,⋯,𝑊𝑡,2𝑛𝑁], and 𝐼. × 𝑊𝑇 is a diagonal matrix. In Eqn (24)

216

219

(23)

Updating

213 214

𝑣𝑡 ― 1,𝑖~𝑁(0,𝑄)

𝜔𝑡,𝑖 =

(26)

𝑖

(2𝑛𝑁 + 1)𝜔𝑡,𝑖 ― 𝑓𝑖𝑥((2𝑛𝑁 + 1)𝜔𝑡,𝑖) 2𝑛𝑁

2𝑛𝑁 + 1 ― ∑𝑖 = 0𝑚𝑖 exp ( ― 0.5/𝑅)(𝑌𝑡 ― 𝐻(𝑋𝑡|𝑡 ― 1, 𝑖))

(25)

(27)

2

2𝑛𝑁

∑𝑖 = 0exp ( ― 0.5/𝑅)(𝑌𝑡 ― 𝐻(𝑋𝑡|𝑡 ― 1, 𝑖))2

(28)

where 𝑚𝑖 = 𝑓𝑖𝑥((2𝑛𝑁 + 1)𝜔𝑡,𝑖). Each predicted state is mapped to an observation using the observation function: 𝑌𝑡, 𝑖 = 𝐻(𝑋𝑡|𝑡 ― 1, 𝑖, 𝑢𝑡,𝑖)

𝑢𝑡,𝑖~𝑁(0,𝑅)

2𝑛𝑁

𝑌𝑡 = ∑𝑖 = 0𝑊𝑡,𝑖𝑌𝑡,

(29) (30)

𝑖

The error covariance matrix of generated measurements is calculated by: 𝑃𝑌𝑡 = 𝑆(𝐼. × 𝑊𝑇)𝑆𝑇

225

𝑆 = [𝑌𝑡, 0 ― 𝑌𝑡,⋯,𝑌𝑡,

226

The error covariance matrix of state variables and the Kalman gain is calculated by:

2𝑛𝑁

― 𝑌𝑡]

(31) (32)

227

𝑃𝑋𝑌𝑡 = 𝐴(𝐼. × 𝑊𝑇)𝑆𝑇

(33)

228

𝐾𝑡 = 𝑃𝑋𝑌𝑡(𝑃𝑌𝑡 + 𝑅) ―1

(34)

229 230

Each ensemble member is updated by: 𝑋𝑡,𝑖 = 𝑋𝑡|𝑡 ― 1, 𝑖 + 𝐾𝑡(𝑌𝑡,𝑖 ― 𝑌𝑡, 𝑖)

(35)

231

where 𝑌𝑡,𝑖 is the generated observation ensemble. The updated value at time t is:

232

𝑋𝑢𝑝 𝑡 = ∑𝑖 = 0𝑊𝑡,𝑖𝑋𝑡,𝑖

2𝑛𝑁

(36)

233

𝑢𝑝 𝐸 = [𝑋𝑡,0 ― 𝑋𝑢𝑝 𝑡 ,⋯,𝑋𝑡,2𝑛𝑁 ― 𝑋𝑡 ]

(37)

234

𝑃𝑋𝑡 = 𝐸(𝐼. × 𝑊𝑇)𝐸𝑇

(38)

235 236 237 238 239 240 241 242 243 244 245 246 247 248

The unscented weighted ensemble Kalman filter procedure is summarized as follows: a. Prediction step (1) Initialization: Draw the ensemble members 𝑋𝑖 (𝑖 = 0,⋯,2𝑛𝑁) at time t-1 based on Eqns (22, 24-26) or Eqns (22, 37-38) (2) Prediction: Predict the state variables 𝑋𝑡,𝑖 (𝑖 = 0,⋯,2𝑛𝑁) at time t based on the above ensemble members using Eqn (23) b. Update step (3) Weights: Calculate the weights 𝑊𝑡,𝑖 (𝑖 = 0,⋯,2𝑛𝑁) for all members 𝑋𝑡,𝑖 based on Eqns (27-30) (4) Resampling: If the weight 𝑊𝑡,𝑖 is too small, resample the ensemble members, and re-calculate the corresponding weights (5) Updating: Update all ensemble members based on Eqns (24-35) and calculate the state estimation using Eqn (36) 3.3.3 The differences between UWEnKF and some new KF-based filters

249

Some new KF-based filters have recently been published, such as the weighted ensemble Kalman

250

filter ---- WEnKF (Papadakis et al., 2010), the weighted ensemble transform Kalman filter ----

251

WETKF (Beyou et al., 2013) and the ensemble unscent Kalman filter ---- EnUKF (Luo and Moroz,

252

2009). However, they differ from UWEnKF proposed in this study.

253

In more detail, the randomly selected ensemble members in WEnKF and WETKF are not

254

symmetric about the mean, which may introduce spurious modes in transformed distribution, similar to

255

EnKF (Luo and Moroz, 2009). Errors are selected randomly in UWEnKF, but the ensemble members

256

(i.e., sigma-points) obtained by Eqn (22) are symmetric about the mean, which reduces the effect of

257

sample error (Luo and Moroz, 2009). For EnUKF, the sigma-points are derived from Eqn (17) with

258

dimension 2𝐿 + 1 and 𝐿 = 𝑛 + 𝐿𝑠 + 𝐿𝑜, where n is the dimension of state variable, 𝐿𝑠 and 𝐿𝑜 are

259

the dimensions of simulation noise and observation noise. Thus, the dimensions of state variable,

260

simulation noise and observation noise determine the dimension of the sigma-points. However, the

261

number of sigma-points in UWEnKF 2𝑛𝑁 + 1 is adjustable, whereby, the number of ensemble

262

members differs between UWEnKF and EnUKF. In addition, the weights are 𝜔0 = 𝐿 + 𝜆 and 𝜔𝑖 =

263

1 2(𝐿 + 𝜆),

264

(27-28). The ensemble members in UWEnKF are therefore random and symmetric while their weights

265

are different. Thus, UWEnKF has the advantages of EnKF, and gives larger weight to the important

266

ensemble member, and reduces the impact of sampling error (Luo and Moroz, 2009).

𝜆

𝑖 = 1,⋯,2𝐿 for each sigma-point in EnUKF, whereas in UWEnKF, they are given by Eqns

267 268

3.4 Data assimilation experiments design

269

To achieve the aim of this study, we set the ensemble size to 161, 1601 and 16001 members and

270

the assimilation interval Δ𝑡𝑎𝑠𝑠 to 6 h, 12 h and 24 h. Thus, nine assimilation experiments (Table 2)

271

were conducted, and the model with each filter was run independently 100 times in each experiment.

272

The time step in the simulation process was set to 1 h, and only the surface (5 cm) observed soil

273

moisture data were assimilated into the Richards equation for each assimilation interval. The

274

assimilation frequencies every 6 h, every 12 h, and daily (every 24 h) corresponded to the assimilation

275

intervals for the data to be assimilated (Δ𝑡𝑎𝑠𝑠= 6 h, 12 h and 24 h). Fig. 2 shows the flowchart of the

276

soil moisture assimilation process using EnKF and UWEnKF.

277

[Insert: Table 2]

278

[Insert: Fig. 2]

279

The nine assimilation experiments were followed by another five brief experiments that were

280

designed and conducted to analyze the impact of initial value, uncertainty of precipitation and soil

281

properties on filter performance (Table 3). The precipitation in assimilation period and soil properties

282

were set to have error ranges of ±5% of observed values. The initial value of soil moisture was set

283

0.2 m3/m3 in all soil layers.

284

[Insert: Table 3]

285 286

3.5 Evaluation criteria

287

The root mean square error (RMSE) and mean absolute error (MAE) were used to evaluate the

288

influence of randomly generated ensemble members on filter performance when assimilating surface

289

soil moisture data: 1

1

293

𝑇𝑇

𝑀𝐴𝐸 = 𝑇𝑇∑𝑡 = 1|𝑋𝑝𝑟𝑒𝑑,𝑡 ― 𝑋𝑜𝑏𝑠,𝑡|

291 292

𝑇𝑇

1/2

𝑅𝑀𝑆𝐸 = (𝑇𝑇∑𝑡 = 1(𝑋𝑝𝑟𝑒𝑑,𝑡 ― 𝑋𝑜𝑏𝑠,𝑡)2)

290

(39) (40)

where 𝑋𝑝𝑟𝑒𝑑,𝑡 is the predicted/assimilated soil moisture at time t, 𝑋𝑜𝑏𝑠,𝑡 is the observed soil moisture at time t, and TT is the total number of time steps.

294 295 296

4. Results and discussion In the study, the model was run for the period 2015-05-17 to 2015-09-21, and the assimilation

297

experiments were conducted from 2015-08-01 (day 213) to 2015-09-21 (day 264) using data from

298

station HY.

299 300

4.1 Performance of EnKF and UWEnKF

301

Figs. 3 and 4 (Cases 1 - 9) show the 100 RMSE and 100 MAE values distributions with 100

302

independent assimilations runs for each of the nine experimental trials (3 assimilation frequencies (6 h,

303

12 h, 24 h) × 3 ensemble sizes (161, 1201, 16001)) using EnKF and UWEnKF at four different soil

304

depths (5 cm, 20 cm, 40 cm, 80 cm). The RMSE and MAE values indicate that random noise affected

305

the performance of both EnKF and UWEnKF.

306

It is importantly to note that the 100 RMSE and MAE values for UWEnKF were less than those

307

for EnKF no matter what the ensemble size or assimilation interval was. The 100 RMSE and MAE

308

values for UWEnKF were more distributed than those for EnKF, especially for the assimilation

309

frequency of every 6 h in soil depths 5 cm (Fig. 3(a1) and Fig. 4(a1)), 40 cm (Fig. 3(a3) and Fig. 4(a3))

310

and 80 cm (Fig. 3(a4) and Fig. 4(a4)). As a result of the distribution, the means of 100 RMSE and

311

MAE values for UWEnKF were smaller than those for EnKF. In other words, the performance of

312

UWEnKF in soil moisture assimilation experiment was better than that of EnKF irrespective of

313

ensemble size and assimilation frequency. The reason for this may be that the more important

314

ensemble members were given larger weights in UWEnKF, whereas each member had the same

315

weight in EnKF, which effectively decreased the weight of important members and increased the

316

weight of unimportant members.

317 318

[Insert: Figs. 3 - 4]

319

4.2 Effects of ensemble size

320

Figs. 3 and 4 also show that the 100 RMSE and MAE values were distributed differently when

321

the ensemble sizes differed. To analyze the effects of ensemble size on EnKF and UWEnKF, the

322

means and standard errors (Std) of 100 independent RMSE values for different ensemble sizes and

323

different assimilation intervals for soil moisture assimilation using EnKF and UWEnKF were shown in

324

Fig. 5 (Cases 1 - 9). Fig. 5 (a1 – a4) and Fig. 5 (c1 – c4) show the means of 100 RMSE values for

325

EnKF and UWEnKF at all soil depths. The mean RMSE for each assimilation frequency showed little

326

change when the ensemble size changed from 161 to 16001. That is to say, the RMSE for the

327

assimilation results was distributed around the same value for a given assimilation frequency whatever

328

the ensemble size changes. Fig. 5 (a4) shows that the mean RMSE for EnKF was larger than the

329

RMSE for simulations (i.e., model results without assimilation) at the fourth soil layer (80 cm)

330

whatever the assimilation frequency and ensemble size was. The reason for this will be discussed later.

331

At other soil depths, the smaller RMSE values show that EnKF gave better soil moisture simulations.

332

UWEnKF improves the accuracy of soil moisture simulations at all soil depths.

333

[Insert: Fig. 5]

334

Figs. 5 (b1 – b4) and 5 (d1 - d4) show the Stds of 100 RMSE values for EnKF and UWEnKF at

335

different soil depths. It is noted that the Std for both EnKF and UWEnKF decreased when the

336

ensemble size changed from 161 to 16001 and were least for ensemble size 16001. That is to say that

337

the range of the 100 RMSE values decreases as ensemble size increases, and the 100 RMSE values

338

were distributed more closely about the mean value for the larger ensemble size. This occurs for the

339

following reasons. Although each member selected from the probability distribution in an ensemble

340

was different, the effect of random noise on the mean value of the members for each ensemble was

341

small and decreased as the ensemble size increased for all 100 independent sampling ensembles. This

342

indicates that increasing the number of ensemble members (ensemble size) selected reduced the effect

343

of random noise on filter performance, which is consistent with the results in Dong et al. (2015). Thus,

344

the filter performance over several independent assimilation runs was more stable with a larger

345

ensemble size.

346

The mean and Std of 100 independent MAE values of soil moisture for EnKF and UWEnKF with

347

different ensemble sizes and different assimilation intervals were shown in Fig. 6 (Cases 1 - 9). The

348

observations made from this figure parallel those made for RMSE from Fig. 5 for different ensemble

349

sizes.

350

[Insert: Fig. 6]

351 352

4.3 Effects of assimilation frequency

353

Figs. 5 - 6 also show that the assimilation frequency influenced filter performance. It is noted that

354

the Std of RMSE (Fig. 5 (b1 – b4)) and MAE (Fig. 6 (b1 – b4)) for EnKF decreased when the

355

assimilation frequency changed from every 6 h to every 24 h except at soil depth 5 cm for ensemble

356

size 161. When the ensemble size increased from 161 to 16001, the differences in Std between the

357

three assimilation frequencies for EnKF became small. Similar conclusions were drawn for UWEnKF

358

(Fig. 5 (d1 – d4) and Fig. 6 (d1 – d4)), except at soil depth 20 cm for ensemble size 161 with

359

assimilation intervals Δ𝑡𝑎𝑠𝑠= 12 h. However, the differences in Std values of RMSE or MAE between

360

the three assimilation frequencies for UWEnKF were larger than those for EnKF at soil depths 5 cm,

361

40 cm, and 80 cm, when the assimilation frequency changed from every 6 h to every 12 h. In other

362

words, UWEnKF was more easily influenced by random noise than EnKF, especially at a higher

363

assimilation frequency. This was due to the unequal weights of the ensemble members in UWEnKF.

364

The important members differed among several sampled ensembles, which led to the large differences

365

in Std of RMSE values between several independent assimilation runs, in comparison to EnKF, in

366

which the members are equally weighted. For a low assimilation frequency (every 24 h) with large

367

ensemble size (16001), the random noise has little effect on the performance of either filter.

368 369

4.4 Computational overhead

370

Computational complexity is also a factor to consider in the evaluation of filter performance. It

371

determined that the calculation complexities of EnKF and UWEnKF are O(M) and O(M2), where M is

372

ensemble size. This implies that UWEnKF has a larger computational cost than EnKF. However, the

373

results obtained above showed that UWEnKF was more effective than EnKF in assimilating soil

374

moisture. Thus, using UWEnKF with low assimilation frequency will reduce the computational burden

375

to some degree but still provide accurate state predictions.

376 377

4.5 Soil moisture assimilation results

378

The Std of RMSE and MAE for assimilation results of several independent assimilation runs

379

decreased as the ensemble size increased. We take ensemble size 16001 as an example to show the

380

change in soil moisture at different assimilation frequencies. Figs. 7 - 9 show the soil moisture

381

assimilations for EnKF and UWEnKF with ensemble size 16001 and Δ𝑡𝑎𝑠𝑠= 6 h (Fig. 7, Case 3), 12 h

382

(Fig. 8, Case 6) and 24 h (Fig. 9, Case 9). The figures show that the simulations (model results without

383

assimilation) well represented the changes in soil moisture when compared to the observation. Soil

384

moisture decreased after precipitation, perhaps because of the sandy loam soil and high evaporation.

385

[Insert: Figs. 7 - 9]

386

The soil moisture model simulations show the greatest improvements from EnKF and UWEnKF

387

when the surface (5 cm depth) in-situ soil moisture observations were assimilated, except at the lowest

388

soil layer (80 cm) for EnKF. For the first three soil layers (5 cm, 20 cm, and 40 cm), the assimilations

389

using UWEnKF were much closer to the observations than those using EnKF, i.e., UWEnKF

390

performed better than EnKF. At soil depth 80 cm, EnKF did not improve the soil moisture simulations,

391

but UWEnKF did, as shown in Figs. 7 - 9. The first reason for this is that the simulations

392

overestimated soil moisture in the upper three soil layers, but underestimated it in the 80 cm layer. The

393

second reason has to do with the ensemble member weighting. The ensemble members were equally

394

weighted in EnKF, and soil moisture was underestimated by EnKF when compared to the simulations

395

for the upper three soil layers. The combination of equal member weighting and underestimation of

396

soil moisture led to the EnKF assimilations underestimating soil moisture compared to the simulation

397

at 80 cm depth because only the surface soil moisture observed data were assimilated (i.e., EnKF did

398

not improve soil moisture model simulations for the lowest soil layer). However, for UWEnKF, the

399

unequal weight for each ensemble member and the members being symmetric distribution about the

400

expectation perhaps led to UWEnKF improving soil moisture simulations at all soil layers.

401

Comparisons between the results shown in Figs. 7 - 9 show that the assimilation results differed

402

most from the observations when the assimilation frequency decreased from every 6 h to every 24 h.

403

This can also be seen in Fig. 10 (Cases 3, 6, 9), which shows the RMSE and MAE values for EnKF

404

and UWEnKF with ensemble size 16001 and different Δ𝑡𝑎𝑠𝑠 values correspond to Figs. 7 - 9. Fig. 10

405

shows that the RMSE and MAE values for EnKF and UWEnKF with different Δ𝑡𝑎𝑠𝑠 frequencies were

406

all less than model simulations (the model results without assimilation) except at 80 cm depth for

407

EnKF, and the RMSE and MAE for UWEnKF were smaller than those for EnKF. This demonstrates

408

that the performance of UWEnKF was better than that of EnKF in soil moisture assimilation. RMSE

409

and MAE values for EnKF and UWEnKF increased as the assimilation frequency changed from every

410

6 h to every 24 h, which also can be found in Figs. 5 - 6, except for EnKF at 80 cm. This is because the

411

number of soil moisture model results updated decreased as the assimilation frequency decreased.

412

[Insert: Fig. 10]

413 414

4.6 Effects of uncertainty of precipitation and soil properties and initial value

415

In cases 10 - 14, ensemble size 16001 is used as an example to analyze the effects of uncertainty

416

in precipitation and soil properties and initial values on filter performance. In the impact experiments,

417

an error range of ±5% was allowed for precipitation and soil properties, and the initial soil moisture

418

0.2 m3/m3 on 2015-05-17 was set in all four soil layers (Table 3). Fig. 11 (Cases 10-11) shows the

419

RMSE and MAE for ±5% errors of observed precipitation during the assimilation period. It is noted

420

that the RMSE and MAE varied for both filters when precipitation error changed from -5% to +5%.

421

That is, they are increased for the topmost three layers and decreased for the lowest layer. For high

422

assimilation frequency (e.g., every 6 h), the variation of RMSE and MAE for two filters was smaller

423

than that for simulations. For low assimilation frequency (e.g., every 24 h), sometimes the variation of

424

RMSE and MAE for two filters was larger than that for simulation, but not apparently. This is because

425

many soil moisture model results were updated with high assimilation frequency. Similar to that of

426

precipitation, the RMSE and MAE for ±5% errors of soil properties were shown in Fig. 12 (Cases 12

427

- 13). The RMSE and MAE decreased for the three topmost soil layers as error in soil properties

428

changed from -5% to +5%, but increased from -5% error to 0 and decreased from 0 to +5% for the

429

lowest layer, while the changes are small.

430

[Insert: Figs. 11 - 12]

431

Ensemble size of 16001 is used as an example to investigate the effects of initial values on filters

432

performance. Fig. 13 shows the soil moisture assimilation results for initial soil moisture of 0.2 m3/m3

433

at all soil layers on 2015-05-17 with assimilation interval 6 h. Comparison of the results shown in Fig.

434

13 to those in Fig. 7, shows that both filters can improve the soil moisture simulations for the first

435

three soil layers, and that UWEnKF performed better than EnKF. However, for the bottom layer, both

436

filters can also improve the soil moisture simulations with initial soil moisture 0.2 m3/m3, but EnKF

437

does not when measured value is used as the initial value, and UWEnKF still performs better than

438

EnKF. This is because when the measured value is used as the initial soil moisture, the simulations

439

underestimate the soil moisture at bottom layer and overestimated it at upper three layers. However,

440

when the initial soil moisture is set to 0.2 m3/m3, the simulations are overestimated at all soil layers,

441

which led to the well performance of EnKF in bottom layer. The results are also can be seen by

442

comparing the RMSE and MAE results in Fig. 10 and Fig. 14 with different initial values.

443

[Insert: Figs. 13 - 14]

444 445

5. Summary and conclusions

446

In this study, based on EnKF and the scaled unscented transformation, we presented a new data

447

assimilation technique, the unscented weighted ensemble Kalman filter (UWEnKF). In UWEnKF, the

448

randomly selected ensemble members are symmetric about the expectation and are unequally weighted

449

(i.e., no two members have equal weights). The performance of UWEnKF was investigated by

450

assimilating surface (5 cm) soil moisture data observed at the HY station in the upper reaches of the

451

Yellow River, China into the Richards equation with different assimilation frequencies and different

452

ensemble sizes. The preceding discussion of the results led to the following conclusions:

453

(1) The RMSE and MAE values for 100 independent assimilation runs show that the performance

454

of the filter was greatly affected by noise perturbations, which was also found by De Lannoy et al.

455

(2006). UWEnKF improves the model simulations better than EnKF at all soil depths no matter what

456

the ensemble size and the assimilation frequency are. At the deepest soil layer, EnKF is influenced by

457

the initial value, but UWEnKF greatly improves soil moisture simulations. Uncertainty in precipitation

458

and soil properties has some impact on filter performance because both impact the model performance.

459

(2) The filter was sensitive to ensemble size and assimilation frequency. Increasing the ensemble

460

size (i.e., the number of randomly selected members) can reduce the effects of random noise on filter

461

performance over several independent assimilation runs. As the ensemble size increases, the

462

differences between the results of independent assimilation runs decrease. Decreasing the assimilation

463

frequency (i.e., increasing the assimilation interval) also reduces the effects of random noise on filter

464

performance.

465

(3) UWEnKF is computationally more cost than EnKF, but more effective. Thus, for future

466

multi-source dataset assimilations (e.g., using satellite soil moisture data) with low assimilation

467

frequency, such as every 24 h, UWEnKF is a better choice if the computational burden is not a

468

concern.

469

Overall, UWEnKF is an effective and practical data assimilation technique that improves soil

470

moisture model simulations, although it requires considerably more computational resources. It is

471

benefit to obtain the high accuracy of catchment or global soil moisture estimations for different soil

472

depths using the sparse in-situ soil observations and remote sensing data, which is helpful for rainfall

473

proportions in watershed flood forecasting, drought warning and management of agriculture

474

production.

475 476

Acknowledgment

477

This study was supported by the National Key R&D Program of China (Grant No.

478

2016YFC0402710); the National Natural Science Foundation of China (Grant No. 51709046,

479

51539003, 41761134090, 41601562); the National Science Funds for Creative Research Groups of

480

China (No. 51421006); the program of Dual Innovative Talents Plan and Innovative Research Team in

481

Jiangsu Province ； the Open Foundation of the State Key Laboratory of Cryospheric Science,

482

Northwest Institute of Eco-Environment and Resources, Chinese Academy of Sciences

483

(SKLCS-OP-2018-03).

484 485

References

486

Allen, R. G. (2002). Using the FAO-56 dual crop coefficient method over an irrigated region as part of

487

an evapotranspiration inter-comparison study. Journal of Hydrology, 229: 27-41.

488

Balsamo, G., Mahfouf, J. F., Bélair, S., and Deblonde, G. (2007). A land data assimilation system for

489

soil moisture and temperature: An information content study. Journal of Hydrometeorology, 8:

490

1225–1242.

491 492

Beyou, S., Cuzol, A., Gorthi, S., and Mémin, E. (2013). Weighted ensemble transform Kalman filter for image assimilation. Tellus A: Dynamic Meteorology and Oceanography, 65 (1): 18803.

493

Brandhorst, N., Erdal., D., and Neuweiler, I. (2017). Soil moisture prediction with the ensemble

494

Kalman fiter: Handling uncertainty of soil hydraulic parameters. Advances in Water Resources,

495

110: 360-370.

496

Chirico, G. B., Medina, H., and Romano, N. (2014). Kalman filters for assimilating near-surface

497

observations into the Richards equation - Part 1: Retrieving state profiles with linear and

498

nonlinear numerical schemes. Hydrology and Earth System Sciences, 18: 2503-2520.

499 500

Clapp, R. B., and Hornberger, G. M. (1978). Empirical equations for some hydraulic properties. Water Resources Research, 14: 601-604.

501

De Lannoy, G. J. M., Houser, P. R., Pauwels, V. R. N., and Verhoest, N. E. C. (2006). Assessment of

502

model uncertainty for soil moisture through ensemble verification. Journal of Geophysical

503

Research, 111, D10101, doi:10.1029/2005JD006367.

504

Dickinson, R. E., Henderson-Sellers, A., and Kennedy, P. J. (1993). Biosphere-Atmosphere Transfer

505

Scheme (BATS) Version 1e as coupled to the NCAR Community Climate Model. Tech. Note,

506

National Center for Atmospheric Research, TN-387+STR, DOI: 10.5065/D67W6959.

507

Dong, J., Steele-Dunne, S. C., Judge, J., and van de Giesen, N. (2015). A particle batch smoother for

508

soil moisture estimation using soil temperature observations. Advances in Water Research, 83:

509

111-122.

510

Dong, J., Steele-Dunne, S. C., Ochsner, T. E., and van de Giesen, N. (2016). Determining soil moisture

511

and soil properties in vegetated areas by assimilating soil temperatures. Water Resources

512

Research, 52: 4280-4300.

513 514

Evensen, G. (1992). Using the extended Kalman filter with a multilayer quasi-geostrophic ocean model. Journal of Geophysical Research, 97: 17905-17924.

515

Evensen, G. (1994). Sequential data assimilation with a nonlinear quasi-geostrophic model using

516

Monte-Carlo methods to forecast error statistics. Journal of Geophysical Research, 99:

517 518 519 520 521

10143-10162. Famiglietti, J. S., and Wood, E. F. (1994a). Multiscale modeling of spatially variable water and energy balance process. Water Resources Research, 30(11): 3061-3078. Famiglietti, J. S., and Wood, E. F. (1994b). Application of multiscale water and energy balance models on a tallgrass prairie. Water Resources Research, 30(11): 3079-3093.

522

Fu, X., Luo, L., Pan, M., Yu, Z., Tang, Y., and Ding, Y. (2018a). Evaluation of Topmodel-based Land

523

Surface-Atmosphere Transfer Scheme (TOPLATS) through a soil moisture simulation. Earth

524

Interactions, 22(15), doi: 10.1175/EI-D-17-0037.1.

525

Fu, X., Yu, Z., Ding, Y., Tang, Y., Lü, H., Jiang, X., and Ju, Q. (2018b). Analysis of influence of

526

observation operator on sequential data assimilation through soil temperature simulation with

527

common land model. Water Science and Engineering, 11(3): 196-204.

528

Fu, X., Yu, Z., Luo, L., Lü, H., Liu, D., Ju, Q., Yang, T., Xu, F., Gu, H., Yang, C., Chen, J., and Wang,

529

T. (2014). Investigating soil moisture sensitivity to precipitation and evapotranspiration errors

530

using SiB2 model and ensemble Kalman filter. Stochastic Environmental Research and Risk

531

Assessment, 28(3): 681-693.

532

Fu, X., Yu, Z., Tang, Y., Ding, Y., Lü, H., Zhang, B., Jiang, X., and Ju, Q. (2019). Evaluating soil

533

moisture predictions based on ensemble Kalman filter and SiB2 model. Journal of Meteorological

534

Research, 33(2): 190-205.

535

Han, X., Franssen, H.J.H., Montzka, C., and Vereecken, H. (2014). Soil moisture and soil properties

536

estimation in the Community Land Model with synthetic brightness temperature observations.

537

Water Resources Research, 50: 6081-6105.

538

Han, X., and Li, X. (2008). An evaluation of the nonlinear/non-Gaussian filters for the sequential data

539

assimilation. Remote Sensing of Environment, 112: 1434-1449.

540

Hanson, J. D., Ahuja, L. R., Shaffer, M. D., Rojas, K. W., DeCoursey, D. G., Farahani, H., and

541

Johnson, K. (1998). RZWQM: simulating the effects of management on water quality and crop

542

production. Agricultural Systems, 57: 161–195.

543 544 545 546 547 548

Heathman, G. C., Starks, P. J., Ahuja, L. R., and Jakson, T. J. (2003). Assimilation of surface soil moisture to estimate profile soil water content. Journal of Hydrology, 279: 1-17. Hoeben, R., and Troch, P. A. (2000). Assimilation of active microwave observation data for soil moisture profile estimation. Water Resources Research, 36: 2805-2819. Houtekamer, P. L., and Mitchell, H. L. (2001). A sequential ensemble Kalman filter for atmospheric data assimilation. Monthly Weather Review, 129(1): 123−137.

549

Huang, C., Li, X., Lu, L., and Gu, J. (2008). Experiments of one-dimensional soil moisture

550

assimilation system based on ensemble Kalman filter. Remote Sensing of Environment, 112:

551

888–900.

552

Jackson, T., Levine, D., Hsu, A., Oldak, A., Starks, P., Swift, C., Isham, J., and Haken, M. (1999). Soil

553

moisture mapping at regional scales using microwave radiometry: The southern great plains

554

hydrology experiment. IEEE Transactions on Geoscience and Remote Sensing, 37: 2136–2151.

555 556 557 558

Julier, S.J., and Uhlmann, J.K. (2004). Unscented filtering and nonlinear estimation. Proceedings of the IEEE Aerospace and Electronic Systems, 92(12): 410-422. Kalman, R. E. (1960). A new approach to linear filtering and prediction problems. Transactions of the ASME-Journal of Basic Engineering, 82(Series D): 35−45.

559

Koster, R. D., Dirmeyer, P. A., Guo, Z., Bonan, G., Chan, E., Cox, P., Gordon, C. T., Kanae, S.,

560

Kowalczyk, E., Lawrence, D., Liu, P., Lu, C-H., Malyshev, S., Mcavaney, B., Mitchell, K.,

561

Mocko, D., Oki, T., Oleson, K., Pitman, A., Sud, Y. C., Taylor, C. M., Verseghy, D., Vasic, R.,

562

Xue, Y., and Yamada, T. (2004). Regions of strong coupling between soil moisture and

563

precipitation. Science, 2004: 1138–1141.

564 565

Koster, R. D., Suarez, M. J., Higgins, R. W., and van den Dool, H. M. (2003). Observational evidence that soil moisture variations affect precipitation. Geophysical Research Letters, 30: 1–4.

566

Kumar, S. V., Reichle, R. H., Peters-Lidard, C. D., Koster, R. D., Zhan, X., Crow, W. T., Eylander, J.

567

B., and Houser, P. R. (2008). A land surface data assimilation framework using the land

568

information system: Description and applications. Advances in Water Resources, 31(11): 1419–

569

1432.

570

Liang, X., D. P. Lettenmaier, E. F. Wood, and S. J. Burges (1994), A Simple Hydrologically Based

571

Model of Land-Surface Water and Energy Fluxes for General-Circulation Models, Journal of

572

Geophysical Research-Atmospheres, 99(D7): 14415-14428.

573

Lievens, H., De Lannoy, G. J. M., Al Bitar, A., Drusch, M., Dumedah, G., Hendricks Franssen, H. J.,

574

Kerr, Y. H., Tomer, S. K., Martens, B., Merlin, O., Pan, M., Roundy, J. K., Vereecken, H.,

575

Walker, J. P., Wood, E. F., Verhoest, N. E. C., and Pauwels, V. R. N. (2016). Assimilation of

576

SMOS soil moisture and brightness temperature products into a land surface model. Remote

577

Sensing of Environment, 180, 292-304, doi: 10.1016/j.rse.2015.10.033.

578

Liu, D., Mishra, A. K., and Yu, Z. (2016). Evaluating uncertainties in multi-layer soil moisture

579

estimation with support vector machines and ensemble Kalman filtering. Journal of Hydrology,

580

538: 243-255.

581

Lü, H., Li, X., Yu, Z., Horton, R., Zhu, Y., Hao, Z., and Xiang, L. (2010). Using a H ∞ filter

582

assimilation procedure to estimate root zone soil water content. Hydrological Processes, 24:

583

3648-3660.

584

Lü, H., Yu, Z., Zhu, Y., Drake, S., Hao, Z., and Sudicky, E. A. (2011). Dual state-parameter estimation

585

of root zone soil moisture by optimal parameter estimation and extended Kalman filter data

586

assimilation. Advances in Water Resources, 34: 395-406.

587 588 589 590

Luo, X., and Moroz, I. M. (2009). Ensemble Kalman filter with the unscented transform. Physica D, 238: 549-562. Ma, L., Scott, H. D., Shaffer, M. D., and Ahuja, L. R. (1998). RZWQM simulations of water and nitrate movement in a tall fescue field. Soil Science, 163: 259–270.

591

Medina, H., Romano, N., and Chirico, G. B. (2014a). Kalman filters for assimilating near-surface

592

observations into the Richards equation - Part 2: A dual filter approach for simultaneous retrieval

593

of states and parameters. Hydrology and Earth System Sciences, 18: 2521-2541.

594

Medina, H., Romano, N., and Chirico, G. B. (2014b). Kalman filters for assimilating near-surface

595

observations into the Richards equation - Part 3: Retrieving states and parameters from laboratory

596

evaporation experiments. Hydrology and Earth System Sciences, 18: 2543-2557.

597 598 599 600

Miller, R. N., Ghil, M., and Ghautiez, F. (1994). Advanced data assimilation in strongly nonlinear dynamic system. Journal of the Atmospheric Sciences, 51: 1037-1055. Moradkhani, H., Sorooshian, S., Gupta, H., and Houser, P. R. (2005). Dual state-parameter estimation of hydrological models using ensemble Kalman filter. Adances in Water Resources, 28: 135-147.

601

Oleson K W, Dai Y, Bonan G, et al. (2004). Technical description of the community land model

602

(CLM). NCAR Technical Note NCAR/TN-461+ STR, National Center for Atmospheric Research,

603

Boulder, CO.

604

Papadakis, N., Mémin, E., Cuzol, A., and Gengembre, N. (2010). Data assimilation with the weighted

605

606

ensemble Kalman filter. Tellus A: Dynamic Meteorology and Oceanography, 62 (5): 673-697. Parinussa, R. M., Holmes, T. R. H., and de Jeu, R. A. M. (2012). Soil moisture retrievals from the

607

WindSat spaceborne polarimetric microwave radiometer. IEEE Transactions on Geoscience and

608

Remote Sensing, 50(7): 2683-2694.

609 610

Richards, L. A. (1931). Capillary conduction of liquids through porous mediums. Journal of Applied Physics, 1(5): 318-333. Doi: 10.1063/1.1745010.

611

Sellers, P. J., Randall, D. A., Collatz, G. J., Berry, J. A., Field, C. B., Dazlich, D. A., Zhang, C.,

612

Collelo, G. D., and Bounoua, L. (1996). A revised land surface parameterization (SiB2) for

613

atmospheric GCMs. part 1: Model formulation. Journal of Climate, 9: 676–705.

614 615 616

Silberstein, R. P., Sivapalan, M., and Wyllie, A. (1999). On the validation of a coupled water and energy balance model at small catchment scales. Journal of Hydrology, 220: 149–168. Van der Merwe, R. (2004). Sigma-Point Kalman filters for probabilistic inference in dynamic

617

state-space

models.

Phd

thesis,

Oregon

Health

and

618

(http://www.cse.ogi.edu/~rudmerwe/pubs/pdf/rvdmerwe_phd_thesis.pdf )

Science

University.

619

Wang, G., Chyi, D., Wang, L., Tan, Y., and Xue, F. (2016). Soil moisture retrieval over Northeast

620

China based on microwave brightness temperature of FY3B satellite and its comparison with

621

other datasets. Chinese Journal of Atmospheric Sciences (in Chinese), 40 (4): 792–804.

622 623 624 625

Western, A. W., and Bloschl, G. (1999). On the spatial scaling of soil moisture. Journal of Hydrology, 217: 203–224. Xie, X., and Zhang, D. (2010). Data assimilation for distributed hydrological catchment modeling via ensemble Kalman filter. Advances in Water Resources, 33(6): 678-690.

626

Xu, S., Yu, Z., Zhang, K., Ji, X., Yang, C., and Sudicky, E.A. (2017). Simulating canopy conductance

627

of the Haloxylon ammodendron shrubland in an arid inland river basin of northwest China.

628

Agricultural and Forest Meteorology, 249: 22–34.

629

Yu, Z., Barron, E. J., Yarnal, B., Lakhtakia, M. N., White, R. A., Pollard, D., and Miller, D. A. (2002).

630

Evaluation of basin-scale hydrologic response to a multi-storm simulation. Journal of Hydrology,

631

257: 212-225.

632

Yu, Z., Carlson, T. N., Barron, E. J., and Schwartz, F. W. (2001). On evaluating the spatial-temporal

633

variation of soil moisture in the Susquehanna River Basin. Water Resources Research, 37:

634

1313-1326.

635

Yu, Z., Fu, X., Lü, H., Luo, L., Liu, D., Ju, Q., Xiang, L, and Wang, Z. (2014a). Evaluating ensemble

636

Kalman, particle and ensemble particle filters through soil temperature prediction. Journal of

637

Hydrologic Engineering, 19, doi:10.1061/(ASCE)HE.1943-5584.0000976.

638

Yu, Z., Fu, X., Luo, L., Lü, H., Ju, Q., Liu, D., Kalin, D. A., Huang, D., Yang, C., and Zhao, L.

639

(2014b). One-dimensional soil temperature simulation with Common Land Model by assimilating

640

in situ observations and MODIS LST with the ensemble particle filter. Water Resources

641

Research, 50: 6950-6965.

642

Yu, Z., Lakhtakia, M. N., Yarnal, B., White, R. A., Miller, D. A., Frakes, B., Barron, E. J., Duffy, C.,

643

and Schwartz, F. (1999). Simulating the river-basin response to atmospheric forcing by linking a

644

mesoscale meteorological model and hydrologic model system. Journal of Hydrology, 218:

645

72-91.

646 647

Yu, Z., Pollard, D., and Cheng, L. (2006). On continental scale hydrologic simulations with a coupled hydrologic model. Journal of Hydrology, 331: 110-124.

648 649

Unscented weighted ensemble Kalman filter for soil moisture assimilation

650

Xiaolei Fu1,2,*, Zhongbo Yu2,*, Yongjian Ding1, Yu Qin1, Lifeng Luo3, Chuancheng Zhao1,

651

Haishen Lü2, Xiaolei Jiang2, Qin Ju2, Chuanguo Yang2

652

1State

653

Resources, Chinese Academy of Sciences, Lanzhou 730000, China

654

2State

655

Nanjing 210098, China

656

3Department

657

Lansing, MI, 48824, USA

658

Corresponding authors: Xiaolei Fu, [email protected]; Zhongbo Yu, [email protected]

Key Laboratory of Cryospheric Science, Northwest Institute of Eco-Environment and

Key Laboratory of Hydrology-Water Resources and Hydraulic Engineering, Hohai University,

of Geography, Environment, and Spatial Sciences, Michigan State University, East

659

Abstract

660

A new data assimilation technique, unscented weighted ensemble Kalman filter (UWEnKF) was

661

developed based on the scaled unscented transformation and ensemble Kalman filter (EnKF). In

662

UWEnKF, the individual members selected are unequally weighted and symmetric about the

663

expectation. To investigate the performance of UWEnKF, nine assimilation experiments with different

664

ensemble sizes (161, 1601, 16001) and different assimilation frequencies (every 6 h, every 12 h, every

665

24 h) were designed to assimilate soil surface (5 cm) moisture data observed at station HY in the upper

666

reaches of the Yellow River, in the northeastern of Tibetan plateau, China into the Richards equation.

667

The results showed that the performance of the filter was greatly affected by random noise, and the

668

filter was sensitive to ensemble size and assimilation frequency. Increasing the ensemble size reduced

669

the effects of random noise on filter performance in several independent assimilation runs (i.e., it

670

decreased the differences between the results of the several independent assimilation runs). Reducing

671

the assimilation frequency also reduced the effects of random noise on filter performance. UWEnKF

672

gave more accurate soil moisture model results than EnKF for all ensemble sizes and assimilation

673

frequencies at all soil depths. Additionally, EnKF may have different performances according to

674

different initial conditions, but not for UWEnKF. Precipitation and soil properties uncertainties had

675

some impact on filter performance. Thus, UWEnKF is a better choice than EnKF, while it is more

676

computationally demanding, for improving soil moisture predictions by assimilating data from many

677

sources, such as satellite-observed soil moisture data, at a low assimilation frequency (e.g., every 24

678

h).

679

Keywords: Soil moisture; Richards equation; ensemble Kalman filter (EnKF); unscented weighted

680

ensemble Kalman filter (UWEnKF)

681 682

Fig. 1. The experiment station HY in Hongyuan, in the upper reaches of the Yellow River, China.

683

Fig. 2. Flowchart of one-dimensional soil moisture assimilation process using (a) EnKF and (b)

684

UWEnKF.

685 686 687 688 689 690 691 692 693 694 695 696 697 698 699

Fig. 3. The 100 RMSE distributions of soil moisture assimilations for EnKF or UWEnKF with different ensemble sizes and different assimilation intervals at different soil depths for 100 independent assimilation runs from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan (Cases 1 - 9). Fig. 4. The 100 MAE distributions of soil moisture assimilations for EnKF or UWEnKF with different ensemble sizes and different assimilation intervals at different soil depths for 100 independent assimilation runs from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan (Cases 1 - 9). Fig. 5. Mean and Std of 100 RMSE values of soil moisture assimilations for EnKF or UWEnKF with different assimilation intervals and different ensemble sizes at different soil depths for 100 independent assimilation runs from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan (Cases 1 - 9). Fig. 6. Mean and Std of 100 MAE values of soil moisture assimilations for EnKF or UWEnKF with different assimilation intervals and different ensemble sizes at different soil depths for 100 independent assimilation runs from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan (Cases 1 - 9).

Fig. 7. Soil moisture assimilation using EnKF and UWEnKF at different soil depths; ensemble size is 16001 and

700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722 723 724 725 726 727 728

assimilation interval is 6 h from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; simulation refers to the model run without assimilation (Case 3). Fig. 8. Soil moisture assimilation using EnKF and UWEnKF at different soil depths; ensemble size is 16001 and assimilation interval is 12 h from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; simulation refers to the model run without assimilation (Case 6). Fig. 9. Soil moisture assimilation using EnKF and UWEnKF at different soil depths; ensemble size is 16001 and assimilation interval is 24 h from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; simulation refers to the model run without assimilation (Case 9). Fig. 10. RMSE and MAE of soil moisture using EnKF and UWEnKF with different assimilation intervals in different soil depths from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; ensemble size is 16001 (Cases 3, 6, 9). Fig. 11. The impact of the error of precipitation on soil moisture predictions using EnKF and UWEnKF with different assimilation intervals in different soil depths from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; ensemble size is 16001 (Cases 10 - 11). Fig. 12. The impact of the error of soil properties on soil moisture predictions using EnKF and UWEnKF with different assimilation intervals in different soil depths from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; ensemble size is 16001 (Cases 12 - 13). Fig. 13, Soil moisture assimilations with initial soil moisture 0.2 m3/m3 using EnKF and UWEnKF at different soil depths; ensemble size is 16001 and assimilation interval is 6 h from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; simulation refers to the model run without assimilation (Case 14). Fig. 14. RMSE and MAE of soil moisture using EnKF and UWEnKF with different assimilation intervals and initial soil moisture 0.2 m3/m3 in different soil depths from 2015-08-01 (day 213) to 2015-09-21 (day 264) of station HY in Hongyuan; ensemble size is 16001 (Case 14).

729 730 731 732

733 734 735 736

Table 1 Soil properties at different soil depths at station HY in Hongyuan, in the upper reaches of the Yellow River, China. Depth(cm)

Sand (%)

Silt(%)

Clay(%)

Bulk density (g/ cm3)

5 20 40 80

52.33 52.24 52.34 70.27

37.96 39.26 37.69 24.02

9.71 8.50 9.97 5.71

0.92 1.06 1.34 1.46

Clay, silt and sand are defined as particles < 0.002, 0.002 - 0.05, and 0.05 - 2 mm in diameter, respectively.

Table 2 Experimental design in the current study. Cases

Initial state values

Assimilation interval (h)

Ensemble size

Error analysis

Soil moisture results

1

161

Figs. 3 - 6

1601

Figs. 3 - 6

3

16001

Figs. 3 - 6, 10

4

161

Figs. 3 - 6

1601

Figs. 3 - 6

6

16001

Figs. 3 - 6, 10

7

161

Figs. 3 - 6

1601

Figs. 3 - 6

16001

Figs. 3 - 6, 10

6

2

5

Measured value

12

24

8 9 737 738 739

Fig. 7

Fig. 8

Fig. 9

Table 3 Experimental design for the impact of initial value, uncertainty of precipitation and soil properties in the current study. Cases

Initial values

Precipitation

Soil properties

Results analysis

10 11 12 13 14

Measured value Measured value Measured value Measured value 0.2 m3/m3

0.95×Measured value 1.05×Measured value Measured value Measured value Measured value

Measured value Measured value 0.95×Measured value (Sand, Clay) 1.05×Measured value (Sand, Clay) Measured value

Fig. 11 Fig. 11 Fig. 12 Fig. 12 Figs. 13 - 14

740 741 742

1. Developed a new filter method--unscented weighted ensemble Kalman filter (UWEnKF).

743

2. UWEnKF is a highly effective and practical assimilation technique

744

3. UWEnKF has better performance than EnKF

745 746 747

Declaration of interests

748 749 750

☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

751 752 753

☐ The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:

754 755 756 757 758