Modelling Hydrological Characteristics Surface Runoff

Published Date: 02 Nov 2017

Short title (<70 characters): Modelling semi-arid catchments hydrological response to wildfires effects

Elias Moussoulis1,2, Giorgos Mallinis3, Nikos Koutsias1 and Ierotheos Zacharias1[����1]*

1 Department of Environmental and Natural Resources Management, University of Western Greece, Seferi 2, GR-30100 Agrinio, Greece

2 DHI Greece, 114 Thiras, Argyroupoli, Athens 16451

3 Department of Forestry & Management of the Environment and Natural Resources, Democritus University of Thrace, Pantazidou 193, GR-68200, Orestiada, Greece

* Author to whom correspondence should be addressed; Tel.: +30 26410 74131;

Fax: +30 26410 74170; E-Mail: [email protected]

Abstract:

Wildfires change the infiltration properties of soils, reduce the amount of interception materials that protect soil and slow runoff. A wildfire at Northeast Attiki, Central Greece, in August 2009, destroyed approximately one third of the study area consisting of a mixture of shrublands, pastures, and pines. The present study can be thought of as a simultaneous modelling assessment of multiple semi-arid, shrubland-dominated Mediterranean catchments� hydrologic response (total, annual and mean monthly runoff) the first few years to following wildfires, during their first few years. The modelling framework that was chosen in this study [EMS2]was MIKE SHE, a physically based, hydrological model. The responses of the model to 26 pre-post land use and fictional escalating rainfall event scenarios were examined for each of the seven catchments. The model behaved well in simulating physical processes of Hortonian overland flowsurface runoff starting between rainfall intensity of 2 mm/ hr-1 and 8 mm/ hr-1, as well as predicting a critical threshold for acceleration of modelled runoff occurring between rainfall intensity of 8 mm/ hr-1 and 26 mm/ hr-1, that could generate potential flooding. Results revealed that total and annual surface runoff would increase by 51% on average as a result of the wildfires. The modelling approach demonstrated here may provide a way of prioritizing catchment selection with respect to post-fire remediation activities. We believe that this modelling assessment methodology would be valuable to other semi-arid areas, while it provides an important means for comprehensively assessing post-fire response over large regions and therefore attempts to bridge some of the gaps in the specific literature field of research.

KEYWORDS: runoff response, hydrological model, MIKE SHE, semi-arid, shrubland, wildfire

INTRODUCTION

Fire is an integral part of many terrestrial biomes including the Mediterranean ones, but is also a major factor of disturbance (Pausas et al., 2008) that affects the interactions between physical, chemical and biological, biotic and abiotic components of terrestrial ecosystems or catchments (Rogers, 1996).

Over the last 20 years, a growing body of literature has been dedicated to wildfire impacts worldwide, notably Australia and South Africa, and particularly in the Mediterranean Basin (Shakesby and Doerr, 2006). During the last decades an increase in the number and size of fires has been observed in Mediterranean areas (Moreno et al., 1998; Pi?ol et al., 1998), which was attributed to land abandonment and afforestation of former agricultural land leading to increased fuel accumulation (Moreira et al., 2001; P?rez et al., 2003), as well as to the influence of climatic changes (Pausas et al., 2008; Pausas, 2004; Pi?ol et al., 1998). Particularly, in response to global warming over the next few decades, fire extent is expected to increase (McCarthy et al., 2001). As a result of above parameters the Mediterranean landscape is more susceptible to fuel accumulation and this may lead to increased fire occurrence (Loepfe et al., 2010; Vega-Garc?a and Chuvieco, 2006, Moreira et al. 2011). Particularly during the 1990s, approximately 50,000 wildfires affected annually 600,000 ha of forest and other wooded land (McCarthy et al., 2001), while burning the vegetation increases the presence of shrubs like Cistus (Margaris, 1980) which are invasive and highly flammable.

One of the most profound effects of wildfires is related to the alterations induced in the watershed processes, altering the normal flow patterns and disrupting the expected hydrologic behaviour several years after fire (Kinoshita and Hogue, 2011). Wildfires change the infiltration properties of soils and reduce the amount of interception materials (canopy, litter, duff and organic debris), that protect soil from raindrop impact and slow runoff (Moody and Martin, 2001a; Ice et al., 2004). The combustion process converts the litter and duff layers into ash and charcoal, while ash and small soil particles can seal soil pores (Morin and Benyamini, 1977; Neary et al., 1999). As a result infiltration generally decreases (Arend, 1941; Anderson, 1949; Fuller et al., 1955), also because chemical (DeBano et al., 1977; Giovannini et al., 1988) and physical (Doehring, 1968; Giovannini and Lucchesi, 1983) properties of soil may have changed. Soil becomes more water-repellent and, thus, affects the hydrologic response (DeBano et al., 1967; Chandler et al., 1983; Moody and Martin, 2001b). Since surface roughness is significantly reduced the subsequent hydrologic response to the normal precipitation regime is often a sudden and dramatic increase in overland flow (or surface runoff) velocities, the size of peak flows and ultimately water discharge at the catchments outlets (e.g. Hoyt and Troxell, 1934; Brown, 1972).

