Statistical Indicators of Groundwater Geochemical Characteristics in a Quaternary Aquifer from the Mala Mlaka Well Field Catchment Area (Zagreb, Croatia)

The Mala Mlaka well field is the largest near Zagreb and is very important for the city’s water supply. The field is located on the right Sava River bank, in the immediate vicinity of the New Zagreb urban complex. The vicinity of the neighbouring industry, agricultural land, illegal gravel pits and dumpsites that developed in abandoned gravel pits threaten the quality of the well field groundwater. The Croatian Groundwater Recording and Management Project clearly determined degradation in the groundwater quality in the Zagreb aquifer system. Continuous and occasional groundwater contamination was detected in a large number of observation wells. Statistical Indicators of Groundwater Geochemical Characteristics in a Quaternary Aquifer from the Mala Mlaka Well Field Catchment Area (Zagreb, Croatia)


INTRODUCTION
The Mala Mlaka well field is the largest near Zagreb and is very important for the city's water supply.The field is located on the right Sava River bank, in the immediate vicinity of the New Zagreb urban complex.The vicinity of the neighbouring industry, agricultural land, illegal gravel pits and dumpsites that developed in abandoned gravel pits threaten the quality of the well field groundwater.The Croatian Groundwater Recording and Management Project clearly determined degradation in the groundwater quality in the Zagreb aquifer system.Continuous and occasional groundwater contamination was detected in a large number of observation wells.

Statistical Indicators of Groundwater Geochemical Characteristics in a Quaternary Aquifer from the Mala Mlaka Well Field Catchment Area
(Zagreb, Croatia) Zoran NAKIĆ, Siniša HORVAT and Andrea BAČANI The article uses results of numerical simulation of the groundwater flow in the Zagreb aquifer system that were employed in determination of the Mala Mlaka well field catchment area under different hydrological conditions.The model results were complemented with geological and hydrogeological interpretation of the southern border of the catchment area (BAČANI & ŠPARICA, 2001).In the model, the southern border was determined to be immediately adjacent to the southern border of the third protection zone, since the number of sites on which groundwater tables are measured is small.The results of the numerical simulation were used as a background for further research into groundwater quality.The pollution sources were recorded in the catchment area (Fig. 1), primarily within the protection zone borders (Fig. 2).Results of the chemical analyses for the period 1994-2000 were provided by the Water Supply and Drainage Utility.
The analysis has revealed that: -the observation well sites sampled for chemical and bacteriological analyses are mostly located within the borders of the Mala Mlaka well field protection zones; -the sampling was generally carried out on a monthly basis.
This article reports on basic hydrogeochemical processes evolving in a shallow water-bearing layer in the Mala Mlaka well field catchment area.The hydrogeochemical facies, i.e. the main water types were determined.The groundwater quality time series were analysed to confirm significant hydrogeochemical processes in the catchment area.The analysis was used to determine the time series for significant chemical elements and compounds from the groundwater quality data for the period 1994-2000.In order to evaluate individual impacts of the pollution sources on the groundwater chemistry, the similarities and interrelationships of the observation well sites hydrochemical data were investigated.The multivariate statistical analyses were used, i.e. factor analysis, cluster analysis and multidimensional scaling (MDS) analysis.Multivariate analysis is based on carrying out simultaneous measurement of several variables on a sample.A primary goal is to determine interrelationships of the measured variables of the tested samples (BROWN, 1998).In the article, the chemical analysis data from 2000 were used.
Statistica for Windows 5.1 software package (STAT-SOFT, 1995) was used for statistical analyses.

