Área de estudio

El área de estudio comprende el sur de los Estados Unidos de América, México y Centroamérica (Figura 1), región que cubre 3 946 177.4 km² y que culturalmente es conocida como las regiones de Mesoamérica y Aridoamérica (^{Kirchhoff, 1943}).

Figura 1 Área de estudio para Phaseolus lunatus.  

Datos de ocurrencia

Se integró una matriz de datos de ocurrencia de Phaseolus lunatus a partir de los datos pasaporte georreferenciados con base en la información de las siguientes entidades: Centro Internacional de Agricultura Tropical (CIAT), Instituto Nacional de Ecología, Universidad de Guadalajara, Instituto Nacional de Investigaciones Forestales, Agrícolas y Pecuarias (INIFAP), ^{Fraitag y Debouck (2002}), así como datos particulares facilitados personalmente por el Dr. Rogelio Lépiz Ildefonso. Profesor-Investigador, CUCBA, Universidad de Guadalajara.

Variables climáticas

Se utilizaron imágenes mensuales normales de temperatura máxima (12), temperatura mínima (12), temperatura media (12) y precipitación (12), además de imágenes de las 19 variables bioclimáticas, todas derivadas del sistema Worldclim y correspondientes al período 1960-2000 (^{Hijmans et al., 2005}). Además, se utilizaron 12 imágenes mensuales normales de evapotranspiración potencial anual (ETP) que se obtuvieron de la base de datos global geoespacial de balance hídrico de suelo del CGIAR (^{Trabuco y Zomer, 2010}). También se agregaron como variables de análisis la longitud, la latitud y la altitud de los sitios de presencia. Para ello se tomó la imagen de altitud del sitio USGS/NASA SRTM (^{Jarvis, 2008}). Con los datos de precipitación y ETP mensual se generaron 12 imágenes mensuales de índice de humedad (IH) mediante la expresión IH= P/ETP. Toda esta información se trabajó con una resolución dada por un tamaño de celda de 30 segundos de arco. En total se analizaron 94 variables.

Caracterización climática de sitios de presencia de Phaseolus lunatus

Se utilizó la información de coordenadas geográficas de los sitios de ocurrencia de frijol lima y las imágenes climáticas y bioclimáticas del sistema Worldclim para caracterizar la climatología de los sitios de distribución de esta especie. Para ello se utilizó el procedimiento de extracción de información del programa ArcGis (ESRI, 2010). La información extraída se incorporó a la matriz de datos de ocurrencia de P. lunatus en la hoja de cálculo Excel de Microsoft. Finalmente se caracterizó el tipo climático de los sitios de presencia mediante el sistema que proponen ^{Ruiz et al. (2018}) y que refieren de manera ilustrativa ^{Sánchez et al. (2018}), en el que son teóricamente posibles 36 ambientes, producto de la combinación de seis zonas hídricas, determinados por intervalos de un índice de humedad (IH= precipitación anual/evapotranspiración potencial anual) e intervalos de temperatura media anual (Tma), respectivamente (Cuadro 1).

Análisis estadístico

A partir de la matriz de sitios de ocurrencia caracterizada climáticamente, se realizó un análisis de correlación múltiple para depurar la lista de variables climáticas. Se utilizó el criterio de coeficiente de correlación ≥ 0.95 para la eliminación de variables. Con la lista depurada de variables se procedió a realizar un análisis de componentes principales e identificar las variables de mayor peso estadístico en la presencia de frijol lima. Para este análisis se empleó el programa NTsys 3.3.0 (^{Rohlf, 2009}).

Descriptores climáticos

A partir, de las variables climáticas seleccionadas para el análisis y considerando la caracterización climática de la matriz de datos de ocurrencia, se dedujeron los descriptores ecológicos de frijol lima, para lo cual se consideraron los valores mínimos, máximos y promedio de cada variable.

Modelación de nichos potenciales de distribución de Phaseolus lunatus

Previo a la modelación de nichos potenciales de distribución, se aplicó una prueba de aleatoriedad para determinar si los sitios de ocurrencia de Phaseolus lunatus están distribuidos al azar (^{Bivand et al., 2014}). Esta prueba se corrió utilizando un programa implementado en R. Se utilizó el modelo Maxent (^{Phillips et al., 2006}) para modelar la distribución potencial de frijol lima bajo el escenario climático actual o de referencia y bajo los escenarios climáticos año 2050 y año 2070.

