In jumanbar/siabox: Datos y Funciones Para Trabajar con el SIA

# Este código define opciones generales para los "chunks" de código ejecutable
# del documento. Las tres opciones elegidas evitan que se muestre el código R en
# el documento final, así como distintos mensajes o advertencias que se puedan
# generar directa o indirectamente.
knitr::opts_chunk$set(echo=FALSE, warning=FALSE, message=FALSE)

# Este chunk es para cargar paquetes necesarios, pero NO para instalarlos. Es
# decir, no incluir comandos del tipo:
#
# > install.packages(...)
#
# > remotes::install_github(...)
#
# Sólo usar llamados a require o library:
library(siabox)
library(tidyverse)
library(patchwork)

## Grupos de parámetros: fisicoquímicos, nutrientes y seleccionados para hacer
## ciertas comparaciones:
#
# Temperatura, OD, Conductividad, pH, Turbidez, ST:
ids_fq <- c(T=2032, OD=2017, Conduc=2009, pH=2018, Turbidez=2035, ST=2028)
# NO2, PT, Fosfato (PO4), NH3L, NO3, NT:
ids_nut <- c(NO2=2102, PT=2098, PO4=2097, NH3L=2091, NO3=2099, NT=2101)
# PT, NT, ST, Conductividad:
ids_comp <- c(PT=2098, NT=2102, ST=2028, Conduc=2009)
# Otros:
ids_otr <- c(NH4 = 2090, Clo_a = 2000, ColifTrm_MF=2111)

# Aquí es donde se importan los datos a usar en el informe. Hay dos opciones:
#
# 1. Traer datos externos. Para la mayoría de las funciones de siabox es
# necesario que los datos tengan cieras características. Bajar los datos de iSIA
# garantiza que así sea:

# d <- readRDS("extraccion_20210205.rds")

# 2. Usar los datos ya incluidos en el paquete (datos_sia o datos_sia_sed). Esta
# es una buena opción para hacer pruebas generales de código, antes de usar
# datos actualizados, descargados desde el iSIA. La función filtrar_datos
# realiza el mismo filtrado que ocurre en la interfaz del iSIA:
d <- datos_sia %>%
  filtrar_datos(id_programa = 10L,
                id_matriz = 6L,
                rango_fechas = c('2019-01-01', '2019-12-31'),
                id_parametro = c(ids_fq, ids_nut, ids_comp, ids_otr),
                tipo_punto_id = 1L)

# La tabla etiqs se define aquí para ajustar el formato de las etiquetas que van
# a aparecer en los ejes de los gráficos. Si no se especifican así las
# etiquetas, las funciones de siabox usarán las columnas param y
# uni_nombre (ver datos_sia) para hacer las etiquetas.
etiqs <- t_eti_add(id_parametro = c(2032L, 2017L, 2009L, 2018L, 2035L, 2028L,
                                    2102L, 2098L, 2097L, 2091L, 2099L, 2101L,
                                    2000L, 2111L),
                   etiqueta = c(expression('Temperatura (ºC)'),
                                expression('OD (mg/L)'),
                                expression('Cond ('*mu*'S/cm)'),
                                expression('pH'),
                                expression('Turbidez (UNT)'),
                                expression('ST (mg/L)'),
                                expression('NT (mg/L)'),
                                expression('PT ('*mu*'g/L)'),
                                expression('PO'[4]*' ('*mu*'g-P/L)'),
                                expression('NH'[3]*' ('*mu*'g-N/L)'),
                                expression('NO'[3]*' (mg-N/L)'),
                                expression('NO'[2]*' (mg-N/L)'),
                                expression('Clo-a ('*mu*'g/L)'),
                                'ColiTrm_MF (UFC/100mL)'))

Gráficos y análisis de datos básicos

En los informes anuales primero se presenta por separando cada subcuenca/zona/cauce; en el caso de laguna Merín, esa división se hace para cada río o arroyo y luego se integra esa iformación con gráficos que comparan las cuencas entre sí.

Análisis de datos

Consideraciones:

• Las estaciones deben ordenarse desde las nacientes hasta la desembocadura

• el PT y PO4 deben estar en la misma unidad; en el caso del IET es necesario el PT en μg/L

• NT y todos los demás, en general están en mg/L

Los cáclulos que deben realizarse antes de empezar son:

• Se calcula NH3 a partir de NH4 según Canadian Water Quality Guidelines for the Protection of Aquatic Life (ver Ammonia):

$$NH_3 (\mu g) = \frac{1000 \times NH_4 (mg)}{(1+10^{(-pH+0.0901821+\frac{2729.92}{T(º C)+273.15} )})}$$

• El valor de IET, se calcula solamente para las estaciones lóticas, con la siguiente ecuación:

$$IET = 10 \times (6-\frac{0.42-0.36 \times ln(\overline{PT (\mu g/L)})}{ln(2)})-20$$

• El valor $\overline{PT (\mu g/L)}$ representa la media geométrica anual (no confundir con la media aritmética) de los valores del Fósforo Total (en microgramos por litro) obtenidos para cada estación, en las distintas campañas.

