Muestreo de material particulado y meteorología

Localizada en una región continental intertropical en el norte-centro de Venezuela, la Gran Caracas es una gran urbe ubicada al sur del parque nacional El Ávila, formación montañosa perteneciente a la Cordillera de la Costa. La Gran Caracas está conformada por Caracas y otras ciudades satélites. Ocupa un área de 4704 km² y tiene 4.25 millones de habitantes (^{INE 2018}). Esta urbe está cruzada por vías rápidas que presentan alto tráfico vehicular en horas críticas y con un número importante de industrias pequeñas y medianas.

Al suroeste de la Gran Caracas se encuentra el Valle de Sartenejas (Fig. 2), a una altitud entre 1200 y 1400 msnm. Este valle ocupa un área de 96 ha, de las cuales 40 ha pertenecen al campus de la Universidad Simón Bolívar. La universidad está rodeada de 255 ha de áreas naturales, de las cuales 45 ha son plantaciones con alta densidad de bosque de pinos; el resto corresponden a bosque nativo, matorrales secundarios y una zona importante de sabana con arbustos diversos y escasos pinos que presentan alta incidencia de incendios (^{GBUSB 2018a}). El Valle de Sartenejas tiene un microclima más fresco que el de la Gran Caracas y habitualmente presenta mayor frecuencia de lluvias.

Fig. 2 Localización de la estación de muestreo AireUSB 

El campus universitario cuenta con una estación de muestreo de material particulado que opera junto con una estación meteorológica, llamadas en conjunto AireUSB. Esta estación está instalada en la azotea de uno de los edificios del campus a una altura de 9 m, para asegurar que el material particulado sea representativo del aire del sector (coordenadas UTM 10.412352 y -66.883558 m). El sector se percibe poco contaminado por estar rodeado de espacios naturales, a pesar de que a una distancia aproximada de 2 km al suroeste se encuentra una autopista con alto tránsito vehicular, pues es la vía de salida de la Gran Caracas hacia el occidente del país, y a una distancia 3.5 km al noreste se emplaza una zona industrial-residencial. La vía de acceso a la universidad presenta alto tráfico vehicular en horas críticas (Fig. 2).

Para la recolección del PM en sus diferentes tamaños se empleó un muestreador de alto volumen (Hi-Vol, por sus siglas en inglés), el cual es recomendado por la Agencia de Protección Ambiental de los Estados Unidos (USEPA, por sus siglas en inglés) para la recolección de PTS cuando se opera con grandes volúmenes de aire (^{EPA 1999a}). El equipo se calibró siguiendo las instrucciones del manual de operación, el cual se basa en procedimientos de la USEPA para medición de PM de aire exterior con equipos de alto volumen.

Al Hi-Vol se le acopló un impactador en cascada (ImpC) modelo Inc. Series 230 (Tisch Environmental). En el ImpC las partículas de PTS se segregan según su diámetro aerodinámico en seis rangos de tamaños, por medio de cinco bandejas de aluminio ranuradas que dividen las PTS en los siguientes rangos: bandeja 1, > 7.2 μm; bandeja 2, 7.2-3.0 μm; bandeja 3, 3.0-1.5 μm; bandeja 4, 1.5-0.95 μm, y bandeja 5, 0.95-0.49 μm. En el fondo del ImpC se coloca como bandeja 6 un filtro < 0.49 μm. Se usaron filtros de fibra de vidrio.

Se destaca que una de las ventajas de usar un impactador en cascada es la posibilidad de combinar las fracciones recolectadas en las diferentes bandejas según el objetivo de la investigación (^{Ehrlich et al. 2007}, ^{Srivastava et al. 2008}, ^{Campos et al. 2015}).

