Resumen

Los líquenes cortícolas responden a las condiciones ambientales y pueden utilizarse como bioindicadores de calidad ambiental. El presente estudio analizó la relación entre el nivel de urbanización y las características de los árboles hospederos con la abundancia, riqueza y composición de las comunidades liquénicas en tres sitios de la provincia de Heredia, Costa Rica. Se realizó un estudio observacional transversal en 45 árboles (15 por sitio), ubicados en la sede Omar Dengo de la Universidad Nacional (alta urbanización), Barva (urbanización media) y Monte de la Cruz (baja urbanización). En un área fija de observación se registraron la abundancia, riqueza y tipos de talos, así como la especie arbórea, el diámetro a la altura del pecho (DAP), el tipo de corteza y la cobertura del dosel. Se aplicaron análisis estadísticos para comparar los sitios y determinar las relaciones entre las variables. La mayor riqueza se registró en la sede Omar Dengo, con diez morfotipos, seguida de Barva con ocho y Monte de la Cruz con seis. La composición varió entre sitios, evidenciando un recambio de morfotipos y la presencia exclusiva de talos gelatinosos en Monte de la Cruz. Barva presentó la mayor abundancia, con predominio de morfotipos foliosos, mientras que Monte de la Cruz registró valores inferiores. El tipo de corteza y la cobertura del dosel no mostraron efectos significativos sobre la abundancia o riqueza, ni se encontró una asociación significativa con el DAP. El análisis multivariable identificó el sitio de muestreo como el principal factor relacionado con la abundancia liquénica. Se concluye que las diferencias entre sitios se manifestaron principalmente en la composición de las comunidades, más que en su abundancia o riqueza. Se recomienda incorporar mediciones directas de contaminantes atmosféricos, controlar las diferencias altitudinales y realizar identificaciones taxonómicas detalladas para fortalecer el uso de los líquenes como bioindicadores de calidad ambiental.

Palabras clave: líquenes cortícolas, urbanización, corteza, morfotipos.


Abstract

Corticolous lichens respond to environmental conditions and can be used as bioindicators of environmental quality. This study analyzed the relationship between the degree of urbanization and the characteristics of host trees, and the abundance, richness, and composition of lichen communities at three sites in the province of Heredia, Costa Rica. A cross-sectional observational study was conducted on 45 trees (15 per site) located at the Omar Dengo Campus of the National University (high urbanization), Barva (moderate urbanization), and Monte de la Cruz (low urbanization). Lichen abundance, richness, and thallus types were recorded within a fixed observation area, along with tree species, diameter at breast height (DBH), bark type, and canopy cover. Statistical analyses were performed to compare the sites and determine relationships among the variables. The highest morphotype richness was recorded at the Omar Dengo Campus, with ten morphotypes, followed by Barva with eight and Monte de la Cruz with six. Community composition varied among sites, revealing morphotype turnover and the exclusive presence of gelatinous thalli at Monte de la Cruz. Barva exhibited the highest lichen abundance, dominated by foliose morphotypes, whereas Monte de la Cruz showed lower abundance. Bark type and canopy cover showed no significant effects on lichen abundance or richness, and no significant association was found with DBH. Multivariate analysis identified the sampling site as the main factor associated with lichen abundance. It is concluded that differences among sites were primarily reflected in community composition rather than in abundance or richness. Future studies should incorporate direct measurements of atmospheric pollutants, control for altitudinal differences, and conduct detailed taxonomic identification to strengthen the use of lichens as bioindicators of environmental quality.

Keywords: corticolous lichens, urbanization, bark, morphotypes.

INTRODUCCIÓN

Los líquenes son organismos formados por una asociación simbiótica entre un hongo y un fotobionte, como una alga o cianobacteria. Esta asociación les permite desarrollarse sobre diferentes tipos de sustratos, como rocas, corteza de árbol y suelos. Debido a que son organismos que obtienen gran parte del agua y otras sustancias necesarias para su desarrollo, directamente del ambiente, pueden evidenciar cambios en las condiciones ambientales, incluida la presencia de contaminantes incluidos los atmosféricos (Hawksworth et al., 2005; Abas, 2021). Por esta razón, son utilizados como bioindicadores de la calidad de aire, ya que factores como la presencia, cobertura y composición pueden aportar información relevante sobre las condiciones ambientales de un sitio (Varela et al., 2018).

Los líquenes presentan distintas formas de crecimiento, conocidas como tipos de talos, entre los cuales se encuentran las formas crustáceas, foliosas, fruticulosos, escamosas y gelatinosas, entre otros. Estas se diferencian por su forma y textura sobre el sustrato. La composición de las comunidades liquénicas también pueden describirse mediante la riqueza, que corresponde al número de especies presentes en un área determinada, y la abundancia, que se refiere a la cantidad de organismos presentes (Cordero-S. et al., 2021; Bartels et al., 2015). Estos parámetros están asociados a las condiciones del entorno y tipos de sustratos, por lo que observar la composición de estas comunidades permite completar el análisis de su cobertura (Hawksworth et al., 2005).

Sin embargo, la presencia y cobertura de los líquenes no dependen únicamente de la contaminación. Características del árbol hospedero, como la especie, el diámetro a la altura del pecho (DAP) y el tipo de corteza también influyen en su desarrollo y distribución (Sebald et al., 2022). Asimismo, la cobertura del dosel puede modificar el acceso de luz a los fotobiontes, y a la humedad disponible, mientras que factores como la altitud y las condiciones del entorno pueden generar diferencias ambientales entre sitios. Por ello, considerar estos factores permite interpretar de manera más completa las diferencias observadas entre las comunidades de líquenes y evitar atribuirse exclusivamente a la contaminación atmosférica (Abas, 2021; Sebald et al., 2022).

En Costa Rica y otros países se ha documentado la relación entre las comunidades de líquenes y la contaminación atmosférica, considerando variables como la cobertura, la diversidad y los niveles de material particulado (Monge-Nájera et al., 2002; Varela et al., 2018). A partir de estos antecedentes, la presente investigación tiene como objetivo evaluar la variación de la abundancia y riqueza de líquenes presentes en árboles de tres zonas de la provincia de Heredia, Costa Rica, con diferentes niveles de urbanización: la Sede Omar Dengo de la Universidad Nacional (alta), Barva (media) y Monte de la Cruz, en San Rafael de Heredia (baja), considerando su relación con las características de los árboles hospederos.

METODOLOGÍA

La investigación se desarrolló bajo un enfoque cuantitativo, con diseño no experimental, observacional y transversal, debido a que se analizó la presencia de líquenes cortícolas en diferentes condiciones de urbanización sin manipular las variables ambientales de los sitios de estudio. El propósito fue comparar las características de las comunidades de líquenes presentes en tres zonas de Heredia, Costa Rica, seleccionadas porque representan un gradiente de urbanización baja, media y alta.

Sitios de estudio

La selección de estos sitios se realizó considerando sus características de urbanización, con el fin de establecer un gradiente espacial de alta, media y baja urbanización. Esta estrategia se fundamenta en estudios de líquenes desarrollados a lo largo de gradientes de urbanos, en los cuales las características del entorno y la cobertura del suelo permiten diferenciar las condiciones ambientales entre los sitios (Koch et al., 2019).

Los tres sitios de estudio se ubicaron en la provincia de Heredia (Figura 1). La zona de alta urbanización correspondió a la Sede Omar Dengo de la Universidad Nacional (UNA) (Sitio 1; 10°00’04.08”N 84°06’35.65”W, 1 163 msnm), ubicada en un área urbana de Heredia, donde se ha documentado una importante dinámica vial, caracterizada por el volumen de tráfico, rutas de autobuses, parqueos y terminales, (Somarribas Chavarría et al., 2011) además de poseer una alta concentración de infraestructura, edificaciones, y tránsito de personas. La zona de urbanización media se ubicó en Barva de Heredia (Sitio 2; 10°00’58.01”N 84°07’45.05”W, 1 156 msnm), donde se combinan áreas residenciales y comerciales con espacios de vegetación y zonas asociadas a actividades agrícolas, charrales y áreas arboladas (Romero-Vargas et al., 2022), lo que representa una condición intermedia entre un ambiente predominantemente urbano y uno con menor intervención. Finalmente, la zona de baja urbanización correspondió al Monte de la Cruz, en San Rafael de Heredia, (Sitio 3; 10°04’03.55”N 84°04’39.46”W, 1 846 msnm) ubicado en un sector de mayor elevación asociado con zonas protegidas de la Cordillera Volcánica Central, donde predominan los bosques secundarios y primarios, presentando una mayor cobertura vegetal (Municipalidad de San Rafael de Heredia, s.f).

