A modelling approach to simulate Chlorophyta and Cyanobacteria biomasses based on historical data of a Brazilian urban reservoir

A process-based model was used to simulate a hindcast based on the worst historical water quality condition of a tropical urban reservoir. Paranoá Lake is located in Brasília-DF, Brazil, and went through intense eutrophication in the 70s and 90s, with an important cyanobacterial bloom event in 1978. The parameters of phytoplankton were calibrated, focusing on the group of Chlorophyta (green algae) and Cyanobacteria (blue-green algae) at four depths (1, 10, 15 and 20m). The results indicated that the model was able to reproduce the Cyanobacteria biomass in comparison with the observations (RMSE=22-29.10 -3 mgC L -1 ). On the other hand, the simulated Chlorophyta biomass showed good agreement with the observed data only in the bottom layer (RMSE=29. 10 -3 mgC L -1 at 20m). In the hindcast simulation, the model was able to predict a significant increase in cyanobacterial biomass facing a water quality deterioration. In the meantime, the simulated Chlorophyta biomass decreased, which may indicate the phytoplankton group succession in response to the environmental conditions.


INTRODUCTION
Water bodies in urban areas play an important role in providing water for multiple uses, especially drinking water supply and recreational activities. The deterioration of water quality in these environments is the result of several factors, mainly the increase in nutrient loadings due to anthropogenic activities and higher sewage release in those watersheds. Therefore, extended eutrophication is frequently observed, as well as the intensification of cyanobacterial blooms (Wells et al., 2015;Reynolds and Walsby, 1975).
Ecological modeling is the result of the demand to comprehend biotic responses to environmental impacts, as abiotic criteria are not sufficient to represent the full complexity of aquatic systems. In recent decades, technological advances, as well as conceptual improvements and increased access for different social actors in the modeling process (Robson, 2014) have allowed the expansion of these tools (Ulańczyk et al., 2018;van der Linden et al., 2015).
Paranoá Lake is located in the urban area of Brasília, Federal District, Brazil. During the 1970s-'90s the reservoir underwent extended eutrophication with water quality improvement only at the end of the 90's. Until 2017, some of the multiple uses of the reservoir were dilution of effluents, recreation and power generation. Nunes et al. (2020) evaluated land use changes in the Paranoá watershed and, together with changes in the uses of the lake, performed calculation of the water balance in the period between 1982 and 2017. The study pointed out the effects on the lake related to the increasing urbanization and the consequent changes in land use in affluent basins. On October 2, 2017, the reservoir also started to provide drinking water supply (700 L s -1 , Lacerda, 2017) for Brazilian Federal District.
Considering the importance of the uses of Paranoá Lake and its history of eutrophication, studies to understand the distribution and dynamics of phytoplankton, with special attention to the group of potentially toxic cyanobacteria, have become important. Thus, this work aimed to simulate a hindcast based on the worst year of water quality in Paranoá Lake using historical monitoring data as input. The initial conditions and meteorological forcing data were not altered. Our hypotheses were: 1) The ecological model would be able to represent the water quality deterioration with focus on the dynamics of two groups of phytoplankton based only on altered inflow concentrations; 2) This simple approach would also be able to predict the responses of phytoplankton groups to water quality deterioration in order to understand reservoir management trends and directions.

Study Area
Paranoá Lake Basin is located in the Brazilian Federal District -DF (15°47′ S, 47°49′ W) and is formed by the hydrographic units: Bananal, Riacho Fundo, Gama and Torto Streams and Paranoá Lake. With a basin drainage area of 1034 km², the reservoir has a surface area of 37.5 km² and total volume of 498.10 6 m³.
According to the Köppen climate classification, a tropical climate with dry winter (Aw) predominates in the DF. The Aw climate is characterized by a wet season in summer, from October to April, and a well-defined dry season in winter, generally from May to September. The historical data from the Brazilian National Institute of Meteorology (INMET) shows that between 2000 and 2018, the annual average rainfall was 1498 mm and the average air temperature was 21.5ºC. The water residence time is ~300 days (Nunes et al., 2020).
Paranoá Lake water quality is susceptible to degradation because its use for the dilution of sewage from two treatment plants (ETE Brasília Norte and ETE Brasília Sul), discharges from the rainwater drainage system of Brasília and also different sources of diffuse pollution (Franz et al., 2014;Machado e Baptista, 2016). Regarding affluent basins, the greatest pollution load comes from the Riacho Fundo tributary. Its sub-basin presents a high rate of urbanization and soil quality degradation and has the highest average daily flow (~4.04 m³ s -1 ) compared to the other tributaries.
Between 1976 and 1998, Paranoá Lake showed a eutrophic pattern, when the maximum chlorophyll concentration reached 166.8 µg.l -1 near the dam. In 1993 and 1994, after the expansion and inclusion of tertiary treatment in the two sewage treatment plants, and the interconnection of the sewage systems of the Riacho Fundo Basin, water quality improved.
Chlorophyll-a concentrations decreased about 10 times between 1999 and 2011 (Souza, 2013). During the worst water-quality period, mainly in the 70's, studies focusing on phytoplankton characteristics and distribution in Paranoá Lake highlighted the cyanobacteria dominance (Giani and Pinto-Coelho, 1986;Branco and Senna, 1994;1996). Recently, Paranoá Lake has shown an increase of trophic status (Machado and Baptista, 2016) and water quality degradation, resulting in new cyanobacterial blooms (Batista and Fonseca, 2018). The monthly average chlorophyll-a concentration was 21 µg. l -1 and the monthly average of total phosphorus and total nitrogen was 25 µg. l -1 and 1.4 mg. l -1 , respectively, from 1976 and 2010, according to monitoring performed by the Environmental Sanitation Company of the Federal District -CAESB.