Antes de cada muestreo los filtros se someten a un tratamiento térmico en una estufa a temperatura de 120 ºC por 24 h. Luego se estabilizan en un desecador por 24 h más, y posteriormente se pesan de forma independiente en condiciones de humedad y temperatura ambiente conocidas. Para determinar la masa de material particulado en cada filtro, éstos se pesan antes y después de cada muestreo (Fig. 3). La concentración de material particulado se calcula como la masa total de partículas recolectadas dividida entre el volumen de aire muestreado (en μg/m³). El volumen total de aire muestreado se calcula a partir del flujo medido y el tiempo de muestreo. Este valor se corrige a condiciones conocidas de presión y temperatura (^{WMO 1978}, ^{Lodge 1988}, ^{EPA 1999b}).

Fig. 3 Esquema para la preparación, colocación y extracción de filtros 

Se tomaron muestras de material particulado cada seis días del 23 de agosto de 2014 al 29 de octubre de 2015 para un total de 85 muestras, durante la época de lluvia y de sequía. Del 23 de agosto al 28 de octubre de 2014 el muestreo se realizó por 24 h cada seis días; del 3 de noviembre de 2014 al 24 de septiembre de 2015 el muestreo se realizó por 48 h cada seis días a partir del encendido del equipo, y del 29 de septiembre al 29 de octubre el muestreo fue de 72 h cada seis días a partir del encendido del equipo. Las fechas reportadas corresponden al día de encendido del equipo. El ajuste del periodo de muestreo se debió a que la masa del PM recolectada por los filtros resultaba estar justo en el umbral mínimo de detección de la balanza disponible, lo cual motivó que se incorporaran los periodos de muestreo de 48 y 72 h, los cuales posteriormente se normalizaron a concentraciones de 24 h, tomando en cuenta tiempo de recolección y volumen de aire.

El cuadro I presenta la concentración promedio de PM de cada bandeja del impactador para cada periodo de muestreo (24, 48 y 72 h), la desviación estándar de las mismas, el valor máximo y mínimo de concentración de cada bandeja, y la cantidad de muestras recolectadas, las descartadas y las utilizadas en esta investigación. En dicho cuadro se observa que las mayores concentraciones tienden a presentarse en las bandejas 6 y 2. De un total de 85 muestras recolectadas se descartaron 21 por manejo inadecuado o contaminación de las muestras y/o problemas operativos con los equipos, resultando un total de 64 muestra útiles: 59 para la construcción del modelo y cinco reservadas para la validación.

En la estación de muestreo AireUSB los datos meteorológicos fueron recolectados con un equipo Davis Instruments 2010, que operó las 24 h del día durante el tiempo de muestreo (agosto de 2014-octubre de 2015). El periodo de toma de datos para la velocidad y dirección del viento fue de 1 h, para la ocurrencia de lluvia fue diariamente y para el resto de las variables meteorológicas cada 10 min. Luego de un análisis se descartaron algunos datos meteorológicos. En el cuadro II se resume la cantidad de datos meteorológicos esperados y los disponibles.

A partir del análisis de los registros meteorológicos puede conocerse que la temperatura promedio mensual en el Valle de Sartenejas durante el periodo de estudio varió de 19.7 a 22.8 ºC. La temperatura mínima se mantuvo en el rango de 15.4-18.5 ºC y la máxima fue de 25.0 a 28.0 ºC. La humedad relativa varió de 73 a 81 % y la radiación solar de 250 a 347 cal/cm². La velocidad del viento promedio fue de 1.57 m/s con 37 % de vientos de calma. Los valores de dirección del viento se descartaron por estar incompletos (Cuadro II).

La estación meteorológica La Mariposa, administrada por el Instituto Nacional de Meteorología e Hidrología (INAMEH), ubicada en las coordenadas UTM 10.40417 y -66.92500 m, a 5.6 km del Valle de Sartenejas, publica registros horarios de velocidad y dirección del viento para 2015. La figura 4 muestra la rosa de vientos en la estación La Mariposa. En ella se observan vientos prevalentes del este con máximos en el este-noreste y sureste, velocidades prevalentes entre 2.5 y 3.5 m/s con 21 % de vientos de calma. En el estudio meteorológico presentado por ^{Goldbrunner (1984)} para el periodo 1951-1970 se indica que el viento dominante en Caracas provenía del este con predominio de vientos del sureste, mostrando un máximo en la dirección este-sureste, velocidad predominante de 2 m/s y 44.3 % de vientos de calma.

