Assessment of Variations in Runoff Due to Landcover Changes Using the SWAT Model in an Urban River in Dublin, Ireland

Abstract

Investigating the impact of land cover change in hydrological modelling is essential for water resources management. This paper investigates the importance of landcover change in the development of a physically-based hydrological model called SWAT. The study area considered is the Dodder River basin located in southern Dublin, Ireland. Runoff at the basin outlet was simulated using SWAT for 1993–2019 using five landcover maps obtained for 1990, 2000, 2006, 2012 and 2018. Results indicate that, in general, the SWAT model-simulated runoff for a chosen time-period are closer to the real-world observations when the landcover data used for simulation was collated as close to the time-period for which the simulations were performed. For 23 (20) years (from 27 years period) the monthly mean (maximum) runoff for the Dodder River generated by the SWAT model had the least error when the nearby landcover data were used. This study indicates the necessity of considering dynamic and time-varying landcover data during the development of hydrological modelling for runoff simulation. Furthermore, two composite quantile functions were generated by using a kappa distribution for monthly mean runoff and GEV distribution for monthly maximum runoff, based on model simulations obtained using different landcover data corresponding to different time-period. Modelling landcover change patterns and development of projected landcover in the future for river basins in Ireland needs to be integrated with SWAT to simulate future runoff.

1. Introduction

Major alterations in land cover indicate the significance of human and anthropogenic activities on land use, which have been deemed responsible for extreme floods and the associated devastations [1,2]. Climate change and urbanization both lead to the generation of excess surface runoff which results in the more frequent occurrence of flood events [3]. Fragmented land cover areas caused by urban development have resulted in alterations to the hydrological cycle [4]. A recent report [5] stated that between one-half to one-third of our planet’s surface has been transformed by human development, leading to the rapid growth of urbanization. There is an increasing concern that urbanization, resulting in increased impervious areas responsible for the reduction in infiltration capacity of the soil and an increase in soil erosion, are causing increased flooding [6,7].

Previous studies reported that high rates of soil erosion are responsible for higher and quicker surface runoff [8,9,10]. The rate of surface runoff is dependent on various factors such as topography, morphology, river network [11,12] and the rainfall-runoff relationship [13]. Karamage et al. [14] reported that rapid urban expansion and agricultural intensification has resulted in an increase in the volume of rainfall-runoff. A geographical information system (GIS) based hydrological model was used in this study to understand the impact of soil texture, land cover and land use, elevation, and rainfall intensity on rainfall-runoff. Other studies [15] have reported that the collection of accurate land cover information is essential for the simulation of runoff at river basin outlets. Detailed information on the nature of land cover obtained in high spatial resolution increases the accuracy of runoff simulations [16] considerably. Another study confirmed that the reduction of vegetation cover has resulted in an increase in annual discharge [17,18]. In many areas the expansion of agricultural activities has resulted in the conversion of native forests and grasslands into agricultural land, leading to faster runoff and hence flooding [19,20]. Yang et al. [21] investigated the impact of combined climate and landcover change at upper Heihe river in China using SWAT model and noted that climate has slightly bigger impact on the runoff. Yang et al. [22] developed another coupled landcover-climate change model at the Yingluoxia and Minxian catchments in northeastern Tibetan Plateau.

Mousavi et al. [23] observed that in the Marand river basin, more than 25.3% of the total land cover is prone to flood risk. The constructed land such as the industrial, urban, and rural vegetated areas are particularly vulnerable to flooding in the river basin [24].

The Corine Land Cover (CLC) database is used to identify the dynamics and rate of change of landcover by many studies in Europe [25,26]. The CLC database aims to characterise landcover by spatial continuity which is non-ambiguous [27]. Previous studies have shown that the CLC data is used for analysing the changes in land use and land cover [28]. For example, Kucsicsa and Dumitrică [29] used the CLC database to estimate the changes in forestry and the impact of deforestation on population dynamics in Europe.

The growth rate of urbanization has been quite significant in Ireland [30] over the past 30 years. Ireland’s annual growth rate of urbanization was found to be 3.1% between 1990 to 2012 [31]. Several studies [31,32,33] have further documented that Ireland experienced an overall higher degree of land conversion relative to other European countries, concluding that the urban land expansion in Ireland has been among the highest in Europe.

The majority of the physically-based hydrological models that are used to simulate the discharge dynamics from river basin outlets require land cover information [34]. One of the major limitations of those models is that they, in general, assume that the land cover information remains the same over time. However, in reality, changes in land cover are evolving over time, which needs to be considered while simulating runoff at a river basin. In situations where the hydrological model is used to simulate runoff for a short period of time, ranging from a few months to a few years, a static land cover information might be sufficient, however, when runoff simulations are required for longer periods of time (more than a decade), it is important to consider the changes in land cover over such a period. This study investigates how changes in land cover impact hydrological runoff simulations using a physically-based rainfall-runoff model called soil water assessment tool (SWAT). The model was run for a 27-year period at a river basin covering a major portion of south Dublin, the capital of Ireland, by considering different land cover data obtained from the CORINE database. A set of SWAT models were developed by considering different land cover information and the corresponding hydrographs were simulated. The performance of the SWAT models was evaluated by comparing the simulated runoff with the observed runoff. In particular, the study focuses on the efficiency of the simulations with respect to when the land cover data used in the simulations had originally been defined. This study, therefore, has investigated the performance of SWAT models developed by considering different land cover data obtained at different time periods.

2. Methodology