As a geomorphic agent, except from increasing runoff (Rowe et al., 1954; Krammes, 1960; Florsheim et al., 1991), wildfires can lower the intrinsic threshold of erosion in a catchment (Schumm, 1973; Wilcox et al., 1997; Davenport et al., 1998), with large increases in the potential for high sedimentation rates, debris flows, and other channel impacts (e.g. Meyer and Wells, 1997; Wilson, 1999; Cannon and Reneau, 2000; Cooper, 1961), and also cause rock weathering (Blackwelder, 1927). In addition, rainfall in semi-arid and mountainous systems with rocky and thin soils is often characterized by convective storms that produce high intensity episodic rainfall which in turn may trigger large floods (Costa, 1987). The Mediterranean catchments are therefore prone to damaging flash floods (Neary et al., 2005). Such events have become more severe during the last decades, especially under the rapid expansion of Wildland-Urban Interface (WUI; the area where structures and other human development meet with undeveloped wildland) observed in Mediterranean countries (Lein and Stump, 2009).

Hydrological responses to wildfire depend upon the precipitation climatology and may be viewed at different temporal scales. Moody and Martin (2001a) maintained that changes in annual runoff or changes in peak discharge can be used to quantify the hydrologic response after wildfire. Effects on the frequency and magnitude of peak discharge events and flooding phenomena resulting from deforestation constitute a significant issue that was not dealt within this study, as these are widely documented in the related literature. As Shakesby and Doerr (2006) state at catchment scales, post-fire accentuated peakflow has received more attention than changes in total flow, reflecting easier measurement and the greater hazard posed by the former.

Also hydrological responses to fire at the catchment scale have received less attention than at smaller spatial scales largely because of the greater difficulties of installing and maintaining instruments (Brunsden and Thornes, 1979), the greater spatial heterogeneity of environmental factors such as geology, topography, vegetation and soils, as well as fire extent over large catchments, fire severity and fire history (Miller et al., 2003). In their review on wildfires as hydrological and geomorphic agent, Shakesby and Doerr (2006) identified research gaps in the related literature that include the need to increase the research effort into past and potential future hydrological and geomorphic changes resulting from wildfire and maintained that past research features mainly post-fire conservation issues or deal more with prescribed fire rather than the impacts of wildfire. We believe that this effort attempts to bridge some of these gaps in the specific literature field of research. Furthermore, Pierson et al. (2001) concluded that the hydrological consequences of fire have been widely examined in forest ecosystems, but few studies have examined wildfire impacts on rangeland hydrology.

Lavabre et al. (1993) mention that even if the general consequences of fire have been largely studied, in particular the ecological and biological impact and the associated soil alteration, there are relatively few studies evaluating the effect on the hydrological response. This relative scarcity of quantitative hydrological studies merely reflects the difficulty of obtaining good quality data to compare hydrological behaviour before and after a fire. Scientific interest in a burnt area generally starts only after the fire, and therefore it is very difficult to find previous studies and pre-fire data to compare (Lavabre et al., 1993). Moreover, while there is an extensive literature on experimental studies of fire impacts, the field of impact modelling is less developed (Bovolo et al., 2009). Finally, any soil aspects, such as erosion and sedimentation rates analyses were not carried out in this study as Shakesby and Doerr (2006) indicated that there is more emphasis in the literature on the impact of wildfire on soil erosion by water than on the timing and quantity of runoff.

Most studies have identified the first one or two post-fire years as the critical period for high runoff and erosion risk, with recovery to pre-fire rates generally within five years (Wright and Bailey 1982). Veenhuis (2002) stated that the hydrologic effects of a wildfire seem to be more pronounced for the three years following the date of the fire. However, in dry areas, such as Mediterranean shrublands much higher runoff and erosion rates are being noticed even five to ten years after the fire (Inbar et al., 1998; Robichaud, 2000; Mayor et al., 2007). In this study, the hydrological effects of wildfires on catchments surface discharge are assessed for the first three years following the fires and the vegetation recovery within or extending over this period are not assessed herein.