Data Availability
The database for this study included measurements from 12 monitoring stations in the Paranoá Lake Basin ( Figure 1, Table 1). The boundary conditions for the simulations include inflows, outflows, meteorological and morphometric data, as well as water quality and phytoplankton time-series. The inflows and bathymetry data were provided by the Environmental Sanitation Company of the Federal District -CAESB. Outflows were made available by the Brasilia Energy Company (CEB) and meteorological forcing data by the INMET for the Brasilia Meteorological Station, located 15 km from the reservoir.
The water quality data corresponded to surface measurements of the four tributaries (Riacho Fundo, Gama, Bananal and Torto), sampling stations in the arms of the reservoir (A; B; D; E) and the effluents of the sewage treatment plants, "ETE Norte" and "ETE Sul". Near to the dam, the monitoring station "C" had available measurements in six depths along the water column: surface, 1, 10, 15 and 20 meters from the surface, and 1 meter from the bottom. The data measured at the surface and 1m from the bottom corresponded to only 40% of the entire period from 1976 to 2010.Therefore, those depths were not considered in the calibration. Water temperature (°C); Salinity (mg L -1 ); Dissolved oxygen (mg L -1 ); P-PO4 (mg L -1 ); Total phosphorus (mg L -1 ); N-NO3 (mg L -1 ); N-NH+4 (mg L -1 ); Chlorophyta biomass (mg L -1 ); Cyanobacteria biomass (mg L -1 );

Model Description
The General Lake Model-(GLM v.2.2) is a one-dimensional hydrodynamic model that calculates vertical profiles of temperature, salinity and density by representing the effects of inflows and outflows, mixing, heating and cooling of water bodies (Hipsey et al., 2019). The Aquatic Ecodynamics Model Library (AED v.1.0) is composed of biogeochemical and ecological modules coupled with GLM through the Framework for Aquatic Biogeochemical Modelling (FABM) (Hipsey et al., 2013). The modules include the components of the carbon, nitrogen and phosphorus cycling and others like oxygen, organic matter and organisms such as the different groups of phytoplankton and zooplankton.
The parameters related to phytoplankton dynamics are based on the processes of photosynthesis and nutrient uptake, respiration, excretion and mortality of the community. The Rev. Ambient. Água vol. 16 n. 5, e2713 -Taubaté 2021 module includes up to three phytoplankton taxa: Bacillariophycea, Chlorophyta and Cyanobacteria. The phytoplankton biomass of a given group is simulated in carbon units (mmol C.m³) and the group can be configured to have a constant C:N:P ratio or variable ratios and uptake dynamic in response to water column changes and intracellular processes. In the present study, we have set up the model to simulate constant C:N:P ratio. It is also possible to configure equations that consider or not photoinhibition as limiting the growth of each group.