Fig. 4 Rosa de vientos en la estación La Mariposa para el año 2015 

Esta estación también ha publicado registros mensuales de precipitación y temperatura para el periodo 1970-2000 (^{MINAMB 2006}). La figura 5 muestra el climograma de la estación La Mariposa para el periodo 1970-2000, en el cual se observa una clara época de lluvia (mayo-noviembre) y otra de sequía (diciembre-abril), mostrando un máximo de lluvias en octubre (125 mm) y un mínimo en febrero (10.6 mm). Este comportamiento corresponde al descrito por ^{Goldbrunner (1984)} para el periodo 1951-1970 en Caracas.

Fig. 5 Climograma de la estación La Mariposa para el periodo 1970-2000 

Información secundaria

La información secundaria consiste en el histórico mensual de precipitación, ocurrencia de incendios, día de inicio del muestreo y días laborables durante el muestreo. El histórico mensual de la precipitación es el promedio de la precipitación mensual acumulada en la estación meteorológica La Mariposa para el periodo 1970-2000, con 328 datos de un total de 360 esperados. La ocurrencia de incendios forestales en los alrededores de la estación AireUSB se obtuvo por medio de los informes técnicos de los bomberos de la Universidad Simón Bolívar y de la agrupación BosquesUSB. Se conoce que durante agosto de 2014 y noviembre de 2015 se reportaron ocho incendios forestales en las zonas verdes aledañas a la universidad (^{GBUSB 2018b}). Las variables día de inicio de muestreo y días laborales se establecieron a partir del calendario; para el caso de días laborales se indica si los días muestreados coinciden con algún día del fin de semana o días festivos.

Selección de variables: dependiente e independientes

La variable dependiente (VP) se define como la concentración de material particulado superfino (PM₁) recolectado por la suma de la masa obtenida en los filtros de las bandejas 5 y 6. Esta masa corresponde al PM cuyo diámetro aerodinámico es menor a 0.95 µm, el cual es referido en esta investigación como PM₁. Esta variable se mide en una escala continua (valores de 0 a infinito). En el cuadro III se describe la variable dependiente.

El marco de referencia para la selección de variables independientes se estableció a partir de publicaciones producto de investigaciones similares (^{Baek et al. 1997}, ^{Chaloulakou et al. 2003}, ^{Silva et al. 2003}, ^{Botero et al. 2004}, ^{Cassidy et al. 2007}, ^{Ramos-Herrera et al. 2010}). La concentración de PM en la atmósfera está relacionada con las concentraciones de otros tamaños de partículas (llamadas variables ambientales), la velocidad y dirección de viento, la radiación solar, la humedad relativa, la temperatura y la precipitación (llamadas variables meteorológicas) y con el tráfico vehicular y la ocurrencia de incendios (llamadas variables de eventos particulares). La selección de variables independientes que se incluye en un modelo es discrecional, pero siempre se basa en los conocimientos físicos y químicos del comportamiento del PM en la atmósfera.

Se seleccionaron once variables independientes: dos ambientales, cinco meteorológicas y cuatro relacionadas con eventos particulares. El detalle de las mismas se presenta en el cuadro IV.

Se utiliza el análisis de correlación de Pearson (r) entre la VP y las VI con la intención de identificar aquellas VI relacionadas de forma significativa (p < 0.1) con la variable dependiente (PM₁). Para corroborar la correlación entre dos variables de diferente escala se aplica adicionalmente la prueba t de Student: VP medida en una escala continua y VI medidas en escala discreta (ocurrencia de incendios, ocurrencia de lluvia, día de inicio del muestreo y día laboral).

Construcción del modelo multivariante

