Holoregmia Nees (Martyniaceae), a monospecific genus is an iconic example of Caatinga’s ancient endemic genera. First collected by Prince Maximilian of Wied in 1817 and formally described by Nees in 1821 (Nees 1821), Holoregmia viscida is a robustly branched shrub attaining a height of up to 3 m (Harley et al. 2003), thus unlikely to be overlooked by botanical collectors in the field. However, due to the incomplete and fragmentary nature of the original material, the correct application of the name Holoregmia and the distinctive nature of this monotypic genus remained in doubt until the early 21^st century. Late 19^th century authors considered Holoregmia doubtfully distinct as a genus and probably a synonym of Craniolaria L. A revision of the Martyniaceae by Van Eseltine in 1929 (Van Eseltine 1929) omitted Holoregmia entirely, reflecting the lack of information about this Caatinga endemic genus (Harley et al. 2003).

Although documentation of the Brazilian flora advanced rapidly in the second half of the 20^th century (Forzza et al. 2010; Morim and Nic Lughadha 2015), specimens of H. viscida collected at that time were often mistakenly identified under other botanical families or left unidentified among the Martyniaceae family specimens, including several collections from the 1980s onward deposited in the HUEFS herbarium (Universidade Estadual de Feira de Santana), which is now reported to have the largest collection of H. viscida species in Brazil (Giulietti and Harley 2013). The lack of accessible, complete, and authoritatively identified reference material in herbaria prevented the recognition of this peculiar genus as a distinct endemic genus. However, new field collections in the early 21^st century were finally recognized as Holoregmia, prompting a review of unidentified plant material in herbaria, including the 1980s collections, resulting in a detailed morphological description of this monospecific genus (Harley et al. 2003).

Holoregmia viscida has long racemes of zygomorphic flowers borne well above the leaves. The calyx is foliaceous, pale green; the tubular corolla is fleshy, dull-yellow mottled with brown and darker at the throat. The two stamens bear cream anthers, and the stigma is prominent (Fig. 1). The fruit is a rounded, two-seeded capsule. The plant becomes dormant during dry periods, sprouting with the first rains. It can be found flowering between October and April, fruiting from January to August (Sinzinando Albuquerque and Daniela Zappi pers. comm.). No direct observations concerning pollination are available, but pollination by large bees has been suggested as probable (Harley et al. 2003). The same authors observed that a trait considered typical of Martyniacae, a pair of horns borne at the apex of the endocarp, is almost imperceptible in Holoregmia, leading to the conclusion that the fruit morphology does not suggest epizoochory, which is characteristic of other family members.

Figure 1. 

Holoregmia viscida Nees showing its inflorescences and details of the flowers. Photos by Sinzinando Albuquerque.

