Research Project on Multi-objective Optimisation: GLOWA

Description
Multiobjective optimization (also known as multiobjective programming, vector optimization, multicriteria optimization, multiattribute optimization or Pareto optimization) is an area of multiple criteria decision making, that is concerned with mathematical optimization problems involving more than one objective function to be optimized simultaneously.

Advances in Geosciences, 5, 19–23, 2005 SRef-ID: 1680-7359/adgeo/2005-5-19 European Geosciences Union © 2005 Author(s). This work is licensed under a Creative Commons License.

Advances in Geosciences

Using multi-objective optimisation to integrate alpine regions in groundwater ?ow models
V. Rojanschi, J. Wolf, R. Barthel, and J. Braun Institut f¨ ur Wasserbau, Universit¨ at Stuttgart, Pfaffenwaldring 61, Stuttgart, Germany Received: 7 January 2005 – Revised: 1 August 2005 – Accepted: 1 September 2005 – Published: 16 December 2005

Abstract. Within the research project GLOWA Danube, a groundwater ?ow model was developed for the Upper Danube basin. This paper reports on a preliminary study to include the alpine part of the catchment in the model. A conceptual model structure was implemented and tested using multi-objective optimisation analysis. The performance of the model and the identi?ability of the parameters were studied. A possible over-parameterisation of the model was also tested using principal component analysis.

2

Modelling structure

1 Introduction Within the framework of the GLOWA-Danube project (Mauser and Barthel, 2004), a groundwater ?ow model was developed for the Upper Danube basin. The model is coupled to a soil water balance model and a hydraulic surface water model. From the former it receives the in?ltration rate through the lower boundary, de?ned at two meters below the land surface at every grid cell. To the latter it delivers a water exchange rate between the aquifers and the surface waters for every river cell. Wolf et al. (2004) reported in detail on the chosen hydrogeological conceptual model, on the dif?culties encountered during the work on the numerical ?ow model (MODFLOW, McDonald and Harbaugh, 1988) and on the solutions found for these dif?culties. The presence of the alpine region in the south of the Upper Danube basin posed a consistency problem between the groundwater and the soil models. Due to their steep, folded and faulted internal structure, the Alps, with the exception of the alluvial aquifers in valleys, are not compatible with the Darcy-Law based MODFLOW approach. A solution had to be found to ?ll the gap between the two models in the alpine part of the catchment. Correspondence to: V. Rojanschi ([email protected])

The task is to develop a model for the subsurface ?ow in the alpine regions which should link the soil water model, concerned with the ?rst two meters of soil and the groundwater model dealing with the ?ow in the alluvial valley aquifers. The absence of deterministic information regarding the fractures dominated subsurface ?ow beneath the mountain slopes obliges the use of a conceptual hydrological approach based on a qualitative description of the involved processes. The proposed modelling structure is presented in Fig. 1. Based on the existing river gauges, alpine subcatchments have been delineated. For every subcatchment the in?ltration computed by the soil water model is split into two parts: the in?ltration above alluvial valleys and above mountains slopes. First, the in?ltration above the alluvial valleys aquifers is injected as vertical groundwater recharge into the MODFLOW model. Second, the in?ltration above the mountainous slopes is aggregated over the subcatchment and again separated in two parts. The ?rst part, named inter?ow in the context of this paper, ex?ltrates along the slope to ?ow directly into the river network. The second part ?ows through the mountain to ex?ltrate in the alluvial aquifer as lateral groundwater recharge. The water ex?ltrating from the valley aquifer into the rivers is named here base?ow. For the water routing through the individual components, conceptual modelling units were used based on the linear storage cascade concept (Nash, 1959). Each of the storage cascades is de?ned by two parameters, namely the number of reservoirs n and the reservoir coef?cient k . The parameter s, determining the separation between inter?ow and base?ow, and the parameters of the storage cascades units ni , ki are being quanti?ed during the calibration process because no direct physically-based information is available for that purpose.

20

Figures
Soil Water Model
Infiltration above valley aquifers
Soil Water Model