El análisis de RLM es una técnica estadística utilizada para estudiar la relación entre un conjunto de más de dos variables independientes (VI), el cual tiene el propósito de averiguar en qué medida la variable dependiente (VP) puede explicarse por las VI. Si se asume una VP (denotada como Yi), su mejor predictor lineal a partir de otras VI (Xi) tendrá la siguiente forma (^{Peña 2002}, ^{Cohen et al. 2003}):

(1)

(2)

donde β ₀ , β ₁ ,... β _k - B ₀ , B ₁ , B ₂ ,…B _k son parámetros fijos desconocidos; X _1i ,...,X _ki son variables independientes cuyos valores son fijados por el investigador, y ∈ i es una variable aleatoria inobservable (error aleatorio).

El modelo A permite identificar qué variable tiene la contribución más importante sobre la varianza total (^{Cohen et al. 2003}). Los coeficientes β están basados en puntuaciones típicas y son directamente comparables entre sí. Indican la cantidad de cambio que se producirá en la VP por cada cambio de una unidad en la correspondiente variable independiente (manteniendo constantes el resto de VI). Una variable tiene tanta más importancia en la ecuación de regresión cuanto mayor sea, en valor absoluto, el coeficiente de regresión estandarizado (β) (^{Cohen et al. 2003}).

El modelo B se construye con coeficientes no estandarizados B de cada variable. El coeficiente B indica el cambio medio que corresponde a la VP por cada unidad de cambio de la variable independiente. Los coeficientes no estandarizados B se leen en las mismas unidades físicas en las cuales se mide cada variable.

La RLM se basa en una serie de supuestos que la base de datos debe cumplir para que el modelo pueda considerarse válido. Estos supuestos son la escala de medida (la VP debe ser una variable en escala continua); el número de variables (debe haber al menos dos VI en la base de datos); el tamaño de muestra o la cantidad de VI (se debe contar con al menos 10 datos de cada VI); la linealidad (debe existir una relación lineal entre la VP y cada una de las VI por separado, y entre la VP y las VI colectivamente, lo cual se prueba con correlaciones); la independencia de los errores (se debe tener independencia de las observaciones o, lo que es igual, independencia de los residuos); la normalidad (los datos deben ser normales, lo cual se prueba analizando el histograma de los residuos o el gráfico de probabilidad P-P acumulada observada/esperada según la normal); la homocedasticidad (significa que la varianza de la distribución de probabilidad del error [€] es constante para todos los valores de la variable independiente, lo cual se prueba usando el diagrama de dispersión de pronósticos tipificados por residuos tipificados); la multicolinealidad (dos o más VI no deben tener una correlación mayor a r > 0.9), y los valores atípicos (la base de datos no debe mostrar más del 5 % de datos atípicos, celdas vacías o datos perdidos) (^{Cohen et al. 2003}).

Validación del modelo

El indicador de error de sesgo promedio (MBE, por sus siglas en inglés) se usa para validar los resultados del modelo contra datos experimentales. El indicador MBE representa el grado de correspondencia entre una predicción media (Pi) y una observación media (Oi) y se utiliza para describir si un modelo sobrestima (valor positivo) o subestima (valor negativo) una observación (ecuación 3). Valores bajos de MBE son aceptados (^{Moustris et al. 2010}). Para la validación del modelo se tomaron las muestras de PM recolectadas para el periodo de 72 h, reservadas para tal fin.

(3)

donde n es el número de puntos de datos, Oi el dato observado y Pi el dato predicho.

Construcción del modelo de regresión lineal múltiple para PM₁