This study considered a conceptual rainfall-runoff model to simulate flow at the outlet of the watershed. The model chosen is the soil water assessment tool (SWAT). The model was run corresponding to different land cover data to simulate runoff using the historically observed meteorological variables. The monthly mean and maximum runoff were extracted from the daily simulated runoff and was compared to historically observed streamflow data for performance evaluation. Subsequently, the monthly mean and maximum runoff were modelled using appropriate frequency distributions to estimate quantiles corresponding to different return periods. Details on data, SWAT modelling, performance measures, parameter estimation of selected distributions and quantile estimations are provided below.

2.1. River Basin and Data Description

Analysis of the SWAT model using different land cover data was performed in the Dodder River Basin for this study. The Dodder River Basin (Figure 1) is located in the southern part of Dublin, the capital of Ireland. The Dodder River is one of the most important rivers in the Dublin area and originating on Kippure mountain to the south and meeting the River Liffey in Dublin city at Ringsend. The length of the Dodder River is approximately 27 km and the overall catchment area is around 113 square km. In general, the flow in the river is low, however, due to high slopes in the catchment area, the upper and middle section of the river is highly susceptible to flooding during periods of extreme rainfall events. The catchment can be broadly classified as lower and upper reaches, where the lower reach belongs to the Dublin city with a milder slope, and the upper reach begins from the boundary of the Dublin city and reaches the Kippure mountains where the Dodder River originated. The average slope in the lower and upper reaches are approximately 2.27° and 9.26° respectively.

