Eficacia del Pennisetum sp. (Poaceae) como fitorremediador en un sistema de humedales artificiales de flujo subsuperficial mediante la medición de parámetros físico-químicos

Ana Abarca Méndez1, Leandro Araya Leitón1, Jeslyn Carranza Chaves1, Addy Echevarría Figueroa1 y Mariana Elizondo Blanco1
1 Escuela de Ciencias Biológicas, Universidad Nacional, Heredia, Costa Rica; , , , ,
Resumen: El uso de humedales artificiales como alternativa para el tratamiento de aguas residuales es una técnica empleada para evaluar la efectividad de distintas plantas como filtradores naturales y reguladores de calidad de agua. Este estudio se enfocó en la evaluación de la eficacia de dos humedales artificiales; uno de ellos con un cultivo de una planta del género Pennisetum (Poales), y otro únicamente con sustrato rocoso para ser utilizado como control. Las variables evaluadas en el agua de ambos humedales fueron: pH, conductividad, oxígeno disuelto, potencial redox y temperatura, las mediciones fueron realizadas en un período de 5 a 8 horas diarias, en tres días distintos, en el Sistema de Tratamiento Alternativo para Aguas Residuales (SATAR-UNA). La conductividad, OD y pH, se vieron afectadas por la fecha de muestreo, por lo tanto presentaron diferencias significativas, asimismo, la hora influyó en los valores de ORP, temperatura y pH. Al realizar la prueba de ANDEVA entre las variables de ambos humedales, los valores de pH, temperatura y OD mostraron leves cambios en relación entre la entrada y salida del agua tratada indiferentemente del humedal, mientras que los valores de conductividad y ORP mostraron un mayor cambio.
Palabras clave: aguas residuales, humedal artificial, contaminantes, fitorremediación.


Introducción

Las aguas residuales en la actualidad son un problema ambiental y sanitario que ha obligado a la sociedad a buscar nuevas soluciones para hallar el tratamiento adecuado para estos desechos, a través de alternativas que sean eficientes y económicamente viables. Una solución que se ha vuelto atractiva consiste en la replicación de mecanismos biológicos que consiguen limpiar las aguas contaminadas. Estos sistemas de filtración, llamados humedales artificiales, se vuelven cada vez más populares debido a que logran efluentes de buena calidad, al mismo tiempo son económicos y fáciles de mantener y operar. Se ha comprobado que estos son efectivos en la reducción de materia orgánica, en asimilación y transformación de nutrientes e incluso, llegan a retener y eliminar sustancias tóxicas (Arias I. & Brix, 2003). Además, cumplen un papel importante en la regulación hídrica y de acuíferos, retención de sedimentos, control de erosión, inclusive, algunos sustentan procesos comerciales, como también proveen servicios de recreación, investigación científica y educación (Chalarca, D, Mejía, R y Aguirre, 2006). Muchas de las aguas residuales no poseen un adecuado manejo y no reciben ningún tipo de tratamiento, es por esta razón, que los sistemas de fitorremediación son una alternativa viable y económica para mejorar la calidad de las aguas.
La filtración es producida por la densa masa de raíces que se forma reteniendo una gran cantidad de sólidos en suspensión, además, actúa como medio de soporte para el crecimiento de organismos aeróbicos que utilizan como alimento buena parte de los contaminantes presentes en el agua residual (Cooper, 1983). Los humedales artificiales (HA), están conformados principalmente por las plantas, el sustrato y la población microbiana. Las plantas suelen ser de diferentes especies y hábitos de enraizamiento y entre sus principales funciones se encuentran la absorción de nutrimentos, la relación simbiótica que se establece con los microorganismos, el suministro de oxígeno y la filtración de partículas (Brix et al. 2001). El sustrato utilizado en el proyecto fue la roca, esta le brinda soporte a la planta y permite la fijación de microorganismos al sistema, que son esenciales en el funcionamiento del humedal, pues se encargan de la remoción de los contaminantes, contribuyen a la degradación de la materia orgánica y a la transformación de compuestos nitrogenados y de fósforo contenidos en las aguas residuales a compuestos más simples (Romero, Colín, Sánchez y Ortiz, 2009).
Las aguas tratadas deben cumplir con una serie de estándares físico-químicos, los cuales se pueden revisar a través de las variables de conductividad, pH, potencial de Oxidación-Reducción (ORP), temperatura y oxígeno disuelto (OD). El pH que se encuentra en el agua debe ser lo más cercano a 7 para evitar que llegue muy básico o muy ácido a los ríos y altere los ecosistemas. El ORP se utiliza para medir la oxidación o reducción que se está llevando a cabo en el medio, donde se puede considerar que el agua cumple los estándares de potabilidad cuando este se encuentra a 650 mV (Paz et al., 2004). Las concentraciones de OD en el agua son importantes debido a que los organismos acuáticos necesitan niveles de oxígeno suficiente para su supervivencia en el agua, además una buena oxigenación del agua asegura que no habitan bacterias en el medio. La conductividad ayuda a medir la actividad eléctrica de los iones en disolución, está ligada con la cantidad de sólidos disueltos, por lo tanto una cantidad grande de conductividad indica que hay contaminación en el agua. Por último medir la temperatura del agua es conveniente debido que esta influye en todas las variables mencionadas anteriormente (Paz et al., 2004).
La adecuada selección de las plantas es fundamental, debido a que la capacidad de filtración y remoción de contaminantes difieren entre especies. Arjona (1987), realizó un estudio acerca de la eficacia del Pennisetum clandestinum, también conocido como pasto kikuyo, en la remoción de sólidos suspendidos en aguas residuales domésticas crudas. Los resultados mostraron que el Pennisetum clandestinum, reduce la concentración de sólidos básicamente por filtración de las raíces, el fraccionamiento de partículas de mayor tamaño, y el consumo de bacterias y detritus por el zooplancton y la fauna asociada al sistema radical.Por lo anterior, el objetivo de esta investigación fue analizar el grado de eficiencia del uso de Pennisetum sp. como fitorremediadora, mediante parámetros físico-químicos en un sistema de Humedales artificiales de flujo subsuperficial, para examinar su potencial fitosanitario.


Materiales y métodos


Sitio de estudio: El proyecto se realizó durante setiembre y octubre del 2019 en el Sistema de Tratamiento Alternativo para Aguas Residuales (SATAR-UNA) (9°58′37″N;84°07′37″O). de las instalaciones de la Escuela de Veterinaria ubicada en el campus Benjamin Nuñez de la Universidad Nacional, localizado en Lagunilla, Heredia. El sistema de tratamiento posee una extensión de 2050 m2 y recibe agua proveniente de los laboratorios, servicios sanitarios y el hospital, por lo tanto, existe un aporte de sustancias químicas (orgánicas e inorgánicas). Está conformado por un conjunto de rejillas, un tanque Inmhoff, cuatro celdas de humedales artificiales de flujo subsuperficial horizontal, un lecho de secado para los lodos y una unidad de desinfección, estos operan en paralelo, por lo tanto, cada uno de los humedales recibe la misma cantidad de agua residual y poseen una capacidad hidráulica para tratar un caudal de 40 m3 d-1, sumando un área superficial de 800 m2. Los humedales de flujo subsuperficial están conformados por plantas macrófitas emergentes, en este caso Pennisetum sp., el cual posee un filtro rellenado de piedra, en donde se siembran las plantas, las cuales lo atraviesan de forma horizontal, por el cual, las aguas residuales fluyen desde una entrada a través de los humedales hasta la salida de recolección del efluente.
Diseño experimental: El sistema de tratamiento se encuentra dividido en 4 celdas, con un sistema de tuberías totalmente aislado una de la otra, cada sección posee un sustrato de piedra donde se plantaron las diferentes especies filtradoras y el Control presentaba únicamente este sustrato. Cada celda del humedal posee 4 filas de pizometros distribuidos en 3 columnas para un total de 12 puntos por celda, a los cuales se les midió cada una de las variables. Se utilizaron únicamente dos de las celdas (Control y Pennisetum sp.) para realizar las mediciones del pH, la temperatura, el potencial redox (ORP), el oxígeno disuelto (OD) y la conductividad que presentaba el agua. Las mediciones se tomaron cada hora, durante un periodo de 5 a 8 horas durante 3 días, empezando desde la zona donde el agua estaba menos contaminada hasta llegar a donde estaba más contaminada. Las variables fisicoquimicas fueron medidas con sensores de pH, conductividad, temperatura, OD y ORP conectados a la interfase Vernier.
Análisis estadístico: Todas las variables (OD, pH, ORP, temperatura y conductividad) se evaluaron mediantes analisis de varianza (ANDEVA) de tres vías para los factores: fecha de muestreo, hora de muestreo y Sistema (Entrada y salida del Control y Pennisetum sp.). No se aplicaron interraciones entre los factores antes mecionado debido a no hubo sufienetes repeticiones en los sistemas. Todos los supuestos estadisticos fueron comprobados (independencia, aleatoridad, normalidad de residuos por pruebas de Shapiro-Wilk, homocedasticidad de varianza por pruebas de Levene). Se aplicaron pruebas a posteriori de LSD de Fisher (Least significant difference) para la comparación de las medias entre los sistemas (Entrada y salida del Control y Pennisetum sp.).Asimismo, se realizaron ANDEVAs de Fisher de una vía o pruebas de Kruskall-Wallis para todos los parámetros físico-químicos, esto para comparar entre las filas o columnas de la posición de los pizometros para cada sistema (Control y Pennisetum sp.), esto para verificar el cambio de los parámetros físico-químicos a lo largo del flujo del agua dentro cada sistema (celda). Igualmente, todos los supuestos estadísticos anteriores fueron comprobados. Pruebas a posteriori de Bonferoni fueron realizadas para todas las Kruskal-Wallis significativos, no se encontró ANDEVAs de Fisher significativas. Finalmente, fueron realizadas correlaciones de Spearman entre todos los parámetros físico-químicos en conjunto y por sistema (Control y Pennisetum sp.). Todos los análisis de datos se realizaron con un nivel de significacia de 95% y en lenguaje de programación R (versión 3.6.1; R Core Team, 2019).