Para la modelación de la distribución potencial de P. lunatus en escenarios futuros o de cambio climático se emplearon los datos derivados de los modelos de circulación general: HadGEM2ES (European Network for Earth System-MetOffice Hadley Centre) y CCSM4 Community Earth System Model (CCSM4) bajo dos rutas de concentración de gases efecto invernadero (RCP) 4.5 y 8.5. El primer modelo, desarrollado en el Reino Unido, incorpora elementos que permiten modelar adecuadamente la vegetación continental pues incluye diversos tipos de vegetación: árboles de hoja ancha, arboles de acículas, hierbas C₃ y C₄, así como arbustos (^{Bellouin et al., 2007}).

El segundo modelo está desarrollado en Estados Unidos de América por el National Center for Atmospheric Research (NACAR) el cual incorpora la cobertura vegetal como uno de sus cuatro componentes acoplados. Estos modelos han sido evaluados estadísticamente por ^{Conde et al. (2011}), concluyen que el Modelo CCSM es uno de los que mejor describe el clima de la región de México y Centroamérica, en tanto que HadGEM ya ha sido utilizado para estudios de la vegetación de México con buenos resultados (^{Instituto de ecología, 2009}).

El modelo Maxent ha sido descrito como especialmente eficiente en el manejo de interacciones complejas entre variables de respuesta y variables predictoras por ser un modelo robusto aún en muestras pequeñas de datos (^{Fourcade et al., 2014}). Maxent utiliza el principio de máxima entropía sobre datos de presencia de la especie para estimar un grupo de funciones que relacionen variables ambientales con idoneidad de hábitat para aproximar la distribución geográfica potencial de la especie (^{Phillips et al., 2006}). Maxent representa la distribución de una especie; a través, de una función de probabilidad P sobre un conjunto X de sitios en el área de estudio. La función de probabilidad para i variables ambientales se expresa de la siguiente manera.

𝑃(𝑥)=𝑒𝜆∙𝑓(𝑥)/𝑍𝜆. Donde: P(x)= función de probabilidad; λ= vector de coeficientes de ponderación; f= vector correspondiente de funciones de variables ambientales; Z= constante de normalización para asegurar que P(x)= 1. Los valores P(x) así obtenidos, representan valores de idoneidad relativa para la presencia de la especie, constituyendo así la base de un modelo de distribución potencial.

En la modelación de la distribución potencial de Phaseolus lunatus con Maxent se utilizaron los parámetros por defecto, mediante la réplica de ejecución crossvalidate con 1 000 repeticiones. Se destinó 50% de los datos para entrenamiento y 50% para prueba. En escenarios de cambio climático, esta configuración se mantuvo. Finalmente se generó un mapa binario presencia/ausencia de Phaseolus lunatus utilizando como método de umbral de corte, aquél que presentara una tasa de omisión cercana a 0% pero con el valor de corte más alto (^{Sánchez et al., 2018}).

Adaptación climática de Phaseolus lunatus

En la Figura 2, se puede observar la distribución de Phaseolus lunatus en las distintas zonas agroclimáticas del área de estudio. Como puede observarse en el mapa de la Figura 2 y Cuadro 2, esta especie está presente en 17 de las 26 zonas agroclimáticas del área de estudio. Sin embargo, de las 17 zonas agroclimáticas en las que ocurre P. lunatus, dos corresponden a zonas áridas.

Las zonas áridas difícilmente pueden albergar el ciclo de desarrollo de una leguminosa como P. lunatus, ya que los niveles de precipitación anual en esos ambientes, es mucho menor que 500 mm, valor señalado como el umbral pluvial mínimo de adaptación de esta especie (^{Baudoin, 1989}). La presencia de P. lunatus en zonas áridas en el área de estudio muy probablemente se debe a que tiene acceso al recurso hídrico, pero no a través de precipitaciones pluviales.

Figura 2 Regiones agroclimáticas y distribución de Phaseolus lunatus. 

En el Cuadro 2 se observa que el mayor número de sitios de ocurrencia de Phaseolus lunatus se produce en la zona agroclimática subhúmeda-húmeda cálida, señalando que este es el ambiente climático que se aproxima más a las condiciones climáticas óptimas para esta planta.

Sin considerar la zona hídrica árida, P. lunatus es capaz de distribuirse en cuatro zonas hídricas (desde semiárida hasta húmeda) y cuatro zonas térmicas (desde templada a muy cálida) aunque parece preferir el intervalo de subhúmeda-húmeda a húmeda y el intervalo de semicálida a muy cálida. La distribución de esta leguminosa en condiciones semiárida y templada es definitivamente muy reducida (Cuadro 2).

Descriptores ecológicos