By capitalizing valuable pre-fire data this paper aimed at analysing and predicting the first three years hydrological response (total, annual and mean monthly runoff) of seven semi-arid, shrubland-dominated catchments to land use changes in the form of deforestation at catchment scale, as a consequence of wildfires. Specifically, the possibility of comparing the hydrological consequences on seven similar catchments affected by the wildfires, lends this study a special relevance. To this end, we considered a three year pre-fire period as the reference period defining the normal behaviour for the catchments. The post assessment of the influences of meteorological characteristics and wildfires was therefore based on this -before the wildfires or past- hydrological regime. A semi-distributed or quasi-lumped, physically-based hydrological model has been calibrated, and used to simulate what the response of catchments surface runoff and other water balance components would have been for the first three years after the wildfires have occurred. Differences between these simulations, scenarios tested and the response observed are analysed herein. Particularly the responses examined were: (i) before-after realistic and fictional land use scenarios, and (ii) before-after fictional escalating rainfall scenarios. Furthermore, a correlation and simple linear regression analysis between burnt land and the increase of surface runoff were carried out. Our analysis illustrates aspects of the dynamics and the interrelationships between deforestation, meteorological characteristics and catchments� surface water resources responses to wildfires.

METHODS

Study area characteristics

The study area concerns a 465 km2 semi-mountainous and lowland area, located east of the Athens periurban area, on the eastern side of Ymittos, Penteli and Parnitha mountains, at the prefecture of Attica, Central Greece. The area was initially divided to 10 catchments, namely: Kato Souli, Lower Marathon (subcatchment downstream of Marathon Lake-dam), Upper Marathon (subcatchment upstream of Marathon Lake-dam), Nea Makri, Grammatiko, Rafina, Vranas, Coastal Catchment 1, Coastal Catchment 2 and Coastal Catchment 3 (Figure 1).

The region is characterized by a Mediterranean semi-arid climate with mild, wet winters and hot, dry summers, with a mean annual temperature around 18 oC, relative humidity 62% and relative sunshine duration 66%. The mean annual number of rain days is 72 (20%) and the mean annual rainfall depth for the three years studied is 465 mm; the reference evapotranspiration rate is more than 1.5 times the rainfall depth while it varies from about 1 mm/day during winter to 3.5 mm/day during summer. Although the hydrographic network of the catchments is particularly dense including fifth-order streams according to the Strahler (1952) method, the low runoff rate in combination with the natural relief did not lead to the formation of significant river networks and river cross-sections. The catchments are mostly drained by streams and ephemeral torrents with little or no water at all for most of the year, dry stream beds, during the summer and quick surface runoff, during storm precipitation events (Diakakis, 2011).

The average elevation of the study area is 204 m, the highest altitude is 1330 m and the average slope is 15% (Figure 1). The western part of the study area is characterized by steep slopes and rocky areas, in the place of which large portions of forest existed in the past. The eastern part presents gentler slopes and is more densely populated. A small-size, but very important, freshwater body, Lake Marathon, divides the upstream from the downstream Marathon subcatchments. Groundwater pumping constitutes the only source for water demand in the study area (30 hm3) of which 94% contributes to irrigation and 5% to water supply, while a very small amount is distributed to stock-breeding (0.5 %) and local industrial users (0.2%).

The study area is of great concern and significance in terms of water resources and environmental management. This is partly because the dry arid conditions that prevail in the area result in higher risk for wildfires occurrence, as well as more frequent fires. At the same time human activities in the area, such as the expansion of residential areas and extension of the Athens periurban zone, have increased the risk of large and high severity fires. Further, shrublands and coniferous forests in the area have higher potential post-fire runoff and erosion rates than densely vegetated broadleaf forests. In addition, a number of human and natural resources in the area are nowadays at risk, especially with respect to water quality and water management issues, such as the Marathon reservoir, that is used as Athens supplementary water supply source.

The area experiences very often wildfires creating a special fire regime characterized by rather small fire return intervals (Pleniou et al., in press). Recently, the study area has been recently affected by two wildfires in August 1995 and in August 2009. Before the first wildfire a large part of the central area was mainly covered by a dense pine forest that was almost totally destroyed. Following the typical regeneration patterns after 14 years just before the second fire in August 2009 the study area presents a complex land use mosaic consisting of shrublands, pastures, cultivated areas and pine forests at their first stages of development (Soulis et al., 2010). In the 21st August 2009, a fire started near Gramatiko village at the north part of the study area and then moved southwards for 5 days burning approximately one third of the study area.

Pre and post-fire land use/land cover mapping

Pre-fire landscape composition was generated through a map-update procedure of a 1:20000 scale Land Use/Land Cover (LULC) map produced originally on 2001. This map was thematically and geometrically updated based on photointerpretation of large scale natural color orthoimagery with 50 cm spatial resolution acquired in 2007 over Greece.