V. Rojanschi et al.: Alpine regions in groundwater ?ow models

Figures

Infiltration above mountain slopes

Surface Runoff

Unit 3
Surface Runoff

s

Unit 1

Infiltration above valley aquifers

Infiltration above mountain slopes

Surface Runoff

Unit 3

s

Unit 1

Surface Runoff

Unit 2

Interflow

measured Discharge River Model

Unit 2

Interflow

Fig. 1: Structure of the proposed conceptual model integrating the Alps in the hydrological Fig. 1. Structure of the proposed conceptual model integrating the Alps in the groundwater model. Fig. 1: Structure of the proposed conceptualcomplex. model integrating the Alps in the hydrological modelling
modelling complex.

Unit 4 Finite-Difference Groundwater Flow Model

Finite-Difference River Model Groundwater Flow Model
Baseflow

measured Unit Discharge

4

Baseflow

a problem solved. Even when using state-of-the-art automatic calibration algorithms one cannot avoid the problems generated by the numerous local optima with very similar performance criteria values (see Duan et al., 1993), by the subjectivity involved in the selection of this criterion and by the “equi?nality” issue (the existence of many parameter sets leading to almost equally good model results, see Beven, 2000). There are several possible answers to these problems, which do not oppose, but rather complement each other. One answer is a generalised sensitivity analysis (Hornberger and Spear, 1981) or the GLUE approach derived from it (Beven and Binley, 1992). Another answer, the one applied here, is the use of a multiobjective calibration as opposed to a simple one objective Fig. 2. The the Pareto in the space of two objective orientated towards upper set rightrepresented corner, which represents the perfect fit. The circles mark calibration (see Gupta et al., 2003). No objective function functions. The Pareto front is orientated towards the upper right the single-objective optima. characterises in an exhaustive manner the quality of the ?t becorner, which represents the perfect ?t. The circles mark the singletween the measured (Qmes ) and the computed (Qsim ) time objective optima. series. For a long time, this was the main argument in favour Fig. 2. The Pareto set represented in the space ofof two objective The Pareto front the manual functions. calibration. Although lessis systematic, less reproducible and much more time consuming, towards upper right corner, which represents the perfect fit. The circles mark the manual calThe task to be orientated solved now is to the determine physicallyibration had the advantage of being able to lead to an opinterpretable sensitive parameters (relative to the available the single-objective optima. timum which takes more than one mathematical expression data) for the model structure presented in Fig. 2. In a prefor the quality of the ?t into consideration. This weakness liminary approach presented here, the MODFLOW model was solved for the automatic approach by the use of multiwas also replaced with a linear storage cascade. The proobjective calibration procedures. The result of the calibration posed testing procedure for the model structure requires a is in this case no longer one single parameter set, but a group large number of model evaluations and would not be appliof parameter sets, termed Pareto sets, which optimise as a cable with the MODFLOW model due to the needed CPU group several prede?ned objective functions. Having a range time. of optimal parameter sets and optimal model results offers an additional advantage. Through the analysis of the spread of these ranges, one can quantify in a more objective manner 3 Methodology the degree of con?dence that one should have in the given model. The dangerous feeling of certainty, which the modThe developments in the last decade in the ?eld of hydrological modelling have made clear that good ?ts between eller has when dealing with one optimal parameter set as a measured and simulated discharge curves, evaluated using ?nal answer of the problem, is thus at least partially eliminated. one performance criteria, are by far not enough to consider
Fig. 2. The Pareto set represented in the space of two objective functions. The Pareto front is

Table 1: Optimal values for the five criteria used in the multi-objective analysis. All five criteria take values in the interval (-?:1], with 1 indicating a perfect fit. The seven subcatchments were sorted according to the quality of the results. 12 cm wide

V. Rojanschi et al.: Alpine regions in groundwater ?ow models

21

Table 1. Optimal values for the ?ve criteria used in the multi-objective analysis. All ?ve criteria take values in the interval (??:1], with 1 indicating a perfect ?t. The seven subcatchments were sorted according to the quality of the results.