Resultados


Comparación de parámetros fisicoquímicos entre entrada y salida del control y Pennisetum sp.: En las ANDEVAs realizadas para los puntos de entrada y de la salida de ambos tratamientos, se encontró que a nivel de sistema (Entrada y Salida del Control y Pennisetum sp.), las variables que tuvieron diferencias estadisticamente significativas fueron el ORP, OD y conductividad (Cuadro 1; p<0.05), mientras que el pH y la temperatura no mostraron diferencia estadisticamente significativa (p>0.05). La fecha del muestreo mostró diferencias estadisticamente significativas en la conductividad, OD y pH (p<0.05), igualmente la hora de medición fue signficativa en ORP, temperatura y pH.


Cuadro 1. Comparación de los parámetros fisicoquímicos en los flujos de entrada y salida (Control y Pennisetum) en un sistema de Humedal artificial de flujo subsuperficial, mediante ANDEVA de tres vías (Factores: Fecha, Hora y Sistema)
Media ± desviación estándar; R2 : coeficiente de determinación; F: valor de Fisher; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre los sistemas (LSD, p>0.05)

Media ± desviación estándar; R2 : coeficiente de determinación; F: valor de Fisher; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre los sistemas (LSD, p>0.05)


Comparación de parámetros fisicoquímicos dentro los sistemas de control y Pennisetum sp.: No se encontró diferencia estadísticamente significativa para todos los parametros físico-químico entre la posición de los pizometros en columnas para ambos sistemas (KW< 4.65; p>0.05, ver Fig. 1-5).Sin embargo, se encontraron diferencias estadísticamente significativa para los parámetros físico-químico con excepción de OD y temperatura, entre la posición de los pizometros en las filas para ambos sistemas (KW>5.07; p<0.05, ver Fig. 1-5).


Fig. 1. Distribución de los niveles de pH dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de pH de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre la posición de los pizometros en filas (Bonferroni, p>0.05).

Fig. 1. Distribución de los niveles de pH dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de pH de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre la posición de los pizometros en filas (Bonferroni, p>0.05).


Fig. 2. Distribución de los niveles de OD dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de OD de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad.

Fig. 2. Distribución de los niveles de OD dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de OD de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad.


Fig. 3. Distribución de los niveles de conductividad eléctrica dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de conductividad de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre la posición de los pizometros en filas (Bonferroni, p>0.05).

Fig. 3. Distribución de los niveles de conductividad eléctrica dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de conductividad de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre la posición de los pizometros en filas (Bonferroni, p>0.05).


Fig. 4. Distribución de los niveles de temperatura (°C) dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de temperatura de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad.

Fig. 4. Distribución de los niveles de temperatura (°C) dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de temperatura de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad.


Fig. 5. Distribucción de los niveles de ORP dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de ORP de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre la posición de los pizometros en filas (Bonferroni, p>0.05).

Fig. 5. Distribucción de los niveles de ORP dentro de los sistemas (Control y Pennisetum sp.) a partir de interpolaciones cuadrática entre los pizometros. Cruces de color blanco indican posición de los pizometros; gradiente de color y curvas de nivel indican la intensidad de ORP de bajo (azul) a alto (rojo); KW: Krukall-Wallis; R2: coeficiente de determinación; p: probabilidad; letras iguales indican no diferencia estadísticamente significativa entre la posición de los pizometros en filas (Bonferroni, p>0.05).


Relación entre los parámetros físico-químicos: El pH y OD no mostró correlación estadisticamente significativa, independientemente del sistema o todos los datos analizados en conjunto (rho<0.1; p>0.05) (Fig. 6). No obstante, todos los demás pares de combinación de parametros fisico-quimicos por sistema o en conjunto mostraron al menos una correlación estadisticamente significativa de Spearman (p<0.05, Fig. 6).


Fig. 6. Correlaciones de Spearman entre los parámetros físico-químicos (pH, temperatura, conductividad, ORP y OD) al nivel de sistema (control o Pennisetum sp.) y ambos sistemas en conjunto.

Fig. 6. Correlaciones de Spearman entre los parámetros físico-químicos (pH, temperatura, conductividad, ORP y OD) al nivel de sistema (control o Pennisetum sp.) y ambos sistemas en conjunto.


Discusión