Figura 1. Sitios donde se realizaron muestreos. Sitio 1: la Sede Omar Dengo de la Universidad Nacional (UNA), Sitio 2: Parque de Patinaje en San Roque (Barva) y Sitio 3: Monte de la Cruz en San Rafael.

Figura 1. Sitios donde se realizaron muestreos. Sitio 1: la Sede Omar Dengo de la Universidad Nacional (UNA), Sitio 2: Parque de Patinaje en San Roque (Barva) y Sitio 3: Monte de la Cruz en San Rafael.

Población y tamaño de muestreo

La población de estudio estuvo constituida por los líquenes cortícolas presentes sobre los árboles de las zonas seleccionadas, ya que estos organismos son ampliamente utilizados como bioindicadores de contaminación ambiental (Hawksworth et al.,2005). La unidad de muestreo fue el árbol (forófito), debido a que los líquenes fueron evaluados directamente sobre la corteza de los individuos seleccionados. Se trabajó con una muestra total de 45 árboles, distribuidos equitativamente entre los tres sitios (15 por sitio) y seleccionados al azar, permitiendo así, realizar comparaciones entre los tres niveles de urbanización bajo un tamaño de muestra equivalente. El muestreo fue de tipo no probabilístico intencional, ya que los sitios se seleccionaron previamente según su grado de urbanización y los árboles debían cumplir el criterio de tener una altura mínima de 4 m.

Variables asociadas al forófito

Además, de registrar las características de los líquenes, se consideraron variables relacionadas con el árbol hospedero que pueden influir en la presencia y cobertura de los líquenes y que permiten considerar posibles diferencias entre los árboles seleccionados (Zárate-Arias et al., 2019): la especie del árbol y el diámetro a la altura del pecho (DAP), medido con cinta diamétrica (cm), el porcentaje de cobertura del dosel, siguiendo la metodología de Suqui et al (2021) y finalmente, se determinó el tipo de corteza de cada forofito, aplicando la clasificación propuesta por Cáceres et al (2007), compuesta por las categorías de: lisa, parcialmente fisurada, moderadamente fisurada y muy fisurada.

Variables registradas en los líquenes

En cada árbol seleccionado se estableció un área de observación con dimensiones iguales de 40x50 cm a una altura de 1,5 m sobre el tronco, según la metodología planteada por Zárate-Arias et al (2019). La superficie de observación se mantuvo constante en todos los individuos y se orientó hacia el norte con respecto al tronco, (Méndez-Estrada & Campos, 2015). El área se dividió en 20 cuadrantes de 10x10 cm, de los cuales se seleccionaron 5 al azar para el conteo. En cada cuadrante se registró el número de talos de cada morfotipo y se asignó cada morfotipo a una forma de crecimiento: crustácea, foliosa, fruticulosa o gelatinosa. Los morfotipos se diferenciaron por sus características morfológicas, sin identificación a nivel de especie. Además, se registró el porcentaje de cobertura que estos representan sobre la superficie de la corteza de cada forofito. A partir de estos registros se definieron las variables de respuesta: la abundancia, como el número total de talos por árbol (suma de los 5 cuadrantes); la riqueza, como el número de morfotipos distintos por árbol; y el tipo de talo (crustáceo, folioso o gelatinoso), que permitió contabilizar el total de los talos de cada tipo por sitio. Los árboles sin líquenes se registraron con el valor de cero.

Análisis de datos

Los datos se organizaron en hojas de cálculo de Excel y se analizaron en el software R, versión 4.6.1 (R Core Team, 2026). La manipulación de bases de datos y gráficos se realizó mediante los paquetes readxl (Wickham & Bryan, 2023), dplyr (Wickham et al., 2026) y ggplot2 (Wickham, 2016). Las pruebas de homogeneidad de varianzas, autocorrelación y valores atípicos se efectuaron con car (Fox & Weisberg, 2019), lmtest (Zeileis et al., 2002) y outlier (Komsta, 2022), mientras que la partición de la diversidad beta se calculó con betapart (Baselga et al., 2025). En todas las pruebas se utilizó un nivel de significancia de α=0,05 y la hipótesis nula fue la ausencia de diferencias, de asociación o de relación. Inicialmente se calcularon el total de talos, riqueza de los morfotipos y la frecuencia de registros de cada tipo de talo.

Características de los forófitos

Para comparar el DAP y el porcentaje de cobertura del dosel entre los tres sitios de estudio de urbanización, se utilizó un tamaño de muestra total de n=45 (15 forófitos por sitio). En cada variable cuantitativa se evaluaron las comparaciones por pares entre los tres sitios (Sitio 1 vs. Sitio 2, Sitio 1 vs. Sitio 3 y Sitio 2 vs. Sitio 3) Inicialmente se planteó un análisis de varianza (ANOVA) de una vía, comprobando los supuestos de normalidad de residuos con la prueba de Shapiro-Wilk y de homogeneidad de varianzas con la prueba de Cochran. Al no cumplirse dichos supuestos, se recurrió a la prueba no paramétrica de Kruskal-Wallis, bajo la hipótesis nula de igualdad de medianas entre los sitios. Tras obtener un resultado significativo en la prueba global, se evaluaron formalmente los contrastes por pares antes descritos mediante la prueba de Wilcoxon. Se aplicó la corrección de Bonferroni a los valores de p justificada por tratarse de un diseño balanceado con menos de seis comparaciones planificadas. Además, para determinar si la distribución de los tipos de corteza (lisa,parcialmente fisurada y muy fisurada) difería entre los sitios, se construyó una tabla de contingencia con el número de forófitos por categoría de corteza y sitios. Sobre esta se aplicó una prueba de Chi-cuadrado, calculando la significancia mediante la prueba exacta de Fisher.

Riqueza y diversidad beta

La riqueza se definió como el número total de morfotipos registrados en cada sitio. Para evaluar las diferencias en la composición de morfotipos entre sitios se construyó una matriz de presencia/ausencia (sitios x morfotipos) y se calculó la diversidad beta mediante la partición del índice de disimilitud de Sørensen ,con la librería betapart (Baselga et al., 2025). Este índice se descompuso en la disimilitud total (βSOR), el componente de recambio (βSIM), que representa el reemplazo de unos morfotipos por otros entre sitios, y el componente de anidamiento (β SNE), que representa la pérdida de morfotipos de un sitio respecto a otro. Los índices se calcularon para los tres sitios en conjunto y para cada par de sitios.

Tipos de talo

Para determinar si la frecuencia de los distintos tipos de talos (crustáceo, folioso y gelatinoso) depende del nivel de urbanización, se construyó una tabla de contingencia con el número de registros de cada tipo de talo por sitio (n=253) y se aplicó una prueba de Chi-cuadrado de independencia. La hipótesis nula fue que la proporción de los tipos de talo es la misma en los tres sitios. Dado que algunas frecuencias esperadas fueron menores a 5, el valor de p se calculó mediante una prueba de Fisher. Cuando la asociación fue significativa, se identificaron las categorías responsables a partir de los residuos estandarizados de la prueba, considerando sobrerrepresentadas y subrepresentadas aquellas con valores absolutos mayores a 1.96. Con el mismo procedimiento se evaluó la asociación entre el tipo de talo y el tipo de corteza del forófito.

Abundancia de líquenes por sitios