El Cuadro 3 muestra los descriptores ecológicos de Phaseolus lunatus a partir, de las variables climáticas más representativas. Como puede verse, aun cuando esta leguminosa tiene preferencia por ambientes semicálidos, cálidos y muy cálidos (Cuadro 2), su presencia en altitudes de 0 a 2 386 m (Cuadro 3) señala una amplia capacidad de adaptación, pues térmicamente puede distribuirse desde zonas templadas (valor mínimo de 13.2 °C de temperatura media anual, Cuadro 3) incluso en ambientes templados (Cuadro 2). El rango altitudinal encontrado coincide de manera aproximada con el reportado por ^{López et al. (2005}) de 0 a 2059 m y difiere un poco del señalado por ^{Baudoin (1989}) de 800 a 2 000 msnm, así como del reportado por ^{Debouck (2020}) de 50-2 800 m. Estas diferencias seguramente se deben a las diferentes áreas geográficas consideradas por estos autores.

De acuerdo con el Cuadro 3, el intervalo de estación de crecimiento en el que ocurre Phaseolus lunatus (83 a 365 días) coincide con el ciclo de desarrollo que ^{Debouck (2020}) reporta de 90 a 365 días para P. lunatus.

Este intervalo de 83 a 365 días de estación de crecimiento, permite deducir que esta planta se adapta a ambientes con temporada de lluvias corta tanto como a ambientes de lluvia prácticamente todo el año, con una precipitación acumulada en el período de crecimiento de 250 a 4 250 mm, por lo que algunas poblaciones de P. lunatus presentes en estas condiciones extremas de humedad podrían estar indicando la existencia de poblaciones portadoras de genes resistentes a estas condiciones, como lo sugiere ^{Maxted et al. (2013}).

De igual manera, para el caso de la temperatura, P. lunatus se adapta a estaciones de crecimiento con una temperatura mínima promedio de 8 °C tanto como a estaciones de crecimiento donde la temperatura máxima puede promediar incluso 37 °C. Las poblaciones representativas de estos valores térmicos también pueden resultar de interés por su posible genética de adaptación a temperaturas extremas. Más aún, de acuerdo con el Cuadro 3, existen poblaciones de Phaseolus lunatus que experimentan el mes más cálido con temperaturas máximas por arriba de 40 °C y otras que presentan el mes más frío con temperaturas por debajo de 0 °C, sobre todo las formas de la especie de ciclo biológico largo, de hábito de crecimiento indeterminado (^{Debouck, 2020}).

Es interesante destacar que los reportes sobre P. lunatus señalan un rango térmico de 16 a 26 °C (^{Debouck, 2020}), en el presente estudio el intervalo obtenido es más amplio, 13.2 a 29.2 °C en escala anual y 13.2 a 29.9 °C a nivel de la estación de crecimiento (Cuadro 3). Estos valores constituyen límites térmicos de adaptación no reportados hasta hoy para P. lunatus, lo cual podría indicar que esta especie ha tenido la capacidad de ir explorando áreas más extremosas térmicamente al paso del tiempo.

De igual forma, la presente investigación obtuvo un límite superior de precipitación acumulada anual de 4 250 mm, cantidad que supera reportes previos como el de ^{Debouck (2020}) de 2 800 mm. Finalmente, de acuerdo con el Cuadro 3, el índice de humedad promedio durante el ciclo de producción señala como límites inferior y superior de adaptación de P. lunatus, 0.51 y 3.09, respectivamente. De estos descriptores no existe referencia antecedente, debido a que variables como evapotranspiración, balance de humedad e índice de humedad no son típicamente parámetros considerados en la descripción de la ecología de cultivos.

Modelación de la distribución potencial de Phaseolus lunatus

La prueba de aleatoriedad reportó que los puntos de ocurrencia de Phaseolus lunatus tienen una distribución aglomerada y no una distribución aleatoria, por lo cual se estimaron las distancias entre los puntos de ocurrencia y se re-muestrearon los datos mediante el programa ArcGis (^{ESRI, 2010}), por medio de una malla construida a partir de las distancias geográficas para obtener muestras aleatoriamente distribuidas. Con estas muestras se realizó la modelación de áreas potenciales de distribución de P. lunatus.

Del total de 94 variables que se consideraron en el presente estudio, el análisis de correlación múltiple y el análisis de componentes principales permitieron reducir la lista a 29 variables con la mayor significancia en la distribución de Phaseolus lunatus. En el Cuadro 4 se muestra el valor de los Eigenvectores de las 29 variables climáticas con mayor significancia en la distribución de P. lunatus en los tres primeros componentes principales.

Se observa entonces que en el primer componente se agrupan variables que tienen que ver con el factor humedad (precipitación), en el segundo componente, variables que tienen que ver con el factor energético (temperatura) y en el tercer componente, otro tipo de variables que se derivan de la interacción de ambos factores