SITE AND HYDROGEOLOGICAL RELATIONS
The Mala Mlaka well field is located in the southern part of Zagreb, on the right Sava River bank, in the immediate vicinity of an urban area.It was designed in 1938, and the first wells were put into operation in 1964 (ČAKARUN et al., 1987).The well field has 10 dug and 6 drilled wells with an average yield of about 1400 l/s.
The geological texture is characterized by younger sediments blended in the Sava alluvium, i.e. the sediments deposited by the Sava River during the last 10,000 years while it cut its course in this area 1 .These are Holocene sediments that reach depths of approximately 35 m in the Mala Mlaka area and are dominated by gravel and sand of the first and the second alluvial terrace.The youngest Holocene sediments are exceptionally important since they form the main aquifer, which is of particular importance for the water supply of the greater Zagreb.At greater depths, deeper aquifers occur in Middle and Upper Pleistocene strata, with frequent lateral and vertical alterations of gravel, sand and clay.
Generally speaking, two regular phenomena have been noticed in the Zagreb aquifer area as regards spreading of the Holocene deposits: (a) a decrease in thickness outwards from the Sava river, and (b) an increase in thickness from West to East (BORČIĆ et al., 1968).
Sediments of the second alluvial terrace consist of alternating coarsely grained gravel and sand.The petrographic composition of the pebbles is different, but carbonate rock pebbles are most frequently encountered.Lithologically, the composition of the first alluvial terrace sediments is dominated by coarsely grained gravel mixed with sand.In the Eastern and South-East areas, the sediment composition includes clay mixed with sand and silt (BAČANI & ŠPARICA, 2001).
The groundwater flow direction during high and low water tables is NW-SE (Fig. 2).An average difference between the groundwater tables during high and low tables at the observation well sites is 1.5 to 2 m.

THE MALA MLAKA WELL FIELD CATCHMENT AREA
The Mala Mlaka well field catchment area was determined within the Croatian Groundwater Recording and Management Project, along with catchment areas of another four Zagreb well fields: Zapruđe, Velika Gorica, Sašnjak and Petruševec.The catchment areas were determined in two ways: -graphically, by using data on the groundwater tables, and -from the numerical simulation results.
The catchment area boundaries were graphically determined by delineating the groundwater flow direction on the contour line maps, during low and high groundwater tables (JONES, 2001).The resulting boundaries of the Zagreb well field sites are for: After verification for these periods, the selected conceptual model was used to observe the catchment area changes depending on the Sava River water levels in 1995 (Fig. 3).

STATISTICAL ANALYSES
In the article, time series analysis and multivariate statistical analyses: factor analysis, cluster analysis and multidimensional scaling analysis (MDS) are used.
The time series analysis is used to determine time trends from the groundwater quality data.The determined trends are used to explain geochemical process- es and potential groundwater contamination risks.The analysis usually consists of determining trend statistics for individual sites.Such an approach demands an adequate continuous data series.Depending on the length of the analyzed series and their homogeneity, a specific statistical analysis is used to accept or reject the statistical trend results.
In this article, a different approach was taken that allows determination of a common trend for a group of sites inside a homogenous area.The analysis comprises several steps: (1) Determination of Spearman rank correlation coefficients for water quality and time variables per observation well.Each well is an identifiable sub-sample of a population inside of homogeneous area.The correlation coefficients are calculated statistics per sub-sample which are used to determine the sample distribution of that statistic.
(2) Determination of the hypothetical null-distribution of parameters for the observation wells group under consideration.This theoretical distribution is derived from a reference population with a given population value and can be formed by drawing an infinite number of samples of definite size from this population (FRAPPORTI et al., 1991;FRAPPORTI, 1994).
(3) Testing a sample distribution against theoretical null hypothesis distribution for equality by χ 2 test.(4) Depending on the critical value of the χ 2 test, the existence of a significant statistical trend is either accepted or rejected.
The above analysis uses Spearman rank correlation coefficient, which is a robust non-parametric equivalent of Pearson correlation coefficient.This correlation statistic is used primarily because of the fact that a time variable, and often a chemical parameter as well, are not normally distributed according to normal or Gauss distribution principles.
The theoretical null hypothesis distribution is defined by all possible combinations of ranks x i and y i for N observations (FRAPPORTI et al., 1994).In the article, a distribution approximation procedure was carried out, using the Student-t test distribution with n-2 degrees of freedom (KENDALL & STUART, 1961;FRAPPORTI et al., 1994): Advantages of this approach to trend determination include: -the analysis can include the maximum possible number of observation well analysis, regardless of the statistical significance of the trend results obtained for individual locations, -the statistical trend distribution study enables the use of low-reliability trends for certain sites for high-reliability determination of trends for a group of observation wells, and -defining multiparametric trends of a group enables significant geochemical processes and pollution/contamination events to be determined over the entire homogenous area.
Factor analysis attempts to simultaneously explain relationships between several variables in order to find more simple relations that ensure an insight into a hidden data structure (SUK & KANG-KUN, 1999).These simple relationships are expressed through a new set of variables called factors, which retain the maximum amount of information and show the interrelationships of the variables.
The most important parameters obtained by factor analysis are communality, factor loading, eigenvalue and the sum of the explained variance.Communality is part of the variance explained by the common factors; factor loadings describe the extent of correlation between the factors and the original variables; eigenvalue is an quantity of variance explained by a particular factor, and the sum of the explained variance is the total quantity of variance obtained from the resulting eigenvalues (BROWN, 1998).
Here, the data on all the variables are standardized as mean value 0 and variance 1 to avoid problems with different variances of variables, which might affect determination of the factor loadings.
The factor scores from factor analysis were used as input data for cluster analysis.The procedure had already been carried out on an example given by SUK & KANG- KUN (1999), and the objective is to reduce those variables and factors that are not significant for explanation of the hydrogeochemical processes.
Cluster analysis is a method used for grouping individual variables or samples so that, depending on the measured characteristics, the distance between individual groups is maximized.Cluster analysis differs from other classification methods in as much as the number and characteristics of individual groups are not known before the analysis, so they are derived directly from the data (BROWN, 1998).
In this article, the hierarchical tree clustering method is used to produce a graphical representation of individual groups using dendograms.Specifically, the Ward method of hierarchical tree clustering is used based on variance analysis to assess the distance between the clusters (WARD, 1963;STATSOFT, 1995).
Multidimensional scaling (MDS) is a form of analysis in which dimensionality is reduced, and it is based on the distance between the points.The result is expressed as a reduced number of dimensions and used to present the data in a way which causes minimum distortion of the original data (BROWN, 1998).The input data in MDS may be matrices of similarity or dissimilarity, and the correlation matrices.The input data used in the article are the data from the distance matrix, obtained by transposition of the nxp matrix of factor scores obtained by factor analysis into a nxn matrix.The procedure is intended to use factor analysis to explain the structure of the observed variables and MDS analysis to explain the differences between individual samples.