We searched digital biological databases for records of Holoregmia viscida both under its current name and under the synonyms Craniolaria unibracteata Nees & Mart., Martynia spathacea Spreng. and Proboscidea unibracteata (Nees & Mart.) Decne. referenced in Flora e Funga do Brasil (Rossini and Gonzaga 2020 – https://floradobrasil.jbrj.gov.br/FB25778, data collected in June 2021) and Plants of the World Online (Govaerts et al. 2021 – https://powo.science.kew.org/taxon/urn:lsid:ipni.org:names:327210-2, data collected in June 2021). We obtained 41 occurrence records for Holoregmia viscida from the Brazilian databases Species Link (https://splink.cria.org.br/, data collected in June 2021), Reflora Virtual Herbarium (http://reflora.jbrj.gov.br/reflora/herbarioVirtual/, data collected in June 2021), and the international database Global Biodiversity Information Facility (http://www.gbif.org/, data collected in June 2021). We checked the taxonomic identification for each of these records by consulting the relevant specimens in the Herbarium of the State University of Feira de Santana (HUEFS), which houses most of these collections, or examining digital images available in the above-named digital collections. Records were cleaned for geographic errors using the R package CoordinateCleaner v.2.0-20 (Zizka et al. 2019) and corrected for spatial sampling bias using the R package spThin v.0.2.0 (Aiello-Lammens et al. 2015). Our final dataset included a total of 34 georeferenced, taxonomically validated records (Supplementary file 1).

To model the potential niche of H. viscida, we used 20 variables as predictors for current and future conditions. Of these, 19 were bioclimatic layers of temperature and precipitation, obtained from WorldClim 2.1 (Fick and Hijmans 2017; http://www.worldclim.org/) and one was the elevation variable, obtained from NASA Shuttle Radar Topographic Mission (SRTM; https://lpdaac.usgs.gov/products/srtmgl1v003/) with 5 arc‐min (0.083° ≈ 10 km) resolution for South America. For the projection of the climatic niche in paleoclimatic scenarios, we used 19 bioclimatic variables obtained from the PaleoClim database (Brown et al. 2018; http://www.paleoclim.org/), and to project future global climate change scenarios, we used the models available in the WorldClim database.

We performed a principal component analysis (PCA) on the current environmental variables based on a correlation matrix to overcome multicollinearity problems and reduce the number of predictor variables (Velazco et al. 2019), avoiding model overfitting that may result in biologically unreliable areas (Jiménez-Valverde et al. 2011). We used the criterion that explains at least 95% of the total variance of the correlation matrix and includes as much environmental information as possible (De Marco and Nóbrega 2018; Destro et al. 2019). Thus, the first six axes derived from the PCA were the new set of variables (Supplementary files 2, 3). Then, we used the PCA’s eigenvectors to calculate scores of each derived PC as new predictors (Velazco et al. 2019).

For paleoclimatic projection, we selected three scenarios: Last Interglacial (ca 130 ka), Last Glacial Maximum (ca 21 ka), and Mid-Holocene (8.326–4.2 ka). We did not include the elevation variable in the paleoclimatic projection. To model future climatic scenarios, we used the climate projection from the Intergovernmental Panel on Climate Change (AR6) as the source of future climate conditions. We evaluated the effect of climate change using three Shared Socioeconomic Pathways (SSPs) based on Coupled Model Intercomparison Project Phase 6 (CMIP6). We selected three scenarios: SSP126 is the updated scenario of RCP2.6, which is an optimistic scenario that would occur in a sustainable world that takes a green path with low greenhouse gas emissions; SSP245 is the updated RCP4.5 scenario that uses an intermediate level of greenhouse gas emissions; and SSP585 is the updated scenario of RCP8.5, which is the most pessimistic and assumes that use of fossil fuels will continue to increase (O’Neill et al. 2014, 2017) forcing a high level of greenhouse gas emissions (Gillett et al. 2016). We used projections for 2050 (mean for the period from 2041 to 2060) and 2090 (mean from 2081 to 2100) using the General Circulation Model IPSL-CM6A-LR.

Ecological niche models (ENMs) are empirical or mathematical approximations of the ecological niche of a species (Barbosa et al. 2012). Based on different statistical approaches, a plethora of algorithms are available to predict the niche of a species (Velazco et al. 2021). Based on the occurrence records of the target species associated with the set of climatic variables, the algorithms create a multidimensional environmental space that infers the niche of the species (Araújo and Guisan 2006; Carvalho et al. 2017). Species occurrence records can present possible errors and biases in the data, so it is important to clean and filter the data (Zurell et al. 2020). However, algorithms are the primary sources of uncertainty in the models. Thus, testing the effect of different algorithms on modelling is essential to detect the best model (Thuiller et al. 2019). We adopted a routine with six different algorithms to evaluate niche modelling (Supplementary file 4): i) Bioclimatic Envelope Method (Booth et al. 2014); ii) Maxent, with default features – MaxNet (Phillips et al. 2006, 2017); iii) Simple Maxent Model (Williams 2010); iv) Support Vector Machine (Tax and Duin 2004); v) Random Forest – RDF (Breiman 2001); and vi) Gaussian Model (Golding and Purse 2016) (see individual models in Supplementary file 5).

