Influence of groundwater quality indicators on nitrate concentrations in the Zagreb aquifer system

Nitrates presents one of the main groups of contaminants in the Zagreb aquifer system. Some natural groundwater quality indicators can have significant influence on their stability and mobility in the saturated zone. Correlation and multivariate statistical analyses were used to test the correlation of average values of NO 3 - with O 2 , ORP, pH, EC and temperature of groundwater, and to allocate observation wells that belong to the same clusters. ORP values didn`t relate to any observed variables, which is probably due to their variability which suggests changes in the oxidation-reduction conditions in the aquifer system. Principal component analysis was used for the determination of variables that are related to the nitrate concentrations and which were then used in cluster analysis. Other variables were excluded from cluster analysis. Three methods were used to perform cluster analysis, where the results calculated with Ward`s method were chosen as the most appropriate. In the end, two clusters were identified, one with smaller, and one with higher NO 3 - , O 2 and EC values. Observation wells from cluster 1 are generally located near the Sava River and have similar nitrate concentrations. Lack of other nitrogen species and moderately aerobic conditions suggest very fast nitrification in the shallow Holocene aquifer.


INTRODUCTION
Nitrates represent one of the most frequent groundwater conta minants in the world (PEÑAHARO et al., 2009) and in Croatia (LARVA et al., 2010;NAKIĆ et al., 2013). Due to their inability to bond in soil by adsorption, they are subject to leaching and in filtration into deeper soil layers and groundwater (CHOWDARY et al., 2005;MKANDAWIRE, 2008). When they reach ground water, they mostly depend on transport and geochemical pro cesses in saturated media. Many studies have shown strong rela tionships between agricultural activity and nitrate concentrations in groundwater (ALMASRI, 2007;PEÑAHARO et al., 2009;LI et al., 2010;HOSONO et al., 2013). Also, nitrate concentrations can be a consequence of the simultaneous influence of other an thropogenic sources, e.g. septic tanks, sewage and landfills. De spite this, there are also different natural factors that have both an influence on their stability and their mobility through the unsatu rated and saturated zones. Nitrate concentrations can depend on the depth to aquifer, recharge, aquifer media, soil media, topo graphy, impact of vadose zone media, hydraulic conductivity, sur face leaching and aquifer type (ALLER et al., 1987;LAKE et al., 2003). Nitrate leaching from the unsaturated zone is a conse quence of complex interactions between many factors, e.g. soil characteristics, nitrogen dynamics, depth to water table, land use practice, onground nitrogen loading and recharge (ALMASRI, 2003;2007). DIMKIĆ et al. (2008) stated that nitrate concentra tions can be predominantly found in aerobic media, where they have high mobility. Other nitrogen species, like nitrite and am monium ions, are generally stable in predominantly anaerobic media. There are different indicators of aerobic and anaerobic conditions, but concentrations of dissolved oxygen (O 2 ) and the oxidationreduction potential (ORP) represent the most important ones. Concentrations of dissolved oxygen between 1 and 2 mg/l, and ORP values between 150 to 200 mV, represent the approxi mate threshold values for defining the boundary between aerobic

Influence of groundwater quality indicators on nitrate concentrations in the Zagreb aquifer system
and anaerobic conditions. KOVAČ et al. (2016) showed that there is a statistically significant positive correlation between nitrate and dissolved oxygen concentrations in the Zagreb aquifer. Corre lation and multivariate statistical analyses are very often used in geochemical research and the identification of potential nitrate contaminant sources (VIDAL et al., 2000;JEONG, 2001;MO RATALLA et al., 2009). Even though nitrates are defined as one of the main group of contaminants in the study area , and are generally the consequence of anthropogenic influence, the main objective of this work was to define and quan tify the correlation of five natural indicators, i.e. O 2 , ORP, pH, electrical conductivity (EC) and temperature of groundwater, on nitrate (NO 3 ) concentrations in the Zagreb aquifer system. For this purpose, correlation and multivariate statistical analyses (principal component and cluster analysis) were used. Principal component analysis (PCA) was used for the reduction of variables that were then used in cluster analysis with the purpose of iden tifying groups of observation wells with similar characteristics.

