An improved Ensemble Kalman Filter for optimizing parameters in a coupled phosphorus model for lowland polders in Lake Taihu Basin, China

An improved Ensemble Kalman Filter for optimizing parameters in a coupled phosphorus model for lowland polders in Lake Taihu Basin, China

Ecological Modelling 357 (2017) 14–22 Contents lists available at ScienceDirect Ecological Modelling journal homepage: www.elsevier.com/locate/ecolm...

3MB Sizes 4 Downloads 3 Views

Ecological Modelling 357 (2017) 14–22

Contents lists available at ScienceDirect

Ecological Modelling journal homepage: www.elsevier.com/locate/ecolmodel

An improved Ensemble Kalman Filter for optimizing parameters in a coupled phosphorus model for lowland polders in Lake Taihu Basin, China Jiacong Huang, Junfeng Gao ∗ Key Laboratory of Watershed Geographic Sciences, Nanjing Institute of Geography and Limnology, Chinese Academy of Sciences, 73 East Beijing Road, Nanjing 210008, China

a r t i c l e

i n f o

Article history: Received 4 September 2016 Received in revised form 28 April 2017 Accepted 28 April 2017 Keywords: Kalman filter Polder Sensitivity analysis Phosphorus Taihu

a b s t r a c t Ensemble Kalman Filter (EnKF) is potential in optimizing parameters of an environmental model, but may lead to a worse performance of the model in case that improper parameters were updated. To overcome this weakness, EnKF was improved by coupling with a dynamic and multi-objective sensitivity analysis. The improved EnKF was applied to update the parameters of a coupled phosphorus model for simulating phosphorus dynamics of Polder Jian located in Lake Taihu Basin, China. Two parameters that were most sensitive to particulate and dissolved phosphorus were identified at each sub-period, and were then updated using EnKF. To evaluate the performance of the improved EnKF, four simulations with different parameter update strategies were implemented, and compared with measured data. The simulation with the improved EnKF well simulated DP dynamics in Polder Jian with a d value of 0.65 and a RMSE value of 0.015 mg/L. This model fit is better than that of other three simulations with different parameter update strategies, implying a success of the improved EnKF in updating parameters of the coupled phosphorus model. This improved EnKF has the advantage to update several parameters simultaneously, and can be applied in other models with minimal changes. © 2017 Elsevier B.V. All rights reserved.