El ORP es una variable que mide el potencial de oxidación-reducción del medio y su capacidad de interactuar como oxidante o reductor con las sustancias que se encuentren en el mismo, es afectado directamente por las concentraciones de iones que lo acidifican o alcalinizan, porque involucra una transferencia de electrones de una sustancia a otra. Según de la Vega., et al (2012), en los procesos de biotratamiento, el ORP refleja las variaciones en la concentración de OD, sustrato orgánico, actividad de microorganismos e inhibición del proceso. Esta inhibición generalmente es causada por una demanda química insuficiente de oxígeno o por la presencia de algún compuesto tóxico en el afluente. El agua con valores elevados de ORP tiende a comportarse como agente oxidante, en caso contrario, si presenta valores bajos, se comporta como agente reductor; los resultados muestran que el valor del flujo de salida de Pennisetum sp., presentó diferencias en comparación al flujo de entrada general del sistema. El flujo de entrada medido para ORP se encontró diferente del flujo de salida de Pennisetum sp., mientras que el flujo de salida de control no mostró significancia (Cuadro 1), esto podría sugerir que en el momento del día cuando se realiza la medición podrían afectar el flujo del agua factores externos que llegan hasta la cañería que escurre al sistema.
El flujo en Pennisetum sp., tiende a entrar en la celda con valores negativos de ORP (reductor) pero conforme avanza por las filas comienza a subir (oxidante) hasta salir(Fig.5). De acuerdo a Lynch y Poole (1979) citado por Fuentes, esto es indicativo de un ambiente altamente reductor, que puede ser ocasionado por las descargas masivas de materia orgánica oxidable del alcantarillado sanitario, lo que aumenta la densidad de bacterias anaerobias facultativas, como las bacterias del grupo Coliformes, quienes se encuentran en grandes cantidades en los ríos estudiados, tal como se muestra en los resultados obtenidos del análisis microbiológico y cuya actividad metabólica altera el potencial Redox. En Control se aprecia como el flujo entra a la celda con valores bajos que comienzan a aumentar pero manteniéndose constantes hasta la salida. Además, el OD también se ve afectado por la presencia de iones en el medio. Un estudio realizado en los ríos Machángara y Monjas en el distrito metropolitano de Quito, demostró la relación entre el OD y potencial Redox, estos son directamente proporcionales, es decir, a medida que el ORP disminuye por la presencia de altas poblaciones microbianas, también disminuye la concentración de oxígeno disuelto, esto causa una reducción de iones y moléculas importantes para la nutrición de microorganismos y otras formas de vida superior (Campaña, Gualoto & Chiluisa, 2017).
El Oxígeno Disuelto (OD), representa un parámetro físico-químico de gran trascendencia, debido a que su reducción indica la presencia de agentes contaminantes de fácil oxidación en el agua. La concentración de oxígeno en el agua proviene principalmente del movimiento de la misma, causado por resaltos hidráulicos naturales y por el intercambio gaseoso, este fenómeno se ve reflejado en la celda que pertenece a control (Fig.2) donde al no haber factores que disminuyen el flujo del agua a través del sustrato, producen que el agua choque libremente contra las rocas y hace que esta se oxigene en el proceso. Contrario a lo que ocurre en Pennisetum sp (Fig.2), donde podemos ver que los niveles de OD tienden a mantenerse constantes desde que entran hasta que salen, esto producto a que las raíces impiden el flujo normal del agua y además las plantas extraen el oxígeno que se encuentra disuelto en el medio para realizar los procesos metabólicos de la planta (Pajoy Muñoz, H. M., 2017), otros factores que afectan en la disminución pueden deberse al aumento en la temperatura, materia orgánica y agentes contaminantes, lo que puede observarse en la correlación que existe entre el OD y las variables de temperatura, ORP y conductividad (Fig.6). (Ramírez, 2015).
El OD en el efluente de los humedales artificiales que contienen plantas, tiende a presentar valores bastante similares y sin mucha fluctuación, esto se debe al aumento en el sistema radicular de las plantas, de esta forma se incrementa la porosidad de la raíz y aumenta el intercambio gaseoso (Fig.2) (Ahmadi et al., 2013),donde Pennisetum sp presenta valores de OD relativamente parecidos entre un rango de 3.7 a 3.9. El flujo de la salida de Pennisetum sp., mostró diferencias significativas en relación a la salida de Control (Cuadro 1); la Fig.2 muestra el flujo del agua en ambas celdas, en Pennisetum sp, se aprecia como este mantiene la tendencia, mencionada anteriormente, a permanecer constante en un rango de valores similares mientras fluye por el sistema, saliendo de este con valores levemente más elevados, es decir, con mayor saturación de oxígeno disuelto pero en menor proporción que en Control, esto podría deberse a un proceso de re-aireación donde el roce del agua contra las piedras del sustrato crea efectos de turbulencia, aumentando la concentración de oxígeno sobre la superficie del agua;. Valores bajos podrían deberse a la acumulación de materia orgánica en las raíces de Pennisetum sp., que puede entrar en descomposición por acción bacteriana, lo que aumenta la absorción de oxígeno. Las plantas transportan desde sus raíces hasta las hojas el oxígeno disuelto en el agua para realizar las funciones metabólicas lo cual implica un consumo del oxígeno que se encuentra en el medio, cosa que no ocurre en Control ya que en este humedal las plantas no utilizan el oxígeno (Lahora, 2003), incluso, otros factores que requieren oxígeno como procesos oxidativos de compuestos, podrían descender los niveles de OD en la parcela. Valores altos de OD muestran mejor calidad de agua.
La conductividad es una medida indirecta de la cantidad de sales o iones disueltos en el medio (Raffo Lecca, 2016), esta medida muestra la capacidad del agua para conducir corriente eléctrica, por lo que también es afectada por la bioacumulación y por ende, por la concentración de solutos principalmente sales disueltas en el agua; la conductividad del agua tenderá a aumentar a medida que aumente la concentración de iones; en este sentido, el flujo de salida de Control y Pennisetum sp., no mostraron diferencia significativa entre sí, pero en comparación con la entrada del flujo, mostraron valores más bajos (Cuadro 1) ; en la Fig. 3 se aprecia que el flujo de Control presenta ligeramente menor conductividad mientras avanza por las filas en la celda, esto podría deberse a que en Pennisetum sp., hay mayor presencia de materia orgánica y por lo tanto de iones y moléculas cargadas que contribuirían a una mejor conductividad. En el flujo de entrada valores elevados de conductividad podrían sugerir la presencia de metales pesados con capacidad de conducción procedentes de detergentes y otros químicos arrastrados desde las tuberías.( Mackinlay Llanza, F.,2017)
Las variables pH y temperatura, al ser parámetros que son afectados por la temperatura ambiente del día y hora de la medición, pueden presentar valores diferentes entre sí; aun así, estas variables se mantuvieron estables con poca fluctuación. La temperatura tiene gran importancia en los humedales debido a que esta surge efecto sobre diversos procesos biológicos, como lo son la absorción de oxígeno por las plantas que se encuentran en el humedal, la precipitación de compuestos, entre otros. La variación de la temperatura que se presentó en los dos humedales, siendo Control el de mayores temperaturas en la salida (Fig 4) es producto a que este se encuentra expuesto a la luz del día, y en Pennisetum sp., las plantas crean una sombra la cual evita la exposición, además que la transpiración de las plantas disminuye la temperatura del humedal, por lo tanto, la temperatura mostró un cambio significativo durante el transcurso del día (Cuadro 1), (Raffo Lecca, 2016).Por último, para la variable de pH, los valores se mantuvieron constantes, no mostrando diferencia entre ambos sistemas (Cuadro 1). En la Fig. 1 se aprecia como ambos flujos tienden a volverse más básicos conforme avanzan por las filas.
Según un estudio realizado por Ramírez, 2015, en el diseño de un humedal artificial con Chrysopogon zizanioides y Pennisetum purpureum, lo ideal es tener un afluente y un efluente con un pH cercano a la neutralidad, el Pennisetum purpureum mejor conocido como “hierba de elefante”, presentó un valor máximo de 7,2. En el presente estudio se da el mismo fenómeno, debido a que tanto en Control como en Pennisetum sp., el máximo valor de pH fue aproximadamente de 7,40. El pH es uno de los parámetros más importantes en los sistemas de tratamiento de aguas residuales ya sea doméstica o industrial, además es uno de los principales indicadores en campo, debido a que permite conocer los niveles de acidez y basicidad del agua residual. Controlar el pH en los humedales artificiales es de suma importancia tanto para las plantas como para los microorganismos, debido a que permite mantener la cinética del crecimiento bacteriano, es decir, estos organismos crecen bajo un rango de pH óptimo, el cual se encuentra entre la neutralidad y medios ligeramente ácidos, por lo tanto, si este rango se altera, y el medio es muy ácido, los microorganismos mueren o disminuyen su tasa de crecimiento y comienzan a proliferar hongos, afectando la remoción de materia orgánica en el afluente y provocando la muerte de las plantas y posteriormente del humedal. (Ramírez, 2015). La variable pH se encuentra correlacionada con la conductividad, el OD y la temperatura (Fig 6), esto es producto a la relación que existe entre los iones disueltos en el agua que son los responsables de aumentar o disminuir el pH y la conductividad, ya que la presencia de estos iones producen que haya una mayor conductividad en el ambiente. (Mackinlay Llanza, F.,2017), lo cual ocurre parecido con temperatura y OD. Rittmann y McCarty (1997), Ávila et al. (2013) y Gikas y Tsihrintzis (2012), demostraron que existe un pH óptimo de funcionamiento en los sistemas de humedales artificiales que beneficia a las plantas y a los microorganismos presentes, este se encuentra en un rango de 6,5 a 7,5, estos datos coinciden con los de la Fig. 1, debido a que el Pennisetum sp. presentó valores dentro de este rango, reafirmando un pH relativamente neutro y óptimo para el desarrollo de microorganismos.


Conclusiones


Para monitorear un sistema de agua residual, indicadores físico-químicos como OD, ORP y conductividad son herramientas viables para evaluar la calidad del agua. Las fluctuaciones observadas en estas variables mostraron diferencias en su comportamiento a través de cada tratamiento independientemente; estas diferencias comprobaron que el Pennisetum sp. produce cambios significativos en la conducta de estas variables. La constancia en los valores de OD en el tratamiento con Pennisetum sp. en comparación con el Control donde se observó un mayor incremento, muestra el consumo de oxígeno llevado a cabo por las plantas, y la disminución en la conductividad representa la eliminación de sólidos contaminantes suspendidos en el agua. Este comportamiento observado en OD y conductividad mostró su relación directa con ORP al elevar paulatinamente los valores de dicha variable a medida que la concentración de OD aumenta levemente y los sólidos disueltos en el agua disminuyen. Los cambios observados respecto al rumbo que toman los valores de las variables a lo largo de un humedal artificial, al fluir las aguas residuales en presencia de Pennisetum sp, mostraron mejoría en la calidad del agua; esto comprueba la efectividad relativa de dicha planta como agente biorremediador. El uso de humedales artificiales muestra ser un recurso viable para el tratamiento de aguas residuales y el evaluar la eficacia de distintas plantas en estos tratamientos, resulta ser de gran relevancia para la identificación de especies recomendables para su utilización en proyectos de biorremediación. Para futuras investigaciones, se recomienda aumentar el número de mediciones y adicionar variables biológicas que permitan evaluar otros agentes indicadores de calidad de agua, además, seleccionar fechas de muestreos aleatoriamente durante las distintas épocas del año, esto para evidenciar la diferenciación de la variables a lo largo de un periodo de tiempo más extenso.