Modelling Setup
The sensitivity analysis and hydrodynamic calibration of the GLM for Paranoá Lake were performed in a previous study (Barbosa et al., 2015). Due to the presence of gaps in the meteorological forcing data, the period in which continuously monitored data were available was from 2007-03-01 to 2009-03-31 (762 days). This was the period selected in this study for model calibration.
The AED boundary conditions include biogeochemical and ecological time-series from the tributaries and the reservoir for the first day of simulation. Linear interpolation was performed through the software MATLAB, Version 2011 to downscale the monthly forcing data into a daily time step.
Calibration was performed using the observed data from Station "C", chosen for data availability along the water column (at 1, 10, 15 and 20 meters from the surface) and for being located in the deepest region of the reservoir (~29m). The data from the other sampling stations and effluents from the sewage treatment plants were considered as inputs in the model. The nutrient concentrations measured in the reservoir arms were also considered as input data for the simulations since the Paranoá arms have independent water quality features in comparison with the main body of the reservoir where monitoring Station C is located. Thus, we used as input variables dissolved oxygen, nitrate, ammonium and soluble reactive phosphorus.
Once the hydrodynamic model was calibrated, the sensitivity analysis of the phytoplankton parameters was started in the same period (762 days). To simulate the phytoplankton module, it was necessary to initialize each biogeochemical module (C, N, P), assuming their original configurations and default parameters. The sensitivity analysis and calibration focused on the biomass and distribution of Chlorophyta and Cyanobacteria along the water column. These two groups were chosen for their dominance in the reservoir. For the simulations, the total biomass data of the two groups were used.
The subgroups of parameters related to growth, light, respiration, nitrogen and phosphorus of the phytoplankton were tested to observe the sensitivity of Chlorophyta and Cyanobacteria biomass to variations of 20% around their default values. The model was configured without considering phytoplankton photoinhibition. Likewise, since data on the zooplankton have not been included, grazing was not considered. Mortality due to grazing was considered in the phytoplankton mass balance equation. Phytoplankton mortality was related to parameters that represented losses by excretion and respiration.
The sensitivity analysis was performed by manual adjustment, searching for parameter values close to the range suggested by Hipsey et al. (2013) and checking their influence on the results of phytoplankton biomass. After the sensitivity analysis, the calibration process was begun, considering the sensitivity parameters individually and also the ensemble for the calculation of the Root Mean Squared Error (RMSE) and Mean Absolute Error (MAE) metrics. The target of calibration was to reduce the difference between observed and simulated biomass of Chlorophyta and Cyanobacteria at the four depths, using these metrics.

Hindcast Simulation
The hindcast simulation was based on water quality data from the most critical period during the eutrophication process of Paranoá Lake, corresponding to 1978. At that moment, there was a cyanobacterial biomass increase with the occurrence of blooms in monitoring Station "A" between September and October 1978. At the same time, an increase in cyanobacterial biomass was registered in Station "C" at the one-meter depth (369. 10 -3 mg L -1 ) at 1m), more than twice the average recorded value between 1976 and 2010 at the same depth of station C (183. 10 -3 mg L -1 ). It should be noticed that in this work the total biomass was analyzed, with the response of the different cyanobacteria species not being considered. No records of Chlorophyta were reported in Station C that year. Another eutrophication consequence, also in 1978, was an anoxia at the bottom of Station C. The average dissolved oxygen concentrations (DO) were 7.71 mg L -1 at 1 meter from the surface and, 3.36 mg L -1 , 1.56 mg L -1 and 0.56 mg L -1 for the depths of 10, 15 and 20 meters from the surface, respectively.
The hindcast simulation was performed changing the input water quality time-series into the historical values recorded in the critical year (1978). The reservoir initial conditions and meteorological forcing data remained unchanged based on the year 2008. Due to the absence of measurements in 1978, the nutrient loads of "ETE Norte" and "ETE Sul" were not considered as inputs in the hindcast simulation. This absence is not considered a major uncertainty factor to the simulation results because of the low historical observed average flow at both stations (approximately 1.5 m³ s -1 between 2000 and 2010). The main drivers of water quality deterioration in this period were non-point source sewage discharges and insufficient treatment of sewage systems in the Riacho Fundo Basin. Input data based on historical inflows concentrations were intended to represent those drivers in the hindcast simulation.

