Phytoregionalisation of the Andean páramo

Vegetation data

This study was based on a dataset of 2,853 vegetation plots fitting the study area, sampled with the phytosociological method (Braun-Blanquet, 1951) and obtained from the VegPáramo database (www.vegparamo.com). All taxa were previously checked for synonymy with the VegPáramo taxon list, which contains over 15,000 plant names from the northern Andes that are frequently updated using the Plant List (www.theplantlist.org) and Tropicos (www.tropicos.org). The initial dataset included plots proceeding from 38 different sources and sampled by various authors over a period from 1981 to 2014. Because these plots were sampled in many different vegetation types, as often stated by their original author, the dataset was rather heterogeneous in terms of vegetation represented, including meadows, bogs, grasslands, giant rosette vegetation, shrublands, bamboo vegetation and forests among others. Moreover, because all plots were sampled with the phytosociological method, which instructs that plot size should vary according to vegetation physiognomy and composition (Ozenda, 1982), the dataset presented a wide range of plot sizes, from 1 to 100 m², which we believed with our data handling (see below) should not affect our analyses. Regarding floristic content, each plot included a list of vascular, and sometimes non-vascular, plant species occurrences with their cover coefficients: + <1%, 1 >5%, 2 <25%, 3 <50%, 4 <75% and 5 >75%. Finally, plots were georeferenced using the UTM coordinate system at a 1 or 10 km² resolution depending on the data source and the precision used in its plots’ geolocation.

In order to homogenize the species list, non-vascular plant species, unidentified taxa and supra-specific taxon names, i.e., genus, family, were removed. In addition, because some authors used sub-specific taxa, i.e., subspecies, variety, in their sampling and others did not, we summed all sub-specific taxa to the species level. Finally, species with low occurrence in the dataset (<2) were removed to reduce statistical noise. Because we aimed at characterising the main páramo biogeographical units under a macroclimatic and broad-scale geographical focus, we removed vegetation plots representing microclimatic or edaphic azonal vegetation, i.e., bogs, marshes, rocky vegetation or Polylepis relict forests, according to the original author’s description. Finally, to take into consideration the geographic sampling bias and varying geolocation precision, we (1) removed plots with >1 km resolution and (2) conducted a preferential geographic stratified resampling at 15 plots per 1 km UTM unit and uniformly distributed with elevation, divided into 200 m altitudinal strata (Knollová et al., 2005). These consecutive reductions lead to a final dataset of 1,647 plots for 1,724 species (Fig. 1B). Finally, to minimise the effect of subjective estimates in species cover, characteristic of phytosociological sampling, the phytosociological scale was transformed into a presence/absence binary scale (Kočí, Chytrý & Tichý, 2003).

Non-hierarchical clustering and floristic characterisation

The following statistical analyses were carried out using the Ginkgo program of the b-VegAna application set (biodiver.bio.ub.es/veganaweb; Bouxin, 2005). The two main clustering methods available, hierarchical and non-hierarchical, apply partitioning to a dataset with an increasing number of clusters, however, the former generates a tree while the latter does not. Although a hierarchical classification on vegetation plot data might seem easier to interpret than a non-hierarchical one, one disadvantage is that plots that are assigned to a cluster at division state n cannot switch to a non-directly descendent cluster at the division n + 1, making the assignment permanent and the clustering rigid (De Cáceres & Wiser, 2012). On the contrary, non-hierarchical clustering is more flexible and might not even yield the same plot assignment result for each iteration of the same partition (Tichý, Chytrý & Smarda, 2011), which makes it more difficult to interpret but more precise and fitted to the data. Because our dataset included certain diversity of vegetation types, whose proportions were probably neither balanced nor representative of real landscapes due to sampling bias, we opted for a non-hierarchical clustering analysis to generate the páramo bioregionalisation, so to ensure independence between clusters of different partitions.