Agradecimientos


Agradecemos al Dr. Junior Pérez Molina y a la Dr. Carola Scholz, ambos profesores de la Universidad Nacional de Costa Rica, por su apoyo y tutoría brindada en esta investigación.


Referencias


Ahmadi, S., Mostajeran, A., & Shokrollahi, S. (2013). Comparing root porosity of sunflower adventitious root segments using cross-sectioning and buoyancy method under hypoxic condition. Plant root, 7, 12-20.


Alfaro, C., Pérez, R. y Solano, M. (2013). Saneamiento de aguas residuales mediante humedales artificiales en el Museo de Cultura Popular de la Universidad Nacional. Revista de Ciencias Ambientales (Trop J Environ Sci). 45(1): 63-71.


Arias I., C., & Brix, H. (2003). Humedales artificiales para el tratamiento de aguas residuales. Ciencia E Ingeniería Neogranadina, 13(1), 17-24. doi: 10.18359/rcin.1321.


Arjona, B., 1987. Evaluación de un cultivo hidropónico de Pennisetum clandestinum Hochst (kikuyo) como tratamiento biológico para aguas residuales domésticas. Tesis Biología. Bogotá. Universidad Nacional de Colombia. Facultad de Ciencias. 152p.


Ávila, C., Garfí, M., y García, J. (2013). Three-stage hybrid constructed wetland system for wastewater treatment and reuse in warm climate regions. Ecological Engineering, 61:43–49. 1, 2, 9.4.1.


Brix H., Arias C. y Bubba M. (2001). Media selection for sustainable phosphorus removal in subsurface flow constructed wetlands. Water Sci. Technol. 44, 47–54.


Castrillón, O., Bedoya, O., & Montoya, D. V. (2006). Effect of pH on the growth of microorganisms during the maturation stage in static compost piles, 1(2), 8-12.


Chalarca Rodríguez, D. A., Mejía Ruiz, R., & Aguirre Ramírez, N. J. (2007). Aproximación a la determinación del impacto de los vertimientos de las aguas residuales domésticas del municipio de Ayapel, sobre la calidad del agua de la ciénaga. Revista Facultad de Ingeniería N.o 40. pp. 41-58.


Cooper, W. J., L. Handley, L. S. Casy, J. L. Lopez, C. H. Fonseca, S.H. Rizvi, J. Richards y A. Vanzella. (1983). Development of the nutrient film technique for renovation of wastewater for non-potable reuse. Report to U S. Bureau of reclamation Drinking Water Research Center Florida International University. 70 p.


De la Vega, P. M., de Salazar, E. M., Jaramillo, M. A., & Cros, J. (2012). New contributions to the ORP & DO time profile characterization to improve biological nutrient removal. Bioresource technology, 114, 160-167.


Gikas, G. D., & Tsihrintzis, V. A. (2012). A small-size vertical flow constructed wetland for on-site treatment of household wastewater. Ecological Engineering, 44, 337-343.


Huamán, F. T. G., Delgado, J. T., & Medrano, S. E. V. (2014). Calidad ecológica del agua del río Utcubamba en relación a parámetros fisicoquímicos y biológicos. Amazonas, Perú. Sciendo, 14(1).


Lahora, A. (2003). Depuración de aguas residuales mediante humedales artificiales: La EDAR de los Gallardos (Almería). Ecología, Manejo Y Conservación De Humedales., 1(1), 99-112.


Mackinlay Llanza, F. (2017) Recuperación de suelos contaminados en campos de tiro por electrorremediación y fitorremediación (Tesis de pregrado). Universidad de Vigo, España.


Morel, A. and Diener, S. (2006) Greywater management in low and middle-income countries, review of different treatment systems for households or neighbourhoods. Swiss Federal Institute of Aquatic Science and Technology (Eawag). Dübendorf, Switzerland.


Pajoy Muñoz, H. M. (2017).Potencial fitorremediador de dos especies ornamentales como alternativa de tratamiento de suelos contaminados con metales pesados( Tesis de Maestría). Universidad Nacional de Colombia, Medellín, Colombia.


Paz, M., Muzio, H., Gemini, V., Magdaleno, A., Rossi, S., Korol, S., & Moretton, J. (2004). Aguas residuales de un Centro Hospitalario de Buenos Aires, Argentina: características químicas, biológicas y toxicológicas. Higiene Y Sanidad Ambiental, 4(1), 83-88.


R Core Team (2019) R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/.


Romero-Aguilar, M., Colín-Cruz, A., Sánchez-Salinas, E., & Ortiz-Hernández, M. A. (2009). Tratamiento de aguas residuales por un sistema piloto de humedales artificiales: evaluación de la remoción de la carga orgánica. Revista internacional de contaminación ambiental, 25(3), 157-167.


Ramírez Cadavid, J. D. (2015). Evaluación de vetiver (Chrysopogon zizanioides) y la elefanta (Pennisetum purpureum) en el diseño de humedales artificiales. Rittman, B. y McCarty, P. (1997). Biotecnología del medio ambiente Principios y aplicaciones. McGraw Hill, España.



ANEXO: Código de R para todos los analisis estadísticos ejecutados

#------------------------------------------------
# Packages
#------------------------------------------------
library("readxl")
library(openxlsx)
library(Rmisc)
library(fields) #https://www.rdocumentation.org/packages/fields/versions/9.8-6/topics/interp.surface
library(plot3D) #https://rpubs.com/yoshio/95844
library(yarrr)
library(broom)
library(car)
library(lsmeans)
library(multcompView)
library(multcomp)
library(dplyr)
library(GGally)
library(factoextra)
library(cowplot)
library(ggplot2)
library(grid)
library(gridExtra)
#------------------------------------------------



#------------------------------------------------
# Loading - Database
#------------------------------------------------
dat<-read_excel("C:/Users/jabar/Downloads/humedal/humedal/Data/Data.xlsx", sheet = "Data")
dat<-data.frame(dat)
dat$Rep<-format(dat$Rep,"%H:%M:%S")
dat$Fecha<-as.Date(dat$Fecha)
pe<-dat[dat$Tratamiento=="Entrada", ]
ps_p<-dat[dat$Tratamiento=="Salida Pennisetum", ]
ps_c<-dat[dat$Tratamiento=="Salida Control", ]
str(dat)
#------------------------------------------------