The burned area extent was delineated applying a GEographical Object Based Image Analysis (GEOBIA) approach to a post-fire satellite image (Figure 1) acquired by the Advanced Land Imager (ALI) instrument on EO-1 satellite on 27th August 2009, 3 days after fire suppression. The ALI instrument provides Landsat type (spatially and spectrally) panchromatic and multispectral bands with three additional bands covering 0.433-0.453, 0.845-0.890, and 1.20-1.30 �m. The satellite data were freely available from United States Geological Survey (USDA). An ortho-rectification process was originally applied using a 10 m digital elevation model (DEM) and ground control points identified on the orthoimagery. The dark object subtraction approach based on equations and rescaling factors provided by Chander et al. (2009), was applied to correct the effects caused by the solar zenith angle, solar radiance and atmospheric scattering. One advantage of this image-based approach is that it is robust, not requiring any in-situ measurements and atmospheric parameters (Chavez Jr, 1996).

Two levels of segments were developed empirically (Mallinis and Koutsias, 2012) for burned area delineation and LULC revision following fire spread based on the multi-scale region-growing segmentation algorithm embedded within the Trimble eCognition software. In the first (lower) level very small segments were produced for classifying likely burned, impervious, and agricultural, vegetated and water segments. A fuzzy rule based classification relied on the use spectral information from the near-infrared and short-infrared part (ALI band6 and ALI band8) of the electromagnetic spectrum. In the upper, more abrupt level, burned patches were classified based on the existence of likely burned sub-objects, and spectral information from the green and near-infrared (ALI band2 and ALI band6) part of the electromagnetic spectrum. To assess the accuracy of the burned area process, true colour orthophotographs of the affected region acquired in 2010 with very high spatial resolution (0.25 meters), were used (www.okxe.gr[EMS3]).

Model calibration and validation

The modelling framework that was chosen in this study for the simulation of the hydrological model and the application of the different land use, meteorological and climatic scenarios was MIKE SHE. MIKE SHE is a derivative of the Systeme Hydrologique Europeen, SHE, developed by Abbott et al. (1986) and is a comprehensive, deterministic, physically based, spatially distributed hydrological model that has been widely used to study a variety of water resource and environmental problems under diverse climatological and hydrological regimes (Refsgaard, 1997). In this case study, a semi-distributed, quasi-lumped approach was followed to describe the processes of overland flowsurface runoff, evapotranspiration (parameters based on land uses distribution), unsaturated vertical flow (two-layer soil column) and saturated subsurface flow (linear reservoirs approach). Analytical solutions were used to describe interception (Rotter model), evapotranspiration (Kristensen and Jensen (1975) model) and snow melt (degree-day approach) (DHI, 2011).

The actual number of parameters in MIKE SHE depends on which modules are included in the model setup and how the catchment is discretized. The guiding principle in the parameterisation was to construct a simple model with as few parameters subject to calibration as possible (Refsgaard, 1997). The catchments model domain cell size was uniformly discretized into a grid of 100 m x 100 m cells, which allowed for the most accurate representation of hydrological variables without, at the same time, placing excessive demand on computational time.

Precipitation, potential evapotranspiration and pumping time series data were used to define temporal variability of the water inflows/outflows. A precipitation and a meteorological network comprising of six and seven gauging stations, respectively, has been utilised, while daily data of recorded precipitation, relative humidity, wind velocity, solar radiation, vapour pressure, and air temperature data have been acquired for the period covering the hydrological years 2006-2009. Precipitation and potential evapotranspiration timeseries were corrected for missing data using bilinear interpolation and correlation methods between adjacent gauging stations. Pumping timeseries were created based on water demands of the municipalities in the study area (data from Ministry of Development, 2006). Thiessen polygons for the corresponding gauging stations were delineated in a GIS environment (ESRI ArcGIS), which allowed the spatial distribution of precipitation and potential evapotranspiration to be considered in the model. Potential evapotranspiration timeseries were calculated using the modified Penman-Monteith equation (Maidment, 1993). These were subsequently used by Kristensen and Jensen (1975) method inside MIKE SHE model in order to estimate the actual evapotranspiration.

Information regarding the spatial variation of the geological and soil substrata within the study area was obtained from hard-copy 1:50000 maps produced from the Greek Institute of Geology and Mineral Exploration (IGME) and National Agricultural Research Foundation (Nakos, 1979). The spatio-temporal variation of each LULC class as concerns leaf area index (LAI), root depth (RD) and crop coefficient (Kc) monthly values were derived from the literature similar to previous studies (V?zquez and Feyen, 2003; V?zquez et al., 2008). Appropriate literature values of Manning number as well as other simple overland parameters (e.g. slope, detention storage) were assigned for the different overland flow zones, the spatial distribution of which was based on the results of the land cover mapping. In order to save computational time and reduce data requirements a simple two-layer unsaturated zone (UZ) model was used with the appropriate types of soils based on the area�s digitised soil map. The upper MIKE SHE UZ layer represents the soil root zone of the existing vegetation cover, while the second layer below extends down to the aquifer�s top level. Soil parameter values of water content at saturated conditions and at field capacity, field wilting point and infiltration rate were defined within the range of literature values and used as calibration parameters both for the experimental subcatchment water balance and for the observed hydrograph. The spatial distribution of UZ and saturated zone (SZ) was implemented using a digitised geological map of the area. A linear reservoir approach was subsequently followed for the SZ geological units that involved slower groundwater movement in the baseflow linear reservoir and a relatively faster groundwater flow in the interflow linear reservoirs. This latter faster groundwater flow component constitutes the quick response to rainfall inputs that subsequently feeds into the stream network or the baseflow reservoir.