RESULTS AND DISCUSSION
The hydrogeochemical facies analysis has shown that the groundwater in the Mala Mlaka well field catchment area mostly belongs to a Ca-Mg-HCO 3 type of water (Fig. 4).The analysis has encompassed the observation well sites with filters installed in the first and second water bearing layer.
Inorganic quality indicators were used in analysis of the time series of the groundwater quality data for the Mala Mlaka well field catchment area.From examination of the data, indicators were determined which were monitored in most observation well sites in 2000.These indicators are: dissolved oxygen, nitrates, chlorides, sulphates, dissolved carbon dioxide, sodium, potassium, iron, manganese and lead.
The trend analysis determined significantly positive trends for chlorides, sodium and potassium and significantly negative trends for dissolved oxygen, nitrates and sulphates (Table 1).Results of χ 2 test also indicate a positive trend for manganese.However, the review of the available data has shown highly pronounced changes in the manganese concentration distribution before and after 1998.This is probably due to the change in heavy metal testing methods from 1998, when atomic absorption spectrophotometry (AAS) was replaced by ICP -Inductively Coupled Plasma-Mass spectrometry.
Agricultural activities in the Mala Mlaka well field catchment area are intensive, and the groundwater con-tamination sources are numerous (Fig. 2).Relatively high concentrations of sodium, potassium, nitrates and sulphates in the groundwater (Fig. 5) confirm the impact of agriculture on the groundwater chemistry.
Reduction in the content of inorganic components highly sensitive to changes in oxidation/reduction conditions in the aquifer has been observed (Table 1).A possible explanation lies in the fact that when oxygen-rich water infiltrates water bearing layers rich in organic matter, the oxygen content is reduced (APPELO & POST-MA, 1994).In low-oxygen or an oxygen-free environment (reduction conditions), biochemical reduction of NO 3 -and SO 4 2-occurs.NO 3 -reduces into NO 2 -and NH 3 , while SO 4 2-reduces into HS -and H 2 S. The time analysis results show a significant positive trend for sodium and potassium within the entire catchment area (Table 1, Fig. 5).Under some circumstances, when the groundwater gradients are high due to pumping of water from the wells, a vertical flow of groundwater is induced through clay and silty intercalations in the aquifer (HEM, 1985), which may cause release of sodium and potassium ions from the crystalline screen, even if irreversible adsorption of potassium ions on the clay minerals is assumed.However, as an increase in sodium and potassium ions, and consequently in chloride ions, has been observed in the groundwater over the entire catchment area, it is more probable that leaking sewage, mineral and organic fertilizers and use of salt to prevent ice formation on roads cause positive trends in sodium, potassium and chlorides in the groundwater.
Similarities and correlations between the hydrochemical data from the observation well sites were checked using the data from two periods: June 2000, during low groundwater table, and October 2000, during the lowest groundwater table in 2000.Data on chemical elements and compounds (Table 2) are standardized as the mean value 0 and variance 1.
The correlation matrices (Table 3) show a sufficiently good correlation and no presence of singularity and multicolinearity between the observed variables.
The factor analysis was performed on 10 variables.The analysis included nitrates, sulphates, chlorides and heavy metals as probable indicators of human impact.The results in Tables 4 and 5 refer to the four most important factors that explain over 77% of variance for the hydrochemical data collected in June 2000 and 82% of variance for hydrochemical data from October 2000.
The cluster analysis results for June 2000 are shown in Fig. 6.The results show that the observation wells are grouped into four clusters.The MDS analysis results (Fig. 7) indicate that only an interrelationship between the observation wells from clusters one and three could be partly reproduced in the MDS analysis.Both analyses give similar results for the observation wells with minimum linkage distance.
The cluster analysis results for October 2000 are shown in Fig. 8.The results indicate grouping of the  observation wells in five clusters.Unlike the June 2000 data, the MDS analysis results (Fig. 9) show a very high similarity compared to the cluster analysis.
Comparison of the cluster and MDS analysis results for June and October 2000 shows an interrelationship between the observation wells in each of these two periods.
PZO-1 and PZO-8 observation wells are grouped in cluster analysis.MDS analysis shows that they depart in June, but there is some similarity so they might be grouped together.In October, the MDS analysis results confirm grouping of these observation wells.Figure 2 shows that these two sites are separated by a larger distance, but the observation well PZO-1 is located down-  stream from the group of the pollution sources which might affect the water chemistry at this site.It may be assumed that similarity in the groundwater chemistry at sites PZO-1 and PZO-8 is probably a consequence of local hydrodynamic and hydrogeochemical conditions in the shallow water bearing layer.
Observation wells MM-321, MM-324 and MM-332 are grouped in cluster analysis during June and October 2000, and the MDS analysis results show their grouping in October 2000.Grouping of these observation wells, which are set at small distances in the second protection zone of the Mala Mlaka well field (Fig. 2), points to the impact of numerous pollution sources (illegal dump sites and gravel pits) located downstream from the observation wells at the sites under consideration.
Observation wells PZO-12 and MM-333 are grouped in cluster analysis for June and October 2000.The MDS analysis results for June 2000 are separated, but there is a certain similarity so they might be grouped together.In October, the MDS analysis results confirmed grouping of these observation wells.The distance between these observation wells is small along the Sava-Odra canal (Fig. 2).It might be presumed that similarity of the groundwater chemistry in the observation wells PZO-12 and MM-333 is probably a consequence of local hydrodynamic and hydrogeochemical conditions in the shallow aquifer.
Observation wells MM-310, MM-311 and PZO-2 are grouped in both analyses, for June and October 2000.The distance between the observation wells MM-310 and MM-311 is small, and they are located downstream from the illegal gravel pits and dump sites (Fig. 2).The results of cluster analysis for October 2000 show the extraordinarily small linkage distance of these sites.There is a high probability that the illegal gravel pits and dump sites will have considerable impact on these observation wells.Observation well PZO-2 is the largest distance from the dump sites and gravel pits, but upstream from this site, a larger group of pollution sources is recorded so the observation well is probably exposed to their impact.
Observation wells MM-32 and MM-319 are grouped in both analyses, in June and October 2000.The distance between these observation wells is larger, and they are located near the first protection zone: MM-32 to the east and MM-319 to the west from the first protection  zone border (Fig. 2).Upstream from MM-32 is a site that accommodates workshops and upstream from MM-319 an illegal gravel pit filled with waste is located.