El IET se normaliza en una escala de 0 a 100, dividida en seis rangos de estados tróficos, desde la Ultraoligotrofía a la Hipereutrofía. Se considera 1 decimal después de la coma para pasar de un estado a otro (por ejemplo: 59,1 corresponde a un estado eutrófico). Los seis estados se indican con colores.

| Clasificación | Criterio | |:------------------|:-------------:| | Ultraoligotrófico | IET ≤ 47 | | Oligotrófico | 47 < IET ≤ 52 | | Mesotrófico | 52 < IET ≤ 59 | | Eutrófico | 59 < IET ≤ 63 | | Supereutrófico | 63 < IET ≤ 67 | | Hipereutrófico | IET > 67 |

La interpretación del IET y el valor final a publicar en el informe debe ser realizado con cautela y tomando en cuenta la globalidad del ecosistema.

d$periodo <- as.factor(stringr::str_to_title(d$periodo))

d <- amoniaco_libre_add(d)

Gráficos

Fisicoquímicos y nutrientes

Para una fase exploratoria se grafican todos los nutrientes, pero luego debe evaluarse cuidadosamente el número de valores cuantificados que hay en cada caso y volver a pensar si corresponde graficar o si es mejor colocar los valores en una tabla.

Olimar

dx <- dplyr::filter(d, nombre_subcuenca_informes == "Olimar")

Fisicoquímico

g_mes_pto_all(dx, ids_fq, ncol = 3, t_eti = etiqs)

Nutrientes

g_mes_pto_all(dx, ids_nut, ncol = 3, t_eti = etiqs)

Comparación entre estaciones

Los puntos representan la media o promedio y las líneas +/- desvío estándar

g_est_dsv_all(dx, ids_comp, ncol = 2, t_eti = etiqs)

Cebollatí

dx <- dplyr::filter(d, nombre_subcuenca_informes == "Cebollatí")

Fisicoquímicos

g_mes_pto_all(dx, ids_fq, ncol = 3, t_eti = etiqs)

Nutrientes

g_mes_pto_all(dx, ids_nut, ncol = 3, t_eti = etiqs)

Comparación entre estaciones

g_est_dsv_all(dx, id_parametro = ids_comp, ncol = 2, t_eti = etiqs)

Comparación entre las subcuencas

Fisicoquímicos y nutrientes

c1 <- g_cue_box(d, id_param = 2009L, t_eti = etiqs)
c2 <- g_cue_box(d, id_param = 2028L, t_eti = etiqs)
c3 <- g_cue_box(d, id_param = 2098L, t_eti = etiqs)
c4 <- g_cue_box(d, id_param = 2101L, t_eti = etiqs)

(c1 + c2) / (c3 + c4) +
  plot_annotation(tag_levels = 'A')

Cloforila-A

En este caso hay un solo valor de clorofila por fecha de muestreo. Las líneas entrecortadas naranja y roja indican los niveles de alerta 1 y 2 de la OMS respectivamente.

d %>%
  dplyr::filter(id_parametro == 2000L) %>%
  ggplot() +
  aes(mes, valor, color = codigo_pto) +
  geom_jitter(alpha = 0.5, width = 0.05) +
  labs(x='Fecha', y = 'Cloa (μg/L)', colour="Estación") +
  geom_hline(yintercept = 10, linetype="twodash", colour='orange', size=.3) +
  geom_hline(yintercept = 50, linetype="twodash", colour='red', size=.3) +
  theme_bw() +
  facet_grid(rows = vars(nombre_subcuenca_informes)) +
  theme(panel.grid.minor = element_blank(),
        text = element_text(size = 10),
        axis.text.x = element_text(angle=0, vjust=1)) +
  guides(colour = guide_legend("Estación"))

Coliformes

d %>%
  dplyr::filter(id_parametro == 2111L) %>%
  ggplot() +
  aes(mes, valor, color = codigo_pto) +
  geom_jitter(alpha = 0.7, width = 0.1) +
  labs(x='Mes', y = 'Colif. Termot. (UFC/100 ml)', colour="Estación")+
  geom_hline(yintercept = 2000, linetype="longdash", colour='red', size=0.3)+
  theme_bw()+
  facet_grid(rows = 'nombre_subcuenca_informes')+
  theme(panel.grid.minor = element_blank(),
        text = element_text(size = 10),
        axis.text.x = element_text(angle=0, vjust=1))+
  guides(colour = guide_legend("Estación"))

NOTA: está pendiente crear una función para estos gráficos (Clo-a y Colif)

IET

# Esta selección tal vez sea incorrecta, pero se basa en la idea de que el IET
# debería usarse solamente para estaciones lóticas (ie: ríos, arroyos):
estaciones_loticas <- c("C0", "C1", "C2", "O2", "SL1", "SM1", "SM2", "T1", "T2",
                        "T3", "Y3", "C3", "O1", "Y1", "Y2")
dfiet <- d %>%
  dplyr::filter(codigo_pto %in% estaciones_loticas) %>%
  iet_tabla()

g_iet_pto(dfiet)

knitr::kable(dfiet, format = "markdown",
             col.names = c("Estación", "IET", "Clasificación"))

jumanbar/siabox documentation built on April 25, 2022, 1:37 p.m.