Subcatch 1 Subcatch 2 Subcatch 6 Subcatch 3 Subcatch 4 Subcatch 5 Subcatch 7

NS 0.952 0.950 0.945 0.916 0.757 0.568 0.406

Calibration values NSdr NSndr SAE 0.893 0.961 0.941 0.920 0.920 0.931 0.939 0.926 0.885 0.888 0.830 0.853 0.767 0.615 0.731 0.565 0.318 0.605 0.231 0.189 0.465

NStr 0.977 0.957 0.930 0.889 0.689 0.516 0.373

NS 0.978 0.965 0.937 0.905 0.615 0.486 0.285

Validation values NSdr NSndr SAE 0.967 0.980 0.938 0.947 0.957 0.923 0.925 0.926 0.868 0.863 0.895 0.855 0.502 0.588 0.728 0.431 0.484 0.588 0.032 0.224 0.404

NStr 0.980 0.967 0.927 0.889 0.671 0.491 0.279

For the case presented in this paper ?ve objective functions were selected for a multi-objective analysis (see Freer et al, 2003): the Nash-Sutcliffe ef?ciency (N S ), Nash-Sutcliffe ef?ciencies computed for the increasing and the decreasing part of the hydrographs (NSdr , NSndr ), SAE =1?S , where S is the sum of the absolute differences between Qmes and Qsim normalised by the sum of Qmes , and the Nash-Sutcliffe ef?ciency computed between Qmes and Qsim after applying a Box-Cox transformation (N Str ). The ?ve functions were computed for four time scales (one, two, seven and thirty days) and the average of the four values was used in the optimisation procedure. The computed Pareto sets were composed of 495 parameter sets, respecting the recommendation given by Gupta et al. (2003) of having around 500 values. 4 Test area: the Ammer catchment

The Ammer catchment, (709 m2 ), located in the southwestern corner of Bavaria upstream of the Ammer lake, was chosen as a test area. Apart from the representativeness of the catchment for the transition zone between the alpine formations and the molasse zone, the choice was also motivated by the very good data availability. Seven subcatchments could be de?ned based on the existing river gauges. The analysed time period was 01.11.1990–01.01.2000. The last seven years of the time series were used for the calibration, the ?rst three years were used for the validation. The necessary input data, the in?ltration rate and the surface runoff, were calculated using the PROMET soil water balance model (Mauser, 1989) and were made available to the authors of this paper by Dr. Ralf Ludwig from the Ludwig Maximilian University in Munich. 5 Results and discussion

The multi-objective analysis was applied on the seven subcatchments of the Ammer catchment. Figure 2 shows the Pareto solution in the criteria space for one pair of objective functions. The Pareto front is clearly de?ned as well as the position of the single – objective optimums at the edge of the Pareto front. It is an indication that the objective func-

tions were chosen correctly. Other performance criteria were tested before selecting the ?ve functions previously mentioned. The root mean square error and the heteroscedastic maximum likelihood estimator proposed by various authors correlated to the Nash-Sutcliffe ef?ciency for this study case, so that the Pareto set was concentrated on the y =x line, thus adding no additional information to the analysis. The optimal values for the ?ve performance criteria, averaged over the four time scales already mentioned, are presented in Table 1 for both the calibration and validation periods. For four of seven subcatchments the results can be quali?ed as very good, with all performance criteria having values between 0.85 and 0.98. The other three subcatchments make a distinct picture, two of them (4 – gauge Oberammergau and 5 – gauge Unternogg) having average results and the third (7 – gauge Obernach) being at the limit between poor and unacceptable. Figure 3, presenting the measured versus the computed time series for the subcatchments with the best (1) and the worst (7) results, is helpful for explaining the poor results for the subcatchments 4, 5 and 7. By comparing the direct results of the soil water model – the sum between the in?ltration rate and the surface runoff – with the measured river discharges, it is noticeable that the soil water model already has attenuated the initial rain signal too much. As the transport model here discussed can only transport the input signal forwards and increase its attenuation, there is no space for it to improve the results of the soil water model, which in this particular case would require a backwards transformation and a de-attenuation. It is interesting to notice that the three subcatchments, whose results are not satisfactory, are situated furthest upstream and are characterised by large altitude differences. The proper parameterisation of the soil layer is an extremely dif?cult process when it comes to very steep mountain slopes. The interpolated rain time series are also affected by a signi?cant degree of uncertainty, although the correlation between rain and elevation was taken into consideration during the interpolation process (Ludwig, 2000). For subcatchment 1, Fig. 3 con?rms in a graphical form the very good ?t between the measured and the computed time series. There is a slight tendency to underestimate the highest peaks, but otherwise the computed Pareto solutions