First, the dataset was converted into a plot-distance matrix using the Bray–Curtis distance, a widely used coefficient for community data in ecology (Bray & Curtis, 1957; Tichý et al., 2010). Second, the matrix was classified using an unsupervised non-hierarchical agglomerative K-means clustering technique, a particularly appropriate approach when working with heterogeneous datasets (MacQueen, 1967; Chytrý et al., 2002; De Cáceres & Wiser, 2012). This clustering technique is based on the random setting of clusters’ seeds and requires for the desired number of clusters to be previously established; in this case, set for partitions of 2–100 clusters with a one-unit pace increase (an exhaustive number given to necessarily include the optimal partition). For each partition, 10,000 iterations were conducted in order to obtain the best clusters configuration available for the páramo matrix. The resulting clusters, could be assimilated to a group of vegetation plots that shared affinities in terms of plant assemblages, dominant species and diagnostic ones (Chytrý et al., 2002). Then, the widely used, context-independent Ochiai indicator index (OI) was calculated to obtain the diagnostic value of each species within each cluster of each partition (Ochiai, 1957; De Cáceres, Font & Oliva, 2008). OI values were thresholded at 0.3 for species to be considered diagnostic of a specific cluster. Furthermore, cluster preference was checked, ensuring that the OI value considered diagnostic at least doubled the next highest OI value for the same species in a different cluster. Finally, particular attention was paid to rare species with low OI value but high exclusivity (Peyre & Font, 2011).

The optimal partition with the best cluster division was selected to represent the main páramo phytogeographical units. To do so, it is usual to rely on statistical criteria, such as the silhouette or pseudo-F, whose values depend on the number and scores of diagnostic species (De Cáceres & Wiser, 2012), even though there is no generally approved criterion to date (Tichý et al., 2010). However, the heterogeneity in our dataset implied that, with increasing partitions, certain units differentiated and separated faster than others, which prevented from using a strictly statistically based criterion that would compare similar speed divisions. Therefore, we built a subjective criterion that ensured sufficient variability, i.e. number of plots—set at 20, and floristic characterisation, i.e. number of diagnostic species—set at five, for every cluster. The optimal partition was hence defined as the most advanced clustering partition that met such criterion. Given the subjective character of the selection criterion employed, we provided five plots for each cluster to be used as seeds in case researchers are interested in reproducing our analyses with a supervised K-means clustering (Supplemental information 1). The optimal partition obtained therefore presented a certain number of natural clusters corresponding to an aggregation of vegetation plots with similar species content. However, some clusters might have not met these characteristics but instead corresponded to an aggregation of plots lacking the floristic similarity required to fit into a natural cluster but representing various vegetation types in too small numbers to create a new cluster. Such clusters should be qualified as artificial and are a common side effect of statistical classifications conducted on heterogeneous vegetation data, because each plot must fit into a cluster (Andrés & Font, 2011).

Bioclimatic data and geo-climatic characterisation