Sensitivity Analysis and Calibration
The most sensitive parameters in the calibration of Chlorophyta and Cyanobacteria were those related to phytoplankton growth, light uptake and respiration: phytoplankton growth rate at 20°C ( ) optimum temperature (Topt), light ½ saturation constant for algal limitation (IK) and phytoplankton respiration rate at 20°C ( ). The influence of growth and respiration parameters on phytoplankton simulations was also highlighted by other authors (Huang et al., 2012;Fragoso Jr., 2008). The possibility of simulating functional groups of phytoplankton has been an opportunity recently explored (Shimoda and Arhonditsis, 2016). Most studies applying one-dimensional ecological models to phytoplankton simulations are still based on chlorophyll-a concentrations to represent phytoplankton biomass (Vinçon-Leite et al., 2017;Fadel et al., 2017;Cavalcanti et al., 2016). Overall, the model underestimated the Chlorophyta and Cyanobacteria biomasses (Table 2, Figure 2). The calibration showed satisfactory performance for Cyanobacteria biomass at the top layers, even though the model overestimated the results at 1m and underestimated at 10, 15 and 20 m. On the other hand, even after the calibration the model performance was unsatisfactory for Chlorophyta along the water column (Table 2), except for the bottom layer. However, the model was able to represent some of the seasonal variations of the Chlorophyta biomass in Paranoá Lake, especially in the wet season.
The simultaneous calibration of the two phytoplankton groups' parameters was a challenge due to their different physiological, morphological and adaptive characteristics. C. raciborski and Monoraphidium sp showed dominance in the reservoir during the period of study.
However, the GLM-AED calibration showed satisfactory results for Chlorophyta biomass at 20 meters depth (RMSE=28.8 10 -3 mgC L -1 ) and Cyanobacteria biomass at four depths (RMSE=22.1-29.210 -3 mgC L -1 ) when compared to the values found by Fadel et al. (2019) for Microcystis aeruginosa (RMSE=21.9-28.48 10 -3 mgC L -1 ) and Chrysosporum ovalisporum (RMSE=14-25.9 10 -3 mgC L -1 ) at Karaoun Reservoir. Even though the model had weak performance to predict the Chlorophyta biomass at the top layers (Table 2), it was able to predict 7 A modelling approach to simulate Chlorophyta and … Rev. Ambient. Água vol. 16 n. 5, e2713 -Taubaté 2021 the biomass increase, generally observed in the summer.  -23.48 -52.11 -49.26 -12.12 7.37 -4.89 -11.37 -13.00 Although we have used the DO and nutrients concentrations of the reservoir as initial conditions, we didn't perform calibration of the aquatic biogeochemistry parameters due to data unavailability at the four depths. Despite the absence of biogeochemical calibration and the likely uncertainties added to the modelling framework, the phytoplankton parameters and processes also consider the uptake of C, N and P to simulate the biomass of the groups. Previous studies performed calibration of DO, nitrogen and phosphorus parameters before the simulation of phytoplankton organisms to improve the simulation results (Bucak et al., 2018;Fenocchi et al., 2019) and that can be an alternative to improve the current modeling approach.

Water Quality Decline
The hindcast simulation was able to represent the higher presence of Cyanobacteria in the top layers of Paranoá Lake in the wet season (Figure 3b), as described by Batista and Fonseca (2018), except for November and December. On the other hand, the model represented a Chlorophyta biomass decrease (Figure 3a). We suggest that a plankton succession may be in response to environmental conditions related to poor water quality, such as increased turbidity. The predominance of Cyanobacteria is attributed to their adaptive characteristics, such as the capacity of migration in the water column, tolerance to low light, ability to fix atmospheric nitrogen and resistance to herbivores (WHO, 2003). Such adaptations reinforce the role of Cyanobacteria as indicators and the dominance tendency under eutrophic conditions, which may, in the future, lead to an increase in their distribution related to the deterioration of water quality and climate warming (Huisman et al., 2018). The present study provides a first attempt to simulate two groups of phytoplankton in Paranoá Lake using a one-dimensional ecological model. Despite the simple approach, the results of the proposed hindcast were satisfactory and allowed predictions of the behavior of phytoplankton groups in the case of water quality deterioration in the Paranoá Lake.
The biotic community's representation is still a challenge, due to the difficulties of systematizing their complex behavior. However, improvements in monitoring efforts to generate long-term ecological data are needed for the use of more complex ecological models. Such tools can generate important answers for water quality management, incorporating climate and socioeconomic scenarios, contributions of the watershed and also lake characteristics and inputs of its main tributaries (Janssen et al., 2019).

CONCLUSIONS
The sensitivity analysis and calibration of phytoplankton identified the parameters related to growth, light, and respiration as being the most important for the simulation of Chlorophyta and Cyanobacteria in Paranoá Lake. Despite the challenge of two phytoplankton group calibrations, the model was able to represent the Cyanobacteria biomass and the Chlorophyta biomass at the bottom layer in comparison with the observed time-series.
Although only prospective, the hindcast showed that the GLM-AED coupled model responded adequately to the changes in the input data for the representation of water quality deterioration, with a predominance of the cyanobacteria.
The model was able to represent the increase in cyanobacteria biomass and the consequent reduction of Chlorophyta in a situation of poor water quality at Paranoá Lake, as occurred in the eutrophication period in 1978. These results may be useful to predict the responses of phytoplankton groups to further processes of water quality deterioration in Lake Paranoá, as well as provide information for the environmental management of the reservoir.
Despite the complexity regarding the required field data, the calibration of each ecological/biogeochemical module is recommended when using the GLM-AED coupled model in phytoplankton simulations, to reduce the uncertainties associated with the modeling process.
We also suggest further studies on the succession of cyanobacteria and other phytoplankton groups in lakes, to better understand the prevalence of cyanobacteria in aquatic environments under eutrophication conditions.