UNED Research Journal (e-ISSN 1659-441X), Vol. 13(1): e3294, June, 2021

Spatial variation of physicochemical parameters in a constructed wetland for wastewater treatment: An example of the use of the R programming language

Junior Pastor Pérez-Molina¹, Carola Scholz², Roy Pérez-Salazar³, Carolina Alfaro-Chinchilla³, Ana Abarca Méndez¹, Leandro Araya Leitón¹, Jeslyn Carranza Chaves¹, Addy Echevarría Figueroa¹, Mariana Elizondo Blanco¹, Rachel Ardón Rivera¹, Sofía Flores Aguilar¹ & Catalina Solís Calderón⁴

1. Laboratorio de Ecología Funcional y Ecosistemas Tropicales (LEFET), Escuela de Ciencias Biológicas, Universidad Nacional, Costa Rica, 86-3000; junior.perez.molina@una.cr; analaura.abarcamndez@gmail.com; arayita88@gmail.com; jeslynfcch@gmail.com; addyef10@gmail.com; marianaelizondo60@gmail.com; rachelardon@hotmail.com; sofia.flores.aguilar@est.una.ac.cr

2. Laboratorio de Fitotecnología (LAFITOTEC), Escuela de Ciencias Biológicas, Universidad Nacional, Costa Rica, 86-3000; carola.scholz@una.cr

3. Laboratorio de Gestión de Desechos y Aguas Residuales (LAGEDE), Escuela de Química, Universidad Nacional, Costa Rica, 86-3000; roy.perez.salazar@una.cr,; carolina.alfaro.chinchilla@una.cr;

4. Laboratorio Nacional de Aguas, Área de Microbiología, Instituto Costarricense de Acueductos y Alcantarillados (AyA), Costa Rica, 30306; csoliscalderon@gmail.com

Received 29-X-2020 – Corrected 18-XI-2020 – Accepted 02-XII-2020

DOI: https://doi.org/10.22458/urj.v13i1.3294

ABSTRACT. Introduction: The implementation of wastewater treatment systems such as constructed wetlands has a growing interest in the last decade due to its low cost and high effectiveness in treating industrial and residential wastewater. Objective: To evaluate the spatial variation of physicochemical parameters in a constructed wetland system of sub-superficial flow of Pennisetum alopecuroides (Pennisetum) and a Control (unplanted). The purpose is to provide an analysis of spatial dynamic of physicochemical parameters using R programming language. Methods: Each of the cells (Pennisetum and Control) had 12 piezometers, organized in three columns and four rows with a separation distance of 3,25m and 4,35m, respectively. The turbidity, biochemical oxygen demand (BOD), chemical oxygen demand (COD), total Kjeldahl nitrogen (TKN), ammoniacal nitrogen (N-NH4), organic nitrogen (N-org.) and phosphorous (P-PO4^-3) were measured in water under in-flow and out-flow of both conditions Control and Pennisetum (n= 8). Additionally, the oxidation-reduction potential (ORP), dissolved oxygen (DO), conductivity, pH and water temperature, were measured (n= 167) in the piezometers. Results: No statistically significant differences between cells for TKN, N-NH4, conductivity, turbidity, BOD, and COD were found; but both Control and Pennisetum cells showed a significant reduction in these parameters (P<0,05). Overall, TKN and N-NH4 removal were from 65,8 to 84,1% and 67,5 to 90,8%, respectively; and decrease in turbidity, conductivity, BOD, and COD, were between 95,1-95,4%; 15-22,4%; 65,2-77,9% and 57,4-60,3% respectively. Both cells showed ORP increasing gradient along the water-flow direction, contrary to conductivity (p<0,05). However, OD, pH and temperature were inconsistent in the direction of the water flow in both cells. Conclusions: Pennisetum demonstrated pollutant removal efficiency, but presented results similar to the control cells, therefore, remains unclear if it is a superior option or not. Spatial variation analysis did not reflect any obstruction of flow along the CWs; but some preferential flow paths can be distinguished. An open-source repository of R was provided. Keywords: Pennisetum alopecuroides, Artificial Wetland System of subsuperficial Flow, oxidation-reduction potential, dissolved oxygen, conductivity.