#------------------------------------------------
# Cuadro 1 - Comparaci?n entrada y salida de los sistemas (Control y Pennisetum)
#------------------------------------------------
g<-merge(pe, ps_c, all=TRUE)
g<-merge(g, ps_p, all=TRUE)
g$Tratamiento_Fecha<-paste(g$Tratamiento, g$Fecha, sep="*")
g$Tratamiento_Fecha_Rep<-paste(g$Tratamiento_Fecha, g$Rep, sep="*")
g$Fecha<-as.factor(g$Fecha)
g$Rep<-as.factor(g$Rep)
g$Tratamiento<-as.factor(g$Tratamiento)
str(g)
#------------------------------------------------
fun.table<-function(model1, model2, Variable){
  shapiro.test(model1$residuals)
  inf<-round(glance(model1),3)
  P.value<-round(as.numeric(inf[5]),3)
  P.value<-ifelse(P.value<0.001, "***",
         ifelse(P.value<0.01, "**",
                ifelse(P.value<0.05, "*","n.s.")))
  
  R2<-round(as.numeric(inf[2]),2)
  F<-round(as.numeric(inf[4]),1)
  coef<-Anova(model1, type="II")
  P<-data.frame(round(coef$`Pr(>F)`,3))
  P<-data.frame(t(P))
  
  sum = summarySE(g, measurevar= Variable, groupvars=c("Tratamiento"), na.rm=TRUE)
  sum<-sum[c(1,2,3,5,6)]
  sum<-data.frame(sum)
  sum<-data.frame(Sistema= sum$Tratamiento, mean=paste(round(sum[,3], 1), round(sum[,4], 2), sep=" ± "))
  sum<-data.frame(t(sum))
  names(sum)<-c("Entrada", "Salida Control","Salida Pennisetum")
  sum<-sum[-1,]
  
  if(P[,3]<0.05){
    lsm = lsmeans(model2, pairwise ~ Tratamiento, adjust="LSD")
    t1<-data.frame(cld(lsm[[1]], alpha=.05, Letters=letters))
    t2<-data.frame(t1[c(1)],t1[c(7)])
    t2<-t2[order(t2$Tratamiento),]
    rm.whitespace <- function (x) gsub("^\\s+|\\s+$", "", x)
    t2$.group<-rm.whitespace(t2$.group)
  } else {
    t2<-data.frame(Tratamiento=c("","",""),".group"=c("","",""))
  }
  
  p<-data.frame(Variable=Variable,Fecha=P[,1],Hora=P[,2],Sistema=P[,3],
                Entrada=paste(sum[1,1], t2[1,2], sep=" "), 
                "Salida Control"=paste(sum[1,2], t2[2,2], sep=" "),
                "Salida Pennisetum"=paste(sum[1,3], t2[3,2], sep=" "),
                R2, F, P.value)
  s<-ifelse(p[,c(2,3,4)]<0.001, "***",
            ifelse(p[,c(2,3,4)]<0.01, "**",
                   ifelse(p[,c(2,3,4)]<0.05, "*","n.s.")))
  p[,2]<-s[,1]
  p[,3]<-s[,2]
  p[,4]<-s[,3]
  return(p)
}
#------------------------------------------------
OD_es<-aov(OD ~  Fecha + Rep + Tratamiento, data=g)
summary(OD_es)
shapiro.test(OD_es$residuals)
leveneTest(g$OD, g$Tratamiento)
OD_es2<-aov(OD ~  Tratamiento, data=g)
summary(OD_es2)
OD_t<-fun.table(model1=OD_es, model2= OD_es2,Variable="OD")
#------------------------------------------------
pH_es<-aov(pH ~  Fecha + Rep + Tratamiento, data=g)
summary(pH_es)
shapiro.test(pH_es$residuals)
leveneTest(g$pH, g$Tratamiento)
pH_es2<-aov(pH ~  Tratamiento, data=g)
summary(pH_es2)
pH_t<-fun.table(model1=pH_es, model2= pH_es2,Variable="pH")
#------------------------------------------------
ORP_es<-aov(ORP ~  Fecha + Rep + Tratamiento, data=g)
summary(ORP_es)
shapiro.test(ORP_es$residuals)
leveneTest(g$ORP, g$Tratamiento)
ORP_es2<-aov(ORP ~  Tratamiento, data=g)
summary(ORP_es2)
ORP_t<-fun.table(model1=ORP_es, model2= ORP_es2,Variable="ORP")
#------------------------------------------------
Temperatura_es<-aov(Temperatura ~  Fecha + Rep + Tratamiento, data=g)
summary(Temperatura_es)
shapiro.test(Temperatura_es$residuals)
leveneTest(g$Temperatura, g$Tratamiento)
Temperatura_es2<-aov(Temperatura ~  Tratamiento, data=g)
summary(Temperatura_es2)
Temperatura_t<-fun.table(model1=Temperatura_es, model2= Temperatura_es2,Variable="Temperatura")
#------------------------------------------------
Conductividad_es<-aov(Conductividad ~  Fecha + Rep + Tratamiento, data=g)
summary(Conductividad_es)
shapiro.test(Conductividad_es$residuals)
leveneTest(g$Conductividad, g$Tratamiento)
Conductividad_es2<-aov(Conductividad ~  Tratamiento, data=g)
summary(Conductividad_es2)
Conductividad_t<-fun.table(model1=Conductividad_es, model2= Conductividad_es2,Variable="Conductividad")
#------------------------------------------------
Cuadro_1<-merge(ORP_t, Conductividad_t, all = TRUE)
Cuadro_1<-merge(Cuadro_1, OD_t, all = TRUE)
Cuadro_1<-merge(Cuadro_1, Temperatura_t, all = TRUE)
Cuadro_1<-merge(Cuadro_1, pH_t, all = TRUE)
Cuadro_1
#------------------------------------------------



#------------------------------------------------
# Cuadro 2 - Comparaci?n entre las filas y columnas dentro de cada sistema
#------------------------------------------------
dat$Columna<-dat$Punto
dat$Fila<-dat$Punto

dat$Columna[dat$Punto==1]<-2.70
dat$Columna[dat$Punto==2]<-5.95
dat$Columna[dat$Punto==3]<-9.2

dat$Columna[dat$Punto==4]<-2.70
dat$Columna[dat$Punto==5]<-5.95
dat$Columna[dat$Punto==6]<-9.2

dat$Columna[dat$Punto==7]<-2.70
dat$Columna[dat$Punto==8]<-5.95
dat$Columna[dat$Punto==9]<-9.2

dat$Columna[dat$Punto==10]<-2.70
dat$Columna[dat$Punto==11]<-5.95
dat$Columna[dat$Punto==12]<-9.2

dat$Fila[dat$Punto==1]<-2.80
dat$Fila[dat$Punto==2]<-2.80
dat$Fila[dat$Punto==3]<-2.80

dat$Fila[dat$Punto==4]<-7.15
dat$Fila[dat$Punto==5]<-7.15
dat$Fila[dat$Punto==6]<-7.15

dat$Fila[dat$Punto==7]<-11.5
dat$Fila[dat$Punto==8]<-11.5
dat$Fila[dat$Punto==9]<-11.5

dat$Fila[dat$Punto==10]<-15.85
dat$Fila[dat$Punto==11]<-15.85
dat$Fila[dat$Punto==12]<-15.85
#------------------------------------------------
e<-as.numeric(c(row.names(dat[dat$Tratamiento=="Entrada",]),
           row.names(dat[dat$Tratamiento=="Salida Control",]),
           row.names(dat[dat$Tratamiento=="Salida Pennisetum",])))
dat2<-dat[-c(e),]
#------------------------------------------------
ORP_aov<-aov( ORP ~ Columna + Fila, data=dat2)
summary.lm(ORP_aov)
shapiro.test(ORP_aov$residuals)
bartlett.test(dat2$ORP, dat2$Fila)
#------------------------------------------------
Conductividad_aov<-aov( Conductividad ~ Columna + Fila, data=dat2)
summary.lm(Conductividad_aov)
shapiro.test(Conductividad_aov$residuals)
bartlett.test(dat2$Conductividad, dat2$Fila)
#------------------------------------------------
OD_aov<-aov( OD ~ Columna + Fila, data=dat2)
summary.lm(OD_aov)
shapiro.test(OD_aov$residuals)
bartlett.test(dat2$OD, dat2$Fila)
#------------------------------------------------
Temperatura_aov<-aov( Temperatura ~ Columna + Fila, data=dat2)
summary.lm(Temperatura_aov)
shapiro.test(Temperatura_aov$residuals)
bartlett.test(dat2$Temperatura, dat2$Fila)
#------------------------------------------------
pH_aov<-aov( pH ~ Fila, data=dat2)
summary.lm(pH_aov)
shapiro.test(pH_aov$residuals)
bartlett.test(dat2$pH, dat2$Fila)
#------------------------------------------------