Frequently, the hydrologic effects of wildfire cannot be studied using the paired-catchment or the calibrated-catchment methods because pre-fire data are usually not available except in a few rare cases, where catchments that were being monitored for other purposes were burned by a wildfire. Such is the case of Pikermi subcatchment where observed values of monthly streamflow (water level gauge at the subcatchment outlet located within the catchment �Rafina�, (Figure ???) 1) and meteorological characteristics were obtained from the National Technical University of Athens experimental subcatchment of �Pikermi� (NTUA, 2010), while a calibrated-catchment approach was used.

The model was calibrated for the hydrological period 2006�2009 with a manual trial-and-error procedure by altering the calibrated parameters, until the simulated surface runoff at the outlet of the experimental subcatchment matched closely the observed values. The range for each calibration parameter in the model was set using previous studies and/or physical reasoning. A standard performance metric, the Nash and Sutcliffe (1970) coefficient of efficiency was chosen as the likelihood measure to evaluate the accuracy of both the magnitude and timing of predicted flows (e.g. Andersen et al., 2001; Beven, 2001; V?squez et al., 2002) and the overall model performance, defined as:

(1)

in which O is the observed flow, O is the mean observed flow, and P is the predicted flow.

In order to minimize the number of parameters used in model calibration a simple comparative screening parameterisation strategy was followed. At the same time, multi-variable checks were carried out during the calibration to ensure total and annual water balance components� values within the literature range for the specific area. Soil hydraulic conductivity, time constants for interflow and baseflow reservoirs and specific yield were the most sensitive calibration parameters. Due to limited stream discharge data the model could not be temporally validated with a split sample test in order to evaluate the robustness of these behavioural sets under a range of hydrological conditions in the experimental subcatchment or any other catchment. It is expected that such data will become available in the next few years, so that the model gets updated, validated and its performance improved.

After model calibration, the group of behavioural parameter sets identified in the calibration period for the Penteli Pikermi experimental subcatchment i.e. hydrological, landuse, soil and hydrogeological parameters were transferred to seven ungauged catchments in the study area and various scenarios were tested in order to simulate and compare the catchments response to runoff before and after the wildfires. The rare opportunity was therefore presented to use the calibrated subcatchment parameters in order to quantify fire effects in the overlapping and neighbouring ungauged catchments before the fires and compare these with the post-fires results.

In total, 26 different scenarios were applied and the runoff response results were examined for each of the seven catchments. The initial model represented the base scenario, LU0, with land uses of the catchments before the wildfires. The following hydrological modelling scenarios were then introduced into the model:

* Scenario, LU1, with land uses of the catchments after the wildfires, compared against scenario LU0;

* Six scenarios with escalating rainfall events, R0i-R0vi for land uses of the catchments before the wildfires,

* Six scenarios with escalating rainfall events, R1i-R1vi for land uses of the catchments after the wildfires, compared against scenarios R0i-R0vi,

* Six scenarios with escalating rainfall events, RBi-RBvi, assuming the entire catchments area as burnt land,

* Six scenarios with escalating rainfall events, RFi-RFvi, assuming the entire catchments area as fully recovered with coniferous forests, compared against scenarios RBi-RBvi.

The six escalating rainfall events, from i to vi, that were used, accordingly, in scenarios R0i-R0vi, R1i-R1vi, RBi-RBvi and RFi-RFvi, were classified by adapting World Meteorological Organization (WMO, 2008), as follows:

* i, 0.3 mm/ hr-1 (drizzle);

* ii, 2 mm/ hr-1 (weak rain);

* iii, 8 mm/ hr-1 (medium intensity rain);

* iv, 26 mm/ hr-1 (strong rain);

* v, 40 mm/ hr-1 (intense rain);

* vi, 55 mm/ hr-1 (very heavy rain).

Daily time steps were considered for the hydrological models used in land use scenarios LU0 and LU1. Model predictions of catchments discharges were made at a daily time step and aggregated to monthly, annual and total values. Time steps of one second were considered for the hydrological models used in escalating rainfall events scenarios R0i-R0vi, R1i-R1vi, RBi-RBvi and RFi-RFvi. For the specific scenarios requirements, the model was simulated for one hour. For post-fire simulations we modified the spatial distribution of input parameters (LAI, Rd, Kc and Manning�s n).

RESULTS AND DISCUSSION

Fire scar mapping