RESUMEN. “Variación espacial de los parámetros fisicoquímicos en un humedal artificial para el tratamiento de aguas residuales: un ejemplo del uso del lenguaje de programación R”. Introducción: La implementación de sistemas de tratamiento de aguas residuales como los humedales artificiales, ha tenido gran interés en la última década por su bajo costo y efectividad en el tratamiento de aguas residuales. Objetivo: Evaluar la variación espacial de algunos parámetros fisicoquímicos en un sistema de humedal construido de flujo sub-superficial de Pennisetum alopecuroides (Pennisetum) y un Control (no plantado). Se proporciona un ejemplo de análisis de dinámica espacial de parámetros fisicoquímicos mediante el uso del lenguaje de programación R. Metodología: Cada una de las celdas (Pennisetum y Control) tenía 12 piezómetros, organizados en tres columnas y cuatro filas con una separación de 3,25 m y 4,35 m, respectivamente. Se midieron en agua del flujo de entrada y salida de ambas condiciones Control y Pennisetum (n = 8): la turbidez, la demanda bioquímica de oxígeno (DBO), la demanda química de oxígeno (DQO), el nitrógeno Kjeldahl total (TKN), el nitrógeno amoniacal (N-NH4), el nitrógeno orgánico (N-org.) y el fósforo (P-PO4^-3). Además, se midieron el potencial de oxidación-reducción (ORP), oxígeno disuelto (DO), conductividad, pH y temperatura del agua (n = 167) en los piezómetros. Resultados: No se hallaron diferencias estadísticamente significativas entre celdas para TKN, N-NH4, conductividad, turbidez, BOD y COD (p>0,05); pero ambas celdas presentaron una reducción significativa en estos parámetros (p<0,05). La eliminación de TKN y N-NH4 fue de 65,8 a 84,1% y de 67,5 a 90,8%, respectivamente. La disminución de turbidez, conductividad, BOD y COD fue entre 95,1-95,4%; 15-22,4%; 65,2-77,9% y 57,4-60,3% respectivamente. Ambas celdas mostraron un gradiente creciente de ORP a lo largo de la dirección del flujo de agua, contrario a la conductividad (p<0,05). Sin embargo, OD, pH y temperatura fueron inconsistentes en la dirección del flujo de agua en ambas. Conclusiones: Pennisetum demostró eficacia en la remoción de contaminantes, pero presentó resultados similares a Control. Hay duda si a largo plazo puede demostrarse como una opción superior. El análisis de variación espacial no reflejó ninguna obstrucción del flujo a lo largo de los CWs; pero se pueden distinguir rutas de flujo preferenciales. Se proporcionó un repositorio de código abierto de R. Palabras clave: Pennisetum alopecuroides, sistema de humedales artificiales de flujo sub-superficial, potencial de oxidación-reducción, oxígeno disuelto, conductividad.

The lack of sanitation has become a worldwide problem. This has exacerbated with the generation of increasing amounts of wastewater from industries and the residential areas (van Puijenbroek, Beusen, & Bouwman, 2019; Ahmad, Majhi, Kothari, & Singh, 2020; Gomes da Silva & Gouveia, 2020). In Costa Rica, the average of wastewater is estimated at 966 455 m³ d^-1 (the average per person is about 0.2 m³ d^-1 (Instituto Nacional de Estadística y Censos [INEC] 2015)). Only 15% of the country's wastewater has some type of treatment (Acueductos y Alcantarillados [AyA], Ministerio de Ambiente y Energía [MINAE] y Ministerio de Salud [MS], 2016). This information is calculated from industrial companies that register the produced wastewater and companies reporting to have some treatment system; which could indicate that wastewater problem sanitation within the country is more serious than is estimated.

At present, there is a progressive explosion in the development of wastewater treatment systems (Xu & Cui, 2019). Among the wastewater treatment technologies used in Costa Rica are aerobic treatment systems such as activated extended aeration sludge and conventional activated sludge. The anaerobic systems used are septic tanks, drainage, up-flow anaerobic sludge blanket digestion, pre-conditioning systems with submarine emissary, and systems by facultative ponds (Ruiz, 2012). Although these treatment technologies are very efficient, their application within the country will largely depend on the balance between their minimum acceptable result and the minimum cost of construction and maintenance.

Constructed wetlands (CW) for wastewater treatment is an emerging ecological technology of low investment, low energy consumption, and convenient management (Tang, Huang, Scholz, & Li, 2011; Cui, Ouyang, Gu, Yang, & Xu, 2013; López-Ocaña et al., 2018). During the 1980s and 1990s, CW with horizontal subsurface flow (CW-HSF) was the most common type of CW used for wastewater treatment, due to the elimination of organic compounds and suspended solids, using filtration beds with anoxic and anaerobic conditions (Vymazal, 2019). In recent years, researchers have recognized the economic and environmental benefits for treating both industrial and residential wastewater effectively (He, Peng, Hua, Zhao, & Xiao, 2018).

CW-HSF with emerging macrophyte plants performs effectively in removing conventional pollutants (e.g.: nitrogen, phosphorus, and organic compounds) due to their excellent sediment filtration (López-Ocaña et al., 2018). Similar systems that have operated for more than 20 years, show removal levels higher than 91,7%; 82,9% and 88,3% for biochemical and chemical oxygen demand, and total suspended solids, respectively during the last five years (Vymazal, 2019). However, despite the diversity of researches conducted in the last decade, the strategies for the spatial monitoring of these wastewater treatment systems (WTS) for proper functioning are still poorly studied. Notwithstanding CWs' main drawbacks are substrate media clogging and high area requirements (Tatoulis, Akratos, Tekerlekopoulou, Vayenas, & Stefanakis, 2017).

Overall, an analysis of spatial variation of total suspended solids (TSS) within CW-HSF could indicate the location of obstruction on the flow and/or possible reduction of the removal. However, in an extensive CW-HSF the measurement of TSS could be laborious and the sample size could be unmanageable. TSS can be estimated from turbidity measurements (Grayson, Finlayson, Gippel, & Hart, 1996; Daphne, Utomo, & Kenneth, 2011); nevertheless, this requires a linear regression model fit, which must be recalibrated at each sample location and time point (Rügner, Schwientek, Beckingham, Kuch, & Grathwohl, 2013). Therefore, we presume that the measurement of physicochemical parameters such as water electrical conductivity could be a promising indicator. Because this is strictly linked to turbidity and total dissolved solids. According to Villa, Fölster, & Kyllmar (2019), the site turbidity and conductivity-total phosphorus relationship was significant at 78% of sites (n= 84; mean r² = 0,62). Other physicochemical parameters such as dissolved oxygen and oxidation-reduction potential could also be tested. It is important to note that all above physicochemical parameters are easy and quick to measure with low-cost sensors (e.g., Sensor-Vernier-LabQuest).