Para comparar la abundancia de líquenes, definida como la cantidad de talos por registro (n = 302), entre los tres niveles de urbanización, se planteó un ANOVA de una vía y se verificaron los supuestos de normalidad ( según Shapiro-Wilk) y de homogeneidad de varianzas ( según Cochran). Debido a que no se cumplieron, se aplicó la prueba de Kruskal-Wallis y, al resultar significativa, se realizó la prueba de comparaciones Wilcoxon por pares con corrección de Bonferroni. Para corroborar el resultado considerando las características del forófito, se ajustó una regresión lineal múltiple (lm()) con un valor por forófito de (n=45), en la que la variable de respuesta fue la raíz cuadrada de la cantidad total de talos, transformación utilizada para datos de conteo, y las variables explicativas fueron el sitio, el DAP y el tipo de corteza. El sitio de baja urbanización (Monte de la Cruz), se utilizó como referencia. La contribución de cada variable se evaluó mediante pruebas F. Posteriormente, se compararon cinco modelos alternativos (completo; solo sitio; sitio y DAP; sitio y corteza; y solo corteza) mediante el criterio de Akaike (AIC), seleccionando como mejor modelo el de menor valor. En el modelo seleccionado se cambió la categoría de referencia para comparar también Barva con la UNA, y se verificaron supuestos de normalidad de los residuos mediante Shapiro y de homocedasticidad con la Prueba de Breusch-Pagan de la librería “car” .

Relación de la abundancia y la riqueza con el DAP y la cobertura del dosel

Para evaluar si el DAP permite predecir la abundancia de líquenes, se ajustó una regresión lineal simple con un valor por forófito donde (n=45), con el fin de evitar la pseudorreplicación. Se utilizó la raíz cuadrada de la cantidad total de talos como variable dependiente y el DAP como variable independiente. Se evaluó la significancia del intercepto (H0 : a=0) y de la pendiente (H0 : b=0), así como el coeficiente de determinación (R2). Los supuestos del modelo se verificaron mediante la prueba de Durbin-Watson para la independencia de los residuos del paquete lmtest (Zeileis et al., 2002), la prueba de Shapiro-Wilk para su normalidad y la prueba de varianza no constante para la homocedasticidad del paquete car. Finalmente, se realizaron predicciones sustituyendo valores del DAP del 1,2 y 3 m en la ecuación obtenida y elevando el resultado al cuadrado para expresarlo en número de talos. Además, para valorar la relación de la abundancia y la riqueza de morfotipos por forófito con la cobertura del dosel y con el DAP se emplearon correlaciones. Se verificó el supuesto de normalidad de cada variable con la prueba de Shapiro-Wilk y, como este no se cumplió, se utilizó el coeficiente de correlación por rangos de Spearman. Se analizaron cuatro pares de variables (abundancia-cobertura del dosel, riqueza-cobertura del dosel, abundancia-DAP y riqueza-DAP), y se interpretaron la dirección y la magnitud del coeficiente, así como su valor p.

RESULTADOS

Panorama general

Se contabilizaron 989 talos de líquenes en 45 árboles, distribuidos en 12 morfotipos. El Parque de Patinaje en San Roque (Barva) presentó el mayor número de talos contabilizando 579, seguido de la Sede Omar Dengo de la Universidad Nacional (UNA) con 280 y el Monte de la Cruz con 130, mientras que la riqueza total de morfotipos fue de 10 en la UNA, 8 en Barva y 6 en el Monte de la Cruz. Del total de 253 registros de morfotipos por cuadrante, 176 (69.6 %) correspondieron a talos foliosos, 68 (26.9%) a crustáceos y 9 (3.6%) a gelatinosos.

# Script utilizado para obtener el Panorama General:
sum(datos$cantidad)                             # Total de talos
## [1] 989
tapply(datos$cantidad, datos$urbanizacion, sum) # Talos por sitio
##  Alta Media  Baja 
##   280   579   130
length(unique(liq$morfotipo))                   # Riqueza total morfotipos
## [1] 12
tapply(liq$morfotipo, liq$urbanizacion, function(x) length(unique(x))) # Riqueza por sitio
##  Alta Media  Baja 
##    10     8     6
nrow(liq)                                       # Registros de morfotipos
## [1] 253
tt <- table(liq$forma_crecimiento); tt          # Frecuencia formas de crecimiento
## 
## cru fol gel 
##  68 176   9
round(100 * prop.table(tt), 1)                  # Porcentajes
## 
##  cru  fol  gel 
## 26.9 69.6  3.6

Características de los forófitos

En cuanto a las características de los forófitos, el DAP no cumplió con los supuestos del ANOVA paramétrico con un nivel de significancia de α=0.05, ya que los residuos no presentaron una distribución normal según la prueba de Shapiro Wilk (W=0.882, p < 0.001) y las varianzas no fueron homogéneas según la prueba de Cochran (C = 0.502, p= 0.165). Por ello, se le aplicó una prueba no paramétrica de Kruskal-Wallis, la cual mostró diferencias significativas en el DAP entre sitios de estudio (H = 9.87, df = 2, p = 0.007). Mediante la prueba de Wilcoxon por pares con la corrección de Bonferroni, se encontró que los forófitos de la UNA, sitio de alta urbanización (mediana = 181 cm), presentaron un DAP mayor que los de Barva (mediana=97 cm; padj = 0.098) y que los del Monte de la Cruz (mediana = 96 cm; padj = 0.005) (Figura 2). Barva y Monte de la Cruz no difirieron entre sí (padj = 1.000).

library(outliers)
# Script utilizado para evaluar el DAP de los forófitos:
m_dap <- aov(dap_cm ~ urbanizacion, data = arb)
shapiro.test(residuals(m_dap))                   
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(m_dap)
## W = 0.88387, p-value = 0.0003096
cochran.test(dap_cm ~ urbanizacion, data = arb)  
## 
##  Cochran test for outlying variance
## 
## data:  dap_cm ~ urbanizacion
## C = 0.5025, df = 15, k = 3, p-value = 0.1654
## alternative hypothesis: Group Alta has outlying variance
## sample estimates:
##      Alta     Media      Baja 
## 10490.829  9959.607   426.981
kruskal.test(dap_cm ~ urbanizacion, data = arb)  
## 
##  Kruskal-Wallis rank sum test
## 
## data:  dap_cm by urbanizacion
## Kruskal-Wallis chi-squared = 9.8777, df = 2, p-value = 0.007163
pairwise.wilcox.test(arb$dap_cm, arb$urbanizacion,
                     p.adjust.method = "bonferroni", exact = FALSE)
## 
##  Pairwise comparisons using Wilcoxon rank sum test with continuity correction 
## 
## data:  arb$dap_cm and arb$urbanizacion 
## 
##       Alta   Media 
## Media 0.0977 -     
## Baja  0.0049 1.0000
## 
## P value adjustment method: bonferroni
tapply(arb$dap_cm, arb$urbanizacion, median)     
##  Alta Media  Baja 
##   181    97    96
# Script utilizado para generar la Figura 2:
boxplot(dap_cm ~ urbanizacion, data = arb,
        names = c("UNA", "Barva", "Monte de la Cruz"),
        col = c("lightgreen", "lightgreen", "lightgreen"),
        xlab = "Sitio", ylab = "DAP (cm)")
Figura 2. DAP de los forófitos por sitio de estudio.

Figura 2. DAP de los forófitos por sitio de estudio.

Por otra parte, no se encontraron diferencias significativas en el porcentaje de cobertura del dosel entre los sitios (H = 3.104, df = 2, p = 0.212) ni en la distribución de los tipos de corteza (X2 = 6.14, p = 0.440). Esto indica que los forófitos de los tres sitios fueron comparables en cuanto a la luminosidad y al tipo de sustrato disponible para los líquenes.

# Script para Cobertura de dosel y Tipo de corteza:
kruskal.test(dosel ~ urbanizacion, data = arb)   
## 
##  Kruskal-Wallis rank sum test
## 
## data:  dosel by urbanizacion
## Kruskal-Wallis chi-squared = 3.1042, df = 2, p-value = 0.2118
tab_cort <- table(arb$urbanizacion, arb$corteza); tab_cort
##        
##         1 2 3 4
##   Alta  1 5 6 3
##   Media 3 5 1 6
##   Baja  2 4 3 6
chisq.test(tab_cort, simulate.p.value = TRUE, B = 10000)   
## 
##  Pearson's Chi-squared test with simulated p-value (based on 10000
##  replicates)
## 
## data:  tab_cort
## X-squared = 6.1429, df = NA, p-value = 0.4394
fisher.test(tab_cort, simulate.p.value = TRUE, B = 10000)
## 
##  Fisher's Exact Test for Count Data with simulated p-value (based on
##  10000 replicates)
## 
## data:  tab_cort
## p-value = 0.4288
## alternative hypothesis: two.sided