1. Introduction Data assimilation techniques, such as Variational Algorithms, Ensemble Kalman Filtering (EnKF) and Particle Filtering, are powerful in assimilating various measured data into hydrological, ecological and environmental models to improve model reliability and reduce model uncertainties (Evensen, 2003; Liu et al., 2012; Mo et al., 2008; Zhao et al., 2013). Many special issues were recently published by international journals to discussed the latest advances and case applications of these data assimilation techniques (Chen, 2013; Franssen and Neuweiler, 2015; Seo et al., 2014). Among these data assimilation techniques, the EnKF showed a good performance in solving non-linear problems, and was one of the most applied data assimilation techniques for parameter and variable update in real time (Chen et al., 2015; Hut et al., 2015; Moradkhani et al., 2005; Nie et al., 2011). The EnKF was proposed by Evensen (1994) to overcome the weaknesses of standard Kalman Filter (e.g., large computational

∗ Corresponding author. E-mail address: [email protected] (J. Gao). http://dx.doi.org/10.1016/j.ecolmodel.2017.04.019 0304-3800/© 2017 Elsevier B.V. All rights reserved.

demands and uselessness in nonlinear model). Two steps (forecast and update steps) are compulsively required to apply the EnKF in a model (Gharamti et al., 2015). In the forecast step, an ensemble of state variables is updated using the model. In the update step, the state variables and parameters are updated by analyzing the covariance matrix of state variables and parameters (Evensen, 2003). From theoretical perspective, in case that all the model parameters were involved in the update step, the EnKF can update the sensitive parameter of the model by calculating the covariance between parameters and target state variables, and ignore the parameters non-sensitive to target state variables. However, in practice, the EnKF is difficult to update all the parameters of a complex model for the following two reasons. (1) The complex model may include a large amount (e.g., tens or hundreds) of parameters, and thus requires intensive computation resources in the update step. (2) Each parameter in the update step requires an empirical noise value, defined by modelers (Nie et al., 2011; Thiboult and Anctil, 2015). Defining these noise values for tens/hundreds of model parameter is challenging and fallible. As a compromise, the EnKF was commonly implemented by identifying a sensitive parameter set (ensemble) of the model, and then updating this sensitive parameter set using EnKF (Chen et al., 2015; Su et al., 2011;

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

Sun et al., 2015). Such implementation steps may be improper in two cases: • The sensitive parameters of the model changes through the simulation periods. For example, the simulation phosphorus (P) from a P model is sensitive to a parameter set (P1 ) during a time period (T1 ), but is sensitive to another parameter set (P2 ) during another time period (T2 ). Such dynamic change of sensitive parameters widely existed in complex environmental models (Huang et al., 2016b; Song et al., 2013). • The model included many state variables that are not sensitive to an identical parameter set. For example, a P model includes two variables: particulate phosphorus (PP) and dissolved phosphorus (DP). PP was sensitive to a parameter set (P1 ), and DP was sensitive to another parameter set (P2 ). In above two cases, EnKF may not achieve high performance by updating all the model parameters (P1 and P2 ) all through the simulation period. Therefore, it would be helpful to identify the sensitive parameters in each time period (T1 , T2 , . . ., Tn ) by sensitivity analysis, and update these sensitive parameters in each time period rather than all parameters. This study aimed to improve EnKF for its proper use in above two cases. For this purpose, a dynamic and multi-objective sensitivity analysis was coupled together with EnKF. A coupled phosphorus (P) model that described the P dynamics in Polder Jian, China was used to evaluate the performance of the improved EnKF. The advantages of the improved EnKF, as well as its proper use for an environmental model, were discussed based on the case study. 2. Material and methods 2.1. Model description Polder systems are closed areas enclosed by dikes to protect lowland areas from flood (Brauer, 2014), and are widely distributed in world-wide lowland areas around aquatic systems, such as seas, lakes and rivers (Vermaat and Hellmann, 2010). Phosphorus (P) export from these polders was considered as high P sources due to the intensive agricultural farming (Tian et al., 2012; Vermaat and Hellmann, 2010), and was one of the main causes of severe eutrophication problems in the aquatic systems near polders (Hellmann and Vermaat, 2012). Simulating P dynamics in a polder system is helpful for water managers to identify the P sources, and thus to take corresponding measures to control the P export from this polder (Huang et al., 2016b). A coupled P model developed by Huang et al. (2016a) was used in this study to describe water and P dynamics in a lowland polder system, and to estimate the P export from this polder. The coupled P model included three modes: a Phosphorus Dynamic model for Polders (PDP), Integrated Catchments model of Phosphorus dynamics (INCA-P) and Universal Soil Loss Equation (USLE) with their interaction in Fig. 1. PDP, originally developed by Huang et al. (2016b), described the polder water balance and P dynamics in the water areas of a polder. The model included six modules: water-area water balance module, residential-area water balance module, paddy-land water balance module, dry-land water balance module, water management module and the phosphorus balance module (Fig. 1). These six modules described manually controlled processes of water exchange (e.g., water export by pumping or through a culvert) and P-related processes based on a daily time step. INCA-P, originally developed by Wade et al. (2002), is a dynamic and process-based model, and was used in this study to estimate DP in the runoff water from dry and paddy lands. The INCA-P model described the

15