This research evaluated the spatial variation of some conventional physicochemical parameters on CW-HSF (1,5 years) of cultivated with Pennisetum alopecuroides and control (not plants). To the best of our knowledge, this is the first study conducted on CW using spatial variability analysis with R programming language. The main R package used was "fields"; this allows to adjust of curves, surfaces, and functions with an emphasis on splines, to later be applied to the interpolation of spatial data and spatial statistics (Nychka, Furrer, Paige, & Sain, 2017). The scope of this research was to provide an example of applications and recommendations for spatial variation analysis using R open-source tools.

MATERIALS AND METHODS

Study area: The research was conducted between May 2018 and October 2019 in a CW-HSF at the Escuela de Medicina Veterinaria located on the Benjamín Núñez Campus, Universidad Nacional [UNA] Costa Rica (9°58′37″N - 84°07′37″W). CW-HSF started to operate on April 2018. According to Instituto Metereológico Nacional [IMN] (2020), average monthly temperature, precipitation, and relative humidity of air during the study period were between 18,5-28°C, 0-186,5mm, and 79-97% (max.-min.), respectively. CW-HSF was built in an area of ~ 800m² (each cell of ~ 200m²) and receives water from toilets, student dining, and research and teaching laboratories of the Escuela de Medicina Veterinaria – UNA. Therefore, there is an important presence of chemical substances (organic and inorganic) from the mentioned facilities. The system in this research consists of a set of grilles, an Inmhoff tank, four cells of constructed wetlands with a horizontal subsurface flow, a sludge drying bed, and a disinfection unit. It is important to note that only two cells out of a total of four cells were used for this study.

All systems are operating in parallel, i.e., each of the wetlands receives the same amount of wastewater with a total hydraulic capacity of 20m³d^-1 flow for each cell. A cell of CW-HSF is formed by emerging macrophyte plants (Pennisetum alopecuroides, hereafter referred to as Pennisetum). This has a filter filled with stone, where the plants are sown to, crossing it horizontally, through which the wastewater flows from its entry through the wetlands and to the effluent collection (see details in Fig. 1). CW has a total of 12 piezometers per cell, evenly distributed in four rows and three columns. Each piezometer has a depth of 60 cm and with lateral perforations at different depths. The Control cells simply consisted of an unplanted system with the same structural characteristics mentioned above.

Fig. 1. Scheme of CW-HSF for wastewater treatment (A: Pennisetum and Control). Images of the substrate and piezometers in both systems (B-C: Pennisetum; and D-E: Control).

Sampling of physicochemical parameters: Three sampling sites were established (Fig. 1A): in-flow and two out-flow (Control and Pennisetum). The physicochemical parameters measured were (n= 8): turbidity (nephelometric turbidity units - NTU), biochemical oxygen demand (BOD, mg L^-1), chemical oxygen demand (COD, mg L^-1), total Kjeldahl nitrogen (TKN, mg L^-1), ammoniacal nitrogen (N-NH₄, mg L^-1), organic nitrogen (N-org., mg L^-1), and phosphorous (P-PO₄^-3, mg L^-1). All the physicochemical parameters were measured strictly following the established by Standard Methods for the Examination of Water and Wastewater (Rice, Baird, & Eaton, 2017).

The electrical conductivity, oxidation-reduction potential (ORP), dissolved oxygen (DO, accuracy ± 0,2 mg L^-1), pH (accuracy ± 0,2 pH) and temperature (accuracy ± 0,02 °C) of water, were measured in each of the piezometers (Control and Pennisetum, 24 piezometers in total). Likewise, measurements at water in-flows and out-flows of the two systems (Control and Pennisetum), were made. All previous these samplings were conducted for three weeks (three days, between August and October 2019) with a frequency hourly between 9:00h and 14:00h (n= 167). Each physicochemical parameter was measured by means of an independent sensor connected and distributed in two portable consoles of the Vernier type (LabQuest^®2, Version 2.8.5, Vernier Software & Technology, Beaverton, OR, USA). Before each sampling, all sensors were calibrated using their respective manufacturer's protocol.

Data analysis: All variables were evaluated using analysis of variance (ANOVA) or Kruskal-Wallis test. Fisher's Least significant difference (LSD) tests were applied to compare the means between the systems (water in- and out-flows of systems). Spearman (r_ρ) correlations were made between all the pairs physicochemical parameters per system and together. Multivariate associations among all physicochemical parameters were analyzed using a principal component analysis (PCA). Piezometers that produced similar responses were clustered using a multivariate technique of grouping analysis according to the method of Tocher, which is based on Euclidean average distances (Rencher, 2005). PCA was performed using the R package “factorextra” (Kassambara & Mundt, 2020).