We tested the performance of our ecological niche modelling using spatial block cross-validation, based on a geographically structured checkerboard to control the potential spatial autocorrelation between training and testing data (Roberts et al. 2017). The test and training data were used to test for spatial autocorrelation (Moran’s I) and environmental similarity (MESS – Multivariate Environmental Similarity Surface) (Destro et al. 2020). Moran’s I ranges from of -1 to 1, distinguishing between negative and positive spatial autocorrelation, respectively (Diniz-Filho et al. 2003). MESS analysis identifies and determines the extent of the environmental differences between model training and model projection (predictions are extrapolations) data (Elith et al. 2010), where the higher the value, the more common the environment of the training point (Camera et al. 2017). Negative values indicate a novel environment, i.e. different environmental locations from the reference region.

Models that use presence-only data are strongly influenced by bias, while models that use presence/pseudo-absence data are not (Phillips et al. 2006). Thus, we selected the method in which pseudo-absences are environmentally constrained to a region with suitability values predicted by a BIOCLIM model (Wisz and Guisan 2009). First, we set the threshold based on the maximization of sensitivity-specificity sum (Liu et al. 2005). Then, an ensemble was created for the species, selecting the best models with True Skill Statistic (TSS) above the average (Marmion et al. 2009; Velazco et al. 2019). We used the Area Under the Curve (AUC) and the True Skill Statistic (TSS) to evaluate the models’ accuracy. The AUC is a threshold-independent index, where values of AUC ≥ 0.8 are considered good to excellent (Komac et al. 2016) to discriminate the places occupied by a species and places where occupation is unknown. The TSS is a threshold-dependent index, where TSS values > 0.6 are considered good to excellent (Komac et al. 2016) and indicate greater agreement between the observed and predicted distribution of the species (Allouche et al. 2006). The models were constructed using the R package ENMTML v.1.0 (Andrade et al. 2020). We established four classes from the threshold set by TSS, ranging from 0 to TSS; TSS to 0.60; 0.61 to 0.80; 0.81 to 1. We also calculate the area suitable for each scenario from the threshold. Our ENM were generated projecting the niche to the Neotropics (all models are available for download in our supplementary data). We calculated suitable areas for the species within the Brazilian territory and for South America. Those interested in the result for the entire continent can recalculate the data using the rasters available for download.

For paleoclimate scenarios, the ensemble model showed a satisfactory performance with AUC 0.95 with a standard deviation of ± 0.06, and TSS values of 0.84 with a standard deviation of ± 0.2 (Supplementary file 6). For the climate change scenarios, the ensemble model showed a satisfactory performance with AUC 0.92 with a standard deviation of ± 0.005 and TSS values 0.81 with a standard deviation of ± 0.06 (Supplementary file 6). Moran’s I values for paleoclimate and climate change scenarios were 0.008 and 0.033 (Supplementary file 7), respectively, indicating low spatial autocorrelation of environmental variables. All values were positive, indicating that the training-test subregions presented similar environments, i.e. they do not suffer from induced low precision due to subsets of data being in climate-discrete regions or questions related to model extrapolation (Peterson et al. 2007; Owens et al. 2013).

The suitable areas for the occurrence of Holoregmia viscida have oscillated in paleoclimatic scenarios. Our results showed that the appropriate environmental conditions during the Last Interglacial were distributed mainly in the Atlantic Forest Domain and to a lesser extent in the south of the Caatinga, resulting in an area of 135,919 km² for Brazil and 212,612 km² for South America (Fig. 3). These suitable areas for the occurrence of H. viscida increased to 198,488 km² for Brazil and 273,621 km² for South America in the Last Glacial Maximum. They were displaced within the boundaries of the Caatinga Domain with small areas in the Cerrado Domain (Fig. 3). In the Mid-Holocene, the geographical extent of suitable areas was reduced to 82,394 km² for Brazil and 153,971 km² for South America, becoming restricted to the south of the Caatinga Domain (Fig. 3). Thus, the adequate areas during the Holocene were likely limited, representing a fraction of Holoregmia’s current total suitable area. This loss explains why Holoregmia is now limited to the southern Caatinga, even though suitable areas are also found in parts of the northern Caatinga (Fig. 3). For the current climate conditions, the models estimated an area of 156,372 km² for Brazil and 215,953 km² for South America that includes the regions with the highest suitability for the occurrence of Holoregmia viscida. The focus of all maps in Fig. 3 is centered in the Caatinga dry forests, but climatically suitable areas for H. viscida can also be found outside this frame. Full data for all models is available as raster files on Figshare scientific repository (see Data Availability Statement).