CONCLUSION
Most of the activities currently carried out in order to protect the groundwater of the Zagreb aquifer system have been focused on the areas immediately adjacent to the well fields or within the well field protection zones.The fact that the well field is only a point where groundwater is tapped for water supply was neglected, and the entire catchment area should be controlled.Efficient protection of the groundwater in the well fields is only possible if the catchment area boundaries were determined from the information on the underground hydrodynamic conditions, and geological and hydrogeological relationships within the entire system.In line with this approach, catchment areas were determined for five Zagreb well fields (Zapruđe, Mala Mlaka, Velika Gorica, Sašnjak and Petruševec) within The Croatian Groundwater Recording and Management Project.The results of the numerical simulation were used as background for further research into groundwater quality and protection.The pollution sources were recorded within the catchment area boundaries, and the chemical analysis data were collected and systematically presented for the period 1994-2000.
Standard methods were used for determination of the hydrogeochemical facies, and it was found that the groundwater in the Mala Mlaka well field catchment area belongs to a Ca-Mg-HCO 3 type of water.
The analysis of time trends has revealed geochemical processes in the groundwater of the Quaternary aquifer, and grouping of the observation wells based on the groundwater chemistry in individual sites confirmed existence of human impact on the groundwater quality.The following was detected: -Relatively high concentrations of sodium, potassium, nitrates and sulphates in the groundwater confirm the impact of agricultural activities on the groundwater chemistry.Positive trends of sodium, potassium and chlorides are a consequence of the application of mineral and organic fertilizers on agricultural land, leaking sewage and use of salt to prevent ice formation on roads during the winter.
-Reduction in the content of inorganic components, nitrates and sulphates, which are sensitive to changes in oxidation/reduction conditions in the water bearing layer is caused by infiltration of oxygen-rich water into the aquifers rich in organic matter, which results in reduction of oxygen due to its consumption for organic matter oxidation.In environments, which are often reducing, nitrate contents decrease by their reduction into nitrites and, depending on the oxygen content, into ammonium.The situation is similar with sulphates which are reduced to sulphite and finally into hydrogen sulphide where oxygen levels are reduced.Both transformation processes are caused by microbial activities.
-Registered pollution sources in the area, illegal dump sites in particular, caused degradation in the groundwater quality in some sites.This impact is traceable even several hundred metres from the registered pollution sources.