Al correr el análisis Maxent con las 29 variables se obtuvo un modelo cuya curva operacional o curva ROC reportó un área bajo la curva (AUC) de 0.94, que de acuerdo con ^{Araújo et al. (2005}) representa un ajuste excelente del modelo y por tanto permite una buena discriminación de áreas adecuadas contra áreas no adecuadas para Phaseolus lunatus (^{Thuiller et al., 2009}) al elaborar un mapa binomial.

El análisis Jackknife permitió identificar las variables de mayor contribución a la modelación de la distribución de Phaseolus lunatus. El Cuadro 5 muestra la contribución de las 11 variables que más aportan a este respecto. Tanto Bio04, como Prec09 y Bio16 se mantuvieron como las tres variables de mayor contribución a la distribución de esta leguminosa en todas las climatologías analizadas, tanto en la de referencia, como en las climatologías de cambio climático.

El análisis de la evolución que tendrían estas tres variables en las climatologías de cambio climático se presenta en la Figura 3. Ahí se puede ver que Bio04 tiende a incrementar en los escenarios de cambio climático. Prec09 en cambio tiende a disminuir, con un porcentaje de -1% en la climatología 2050 RCP4.5 modelo CCMS4, a 15% en la climatología 2070 RCP8.5 modelo HADGEM-ES.

Con respecto a Bio16 el modelo CCMS4 estima una reducción tanto al pasar del año 2050 a 2070 como al cambiar del RCP 4.5 a 8.5, mientras que el modelo HADGEM-ES estima un valor estable de Bio16.

Figura 3 Evolución esperada de tres variables en diversos escenarios de cambio climático. 

El mapa binario para el escenario de referencia reveló que, a gran escala, la distribución actual de Phaseolus lunatus coincide con la distribución potencial simulada por el modelo Maxent (Figura 4). Sin embargo, a menor escala son visibles áreas de distribución potencial que carecen de sitios de ocurrencia actual de esta leguminosa. Esto indica la presencia de sitios adecuados para la futura exploración y posible colecta de nuevas poblaciones de P. lunatus.

Figura 4 Distribución actual (verde) y potencial (naranja) de Phaseolus lunatus en diversas climatologías: A) 1961-2010; B) HADGEM-ES RCP4.5 2050; C) HADGEM-ES RCP4.5 2070; D) HADGEM-ES RCP8.5 2050; E) HADGEM-ES RCP8.5 2070; F) CCMS4 RCP4.5 2050; G) CCMS4 RCP4.5 2070; H) CCMS4 RCP8.5 2050; I) CCMS4 RCP8.5 2070. 

En las climatologías de cambio climático, el modelo HADGEM-ES en los dos RCPs señala un descenso en la superficie de distribución potencial de P. lunatus principalmente en la Península de Yucatán y Centroamérica (Figuras 4B, 4C, 4D y 4E) que va desde -8 a -10% en el RCP4.5 a un intervalo de -13 a -37% (Figura 5) en el RCP8.5. Estos efectos sobre la Península de Yucatán son particularmente importantes, debido a que esta zona representa uno de los centros de diversidad más importantes del grupo genético MII de frijol lima.

Figura 5 Efecto del cambio climático en la superficie potencial de distribución de P. lunatus. 

Desde el escenario 2050 RCP4.5 (Figura 4B) se advierte una disminución de áreas potenciales para esta planta, en dirección de sur a norte de la Península de Yucatán, esto es norte de Guatemala y Campeche, lo cual se va acentuando en los siguientes escenarios (Figuras 4C, 4D y 4E) abarcando cada vez mayor superficie de este territorio, hasta reducir las áreas potenciales a una pequeña zona en el extremo noreste de la Península, en el estado de Yucatán (Figura 4E).

En el caso del modelo CCMS4, las Figuras 4F, 4G, 4H y 4I muestran la secuencia de mapas de áreas potenciales de distribución del frijol lima. En ella se puede ver que el RCP4.5 representa condiciones ambientales favorables, que causarían un incremento de 7% en el año 2050 y +2% en 2070 (Figura 5), ocurriendo principalmente en la región norte-centro y occidente de México, además de El Caribe (Figuras 4F y 4G), lo cual podría en un futuro dar lugar a la presencia de nuevas poblaciones de frijol lima en estas regiones y con ello una nueva configuración de la distribución de esta especie en el área de estudio. Mientras que el RCP8.5 propiciaría una reducción en la superficie potencial de frijol lima en la Península de Yucatán y Centroamérica (Figura 4) de -0.17% en 2050 a -22% en 2070 (Figura 5). Estos resultados positivos y negativos para P. lunatus concuerdan con lo reportado para Tithonia diversifolia en México por ^{Durán et al. (2020}).