Figure 3. 

Map zooming in the Caatinga Domain showing the climatic suitability for Holoregmia viscida in current conditions (centre of the figure), paleoclimatic scenarios (top of the figure) and three future climate change scenarios (bottom of the figure). Brazilian phytogeographical domains: Am – Amazon Rainforest; At – Atlantic Rainforest; Ca – Caatinga; Ce – Cerrado. TSS represents the True Skill Statistic, used to define the best model and the threshold value in the model to decide whether a pixel has suitable conditions for the occurrence of the species. SSPs are the Shared Socioeconomic Pathways, which are possible scenarios regarding the future emission of greenhouse gases. The SSP126 is an optimistic scenario that would occur in a sustainable world with low greenhouse gas emissions; SSP245 considers an intermediate level of greenhouse gas emissions; and SSP585 is the most pessimistic scenario and assumes that the use of fossil fuels will continue to increase.

The species' niche will be affected by future climatic conditions, even under the most optimistic scenarios (Fig. 3). The SSP126 scenario predicts that the niche will be reduced to 84,928 km² for Brazil (i.e. 54.31% of the current range) and 109,124 km² for South America (i.e 50.53% of the current range) in 2050 and 79,418 km² (i.e 50.78% of the current range) for Brazil and 99.909 km² for South America (i.e 46.26% of the current range) in 2090. Niche losses tend to increase under the SSP245 scenario, with an average reduction to 67,248 km² (i.e 43% of the current range) for Brazil and 83,522 km² (i.e 38.67% of the current range) for South America in 2050 and 30,309 km² (i.e 19.38% of the current range) for Brazil and 37,262 km² (i.e 17.25% of the current range) for South America in 2090. The SSP585 scenario estimates more significant niche losses than the other scenarios, reducing to 37,274 km² (i.e 23.83% of the current range) for Brazil and 51,151 km² (23.68% of the current range) for South America in 2050 and a total loss in 2090 for Brazil and 1,341 km² (i.e 0.62% of the current range) for South America . Under this latter scenario, it is possible that H. viscida will become extinct, as no suitable areas will remain (Fig. 3).

Martyniaceae is a New World family primarily associated with dry ecosystems (Gormley et al. 2015), suggesting the lineages conserve their niche adapted to dry areas. Considering that H. viscida is estimated to have diverged from its sister taxa 9.4 million years ago (Queiroz et al. 2017), this lineage has undoubtedly been subjected to multiple climatic fluctuations throughout its evolutionary history. Currently, the genus is limited to the southern part of the Caatinga. This restricted occurrence is likely the result of the paleoclimatic fluctuations, as discussed below. The potential niche of H. viscida appears to have expanded from the Last Interglacial to the Last Glacial Maximum, however, models show a considerable niche loss in the Holocene. In addition, the potential niche in the Holocene was much smaller and restricted compared to the current niche.

In the current scenarios, the potential niche goes beyond the southern Caatinga, since suitable areas are also found in the northern Caatinga today. For future conditions, the effects of climate change on the ecological niche of H. viscida may be irreversible. The models showed that the species tends to lose suitable areas in optimistic scenarios and especially in pessimistic scenarios, with the possible extinction of the species in SSP585 by the year 2100.

Our paleoecological reconstructions show a major displacement of H. viscida throughout the Quaternary. While in the Last Interglacial (Fig. 3), suitable niches for the species occurred in the northeastern part of the Caatinga, and also in areas within the Cerrado and the Atlantic Forest Domains; during the Last Glacial Maximum (LGM) (Fig. 3), suitable niches appear to have expanded to almost the entire extent of the Caatinga. The model showed that during the LGM, the suitable niches were within the current limits of the Caatinga. The Mid-Holocene had the smallest projected suitable niches of the paleoclimatic scenarios analysed here. In this case, the potentially suitable niches for the species remained restricted to small patchy areas in the southern Caatinga and ecotonal areas between Caatinga and Atlantic Forest and drier areas of the Atlantic Forest.