#------------------------------------------------
# Figure - Analisis espacial
#------------------------------------------------
fun.plot3d<-function(data, var1, var2, tratamiento1, tratamiento2, Variable){
  # Paso 1
  sum = summarySE(data, measurevar= Variable, groupvars=c("Tratamiento", "Punto"), na.rm=TRUE)
  sum<-sum[,c(1,2,3,4,6,7)]
  sum<-data.frame(Variable, sum)
  names(sum)<-c("Variable","Tratamiento","Punto","N","Mean","S.E.","C.I.95")
  sum
  sum1<-matrix(sum$Mean[sum$Tratamiento=="Control"],nrow = 3, ncol = 4)
  sum2<-matrix(sum$Mean[sum$Tratamiento=="Pennisetum"],nrow = 3, ncol = 4)
  # Paso 2A
  layout(matrix(c(1,1, 2,2, 3,3, 4,4,
                  1,1, 2,2, 3,3, 4,4, 
                  7,7, 5,5, 8,8, 6,6,
                  0,0, 0,0, 0,0, 0,0,
                  0,0, 0,0, 0,0, 0,0), nrow = 5, byrow=T))
  pm <- par("mfrow")
  
  par(xpd = FALSE, mgp = c(1.5,0.5,0), mar = c(1.5,4,1.5,1))
  boxplot(var1 ~ dat$Fila[dat$Tratamiento=="Control"], xlab=Variable, ylab= "Fila: Distancia (m)",horizontal=TRUE, col="gray45")
  
  #-------------------
  x=c(2.70, 5.95, 9.2)
  y=c(2.8, 7.15, 11.5, 15.85)
  #-------------------
  
  par(xpd = TRUE, mgp = c(1.5,0.5,0), mar = c(1.5,0.5,2,2.5)) #contour = list(lwd = 2, col = jet.col(11))
  obj<- list( x= x, y=y, z= sum1)
  set.seed(123)
  grid.list<- list( x= seq( min(x),max(x),,100), y=  seq( min(y),max(y),,100))
  m<-interp.surface.grid(obj, grid.list)
  image2D(z = m, lwd = 3, shade = 0.2, rasterImage = TRUE, contour=TRUE, main = tratamiento1, clab = sum$Variable[1], xlab="", ylab="")
  grid <- mesh(dat$Columna, dat$Fila)
  points(grid, pch=3, lwd=2, cex=1, col="White")
  
  par(xpd = FALSE, mgp = c(1.5,0.5,0), mar = c(1.5,4,1.5,1))
  boxplot(var2 ~ dat$Fila[dat$Tratamiento=="Pennisetum"],  xlab=Variable, ylab= "Fila: Distancia (m)",horizontal=TRUE, col="gray45")
  
  par(xpd = TRUE, mgp = c(1.5,0.5,0), mar = c(1.5,0.5,2,2.5)) #contour = list(lwd = 2, col = jet.col(11))
  obj<- list( x= x, y=y, z= sum2)
  set.seed(123)
  grid.list<- list( x= seq( min(x),max(x),,100), y=  seq( min(y),max(y),,100))
  m<-interp.surface.grid(obj, grid.list)
  image2D(z = m, lwd = 3, shade = 0.2, rasterImage = TRUE, contour=TRUE, main = expression(paste(italic("Pennisetum"), " sp.")), clab = sum$Variable[1], xlab="", ylab="")
  grid <- mesh(dat$Columna, dat$Fila)
  points(grid, pch=3, lwd=2, cex=1, col="White")
  
  par(xpd = FALSE, mgp = c(1.5,0.5,0), mar = c(3,1,1,3))
  boxplot(var1 ~ dat$Columna[dat$Tratamiento=="Control"],  ylab=Variable, xlab= "Columna: Distancia (m)", horizontal=FALSE, col="gray45")
  
  
  par(xpd = FALSE, mgp = c(1.5,0.5,0), mar = c(3,1,1,3))
  boxplot(var2 ~ dat$Columna[dat$Tratamiento=="Pennisetum"],  ylab=Variable, xlab= "Columna: Distancia (m)", horizontal=FALSE, col="gray45")
  
  #------------------------------------------------
  par(xpd = FALSE, mgp = c(1.5,0.5,0), mar = c(0,0,0,0)) #contour = list(lwd = 2, col = jet.col(11))
  obj<- list( x= x, y=y, z= sum1)
  set.seed(123)
  grid.list<- list( x= seq( min(x),max(x),,100), y=  seq( min(y),max(y),,100))
  m<-interp.surface.grid(obj, grid.list)
  x <- 1 : nrow(m$z)
  y <- 1 : ncol(m$z)
  panelfirst <- function(pmat) {
    XY <- trans3D(x = rep(1, ncol(m$z)), y = y,
                  z = m$z[50,], pmat = pmat)
    scatter2D(XY$x, XY$y, colvar = m$z[50,],
              type = "l", lwd = 3, add = TRUE, colkey = FALSE)
    XY <- trans3D(x = x, y = rep(ncol(m$z), nrow(m$z)),
                  z = m$z[,50], pmat = pmat)
    scatter2D(XY$x, XY$y, colvar = m$z[,50],
              type = "l", lwd = 3, add = TRUE, colkey = FALSE)
  }
  pmat <- persp3D(z = m$z, x = x, y = y, scale = FALSE, theta = 30,
                  expand = 0.1, panel.first = panelfirst, colkey = FALSE,contour=FALSE)
  XY <- trans3D(x = rep(50, ncol(m$z)), y = y, z = m$z[50,],
                pmat = pmat)
  lines(XY, lwd = 1, lty = 3)
  XY <- trans3D(x = x, y = rep(50, nrow(m$z)), z = m$z[,50],
                pmat = pmat)
  lines(XY, lwd = 1, lty = 3)
  
  #-------------------
  x=c(2.70, 5.95, 9.2)
  y=c(2.8, 7.15, 11.5, 15.85)
  #-------------------
  
  obj<- list( x= x, y=y, z= sum2)
  set.seed(123)
  grid.list<- list( x= seq( min(x),max(x),,100), y=  seq( min(y),max(y),,100))
  m<-interp.surface.grid(obj, grid.list)
  x <- 1 : nrow(m$z)
  y <- 1 : ncol(m$z)
  panelfirst <- function(pmat) {
    XY <- trans3D(x = rep(1, ncol(m$z)), y = y,
                  z = m$z[50,], pmat = pmat)
    scatter2D(XY$x, XY$y, colvar = m$z[50,],
              type = "l", lwd = 3, add = TRUE, colkey = FALSE)
    XY <- trans3D(x = x, y = rep(ncol(m$z), nrow(m$z)),
                  z = m$z[,50], pmat = pmat)
    scatter2D(XY$x, XY$y, colvar = m$z[,50],
              type = "l", lwd = 3, add = TRUE, colkey = FALSE)
  }
  pmat <- persp3D(z = m$z, x = x, y = y, scale = FALSE, theta = 30,
                  expand = 0.1, panel.first = panelfirst, colkey = FALSE,contour=FALSE)
  XY <- trans3D(x = rep(50, ncol(m$z)), y = y, z = m$z[50,],
                pmat = pmat)
  lines(XY, lwd = 1, lty = 3)
  XY <- trans3D(x = x, y = rep(50, nrow(m$z)), z = m$z[,50],
                pmat = pmat)
  lines(XY, lwd = 1, lty = 3)
  #------------------------------------------------
  return(sum)
}
#------------------------------------------------
ORP<-fun.plot3d(data= dat, var1=dat$ORP[dat$Tratamiento=='Control'],var2=dat$ORP[dat$Tratamiento=='Pennisetum'],tratamiento1= "Control", tratamiento2= "Pennisetum", Variable="ORP")
ORP

pH<-fun.plot3d(data= dat, var1=dat$pH[dat$Tratamiento=='Control'],var2=dat$pH[dat$Tratamiento=='Pennisetum'],tratamiento1= "Control", tratamiento2= "Pennisetum", Variable="pH")
pH

OD<-fun.plot3d(data= dat, var1=dat$OD[dat$Tratamiento=='Control'],var2=dat$OD[dat$Tratamiento=='Pennisetum'],tratamiento1= "Control", tratamiento2= "Pennisetum", Variable="OD")
OD

Temperatura<-fun.plot3d(data= dat, var1=dat$Temperatura[dat$Tratamiento=='Control'],var2=dat$Temperatura[dat$Tratamiento=='Pennisetum'],tratamiento1= "Control", tratamiento2= "Pennisetum", Variable="Temperatura")
Temperatura

Conductividad<-fun.plot3d(data= dat, var1=dat$Conductividad[dat$Tratamiento=='Control'],var2=dat$Conductividad[dat$Tratamiento=='Pennisetum'],tratamiento1= "Control", tratamiento2= "Pennisetum", Variable="Conductividad")
Conductividad
#------------------------------------------------


 
#------------------------------------------------
# Fig. Spearman - Correlations.
#------------------------------------------------
shapiro.test(dat2$pH)
shapiro.test(dat2$OD)
shapiro.test(dat2$Temperatura)
shapiro.test(dat2$ORP)
shapiro.test(dat2$Conductividad)

shapiro.test(dat2$pH[dat2$Tratamiento=='Control'])
shapiro.test(dat2$OD[dat2$Tratamiento=='Control'])
shapiro.test(dat2$Temperatura[dat2$Tratamiento=='Control'])
shapiro.test(dat2$ORP[dat2$Tratamiento=='Control'])
shapiro.test(dat2$Conductividad[dat2$Tratamiento=='Control'])