RESEARCH AREA
The Zagreb aquifer system is located in the NW part of the Re public of Croatia ( Fig. 1) and it represents the only source of po table water for the inhabitants of the City of Zagreb and Zagreb County. The Zagreb aquifer system is designated as a part of country`s strategic water reserves. It is located between Medved nica Mountain in the north and the Vukomeričke Gorice hills in the south and it covers an area of approximately 350 km 2 . The wider region is characterized by great variability in lithology, land use, pedological features and the hydraulic properties of the aquifers.
The Zagreb aquifer system is composed of sediments of Quaternary age, deposited during the Middle and Upper Pleis tocene and Holocene. Pleistocene deposits are represented by la custrine-marshy deposits, while Holocene ones are alluvial de posits. Microfaunal and microfloral analysis were used for the definition of the main stratigraphic units (SOKAČ, 1978;HER NITZ et al., 1981). Changes in the sedimentary environment and petrographical origin of the clasts in gravels and sands was used for the determination of the boundary between Pleistocene and Holocene deposits. In the beginning of the Holocene, the Sava River started to flow and transport material from the Alps which was mainly carbonate, in contrast to the Pleistocene deposits, which were generally siliciclastic (VELIĆ & SAFTIĆ, 1991). Overall, Lower Pleistocene deposits are mostly composed of clayey silts and silty clays with sporadic interbeds and lenses of gravelly sands, while the lower and middle part of the Middle Pleistocene is composed of sands, with silts and clays discovered in the upper part. The Late Pleistocene is characterized by fre quent lateral changes in gravels, sands, silts and clays, while the Holocene is generally composed of gravels and sands (VELIĆ & DURN, 1993). Hydrogeologically, the Zagreb aquifer system is divided into three main units (Figs. 2 and 3). The first unit is over burden, which is in the most part disintegrated by anthropogenic influences. The thickness of this unit generally varies from 2 to 8 m (RUŽIČIĆ et al., 2012). Generally, Fluvisols, Stagnic Pod zoluvisols and Eutric Cambisols are developed in the research area (SOLLITTO et al., 2010). The second unit is represented by the shallowest Holocene aquifer composed mostly of gravels and sands. The third unit is the deeper Pleistocene aquifer characteri zed by frequent lateral and vertical alternations of sand, gravel and clays . The thickness of the Holocene aquifer varies from 5 to 40 m, while the deeper aquifer extends up to 60 m depth in the eastern area (NAKIĆ et al., 2011). Even though these aquifer layers are hydraulically connected, geo chemical stratification with depth is recognized. Groundwater from the Holocene aquifer generally belongs to a CaMg-HCO 3 type while the Pleistocene aquifer can be additionally characterized by CaMgNa-HCO 3 hydrogeochemical facies. Higher sodium concentrations in groundwater can also be a consequence of an thropogenic influence in some areas (VLAHOVIĆ et al., 2009;MARKOVIĆ et al., 2013).
The Holocene aquifer is an unconfined aquifer. It is in direct contact with the Sava River, while the general groundwater flow is from W/NW to E/SE. The Sava River represents the main source of recharge and the main boundary condition. POSAVEC (2006) showed that during medium and low water levels, the Sava River drains the aquifer in some areas, while during high water levels it gives water to the aquifer. This means that groundwater levels and the thickness of the unsaturated zone mainly depend on the Sava River, which significantly contributes to groundwa ter recharge in the study area (MILETIĆ & BAČANI, 1999).
Industrial development and fast growth of the City of Zagreb have affected the groundwater quality in the research area. NAKIĆ et al. (2013) have identified pesticides, nitrates, poten tially toxic metals, pharmaceuticals and chlorinated aliphatics as the main contaminants. Leakage from septic tanks, sewage sys tems and agricultural activity present the main potential sources of nitrate contamination in the Zagreb aquifer area. Nitrate trends are generally decreasing except in the eastern part on the left bank of the Sava River . It has been noted that threshold values, calculated with a VB macro BACK GROUND (NAKIĆ et al., 2007), of nitrate concentrations in groundwater of the wider Zagreb area generally range from 7.6 to 18.9 mg/l (NAKIĆ et al., 2010;KOVAČ et al., 2013;NAKIĆ et al., 2016), depending on the hydrological conditions in which they were calculated. Those concentrations suggest existence of very high ambient background nitrate concentrations in the study area which represents direct evidence of anthropogenic influence on the groundwater quality of the Zagreb aquifer. Also, ground water levels are declining, on average, for 1-2 metres every ten years, while the permanent groundwater reserves have decreased for about 4% from 19764% from to 20064% from (BAČANI et al., 2010. The main reasons for groundwater decrease are associated with the deep ening of the Sava riverbed, increased groundwater abstraction and construction of dykes along the Sava River (POSAvEC, 2006).

DATA AND METHODS
Groundwater quality data from 1991 to 2015 were used for this research. They originate from 153 observation wells of the Na tional Monitoring Programme of Croatian Waters and the moni toring programme of the Jakuševec landfill. It has been noted that some observation wells were used in the monitoring programme for only a few years, after which they were excluded. Also, dif ferent sampling intervals have been observed, from monthly to yearly. Therefore, given the sampling interval variation and in consistent exclusion and inclusion of some observation wells from the monitoring network through observed time period, the data were aggregated at the level of a given observation well. Av erage values were calculated for all the observed parameters (NO 3 , O 2 , ORP, pH, EC and temperature of groundwater) for each observation well. All values below the limit of quantifica tion for nitrate were not taken into account (237 values of NO 3 in    about 16000 analysis, i.e. ~1.5% of values). Due to the lack of ORP data for some observation wells, statistical analyses were done using data from just 126 observation wells. A normal dis tribution of selected variables was tested by the Kolmogorov-Smirnov test (D-calculated value; D 0 -critical value based on number of cases; α=0.05), using Statistica 64 (version 13.1) soft ware. All data were standardized to Z-scores. Nitrate concentra tions from three river stations (Jankomir, Petruševec and Rugvica, Fig. 1) were compared to the nitrate groundwater concentrations.
Pearson r and Spearman ρ correlation coefficients were cal culated for all the observed parameters. Parametric and nonpara metric correlation coefficients were used due to differences in the normal distribution of the selected parameters. There are different classifications concerning the interpretation of correlation co efficients. For example, UDOVIČIĆ et al. (2007) (HAIR et al., 2010). Cluster analysis was then performed using only those variables that had loadings >±0.5 (which is consistent with guidelines provided in HAIR et al. (2010)), and were in the same principal component with nitrates. The purpose of cluster analysis was to group obser vation wells based on variables obtained from principal compo nent analysis. For those purposes cluster analysis was tested us ing Ward`s method, single linkage and complete linkage rules, while Euclidean and squared Euclidean distances were used as distance measures. Even though there are different rules of thumb that prescribe the minimal requirement of sample size for factor analysis, a minimum of 100 cases (MACCALLUM et al., 1999) and a case/variable ratio of 5:1 (BRYANT & YARNOLD, 1995;HAIR et al., 2010) rule was satisfied. Regarding usage of mini mal sample size in cluster analysis, MOOI & SARSTEDT (2011) stated that there is no general rule of thumb that provides mini mum sample size, or the relationship between the cases and number of clustering variables. In the end, all available pH and ORP data were placed on a Pourbaix diagram for nitrogen com pounds, which shows the possible stable equilibrium phases of an aqueous electrochemical system. For this purpose 3202 pairs of pH and ORP values available from groundwater chemical anal yses were used. Calculations and figure construction was done using Microsoft © Excel, Statistica 64 (version 13.1) and ArcMap 10.1, while a geocoded terrain (georeferenced orthophoto) image was obtained from the geoportal of the Croatian Geodetic Ad ministration. All maps are presented using the official coordinate system of the Republic of Croatia (HTRS96/TM).