Riqueza y diversidad beta

Respecto a la composición de las comunidades de líquenes, la riqueza de los morfotipos cambió a lo largo del gradiente de urbanización. La riqueza más alta se registró en la UNA con 10 morfotipos, seguida por Barva con 8 y por último el Monte de la Cruz con 6. La partición del índice de disimilitud de Sørensen (βSOR) mostró una diversidad beta total de βSOR = 0.429. Al descomponer este valor, se encontró que el componente principal de la disimilitud fue el recambio de morfotipos (βSIM = 0.294), el cual representa el reemplazo de unos morfotipos por otros entre sitios. El componente de anidamiento fue menor (βSNE = 0.134), lo que indica que las diferencias en riqueza entre sitios contribuyen en menor medida a la diversidad beta total. Las comparaciones por pares de sitios mostraron que la mayor similitud se presentó entre la UNA y Barva (βSOR = 0.111). En este par la disimilitud se debió por completo al anidamiento (βSNE = 0.111), es decir, los morfotipos de Barva constituyeron un subconjunto de los registros en la UNA. La mayor disimilitud se encontró entre Barva y Monte de la Cruz (βSOR = 0.571), atribuida principalmente a recambio (βSIM = 0.500), seguida por la UNA y Monte de la Cruz (βSOR = 0.500); (βSIM = 0.333). Lo anterior sugiere que Monte de la Cruz posee una composición de morfotipos notablemente distinta a la de los dos sitios urbanizados. Además, fue el único con morfotipos exclusivos de tipo gelatinoso (M11) y crustáceo (M12).

library(betapart)
# Script para Matriz de presencia/ausencia e Índices de Sørensen:
pa_sitio <- as.data.frame.matrix(table(liq$urbanizacion, liq$morfotipo))
pa_sitio[pa_sitio > 0] <- 1
rowSums(pa_sitio)                               
##  Alta Media  Baja 
##    10     8     6
beta.multi(pa_sitio, index.family = "sorensen")  
## $beta.SIM
## [1] 0.2941176
## 
## $beta.SNE
## [1] 0.1344538
## 
## $beta.SOR
## [1] 0.4285714
bp <- beta.pair(pa_sitio, index.family = "sorensen")
bp$beta.sor; bp$beta.sim; bp$beta.sne            
##            Alta     Media
## Media 0.1111111          
## Baja  0.5000000 0.5714286
##            Alta     Media
## Media 0.0000000          
## Baja  0.3333333 0.5000000
##             Alta      Media
## Media 0.11111111           
## Baja  0.16666667 0.07142857

Tipos de talo

En cuanto a los tipos de talo, se encontró una asociación estadísticamente significativa entre el tipo de talo y el sitio de estudio (X2 = 35.33, df = 4, p < 0.001; Tabla 1). Este resultado se corroboró mediante la prueba de Fisher (p < 0.001), debido a la baja frecuencia de los talos gelatinosos, y confirma que la proporción de los diferentes tipos de talo no es homogénea entre los sitios de estudio. El análisis de residuos estandarizados respaldó estos resultados: en la UNA estuvieron sobrerrepresentados los talos de tipo crustáceos (residuo = 2.27), en Barva los de tipo folioso (residuo = 3.35) y en Monte de la Cruz los de tipo gelatinoso (residuo = 5.25).

# Script para generar la Tabla 1 y su prueba estadística:
tab_talo <- table(liq$urbanizacion, liq$forma_crecimiento)[, c("cru","fol","gel")]
colnames(tab_talo) <- c("Crustaceo","Folioso","Gelatinoso")

ct <- chisq.test(tab_talo)                      
## Warning in chisq.test(tab_talo): Chi-squared approximation may be incorrect
ct
## 
##  Pearson's Chi-squared test
## 
## data:  tab_talo
## X-squared = 35.332, df = 4, p-value = 3.97e-07
fisher.test(tab_talo, simulate.p.value = TRUE, B = 10000)   
## 
##  Fisher's Exact Test for Count Data with simulated p-value (based on
##  10000 replicates)
## 
## data:  tab_talo
## p-value = 9.999e-05
## alternative hypothesis: two.sided
round(ct$stdres, 2)                             
##        
##         Crustaceo Folioso Gelatinoso
##   Alta       2.27   -1.30      -2.21
##   Media     -2.43    3.35      -2.51
##   Baja       0.26   -2.36       5.25
knitr::kable(tab_talo, caption = "Tabla 1. Frecuencia de los tipos de talo de líquenes registrados en los forófitos de cada sitio de estudio.")
Tabla 1. Frecuencia de los tipos de talo de líquenes registrados en los forófitos de cada sitio de estudio.
Crustaceo Folioso Gelatinoso
Alta 31 56 0
Media 19 83 0
Baja 18 37 9

Como se muestra en la Figura 3, el talo gelatinoso presentó una distribución totalmente sesgada hacia el Monte de la Cruz, ya que fue el único sitio registrado (100%). Los talos crustáceos se registraron principalmente en la UNA (45.6%) y en menor medida en Barva (27.9%) y Monte de la Cruz (26.5%). Los foliosos se concentraron en Barva (47.2%), seguida por la UNA (31.8%) y Monte de la Cruz (21.0%).

# Script para generar la Figura 3:
library(ggplot2) 
library(dplyr)
g3 <- as.data.frame(tab_talo)
names(g3) <- c("Urbanizacion", "Tipo_talo", "Conteo")
g3 <- g3 %>%
  group_by(Tipo_talo) %>%
  mutate(Porcentaje = Conteo / sum(Conteo) * 100) %>%
  ungroup()

g3$Urbanizacion <- factor(g3$Urbanizacion, levels = c("Baja", "Media", "Alta"))
g3$Tipo_talo    <- factor(g3$Tipo_talo, levels = c("Crustaceo", "Gelatinoso", "Folioso"))
col_sitio       <- c("Baja" = "lightgreen", "Media" = "skyblue", "Alta" = "pink")

ggplot(g3, aes(x = Tipo_talo, y = Porcentaje, fill = Urbanizacion)) +
  geom_bar(stat = "identity", position = "fill", color = "black", width = 0.7) +
  geom_text(aes(label = ifelse(Porcentaje > 0, paste0(round(Porcentaje, 1), "%"), "")),
            position = position_fill(vjust = 0.5), size = 3.5) +
  scale_y_continuous(labels = function(x) x * 100, name = "Frecuencia relativa por tipo de talo") +
  scale_x_discrete(name = "") +
  scale_fill_manual(name = "", values = col_sitio) +
  theme_classic() +
  theme(axis.title.y = element_text(size = 12), axis.text = element_text(size = 11), legend.position = "right")
Figura 3. Distribución de frecuencias relativas de cada tipo de talo entre los sitios de estudio.

Figura 3. Distribución de frecuencias relativas de cada tipo de talo entre los sitios de estudio.

Por otra parte, no se encontró asociación significativa entre el tipo de talo y el tipo de corteza (X2 = 10.49, p = 0.106).

# Script para Tipo de talo x Tipo de corteza:
tab_tc <- table(liq$tipo_corteza, liq$forma_crecimiento)[, c("fol","gel","cru")]
tab_tc
##    
##     fol gel cru
##   1  28   1   4
##   2  56   2  32
##   3  38   1  15
##   4  54   5  17
chisq.test(tab_tc, simulate.p.value = TRUE, B = 10000)      
## 
##  Pearson's Chi-squared test with simulated p-value (based on 10000
##  replicates)
## 
## data:  tab_tc
## X-squared = 10.49, df = NA, p-value = 0.1012

Abundancia de líquenes por sitios