Discharge at the Dodder River has been measured at the Waldron Bridge gauging station from 3rd January 1986 at a daily scale, and are available from the Irish Environmental Protection Agency HydroNet database (https://www.epa.ie/hydronet/#09010 accessed on 26 November 2021). However, since the majority of the streamflow data are missing till 1992, data ranging from 1 January 1993 to 31 December 2019 were used for further analysis in this study. The average discharge in the river over this monitoring period was 2.07 m³/s. However, during flood events, the maximum flow in the river was noted to reach up to 70 m³/s. As the size of the river basin is small with steep slopes, it usually takes 2–3 h to reach the water from upstream at the mountains to the downstream of river Dodder, which can create flash flooding the area. Time series plots of the daily, monthly mean and monthly maximum data are shown in Figure 2.

Three watershed related types of data (elevation, land cover and soil type) and five meteorological variables data (rainfall, maximum and minimum temperature, wind speed, relative humidity, and solar radiation) are needed for the river basin for running the SWAT model. A digital elevation model (DEM) was used to delineate the watershed/basin boundary, develop the stream network in the watershed and estimate the slope of the watershed. The DEM data used in this study is the 25 m resolution EU-DEM data available from “https://www.eea.europa.eu/data-and-maps/data/copernicus-land-monitoring-service-eu-dem accessed on 26 November 2021”. The EU-DEM was generated for 39 countries across Europe by DHI GRAS Geocenter Denmark [35] as a hybrid DEM based on 90 m resolution shuttle radar topographic mission DEM, 30 m resolution advanced spaceborne thermal emission and reflection radiometer DEM and publicly available topographic maps. The DEM map for the Dodder River Basin is shown in Figure 3a.

The land cover map for the Dodder River Basin was obtained from the COPERNICUS CORINE land cover inventory (https://land.copernicus.eu/pan-european/corine-land-cover accessed on 26 November 2021). The land cover maps were prepared across European countries using visual interpretation of high-resolution satellite imageries. The land cover maps are available for the years 1990, 2000, 2006, 2012 and 2018 at a spatial resolution of 100 m. Landsat-5 satellite data obtained from 1986–1998 were used to create a land cover map for 1990 (LC1990); Landsat-7 satellite data collected for 1999–2001 were used to develop the 2000 land cover maps (LC2000); land cover maps for the year 2006 (LC2006) was generated using SPOT-4, SPOT-5 and IRS P6 satellite images collected for the period 2005–2007; IRS P6 and RapidEye satellite data obtained from 2011–2012 were used to create the 2012 land cover maps (LC2012); whereas Sentinel-2 and Landsat-8 satellite data collected for 2017–2018 were used to generate the 2018 land cover images (LC2018). The CORINE land cover maps consist of 44 different land cover classes, shown in

Table 1. Linkage of land cover classes available in CORINE data and SWAT model. Details of the SWAT classes are available in SWAT Model Database Appendix A (https://swat.tamu.edu/media/69419/Appendix-A.pdf accessed on 26 November 2021).

Sl.	Land Cover CODE	Land Cover Description	SWAT Class	Sl.	Land Cover CODE	Land Cover Description	SWAT Class
1	111	Continuous urban fabric	URBN	23	311	Broad-leaved forest	FRSD
2	112	Discontinuous urban fabric	URML	24	312	Coniferous forest	FRSE
3	121	Industrial or commercial units	UCOM	25	313	Mixed forest	FRST
4	122	Road and rail networks and associated land	UTRN	26	321	Natural grasslands	RNGE
5	123	Port areas	UTRN	27	322	Moors and heathland	RNGB
6	124	Airports	UTRN	28	323	Sclerophyllous vegetation	RNGB
7	131	Mineral extraction sites	SWRN	29	324	Transitional woodland-shrub	RNGB
8	132	Dump sites	UIDU	30	331	Beaches, dunes, sands	WETN
9	133	Construction sites	UIDU	31	332	Bare rocks	BARR
10	141	Green urban areas	RNGE	32	333	Sparsely vegetated areas	BARR
11	142	Sport and leisure facilities	UCOM	33	334	Burnt areas	BARR
12	211	Non-irrigated arable land	AGRL	34	335	Glaciers and perpetual snow	WATR
13	212	Permanently irrigated land	AGRL	35	411	Inland marshes	WETN
14	213	Rice fields	RICE	36	412	Peat bogs	WETL
15	221	Vineyards	GRAP	37	421	Salt marshes	WETL
16	222	Fruit trees and berry plantations	AGRL	38	422	Salines	WETL
17	223	Olive groves	OLIV	39	423	Intertidal flats	WETL
18	231	Pastures	PAST	40	511	Water courses	WATR
19	241	Annual crops associated with permanent crops	AGRL	41	512	Water bodies	WATR
20	242	Complex cultivation patterns	AGRL	42	521	Coastal lagoons	WATR
21	243	Land principally occupied by agriculture, with significant areas of natural vegetation	AGRL	43	522	Estuaries	WATR
22	244	Agro-forestry areas	FRST	44	523	Sea and ocean	WATR

. Since the SWAT model was developed mainly for simulation of runoff for United States river basins, the land cover classes present in the SWAT database do not match with the land cover classes in the CORINE maps. Hence, each of the CORINE land cover classes was reclassified to match them to the SWAT database. Details of the reclassified land cover classes are shown in Table 1 and the reclassified land cover maps for the Dodder River basin corresponding to years 1990, 2000, 2006, 2012 and 2018 are shown in Figure 4. It can be noted from the figure that a majority portion of the upper Dodder River catchment had been transformed from wetlands (peat bogs) to a combination of natural grasslands, moors and heathlands, and burnt/barren land from 2006 to 2018. One of the reasons is the rise in temperature during the summer months, scarcity of water, and increase in the forest fire in Ireland.

The soil maps for Ireland were obtained from the General Soil Map of Ireland, available at http://gis.teagasc.ie/soils/downloads.php (accessed on 26 November 2021). Similar to the land cover classes, the soil classes also needed to be reclassified for the SWAT model. The reclassified soil map is shown in Figure 5 and details on the soil categories is provided in Basu [36]. It should be noted that the land cover and soil maps for Ireland needed to be integrated with the SWAT model’s default database. This was performed using a lookup table for land cover and soil. Furthermore, the soil database is needed to be integrated into the SWAT geodatabase. The newly developed geodatabase (SWAT2012.mdb), the lookup table for land cover (lulc_classes.txt) and the lookup table for soil (soillookup.txt) for Ireland is provided in the following link https://doi.org/10.5281/zenodo.4767926 (accessed on 26 November 2021) [37] as open-source data for public use. Details on soil type in the figure can be obtained from the Teagasc Soil Guide (http://gis.teagasc.ie/soils/soilguide.php accessed on 26 November 2021).

The five meteorological types of data needed to run the SWAT model were obtained either from on-site weather stations data or from satellite-based data. In general, the performance of the SWAT model has been shown to be better in situations where multiple point-source weather stations meteorological data are used for analysis [34]. Figure 1 shows the locations of the Met Eireann rainfall gauges, temperature sensors, and weather stations that measure relative humidity (RH), wind speed (WS) and solar radiation (SR) in the surrounding vicinity of the Dodder River Basin, while

Table 2. Location of Met Eireann rain (R) gauges, temperature (T) gauges and weather stations measuring relative humidity (RH), wind speed (WS) and solar radiation (SR) in and around Dodder River Basin. Stations located inside the Dodder River Basin are highlighted in bold.

Gauge Number	Location Details		Data Range	Available Variables
	Lat	Lon
S175	53.364	−6.350	16/08/2003–31/12/2020	R, T, RH, WS, SR
S1823	53.370	−6.270	01/01/1984–31/10/2020	R, T
S1923	53.239	−6.367	01/01/1984–30/09/2014	R
S3723	53.306	−6.439	01/01/1964–31/10/2020	R, T, RH, WS, SR
S3923	53.341	−6.253	01/01/1948–31/10/2020	R, T
S5623	53.239	−6.364	01/10/1959–31/10/2020	R
S6623	53.276	−6.304	01/02/1967–31/10/2020	R
S7523	53.325	−6.225	01/08/1972–28/02/2015	R

provide details of those stations. Only three rain gauges (S1923, S5623, S6623) are located inside the Dodder River Basin, and no Met Eireann stations are present inside the basin measuring temperature, relative humidity, wind speed or solar radiation. Furthermore, data from station S1923 were unavailable from October 2014.

Due to the scarcity of on-site weather data for all meteorological variables except rainfall, only rainfall data from observed rain gauges were considered in the development of the SWAT model. On-site daily rainfall data from three rain gauges S3723, S5623 and S6623 were used for the analysis. Those three gauges were considered as there are very little rainfall data missing and they are located either inside the Dodder River Basin or at close proximity to the basin. Based on the analysis performed by Met Eireann, 60 mm increase in the annual average rainfall were noted from the period 1981–2010 [38]. Based on analysis performed using several regional climate models, Nolan et al. [38] documented that for the period 2041–2060, annual, spring and summer rainfall in Ireland are expected to decrease by 0–13%, while heavy rainfall events are expected to increase in the winter and autumn months by ~20%.

The other four meteorological variables temperature, relative humidity, wind speed and solar radiation required for the development of the SWAT model were obtained at a gridded scale from COPERNICUS E-OBS and ERA-5 land databases.

A readily available data source of the meteorological variables is the COPERNICUS E-OBS dataset. The E-OBS dataset was generated based on a network of station data observations from European National Meteorological and Hydrological Services and other European database centres. The dataset covers the majority of Europe and provides rainfall, maximum and minimum temperature, relative humidity, and shortwave solar radiation at 0.1° grids spacing at a daily time scale from the year 1950 till the present time. The data has been obtained from the following link: https://cds.climate.copernicus.eu/cdsapp#!/dataset/insitu-gridded-observations-europe?tab=form (accessed on 26 November 2021). It should be noted that wind speed data are not available in this dataset. The ERA5-Land data from the COPERNICUS database has been used to obtain the wind speed data in this study. The ERA-5 Land data are available from the year 1981 to the present date at an hourly time scale at 0.1° spatial interval covering the entire globe. The data is available from the link https://cds.climate.copernicus.eu/cdsapp#!/dataset/reanalysis-era5-land?tab=form (accessed on 26 November 2021). Locations of the E-OBS grids and the ERA5 grids used for this study covering the Dodder River Basin are shown in Figure 1.

2.2. Soil Water Assessment Tool (SWAT) Model and Performance Evaluation

The soil water assessment tool model (SWAT) is a conceptual semi-distributed continuous-time hydrological model developed by the United States Department of Agriculture to simulate the water transport at a basin/watershed scale on a daily or sub-daily time scale [34,39]. A flowchart of the SWAT model is shown in Figure 6 which shows how watershed-related data (elevation, land cover and soil type), and meteorological data (rainfall, temperature (daily maximum and daily minimum), wind speed, relative humidity and solar radiation) are used in the different steps of the SWAT model.

The SWAT model uses the basin DEM data to delineate the watershed boundary, the river network inside the watershed, and the slope of the entire watershed. Once the watershed and the river network were generated, the land cover and the soil data for the watershed are used by the SWAT model along with the estimated slope to divide the entire watershed into multiple sub-basins. Each of those sub-basins, termed as hydrological response units (HRUs), contains the same land cover, soil, and similar slope characteristics. Once the HRUs for the watershed was formed, the flow dynamics in the watershed are modelled in two phases: the land phase and routing phase. In the land phase, SWAT simulates the hydrologic cycle for each HRU based on the following water balance equation:

$$S{W}_t=S{W}_0+{\displaystyle \sum }_{i=1}^n{\left(R_{day}-Q_{sur}-E_a-W_{seep}-Q_{gw}\right)}_i$$

(1)

where $t$ is the time in days, $S{W}_t$ is the final soil water content (mm), $S{W}_0$ is the initial soil water content (mm), $R_{day}$ is the amount of precipitation on day $i$ (mm), $Q_{sur}$ is the amount of surface runoff on day $i$ (mm), $E_a$ is the amount of evapotranspiration on day $i$ (mm), $W_{seep}$ is the amount of water entering the vadose zone from the soil (percolation) on day $i$ (mm), and $Q_{gw}$ is the amount of return flow on day $i$ (mm) [40], and $n$ is the total number of days considered in the simulation study. Once the land phase is modelled, the water is routed through reservoirs and channels in the routing phase by using the variable storage coefficient method developed by Williams [41]. The water balance of each HRU is represented by four storage volumes including: snow, soil profile (0–2 m), shallow aquifer (typically 2–20 m) and deep aquifer (>20 m).

SWAT uses the SCS curve number approach to estimate the surface runoff in the water balance equation, and the Penman–Monteith equation [42] to estimate the evapotranspiration. To estimate evapotranspiration at a daily time scale, daily values of four meteorological variables namely, mean temperature, wind speed, relative humidity and solar radiation are required. The daily mean temperature was estimated as the average of the daily maximum and daily minimum temperature. In the SWAT model, lateral subsurface flow in the soil profile (0–2 m) is calculated simultaneously with percolation. Using the degree of slope, slope length and saturated hydraulic conductivity of the soil layer, SWAT performs a kinematic storage routing to predict the lateral flow in each soil layer. Lateral flow occurs in situations where the storage in any layer exceeds field capacity after percolation. SWAT creates a shallow aquifer storage to simulate the groundwater flow contribution to the total streamflow [43]. Percolation from the bottom of the root zone is considered as recharge to the shallow aquifer. The SWAT model simulates surface runoff, sub-surface runoff and groundwater recharge as its model output.

In this study, the SWAT model is used to simulate surface runoff at the basin outlet for the historical period by considering different land cover data as it changes over time. This study uses land cover data obtained at different time periods and uses them separately to simulate the surface runoff for the entire historical time period. In general, it can be expected that the simulated runoff should be closer to the observed runoff/streamflow at the basin outlet when appropriate land cover data has been used to run the SWAT model. To verify this hypothesis, this study quantifies the error in runoff simulations at different time periods for the entire historical period using different land cover data. In situations where the hypothesis is correct, for a land cover data collected at a chosen year, the error in runoff simulation using SWAT model should be at a minimum; conversely the error in runoff simulation should be higher for the simulated years that are further away from the year when the land cover data was collected. The land cover data were collected on multiple different years in this study, and the error in monthly mean and monthly maximum runoff obtained by comparing the model simulation with observed runoff for each land cover data was estimated. using the Nash–Sutcliffe efficiency (NSE), as follows:

$$NSE=1-\frac{{{\displaystyle \sum }}_{t=1}^T{\left[\widehat{y}\left(t\right)-y\left(t\right)\right]}^2}{{{\displaystyle \sum }}_{t=1}^T{\left[y\left(t\right)-\overline{y}\right]}^2}$$

(2)

where $y\left(t\right)$ is the observed runoff/streamflow at time $t$, $\widehat{y}\left(t\right)$ is the SWAT model predicted runoff at time $t$, $T$ is the number of data points (monthly mean or monthly maximum), and $\overline{y}={\displaystyle \sum }_{t=1}^{T}y\left(t\right)/T$ is the mean observed flow for the chosen period.

The value of NSE ranges from $\left(-\infty ,\left.1\right]\right.$, where the performance measure values tend to unity as model predictions become closer to the observed values.

2.3. Kappa and Generalized Extreme Value Distribution

The simulated monthly mean and maximum runoffs obtained were then used to estimate the runoff quantile corresponding to different return periods. According to the extreme value theory [44], the monthly maximum runoff can be modelled using the extreme value distribution. This study considers the generalised extreme value (GEV) distribution to model the monthly maximum runoff obtained from the SWAT model. There appear to be no established distributions that are being used to model the monthly mean runoff. Only a limited number of studies [45,46,47] have performed distribution fitting to the monthly mean streamflow/runoff. Blum et al. [45] used a five parameter Wakeby (WAK5) distribution, four-parameter kappa (KAP4) distribution, and three-parameter generalized Pareto (GPA) and log-normal (LN3) distributions to fit the mean monthly streamflow from 398 stream gauges in the conterminous United States and noted that kappa is the most suited distribution. Svensson et al. [46] analysed mean monthly streamflow obtained from 121 catchments in the United Kingdom using 12 frequency distributions: KAP4, three-parameter distributions kappa (KAP3), generalised logistic (GLO3), Pearson type 3 (PE3), generalised extreme value (GEV), Tweedie, and two-parameter distributions generalized logistic (GLO2), gamma (GAM), Weibull, Gumbel, normal and log-normal (LN2). Based on the analysis they noted that both three and four-parameter kappa distributions fit the data well. The other study [47] used LN3, GEV, Weibull, Gamma, Gumbel and normal distributions to fit annual mean streamflow data at the Tana River in Kenya and found LN3 to be the best-fit distribution. In this study, four-parameter kappa distribution has been used to fit the monthly mean runoff. The parameters for KAP4 and GEV (in the case of monthly maximum runoff) were estimated using the L-moments. Details of parameter estimation and quantile estimation corresponding to different return periods using KAP4 and GEV are provided below.

The four-parameter kappa distribution has been used to fit the monthly mean runoff generated using the SWAT model. The cumulative distribution function (CDF) of kappa is given as,

$$F\left(x\right)={\left[1-h{\left\{1-\frac{k}{\alpha }\left(x-\xi \right)\right\}}^{1/k}\right]}^{1/h}$$

(3)

where $\xi$ is the location parameter, $\alpha$ is the scale parameter and $k,h$ are the shape parameters respectively. The distribution parameters can be estimated using the first two L-moments ${\lambda }_1$, ${\lambda }_2$, L-skewness ${\tau }_3$ and L-kurtosis ${\tau }_4$ based on the following relationship,

$$\begin{gathered} {\lambda }_1=\xi +\alpha \left(1-g_1\right)/k \\ {\lambda }_2=\alpha \left(g_1-g_2\right)/k \\ {\tau }_3=\left(-g_1+3g_2-2g_3\right)/\left(g_1-g_2\right) \\ {\tau }_4=\left(-g_1+6g_2-10g_3+5g_4\right)/\left(g_1-g_2\right) \\ {g}_r=\left\{\begin{array}{c}\frac{r\mathsf{\Gamma }\left(1+k\right)\mathsf{\Gamma }\left(r/h\right)}{h^{1+k}\mathsf{\Gamma }\left(1+k+r/h\right)}, h>0\\ \frac{r\mathsf{\Gamma }\left(1+k\right)\mathsf{\Gamma }\left(-k-r/h\right)}{{\left(-h\right)}^{1+k}\mathsf{\Gamma }\left(1-r/h\right)}, h<0\end{array}\right. \end{gathered}$$

(4)

where $\mathsf{\Gamma }\left(z\right)={{\displaystyle \int }}_0^{\infty }{t}^{z-1}{e}^{-t}dt$ is the gamma function.

Considering there are $n$ number of data points in the monthly mean time series $\left\{x_1,x_2,\cdots ,x_n\right\}$, the series is first arranged in an ascending order $\left\{x_{1:n},x_{2:n},\cdots ,x_{n:n}\right\}$ such that $x_{1:n}\le x_{2:n}\le \cdots \le x_{n:n}$. Subsequently, the probability-weighted moments $b_r$ can be estimated as,

$$b_r=\frac{1}{n}{\left(\begin{array}{c}n-1\\ r\end{array}\right)}^{-1}{\displaystyle \sum }_{j=r+1}^n\left(\begin{array}{c}j-1\\ r\end{array}\right)x_{j:n}, r=0,1,2,3,4,\cdots$$

(5)

The first four sample L-moments can be estimated as,

$$\begin{gathered} {\lambda }_1=b_0 \\ {\lambda }_2=2b_1-b_0 \\ {\lambda }_3=6b_2-6b_1+b_0 \\ {\lambda }_4=20b_3-30b_2+12b_1-b_0 \end{gathered}$$

(6)

where ${\tau }_3={\lambda }_3/{\lambda }_2$ and ${\tau }_4={\lambda }_4/{\lambda }_2$.

The parameters of the kappa distribution can be solved numerically. Further details on the parameter estimation can be found in Hosking and Wallis [48].

Once the parameters of the kappa distribution are obtained, the quantile corresponding to a non-exceedance probability $F$ or a return period $T=1/\left(1-F\right)$ can be obtained as,

$$x\left(F\right)=\xi +\frac{\alpha }{k}\left\{1-{\left(\frac{1-F^h}{h}\right)}^k\right\}$$

(7)

The generalized extreme value (GEV) distribution used for fitting the monthly maximum runoff data are defined by its CDF as,

$$F\left(x\right)=exp\left\{-{\left[1-\frac{k}{\alpha }\left(x-\xi \right)\right]}^{\frac{1}{k}}\right\}$$

(8)

where $\xi$, $\alpha$ and $k$ are respectively the location, scale, and shape parameters of the distribution. Once the L-moments ${\lambda }_1$, ${\lambda }_2$ and L-skewness ${\tau }_3$ were obtained from the monthly maximum runoff data, the GEV distribution parameters can be estimated as [48],

$$\begin{gathered} k=7.8590c+2.9554c^2 \mathrm{where} c=\frac{2}{3+{\tau }_3}-\frac{log\left(2\right)}{log\left(3\right)} \\ \alpha =\frac{{\lambda }_{2}k}{\left(1-2^{-k}\right)\mathsf{\Gamma }\left(1+k\right)} \\ \xi ={\lambda }_1-\alpha \left\{1-\mathsf{\Gamma }\left(1+k\right)\right\}/k \end{gathered}$$

(9)

where $\mathsf{\Gamma }\left(·\right)$ is the gamma function.

The quantile function for GEV can be estimated as,

$$x\left(F\right)=\left\{\begin{array}{ll}\xi +\frac{\alpha }{k}\left[1-{\left(-lnF\right)}^k\right],& if k\ne 0\\ \xi -\alpha \left[ln\left(-lnF\right)\right],& if k=0\end{array}\right.$$

(10)

The non-exceedance probability $F$ can be replaced with the return period $T$ where $T=1/\left(1-F\right)$. The flood quantile functions are generally plotted against the return period $T$.

The CDF/quantile function obtained based on fitting the SWAT simulated runoff can be compared to the CDF/quantile function generated using the observed monthly mean/monthly maximum streamflow. This study uses the Weibull plotting position formula to generate the empirical CDF for the observed streamflow. The Weibull plotting position formula is a commonly used approach to develop empirical CDF for daily/monthly streamflow [49]. The Weibull plotting position has been found to have unbiased exceedance probabilities for all distributions.

3. Results and Discussion

The watershed upstream of Waldron bridge river gauging station was delineated and gives a basin area of 92.266 km². For comparison purposes, the watershed area estimated by the Office of Public Works (OPW), Ireland based on the topographic maps of Ireland is documented as 94.3 km², indicating that the difference in the basin area obtained based on topographic maps and EU-DEM based SWAT model is around 2%, which as deemed acceptable for practical purposes.

Based on the DEM, the slope of the watershed has been estimated using the SWAT model and classified into five classes: 0–2°, 2–5°, 5–10°, 10–25° and >25° (see Figure 3b). The SWAT model combines the land cover, soil and slope data to create the hydrological response units (HRUs). In this analysis, five different SWAT models had been generated by using five land cover maps obtained for the years 1990, 2000, 2006, 2012 and 2018. Due to changes in the land cover maps for each model, the HRUs generated by the model also varied for each case. The total number of HRUs generated while using 1990, 2000, 2006, 2012 and 2018 land cover maps are respectively 1446, 1481, 1497, 1466 and 1504.

Once the HRUs were developed for each SWAT model, the five meteorological variables at a daily time scale were fed into each SWAT model to simulate surface runoff for the period of 1 January 1993 to 31 December 2019. In general, the SWAT model parameters are being calibrated to improve the performance of the SWAT model. However, this study primarily focuses on the effectiveness of the SWAT model while considering different land cover data. Hence, for every SWAT model, the DEM, soil, and weather data as well as the model parameters were retained to be the same, while varying only the land cover data. This ensures to compare the performance of different land cover data in the simulation of runoff at watershed outlet using SWAT model.

The monthly mean and the monthly maximum runoff were then extracted from the daily simulated runoff value for each of the five SWAT models and compared to the observed monthly mean and monthly maximum streamflow values obtained based on the daily observed streamflow data available at Waldron bridge. The error in monthly mean and monthly maximum runoff was quantified in terms of NSE for each year and the values are shown in

Table 3. Nash–Sutcliffe efficiency obtained for the period of 1993–2019 estimated based on the observed and simulated runoff using five different land cover data corresponding to (a) monthly mean runoff and (b) monthly maximum runoff at Waldron bridge. The best-fit SWAT models are shown in bold.

Year	(a) Monthly Mean Runoff					(b) Monthly Maximum Runoff
	LC1990	LC2000	LC2006	LC2012	LC2018	LC1990	LC2000	LC2006	LC2012	LC2018
1993	0.55	0.51	0.47	0.45	0.41	0.56	0.53	0.47	0.41	0.44
1994	0.85	0.77	0.80	0.82	0.74	0.58	0.46	0.49	0.55	0.52
1995	0.65	0.60	0.55	0.57	0.62	0.56	0.47	0.50	0.40	0.43
1996	0.52	0.43	0.49	0.54	0.46	0.52	0.43	0.47	0.55	0.49
1997	0.59	0.52	0.50	0.55	0.46	0.57	0.51	0.55	0.48	0.45
1998	0.54	0.47	0.44	0.51	0.41	0.52	0.48	0.42	0.45	0.38
1999	0.50	0.57	0.45	0.53	0.42	0.57	0.46	0.51	0.48	0.39
2000	0.80	0.83	0.77	0.71	0.75	0.82	0.91	0.85	0.79	0.88
2001	0.45	0.53	0.34	0.48	0.39	0.26	0.44	0.51	0.23	0.35
2002	0.77	0.74	0.71	0.79	0.82	0.43	0.55	0.50	0.39	0.47
2003	0.46	0.60	0.43	0.55	0.51	0.54	0.57	0.48	0.51	0.45
2004	0.54	0.46	0.58	0.51	0.43	0.57	0.53	0.66	0.60	0.63
2005	0.51	0.46	0.59	0.54	0.43	0.48	0.51	0.43	0.36	0.40
2006	0.47	0.39	0.52	0.44	0.36	0.53	0.40	0.56	0.51	0.45
2007	0.36	0.41	0.50	0.46	0.29	0.45	0.49	0.57	0.53	0.40
2008	0.58	0.43	0.52	0.47	0.39	0.48	0.43	0.51	0.45	0.39
2009	0.41	0.44	0.52	0.47	0.37	0.47	0.37	0.52	0.34	0.31
2010	0.48	0.44	0.41	0.50	0.38	0.44	0.49	0.55	0.52	0.46
2011	0.51	0.54	0.57	0.45	0.48	0.78	0.83	0.86	0.89	0.81
2012	0.47	0.44	0.41	0.54	0.35	0.51	0.47	0.44	0.54	0.41
2013	0.47	0.41	0.38	0.51	0.45	0.44	0.41	0.47	0.52	0.49
2014	0.69	0.66	0.63	0.71	0.60	0.51	0.42	0.45	0.54	0.48
2015	0.73	0.70	0.67	0.77	0.64	0.51	0.54	0.46	0.40	0.43
2016	0.87	0.82	0.85	0.67	0.90	0.48	0.52	0.42	0.34	0.40
2017	0.86	0.81	0.80	0.84	0.91	0.47	0.38	0.42	0.35	0.52
2018	0.61	0.58	0.63	0.55	0.66	0.88	0.79	0.82	0.85	0.90
2019	0.71	0.66	0.68	0.63	0.74	0.49	0.52	0.55	0.46	0.58

. The land cover map for the year 1990 (LC1990) was developed based on data collected for the period of 1986–1998, while the same for LC2000 was 1999–2001, LC2006 was 2005–2007, LC2012 was 2011–2012 and LC2018 was 2017–2018, respectively. Hence, each of the land cover maps was considered to represent the following years in this study: LC1990 for the years 1993–1998, LC2000 for the years 1999–2003, LC2006 for the years 2004–2009, LC2012 for the years 2010–2015 and LC2018 for the years 2016–2019. The SWAT simulated month mean/monthly maximum runoff that is closer to the observed runoff for a chosen year obtained using a chosen land cover data (out of the available five maps) was identified based on the three performance measures. The best-fit land cover map based on the NSE values for each year ranging from 1993 to 2019 are shown in Figure 7. Table 3 and Figure 7 indicates that, in general, the performance of the SWAT model is superior when the land cover map representing a chosen year is used to simulate the runoff for that year. The difference in NSE value for a chosen year corresponding to different land cover data ranged from 0.1 to 0.23 for monthly mean runoff and 0.11 to 0.28 for monthly maximum runoff, indicating the importance of selecting appropriate land cover data while developing the SWAT model. Out of the simulated 27 years data, the best-fit SWAT model for the monthly mean data simulation matched with the expected model for 23 years (85%), while that for the monthly maximum data simulation matched for 20 years (74%). The results indicate the advantage of using the most appropriate land cover data while simulating the runoff using a hydrological model.

The NSE values from Table 3 indicate that the NSE values corresponding to the best-fit model range from 0.50–0.91 for monthly mean runoff, where the NSE value is less than 0.6 for 16 years (59%), between 0.6 to 0.8 for 6 years (22%), and greater than 0.8 for 5 years (19%). In the case of the monthly maximum runoff, the range of NSE value corresponding to the best-fit model is 0.51–0.91, where the NSE is less than 0.6 for 23 years (85%), between 0.6 to 0.6 for 1 year (4%) and greater than 0.8 for 3 years (11%). This result indicates that even though the use of effective land cover data considerably improves the accuracy of the SWAT model in runoff simulation, the performance of the developed SWAT models can be termed fair at the best for the majority of the years. One of the primary reasons is that the performance of the SWAT model depends not only on the land cover data but also on the DEM, the soil, and most importantly on the five weather variables [21,34,43]. Since no in-situ measurements were available inside the Dodder River catchment for temperature, relative humidity, wind speed and solar radiation, gridded data were considered for the development of the SWAT models obtained from the COPERNICUS E-OBS and ERA-5 land database. Out of the two rain gauges that were located inside the catchment, rainfall station S5623 had 3% data missing for the period 1993–2019, while rainfall station S6623 had 40% missing data for the considered period. The use of appropriate weather data is essential for the improvement of the model performances. Furthermore, it needs to be noted that this study considered 25 m resolution EU-DEM to identify the watershed boundary, generate the stream network, and estimate the slope of the watershed. A finer resolution DEM might also be useful in improving the performance of the SWAT models. However, it needs to be noted that the objective of this manuscript is to investigate the relative performance of the SWAT models corresponding to different land cover data while retaining all the other input data (DEM, soil, weather variables) to be the same.

A composite runoff time series for monthly mean and monthly maximum time period were constructed by considering different land cover data across different time periods. The land cover data selected to simulate runoff for a chosen year are described in the previous paragraph. A scatter plot between the observed monthly mean and composite monthly mean, as well as observed and composite monthly maximum runoff were plotted in Figure 8. A point belonging to the 45° line (shown in green) indicate that the observed and model simulated composite runoff are the same, while a point located away from the line indicates that model-simulation is different from the observed value. The distance of the point from the 45° line measures the error in the model simulation (under-prediction or over-prediction depending on the side). It should be noted that use of real-world weather data is expected to improve the runoff obtained from the composite SWAT model.

Subsequently, quantile curves were generated using the monthly mean and maximum time series runoff simulations. The kappa distribution was considered for the monthly mean series while the generalized extreme value (GEV) distribution was fitted for the monthly maximum data. The parameters of the kappa distribution and the GEV distribution are provided in

Table 4. Parameters of (a) kappa distribution estimated using the monthly mean flow and (b) generalized extreme value distribution estimated using monthly maximum flow data simulated using the SWAT models by considering composite time series, and land cover maps obtained for the year 1990, 2000, 2006, 2012 and 2018.

Data Type	Location ( $\mathit{\xi }$ )	Scale ( $\mathit{\alpha }$ )	Shape
			$\mathit{k}$	$\mathit{h}$
Composite	0.345	1.668	0.126	0.781
LC1990	−0.137	2.063	0.187	1.072
LC2000	−0.097	2.024	0.178	1.052
LC2006	−0.064	1.983	0.171	1.036
LC2012	−0.116	2.030	0.179	1.066
LC2018	0.261	1.727	0.138	0.830
(a) Kappa distribution parameters
Data Type	Location ($\mathit{\xi }$)	Scale ($\mathit{\alpha }$)	Shape ($\mathit{k}$)
Composite	2.313	2.487	−0.395
LC1990	2.344	2.470	−0.396
LC2000	2.369	2.511	−0.392
LC2006	2.387	2.503	−0.393
LC2012	2.330	2.442	−0.396
LC2018	2.410	2.488	−0.390
(b) Generalized extreme value distribution parameters

and the quantile curves are shown in Figure 9. The quantile curves were plotted corresponding to different return periods ranging from 100 years to 1000 years.

The quantile curves obtained based on the SWAT simulations were compared to the quantile curves generated using the historical observations at Waldron bridge river gauge, located at the basin outlet. The quantile curve from the historical data was generated using the Weibull plotting position. For comparison purposes, the quantile curves obtained based on simulated runoff corresponding to each of the five SWAT models, where fixed land cover data were only considered for each simulation, are also plotted. The comparison of the quantile curves indicates that the composite runoff is closest to the quantile curve obtained from historical data. However, it can be noted that all of the different model simulations underpredict the flow quantile when compared to the historical data. However, the quantile values obtained using the composite SWAT simulated time series were closer to the historical one. It needs to be noted that the quantile curves were developed based on 27 years of historical data, while the curve estimate runoff/streamflow quantiles up to 1000 years. Based on the short period of data, it is unwise to claim that the historical quantile curve represents the real-world scenario. The availability of longer historical records can improve the development of the true quantile curve. Furthermore, the use of real-world weather data in the development of the SWAT model can improve the quantile functions generated based on the model outputs.

Based on the overall results it can therefore be concluded that a dynamic approach is necessary to simulate runoff using the SWAT hydrological model where the land cover data needs to be updated over time to obtain realistic runoff.

4. Conclusions

The impact of changes in land cover over time with respect to runoff simulations using a physically-based hydrological model called soil water assessment tool has been investigated in this study. The majority of the hydrological models assume that land cover remains fixed across the periods that they are used to simulate, whereas this study has developed multiple SWAT models that compare different land cover data that has been obtained at different time periods. The comparison of the simulated runoff with real-world observations revealed that the model performance is highest when the land cover data used in the model has been collated as close to the time period across which the simulations need to be performed. Based on the results, it has been concluded that a dynamic land cover database should be considered for modelling runoff in order to obtain better performance. This study, therefore, proposes a method whereby the entire time period of simulation is subdivided into multiple time segments. The subdivision should be based on the land cover data that are available for the river basin at several temporal scales. Subsequently, unique land cover data needs to be identified for each of those subdivisions and the SWAT model is then be run using different land cover data for different time periods. Following this, a composite time series of runoff simulations can be constructed from which monthly mean and monthly maximum runoff for the entire period of simulation can be extracted. This methodology was carried out on the studied watershed (the River Dodder) and the monthly mean data were fitted using four-parameter kappa distribution and the monthly maximum data were fitted using generalized extreme value distribution to obtain quantile functions corresponding to different return periods. The performance of the composite approach was found to be superior in obtaining effective flood quantiles. Furthermore, the performance of the composite SWAT model can be improved by considering real-world weather data. Overall, this research shows that in order to simulate runoff for future periods in a river basin, it is important to understand the land cover dynamics and use it to project future scenarios of land cover patterns before simulating future runoff. Ongoing research focuses on the development of land cover change model for river basins in Ireland.