Average values of observed variables
Average nitrate concentrations from 153 observation wells are shown in Fig. 4. They vary from 0.8 to 44.20 mgl 1 NO 3 . The highest concentrations are registered in groundwater of the urban part of the City of Zagreb (the left bank of the Sava River) and in the predominantly agricultural area between Mala Mlaka and Velika Gorica City (the right bank of the Sava River). Dissolved oxygen concentrations are shown in Fig. 5 and are divided into 3 groups (<1, 12 and >2 mgl 1 ). Dissolved oxygen concentrations were divided into 3 groups to identify observation wells with pre dominantly aerobic conditions (>2 mgl 1 ), predominantly anaero bic conditions (<1 mgl 1 ), and those that have threshold values between the two types of conditions (1-2 mgl 1 ). It can be seen that predominantly aerobic conditions prevail in the western part of the aquifer system, while more anaerobic conditions occur in the eastern part of the aquifer system and in some observation wells that are located in the vicinity of the Sava River. Observa tion wells that have dissolved oxygen concentrations between 1 and 2 mgl 1 only occur in the eastern part, near the Sava River. Average ORP values vary from 48.82 to 524.89 mV, while pH values vary from 6.98 to 7.63. EC ranges from 425.28 to 1241.4 μS/cm, where the higher values are mostly located on the left bank of the Sava River, in the urban part of the City of Zagreb. Water temperature varies on average from 11.4 to 17 °C.