Se encontraron diferencias significativas en la cantidad de líquenes (abundancia) entre los sitios de estudio (H = 36.293, df = 2, p < 0.001). Los datos no cumplieron los supuestos del ANOVA paramétrico, ya que los residuos no presentaron distribución normal (p < 0.001) y las varianzas no fueron homogéneas (C = 0.836, p < 0.001). La prueba de Wilcoxon por pares con corrección de Bonferroni, indicó diferencias significativas entre los sitios (Figura 4). Barva presentó la mayor cantidad promedio de talos por registro (5.3), seguida de la UNA (2.6; padj = 0.006) y del Monte de la Cruz (1.5). El Monte de la Cruz presentó diferencias tanto de la UNA (padj = 0.002) como de Barva (padj < 0.001)

# Script para evaluar diferencias de abundancia entre sitios:
anova_ab <- aov(cantidad ~ urbanizacion, data = datos)
shapiro.test(residuals(anova_ab))               
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(anova_ab)
## W = 0.80315, p-value < 2.2e-16
cochran.test(cantidad ~ urbanizacion, data = datos)         
## 
##  Cochran test for outlying variance
## 
## data:  cantidad ~ urbanizacion
## C = 0.83629, df = 100.67, k = 3.00, p-value < 2.2e-16
## alternative hypothesis: Group Media has outlying variance
## sample estimates:
##      Alta     Media      Baja 
##  4.860737 34.416097  1.876334
kruskal.test(cantidad ~ urbanizacion, data = datos)        
## 
##  Kruskal-Wallis rank sum test
## 
## data:  cantidad by urbanizacion
## Kruskal-Wallis chi-squared = 36.293, df = 2, p-value = 1.315e-08
tapply(datos$cantidad, datos$urbanizacion, mean)            
##     Alta    Media     Baja 
## 2.641509 5.263636 1.511628
pairwise.wilcox.test(datos$cantidad, datos$urbanizacion,
                     p.adjust.method = "bonferroni", exact = FALSE)
## 
##  Pairwise comparisons using Wilcoxon rank sum test with continuity correction 
## 
## data:  datos$cantidad and datos$urbanizacion 
## 
##       Alta   Media  
## Media 0.0056 -      
## Baja  0.0022 1.6e-08
## 
## P value adjustment method: bonferroni
# Script para generar la Figura 4:
boxplot(cantidad ~ urbanizacion, data = datos,
        names = c("UNA", "Barva", "Monte de la Cruz"),
        col = c("pink", "skyblue", "lightgreen"),
        xlab = "Urbanización", ylab = "Cantidad de líquenes")
Figura 4. Cantidad de talos de líquenes totales en los tres sitios de estudio: alta (UNA), media (Barva) y baja (Monte de la Cruz)

Figura 4. Cantidad de talos de líquenes totales en los tres sitios de estudio: alta (UNA), media (Barva) y baja (Monte de la Cruz)

Estos resultados se corroboraron mediante una regresión lineal múltiple, en la que se utilizó la raíz cuadrada de la cantidad de talos como variable respuesta. El modelo completo, que incluyó el sitio, el DAP y el tipo de corteza, fue significativo (F(6,38) = 3.12, p = 0.014,R2 ajustado = 0.22). Sin embargo, únicamente el sitio aportó de forma significativa a la explicación de la abundancia (p = 0.003); ni el DAP (p = 0.925) ni el tipo de corteza (p = 0.286) lo hicieron.

A pesar de la falta de significancia estadística del tipo de corteza, la Figura 5, permite observar que los forófitos de corteza parcialmente fisurada presentaron la mediana de abundancia más alta (22 talos), seguidos por los de corteza moderadamente fisurada (17 talos). Las cortezas lisas (9.5 talos) y muy fisuradas (9 talos) mostraron las medianas más bajas, aunque esta última tuvo la mayor dispersión. De igual manera, en la Figura 6 no se observa una tendencia clara entre el DAP y la abundancia de talos. Los forófitos de mayor DAP se concentraron en la UNA y Barva, mientras que en el Monte de la Cruz los DAP fueron más homogéneos (57-125 cm) y la abundancia se mantuvo baja, independientemente del tamaño del árbol.

# Script para Regresión Lineal Múltiple y selección de modelos (AIC):
library(lmtest) # <-- Carga la librería para dwtest()
## Cargando paquete requerido: zoo
## 
## Adjuntando el paquete: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(car)
## Cargando paquete requerido: carData
## 
## Adjuntando el paquete: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
reg <- arb
reg$sitio <- relevel(reg$urbanizacion, ref = "Baja")        
m1 <- lm(raiz ~ sitio + dap_m + corteza, data = reg)        # completo
m2 <- lm(raiz ~ sitio, data = reg)                          # solo sitio
m3 <- lm(raiz ~ sitio + dap_m, data = reg)                  # sitio + DAP
m4 <- lm(raiz ~ sitio + corteza, data = reg)                # sitio + corteza
m5 <- lm(raiz ~ corteza, data = reg)                        # solo corteza

summary(m1)            
## 
## Call:
## lm(formula = raiz ~ sitio + dap_m + corteza, data = reg)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6938 -0.6766  0.0953  0.9996  4.6406 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.33182    0.96392   1.382  0.17514    
## sitioAlta    0.79366    0.85641   0.927  0.35991    
## sitioMedia   2.79052    0.77116   3.619  0.00086 ***
## dap_m        0.03609    0.41252   0.087  0.93075    
## corteza2     2.13446    1.08024   1.976  0.05546 .  
## corteza3     1.50081    1.13359   1.324  0.19343    
## corteza4     1.54616    0.99699   1.551  0.12923    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.04 on 38 degrees of freedom
## Multiple R-squared:  0.3297, Adjusted R-squared:  0.2239 
## F-statistic: 3.115 on 6 and 38 DF,  p-value: 0.01396
drop1(m1, test = "F")  
## Single term deletions
## 
## Model:
## raiz ~ sitio + dap_m + corteza
##         Df Sum of Sq    RSS    AIC F value  Pr(>F)   
## <none>               158.07 70.537                   
## sitio    2    58.421 216.49 80.690  7.0221 0.00254 **
## dap_m    1     0.032 158.10 68.546  0.0077 0.93075   
## corteza  3    16.544 174.62 69.017  1.3257 0.28026   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AIC(m1, m2, m3, m4, m5)                         
##    df      AIC
## m1  8 200.2418
## m2  4 197.4779
## m3  5 198.7211
## m4  7 198.2509
## m5  5 208.6440
summary(m2)            
## 
## Call:
## lm(formula = raiz ~ sitio, data = reg)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6007 -0.8542  0.1458  1.1075  4.7434 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   2.8542     0.5309   5.376 3.12e-06 ***
## sitioAlta     0.9643     0.7508   1.284 0.206090    
## sitioMedia    2.7464     0.7508   3.658 0.000703 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.056 on 42 degrees of freedom
## Multiple R-squared:  0.247,  Adjusted R-squared:  0.2111 
## F-statistic: 6.888 on 2 and 42 DF,  p-value: 0.002588
shapiro.test(residuals(m2))                      
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(m2)
## W = 0.96851, p-value = 0.2554
ncvTest(m2)                                      
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 10.40885, Df = 1, p = 0.0012541
reg$sitio_alta <- relevel(reg$sitio, ref = "Alta")
summary(lm(raiz ~ sitio_alta, data = reg))       
## 
## Call:
## lm(formula = raiz ~ sitio_alta, data = reg)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6007 -0.8542  0.1458  1.1075  4.7434 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       3.8185     0.5309   7.192 7.73e-09 ***
## sitio_altaBaja   -0.9643     0.7508  -1.284   0.2061    
## sitio_altaMedia   1.7822     0.7508   2.374   0.0223 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.056 on 42 degrees of freedom
## Multiple R-squared:  0.247,  Adjusted R-squared:  0.2111 
## F-statistic: 6.888 on 2 and 42 DF,  p-value: 0.002588
# Script para generar la Figura 5:
boxplot(raiz ~ corteza, data = reg,
        names = c("Lisa", "Parc. fisurada", "Mod. fisurada", "Muy fisurada"),
        col = c("khaki", "tan", "tan4", "saddlebrown"),
        xlab = "Tipo de corteza", ylab = "Raíz de la cantidad de talos")    