The estimated surface of the area burned was 14919.72 ha with an overall accuracy of 92% and a Kappa index of agreement of about 0.84. The accuracy is very similar to the one reported in a recent study over the same area, using an EO-1 ALI image acquired 3 days later (i.e. 30th of August) and machine learning algorithms (Petropoulos et al., 2012) and within the observed range of accuracy values achieved from studies using similar satellite data (Mallinis and Koutsias, 2012).

Table 1 presents the surface extent for each catchment as well as the surface extent of burnt land within each catchment. It is evident that burnt areas represent more than 60% within the corresponding catchments in Vranas, Lower Marathon and Grammatiko catchments, which were most affected by the wildfires. Burnt areas represent between 24-25% within the corresponding catchments, of Rafina and Kato Souli catchments, while smaller percentages of burnt areas were observed at Nea Makri (18%) and Upper Marathon (14%) catchments. All three Coastal catchments are isolated small coastal catchments at the eastern part of the study area that were to a great extent unaffected by the wildfires (~0% burnt land) and are therefore not analysed herein.

Model calibration

Overland and land use values for shrublands were primarily modified, since according to Table 2 more than 91% of the experimental subcatchment area was classified as shrublands (33% of the study area catchments) while it was the main vegetation type affected by the fire (decrease of 20% of surface extent). Unsaturated and saturated zone values for low permeability/impermeable soils/geological units were also primarily modified, since more than 67% of the experimental subcatchment was classified as geological units of impermeable flysch, schist and marls (as well as 54% of the study area surface extent) and it was the main soil/geological type affected by the fire (88% of study area surface extent) (Table 2).

Calibration results are presented at Figure ??? 3 to Figure ???. 5. The results show satisfactory calibration of the model with a coefficient of efficiency (Nash-Sutcliffe) of 38% and statistical error criteria MAE, RRMSE and BIAS of 0.01, 0.52 and 0.005, respectively. As it is evident from Figure ??? 3 deviations of simulated discharge occur in March 2007, where the model overestimates the observed discharge by ~0.025 m3/s, and from December 2007 to April 2008 where the model underestimates the observed discharge by 0.016 m3/s. This decrease in performance may be attributed to uncertainties associated with precipitation and water level/discharge gauges, as well as other data inputs to the model. The model behaves quite well during the rest of the simulation and particularly during the last hydrological year 2008-2009, when the coefficient of efficiency reaches 68%. Discharge generally peaks between the months of December-April, while the minimum is observed during the summer months, June to September, when baseflow contributes entirely to stream discharge.

The scatterplot (Figure ???) 4) of the observed vs simulated mean monthly discharge shows that simulated discharge is generally underestimated in the low range of observed discharge, particularly between 0.01 and 0.04 m3/s. Prediction errors (i.e. observations falling outside the uncertainty bounds) are indicative of deficiencies in model structure, uncertainty in input data, and/or errors in observed streamflow.

The fluctuation of observed accumulated vs model accumulated mean monthly discharge shows are presented at Figure ??? 5 where it is evident that the model initially overestimates the discharge compared to the observed values. Between December 2006 and December 2007 there is generally an agreement between observed and predicted discharge, then for the next hydrological year 2007-2008 a deviation of simulated discharge occurs compared to the observed values, and finally, during the last hydrological year 2008-2009, there is again an agreement between the observed accumulated and model accumulated values.

Figure ??? 6 presents the differenced (or deviation of) simulated (Post-Pre) mean monthly surface runoff in each of the seven catchments. This parameter was calculated when subtracting mean monthly surface runoff on land uses before wildfires from the simulated surface runoff on land uses after the wildfires. As it is evident from Figure ???, 6, monthly surface runoff is higher after introducing into the model the changes in land cover due to fire, since the deviations are positive. Also the differences in mean monthly surface runoff are higher and more obvious with higher precipitation occurring in winter and spring, while the differences are smaller at lower surface runoff values and with lower precipitation, particularly during the summer months. The absolute and proportional increase in surface runoff is more noticeable at the catchments of Lower Marathon, Grammatiko, Vranas and Nea Makri, whereas it is found to be lower for Upper Marathon, Rafina and Kato Souli catchments (Figure ??? 6 and Table ???).3).

Correlation, regression and hydrological Model model validationresults[EMS4]

Table ??? 3 presents the average increase of surface runoff expressed as absolute volume (hm3) and as percentage increase; and burnt area as surface extent (km2) and as percentage of catchment area. It is evident that as the percentage coverage of burnt areas increases, so does the surface runoff. Another observation is that there are two different groups of catchments with burnt land percentage ranges: the first group reflects to slightly affected catchments Rafina, Upper Marathon, Kato Souli, Nea Makri with low proportion of burnt land between 13-25 % and the second group concerns heavily affected catchments Lower Marathon, Vranas, Grammatiko with high proportion of burnt land between 60-63%.