22

V. Rojanschi et al.: Alpine regions in groundwater ?ow models

Figure 3. Computed versus measured time series for best the subcatchment with best To (1)explain and the poor results in Fig. 3. Computed versus measured time series for the subcatchment with the (1) and the worst (7)the results. subcatchment 7 the model input (the results of the soil water model) was also represented. the worst (7) results. To explain the poor results in subcatchment 7 the model input (the results of the soil water model) was also represented.

are able to reproduce the dynamics of the measured discharge curve well. Notice should be also given to the thin range of values characterising the Pareto solutions. Although the Nash-Sutcliffe ef?ciency for this subcatchment is 0.95, the measured discharge is located inside the interval de?ned by the Pareto set in only 53% of the days from the calibration period and 48% of the days from the validation period. One goal of the multi-objective analysis, namely to create a solution set fully “including” the measured points (Gupta, 2003), could thus not be achieved. This result was also noticed during other studies and was used to criticise the overinterpretation of the Pareto set as a measure of the uncertainty of a model (Freer et al., 2003). In addition to testing model performance, it is important to test whether the calibrated parameters are more than the results of a mathematical optimisation and can be interpreted in a physical way. Figure 4a shows the distribution of the 495 normalised Pareto solutions for the eleven model parameters for subcatchment 1. During the calibration, the parameters were restricted to positive values. The upper bound was imposed by restricting the two coef?cients for every storage cascade unit and also their product (which is the time distance with which the centre of gravity of the input signal is translated into the output signal) to values predetermined on the basis of hydrograph separation methods (Schwarze et al., 1991). The normalised values of the Pareto parameters in Fig. 4a show no trend, as it seems that combinations of values

throughout the whole allowed spectrum were obtained after the optimisation process. The storage cascades’ coef?cients are agglomerated into the lower range only due to the forced restriction for every ni ? ki product to an upper bound. One possible explanation of the poor identi?ability of the results is the over-parameterisation of the model and the parameter interdependence that comes with it. To test this hypothesis, the correlation matrix inside the Pareto set was computed and its eigenvalues and eigenvectors determined (principal component analysis, see Bishop, 1995). The analysis lead to relatively high correlation coef?cients and to few dominant eigenvalues, strongly suggesting that the parameters are compensating each other in the optimisation process. Figure 4b shows the Pareto set tranformed into the eigenvectors. Although a certain degree of variability remains, the values are clearly de?ned, proving the over-parameterisation hypothesis. It is also worth mentioning that this analysis “catches” the linear interdependences only, which means that the computed number of needed independent parameters (the number of dominant eigenvalues) is certainly overestimated. Additional studies are needed to determine the non-linear dimensionality of the system. 6 Conclusions

For the integration of the alpine part of the catchment into a regional groundwater ?ow model, a conceptual model structure was tested using a multi-objective optimisation analysis.

V. Rojanschi et al.: Alpine regions in groundwater ?ow models

23

Fig. 4. (a) Normalised values of the Pareto set’s parameters; (b) Normalised values of the Pareto set’s eigenvectors. v 11 is the eigenvector Pareto set‘s eigenvectors. v11 is the explaining the smallest amount of parameter variability, v 1 eigenvector the largest. explaining the smallest amount of parameter

Figure 4. a. Normalised values of the Pareto set‘s parameters; b. Normalised values of the