Figura 5. Cantidad de talos de líquenes por forófito transformado por su raíz cuadrada según el tipo de corteza.

Figura 5. Cantidad de talos de líquenes por forófito transformado por su raíz cuadrada según el tipo de corteza.

# Script para generar la Figura 6:
plot(reg$dap_cm, reg$raiz, pch = 21, cex = 1.5,
     bg = c("pink", "skyblue", "lightgreen")[as.integer(reg$urbanizacion)], col = "black",
     xlab = "DAP (cm)", ylab = "Raíz de la cantidad de talos")
legend("topright", legend = c("Monte de la Cruz", "Barva", "UNA"),
       pt.bg = c("lightgreen", "skyblue", "pink"), pch = 21, pt.cex = 1.5, bty = "n")
Figura 6. Relación entre el DAP de los forófitos y la cantidad de talos de líquenes, transformada mediante raíz cuadrada, en los tres sitios de estudio.

Figura 6. Relación entre el DAP de los forófitos y la cantidad de talos de líquenes, transformada mediante raíz cuadrada, en los tres sitios de estudio.

De acuerdo con el criterio de información de Akaike, el mejor modelo fue el que incluyó sólo el sitio (AIC = 197.5), y fue estadísticamente significativo (F(2,42) = 6.94, p = 0.003, R2 ajustado = 0.21). Según este modelo, Barva presentó una mayor abundancia de talos que Monte de la Cruz (= 2.76, p < 0.001) y que la UNA ( = 1.78, p = 0.022), mientras que la UNA no difirió de Monte de la Cruz ( = 0.98, p = 0.201). Los residuos del modelo cumplieron el supuesto de normalidad (p = 0.259), pero no el de homocedasticidad (p = 0.001), debido a la alta variabilidad en la abundancia de los talos entre los forófitos de Barva.

Relación de la abundancia y la riqueza con el DAP y la cobertura del dosel

La regresión lineal simple también se utilizó la raíz cuadrada de la cantidad total de talos como variable dependiente y el DAP de los forófitos. El modelo generó la ecuación (cantidad de talos= 3.57 + 0.37 * DAP) como se muestra en la Figura 7. El intercepto fue significativo (a = 3.57, p < 0.001), pero la pendiente no difirió de cero (b = 0.37, p = 0.340). El modelo no fue significativo (F(1,43) = 0.93, p = 0.340) y el DAP explicó menos del 1% de la variación en la cantidad de talos (R2 = 0.021). En consecuencia, las predicciones obtenidas a partir de la ecuación fueron prácticamente constantes: 15.9 talos para un DAP de 1m, 18.5 para un DAP de 2m y 21.8 para un DAP de 3m. Esto indica que el tamaño del forófito no permite predecir la abundancia de líquenes. El modelo cumplió los supuestos de la regresión lineal. Los residuos fueron independientes (DW = 1.75, p = 0.19). Además, los residuos presentaron distribución normal (W = 0.966, p = 0.205) y varianza constante (X2 = 1.22, p = 0.270).

# Script para Regresión Lineal Simple y comprobación de supuestos:
reg_dap <- lm(raiz ~ dap_m, data = arb)
summary(reg_dap)      
## 
## Call:
## lm(formula = raiz ~ dap_m, data = arb)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.0216 -1.3347 -0.2701  1.4272  6.4211 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   3.5687     0.6423   5.556 1.61e-06 ***
## dap_m         0.3652     0.3786   0.965     0.34    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.317 on 43 degrees of freedom
## Multiple R-squared:  0.02118,    Adjusted R-squared:  -0.001583 
## F-statistic: 0.9305 on 1 and 43 DF,  p-value: 0.3401
dwtest(reg_dap)        
## 
##  Durbin-Watson test
## 
## data:  reg_dap
## DW = 1.7487, p-value = 0.1863
## alternative hypothesis: true autocorrelation is greater than 0
shapiro.test(residuals(reg_dap))                
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(reg_dap)
## W = 0.96596, p-value = 0.2055
ncvTest(reg_dap)                                 
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 1.218798, Df = 1, p = 0.2696
# Predicciones a 1m, 2m y 3m:
b <- coef(reg_dap)
(b[1] + b[2] * c(1, 2, 3))^2                     
## [1] 15.47569 18.48214 21.75529
# Script para generar la Figura 7:
plot(arb$dap_m, arb$raiz, pch = 21, cex = 1.5,
     bg = c("pink", "skyblue", "lightgreen")[as.integer(arb$urbanizacion)], col = "black",
     xlab = "DAP (m)", ylab = "Raíz cuadrada de la cantidad de talos")
abline(reg_dap, col = "red", lwd = 2)
legend("topright", legend = c("UNA", "Barva", "Monte de la Cruz"),
       pt.bg = c("pink", "skyblue", "lightgreen"), pch = 21, pt.cex = 1.5,
       title = "Urbanización", bty = "n")
Figura 7. Relación entre el DAP de los forófitos y la raíz cuadrada de la cantidad de talos de líquenes por forófito. La línea roja representa la recta de regresión lineal.

Figura 7. Relación entre el DAP de los forófitos y la raíz cuadrada de la cantidad de talos de líquenes por forófito. La línea roja representa la recta de regresión lineal.

Finalmente, no se encontraron correlaciones significativas entre la abundancia y riqueza de líquenes respecto al porcentaje de cobertura del dosel (p > 0.05). En el caso de la abundancia, el coeficiente de Spearman fue de ρ = -0.067, lo que indica que la cobertura del dosel no se relaciona de manera directa con la abundancia de líquenes en los forófitos muestreados. En cuanto a la riqueza, el coeficiente fue de ρ = -0.051, una correlación negativa muy débil. De igual manera, no se hallaron correlaciones significativas entre el DAP de los forófitos y la abundancia (ρ = 0.175, p = 0.249) o la riqueza de líquenes (ρ = -0.031, p = 0.983), lo que concuerda con los resultados de regresión, donde el tamaño del forófito no se relacionó con la cantidad de talos ni con el número de morfotipos.

# Script para Correlaciones de Spearman:
shapiro.test(arb$talos); shapiro.test(arb$riqueza)          
## 
##  Shapiro-Wilk normality test
## 
## data:  arb$talos
## W = 0.80943, p-value = 3.859e-06
## 
##  Shapiro-Wilk normality test
## 
## data:  arb$riqueza
## W = 0.91117, p-value = 0.00215
shapiro.test(arb$dosel);  shapiro.test(arb$dap_cm)          
## 
##  Shapiro-Wilk normality test
## 
## data:  arb$dosel
## W = 0.77202, p-value = 6.088e-07
## 
##  Shapiro-Wilk normality test
## 
## data:  arb$dap_cm
## W = 0.8386, p-value = 1.888e-05
cor.test(arb$talos,   arb$dosel,  method = "spearman", exact = FALSE)  
## 
##  Spearman's rank correlation rho
## 
## data:  arb$talos and arb$dosel
## S = 16134, p-value = 0.6817
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##         rho 
## -0.06285055
cor.test(arb$riqueza, arb$dosel,  method = "spearman", exact = FALSE) 
## 
##  Spearman's rank correlation rho
## 
## data:  arb$riqueza and arb$dosel
## S = 15955, p-value = 0.7391
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##         rho 
## -0.05105382
cor.test(arb$talos,   arb$dap_cm, method = "spearman", exact = FALSE) 
## 
##  Spearman's rank correlation rho
## 
## data:  arb$talos and arb$dap_cm
## S = 12520, p-value = 0.2496
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##       rho 
## 0.1752292
cor.test(arb$riqueza, arb$dap_cm, method = "spearman", exact = FALSE)  
## 
##  Spearman's rank correlation rho
## 
## data:  arb$riqueza and arb$dap_cm
## S = 15227, p-value = 0.9838
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
##          rho 
## -0.003121972

DISCUSIÓN