Correlation
Results of the Kolmogorov-Smirnov test are shown in Table 1. Results show that ORP and temperature of groundwater are not normally distributed, while NO 3 , O 2 , pH and EC show normal distribution patterns.
Correlation coefficient matrices for r and ρ are shown in Ta ble 2. All statistically significant results are marked in red. Generally, the results show that nitrate concentrations are more related to pH, O 2 and EC, where correlation with pH is negative, while correlation with the other two parameters is positive. Positive cor relations of NO 3 with O 2 and EC seem very logical because higher O 2 concentrations should provide a more stable geo chemical

Principal component analysis
Principal component analysis generated 2 principal components explaining 76% of the total variance, where 54% was explained by the first component. Varimax raw rotated loadings are shown in Table 3, together with its communalities. PC 1 is presented with pH, O 2 , EC and NO 3 , while PC 2 is presented with groundwater temperature and ORP. It is evident that an oxygenated environ ment generates higher nitrate concentrations and EC values which results in a more acidic environment. Also, results indicate that the ORP and temperature of groundwater are not related to nitrate concentrations, which coincides with previous results. Further more, results suggest that ORP and groundwater temperature present variables that are less under human influence than the vari ables from PC 1. Due to very low communalities of ORP (<0.5) and the affiliation of the temperature of groundwater to different principal components, these two variables were excluded from the cluster analysis. Varimax rotation loadings are shown in Fig 6.

Cluster analysis
Cluster analysis was performed using three different methods, where Ward`s method generated the two most distinctive clusters ( Fig. 7 and 8). In all cases squared Euclidean distances and Eu clidean distances gave similar results. However, the squared dis tances facilitated the drawing of conclusions regarding the num ber of clusters.
The results of cluster analysis obtained using Ward`s method have been taken as the most representative one, which is probably a consequence of Ward`s method of calculation. Single and link age methods define the similarity of clusters using minimum and maximum distances between objects, while Ward`s methods maximizes the homogeneity between clusters using the sum of squares within the cluster. Evaluation of projected clusters have shown two main results. First, clusters have very distinctive con centrations of NO 3 , O 2 and EC (Table 4). Cluster 1 has generally lower values of variables, while cluster 2 has higher values. In cluster 1 the average value of O 2 from 61 observation wells is 1.77 mgl 1 , of EC is 573.31 μS/cm, and of NO 3 is 7.84 mgl 1 . In cluster 2 the average value of O 2 from 65 observation wells is 4.8 mgl 1 , of EC is 859.03 μS/cm, and of NO 3 is 25.93 mgl 1 . These results indicate that NO 3 concentrations are generally controlled by O 2 concentrations. Also, where higher concentrations of EC are re corded, higher nitrate concentrations can also be expected, which is also confirmed by correlation analysis. Secondly, if clusters are evaluated spatially (Fig. 9), observation wells in different clusters generally coincide with the average values of dissolved oxygen shown in Fig. 3. Most of the observation wells from cluster 1 are near the Sava River, in the area with the lower oxygen content. Nitrate concentrations in the Sava River have average concentra tions of 7.33 mgl 1 at the Jankomir station, 7.5 mgl 1 at the Petruševec, and 6.67 mgl 1 at Rugvica (Table 5), while the ave-  rage NO 3 concentrations in cluster 1 are similar, i.e. 7.84 mgl 1 .
Lower values of O 2 in the eastern part of the Zagreb aquifer sys tem are probably a consequence of greater aquifer depth and mix ing of different groundwaters from Holocene (more oxidative condition) and Pleistocene aquifers (more reductive condition). Also, results indicate that the Sava River has an influence on the dissolved oxygen concentrations. The reasons for this could be varied. They can be associated with too many bacteria and an excess amount of biological oxygen demand, and maybe with fer tilizer runoff from farm fields. The first reason is the more likely in this case because it can be associated with some kind of or ganic discharge, probably sewage. These results indicate that the quality of the Sava River is under very significant anthropogenic influences.