variability, v1 the largest. For most of the subcatchments of the test area, a good performance was achieved. The optimised parameters were poorly de?ned, and a clear over-parameterisation was identi?ed which lead to strong correlations and compensation in the parameter space. Further studies are planed to test whether the inclusion of the groundwater model resolves this issue or whether a rethinking of the structure is needed.
Acknowledgement. The authors acknowledge the help of R. Ludwig, Ludwig Maximilian University in Munich, he provided us with measured data and model results, without which the work presented here would not have been possible. Edited by: P. Krause, K. Bongartz, and W.-A. Fl¨ ugel Reviewed by: anonymous referees

References
Beven, K. J.: Rainfall-Runoff Modelling: the Primer, John Willey & Sons Ltd, Chichester, England, 2000. Beven, K. and Binley, A.: The Future of Distributed Models: Model Calibration and Uncertainty Prediction, Hydrol. Processes, 6, 279–298, 1992. Bishop, C. M.: Neural Networks for Pattern Recognition, Oxford University Press, 1995. Duan, Q., Gupta, V. K., and Sooroshian, S.: A Shuf?ed Complex Evolution Approach for Effective and Ef?cient Global Minimization, Journal of Optimization Theory and its Applications, 76, 501–521, 1993. Freer, J., Beven, K., and Peters, N.: Multivariate Seasonal Period Model Rejection within the Generalised Likelihood Uncertainty Estimation Procedure, Water Science and Application, 6, 9–29, 2003.

Gupta, H. V., Bastidas, L. A., Vrugt, J. A., and Sorooshian, S.: Multiple Criteria Global Optimization for Watershed Model Calibration, Water Science and Application, 6, 125–132, 2003 Hornberger, G. M. and Spear, R. C.: An Approach to the Preliminary Analysis of Environmental Systems, J. Environ. Manag., 12, 7–18, 1981. Ludwig, R.: Die ?¨ achenverteilten Modellierung von Wasserhaushalt und Ab?ussbildung im Einzugsgebiet der Ammer, M¨ unchener Geographische Abhandlungen, Reihe B 32, M¨ unchen, 2000. Mauser, W.: Die Verwendung hochau?¨ osender Satellitendaten in einem Geographischen Informationssystem zur Modellierung von Fl¨ achenverdunstung und Bodenfeuchte, Habilitationsschrift, Albert-Ludwigs-Universit¨ at, Freiburg i. Br., 1989. Mauser, W. and Barthel, R.: Integrative Hydrologic Modeling Techniques for Sustainable Water Management regarding Global Environmental Changes in the Upper Danube River Basin, in: Research Basins and Hydrological Planning, edited by: Xi, R.-Z., Gu, W.-Z., and Seiler, K.-P., A. A. Balkema Publishers, Rotterdam, The Netherlands, Brook?eld, USA, 53–61, 2004. McDonald, M. G. and Harbaugh, A. W.: A Modular ThreeDimensional Finite-Difference Ground-Water Flow Model: U.S. Geological Survey Techniques of Water-Resources Investigations, book 6, chap. A1, Washington, USA, 1988. Nash, J. E.: Systematic Determination of Unit Hydrograph Parameters, J. Geophys. Res., 64, 111–115, 1959. Schwarze, R., Herrmann, A., M¨ unch, A., Gr¨ unewald, U., and Sch¨ oniger, M.: Rechnergest¨ utzte Analyse von Ab?ußkomponenten und Verweilzeiten in kleinen Einzugsgebieten, Acta hydrophis., 35/2, 143–184, 1991. Wolf, J., Rojanschi, V., Barthel, R., and Braun, J.: Modellierung der Grundwasserstr¨ omung auf der Mesoskala in geologisch und geomorphologisch komplexen Einzugsgebieten, in: 7. Workshop zur großskaligen Modellierung in der Hydrologie, edited by: Ludwig, R., Reichert, D., and Mauser, W., Kassel University Press GmbH, 155–162, 2004.



doc_420915040.pdf
 

Attachments

Back
Top