To correlate the natural cluster previously obtained with geo-climatic data, we relied on several climatic variables: annual mean temperature, temperature annual seasonality, annual mean precipitation, precipitation annual seasonality and cloud annual cover. The temperature and precipitation variables were obtained at a resolution of 30 arc-seconds (∼1 km) from the CHELSA project 1.2, Climatologies at High resolution for the Earth’s Land Surface Areas (http://chelsa-climate.org), which enhances bioclimatic data quality in tropical areas (Karger et al., 2017). The cloud variable was also obtained at a resolution of 30-arc-seconds from the EarthEnv Project 1 (http://www.earthenv.org; Wilson & Jetz, 2016). Using the R 3.4.3 software, all variables were cropped at the vegetation data extent and elevation (Supplemental Information 2). Each variable was then extracted for the latitude and longitude decimal coordinates corresponding to each UTM central point of the plot data. Finally, a data-frame containing for each plot its UTM central decimal coordinates, elevation and a value for each of the five bioclimatic variables considered was created.

Then, a geo-climatic principal component analysis (PCA) of the plot data was conducted, analyzed and visualised after applying the optimal clustering partition. Therefore, every cluster could be assimilated to a páramo phytogeographical unit, characterised by its list of dominant and diagnostic species as well as its climate and geography.

Hierarchical classification—similarity between phytogeographical units

One advantage of hierarchical clustering compared to non-hierarchical clustering is that it allows to understand the relationships within and between partitions thanks to the order of successive division (Tichý, Chytrý & Smarda, 2011). Since one of this study’s goals was to statistically quantify the relationships between the final clusters, the previously obtained K-means clustering results from partitions 2 to the optimal partition were transformed into a distance hierarchical classification. To do so, the Bray–Curtis distance was calculated, based on plot composition, for each pair of clusters of two successive K-means partitions, and then transformed into a similarity measure (S = 1 − D). All obtained values, corresponding to plot shifts and cluster forming with progressing divisions of the dataset, i.e., increasing partitions, were merged into a similarity matrix that quantified the relationships between the clusters of the optimal partition and all clusters from previous partitions. Then, a hierarchical divisive classification was built, using all similarity values superior to 0.1 (max. 1). Finally, this similarity tree up to the optimal partition was interpreted to (1) understand how easily the final clusters formed and segregated from the rest of the dataset, i.e., when no important plot exchange occurred and S values remained constantly high, and (2) evaluate their closeness to the other phytogeographical units.

PCA results

The first three components (PC) of the geo-climatic PCA represented a cumulative 76.9% of the variance. According to our results (Figs. 4A and 4B), PC1 (38.3%) was mostly a geographical axis positively marked by increasing latitude and longitude, as shown by the northern clusters Sub-1 and Mid-1, and negatively marked by increasing elevation with Sup-3 and Sup-4. PC2 (22.8%) was more of a humidity axis, positively marked by increasing precipitation and cloud cover, as shown by clusters Mid-3 and Mid-4, and negatively marked with increasing precipitation seasonality, as shown by clusters Mid-9 and Sup-4. Finally, PC3 (15.8%) was more of a temperature-related axis, positively marked by increasing annual temperature and temperature seasonality, as shown by clusters Mid-1 and Mid-6, and negatively marked by increasing elevation and decreasing latitude with clusters Sup-2 and Sup-3.

Description of the phytogeographical units

A description of the phytogeographical units based on their floristics, i.e. diagnostic species and dominant species (for further detail, see synoptic table in Supplemental Information 4), geo-climate, and complemented by bibliographical research on ecology is given below following an elevational and North–South order.

The sub-páramo

Two main sub-páramo clusters were identified.

Sub-1 (38 plots)—Guaramacal sub-páramo—The plots contained in this cluster were mostly distributed in the low-elevation páramos of Guaramacal (∼9°N–70°W) and punctually Táchira (∼8°N–72°W) in Venezuela. Among the diagnostic species were the small tree Libanothamnus griffinii (OI: 0.43) and shrub-tall Ruilopezia lopez-palacii (OI: 0.78) giant rosettes, bamboos including Chusquea angustifolia (OI: 0.78) and Chusquea steyermarkii (OI: 0.43), and shrubs such as Hypericum paramitanum (OI: 0.65). Dominant species present in most plots included the shrub Pernettya prostrata, tall grass Cortaderia hapalotricha and clubmoss Lycopodium clavatum. Considering such characteristics, the cluster was considered representative of the isolated and rich in endemics sub-páramo from Guaramacal dominated by mixed sub-páramo dwarf forests and shrublands with giant rosettes and bamboos (Cuello & Cleef, 2009).

Sub-2 (190 plots)—Widespread sub-páramo—This cluster was unresolved and considered artificial, because it included widespread vegetation plots that lacked floristic similarity. It was considered mostly a sub-páramo cluster due to its general low-elevation (<3,100 m), however some plots included here came from high-elevation. At this stage, no significant floristic coherence could be detected at species level, hence no valid list of diagnostic species could be provided. The dominant species, Pernettya prostrata and Hesperomeles obtusifolia were very common páramo shrubs present in respectively 22 and 20% of this cluster’s plots and also important in other clusters. Because of the cluster’s heterogeneity, its large amounts of plots and extensive species list (1,186), we conducted a complementary clustering analysis to intend revealing sub-clusters. To do so, we carried out a K-means clustering analysis on the cluster’s dataset at genus level, followed by the calculus of Ochiai Index values. In this case, the optimal partition could be identified using the silhouette statistical criterion, which showed a peak value of 0.08326 at the division into three clusters (Rousseeuw, 1987). Three main sub-clusters were identified: (1) sub-cluster 1—semi-dry grassland, (2) sub-cluster 2—shrubland and dwarf forests and (3) sub-cluster 3—secondary succession vegetation (Table 1). However, because the OI calculated at genus level were only considering this particular cluster, this analysis provided information on tendencies and not diagnostic taxa. Additional data would be required to segregate better this cluster.

The mid-páramo

Nine clusters of mid-páramo phytogeographical units spread over a large geographic gradient could be distinguished, ranging from giant rosette dominated communities in Venezuela to the North–East to mixed grasslands with giant rosettes in Colombia and tussock grasslands in Ecuador and Peru to the South–West. Altitudinal divisions between lower and upper mid-páramo were also perceived in some Ecuadorian and Colombian páramo areas.

Mid-1 (119 plots)—Périja-Santa Marta mid-páramo—Most plots from this cluster came from the Sierra de Périja and Sierra Nevada de Santa Marta close mountain ranges from northern Colombia. Few plots from the Colombian eastern cordillera and Venezuela that shared a similar semi-dry and seasonal climate were also included in this cluster. Its diagnostic species were the shrubs Pentacalia albotecta (OI: 0.41) and several Hypericum species, including Hypericum magdalenicum (OI: 0.41), Hypericum stenopetalum (OI: 0.36) and Hypericum baccharoides (OI: 0.33). Diagnostic herbs included Ranunculus spaniophyllus (OI: 0.38) and Lupinus carrikeri (OI: 0.32). Given that the tussock grass Calamagrostis effusa was dominant, this cluster mostly represented the mid-páramo mixed Calamagrostis grasslands with many locally diversified and endemic shrubs from this northern Colombian biogeographic complex (Rivera-Díaz & Fernández-Alonso, 2003; Pinto-Zárate & Rangel-Churio, 2010a).

Mid-2 (137 plots)—Eastern cordillera mid-páramo—This cluster was mainly distributed in the Colombian eastern cordillera. Among its diagnostic species were the dominant giant rosette Espeletia grandiflora (OI: 0.62), shrubs such as Diplostephium phylicoides (OI: 0.73) and Arcytophyllum nitidum (OI: 0.51), the grass Calamagrostis bogotensis (OI: 0.50), and herbs such as Bartsia santolinifolia (OI: 0.45) and Castratella piloselloides (OI: 0.54). The bamboo Chusquea tessellata and tussock grass Calamagrostis effusa were also dominant species, which suggested that both the mid-páramo from the drier eastern slope dominated by grasslands and the mid-páramo from the wetter western slope dominated by mixed-bamboo vegetation (Cleef, 1981) were represented in this cluster.

Mid-3 (164 plots)—Central and western cordilleras mid-páramo—This cluster’s plots were mostly distributed in the central and western cordilleras as well as the southern Andes in Colombia. The climatic conditions associated with this cluster informed of certain humidity, higher in the western cordillera and lower in the central cordillera. Diagnostic species included the shrubs Diplostephium schultzii (OI: 0.43), Monnina revoluta (OI: 0.32) and Baccharis macrantha (OI: 0.31) as well as the herbs Niphogeton ternata (OI: 0.40) and Bartsia orthocarpiflora (OI: 0.29). Among the dominant species were found the tussock grass Calamagrostis effusa, the giant rosette Espeletia hartwegiana, the shrub Pentacalia vaccinioides and fern Blechnum loxense. Considering such species composition, this cluster was therefore considered representative of the mixed grasslands with shrubs from the semi-humid and humid mid-páramo of the central and western Colombian cordilleras (Pinto-Zárate & Rangel-Churio, 2010b).

Mid-4 (161 plots)—Mixed group of humid mid-páramo—This cluster did not represent a fully coherent biogeographical unit and was considered artificial, because it included geographically widespread plots in Colombia and showed no significant diagnostic species, while its dominant species were common páramo plants. However, the geo-climatic PCA suggested that this cluster had a very strong humidity component, which was also sustained by the presence of common species such as the bamboo Chusquea tessellata and herbs like Arcytophyllum muticum and Carex bonplandii (Luteyn, 1999). An often repeated species combination within this cluster’s plots was the assemblage of the dominant tussock grass Calamagrostis effusa with the low shrub Pernettya prostrata, both widespread páramo species, and the prostrate plant Arcytophyllum muticum. As a result, this cluster undoubtedly represented humid mid-páramo from Colombia and showed certain floristic affinities with Mid-2 and Mid-3, but it was probably generated by grouping plots that lacked the diagnostic species of the other clusters but could not create a new valid one at this stage.

Mid-5 (55 plots)—Carchi mid-páramo—The plots included in this cluster came from the Ecuador-Colombia Andean border. Among the diagnostic species encountered were the local giant rosette Espeletia pycnophylla (OI: 0.90), shrubs such as Brachyotum lindenii (OI: 0.52) and Diplostephium rhododendroides (OI: 0.75), as well as the herbs Chaptalia cordata (OI: 0.40) and Lupinus pubescens (OI: 0.66). Consequently, and considering the dominance of the tussock grass Calamagrostis effusa, this cluster was revealed as the particular mid-páramo of mixed grasslands with the only Espeletia giant rosette species known to Ecuador (Moscol-Olivera & Cleef, 2009), and where transitionally occurs the southern-ending distribution for both Espeletia spp. and Calamagrostis effusa.

Mid-6 (139 plots)—Ecuadorian mid-páramo—This cluster’s plots were mostly located in Ecuador and to a lesser extent in Peru, and showed a certain warm and relatively seasonal climatic character. Diagnostic species for this cluster were mostly herbs, including Ranunculus peruvianus (OI: 0.31), Carex pygmaea (OI: 0.32), Galium corymbosum (OI: 0.40), Dorobaea pimpinellifolia (OI: 0.23) and Senecio chionogeton (OI: 0.35). Among the most dominant species was Calamagrostis intermedia, which confirmed this cluster’s representation of the common semi-dry grasslands mixed with shrubs of Pentacalia spp., Diplostephium spp. and Monnina spp. (Ramsay, 1992) of the Ecuadorian and Peruvian mid-páramo.

Mid-7 (87 plots)—Venezuelan mid-páramo and lower super-páramo—These plots were distributed in most Venezuelan páramos, especially in the Cordillera de Mérida (∼8°N–71°W), and covered a wide altitudinal range over the mid-páramo and lower super-páramo altitudinal belts. Among the most important diagnostic species were the giant rosette Espeletia schultzii (OI: 0.73), shrubs such as Baccharis prunifolia (OI: 0.44) and Oxylobus glanduliferus (OI: 0.44) and herbs like Azorella julianii (OI: 0.30) and Poa petrosa (OI: 0.48). In addition, the shrub Hypericum laricifolium and prostrate herb Acaena cylindrostachya were common. The lack of tussock grass dominance in this cluster contrasted to the other mid-páramo clusters. It therefore represented the dominant semi-dry páramos of Venezuela, where the diversified giant rosettes Espeletia spp., Coespeletia spp. and Ruilopezia spp., co-occur with shrubs such as Baccharis spp. and Chaetolepis spp. (Monasterio & Reyes, 1980; Diazgranados, 2012). However, the lower-super-páramo from Venezuela was also represented in this cluster as shown by the diagnostic herbs Hinterhubera imbricata (OI: 0.50) and Draba pulvinata (OI: 0.36) of high-elevation deserts (Berg, 1998).

Mid-8 (72 plots)—The Nevados upper mid-páramo—The plots from this cluster were mostly distributed in the Nevados páramo (∼4.8°N–75.3°W) and punctually in the Sumapaz páramo (∼4°N–74.2°W) in Colombia around 4,000 m elevation. This distinctive cluster was characterised by the diagnostic tussock grass Calamagrostis recta (OI: 0.69), shrubs such as Pentacalia vernicosa (OI: 0.44) and Baccharis rupicola (OI: 0.39), and herbs including Gentianella dasyantha (OI: 0.42) and Aa colombiana (OI: 0.29). Among the dominant species were the giant rosette Espeletia hartwegiana and herbs such as the common Oreomyrrhis andicola and Hypochaeris sessiliflora. With such floristic characteristics, this cluster was therefore distinguished from the common mixed grassland of Calamagrostis effusa and Espeletia spp. of the Colombian mid-páramo (Salamanca, Cleef & Rangel-Churio, 2003), and instead represented the ecotone between mid-páramo and super-páramo in these selected mountain ranges.

Mid-9 (93 plots)—The Ecuadorian upper mid-páramo—The plots contained in this cluster were distributed in Ecuador around 4,000 m in relatively seasonal páramos. Among the diagnostic species were the grasses Calamagrostis fibrovaginata (OI: 0.33) and Festuca andicola (OI: 0.47), as well as the herbs Gentianella cerastioides (OI: 0.55), Cerastium imbricatum (OI: 0.42) and the acaulescent rosette Valeriana rigida (OI: 0.35). Common species included other grasses, in particular Calamagrostis intermedia, and cushion forming plants, mostly from the genus Azorella, such as Azorella pedunculata and Azorella aretioides. This cluster therefore represented the upper mid-páramo transition from Ecuador dominated by mixed grasslands with cushion plants (Ramsay, 1992).

The super-páramo

Four clusters of Colombian and Ecuadorian super-páramos were revealed.

Sup-1 (156 plots)—Lower humid super-páramo—The plots included in this cluster were distributed in relatively humid environments around 4,200 m in Ecuador and southern Colombia. Among the diagnostic species were the shrub Diplostephium rupestre (OI: 0.50), grasses such as Festuca asplundii (OI: 0.46) and Calamagrostis ecuadoriensis (OI: 0.35), and herbs including Gentianella nummulariifolia (OI: 0.39) and Valeriana bracteata (OI: 0.36). Commonly found species included cushion plants such as Xenophyllum humile, Azorella aretioides and Plantago rigida. Therefore, this cluster representated of the transitional cushion plant communities with small shrubs from the lower-sub-páramo of semi-humid and humid mountains in Ecuador and southern Colombia. In contrast to the Mid-9 cluster, which contained mix grass-cushion communities dominated by grasses of the lower ecotone, Sup-1 represented the cushion-dominated vegetation of the upper ecotone where environmental humidity is constant, soils are deep and frost is limited (Sklenář, 2009). This zonal cluster resembled azonal bogs and mire vegetation due to shared dominant species, but its accompanying diagnostic species were key to differentiate them.

Sup-2 (81 plots)—The Nevados super-páramo—Most of the plots included here came from the Nevados páramo and secondarily from the Sumapaz páramo. The grass Bromus lanatus and the small herbs Hypochaeries sessiliflora and Valeriana pilosa were dominant species. Among the diagnostic species were the small shrubs Loricaria colombiana (OI: 0.20) and Pentacalia gelida (OI: 0.47), the grass Agrostis araucana (OI: 0.44) and several herbs including Lupinus alopecuroides (OI: 0.36), Erigeron chionophilus (OI: 0.57), Senecio isabelis (OI: 0.34) and Draba pennell-hazenii (OI: 0.27). This cluster therefore represented the super-páramo from the Nevados páramo, including the lower super-páramo, with Loricaria shrublands and blue meadows, and upper super-páramo, characterised by high-elevations deserts (Salamanca, Cleef & Rangel-Churio, 2003).

Sup-3 (94 plots)—Upper humid super-páramo—This cluster was mainly distributed around 4,400 m in central and northern Ecuador, and also in the Nevados páramo. Diagnostic species included the herbs Senecio nivalis (OI: 0.70), Erigeron ecuadoriensis (OI: 0.46) and Draba aretioides (OI: 0.43), as well as two indicators of certain humidity, the grass Calamagrostis ligulata (OI: 0.45) and the prostrate herb Ourisia muscosa (OI: 0.40). Other associated common species were Agrostis foliata, Xenophyllum humile, Cerastium floccosum and Luzula racemosa. As a result, this cluster represented the cold and semi-humid to humid upper-super-páramo found mostly in Ecuador. At this elevation, climatic conditions usually are very drastic with permanent night frost and high solifluction that confines the vegetation to few available microsites, making it almost desertic. In the case of humid upper super-páramo, the vegetation is organised in small patches with an overall ground-cover of 15–20%, which contrasts with dry upper super-páramos (Sklenář, 2000).

Sup-4 (61 plots)—Upper dry Ecuadorian super-páramo—The plots included here are restricted to the dry inter-Andean valley slopes of Ecuador, at elevations around 4,500 m. The grasses Calamagrostis mollis (OI: 0.49) and Festuca vaginalis (OI: 0.39) as well as small herbs such as Astragalus geminiflorus (OI: 0.72), Nototriche jamesonii (OI: 0.59) and Draba depressa (OI: 0.50) were diagnostic. Other common species were the basal rosette Hypochaeris sessiliflora and prostrate plant Baccharis caespitosa. This cluster represented the dry upper super-páramo, including rain-shadow deserts (e.g., the Chimborazo mountain) in Ecuador, which contrarily to its humid counterparts (Sup-3), shows 10–15% vegetation ground cover and includes highly specialised and endemic species (Sklenář, 2000).

Similarities between phytogeographical units

The results of the transformation of the K-means clustering into a hierarchical classification are shown in Fig. 5 and described below following the increasing order of division observed. Already at the division into four clusters, the dataset broadly divided latitudinally into four main units: one Ecuadorian mid-páramo unit (P4-b), one sub-páramo/Venezuelan mid-páramo unit (P4-c), one Colombian mid-páramo unit (P4-d), and one altitudinally discriminated unit, the super-páramo unit (P4-a). In further divisions, the Ecuadorian mid-páramo and super-páramo branches continued to interchange plots, as did the sub-páramo/Venezuelan mid-páramo and Colombian mid-páramo branches, but they remained mostly separated from one another.

First, early at the fifth clusters partition, the super-páramo group (P4-a) divided, first separating the transitional lower super-páramo from Ecuador (P5-a) from the rest, and then at seventh clusters partition, isolating the Nevados high-elevation páramo (P7-b) from the Ecuadorian high elevation super-páramos (P7-c). The latter group redivided at the 10th clusters partition by combining plots with the Ecuadorian mid-páramo cluster (P9-d) to create the humid upper super-páramo cluster (P10-d) and isolate the dry upper super-páramo cluster (P10-c). At the 15th clusters partition, the Nevados cluster P7-b divided by elevation into Mid-8, which was a transitional páramo/super-páramo ecotone, from its directly above super-páramo cluster Sup-2.

Second, the Ecuadorian mid-páramo group (P4-b) separated first at the 10th clusters partition based on elevation between lower mid-páramo (P10-e) and humid upper mid-páramo (P10-d). This latter cluster parted later at the 14th clusters partition, and divided by elevation into upper mid-páramo (P14-e) and upper dry super-páramo (P14-d).

Third, the Colombian mid-páramo group (P4-d) separated quickly at the sixth clusters partition into western (P6-f) and eastern clusters (P6-e). The eastern cluster divided progressively between the eighth and 13th clusters partitions, isolating clusters of the Colombian eastern cordillera (P13-j), the Santa Marta/Périja complex (P13-i) and the mixed humid Colombian grasslands (P13-k). The western cluster divided later at the 12th clusters partition into the Carchi mid-páramo cluster (P12-l) and the Colombian central and western cordilleras cluster (P12-k).

Lastly, the sub-páramo/Venezuelan mid-páramo group (P4-c) divided at the eighth clusters partition into a mixed cluster with the Colombian mid-páramo group, which eventually lead to the particular Guaramacal sub-páramo cluster (P9-f). Later, at the 12th clusters partition, the other group divided into the general sub-páramo cluster (P11-f) and Venezuelan páramo cluster (P11-g).