P dynamics in soils and groundwater and incorporated different P sources (e.g., fertilization and manure) with a daily time step (Whitehead et al., 2014). The USLE model was used in this study to estimate the annual sediment yield. PP concentration in the runoff water can then be estimated based on the calculated annual sediment yield. Calculation details of the USLE model can be found in Neitsch et al. (2005). The coupled P model included 39 parameters, with the sensitive parameters calibrated by Huang et al. (2016a) based on a dataset collected in 2014 from Polder Jian, China. The validation results showed that the model performance was acceptable with a coefficient of determination value of 0.72 for total phosphorus (TP, TP = PP + DP). The trend of phosphorus dynamics in Polder Jian was well simulated. Such performance was acceptable, because P dynamics in Polder Jian was strongly influenced by human activity, and was more difficult to simulate compared with other aquatic ecosystems without strong human activity. Detailed descriptions and calculation equations can be found in Huang et al. (2016a) and Huang et al. (2016b). 2.2. EnKF improvement To apply in a model, standard EnKF was generally implemented in three steps, i.e., parameter identification, forecast and update steps. The parameter identification step aimed to identify the important (sensitive) parameter set in the model, and can be implemented either empirically or based on sensitivity analysis (Huang et al., 2013; Su et al., 2011; Sun et al., 2015). The forecast step updated the state variables in the model, while the update step updated the sensitive parameter set, identified in the parameter identification step (Su et al., 2011; Sun et al., 2015). Compared with the standard EnKF, the EnKF was improved for the coupled P model (Section 2.1) by coupling a dynamic and multi-objective sensitivity analysis. • The dynamic sensitivity analysis. A variance-based global sensitivity analysis approach was used to identify the sensitive parameters in each step when observed data were available for updating parameters. This aimed to guarantee that the update parameter set by EnKF was the most sensitive parameter set in the coupled P model. • The multi-objective sensitivity analysis. The most sensitive parameter to the state variables of PP and DP was identified, respectively. The sensitive parameters to PP and DP were updated by assimilating observed PP and DP, respectively. This aimed to achieve a better simulation of both PP and DP. A schematic illustration of the improved EnKF for updating model parameters in the coupled P model was given in Fig. 2. Four primary steps were required as follows. Step 1. Sensitivity analysis. Two parameters that were most sensitive to PP and DP, respectively, were firstly identified. Parameter sensitivity was compared based on a sensitivity index (STi ) computed by the variance-based sensitivity analysis approach. STi =

EX∼i (VXi (Y |X∼i )) V (Y )

(1)

where Xi is the i-th parameter and X∼i denotes the matrix of all parameters but Xi . Y is the model output. V(Y) is the total unconditional variance, and can be simply estimated by calculating the expected variance of n simulations with random parameter sets that all parameters are not fixed. EX∼i (VXi (Y |X∼i )) is the expected variance that would be left if all factors but Xi could be fixed, and can be estimated based on k*n simulations. k is the number of input parameters. A higher value of STi index implies that the i-th parameter was more sensitive to the target state variable (PP or DP).

16

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

Fig. 1. Conceptual diagram for the tight coupling of the Phosphorus Dynamic model for Polders (PDP), Integrated Catchments model of Phosphorus dynamics (INCA-P) and Universal Soil Loss Equation (USLE) (Huang et al., 2016a). DOP, dissolved organic phosphorus; DIP, dissolved inorganic phosphorus; DP, dissolved phosphorus; PP, particulate phosphorus.

This sensitivity index has been widely used and greatly encouraged for its potential in screening the sensitive parameters (Dai and Ye, 2015; Saltelli and Annoni, 2010; Shahsavani and Grimvall, 2011). Compared with a local sensitivity analysis approach, the global sensitivity analysis approach (variance-based sensitivity analysis) needs more computation resources. However, it was necessary for this study, because the coupled P model (Section 2.1) used in this study was a non-linear model, while the local sensitivity analysis approach was encouraged to applied in a linear model (Saltelli and Annoni, 2010). Further details about this variance-based sensitivity analysis approach can be found in Saltelli and Annoni (2010). Step 2. Ensemble simulations. In the forecast step, the PP and DP for each ensemble member at time T (PiT ) was forecasted by: i+

PiT = f (PiT −T , uiT , ␪T −T ) + Ti , i = 1, ..., n

(2)

where f represents the coupled P model, PiT is the simulation PP and DP at time T, n is the ensemble size, uiT represents the forcing data with a i noise (a zero-mean random variable with a nor-

mal distribution and covariance ), Ti represents model structure i+

uncertainty. ␪T −T represents two parameter ensembles for the i+

most sensitive parameters of PP and DP at time T − T . ␪0 is generated by adding a i noise to the initial values of ␪0 including two parameters that were most sensitive to PP and DP, respectively.

Step 3. Parameter update. When observed data were available, i+ the parameter ensemble (␪T ) was updated by

i+

i+

␪T = ␪T −T + KT

i



i

PiT − Pˆ T

 (3)

where Pˆ T represents two ensembles generated by adding a noise to the observed PP and PP, respectively. KT is Kalman gain for correcting the parameter trajectories. Further calculation details on parameter update using EnKF can be found in many previous studies (Burgers et al., 1998; Evensen, 2009). Step 4. Sensitivity analysis for next simulation. After the parameter update at time T, two parameters that were most sensitive to PP and DP at next time step T + T were identified. If these two sensitive parameters were identical to the updated parami+ eters at time T, the parameter ensemble (␪T ) was used directly at the time step T + T . Otherwise, two new parameter ensembles were generated by adding a noise to their initial values. All parameters except these two sensitive parameters (non-sensitive parameters) were set to their initial values. Such setting for nonsensitive parameters was acceptable because some changes in the non-sensitive parameters would not resulted in large change of simulation results.

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