However, average increase of surface runoff, as a response to wildfires, is better correlated with burnt area expressed in absolute values (km2) rather than proportion (%) of burnt land extent within a catchment (Figure ???7b and ???7a). In the former case, the relationship is very strong with a correlation coefficient of 96.1 %, whereas in the latter case the relationship is less strong (47.9%, i.e. half of what it was with absolute value of burnt land extent).

The latter correlation coefficient value improves to the significantly high levels (>95%) when two sub-groups of catchments are differentiated (Figure ???7a). Particularly, Group I with catchments Lower Marathon, Rafina and Upper Marathon presents a fitted trend line with correlation coefficient of 97.2%, while Group II with catchments Vranas, Grammatiko, Kato Souli and Nea Makri a correlation coefficient of 99%. The strong correlation for the specific groups of catchments may give indication that surface runoff behaves differently, in a way that may be attributed to catchment size. Particularly the catchments of Group II are comparatively small-sized (<50 km2), whereas the catchments of Group I are medium-sized (Group I: >70 km2). Such simple linear regression relationships may be used in semi-arid areas that present similar land use, soil and geological characteristics, i.e. with significant areas of shrublands and large extent of impermeable soils and geological units (flysch, schists and marls).

Figure ??? 8 presents the results of simple linear regression and correlation analysis between total surface runoff (hm3) before vs after the wildfires for the simulation period 2006-2009. The graph on the left shows the linear and multinomial relationships between total surface runoff before and after the wildfires and includes separately the two Upper and Lower Marathon subcatchments. As it is shown the second order multinomial relationship presents a very strong correlation of 98.6%, while the linear relationship also presents a strong correlation of 88%. The linear relationship improves even more, with a correlation coefficient of 99.4%, when the values of surface runoff for the two Marathon subcatchments are added together (right graph below). The strong linear relationship is given by:

Total annual surface runoff after (hm3) = 1.32 x Total annual surface runoff before + 2.58

The results revealed that total and annual surface runoff would increase by 51% on average as a result of the wildfires. This response of simulated annual runoff to the wildfires is comparable with various similar cases of monitored sites where annual runoff after the wildfires would increase significantly from 0.09- to 21-fold (Hoyt and Troxell, 1934; Helvey, 1980; Langford, 1976; Campbell et al., 1977; Troendle and Bevenger, 1996; Mackay and Cornish, 1982; Prosser, 1990; Moody and Martin, 2001a; Soler et al., 1994; Scott, 1997; Lavabre et al., 1993; Scott and Van Wyk, 1990 Badia and Marti, 2008).

The annual surface runoff coefficients results presented comparatively high runoff coefficients for the catchments Lower Marathon, Grammatiko and Vranas (between 19% and 43% on average), while lower runoff coefficients are observed at catchments Kato Souli, Upper Marathon, Nea Makri and Rafina (between 0.07 and 0.20 on average). This finding is consistent with the occurrence of higher runoff coefficients at catchments where a higher percentage of burnt areas is observed. Figure ??? 9 presents the deviation of (or differenced) simulated annual surface runoff coefficients at seven catchments for the period 2006-2009, between the actual land uses (Pre wildfires) and the land uses that resulted after Attiki summer 2009 wildfires (Post wildfires). Higher increase of surface runoff coefficients is also observed at Grammatiko, Lower Marathon and Vranas catchments, which is attributed to the higher percentage of burnt land within each catchment. On average an increase of 1.5-fold of annual surface runoff coefficients occurred across the seven catchments from 16.1% to 25.2%.

The increase of annual surface runoff coefficient is higher for the hydrological year 2007-2008 due to lower annual precipitation values for the specific year, rather than the increase of absolute surface runoff values. This finding is in accordance with a number of authors who have reported higher overland flowsurface runoff coefficients following dry than wet periods, which were mainly attributed to the enhancement or development of soil water repellency (e.g. Shahlaee et al., 1991; Walsh et al., 1994), reflecting Hortonian rather than saturation overland flowsurface runoff (Soto and D?az-Fierros, 1998) and supporting the validity of our modelling approach.

The results of differenced total surface runoff i.e. the increase of surface runoff as a response to escalating rainfall scenarios for land uses before the wildfires (R0i-R0vi) compared to the results of the respective scenarios for land uses after the wildfires (R1i-R1vi) are presented in Figure ???. 10 Differenced (or increase of) total surface runoff is larger for catchments Lower Marathon, Rafina and Vranas, which are the most affected by the fire in terms of absolute extent of burnt land. The differenced total surface runoff is smaller for Nea Makri, Kato Souli, Upper Marathon and Grammatiko catchments.