Se construyó un modelo de regresión lineal usando solamente las seis variables que se relacionan de forma significativa con PM₁ (p < 0.1) y 52 observaciones (se descartaron siete observaciones para un mejor cumplimiento de los supuestos de la RLM). El cuadro VII muestra el modelo que mejor explica la varianza de PM₁ con el menor número de VI. Al modelo se le comprobaron los supuestos de homocedasticidad, independencia de los errores y normalidad para verificar que los datos cumplieran con los criterios mínimos para ser utilizados en las predicciones. La varianza de PM₁ se explica en más de un 75 % (R² ajustado = 0.759; p < 0.000) a partir de tres VI: PMcuasifino (PM_3.0-0.95), precipitación histórica mensual (histórico del promedio mensual de precipitación acumulada) e IncFor. El análisis mostró que la RLM tenía valores R² significativos, siendo este resultado similar al encontrado por ^{Thomas y Jacko (2007)}), quienes reportan un R² no ajustado de 80 % para una regresión lineal relacionada con predicción de PM_2.5 a través de cuatro VI. Asimismo, resultó ser mayor al R² ajustado publicado por ^{Nazif et al. 2016}) para un análisis de regresión lineal que predice PM₁₀ a partir de diferentes combinaciones de VI.

Con la información del cuadro VII se construyen el modelo A (ecuación 4) y el modelo B (ecuación 5). El modelo A sólo indica la importancia relativa de cada variable.

(4)

Respecto de este modelo, se tiene que:

Al aumentar la concentración de PMcuasifino (PM_3.0-0.95) aumenta la concentración de PM₁, lo que coincide con los resultados obtenidos en el análisis correlacional. Esta variable es la de mayor importancia al momento de predecir la concentración de PM₁ (véase β en el cuadro VII).

La asociación entre IncFor y PM₁ señala que al ocurrir los primeros se agrega un valor de 0.244 unidades al valor predicho de la PM₁.

El signo positivo de la variable que mide el comportamiento histórico del promedio mensual de precipitación acumulada, muestra que los meses con mayor registro histórico de precipitaciones se asocian con mediciones mayores de PM₁. Este resultado lo explica la variación en la distribución de la intensidad de la precipitación mensual de 2015, como ya se mencionó en el análisis correlacional.

El modelo B o modelo predictivo se usa cuando la intención es predecir la concentración de PM₁:

(5)

En este modelo, el coeficiente que acompaña a la variable PM_3.0-0.95 indica que si el resto de las variables se mantuvieran constantes, a un aumento de una unidad en la variable PM_3.0-0.95 le correspondería, en promedio, un aumento de 0.572 veces la PM₁ (en μg/m³). El mismo análisis se aplica para las otras dos variables.

Validación del modelo predictivo

En el cuadro VIII se presentan los valores medidos y predichos con el modelo B de PM₁ usando datos reservados para la validación, así como el valor de MBE para el modelo predictivo. El valor de MBE igual a -1.5 µg/m³ revela que el modelo subestima la concentración de PM₁, por lo cual la información que arroja el modelo sobre la concentración promedio de PM₁ para una región intertropical en el centro-norte de Venezuela resulta menor que el valor que reportaría el equipo (cuando el modelo prediga una concentración de PM₁, el valor medido será siempre mayor). Esta información es de gran utilidad en los casos en que valores obtenidos por el modelo resulten cercanos a valores considerados no seguros para la salud.

Aún no hay un consenso unificado sobre los valores máximos aceptables para el MBE. La literatura indica que entre más cercano a cero sea el valor de MBE mejor será el modelo predictivo. Se comparó el resultado del MBE obtenido (-1.5 µg/m³) con el de estudios similares para predecir contaminantes atmosféricos: los valores de MBE calculados por ^{Corani y Barazzetta (2004)} para modelos de predicción de PM₁₀ en Milán variaron de -2.46 a 2.34 µg/m³, los cuales son mayores que el obtenido en esta investigación. Los resultados mostrados por ^{Moustris (2010)} para modelos predictivos con NO₂, CO, SO₂, y O₃ resultaron menores a la unidad (-0.771 y 0.615 µg/m³). Los resultados de ^{Dutot et al. (2007)} para O₃ resultaron comparativamente elevados (-13 y 11 µg/m³). Los modelos diseñados en esta investigación representan los primeros resultados obtenidos en predicciones de PM₁ en Venezuela para un área continental intertropical.