During the Mid-Holocene in the Caatinga, climatic conditions were characterized by increased precipitation, which likely affected the range of optimal climatic conditions for the species. For example, pollen records of the genus Ziziphus indicate that in the period from 10,000 to 6,000 years before present (BP), there were taxa adapted to high humidity conditions. However, drier conditions are assumed to have occurred between 6,400 and 1,800 BP based on the lack of sediment deposits in this timeframe (Medeiros et al. 2018). Similarly, the vegetation history of the Caatinga based on peat-bog sediments deposited over > 10,000 years showed that during 10,540–6790 BP, there was a warmer and wetter period evidenced by the prevalence of Mauritia pollen, with the establishment of current climatic conditions and vegetation patterns after 4535 BP (Oliveira et al. 1999).

Overall, all results show a great variation in the climatic suitability for H. viscida since the Last Interglacial. Its niche exhibits alternation between periods of contraction and expansion and shifts in its geographical placement. Events such as this are likely to leave a signature in the phylogenetic and/or phylogeographic patterns of individual taxa (see Hewitt (2004) for a review). There are many determining factors for the ecological niche of a species (Grinnell 1917; Hutchinson 1957), but the scenario presented here could be helpful as an initial hypothesis explaining the current distribution of H. viscida. Given the reduction and subsequent expansion of populations from the LGM onwards, we expect that phylogeographical studies will likely show a signature of a recent population bottleneck in H. viscida. Such a bottleneck could have severe implications for the adaptation of the species in the face of future climate change scenarios since such events generally reduce species’ genetic diversity (Bellard et al. 2012), making them less likely to survive disturbances.

Our results indicate that all projected future scenarios involving an increase in temperature and a reduction in precipitation would lead to a major decrease in areas that are climatically suitable for the survival of Holoregmia viscida. Such a decrease in the species niche is generally accompanied by a shift in its distribution towards more southern latitudes in areas that currently harbour the Atlantic Rainforest. The worst-case scenario (SSP585 2080–2100) indicates the complete disappearance of suitable climatic conditions for H. viscida, suggesting a possible extinction of the species within 80 years.

A shift of Caatinga towards the Atlantic Forest suggests that climatic conditions related to a combination of temperature and precipitation will be even harsher in the future. Similar results, where a shift in the distribution of Caatinga species towards areas that currently encompass wetter habitats have already been reported for other endemic plants (e.g. Silva et al. 2019) and for the iconic and endemic Caatinga cacti (e.g. Simões et al. 2020).

Changes in the limits between Caatinga, Cerrado, and humid forests were detected during past climatic fluctuations (Costa et al. 2017). We can infer a similar situation by observing Holoregmia’s potential niche shifting to areas further south than today, reaching areas now occupied by the Atlantic Forest. This suggests that there will be significant changes in the wetter vegetation adjacent to what is now the Caatinga Domain. This pattern is consistent with Wang et al. (2019), who demonstrated that Tropical and Subtropical Moist Broadleaf Forests in neotropical realms would be more vulnerable to climate change in terms of loss in niches. For example, models projecting suitable areas for Eschweilera tetrapetala, a tree species endemic to humid forests of the highlands of Bahia, Brazil, predict a strong reduction in suitable areas for the species by 2070 (Menezes et al. 2021).

Therefore, Holoregmia viscida models under global warming scenarios indirectly indicate that areas within the Caatinga Domain may become even drier than they currently are. In contrast, wetter areas within the Atlantic Forest may be under drier and seasonal climates, suitable for species from the Caatinga. While shifts in the biome limits and biotas occurred during Pleistocene climatic fluctuations (Costa et al. 2017; Ledo and Colli 2017), anthropogenic climate change is occurring more rapidly than past natural variations. In the past, slow climate changes happened over thousands of years and allowed the migration of many taxa. Now, with rapid changes expected to occur over decades, associated with the severe deforestation and fragmentation to which biomes are already subjected, many species might not be able to migrate. The Caatinga, for example, has already lost half of its area (Brazil 2015; Antongiovanni et al. 2018), and the majority of the Atlantic Forest has been removed (Brazil 2015).