## 2 Material and Methods
### Database Design  
In Italy the Regional Environmental Protection Agencies (ARPA) are the institutions in charge of monitoring of the quantity and quality of surface and groundwater bodies, according to European (Directive 2000/60/EC, Water Framework Directive WFD, [6]), and national directives [5; 49; 50]. The monitoring of the quantitative state is carried out by measuring the piezometric levels with reference to the mean sea level; while the monitoring of the chemical status is based on the analysis of: pollutants subject to quality standards, identified at EU level (Annex III of D.Lgs 30/09, Table 2, [49]); and element subject to threshold values, identified at national level (Annex III of D.Lgs 30/09, Table 3, [49]). The decree provides indications on minimum sampling frequencies depending on the characteristics of the water body, the base parameters, and some additional ones to be controlled. In this study we used the 2018 databases of qualitative and quantitative groundwater monitoring provided by the ARPAs of Piedmont (hereinafter referred to as PIE), Lombardy (LOM), Veneto (VEN), Friuli-Venezia Giulia (FVG), and Emilia-Romagna (EMR) through their OpenData portals. The databases were integrated with information of sampling stations (well or spring, coordinates, etc.,), related aquifer (phreatic, semi-confined, or confined), and monitoring depths. Since each agency provided data according to its own standards (coordinate system, detection limits of analytical protocols, concentration units, and monitoring time interval), a great effort was dedicated to data homogenization and management to create an operational database describing the hydrogeochemical characterization of the 2018 groundwater of the Po and Venetian-Friulian plains. The inhomogeneities in the dataset depend on the choices made by each individual ARPAs, especially as regards metals and trace elements, as established in part III D.Lgs 30/09 [49]. Table S1 in Supplementary Material shows the parameters analyzed by each Region with the relative measurement units. After the verification and cleaning processes, consisting in the eliminations of samples without clear aquifer identification or located in aquifers other than semi-confined, confined, and phreatic (i.e., colluvial and mountainous aquifers), 3671 validated samples, out of the initial 4049, were used for data analysis and map elaborations. The water samples were analyzed by different laboratories, and therefore had different detection limits for the same analyte, both at regional and in some cases even at the provincial level, as shown in Tables S1 and S2 in Supplementary Material. Another factor that caused low homogeneity of the final database concerns the different analytical pool analyzed between regions. Specific electrical conductivity at 20 \({}^{\circ}\)C (EC), hardness, pH, and temperature (T) were provided by at least 4 out of 5 agencies. For the other chemical and physical characteristics, instead, data obtained were fragmentary and incomplete over the whole considered time period. Concentrations of metals and compounds of environmental interest were not homogeneously collected by the agencies; the divergence is justified by the different operational choices adopted by the agencies, as allowed by the relevant legislation (Annex III of the D.Lgs. 30/09 [49]) as shown in Table S2 of Supplementary Material. Several authors established strategies for using values below detection limit (<D.L.) as valid data [51; 52; 53]. However, in the present study the management of these data was not done using the proposed methods, because the characteristics of the database were not in line with the conditions set by the authors. The first difference concerns the presence of many different D.L. for each analyte, and secondly data population of several analytes was predominantly composed by values <D.L. The adopted solutions varied according to the data considered.