Kruskal-Wallis tests were carried out with Bonferoni tests to compare between the rows of the piezometers positions for all the physicochemical parameters for each system (Control or Pennisetum), this to verify the change of the parameters along the flow of water within the systems. All statistical assumptions were checked. Likewise, bilinear interpolators between the piezometers were done, to obtain a distribution map of the levels of each physicochemical parameter to understand the dynamics within the systems to due to the direction of the water flow (see “fields” R package; Nychka et al. (2017)). All statistical analyses were performed using R programming language, version 3.6.1 (RCoreTeam, 2020) with a significance level of α=0,05. Finally, a GitHub repository was made to the access to R script, database and results of statistical analyses (see DOI: http://doi.org/10.5281/zenodo.3726063).

RESULTS

CW-HSF efficiency: Water out-flow of Pennisetum and Control showed a significant reduction with respect to water in-flow, for turbidity, BOD, COD, TKN, and N-NH₄ (P<0,05, Table 1). No significant difference was found between water in- and out-flow for N-org. and P-PO₄^-3 (P>0,05), a high variability was registered in these parameters. In addition, statistically significant differences were found between the sampling dates for conductivity, DO and pH, and the time of measurement for ORP, pH, and temperature (P<0,05, Table 2). These parameters in the time-course did not followed a determined pattern concerning the time of day except for temperature (Fig. 2). Nevertheless, statistically significant differences between water in- and out-flows cells (Control and Pennisetum) for conductivity, ORP and DO parameters, were found (P<0,05; Table 2).

Both systems showed a significant reduction in electrical conductivity and water turbidity (Table 3). Also, significant removals between 65,8% and 84,1% of TKN were found for the Control and Pennisetum, respectively. Both systems reported significant reductions in both BOD and COD (P<0,05). Control evidenced an increase in DO (20%, P<0,05), contrary to Pennisetum (1,1%). In contrast, ORP showed an increase in Pennisetum (114,8%, P<0,05).

The interplay between physicochemical parameters: Conductivity of water was correlated negatively with all independent parameters of the system and for all the data together, except for temperature and ORP in Pennisetum and Control, respectively (Table 4). The variability explained by the first two axes of the principal component (PC) analysis for all combinations of physicochemical parameters was 67,3% (Fig. 3A), which, it was not possible to separate the systems (Control vs Pennisetum). A more detailed analysis by distance in the direction of the out-flow (row of piezometers position) showed a variability explained by the first two axes of PC between 65,6% and 71,7% (Fig. 3B-E), additionally, the separation between the systems from the measured parameters, was not possible.

TABLE 1

Physicochemical parameters comparison in water in- and out-flows (Control and Pennisetum) in the CW-HSF, using one-way ANOVA

Parameters (abbreviation, unit)	In-flow	Out-flow		R2	F or KW (†)	P
		Control	Pennisetum
Turbidity (NTU)	32,8 ± 4,3 a	1,6 ± 0,4 b	1,5 ± 0,4 b	---	15,4 †	***
Biochemical oxygen demand (BOD, mg L-1)	68,7 ± 6,7 a	23,9 ± 5,1 b	15,2 ± 4,2 b	0,73	27,9	***
Chemical oxygen demand (COD, mg L-1)	111,2 ± 19,6 a	44,1 ± 6,9 b	47,4 ± 6,6 b	---	13,8 †	**
Total Kjeldahl nitrogen (TKN, mg L-1)	19,3 ± 2,6 a	6,6 ± 1,4 b	3,1 ± 0,7 b	0,69	23,4	***
Ammoniacal nitrogen (N-NH4, mg L-1)	16,0 ± 2,7 a	5,2 ± 1,3 b	1,5 ± 0,7 b	0,63	17,5	***
Organic nitrogen (N-org., mg L-1)	3,3 ± 0,7	1,4 ± 0,7	1,6 ± 0,4	0,21	2,9	n.s.
Phosphorous (P-PO4-3, mg L-1)	4,6 ± 2,3	1,4 ± 0,5	1,7 ± 0,8	---	1,7 †	n.s.

Mean ± Standard Error; R2: coefficient of determination; F: Fisher value; KW (†): Kruskal-Wallis test; P: probability; Equal letters indicate no statistically significant difference between water out-flows and in-flows (Control and Pennisetum; LSD, P>0,05); n.s.: non-significative; *: P<0,05; **: P<0,01; ***: P<0,001.

TABLE 2

Physicochemical parameters comparison in water in- and out-flows (Control and Pennisetum) in the CW-HSF, using three-way ANOVA (factors: date, time, and system)

Parameters (abbreviation, unit)	Date	Time	Systems	In-flows	Out-flows		R2	F	P
					Control	Pennisetum
Conductivity (µS cm-1)	***	n.s.	***	353,4 ± 11,9 b	274,3 ± 8,3 a	300,3 ± 10,4 a	0,61	6,3	***
Oxidation-reduction potential (ORP, mV)	n.s.	***	**	-12,8 ± 31,0 b	57,1 ± 21,2 ab	75,6 ± 17,0 a	0,64	6,7	***
Dissolved oxygen (DO, mg L-1)	**	n.s.	***	4,3 ± 0,2 b	5,4 ± 0,2 a	4,7 ± 0,2 b	0,43	3,6	**
pH	**	*	n.s.	7,1 ± 0,1	7,3 ± 0,1	7,4 ± 0,1	0,42	3,5	**
Temperature (°C)	n.s.	***	n.s.	24,5 ± 0,2	24,8 ± 0,5	24,7 ± 0,4	0,61	6,3	***

Mean ± Standard Error; R2: coefficient of determination; F: Fisher value; P: probability; Equal letters indicate no statistically significant difference between water out-flows and in-flows (Control and Pennisetum; LSD, P> 0,05); n.s.: non-significative; *: P<0,05; **: P<0,01; ***: P<0,001.

Fig. 2. Time-course of ORP (A), pH (B), and temperature (C) for water in- and out-flows of two systems (Control and Pennisetum). Standard Error (SE). Statistical analysis sees Table 2.

TABLE 3

Statistically significant difference and percentage`s removal (LSD, P<0,05) of physicochemical parameters between water out-flows and in-flows (Control and Pennisetum) in the CW-HSF

Parameters (abbreviation, unit)	Control		Pennisetum
	Difference	Removal (%)	Difference	Removal (%)
Turbidity (NTU)	-31,2	95,1	-31,3	95,4
Conductivity (µS cm-1)	-79,1	22,4	-53,1	15,0
Nitrogen Total Kjeldahl (NTK, mg L-1)	-12,7	65,8	-16,3	84,1
Ammoniacal nitrogen (N-NH4, mg L-1)	-10,8	67,5	-14,6	90,8
Biochemical oxygen demand (BOD, mg L-1)	-44,8	65,2	-53,5	77,9
Chemical oxygen demand (COD, mg L-1)	-67,1	60,3	-63,8	57,4

Spatial variation of physicochemical parameters: Conductivity and ORP showed statistically significant differences along the direction of water flow in both systems (KW> 5,07; p<0,05; see Fig. 4-5). The conductivity showed a gradient reduction along the flow water direction, contrary to ORP. Nevertheless, DO, pH, and temperature were inconsistent along the direction of the water flow, i.e.: DO in the systems were irregular along the flow, however, it was noted that the Control showed a slight increase in DO at the end of the system (out-flows; Fig. 6); pH and temperature showed non-significant changes throughout the systems with average difference of 0,2 pH and 1°C (see figures of temperature and pH in repository, DOI: http://doi.org/10.5281/zenodo.3726063).

TABLE 4

Spearman's (r_ρ) correlation coefficient matrix between pairs of the physicochemical parameters at the system level (C: Control; or P: Pennisetum; bottom diagonal) and both systems together (top diagonal)

Parameters	pH	DO	Temperature	ORP	Conductivity
pH		0,07 ^n.s^.	-0,33 ^***	-0,10 ^n.s^.	-0,31^***
DO	C= 0,10 ^n.s^. P= 0,07 ^n.s^.		0,18^***	0,17^***	-0,31^***
Temperature	C= 0,02 ^n.s^. P= -0,49^***	C= 0,36^*** P= -0,03 ^n.s^.		0,47^***	-0,27^***
ORP	C= 0,15^* P= -0,19^*	C= 0,18^* P= 0,09 ^n.s^.	C= 0,21^** P= 0,12 ^n.s^.		-0,39^***
Conductivity	C= -0,46^*** P= -0,38^***	C= -0,31^*** P= -0,27^***	C= -0,14^* P= 0,04 ^n.s^.	C= -0,15 ^n.s^. P= -0,29^***

n.s.: non-significative; *: p<0,05; **: p<0,01; ***: p<0,001.

Fig. 3. The first two axes of a principal component (PC) analysis for all combinations of physicochemical parameters for all system (A) and by each row of the position of the piezometers (B, C, D, and E).

Fig. 4. Spatial variation of electrical conductivity of water within the systems (Control and Pennisetum) based on bilinear interpolations between piezometers. Colour gradient and contour lines indicate parameter intensity from low (blue) to high (red), and white crosses indicate the position of piezometers. Boxplots show the comparison between rows and columns of the piezometers position where box marks Q1 and Q3, the black line is median (Q2), lines shown maximum and minimum values, and circles are values outliers with three times greater than the mean. Subplot with colour gradient indicate flow direction in y axis. KW: Krukall-Wallis; p: probability; equal letters indicate no statistically significant difference between the position of the piezometers in rows (Bonferroni, p>0,05).

Fig. 5. Spatial variation of ORP within the systems (Control and Pennisetum) based on bilinear interpolations between the piezometers. See legend explanation in Fig. 4.

Fig. 6. Spatial variation of the DO within the systems (Control and Pennisetum) based on bilinear interpolations between the piezometers. See legend explanation in Fig. 4.

DISCUSSION

Flow conditions and other practical issues: A challenge in studying these systems is to keep in-flow disturbances to a minimum to reduce variation in BOD and hydraulic load rates (Haddis, Van der Bruggen, & Smets, 2020). Two problems in this investigation caused alterations in the in-flow in the CWs: (i) large waste pieces (clothing, bags, plastic material, laboratory waste, etc.) caused the obstruction of the flow in the set of grilles in the Inmhoff tank; and (ii) around a year of Pennisetum growth, the root system clogged the pipes. Both problems were solved by periodic cleaning maintenance. The control system was the one that presented the minimum problems of changes in in-flow. The only problem observed was the natural growth of herbaceous plants within the CW, but it was solved by their removal during the maintenance.

Despite all the problems, the preliminary advances of this research evidence for both systems (Control and Pennisetum) acceptable performances in the general removal rates of organic pollutants. Even both CW-HSF from the first month of construction complied with the provisions of Costa Rican national legislation (Perez-Salazar, Alfaro-Chinchilla, & Scholz, 2019). Furthermore, preliminary data (unpublished data) shows that both systems report a reduction of ~ 29% of faecal coliform (log) and ~ 44% of Escherichia coli (log), where Pennisetum highlighted with 34% microorganism removal compared to Control. Overall, we do not have enough information to ensure that the use of Pennisetum emerging plants could have superior efficiency on an unplanted system, however, the results that will be discussed below show that both CWs have promising performances in wastewater treatment.

CW-HSF efficiency and the interplay between physicochemical parameters: The differences found between sampling dates for conductivity, DO and pH, the time of measurement for ORP, pH, and temperature, possibly indicate an irregularity in water in-flow quality (Haddis et al., 2020). The temperature (21,1min. – 28,2max.) and pH (6,1min. – 7,7max.). Recorded between 8:00 and 14:00 were within optimal conditions for the development of all heterotrophic microbial activities (Viessman & Hammer, 2005). CW performance is temperature dependent, as it has a strong effect on nitrogen removal, but may have less effect on BOD and phosphorous removal (Kadlec & Reddy, 2001; Tatoulis et al., 2017).

Both systems showed a significant reduction in electrical conductivity and water turbidity, explained by the rocky bed that filters most of the sediments during the horizontal flow through the system (Vymazal, 2019). In contrast, the ORP shows an increase in Pennisetum; it is a slight improvement in the ability to reduce some of the organic compounds in the medium; similar results found in Srivastava, Abbassi, Garaniya, Lewis, & Yadav (2020). Evapotranspiration is a factor influencing the redox condition, as are the flow regimes of the system (García et al., 2010). Planted wetlands produce up to 20% higher redox gradient as compared to unplanted ones (Doherty et al., 2015); where, it is in the CW depth profile, ORP less negative values near the surface have lower anaerobic conditions, while ORP more negative values in the lower layer indicate more anaerobic conditions (Srivastava et al., 2020).

It is well documented that CW-HSFs have organic compound removal rates, with average rates from 80% to 90% of nitrogen, phosphorus, and suspended solids (López-Ocaña et al., 2018; Tatoulis et al., 2017; Vymazal, 2019). Significant reductions of TKN, BOD, and COD were found for the Control and Pennisetum. These removal efficiencies of vegetated beds in CW-HSF has been reported of 75% (Vymazal, 2010); it is an unanimously agreed outcome on CW-HSFs with emerging macrophyte plants for wastewaters treatment from municipal, agriculture, industry and landfill leachate (Vymazal, 2010; Vymazal & Kröpfelová, 2009). All these results are explained by the prevailing anoxic and anaerobic conditions in the filter bed of these systems (Vymazal, 2019).

The weak but significant correlations of water electrical conductivity with the analysed variables reflect the plasticity of the conductivity and the potential use that could be made as an indicator variable of important changes within the system. Among other uses, the total dissolved solids can be estimated from conductivity (Villa et al., 2019), since it is a measure of the volume of ionized solids, such as salts and minerals. Finally, the variability explained by the PCs of all the variables measured in each of the piezometers did not allow to separate the two treatments (Control versus Pennisetum), reflecting a similar consistency of all the variables throughout the flow of water within CW.

The bilinear interpolations of the physicochemical parameters between the piezometers allowed the establishment of a gradient of the measured levels of conductivity, ORP and OD. The spatial variation of conductivity in the control showed a slight increase in the position of the second row of the piezometers, possibly an irregularity mediated by the retention times of the water within the systems. Meanwhile, in Pennisetum a gradual reduction of the average conductivity between the rows of the piezometer positions was recorded; at the same time, a reduction in its variability was observed. A change in water electrical conductivity throughout the treatment system indicates that plants` salts extraction can overcome the increase in conductivity caused by their evapotranspiration (Teixeira, Matos, Matos, Hamakawa, & Teixeira, 2020). In conclusion, the preliminary results of this research showed no spatial obstruction of the flow and/or possible reduction of the removal. However, some preferential flow paths can be distinguished throughout the system.

Applications and recommendations for spatial variation analysis with R script: This research provides the GitHub repository where it is explained step by step how to use the R script (see DOI: http://doi.org/10.5281/zenodo.3726063). Likewise, the references of the main packages used are provided; all packages used are free and constantly updated. For this reason, is recommended to install the R packages directly from Github and verify updates before running the script. Other types of CWs, plant species, a greater number of variables, no restrictions on the repetition number, and treatments can be applied to the R scripts and packages, with a simple modification of the original R script. There are no restrictions on the repetition number (number of piezometers).

This research is a preliminary advance of spatial variability monitoring results. Therefore, it is recommended to establish a long-term analysis of spatial-temporal variability of conventional pollutants on CWs (e.g.: nitrogen, phosphorous, sediments, etc.) to integrate and model the long-term removal efficiency. Also, the seasonality (dry and rainy) of tropical climate could show the performance of systems subjected to flow change effect. Finally, the main recommendation of this research for CW-HSF is the periodic monitoring of water in-flow rate, to detect obstructions in the pipelines and avoid changes in the hydraulic load rate.

Conclusion: Pennisetum demonstrated pollutant removal efficiency, but presented results like to the Control, therefore, it remains in doubt if in long-term it can be demonstrated as superior option respect to the control. Other emerging macrophyte plants should be tested. However, we can conclude that both systems had a significant effect on the removal of organic pollutants within the acceptable range (80-90%). In perspective, spatial variation analysis did not reflect any obstruction of flow along the CWs; but some preferential flow paths can be distinguished. Spatial variation is a tool that can be used in developing countries which wish to carry out monitoring campaigns to their treatment systems at a very low cost. Which can be used for the possibility of making decisions at the maintenance level and improvements in the systems. All R scripts used can be applied to other plants in CWs or similar treatments.

ACKNOWLEDGEMENTS

We gratefully give thanks to the Escuela de Ciencias Biológicas and Escuela de Química of the UNA, for allowing us the development and conclusion of this investigation. This work was supported by Funds-UNA-Costa Rica, within the framework of the project: SIA 0614-17. Finally, but no less important to Deanna Sekulich and Stefany Solano González for revising the English version of this manuscript.

ETHICAL, CONFLICT OF INTEREST AND FINANCIAL STATEMENTS

The authors declare that they have fully complied with all pertinent ethical and legal requirements, both during the study and in the production of the manuscript; that there are no conflicts of interest of any kind; that all financial sources are fully and clearly stated in the acknowledgements section; and that they fully agree with the final edited version of the article. A signed document has been filed in the journal archives.

The declaration of the contribution of each author to the manuscript is as follows: J. P. P. M.: Experiment design, data collection, statistical data analysis, Writing  -  original  draft,  Writing  -  review  &  editing; C. S.: Experiment design, data collection, Writing - review  &  editing; R. P. S.: Experiment design, data collection, Writing - review  &  editing; C. A. C.: Experiment design, data collection, Writing - review  &  editing; A. A. M.: Data collection; L. A. L.: Data collection; J. C. C.: Data collection; A. E. F.: Data collection; M. E. B.: Data collection; R. A. R.: Data collection; S. F. A.: Data collection & C. S. C.: Data collection.

REFERENCES

Ahmad, S., Majhi, P. K., Kothari, R., & Singh, R. P. (2020). Industrial Wastewater Footprinting: A Need for Water Security in Indian Context. Environmental Concerns and Sustainable Development, 1, 197–212. https://doi.org/10.1007/978-981-13-5889-0_10

AyA, MINAE & MS. (2016). Política Nacional de Saneamiento en Aguas Residuales 2016-2045 (1^st ed.). AyA-MINAE-MS. Retrieved from https://www.aya.go.cr/Noticias/Documents/Politica Nacional de Saneamiento en Aguas Residuales marzo 2017.pdf

Cui, L., Ouyang, Y., Gu, W., Yang, W., & Xu, Q. (2013). Evaluation of nutrient removal efficiency and microbial enzyme activity in a baffled subsurface-flow constructed wetland system. Bioresource Technology, 146, 656–662. https://doi.org/10.1016/j.biortech.2013.07.105

Daphne, L. H. X., Utomo, H. D., & Kenneth, L. Z. H. (2011). Correlation between Turbidity and Total Suspended Solids in Singapore Rivers. Journal of Water Sustainability, 1(3), 313–322. https://doi.org/10.11912/jws.1.3.313-322

Doherty, L., Zhao, Y., Zhao, X., Hu, Y., Hao, X., Xu, L., & Liu, R. (2015). A review of a recently emerged technology: Constructed wetland – Microbial fuel cells. Water Research, 85, 38–45. https://doi.org/10.1016/j.watres.2015.08.016

García, J., Rousseau, D. P. L., Morató, J., Lesage, E., Matamoros, V., & Bayona, J. . (2010). Contaminant Removal Processes in Subsurface-Flow Constructed Wetlands: A Review. Critical Reviews in Environmental Science and Technology, 40(7), 561–661. https://doi.org/10.1080/10643380802471076

Gomes da Silva, F. J. & Gouveia, R. M. (2020). Global Population Growth and Industrial Impact on the Environment. Cleaner Production, 33–75. https://doi.org/10.1007/978-3-030-23165-1_3

Grayson, R. B., Finlayson, B. L., Gippel, C. J., & Hart, B. T. (1996). The potential of field turbidity measurements for the computation of total phosphorus and suspended solids loads. Journal of Environmental Management, 47(3), 257–267.

Haddis, A., Van der Bruggen, B., & Smets, I. (2020). Constructed wetlands as nature based solutions in removing organic pollutants from wastewater under irregular flow conditions in a tropical climate. Ecohydrology & Hydrobiology, 20(1), 38–47. https://doi.org/10.1016/j.ecohyd.2019.03.001

He, Y., Peng, L., Hua, Y., Zhao, J., & Xiao, N. (2018). Treatment for domestic wastewater from university dorms using a hybrid constructed wetland at pilot scale. Environmental Science and Pollution Research, 25(9), 8532–8541. https://doi.org/10.1007/s11356-017-1168-7

IMN. (2020). Instituto Meterológico Nacional de Costa Rica. Retrieved January 9, 2020, from https://www.imn.ac.cr/

Kadlec, R. H. & Reddy, K. R. (2001). Temperature Effects in Treatment Wetlands. Water Environment Research, 73(5), 543–557. https://doi.org/10.2175/106143001X139614

Kassambara, A. & Mundt, F. (2020). Factoextra: extract and visualize the results of multivariate data analyses. R Package version 1.0.7. Retrieved January 9, 2020, from https://cran.r-project.org/web/packages/factoextra/index.html

López-Ocaña, G., Bautista-Margulis, R. G., Ramos-Herrera, S., Torres-Balcazar, C. A., López-Vidal, R., & Pampillón-González, L. (2018). Phytoremediation of wastewater with thalia geniculata in constructed wetlands: basic pollutants distribution. Ecology and the Environment, 228, 53–63. https://doi.org/10.2495/WP180071

Nychka, D., Furrer, R., Paige, J., & Sain, S. (2017). Package “fields” Tools for Spatial Data. R package version 11.6. January 9, 2020. Retrieved from https://doi.org/10.5065/D6W957CT

Perez-Salazar, R., Alfaro-Chinchilla, C., & Scholz, C. (2019). Experiencia en la fase de estabilización del nuevo sistema de tratamiento de aguas residuales en la Escuela de Medicina Veterinaria, Campus Benjamín Núñez, UNA. In Y. Morales-López (Ed.), Memorias del I Congreso Internacional de Ciencias Exactas y Naturales (pp. 1–8). Universidad Nacional. https://doi.org/10.15359/cicen.1.68

RCoreTeam. (2020). R: A language and environment for statistical computing (R version 3.6.1). Retrieved January 9, 2020, from https://www.r-project.org/

Rencher, A. C. (2005). A Review Of “Methods of Multivariate Analysis, Second Edition.” IIE Transactions, 37(11), 1083–1085. https://doi.org/10.1080/07408170500232784

Rice, E. W., Baird, R. B., & Eaton, A. D. (2017). Standard Methods for the Examination of Water and Wastewater, 23^rd Edition. American Public Health Association.

Rügner, H., Schwientek, M., Beckingham, B., Kuch, B., & Grathwohl, P. (2013). Turbidity as a proxy for total suspended solids (TSS) and particle facilitated pollutant transport in catchments. Environmental Earth Sciences, 69(2), 373–380. https://doi.org/10.1007/s12665-013-2307-1

Ruiz, F. F. (2012). Gestión de las Excretas y Aguas Residuales en Costa Rica: Situación Actual y Perspectiva. Retrieved January 9, 2020, from https://www.aya.go.cr/centroDocumetacion/catalogoGeneral/Gesti%C3%B3n%20de%20las%20Excretas%20y%20Aguas%20Residuales%20en%20Costa%20Rica%20%20Situaci%C3%B3n%20Actual%20y%20Perspectiva.pdf

Srivastava, P., Abbassi, R., Garaniya, V., Lewis, T., & Yadav, A. K. (2020). Performance of pilot-scale horizontal subsurface flow constructed wetland coupled with a microbial fuel cell for treating wastewater. Journal of Water Process Engineering, 33, 100994. https://doi.org/10.1016/j.jwpe.2019.100994

Tang, X., Huang, S., Scholz, M., & Li, J. (2011). Nutrient removal in vertical subsurface flow constructed wetlands treating eutrophic river water. International Journal of Environmental Analytical Chemistry, 91(7–8), 727–739. https://doi.org/10.1080/03067311003782674

Tatoulis, T., Akratos, C. S., Tekerlekopoulou, A. G., Vayenas, D. V., & Stefanakis, A. I. (2017). A novel horizontal subsurface flow constructed wetland: Reducing area requirements and clogging risk. Chemosphere, 186, 257–268. https://doi.org/10.1016/j.chemosphere.2017.07.151

Teixeira, D. L., Matos, A. T., Matos, M. P., Hamakawa, P. J., & Teixeira, D. V. (2020). Evapotranspiration of the Vetiver and Tifton 85 grasses grown in horizontal subsurface flow constructed wetlands. Journal of Environmental Science and Health, Part A, 1–8. https://doi.org/10.1080/10934529.2020.1727703

van Puijenbroek, P. J. T. M., Beusen, A. H. W., & Bouwman, A. F. (2019). Global nitrogen and phosphorus in urban waste water based on the Shared Socio-economic pathways. Journal of Environmental Management, 231, 446–456. https://doi.org/10.1016/j.jenvman.2018.10.048

Viessman, W. & Hammer, M. J. (2005). Water supply and pollution control. Michigan, USA: Pearson Prentice Hall.

Villa, A., Fölster, J., & Kyllmar, K. (2019). Determining suspended solids and total phosphorus from turbidity: comparison of high-frequency sampling with conventional monitoring methods. Environmental Monitoring and Assessment, 191(10), 605. https://doi.org/10.1007/s10661-019-7775-7

Vymazal, J. (2010). Constructed Wetlands for Wastewater Treatment. Water, 2(3), 530–549. https://doi.org/10.3390/w2030530

Vymazal, J. (2019). Is removal of organics and suspended solids in horizontal sub-surface flow constructed wetlands sustainable for twenty and more years? Chemical Engineering Journal, 378, 122117. https://doi.org/10.1016/j.cej.2019.122117

Vymazal, J. & Kröpfelová, L. (2009). Removal of organics in constructed wetlands with horizontal sub-surface flow: A review of the field experience. Science of the Total Environment, 407(13), 3911–3922. https://doi.org/10.1016/j.scitotenv.2008.08.032

Xu, Q. & Cui, L. (2019). Removal of COD from synthetic wastewater in vertical flow constructed wetland. Water Environment Research, 91(12), 1661–1668. https://doi.org/10.1002/wer.1168

Edited by Melissa Garro Garita.