La mayor cantidad de líquenes en los sitios con más urbanización, junto con el predominio de los talos foliosos, no necesariamente refleja una mejor calidad ambiental en estos sitios. En los ambientes urbanos tienden a predominar especies tolerantes a la sequía, a la alta exposición de luz y al enriquecimiento con nitrógeno, las cuales pueden alcanzar abundancias elevadas cuando las especies más sensibles desaparecen (Claerhout et al., 2026). Esto explicaría por qué la comunidad de Barva estuvo dominada casi por completo por el tipo de morfotipo folioso, que posiblemente corresponda a un grupo tolerante a las condiciones urbanas (Gauslaa, 2024) y están bien adaptados a funcionar bajo la contaminación atmosférica (Styburski & Skubała, 2023). Por lo tanto, la abundancia por sí sola no parece ser un buen indicador del efecto que podría tener la urbanización, y sería necesario considerar también la composición de la comunidad. Por otro lado, que la cobertura del dosel y el tipo de corteza fueran similares entre los sitios, permitió descartar que las diferencias en luz y sustratos expliquen los patrones observados en las comunidades de líquenes. La única característica que difirió fue el DAP, probablemente porque los árboles de la UNA corresponden a individuos de mayor edad dentro del campus y el DAP de los árboles aumenta con la edad (Shoda et al., 2020). Sin embargo, como se discute más adelante, esta diferencia no se reflejó en la abundancia de líquenes, lo que refuerza la idea de que las diferencias entre las comunidades responden más a las condiciones de cada sitio que de las características de los árboles hospederos (Sebald et al., 2022). La mayor riqueza de morfotipos en el sitio de alta urbanización puede explicarse porque la diversidad de líquenes tiende a ser mayores en niveles intermedios de alteración urbana, ya que tanto las condiciones muy alteradas, como las muy húmedas y sombreadas limita el número de especies que pueden establecerse (Claerhout et al., 2026). Este patrón contrasta con lo reportado por Correa-Ochoa et al. (2021), quienes encontraron una menor riqueza de líquenes en las zonas con más material particulado. Además, esta diferencia podría deberse a que en el presente estudio no se midieron directamente los contaminantes y a que el sitio de baja urbanización presenta una mayor altitud y humedad (Saiz et al., 2021). Por otro lado, que la disimilitud entre sitios se debiera principalmente al recambio, indica que la urbanización no sólo reduce el número de morfotipos sino que reemplaza a la comunidad por una diferente (Sebald et al., 2022). Esto también explica por qué los dos sitios urbanos compartieron casi la misma comunidad mientras que el Monte de la Cruz presentó morfotipos propios. La asociación entre el tipo de talo y el sitio refleja que cada forma de crecimiento responde de manera distinta a las condiciones del ambiente. La mayor presencia de talos crustáceos en la UNA podría deberse a que estas formas crecen fuertemente adheridas al sustrato, lo que reduce su superficie expuesta y les permitiría resistir mejor la desecación propia de los ambientes urbanos (Claerhout et al., 2026). En Barva, la asociación con los talos foliosos se explica por la dominancia del morfotipo tolerante discutida anteriormente. Por su parte, la presencia exclusiva de talos gelatinosos en el Monte de la Cruz es especialmente relevante, ya que estos suelen tener cianobacterias como fotobionte, y los cianolíquenes son muy sensibles a la contaminación atmosférica, en particular al dióxido de azufre (Gunawardana et al., 2023), por lo que su presencia sugiere condiciones ambientales más favorables en este sitio. En cambio, la falta de asociación entre el tipo de talo y el tipo de corteza contrasta con lo reportado por De Freitas y Ribeiro (2024), quienes encontraron líquenes foliosos en cortezas rugosas. Esto podría deberse a que los tipos de cortezas se distribuyeron de forma similar entre los sitios, de modo que el efecto del sitio predominó sobre el sustrato. La mayor abundancia en Barva puede explicarse por la dominancia de un morfotipo tolerante, mientras que la menor abundancia en el Monte de la Cruz, aun siendo el sitio con condiciones ambientales más favorables, podría deberse a la competencia con briófitos. En este sitio se observó musgo con frecuencia sobre la corteza, y en ambientes más húmedos y sombreados la menor luminosidad y la competencia con briófitos pueden reducir la presencia de líquenes en los árboles (Cordero-S. et al., 2021). Que en la regresión múltiple solo el sitio explica la abundancia, y que el mejor modelo incluyera únicamente esta variable, confirma que las condiciones propias de cada lugar tuvieron más peso que las características del árbol (Sebald et al., 2022). Por otro lado, la falta de relación entre el DAP y la abundancia o la riqueza de líquenes indica que un árbol más grande no necesariamente alberga más líquenes. De Freitas y Ribeiro (2024) tampoco encontraron una relación consistente entre el DAP y la abundancia de líquenes en un bosque tropical, lo que sugiere que las condiciones del entorno pueden ser más determinantes que la superficie de corteza disponible. Por ello, la regresión con el DAP no permitió realizar predicciones útiles sobre la abundancia de líquenes. En el caso de la cobertura del dosel, la ausencia de correlación podría explicarse porque la mayoría de los forófitos presentaban coberturas altas y similares, de modo que la variación en luminosidad entre árboles fue insuficiente para influir en las comunidades líquenes. Estos resultados se pueden interpretar considerando algunas limitaciones. En primer lugar, no se midieron directamente los contaminantes atmosféricos, por lo que las diferencias observadas no pueden atribuirse exclusivamente de la contaminación, como sí lo hicieron otros estudios que sí relacionaron las comunidades de líquenes con el material particulado (Correa-Ochoa et al., 2021). En segundo lugar, el Monte de la Cruz se ubica a una elevación considerablemente mayor que los otros sitios, por lo que las diferencias de clima asociadas a la altitud podrían influir en los resultados junto con la urbanización (Saiz et al., 2021).

CONCLUSIONES

El análisis de las comunidades de líquenes en los tres sitios de estudio demostró que el grado de urbanización influyó principalmente en su composición, más que en su abundancia o riqueza. Los sitios con mayor urbanización, la Sede Omar Dengo de la UNA y Barva, presentaron más talos y morfotipos, aunque la mayor parte de los talos se concentró en uno o dos morfotipos, posiblemente tolerantes a las condiciones urbanas. Por ello, una mayor abundancia o riqueza no reflejó necesariamente una mejor calidad ambiental. Además, la disimilitud en la composición entre el Monte de la Cruz y los dos sitios urbanos se debió principalmente al recambio de morfotipos, mientras que la UNA y Barva compartieron prácticamente la misma comunidad, lo que sugiere que la urbanización impulsa el reemplazo de unos morfotipos por otros y no solo a la pérdida de morfotipos. En cuanto a las características de los forófitos, la cobertura del dosel y el tipo de corteza fueron similares entre los sitios, por lo que los árboles fueron tanto comparables en cuanto a la luminosidad y sustrato. Aunque el DAP fue mayor en la UNA que en el Monte de la Cruz, no se pudo establecer una relación significativa entre el DAP o la cobertura del dosel y la abundancia o riqueza de líquenes, y el DAP no permitió predecir la cantidad de talos. Asimismo, el tipo de corteza no se asoció con los tipos de talo, y al considerar en conjunto el sitio, el DAP y el tipo de corteza, únicamente el sitio explicó la abundancia de líquenes. Esto indica que las condiciones propias de cada sitio tuvieron mayor peso que las características individuales de los árboles hospederos. Del mismo modo, se encontró una asociación entre los tipos de talo y el grado de urbanización, lo que refleja que cada forma de crecimiento responde de manera distinta a las condiciones ambientales de cada sitio. Los líquenes crustáceos estuvieron sobrerrepresentados en la UNA y los foliosos en Barva, mientras que los gelatinosos únicamente en el Monte de la Cruz, el sitio menos urbanizado. Finalmente, estos resultados respaldan el uso de los líquenes cortícolas como indicadores de las condiciones ambientales, siempre que se considere la composición de la comunidad y no únicamente su abundancia. Para futuros estudios se recomienda medir directamente los contaminantes atmosféricos, seleccionar sitios con elevaciones similares para separar el efecto de la urbanización de la altitud, incluir más de un sitio por nivel de urbanización e identificar los morfotipos a nivel de especie.