shapiro.test(dat2$pH[dat2$Tratamiento=='Pennisetum'])
shapiro.test(dat2$OD[dat2$Tratamiento=='Pennisetum'])
shapiro.test(dat2$Temperatura[dat2$Tratamiento=='Pennisetum'])
shapiro.test(dat2$ORP[dat2$Tratamiento=='Pennisetum'])
shapiro.test(dat2$Conductividad[dat2$Tratamiento=='Pennisetum'])

d <- dat2[, c(3, 5:9)]
str(d)
d<-na.omit(d)
names(d)<-c("Tratamiento","pH","OD","Temp","ORP","Cond")
d<-d[order(d$Tratamiento),]
d2<-d[,1]
d<-d[,-1]

hc <- hclust(as.dist(1-cor(d, method='spearman', use='pairwise.complete.obs')))
hc.order <- order.dendrogram(as.dendrogram(hc))
#d <- d[ ,hc]
d[ ,hc.order]
gr <- as.factor(d2)

cols.key <- scales::muted(c('black', 'black'))
cols.key <- adjustcolor(cols.key, alpha.f=1)
pchs.key <- c(19,17)

panel.hist <- function(x, ...) {
  usr <- par('usr'); on.exit(par(usr))
  par(usr=c(usr[1:2], 0, 1.5))
  h <- hist(x, plot=FALSE)
  breaks <- h$breaks
  nB <- length(breaks)
  y <- h$counts; y <- y/max(y)
  rect(breaks[-nB], 0, breaks[-1], y, col='gray', ...)
}
panel.cor <- function(x, y, ...){
  usr <- par('usr'); on.exit(par(usr))
  par(usr=c(0,1,0,1))
  r <- cor(x, y, method='spearman', use='pairwise.complete.obs')
  zcol <- lattice::level.colors(r, at=seq(-1, 1, length=81), col.regions=colorRampPalette(c(scales::muted('red'),'white',scales::muted('blue')), space='rgb')(81))
  ell <- ellipse::ellipse(r, level=0.95, type='l', npoints=50, scale=c(.2, .2), centre=c(.5, .5))
  polygon(ell, col=zcol, border=zcol, ...)
  text(x=.5, y=.5, lab=100*round(r, 2), cex=2, col='black')
  pval <- cor.test(x, y, method='spearman', use='pairwise.complete.obs')$p.value
  sig <- ifelse(pval<0.001,"***",ifelse(pval<0.01,"**", ifelse(pval<0.05,"*","n.s.")))#symnum(pval, corr=FALSE, na=FALSE, cutpoints = c(0, 0.001, 0.01, 0.05, 0.1, 1), symbols = c('***', '**', '*', '.', ' '))
  text(.6, .8, sig, cex=2, col='gray20')
}
panel.scatter <- function(x, y){
  cols<-ifelse(as.numeric(gr)==1,'gray25','gray50')
  points(x, y, col=cols, pch=pchs.key[gr], cex=1.15)
  lines(lowess(x, y))
}
#------------------------------------------------
dat3<-dat2[,c(3, 5:9)]
dat3$Grupo<-dat3$Tratamiento
dat3$Grupo[dat3$Tratamiento=="Control"]<-"C"
dat3$Grupo[dat3$Tratamiento=="Pennisetum"]<-"P"
dat3 %>% ggpairs(.,legend = 1,columns = 2:5,mapping = ggplot2::aes(colour=Grupo),upper = list(continuous = wrap('cor', method = "spearman")),
  lower = list(continuous = wrap("smooth", alpha = 0.5, size=2, pch=c(19)))) +
  theme(legend.position = "bottom") +
  theme_bw()
#------------------------------------------------


#------------------------------------------------
# Fig. Efecto tiempo (ver Cuadro 1)
#------------------------------------------------
str(g$Rep)
g$Rep<-as.character(g$Rep)
g$Time<-as.character(g$Rep)
g$Time[g$Rep=="08:00:00"]<-"8"
g$Time[g$Rep=="09:00:00"]<-"9"
g$Time[g$Rep=="10:00:00"]<-"10"
g$Time[g$Rep=="11:00:00"]<-"11"
g$Time[g$Rep=="12:00:00"]<-"12"
g$Time[g$Rep=="01:00:00"]<-"13"
g$Time[g$Rep=="02:00:00"]<-"14"
g$Time[g$Rep=="12:40:00"]<-"12"
g$Time[g$Rep=="10:40:00"]<-"11"
g$Time<- as.numeric(g$Time)
#------------------------------------------------
layout_matrix = rbind(c(1, 1, 1),
                      c(2, 2, 2))
layout(layout_matrix)
#------------------------------------------------
par(xpd = FALSE,mgp = c(2.5,0.5,0), mar = c(5,7,1,1))
pirateplot(formula = ORP ~  Time + Tratamiento, data = g, 
           main = "", xlab = "", ylab = "ORP (mV, ? SE)",
           ylim=c(-180,160),point.pch=NA,
           bar.f.col=c("white", "gray20"),avg.line.fun=mean,plot=TRUE,
           theme = 4, gl.col = NA,
           avg.line.col = "gray10",
           inf.f.col = "gray10",
           jitter.val = 0.09,
           avg.line.o = 1,bar.b.o=1,
           avg.line.lwd=3,
           inf.lwd=1,
           inf.method="se")
#------------------------------------------------
par(xpd = FALSE,mgp = c(2.5,0.5,0), mar = c(5,7,1,1))
pirateplot(formula = Temperatura ~  Time + Tratamiento, data = g, 
           main = "", xlab = "", ylab = "Temperatura (?C, ? SE)",
           ylim=c(20,30),point.pch=NA,
           bar.f.col=c("white", "gray20"),avg.line.fun=mean,plot=TRUE,
           theme = 4, gl.col = NA,
           avg.line.col = "gray10",
           inf.f.col = "gray10",
           jitter.val = 0.09,
           avg.line.o = 1,bar.b.o=1,
           avg.line.lwd=3,
           inf.lwd=1,
           inf.method="se")
#------------------------------------------------

datc<-dat2[dat2$Tratamiento=="Control",]
datp<-dat2[dat2$Tratamiento=="Pennisetum",]

#Correlación
#OD y pH

ODc<-(datc$OD)
shapiro.test(ODc)
pHc<-(datc$pH)
shapiro.test(pHc) 
cor.test(ODc, pHc, method = "spearman") 

ODp<-(datp$OD)
shapiro.test(ODp)
pHp<-(datp$pH)
shapiro.test(pHp) 
cor.test(ODp, pHp, method = "s") 

ODg<-(dat2$OD)
shapiro.test(ODg)
pHg<-(dat2$pH)
shapiro.test(pHg) 
cor.test(ODg, pHg, method = "s") 

#ORP y Temperatura

ORPc<-(datc$ORP)
shapiro.test(ORPc) 
Tempc<-(datc$Temperatura)
shapiro.test(Tempc)
cor.test(ORPc, Tempc, method = "s") 

ORPp<-(datp$ORP)
shapiro.test(ORPp) 
Tempp<-(datp$Temperatura)
shapiro.test(Tempp) 
cor.test(ORPp, Tempp, method = "s") 

ORPg<-(dat2$ORP)
shapiro.test(ORPg)
Tempg<-(dat2$Temperatura)
shapiro.test(Tempg) 
cor.test(ORPg, Tempg, method = "s") 

#ORP + pH

ORPC <- c(datc$ORP)
shapiro.test(ORPC)
pHC <- c(datc$pH)
shapiro.test(pHC) 
cor.test(ORPC, pHC, method = "spearman") 

ORPP <- (datp$ORP)
shapiro.test(ORPP)
pHP <- (datp$pH)
shapiro.test(pHP)
cor.test(ORPP, pHP, method = "spearman")

cor.test(ORPg, pHg , method = "spearman") 

#Temperatura con OD

Tempc<-(datc$Temperatura)
shapiro.test(Tempc) 
ODc<-(datc$OD)
shapiro.test(ODc) 
cor.test(Tempc, ODc, method = "s")

Tempp<-(datp$Temperatura)
shapiro.test(Tempp) 
ODp<-(datp$OD)
shapiro.test(ODp) 
cor.test(Tempp, ODp, method = "s") 

Tempg<-(dat2$Temperatura)
shapiro.test(Tempg) 
ODg<-(dat2$OD)
shapiro.test(ODg) 
cor.test(Tempg, ODg, method = "s") 

#TEMPERATURA + pH