The results of differenced total surface runoff i.e. the increase of surface runoff as a response to escalating rainfall scenarios for uniform land use assuming the entire catchments as burnt area (RBi-RBvi) compared to the results of the respective scenarios for uniform land use, assuming the entire catchments as area fully recovered with coniferous forests (RFi-RFvi) are presented at Figure ???. 11 Differenced (or increase of) total surface runoff is higher for catchments Upper Marathon, Rafina and Lower Marathon. On the contrary, differenced total surface runoff is lower for smaller catchments of Nea Makri, Grammatiko, Vranas and Kato Souli. This finding may be attributed to the fact that the former catchments are the largest in surface extent, which in turn may indicate that total surface runoff for the scenarios transition from uniform forest land use to uniform burnt land is influenced to a large extent by the geometrical parameters of catchments size and perimeter.

As it is evident from Figures ??? 10 and ??? 11 total surface runoff for the first and second rainfall events of drizzle and weak rain (0.3 mm/hr hr-1 and 2 mm/hr hr-1, respectively) result in zero or close to zero values for all catchments and scenarios, as a result of the soil�s detention storage trapping all the raindrops in small depressions, as well as of the evapotranspiration and infiltration processes limiting overland flowsurface runoff. According to the hydrological model output, Hortonian overland flowsurface runoff resulting into streamflow and finally discharge at the catchments outlets starts between 2 mm/hr hr-1 and 8 mm/hr hr-1.

Another observation is that the larger differences of increase of surface runoff for the two pairs of tested scenarios occur between the third and the fourth rainfall events (8 mm/hr hr-1 and 26 mm/hr hr-1), compared to the results of surface runoff increase for the other three rainfall events (26 mm/hr hr-1, 40 mm/hr hr-1 and 55 mm/hr hr-1). This indicates some acceleration of runoff occurring between rainfall of 8 mm/hr hr-1 and 26 mm/hr hr-1.

These observations are in accordance with other authors reporting the existence of a critical threshold intensity. Data from wildfires in Colorado, New Mexico, and western South Dakota indicate that storms with a maximum 30-minute (I30) rainfall intensity of only 7-10 mm h-1 can induce Horton overland flowsurface runoff (Cannon et al. 2001; Moody and Martin 2001b; Benavides-Solorio 2003; Kunze and Stednick 2006; Wagenbrenner et al. 2006). MacDonald and Larsen (2009) reported that rainfall intensities of only 8-10 mm h-1 can generate substantial amounts of runoff and surface erosion. Kinner and Moody (2008) suggested that in their experimental case using multiple rainfall intensities a rainfall intensity around 20 mm/ hr-1 might have been the threshold at the 1-square-meter scale. Moody and Martin (2001a) reported a threshold of rainfall intensity of 10 mm h-1 below which the majority of the rainfall infiltrates and above which the rainfall intensity may exceed the average catchment infiltration rate, such that the runoff is dominated by sheetflow, which produces flash floods. In their study this threshold changed with time (from 10 mm h-1 the first one-two years after the wildfire to 16 mm h-1 the third year) reflecting the recovery of the hillslopes as infiltration and the canopy density of the fire-adapted, under-story vegetation increase. However, at large spatial scales (> 100 km2), runoff would be less variable as the variability of rainfall intensity would be averaged over a larger surface area (Moody et al., 2008). Although flash floods and stormflows were not studied in this paper, an indication of flashy behaviour occurring in the model is provided by the finding of such a threshold existing between 8 mm/ hr-1 and 26 mm/ hr-1 of rainfall. Nevertheless, this does not necessarily imply that our model could potentially be used as a forecasting tool e.g. to set threshold limits in rain gauges that are part of an early-warning flood system after wildfires.

CONCLUSIONS[����5]

Hydrological responses to fire at the catchment scale have received less attention than at smaller scales (Brunsden and Thornes, 1979), as well as examining wildfires effects on multiple catchments simultaneously rather than calibrated catchment or paired catchment approaches. The possibility of comparing the hydrological consequences on seven similar catchments affected by the wildfires, lends this study a special relevance. We believe that this effort attempts to bridge some gaps in the literature field of wildfires as hydrological agent and their impacts on semi-arid, shrubland-dominated hydrology.

The probability of landscape-scale disturbances such as fire is expected to increase in the future due to anticipated climate changes and past land management practices. Particularly, the frequency and intensity of disturbance events such as fires are expected to increase in future (Lenihan et al., 1998; Swetnam et al., 1999; Easterling et al., 2000). These disturbances can produce dramatic changes in hydrologic responses (e.g. surface runoff) that can pose risks to human life, infrastructure, and the environment (Beeson et al., 2001). Hence, it will be increasingly important to be able to assess rapidly and effectively pre-post wildfires changes and relationships in hydrology and to apply these assessments to evaluate such risks.

ACKNOWLEDGEMENTS

The ALI satellite data used in the study were freely obtained from the U.S. Geological Survey (USGS). We would also like to thank DHI for the kind provision of the hydrological software MIKE SHE for use within the scope of this research effort.