Fig. 1
Fig.1The catchment area and the sanitary protection zones of the Mala Mlaka well field.

Fig. 2
Fig.2Groundwater sampling sites and impacts in the Mala Mlaka well field catchment area.
(a) high water tables, on March 9, 1995, (b) low water tables, on April 2, 1998.The catchment areas of the Zagreb well fields were numerically determined on the basis of the results obtained by numerical simulation conducted using Visual Modflow software (by Waterloo Hydrogeologic Inc.).

Fig. 3
Fig.3The Mala Mlaka well field catchment area at the time of low and high groundwater tables.

Fig. 4
Fig. 4 The Piper diagram of the groundwater ionic composition in the Mala Mlaka well field catchment area.

Fig. 6
Fig. 6 Results of cluster analysis shown by a dendogram of June 2000 data.Fig. 7 Results of multidimensional scaling analysis shown by the two most important dimensions (June 2000).

Fig. 8
Fig. 8 Results of cluster analysis shown by a dendogram of October 2000 data.

Fig. 9
Fig. 9 Results of multidimensional scaling analysis shown by the two most important dimensions (October 2000).

Table 1 -
continued on the next page.

Table 2
Chemical analysis data for (a) June 2000, and (b) October 2000.

Table 3
Correlation matrix of 10 variables computed from chemical analysis data for (a) June 2000, and (b) October 2000.
Table 4 Eigenvalues derived by factor analysis, percentages of total variance, cumulative eigenvalues and cumulative percentages explained by factors: (a) June 2000, (b) October 2000.