Pourbaix diagram and dominant nitrogen species
In general, all statistical analyses showed the close relationship of NO 3 , O 2 , EC and pH. ORP and groundwater temperature showed a very poor relationship with the other variables. When evaluating the Pourbaix diagram for nitrogen species and mea sured data (Fig. 9) in the Zagreb aquifer, it can be seen that all values are oriented more in a vertical than in a horizontal direc tion. This indicates constant change in oxidative and reductive conditions in the Zagreb aquifer, which is also probably the rea son why the correlation between NO 3 and ORP can`t be ob   served. In Fig. 10 it can be seen that the most frequent average ORP class is that between 150 and 200 mV, which presents the threshold class between oxidative and reductive conditions in the aquifer (DIMKIĆ et al., 2008). Also, Fig. 10 shows that about 50% of observation wells have ORP values higher than 200 mV, which suggests that in the Zagreb aquifer system moderately oxi-   dative conditions prevail. Furthermore, from Fig. 11 it can be seen that many samples fall in the area of high NH 4 + stability. Due to the very rare occurrence of NH 4 + (367 values at 49 observation wells) and NO 2 concentrations (249 values at 69 observation wells) observed in 25 years of groundwater quality monitoring of the Zagreb aquifer system (Fig. 12), it can be assumed that very rapid nitrification of NH 4 + and NO 2 to NO 3 occurs in the aquifer system. Moreover, Fig. 12 shows that concentrations of NH 4 + and NO 2 on most observation wells occur less than five times in the observed period.

CONCLUSION
Nitrates represent one of the main groups of contaminants in the Zagreb aquifer system. Their concentrations are mainly the result of anthropogenic influences. Despite this, some natural ground water quality indicators, more or less influenced by humans, may affect their stability and mobility in groundwater. Average values of NO 3 , O 2 , ORP, pH, EC and the temperature of groundwater were used to test their relationship using correlation and multi variate statistical analysis. Pearson and Spearman correlation sta tistics have produced similar results in that NO 3 is generally positively correlated with O 2 and EC, and negatively with pH. ORP values showed no to very poor correlation with all other variables, as well as with the temperature of groundwater. Cor relation analysis indicated that only O 2 , EC and pH were variables related to nitrate, which was confirmed with multivariate statis tical analysis results. PCA was used to identify principal compo nents and variables which were then used as variables in the clus ter analysis. It generated 2 PCs where only pH, O 2 , EC and NO 3 from PC 1, with factor loadings >±0.5, were used in the cluster analysis. Cluster analysis was tested with Ward`s method, single linkage and complete linkage rules, while Euclidean and squared Euclidean distances were used as distance measures. Different distance measures did not provide any significantly dif ferent results of the cluster analysis, but linkage rule methods did. Ward`s method generated two most distinctive clusters, probably because of the difference in calculation of the difference between clusters, and its results were used for the evaluation of two clus ters. Clusters mostly differ in O 2 , EC and NO 3 concentrations.
In cluster 1 lower values of NO 3 are probably a consequence of the influence of the Sava River and greater aquifer depth in the eastern part of the aquifer system, where mixing of anaerobic and aerobic water occurs. Due to the great variability of ORP in the study area and lack of NH 4 + and NO 2 species in groundwater, it can be assumed that very rapid nitrification of NH 4 + and NO 2 to NO 3 occurs in the Zagreb aquifer system, particularly in the shal low Holocene aquifer.

ACKNOWLEDGMENT
Authors would like to thank Croatian Waters for the data obtained within the project "Groundwater trend and status assessment in the Pannonian part of Croatia" and to Zagreb holding (branch ZGOS) for the data obtained from monitoring program of landfill Jakuševec.