temperaturac <- (datc$Temperatura)
shapiro.test(temperaturac) 
pHC <- (datc$pH)
shapiro.test(pHC) 
cor.test(temperaturac, pHC, method = "spearman")  

temperaturap <- (datp$Temperatura)
shapiro.test(temperaturap) 
pHP <- (datp$pH)
shapiro.test(pHP) 
cor.test(temperaturap, pHP, method = "s") 

temperaturag <- (dat2$Temperatura)
shapiro.test(temperaturag) 
pHG <- (dat2$pH)
shapiro.test(pHG) 
cor.test(temperaturag, pHG, method = "spearman") 

#Conductividad y pH

condg<-(dat2$Conductividad)
shapiro.test(condg) #no normal
cor.test(condg, pHG, method = "s") 

condc<-(datc$Conductividad)
shapiro.test(condc) #no normal
cor.test(condc, pHC, method = "s") 

condp<-(datp$Conductividad)
shapiro.test(condp) #no normal
cor.test(condp, pHP, method = "s") 

#Conductividad y OD

cor.test(condg, ODG, method = "s") 
cor.test(condc, ODC, method = "s") 
cor.test(condp, ODp, method = "s") 

#Conductividad y temperatura

cor.test(condg, temperaturag, method = "s") 
cor.test(condc, temperaturac, method = "s") 
cor.test(condp, temperaturap, method = "s")

#Conductividad y ORP

cor.test(condg, ORPG, method = "s") 
cor.test(condc, ORPC, method = "s") 
cor.test(condp, ORPP, method = "s")

#OD + ORP

ODC <-(datc$OD)
shapiro.test(ODC) 
ORPC <-(datc$ORP)
shapiro.test(ORPC) 
cor.test(ODC, ORPC, method = "spearman") 

ODP <-(datp$OD)
shapiro.test(ODP) 
ORPP <-(datp$ORP)
shapiro.test(ORPP) 
cor.test(ODP, ORPP, method = "spearman") 

ODG <- (dat2$OD)
shapiro.test(ODG) 
ORPG <- (dat2$ORP)
shapiro.test(ORPG) 
cor.test(ODG, ORPG, method = "spearman") 

#ANDEVAS COLUMNAS

odcontrol<- aov(OD~Columna, data = datc)
shapiro.test(odcontrol$residuals)
tapply(datc$OD, datc$Columna, length)
fligner.test(datc$OD~datc$Columna)
kruskal.test(datc$OD~datc$Columna)

odpeni<- aov(OD~Columna, data=datp)
shapiro.test(odpeni$residuals)
tapply(datp$OD, datp$Columna, length)
bartlett.test(datp$OD~datp$Columna)
summary(odpeni)

orpcontrol <- aov(ORP~Columna, data = datc)
shapiro.test(orpcontrol$residuals)             
tapply(datc$ORP, datc$Columna, length)
fligner.test(datc$ORP ~ datc$Columna)
kruskal.test(datc$ORP ~ datc$Columna)


orppennisetum <- aov(ORP~Columna, data = datp)
shapiro.test(orppennisetum$residuals)  
tapply(datp$ORP, datp$Columna, length)
fligner.test(datp$ORP ~ datp$Columna)
kruskal.test(datp$ORP ~ datp$Columna)


pHc<-aov(pH~Columna,data = datc)
shapiro.test(pHc$residuals)
tapply(datc$pH,datc$Columna,length)
fligner.test(pH~Columna, data = datc)
kruskal.test(pH~Columna, data = datc)

pHp<-aov(pH~Columna,data = datp)
shapiro.test(pHp$residuals)
tapply(datp$pH,datp$Columna,length)
fligner.test(pH~Columna, data = datp)
kruskal.test(pH~Columna,data = datp)

Tc<-aov(datc$Temperatura~datc$Columna)
shapiro.test(Tc$residuals)
tapply(datc$Temperatura,datc$Columna,length)
fligner.test(Temperatura~Columna,data = datc)
kruskal.test(datc$Temperatura~datc$Columna)

Tp<-aov(datp$Temperatura~datp$Columna)
shapiro.test(Tp$residuals)
tapply(datp$Temperatura,datp$Columna,length)
fligner.test(Temperatura~Columna,data = datp)
kruskal.test(datp$Temperatura~datp$Columna)

conduc<-aov(datc$Conductividad~datc$Columna)
shapiro.test(conduc$residuals)
tapply(datc$Conductividad,datc$Columna,length)
fligner.test(datc$Conductividad~datc$Columna)
kruskal.test(datc$Conductividad~datc$Columna)

condup<-aov(datp$Conductividad~datp$Columna)
shapiro.test(condup$residuals)
tapply(datp$Conductividad,datp$Columna,length)
fligner.test(datp$Conductividad~datp$Columna)
kruskal.test(datp$Conductividad~datp$Columna)

#ANDEVAS FILAS 

tempfilacontrol<- aov(Temperatura~Fila, data = datc)
shapiro.test(tempfilacontrol$residuals) 
tapply(datc$Temperatura, datc$Fila, length)
fligner.test(datc$Temperatura~datc$Fila)
kruskal.test(datc$Temperatura~datc$Fila)


tempfilapen<- aov(Temperatura~Fila, data=datp)
shapiro.test(tempfilapen$residuals)
tapply(datp$Temperatura, datp$Fila, length)
fligner.test(datp$Temperatura~datp$Fila)
kruskal.test(datp$Temperatura~datp$Fila)
pairwise.wilcox.test(datp$Temperatura, datp$Fila, p.adjust.method = "b", exact=F)


orpcontrol <- aov(ORP~Fila, data = datc)
shapiro.test(orpcontrol$residuals)
tapply(datc$ORP, datc$Fila, length)
bartlett.test (datc$ORP ~ datc$Fila)
kruskal.test(datc$ORP ~ datc$Fila)
pairwise.wilcox.test(datc$ORP, datc$Fila, p.adjust.method = "b", exac=F)


orppenni <- aov(ORP~Fila, data = datp)
shapiro.test(orppenni$residuals)
tapply(datp$ORP, datp$Fila, length)
fligner.test(datp$ORP ~ datp$Fila)
kruskal.test(datp$ORP ~ datp$Fila)
pairwise.wilcox.test(datp$ORP, datp$Fila, p.adjust.method = "b", exac=F)


condc<-aov(Conductividad~Fila, data = datc)
shapiro.test(condc$residuals)
tapply(datc$Conductividad, datc$Fila, length)
leveneTest(datc$Conductividad, datc$Fila)
kruskal.test(datc$Conductividad, datc$Fila)
pairwise.wilcox.test(datp$Conductividad, datp$Fila, p.adjust.method = "b", exact=F)

condp<-aov(Conductividad~Fila, data = datp)
shapiro.test(condp$residuals)
tapply(datp$Conductividad, datp$Fila, length)
bartlett.test(Conductividad~Fila, data = datp)
kruskal.test(Conductividad~Fila, data = datp)
pairwise.wilcox.test(datp$Conductividad, datp$Fila, p.adjust.method = "b", exact=F)


pHF<-aov(pH~Fila, data=datp)
shapiro.test(pHF$residuals)  
tapply(datp$pH, datp$Fila, length)  
fligner.test(pH~Fila, data=datp)
kruskal.test(datp$pH~datp$Fila)
pairwise.wilcox.test(datp$pH, datp$Fila, p.adjust.method = "b", exact=F)


pHFC<- aov(pH~Fila, data=datc)
shapiro.test(pHFC$residuals)  
tapply(datc$pH, datc$Fila, length)  
fligner.test(pH~Fila, data=datc)
kruskal.test(datc$pH~datc$Fila)
pairwise.wilcox.test(datc$pH, datc$Fila, p.adjust.method = "b", exact=F)


ODF<-aov(Conductividad~Fila, data = datc)
shapiro.test(ODF$residuals)
tapply(datc$Conductividad, datc$Fila, length)
fligner.test(Conductividad~Fila, data = datc)
kruskal.test(datc$Conductividad~datc$Fila)
pairwise.wilcox.test(datp$Conductividad, datp$Fila, p.adjust.method = "b", exact=F)

ODPF<-aov(Conductividad~Fila, data = datp)
shapiro.test(ODPF$residuals)
tapply(datp$Conductividad, datp$Fila, length)
fligner.test(Conductividad~Fila, data = datp)
kruskal.test(datp$Conductividad~datp$Fila)
pairwise.wilcox.test(datp$Conductividad, datp$Fila, p.adjust.method = "b", exact=F)