17

Fig. 2. Schematic illustration of the improved Ensemble Kalman Filter for updating model parameters in the coupled P model.

2.3. Data sources

and water surface). The vegetation coverage data were obtained by monthly survey in the polder. Among the meteorological data, precipitation was measured using an automatic rain gauge located at Polder Jian, while other meteorological data were collected from a national weather station (Liyang Station; No.58345), located 8.5 km southeast of Polder Jian. The water quality data were collected by water and sediment sampling. Further details about this dataset can be found in Huang et al. (2016b).

To evaluate the performance of the improved EnKF, a dataset was collected from Polder Jian, a typical lowland polder located in Lake Taihu Basin, China. The dataset, including land use, vegetation coverage, meteorological and water quality data during Jan. and Dec. 2014 (Table 1), was collected for model inputs. The land use data were obtained from satellite image interpretation, and included four land use types (paddy land, dry land, residential area

Table 1 Data collected for data assimilation of the coupled P model. Type

Indicator

Source

Time period

Temporal resolution

Land use

Land use type Vegetation coverage

Satellite image and surveying Surveying

2014 2014

– Monthly

Meteorology

TMax , TMin , TAve , Wet, WS and HSun Pr

Weather station Automatic rain gauge

2014 2014

Daily Daily

Water quality

DP, PP

Sampling

2014

Twice a month/Monthly



Note: TMax , TMin and TAve : daily maximum, minimum and average of air temperature ( C); Wet: daily average humidity (%); WS: daily average wind speed (m/s); HSun : daily sunshine hours (h); Pr: daily precipitation (mm); DP: dissolved phosphorus concentration (mg/L); PP: particulate phosphorus concentration (mg/L).

18

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22 Table 2 Correlation matrix of simulation PP, DP and TP from four simulations (Sim None, Sim EnKF, Sim EnKF+ and Sim EnKF++). Definitions of these four simulations can be found in Section 2.4. Sim None

Sim EnKF

Sim EnKF+

Sim EnKF++

PP

Sim Sim Sim Sim

None EnKF EnKF+ EnKF++

1.000 0.976 0.963 0.993

0.976 1.000 0.998 0.977

0.963 0.998 1.000 0.966

0.993 0.977 0.966 1.000

DP

Sim Sim Sim Sim

None EnKF EnKF+ EnKF++

1.000 0.979 0.365 0.677

0.979 1.000 0.230 0.587

0.365 0.230 1.000 0.829

0.677 0.587 0.829 1.000

TP

Sim Sim Sim Sim

None EnKF EnKF+ EnKF++

1.000 0.972 0.782 0.950

0.972 1.000 0.647 0.874

0.782 0.647 1.000 0.891

0.950 0.874 0.891 1.000

(d). Lower values of MAE, MAPE and RMSE indicate less error variance, while higher values of R2 and d indicate less error variance. Calculation details of these five measures can be found in Hauduc et al. (2015). 3. Results Fig. 3. Parameter update strategies for four simulations (Sim None, Sim EnKF, Sim EnKF+ and Sim EnKF++).

2.4. Simulations to evaluate the potential of the improve EnKF Four simulations (Sim None, Sim EnKF, Sim EnKF+ and Sim EnKF++) using the coupled P model were implemented and compared to evaluate the performance of the improved EnKF (Fig. 3). The simulation period (T) was divided into 24 sub-periods (T1 , T2 , . . ., T24 ). Parameter settings of these four simulations were described as follows. • Sim None: No parameter was updated during the simulation period. At each sub-period, all the parameter values kept fixed (Fig. 3). Its simulation results were compared with other simulations to evaluate the potential of different parameter update strategies using EnKF. • Sim EnKF: A standard EnKF was used to update the parameter (TTP ) that was most sensitive to TP during the simulation period (T). All parameter values kept fixed at each sub-period except the parameter TTP (Fig. 3). Its simulation results were compared with Sim None to evaluate the potential of the standard EnKF in parameter optimization. • Sim EnKF+: EnKF was used to update two parameters ( PP and T TDP ). These two parameters were most sensitive to PP and DP, respectively, during the simulation period (T). All parameter values kept fixed at each sub-period except the parameter TPP and TDP (Fig. 3). Its simulation results were compared with Sim EnKF to evaluate the value of update more parameters using EnKF. • Sim EnKF++: Different from Sim EnKF+, two parameters ( PP and T

3.1. Evidences in the necessity of the improved EnKF The results of dynamic sensitivity analysis based on Eq. (1) were shown in Fig. 4. Although the coupled P model includes 39 parameters, 23 parameters have a sensitivity index for both PP and DP less than 0.1 at all these 24 sub-periods, and were not included in Fig. 4. The sensitivity analysis results (Fig. 4) showed two direct evidences to couple the dynamic and multi-objective sensitivity analysis with EnKF. • The most influential model parameters for both PP and DP changed through the simulation periods T. For example, p6 was the most influential model parameters for PP during the simulation periods of T1 -T10 , while p5 was the most influential model parameters for PP during the simulation periods of T20 -T24 . p6 had a sensitivity index of for PP larger than 0.9 during T2 and T3 . This revealed that PP was mostly determined by p6 , and was nearly not affected by other parameters during T2 and T3 . This encouraged us to update the most influential model parameters at each step, rather than update an identical parameter through all the simulation period. Dynamic sensitivity analysis can be used to identify the most influential model parameters at each step. • PP and DP had different sensitive parameters during an identical time period. Fig. 4 showed clearly different patterns of the sensitivity indexes for PP and DP. PP dynamics were most strongly affected by two parameters (p5 and p6 ), while DP dynamics were strongly affected by another two parameters (p2 and p12 ). Multiobjective sensitivity analysis was thus needed to identify the most influential model parameters for PP and DP.

i

TDP ) that were most sensitive to PP and DP were identified at i

each sub-period (Ti ). These two parameters were then updated using EnKF. I.e., the update parameters may change through the simulation period T (Fig. 3). Its simulation results were compared with Sim EnKF to evaluate the potential of the improved EnKF in parameter optimization. To quantify the model fits of above four simulations, five widely used measures were used including mean absolute error (MAE), mean absolute percent error (MAPE), root mean square error (RMSE), coefficient of determination (R2 ) and index of agreement

3.2. Performance of the improved EnKF These four simulations (Sim None, Sim EnKF, Sim EnKF+ and Sim EnKF++) provided an acceptable description of PP dynamics with a R2 value ranging from 0.45 to 0.57. Its trend of increase from Jan. to May was well simulated (Fig. 4(a)). Simulation PP from these four simulations showed very slight difference with the correlation coefficient among any two simulations higher than 0.96 (Table 2). Compare with Sim None, Sim EnKF and Sim EnKF+, Sim EnKF++ better described the DP dynamics with a highest d value (0.65) and lowest RMSE value (0.015 mg/L) (Table 3). Sim-

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

19

Fig. 4. Sensitivity index of 16 parameters in the coupled P model during 2014. The simulation period T (the year of 2014) was divided into 24 sub-periods (Ti , i = 1, ..., 24) based on the observation data. 23 parameters with sensitivity index for PP and DP less than 0.1 at these 24 sub-periods were not included. Table 3 Model fits of four simulations (Sim None, Sim EnKF, Sim EnKF+ and Sim EnKF++) in 2014. Definitions of these four simulations can be found in Section 2.4. Model

Items

MAE (mg/L)

MAPE (%)

RMSE (mg/L)

R2

d

Sim None

PP DP TP

0.021 0.017 0.018

30.1 66.2 17.7

0.026 0.020 0.025

0.57 0.00 0.67

0.85 0.41 0.89

Sim EnKF

PP DP TP

0.019 0.017 0.023

27.5 63.4 20.2

0.028 0.020 0.030

0.48 0.00 0.54

0.80 0.42 0.81

Sim EnKF+

PP DP TP

0.020 0.035 0.035

27.9 142.1 31.7

0.028 0.051 0.048

0.45 0.22 0.55

0.79 0.40 0.78

Sim EnKF++

PP DP TP

0.020 0.013 0.017

28.5 54.2 17.3

0.027 0.015 0.023

0.53 0.19 0.70

0.84 0.65 0.91

ulation DP from Sim EnKF and Sim None were similar with a high correlation coefficient of 0.979 (Table 2). Sim EnKF+ failed to simulate the DP dynamics from Sep. to Dec. 2014 with a RMSE value as high as 0.035 mg/L (Table 3). All the five model fit value in Table 3 showed that Sim EnKF++ provided a slightly better description of TP. Parameter dynamics in Sim EnKF++ showed that the parameter value varied significantly during the simulation period (T1 -T24 ) with a ␪Ti /␪T1 value ranging from 0.5-1.5. Two parameters (p1 and p2 ) in Sim EnKF++ using the improved EnKF were significantly updated during T21 to T24 (Fig. 6). This update resulted in a better description of DP dynamics during this period (Fig. 5(b)).

4. Discussion 4.1. Advantages of the improved EnKF The improved EnKF was successfully applied in the coupled P model for parameter optimization, and demonstrated several advantages in the process of parameter optimization.

(1) The dynamic and multi-objective sensitivity analysis guaranteed that proper parameters can be updated by the improved EnKF in each sub-period. This advantage is particularly useful for a complex model with tens/hundreds parameters, because the sensitive parameters for such complex model are most likely to change through time, and to be sensitive to different parameters for various state variables. (2) Based on the multi-objective sensitivity analysis strategy, several sensitive parameters in the coupled P model can be updated simultaneously. Although the coupled P model targeted at a better description of TP, the sensitivity analysis results (Fig. 4) revealed that it is necessary to screen the model parameters for the detailed components of TP (i.e., DP and PP). The case study updated only two parameters by assimilating the measured PP and DP data, more parameters can be updated using this improve EnKF with more measured data. For example, the measured water level can be used to update the sensitive parameter in the water-area water balance module in PDP (Fig. 1). This advantage allows the improved EnKF to make better use of various measured data in a complex model. (3) The improved EnKF was loosely coupled with the coupled P model. I.e., the information exchange between the improved

20

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

Fig. 5. Phosphorus from measurement data and simulation results from four simulations (Sim None, Sim EnKF, Sim EnKF+ and Sim EnKF++) in 2014. Definitions of these four simulations can be found in Section 2.4.

EnKF and coupled P model was based on files. Such loose coupling allowed the improved EnKF to be applied in other environmental models with minimal changes in model structure.

4.2. Proper use of the improved EnKF The case study demonstrated the potential of the improved EnKF in parameter optimization of the coupled P model, and provided an implementation details for environmental modelers. However,

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

21

(3) The improved EnKF cannot solve all sources of model uncertainty, but only solve the model uncertainty from parameters. The improved EnKF (Sim EnKF++) did not result in a significantly better simulation of PP and TP (Fig. 5 and Table 3). This is probably because that other uncertainty sources than parameter uncertainty, such as input data and model structure, existed in the coupled P model, and were not solved in the improved EnKF. The DP dynamics were not well simulated compared with PP dynamics (Table 3), implying that the processes (e.g., P uptake by phytoplankton and transformation between inorganic and organic P) related to DP dynamics were not adequately described. More detailed descriptions of these processes can be further implemented based on the existing models for aquatic ecosystems, such as Environmental Fluid Dynamics Code and PCLake (Janse et al., 2010; Tetra Tech, 2007). 5. Conclusions EnKF was improved by coupling with a dynamic and multiobjective sensitivity analysis, and was applied for updating parameter of a coupled P model for simulating P dynamics in Polder Jian, China during 2014. Compared with the simulation using the standard EnKF, the simulation using the improved EnKF better simulated the phosphorus dynamics, especially the DP dynamics, implying the potential of the improved EnKF in parameter optimization by assimilating observed data. This improved EnKF coupled a dynamic and multi-objective sensitivity analysis to identify the most sensitive parameters in each time step, and thus guaranteed that proper parameters can be updated in real time. It can be easily adapted for other environmental models due to its loose coupling with the coupled P model, and can be potentially used to improve model performance by investigating parameter trajectories through time. Acknowledgments

Fig. 6. Parameter dynamics of three simulations (Sim EnKF, Sim EnKF+ and Sim EnKF++) during 2014 represented by ␪T /␪T1 value. ␪T is the parameter value i

i

at Ti . ␪T1 is the initial value of the parameter. The simulation period T (the year of 2014) was divided into 24 sub-periods (Ti , i = 1,..., 24) based on the observation data. 31 parameters without any update were not included.

The project was financially supported by Natural Science Foundation of Jiangsu, China (BK20161614) and National Natural Science Foundation of China (41301574). The authors would like to thank China Meteorological Data Sharing Service System for providing the measured data for model development. Special thanks to Dr. Hongbin Yin, Dr. Wei Huang, Dr. Qi Huang and Mr. Desheng Zhu for their helps on the samplings in Polder Jian. References

several issues should be emphasized to achieve a better application of the improved EnKF. (1) Compared with the standard EnKF, the improved EnKF need more computation resources due to its additional sensitivity analysis in each update step. It was not necessarily needed for all environmental models, but limited to complex models with dynamic parameter sensitivity as shown in Fig. 4. For a simple model with few parameters, the standard EnKF can update all parameters, and does not required the additional sensitivity analysis to identify the sensitive parameters in each update step. Necessity of the improved EnKF should be adequately evaluated based on dynamic sensitivity analysis of the model as shown in Fig. 4. (2) Proper parameters should be well identified and then updated by the improved EnKF. The potential of the improved EnKF attributed to its power in identifying the sensitive parameter in each sub-period. An improper update of model parameters may even result in a worse fit as shown in the failure of Sim EnKF+ in predicting the DP dynamics (Fig. 5(b)). This is because the update parameters were not sensitive to DP, but was updated based on the error covariance of DP in EnKF.

Brauer, C.C., 2014. Modelling Rainfall-runoff Processes in Lowland Catchments. Ph.D. Thesis. Wageningen University, the Netherlands. Burgers, G., Jan van Leeuwen, P., Evensen, G., 1998. Analysis scheme in the Ensemble Kalman Filter. Mon. Weather Rev. 126 (6), 1719–1724. Chen, W., Huang, C., Shen, H., Li, X., 2015. Comparison of ensemble-based state and parameter estimation methods for soil moisture data assimilation. Adv. Water Resour. 86, 425–438. Chen, Y., 2013. Ensemble Kalman Filter for data assimilation. Comput. Geosci. 55, 1–2. Dai, H., Ye, M., 2015. Variance-based global sensitivity analysis for multiple scenarios and models with implementation using sparse grid collocation. J. Hydrol. 528, 286–300. Evensen, G., 1994. Sequential data assimilation with a nonlinear quasi-geostrophic model using Monte Carlo methods to forecast error statistics. J. Geophys. Res. 99 (C5), 10143–10162. Evensen, G., 2003. The Ensemble Kalman Filter: theoretical formulation and practical implementation. Ocean Dyn. 53 (4), 343–367. Evensen, G., 2009. Data Assimilation: The Ensemble Kalman Filter. Springer Verlag, Bergen. Franssen, H.J.H., Neuweiler, I., 2015. Data assimilation for improved predictions of integrated terrestrial systems. Adv. Water Resour. 86, 257–259. Gharamti, M.E., AitElFquih, B., Hoteit, I., 2015. An iterative Ensemble Kalman Filter with one-step-ahead smoothing for state-parameters estimation of contaminant transport models. J. Hydrol. 527, 442–457. Hauduc, H., Neumann, M.B., Muschalla, D., Gamerith, V., Gillot, S., Vanrolleghem, P.A., 2015. Efficiency criteria for environmental model quality assessment: a

22

J. Huang, J. Gao / Ecological Modelling 357 (2017) 14–22

review and its application to wastewater treatment. Environ. Modell. Softw. 68, 196–204. Hellmann, F., Vermaat, J.E., 2012. Impact of climate change on water management in Dutch peat polders. Ecol. Model. 240, 74–83. Huang, J., Gao, J., Liu, J., Zhang, Y., 2013. State and parameter update of a hydrodynamic-phytoplankton model using Ensemble Kalman Filter. Ecol. Model. 263, 81–91. Huang, J., Gao, J., Yan, R., 2016a. How can we reduce phosphorus export from lowland polders?: Implications from a sensitivity analysis of a coupled model. Sci. Total Environ. 562, 946–952. Huang, J., Gao, J., Yan, R., 2016b. A phosphorus dynamic model for lowland polder systems (PDP). Ecol. Eng. 88, 242–255. Hut, R., Amisigo, B.A., Steele-Dunne, S., van de Giesen, N., 2015. Reduction of Used Memory Ensemble Kalman Filtering (RumEnKF) A data assimilation scheme for memory intensive, high performance computing. Adv. Water Resour. 86, 273–283. Janse, J.H., Scheffer, M., Lijklema, L., Van Liere, L., Sloot, J.S., Mooij, W.M., 2010. Estimating the critical phosphorus loading of shallow lakes with the ecosystem model PCLake: sensitivity, calibration and uncertainty. Ecol. Model. 221 (4), 654–665. 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., Restrepo, P., 2012. Advancing data assimilation in operational hydrologic forecasting: progresses, challenges, and emerging opportunities. Hydrol. Earth Syst. Sci. 16 (10), 3863–3887. Mo, X., Chen, J.M., Ju, W., Black, T.A., 2008. Optimization of ecosystem model parameters through assimilating eddy covariance flux data with an Ensemble Kalman Filter. Ecol. Model. 217 (1–2), 157–173. Moradkhani, H., Sorooshian, S., Gupta, H.V., Houser, P.R., 2005. Dual state-parameter estimation of hydrological models using Ensemble Kalman Filter. Adv. Water Resour. 28 (2), 135–147. Neitsch, S., Arnold, J., Kiniry, J., Williams, J., King, K., 2005. SWAT theoretical documentation. Soil Water Res. Lab.: Grassl. 494, 234–235. Nie, S., Zhu, J., Luo, Y., 2011. Simultaneous estimation of land surface scheme states and parameters using the Ensemble Kalman Filter: identical twin experiments. Hydrol. Earth Syst. Sci. 15 (8), 2437–2457.

Saltelli, A., Annoni, P., 2010. How to avoid a perfunctory sensitivity analysis. Environ. Modell. Softw. 25 (12), 1508–1517. Seo, D.J., Liu, Y., Moradkhani, H., Weerts, A., 2014. Ensemble prediction and data assimilation for operational hydrology. J. Hydrol. 519, 2661–2662. Shahsavani, D., Grimvall, A., 2011. Variance-based sensitivity analysis of model outputs using surrogate models. Environ. Modell. Softw. 26 (6), 723–730. Song, X., Bryan, B.A., Almeida, A.C., Paul, K.I., Zhao, G., Ren, Y., 2013. Time-dependent sensitivity of a process-based ecological model. Ecol. Model. 265, 114–123. Su, H., Yang, Z.L., Niu, G.Y., Wilson, C.R., 2011. Parameter estimation in ensemble based snow data assimilation: a synthetic study. Adv. Water Resour. 34 (3), 407–416. Sun, L., Nistor, I., Seidou, O., 2015. Streamflow data assimilation in SWAT model using extended kalman filter. J. Hydrol. 531, 671–684. Tetra Tech, Inc., 2007. The Environmental Fluid Dynamics Code: Theory and Computation. US EPA, Fairfax, VA. Thiboult, A., Anctil, F., 2015. On the difficulty to optimally implement the Ensemble Kalman Filter: an experiment based on many hydrological models and catchments. J. Hydrol. 529 (Part 3), 1147–1160. Tian, P., Zhao, G., Li, J., Gao, J., Zhang, Z., 2012. Integration of monthly water balance modeling and nutrient load estimation in an agricultural catchment. Int. J. Environ. Sci. Technol. 9 (1), 163–172. Vermaat, J.E., Hellmann, F., 2010. Covariance in water-and nutrient budgets of Dutch peat polders: what governs nutrient retention? Biogeochemistry 99 (1–3), 109–126. Wade, A., Whitehead, P., Butterfield, D., 2002. The integrated catchments model of phosphorus dynamics (INCA-P), a new approach for multiple source assessment in heterogeneous river systems: model structure and equations. Hydrol. Earth Syst. Sci. 6 (3), 583–606. Whitehead, P.G., Jin, L., Crossman, J., Comber, S., Johnes, P.J., Daldorph, P., Flynn, N., Collins, A.L., Butterfield, D., Mistry, R., Bardon, R., Pope, L., Willows, R., 2014. Distributed and dynamic modelling of hydrology, phosphorus and ecology in the Hampshire Avon and Blashford Lakes: Evaluating alternative strategies to meet WFD standards. Sci. Total Environ. 481, 157–166. Zhao, Y., Chen, S., Shen, S., 2013. Assimilating remote sensing information with crop model using Ensemble Kalman Filter for improving LAI monitoring and yield estimation. Ecol. Model. 270, 30–42.