REFERENCIAS

Abas, A. (2021). A systematic review on biomonitoring using lichen as the biological indicator: A decade of practices, progress and challenges. Ecological Indicators, 121, 107197. https://doi.org/10.1016/j.ecolind.2020.107197

Bartels, S. F., Chen, H. Y. H., Wiersma, Y. F., & Macdonald, S. E. (2015). Dynamics of epiphytic macrolichen abundance, diversity and composition in boreal forest. Journal of Applied Ecology. https://doi.org/10.1111/1365-2664.12360

Cáceres, M. E. S., Lücking, R., & Rambold, G. (2007). Phorophyte specificity and environmental parameters versus stochasticity as determinants for species composition of corticolous crustose lichen communities in the Atlantic rain forest of northeastern Brazil. Mycological Progress, 6(3), 117-136. https://doi.org/10.1007/s11557-007-0532-2

Claerhout, T., Sparrius, L. B., Keßler, P. J. A., & Stech, M. (2026). Urban heat Islands shape epiphytic communities of lichens and bryophytes. Urban Ecosystems, 29(2), 63. https://doi.org/10.1007/s11252-026-01930-8

Cordero S., R. A., Garrido, A., Pérez-Molina, J. P., Ramírez-Alán, Ó., & Chávez, J. L. (2021). Lichen community structure and richness in three mid-elevation secondary forests in Costa Rica. Revista de Biología Tropical, 69(2), 688-699. https://doi.org/10.15517/rbt.v69i2.46162

Correa-Ochoa, M. A., Vélez-Monsalve, L. C., & Saldarriaga-Molina, J. C. (2021). Spatial distribution of lichen communities and air pollution mapping in a tropical city: Medellín, Colombia. Revista de Biología Tropical, 69(3), 1107-1123. https://doi.org/10.15517/rbt.v69i3.46934

De Freitas, G. V., & Ribeiro, D. R. (2024). Effect of host tree traits on the abundance of epiphytic lichens in a tropical Atlantic rainforest. Community Ecology, 25(2), 237-247. https://doi.org/10.1007/s42974-024-00198-2

Fox, J., & Weisberg, S. (2019). An R companion to applied regression (3rd ed.). Sage Publications. https://www.john-fox.ca/Companion/

Gauslaa, Y. (2024). Changes in epiphytic lichen diversity along the urban-rural gradient before, during, and after the acid rain period. Biodiversity and Conservation, 33(8-9), 2247-2263. https://doi.org/10.1007/s10531-024-02871-4

Gunawardana, D., Wanigatunge, R. P., Wewalwela, J. J., Vithanage, M., & Wijeyaratne, C. (2023). Sulfur is in the Air: Cyanolichen Marriages and Pollution. Acta Biotheoretica, 71(3), 14. https://doi.org/10.1007/s10441-023-09465-7

Hawksworth, D. L., Iturriaga, T., & Crespo, A. (2005). Líquenes como bioindicadores inmediatos de contaminación y cambios medio-ambientales en los trópicos. Revista Iberoamericana de micología, 22(2), 71-82.

Koch, N. M., Matos, P., Branquinho, C., Pinho, P., Lucheta, F., Martins, S. Ma. D. A., & Vargas, V. Ma. F. (2019). Selecting lichen functional traits as ecological indicators of the effects of urban environment. Science of The Total Environment, 654, 705-713. https://doi.org/10.1016/j.scitotenv.2018.11.107

Komsta, L. (2022). outliers: Test for outliers (Version 0.15) [Software], CRAN. 10.32614/CRAN.package.outliers

Méndez-Estrada, V. H., & Campos, C. A. (2015). Cobertura de líquenes arborícolas y su relación con la orientación cardinal en parques municipales de la Gran Área Metropolitana de Costa Rica. 7(2), 313-317. UNED Reserch Journal https://doi.org/10.22458/urj.v7i2.1159

Municipalidad de San Rafael de Heredia. (s.f). Estudios Generales. https://www.munisrh.go.cr/estudios-generales/

Monge-Nájera, J., González, M. I., Rivas Rossi, M., & Méndez-Estrada, V. H. (2002). Twenty years of lichen cover change in a tropical habitat (Costa Rica) and its relation with air pollution. Revista de Biología Tropical, 50(1), 309–319.

R Core Team (2026). R: A language and environment for statistical computing. R Foundation for Statistical Computing. https://www.R-project.org

Romero-Vargas, Marilyn, Bermúdez-Rojas, Tania, Durán-Apuy, Alejandro, Alfaro-Sánchez, Marvin, & Bonilla-Soto, Sebastián. (2022). Áreas verdes urbanas, una caracterización paisajística y biológica aplicada a una microcuenca de la Gran Área Metropolitana de Costa Rica. Revista Geográfica de América Central, (69), 23-74. https://dx.doi.org/10.15359/rgac.69/2.1

Saiz, H., Dainese, M., Chiarucci, A., & Nascimbene, J. (2021). Networks of epiphytic lichens and host trees along elevation gradients: Climate change implications in mountain ranges. Journal of Ecology, 109(3), 1122-1132. https://doi.org/10.1111/1365-2745.13538

Sebald, V., Goss, A., Ramm, E., Gerasimova, J. V., & Werth, S. (2022). NO₂ air pollution drives species composition, but tree traits drive species diversity of urban epiphytic lichen communities. Environmental Pollution, 308, 119678. https://doi.org/10.1016/j.envpol.2022.119678

Shoda, T., Imanishi, J., & Shibata, S. (2020). Growth characteristics and growth equations of the diameter at breast height using tree ring measurements of street trees in Kyoto City, Japan. Urban Forestry & Urban Greening, 49, 126627. https://doi.org/10.1016/j.ufug.2020.126627

Somarribas Chavarría, L., Rodríguez Soto, F., & León Arce, R. (2011). Reordenamiento vial e impacto ambiental del transporte en la ciudad de Heredia. Revista Geográfica De América Central., 1(38), 27-41. https://www.revistas.una.ac.cr/index.php/geografica/article/view/1721

Styburski, J., & Skubała, K. (2023). Do urban air pollutants induce changes in the thallus anatomy and affect the photosynthetic efficiency of the nitrophilous lichen Physcia adscendens? Environmental Science and Pollution Research, 30(52), 112336-112346. https://doi.org/10.1007/s11356-023-30194-4

Suqui, A., Célleri, R., Crespo, P., & Carrillo-Rojas, G. (2021). Interacciones entre índice de área foliar, densidad del dosel y precipitación efectiva de un bosque de Polylepis reticulata ubicado en un ecosistema de páramo. La Granja, 34(2), 63-79. https://doi.org/10.17163/lgr.n34.2021.04

Varela, Z., López-Sánchez, G., Yáñez, M., Pérez, C., Fernández, J. A., Matos, P., … & Aboal, J. R. (2018). Changes in epiphytic lichen diversity are associated with air particulate matter levels: The case study of urban areas in Chile. Ecological Indicators, 91, 307-314.

Wickham, H. (2016). ggplot2: Elegant graphics for data analysis. Springer. https://doi.org/10.1007/978-3-319-24277-4

Wickham, H., François, R., Henry, L., Müller, K., & Vaughan, D. (2026). dolyr: A grammar of data manipulation (Versión 1.2.1) [Paquete de R]. CRAN. https://doi.org/10.32614/CRAN.package.dplyr

Wickham, H., & Bryan, J. (2023). readxl: Read Excel files (Version 1.4.3)[Software]. CRAN. 10.32614/CRAN.package.readxl

Zárate-Arias, S., Moreno-Palacios, M., & Torres-Benítez, A. (2019). Diversidad, especificidad de forófito y preferencias microambientales de líquenes cortícolas de un bosque subandino en la región Centro de Colombia. Revista de la Academia Colombiana de Ciencias Exactas, Físicas y Naturales, 43(169), 737-745. https://doi.org/10.18257/raccefyn.886

Zeileis, A., & Hothorn, T. (2002). Diagnostic checking in regression relationships. R News, 2(3), 7-10.https://CRAN.R-project.org/doc/Rnews/