Herb stratum diversity and community structure in Gurez valley of Kashmir Himalaya: application of multivariate techniques in community analyses

ABSTRACT

In temperate forests, herbs comprise up to 90% of species richness yet received less attention, and their community structure, diversity patterns, and driving mechanisms still remain understudied. We provide a phytosociological overview of herb stratum to develop a vegetation model to characterize and designate vegetation groups and interpret their diversity and species distribution patterns. Multivariate analyses – classification (TWINSPAN) and ordination (indicator species analysis [ISA], non-metric multidimensional scaling [NMDS], and Canonical Correspondence Analysis [CCA]) were applied to data set collected from 32 plots across Gurez valley, Kashmir Himalaya to recognize the understory plant communities, determine the environmental predictors and highlight their significance. Altogether, 131 herbs from 38 families were recorded. TWINSPAN classified understory vegetation into three communities: Tanacetum multicule – Pedicularis pectinata – Aconitum heterophyllum (TPA), Taraxacum officinale – Trifolium repens – Plantago major (TTP), and Impatiens brachycentra – Tussilago farfara – Galium boreale (ITG). Diversity indices and species richness vary significantly, following the trend TTP>TPA>ITG. CCA revealed that disturbances, canopy, altitude, and moisture were the strongest parameters determining species differentiation. The study outlines a methodological workflow based on analytical methods and vegetation-plot data that describe vegetation groups and might be helpful in further in-depth vegetation classification syntheses and decision-making in conservation, global change issues, and management.

KEYWORDS:

• Temperate forests
• Herbs
• Altitude
• Conservation
• Disturbance
• Classification

1. Introduction

Mountain ecosystems support one-third of all terrestrial biodiversity characterized by topographic and habitat heterogeneity (Spehn et al., 2011). To ascertain the extent of species density and diversity for conservation interests, phytosociological explorations in mountain ecosystems are highly recommended (Gómez-Díaz et al., 2017; López-Angulo et al., 2018). Understory herbaceous stratum or ground herbs constitute majority of the cumulative floristic diversity in temperate forests, accounting for 90% of plant species (Gilliam, 2007; Spicer et al., 2020). Herb stratum, a taxonomically and ecologically heterogeneous group, governs various ecological processes in forests (Gilliam, 2007; Landuyt et al., 2019; Whigham, 2004). Consequently, several herbs are well-acknowledged indicators of environmental conditions, forest management, climate change, sustainability, health, and conservation status in forest ecosystems (Chávez & Macdonald, 2012; Lencinas et al., 2008). Exploring the mechanisms shaping the formation, development, and structuring of ecological communities has drawn much attention in ecological research. Such cognizance not only aids in comprehending the geographical heterogeneity in species distribution and plant community structure but also environmental requirements needed for effective ecological restoration and biodiversity conservation (Solefack et al., 2018). In the backdrop of changing environmental conditions and declining biodiversity, characterizing and monitoring herb-layer diversity remains challenging.

Experimental and theoretical investigations in mountain systems revealed that complex vegetation, climate, and soil characteristics (A. A. Dar & Parthasarathy, 2022a; Mestre et al., 2017; X. Zhang et al., 2021) drive species diversity patterns and community structure. Thus, it is evident that complex interrelationships exist among environmental covariates over species distribution; however, responses to environmental variables vary across plant groups (Kreft et al., 2010; Ricklefs & Renner, 2012). Topographic complexity alters solar radiation distribution, hydrology, and nutrient availability, resulting in significant fluctuations in climate, illumination hour, seasonal duration, geology, and anthropogenic disturbances; consequently, physiognomic differentiation of vegetation (Körner, 2007; Moeslund et al., 2013).

Over the past three decades, biodiversity loss triggered by anthropogenic activities has emerged as a major concern for field ecologists, notably following the Rio de Janeiro Convention on Biological Diversity summit on global biodiversity. Research in Himalaya has shown the impact of human activities on edaphic and environmental conditions by altering the vegetation diversity and distribution patterns (Wang et al., 2019). The impact of grazing on vegetation diversity is obscure since it depends on several temporal, geographical, and ecological constraints (Olsvig-Whittaker et al., 2006). Excessive degradation or the establishment of roads, powerlines, or other anthropogenic pathways may promote the invasion of non-native species (Le et al., 2018; Wedegärtner et al., 2022). Non-native taxa may change the structure and inherent functional aspects of the vegetative community without altering diversity, which might, in turn, affect ecological processes (Tilman et al., 1997). The disruptions caused by human activities, such as ecotourism, harvesting wild herbs and vegetables, and deforestation, are other key aspects determining the community structure (J. -T. Zhang et al., 2016). Considering the aforementioned threats and challenges, it is essential to investigate the interactions among forests, environment, and disturbance as well as comprehend the status of biodiversity in hitherto unexplored and species-diverse mountain areas.

The Himalaya represents one of the planet’s most intricate and diverse highland habitats, distinguished by extreme climatic variations and strong seasonality, harboring a wide range of ecosystems and species (Kala & Mathur, 2002; Oommen & Shanker, 2005). With a poorly protected geographical extent of over 600,000 km², the Himalaya consists of a network of mountains that span the borders of five countries – China, India, Bhutan, Nepal, and Pakistan (Wambulwa et al., 2021); nonetheless, recent data suggest that not more than 42% of the forest remains. Forests cover 7,13,789 km² of the land surface of India (21.71% of the total land area), of which more than 28% correspond to Himalayan landscape (Forest Survey of India (FSI), 2021). However, owing to remote access, uneven terrain, and unsettled geopolitical status, a significant area of the landscape awaits ecological exploration (Liang et al., 2020; Oommen & Shanker, 2005; S. M. Khan et al., 2013). The Kashmir Himalaya is no exception, exhibiting heterogeneous relief patterns. These differences have generated numerous habitats with distinct climates, plant groups, and species richness. Research characterizing species richness in these mountains covers specific areas, leaving several areas unstudied. Gurez valley, a part of Kashmir Himalayan landscape, is one such example. Prior investigations of the montane vegetation in Gurez valley have tended to be descriptive with limited quantitative approaches (Dad & Khan, 2011a, 2011b; Dad & Reshi, 2015), especially computer-based vegetation data analysis. No attempts have been undertaken to quantify the plant communities along with environmental predictors to determine the primary drivers of vegetation patterns. Such quantitative efforts may provide new perspectives on the processes governing species distribution, which are crucial for management as well as conservation initiatives in the context of global environmental change.

Application of computer-aided methods to analyze vegetation, especially in connection to species heterogeneity across environmental gradients, is a well-established field in ecology (S. M. Khan et al., 2013). Multivariate analyses are robust procedures intended to investigate multivariate interacting systems. Such analyses have been extensively implemented in ecological research and can efficiently characterize plant communities (Gauch, 2010). The classification (cluster analysis) and ordination tools (indicator species analysis [ISA], non-metric multidimensional scaling [NMDS], and Canonical Correspondence Analysis [CCA]) are used to recognize plant communities by creating a low-dimensional space in which related sample plots/species aggregate together, determine the environmental predictors and highlight the significance of interactions between vegetation and environmental matrices (Dufrêne & Legendre, 1997; Ter Braak, 1987). In this study, we applied multivariate analysis to an observational data set collected from 32 plots across the Gurez valley. The data set contains temporal measurements of understorey herb diversity and abundance with precise measurements and records of climatic, topographic, disturbance, and canopy cover variables. Combined with multivariate analysis, the data set renders a good opportunity to examine relative contribution of predictor variables to the observed species distribution patterns. The present study intends to develop an empirical vegetation model to characterize and designate the vegetation groups and present them in an ecological and vegetation paradigm to recognize indicators for specific micro-environmental conditions. Specifically, we pose and attempt to address: (Q1) Do sampled data exhibit clustering, i.e., presence of different communities; (Q2) Do species composition and diversity and density of the herb layer differ among communities? and (Q3) What are the principal environmental predictors of understory herb stratum, and are these drivers equally important? To address the first and third questions, we used multivariate analytical techniques to tease apart the communities and check for the drivers of species distribution patterns. We hypothesized – First (H1), species composition, density, and diversity vary among the plant communities due to micro-environmental variation; Second (H2), we expected that understory vegetation is influenced by topography, climate, disturbance, and canopy cover to varying extent. Altitude is likely the most significant determinant because it affects the microenvironment that supports the survival of understory vegetation (Rahman et al., 2020).

2. Materials and methods

2.1. Study area

The study was conducted in Gurez valley, located on the banks of river Neelam extending between 34°30′–74°41′N and 74°37′–74°46′E in Kashmir Himalaya, India (Figure 1). The area shares an international boundary with Pakistan and is, therefore, largely remained unexplored due to security concerns. The valley lies at an average altitude of about 2370 m, experiencing a mean annual temperature and precipitation of 4.6°C and 742 mm, respectively (https://www.worldclim.org/). The climate is temperate, characterized by brief, cool summer, and prolonged, harsh winter with frequent snow and frost. Lofty snow-covered mountains bound the valley with the highest peak at Razdan-top. The understory begins to emerge when the snow melts in the spring and summer, blooming fully from June to September and drying up by the end of October. The soil is shallow, and typically, only bed rock is visible at the summits of steep, bare cliffs and protruding ridges; however, at lower altitudes, the soil is deeper and has a medium to fine texture (Dad & Reshi, 2015). The study area is an important ecological landscape characterized by dense coniferous forests dominated by Abies pindrow and Pinus wallichiana with intermittent patches of Betula utilis at higher altitudes followed by alpine vegetation. Microclimatic variations due to complex topography support highly diverse understory vegetation (Dad & Reshi, 2015). Grazing by livestock is common during the short summer.

Figure 1. Study area location and landscape overview of Gurez valley, Kashmir Himalaya, India.

2.2. Data collection and analyses

A floristic survey was carried out in summer of 2020 to investigate the herb-layer species diversity of hitherto unstudied coniferous forest of Gurez valley in Kashmir Himalaya. Phyto-sociological approach (F. Haq et al., 2017; Laurance et al., 2004) was adopted to record quantitative features of understory herb stratum in five 1 m × 1 m sub-quadrats (four at corners and one at the centre) nested within 0.25 ha ([1 m²] × [5] × [32]) plots established randomly along topographic, climatic, canopy and disturbance gradients (A. A. Dar & Parthasarathy, 2022a; Cochran, 1977) with a minimum distance of 500 m between plots. Attempts were made to replicate plots evenly across the landscape, but this was not always possible because of logistic constraints and steep terrain. Species composition, richness (total number of species), and abundance encountered per plot were recorded. Botanical specimens were collected, poisoned, pressed, dried, and layered between newspapers. Specimens were provisionally identified by means of available web resources – Flora of Pakistan (http://legacy.tropicos.org/Project/Pakistan), Flowers of India (http://www.flowersofindia.net/), Indian Biodiversity Portal (https://indiabiodiversity.org/) and Plants of the World Online (https://powo.science.kew.org/) and the book – Flowers of Himalaya (Polunin & Stainton, 1984) before being validated at Centre for Biodiversity and Taxonomy, University of Kashmir. Nomenclature of families is based on the Angiosperm Phylogeny Group classification (Angiosperm Phylogeny Group et al., 2016), and species follow The World Flora Online (http://www.worldfloraonline.org/) to attain global homogeneity and to remain in touch with revised taxonomic nomenclature. Coordinates (latitude and longitude) and altitude of the sampling plots were recorded using a JUNO 3E (accuracy 2–5 m) handheld GPS, and other corresponding topographic variables – slope (°) and aspect (°), are determined through modelling 30 m digital elevation model (DEM) using ArcGIS (10.3; ESRI, Redlands, CA, U.S.A). Climate data (temperature and precipitation) of respective vegetation plots retrieved from WorldClim global database (https://www.worldclim.org/) using “raster” library (v. 3.5–11; Hijmans, 2021) are downscaled to 1-arc second (30 m resolution) to maintain spatial homogeneity with DEM data. Canopy (%) is assessed by observing the tree canopy and assigned to one of five classes (1–5, open to closed). To account for disturbance, grazing was determined by observing recent signs and assigning “Yes” or “No”, number of trails crossing through the plot, and number of invasive alien species per plot (A. A. Dar & Parthasarathy, 2022a; Ebrahimi et al., 2014; Wangchuk et al., 2014). IUCN criteria (https://www.iucnredlist.org/) is applied to identify the Red List status (critically endangered, endangered, rare, and vulnerable) of the species. Endemics (restricted to Himalayan landscape) are categorized based on literature (Chawla et al., 2008) and online web sources (https://powo.science.kew.org/). Density, and frequency of each species are computed and changed to their relative values to evaluate modified importance value index (IVI; Curtis & McIntosh, 1950; Rasingam & Parthasarathy, 2009) in “BiodiversityR” library (v. 2.14–4; Kindt, 2023) in R software (R Core Team, 2022). To check sampling completeness and to visualize the richness differences among plant associations, the sampling coverages are computed and plotted using 50 bootstrap replicates with 95% confidence interval in “iNEXT” library (v. 2.0.20; Hsieh et al., 2022). The diversity function of “vegan” library (v. 2.6–4; Oksanen et al., 2022) is used to calculate Shannon diversity index and Pielou evenness following Magurran (2004). For β-diversity analysis, we used the spatial turnover and nestedness-resultant components, applying “Jaccard” as family of dissimilarity index using “adespatial” (v. 0.3–21; Dray, Bauman, et al., 2023) and “ade4” (v. 1.7–22; Dray, Dufour, et al., 2023) libraries. The species are divided into five frequency classes – Rare (1 − 20%), Occasional (21 − 40%), Frequent (41 − 60%), Abundant (61 − 80%), and very abundant (81 − 100%), following Tansley and Chipp (1926) and Tansley (1946). Kruskal-Wallis test/one-way ANOVA is used to test the statistical differences (p < 0.05) of different variables among plant communities, followed by post hoc pairwise comparisons given that the χ²/F is significant (p ≤0.05). Descriptive statistical tools are used to calculate the mean and standard error (SE) of the data.

2.3. Multivariate analyses

IVI of all the recorded plant species, along with climatic, disturbance, and topographic data, was used for multivariate analytical techniques (Rahman et al., 2020). The Hopkins statistics (Lawson & Jurs, 1990) analyze the likelihood that a particular dataset was produced by a uniform data distribution, which enables it to determine the spatial randomness (clustering tendency) of the dataset. The null hypothesis (H0) asserts that the data follow a continuous distribution (implying no meaningful clusters). H > 0.5 allows rejecting the null hypothesis (Hopkins & Skellam, 1954). We performed a hierarchical unsupervised divisible classification – two-way indicator species analysis (TWINSPAN; Hill, 1979) and a non-parametric multiscale ordination technique – NMDS (Clarke, 1993) using the Bray – Curtis distance of the whole data to identify plant communities using “twinspan” (v. 0.9–0; Oksanen & Hill, 2019) and “vegan” libraries. TWINSPAN classified sample plots by reciprocal averaging (Hill, 1979), and a dendrogram was constructed based on eigenvalues. NMDS aims to depict the relative position of plant associations in multidimensional space using minimal recognizable dimensions (Thakur et al., 2022). Non-parametric Permutational Multivariate Analysis of Variance (PERMANOVA) post hoc test utilizing Bray – Curtis dissimilarity is carried out with ”adonis2()” function of ‘vegan’ library to analyze the significant differences (after 999 permutations) among plant associations. Once the number of clusters is determined, the affinity of species (specificity and fidelity) in each habitat type is evaluated by indicator species analysis (ISA; sensu Dufrêne & Legendre, 1997) using ‘indicspecies’ library (v. 1.7.12; De Caceres & Legendre, 2009). Indicator values of species within each determined association are obtained and tested for statistical significance using a Monte Carlo test (n = 999 random permutations). The cut-off for identifying indicator species is set at 30% with a 95% confidence level (Dufrêne & Legendre, 1997; Ter Braak & Prentice, 1988), and the selected indicator species are then utilized to name the communities. To link species distribution patterns (species matrix) to changes in environmental conditions (environmental matrix) using ”cca()” (canonical correspondence analysis) function in “vegan” library, we measured 11 variables related to topography, climate, canopy, and disturbance characterizing the environment of each plot. Variance inflation factor (VIF) of variables is used to measure how much other explanatory variables inflate the regression coefficient; therefore, explanatory variables with VIF > 5 are removed through the backward selection. Finally, we used the Monte Carlo test (n = 999 permutations) to assess order (ranking), significance (p < 0.05) of environmental variables on the response variables (species and plots) and to test the significance of each environmental variable in the canonical axes (Hill & Gauch, 1980; Hill, 1979; Ter Braak, 1986).

3. Results

3.1. Floristic composition

We recorded 131 understory herbs representing 38 families in 32 sampling plots with a mean species richness of 13.97 (±1.11) species m⁻² across the study area (Table 1 and S1). Asteraceae (17 genera, 23 species), Fabaceae (7, 12), Apiaceae (7, 10), Lamiaceae (7, 9), and Rosaceae (5, 8) exhibit highest diversity in both genera and species, respectively (Figure 2). Angiosperms represent 37 families and 130 (90.1%) species, and pteridophytes by one family harboring just single species. Eight families are represented by two species, whereas four families are by three, and 15 are monotypic, represented by single species. Viola (4 spp.), Artemisia (3), Cynoglossum (3), Gentiana (3), and Saussurea (3) constitute the most diversified genera. The categorization of species into major plant groups showed that the majority of them are allied to the Magnoliopsida (dicotyledonous; 116 spp., 88.55%), followed by Liliopsida (monocots; 16 spp., 10.69%), and Polypodiopsida (ferns; one sp.). List of the species, along with letter codes are presented in Table S1.

Figure 2. Differentiation of species by family across three understory herb associations in Gurez valley of Kashmir Himalayan region.

Table 1. Descriptive statistics (mean ± SE) of community attributes across three understory herb communities in Gurez valley of Kashmir Himalayan region.

Community parameter	ITG community		TTP community		TPA community		Total
Sampled plots (1m^2)	50 (10×5)		70 (14×5)		40 (8×5)		160 (32×5)
Species richness	55		87		57		131
Species (m^−2)	12.10±1.63		16.43±1.76		12.00±2.22		13.97±1.11
Density (m^−2)	37.58±5.65		43.26±3.36		23.00±5.01		34.58±2.79
Abundance	1492		3028		1012		5532
Geographical parameter
Altitude (m)	3066.10±55.11		2758.71±66.45		3208.00±190.50		2967.09±65.68
Slope (°)	25.00±2.25		26.17±2.83		25.02±2.48		25.52±1.51
Aspect (°)	29 (18.02–351.87)		13 (28.30–350.09)		21 (3.81–359.59)		34 (3.81–359.59)
Species diversity
Shannon diversity (Hˈ)	1.53±0.16		2.07±0.14		1.83±0.15		1.84±0.09
Pielou evenness (J)	0.63±0.05		0.75±0.03		0.77±0.04		0.72±0.03
Beta (β)-diversity	0.36		0.35		0.42		0.40
Disturbance
Trails	0.40±0.16		0.79±0.24		0.38±0.18		0.56±0.13
Grazed plots	0.30±0.15		0.71±0.13		0.75±0.16		0.59±0.09
Invasive species	0.40±0.22		3.64±0.63		1.25±0.53		2.03±0.40

3.2. Community associations

Hopkins’ statistics show that the null hypothesis is rejected (H > 0.5), indicating that dataset very likely contains clusters. TWINSPAN effectively discriminated the plots into three interpretable plant communities, which can be clearly seen in three main branches of the dendrogram (Figure 3). The variation in forest communities is evident and related to intricate ecological gradients: (i) a relatively lower altitude (2758.71 m) association within mixed evergreen forest; (ii) a mid-elevation (3066.10 m) community dominated by moist temperate and sub-alpine elements and (iii) a high altitude (3208 m) plant community segregating the alpine vegetation (Table 1).

Figure 3. Classification (TWINSPAN) of understory herb vegetation presenting three communities structured in Gurez valley of Kashmir Himalayan region. The numerical numbers within square at each node represent the division number. N is the number of plots. Communities are named based on the indicator species.

3.2.1. Impatiens brachycentra–Tussilago farfara–Galium boreale (ITG)

The ITG community is comprised of 37.58 (±5.65) individuals per plot representing 55 (60%) plant species with a mean of 12.10 (±1.63) species m⁻² (Table 1). This plant community is recognized in 10 sampling plots occupying north-to-northeast (steep) slopes at 3066.10 (±55.11) m. The leading indicators are Impatiens brachycentra (stat = 0.87; p = 0.001), Tussilago farfara (stat= 0.62; p = 0.025), and Galium boreale (stat = 0.61; p = 0.025), thus referred as ITG community. Besides, other associated significant indicators include Epitalem elatum and Trillium govanianum. Fragaria nubicola is the most abundant, while Galium boreale, Geranium wallichianum, Impatiens brachycentra, and Tussilago farfara are abundant. Moreover, there are two frequent, 19 occasional, and 29 rare species. Despite being least affected by anthropogenic factors (Table 1), ITG harbours least number of indicator species (25). The community supports three critically endangered (Saussurea costus), two endangered (Trillium govanianum and Aconitum heterophyllum), two vulnerable (Cypripedium cordigerum and Polygonatum verticillatum), and two endemic species (Figure 4).

Figure 4. Distribution of endemic and IUCN threatened species across three understory herb communities in Gurez valley of Kashmir Himalayan region. The number on the scale represents the species abundance. Letter codes represent the species as: Aco_het: Aconitum heterophyllum; Del_ves: Delphinium vestitum; Geu_roy: Geum roylei; Sau_cos: Saussurea costus; Bup_lon: Bupleurum longicaule; Dip_ine: Dipsacus inermis; Gen_car: Gentiana carinata; Lom_cae: Lomatogonium caeruleum; Oxy_cae: Oxytropis cachemiriana; Ped_pec: Pedicularis pectinata; Rhe_web: Rheum webbianum; Sel_veg: Selinum vaginatum; Lac_dol: Lactuca dolichophylla; Pri_mac: Primula macrophylla; Sal_mor: Salvia moorcroftiana; and Sca_spe: Scabiosa speciosa.

3.2.2. Taraxacum officinale–Trifolium repens–Plantago major (TTP)

Collectively, 14 sample plots hosted this plant association exhibiting diverse and heterogeneous floristic composition (87 [66.41%]) with the highest mean species count (16.43 ± 1.76 m⁻²) distributed at relatively lower altitude (2758.71 ± 66.45 m). This community was recorded mostly on northern aspect except for two plots (TTP_16 and TTP_19) which lie on southern side. TTP comprises the highest number of significant indicators (37 spp.) with Taraxacum officinale (stat = 0.82; p = 0.002), Trifolium repens (0.82; 0.001), and Plantago major (0.76; 0.004) being the topmost species. Pimpinella diversifolia, and Ranunculus laetus are other associated significant indicator herbs. Maximum density (43.26 ± 3.36 ind. m⁻²) with intensive grazing has led to an increase in number of invasive aliens. Among the recorded species, most abundant, abundant, frequent, occasional, and rare species account for 3.45%, 1.15%, 5.756%, 25.29%, and 64.37%, respectively. One critically endangered (Aconitum chasmanthum), one endangered (Aconitum heterophyllum), one vulnerable (Silene vulgaris), and four endemic species are supported by this community (Figure 4).

3.2.3. Tanacetum multicule–Pedicularis pectinata– Aconitum heterophyllum (TPA)

TPA is a high-altitude (3208 ± 190.50 m) plant community recognized in eight sampling plots on steep slopes above the tree line in alpine Himalaya, contributing 57 (43.51%) plant species (Table 1). Tanacetum multicule, Pedicularis pectinata, and Aconitum heterophyllum characterize the top three indicators with an indicator value of 0.87 (p = 0.00), 0.63 (p = 0.00), and 0.62 (p = 0.00), respectively. Distributed on the northeastern slopes with lowest density (23.00 ± 5.01 ind. m⁻²) and species richness (12.00 ± 2.22 spp. m⁻²), this association prefers low temperature and precipitation. The community suffers intense grazing pressure (0.75 ± 0.16) during the brief growing season in summer. The abundant species include Tanacetum multicaule and Koeleria macrantha, whereas 15 species occurred occasionally, and 36 are rare species. The community exhibits two critically endangered (Aconitum chasmanthum and Fritillaria roylei), one endangered (Aconitum heterophyllum), one vulnerable (Rheum webbianum), and 10 endemic species.

3.3. Species diversity and significance testing

NMDS generated a two-dimensional solution (stress = 0.19) with clearly distinguished three plant associations, as indicated by the distinct separation of the communities along the first two axes (Figure 5). Plots in close proximity indicate that they share similar species composition. The species distribution in ordination space demonstrates a consistent variation in species composition across communities. Species congregated close to the center of the biplot indicate their lack of affinity for a particular environment and are omnipresent. The post hoc PERMANOVA characterizes a strong compositional difference (F = 10.77, p < 0.001) among communities. TTP community has a significantly higher H′ (2.07 ± 0.14), followed by TPA and ITG communities (χ² = 6.61; p = 0.03; Figure 6; Table 1). The J ranged between 0.63 and 0.77 and followed the trend TPA > TTP > ITG. The three communities have nearly identical sample coverage (ITG = 99.60 [±0.005]%; TPA = 99.61 [±0.008]% and TTP = 99.70 [±0.002]%), implying they are equally complete (Figure 7). The extrapolation across the communities increases the sample coverage by mere 0.002%, 0.004%, and 0.002% in ITG, TPA, and TTP, respectively. The estimated ITG and TPA richness is almost equal; however, TTP is 1.53 and 1.58 times richer than ITG and TPA communities (Figure 7). Across the three communities, the beta diversity ranges between 0.35 and 0.42, with highest turnover and nestedness in TPA (0.26) and TTP (0.08). The maximum LCBD corresponds to TPA, followed by ITG and TTP communities (F = 3.47; p = 0.04), whereas SCBD follows ITG > TPA > TTP (Figure 6e,f). Pairwise comparison results show that altitude of TTP is significantly lower (p = 0.003) than ITG, whereas densities in ITG (p = 0.04) and TPA (p = 0.01) are significantly lower than TTP. Likewise, the number of invasive aliens is statistically higher in TTP compared to ITG (p = 0.000) and TPA (p = 0.014; Figure 6).

Figure 5. Two-dimensional ordination plot of the distance matrices from a NMDS (non-metric multidimensional scaling; stress = 0.19) segregating 32 sampling sites representing 131 plant species into three plant communities (ITG: Impatiens brachycentra–Tussilago farfara–Galium boreale; TPA: Taraxacum officinale–Trifolium repens–Plantago major; TTP: Tanacetum multicule–Pedicularis pectinata–Aconitum heterophylum) in Gurez valley of Kashmir Himalayan region. Refer to the supplementary (Table S1) for scientific names of the species letter codes.

Figure 6. Statistical variations in community attributes in Gurez valley of Kashmir Himalayan region.

Figure 7. Coverage-based rarefaction and extrapolation curves with 95% confidence intervals comparing plant species richness across three understory herb communities (ITG, TPA, and TTP) in Gurez valley of Kashmir Himalayan region. Note solid lines represent the reference samples, and the extrapolation part is almost invisible.

3.4. Species environment relationships

Results of CCA are summarized in Table 2 and Figure 8. Temperature is removed from the list of potential predictors owing to its co-linearity (VIF >10), enabling other variables to account for 41.92% variation (eigenvalues: 1.20). Together, the first two canonical axes account for 56% inertia representing an eigenvalue of 0.28 for each (CCA1: F = 3.5; p = 0.001 and CCA2: F = 2.67; p = 0.034). Overall, the CCA model is significant (F = 1.76; p = 0.001), suggesting further testing (Monte Carlo permutation) of the constrained axes and variables. Significance testing depicts that key factors influencing species distribution are canopy, altitude, number of trails, and invasive aliens (p < 0.05), followed by (in order) grazing, aspect, slope, and precipitation (though insignificant). CCA1 exhibits strong positive and negative correlations with disturbance and precipitation, and altitude, suggesting a complex disturbance-topographic-precipitation gradient from lower altitude with high precipitation and grazing to high mountains with low precipitation. The CCA2 is positively and negatively correlated with altitude and canopy, reflecting altitude-canopy gradient, separating canopy-dense plots at lower altitudes from open plots at higher altitudes. Biplot demonstrates that key diagnostic species of ITG community strongly correlate with precipitation-disturbance gradient. The indicator species of TPA are strongly associated with open slopes along CCA2 at higher altitudes, whereas TTP community prefers precipitation-altitude gradient.

Figure 8. Species distributions on topographic, climatic (PPT: precipitation; MI: moisture index), and disturbance (IA: invasive aliens) gradients across the three understory plant communities in Gurez valley of Kashmir Himalayan region. Green circles represent the sampled plots. The six-letter code represents the name of the plant species, with first three letters corresponding to the genus name and the subsequent letters being initials of species name. Refer to the supplementary (Table S1) for scientific names of the species letter codes.

Table 2. Summary of canonical correspondence analysis (CCA).

Variable	VIF	Significance				Axis	Eigenvalue (% explained inertia)	Significance
		χ^2	F	p				χ^2	F	p
Altitude	5.48	0.18	2.33	0.002		CCA1	0.28 (24)	0.28	3.5	0.001
Slope	1.89	0.07	0.87	0.66		CCA2	0.28 (23)	0.21	2.67	0.034
Canopy	1.59	0.18	2.39	0.001		CCA1+CCA2	0.56 (47)
Aspect	1.54	0.09	1.17	0.24
Trails	1.54	0.12	1.61	0.02
Grazing	1.95	0.1	1.31	0.09
Moisture index	1.94	0.12	1.62	0.02			Significance of CCA model
Precipitation	3.71	0.09	1.18	0.183			χ^2	F	p
Invasive aliens	3.26	0.13	1.66	0.016			1.20	1.76	0.001

4. Discussion

4.1. Species diversity

Despite being described as structurally diverse and species-poor ecosystems with oligarchic patterns, temperate forests may significantly enhance the diversity of distinct plant groups (e.g., herbaceous flora) or even promote regional diversity (Draper et al., 2018). The understory vegetation has major ecological significance because it affects tree regeneration and inhibits the establishment of characteristic seedlings (Thrippleton et al., 2016). The understory flora of Gurez valley is fairly heterogeneous across the complex mountain landscape (Table 1). Regardless of the area sampled, the number of species documented in our study was greater than many prior Himalayan studies (Altaf et al., 2021; Tiwari et al., 2020; Ullah et al., 2022) and studies undertaken elsewhere (Ebrahimi et al., 2014; Wang et al., 2019). Compared to ITG and TPA communities, the abundance and richness in TTP community are greater, presumably due to diverse microhabitats, greater resource availability, and disturbance mechanisms that led to the localization of non-native species too supporting the H1. The enhancement of species richness by intermediate disturbance levels, as proposed by Connell (1978), is well supported by our results. Since the survival of several understory taxa depends on trees, “sheltering effect” might be another potential mechanism for greater herbaceous diversity (Sánchez-Jardón et al., 2010; Treydte et al., 2007). Decreasing trend in species richness from TTP to TPA along the altitudinal gradient may have resulted from eco-physiological limitations due to persistent snow cover and freezing temperature, thereby inhibiting plant growth and ecosystem productivity and, consequently, reducing diversity (Gómez-Díaz et al., 2017; Körner et al., 2011). The alpine flora acquired various adaptations, viz., enhanced resource investment in reproductive and underground components (e.g., rhizomes, bulbs, etc.), extended pollinator-attractive parts, and improved dispersal methods to endure severe conditions (Fabbro & Körner, 2004). The Shannon diversity of the observed plant associations is consistent with prior vegetation studies in the Himalayan landscape (Altaf et al., 2021; J. A. Dar & Sundarapandian, 2015; Kumar et al., 2020). However, compared to western Himalaya (Joshi et al., 2022; Rana et al., 2022), coniferous forests in Bhutan (Wangchuk et al., 2014), and Bogd Khan Mountain, Magnolia (Bazarragchaa et al., 2022), our study area is less diverse. This could be attributed to plot type and size variations, area inventoried, and vegetation type since distinct habitats encourage heterogeneous species richness (A. A. Dar & Parthasarathy, 2022b). Lesser evenness observed in ITG community indicates the uneven distribution of abundance among species resulting in the dominance of certain species compared to TTP and TPA communities. Majority of the species exhibit patchy distribution specialized for a particular community; over half the species occurred in either one of the communities. For example, the least overlap between TPA and ITG communities suggests that many alpine species are not prevalent and/or susceptible to competition, characteristic environmental conditions, and dispersal restrictions because of their smaller fruit size and short stature (Leal et al., 2014; Lima et al., 2015). The abundance of these characteristic species might be used as effective markers for gauging the efficacy of management strategies (D’Amato et al., 2009). However, several species (ubiquitous) possess characteristics such as effective seed dispersal, efficient water and mineral uptake, vegetative reproduction, etc., that enhance their competitiveness and capacity for colonization (Lima et al., 2015). Asteraceae, Fabaceae, and Apiaceae were in the order of dominance. Asteraceae exhibits a wide ecological range and has a firm commitment to finite resources, enabling it to dominate in a wide environmental setting. The observed family dominance resembles other vegetation studies conducted across Himalayan landscape (Joshi et al., 2022; S. M. Haq et al., 2021; Tiwari et al., 2020) and elsewhere in the world (Tserendulam et al., 2018; Wang et al., 2019).

4.2. Multivariate techniques in community analysis

With the development of new analytical tools and approaches, several attempts in recent vegetation ecology have uncovered the intricate pattern of vegetative structure (Mayor et al., 2017). Such knowledge is essential for devising future conservation and plan management strategies. Forest policy and management efficiency may increase by using conservation techniques based on spatial patterns of species richness and diversity, including the recognition of key floristic groups and target locations for conservation (W. Khan et al., 2016). The multivariate analyses carried out as part of this study established three distinct plant communities in the Gurez valley of Kashmir Himalaya. The number of plant associations varies depending on the area studied, sampling method, and analytical techniques used. For instance, Bano et al. (2018) identified four vegetation communities across Beer Hills, Pakistan. Likewise, the number of associations can be three (Rahman et al., 2020), four (Ullah et al., 2022), five (W. Khan et al., 2016), and can go all along up to 24 (A. M. Khan et al., 2019; Ahmed et al., 2006; F. Haq et al., 2017). The NMDS strongly validated the cluster analysis as communities identified can easily be projected on the two-dimensional NMDS ordination space.

In recent years, the application of indicator species analysis for monitoring environmental conditions or identifying habitat or community types has become a well-established approach in vegetative ecology for both theoretical and practical reasons (W. Khan et al., 2016). One or more indicator species typically exist for each kind of environment and plant group (Vilches et al., 2013), which is clearly elucidated in our findings. These indicators reflect a particular microclimatic state or environmental change. Indicator species have been utilized as a powerful management tool, with the flexibility to characterize any ecoregion for the purpose of analyzing biodiversity and overall ecosystem health (Legendre, 2013; S. M. Khan et al., 2011). Employing the indicator species framework, our results highlight realistic targets for designing a long-term and pragmatic biodiversity conservation plan for the Gurez valley.

4.3. Environmental drivers facilitating species distribution

The role of ecological gradients exhibiting heterogeneity over space and time in structuring plant communities is well acknowledged (Bueno et al., 2018; J. -T. Zhang et al., 2015; Radhamoni et al., 2023). This was proven true in Gurez valley (Table 2 and Figure 8). Each community displayed a distribution area corresponding to a particular set of environmental factors. The findings are congruent with the environmental heterogeneity hypothesis (Huston, 1994), which asserts that subtle environmental variations govern species diversity. This could be attributed to the synchronous existence of several plant groups, such as sciophytes, together with resource-demanding heliophytes and generalist species (Chávez & Macdonald, 2012). Our findings emphasize the significant contributions of the aforementioned drivers to species distribution, particularly in relation to the distribution of the indicators.

4.3.1. Role of altitude

CCA ordination of herb stratum revealed a vegetation continuum that emerges to be a function of disturbance and altitude since both are strongly linked with the first axis (CCA1). Consequently, it can be claimed that the disturbance-altitude gradient is the main factor regulating understorey vegetation, which is consistent with the findings of prior research on mountain vegetation (e.g., Bisht et al., 2022; Gómez-Díaz et al., 2017). Topographic limitations are principal contributors to regional habitat heterogeneity regulating the vegetation distribution patterns in mountain landscapes (Titshall et al., 2000). In the present study, the relative contribution of topography, in particular the altitude, is significant (Table 2), probably because altitude, together with slope and aspect, shapes the solar azimuth, hydrology, climate, and edaphic characteristics, indirectly affecting the vegetation by creating microenvironments (A. A. Dar & Parthasarathy, 2023; S. M. Khan et al., 2013). This is additionally backed up by the fact that the majority of the vegetation in our area is distributed between 2300 and 3800 m. These findings are analogous to vegetation studies carried out in the western (Ilyas et al., 2018; S. M. Khan et al., 2013), eastern (Tambe & Rawat, 2010), and northwestern Himalayas (Dvorský et al., 2011; Ullah et al., 2022), which likewise noted the significant influence of altitude.

4.3.2. Role of climate

Altitudinal increase in mountain systems is often accompanied by reduced temperatures and increased humidity. This climate gradient which serves as a stressor for plants governs the regional vegetation trends by encouraging species that can survive at various altitudes (Klimes, 2003). The ideal limits of resilience to temperature, precipitation, and moisture levels vary across understory herbs. Temperature limitations cause a reduction in species richness at higher altitudes (Brinkmann et al., 2009). Alpine species are susceptible to temperature increases since they already occur at the mountain peaks and are unable to ascend uphill. The higher species richness of specialized plant groups at different altitudes reflects an evolutionary inheritance or greater speciation or reduced extinction under certain conditions, as reported for the family richness of herbaceous angiosperms across the latitudinal gradient (Hawkins et al., 2011; Jablonski et al., 2006). Even though precipitation generally decreases with altitude, humidity tends to rise because lower temperatures induce less evaporation, and clouds and mist generate additional water input that is not taken into account by conventional precipitation estimates (Gehrig-Downie et al., 2011; Karger et al., 2012).

4.3.3. Role of disturbance

Disturbance intensity directly (e.g., trampling, flexing shoots, chopping stems, harvesting plants and plant components, and grazing) as well as indirectly (e.g., soil compactness, accumulation, and decay of litter) affects plant growth, development, and regeneration (Austrheim, 2002; Sproull et al., 2015). Multiple studies (Jamtsho & Sridith, 2015; S. M. Khan et al., 2013; Shaheen et al., 2015) have revealed the prominent role of human disturbances in unsettling the sensitive Himalayan forest ecosystems. Our findings demonstrate the essential role of disturbances in framing the plant communities (Table 2). Overgrazing may exacerbate the problem by enhancing aridity and insolation, making mountainous terrain less suitable for plant growth. Additionally, grazing may alter the pattern of species dominance (Li et al., 2021; Múgica et al., 2021) and modify communities based on species’ sensitivity to grazing. Grazing eliminates many opportunistic rivals (mat- and rosette-forming species that grow near the surface to dodge grazing) of the showy blooming plants that distinguish the high alpine meadows. Several species, including Podophyllum hexandrum and Aconitum hereophyllum, designated by the IUCN as Red-List species, were collected for medicinal purposes. Additionally, these plant biosystems are threatened by rapid economic transformation, tourist boom, and societal shifts. Proactive efforts can assist in safeguarding the biota of these mountains and mitigate multiple problems that other Himalayan highlands are encountering (S. M. Khan et al., 2013).

Grazing can influence ecological conditions, notably competition among species (Anderson & Hoffman, 2007; Vera, 2000), alter vegetation structure and composition, and facilitate plant invasion. The rising temperature and retreating snow cover are other challenges that are further attributed to enhanced weed invasion in Himalayan landscape (A. M. Khan et al., 2018). Invasive species, distinguished by high reproductive success and strong competitive ability, vigorously compete with native species and have the potential to profoundly alter vegetation characteristics and ecological processes (Godefroid et al., 2005; Martins da Silva et al., 2009). Although the migration of invasive species into natural habitats often takes time, human involvement and site deterioration speed up the process (Martins da Silva et al., 2009). The higher disturbance intensity and larger numbers of landscape passages, such as trails, facilitate the dispersion of alien species, as reported in the previous studies (Liedtke et al., 2020; Wedegärtner et al., 2022). Since onsite disturbances directly impact species diversity and richness, conservation initiatives should be designed in conjunction with the local community to encourage education and proper land use in order to safeguard natural resources and diversity.

4.3.4. Role of canopy

CCA also highlights the important role of forest canopy in composition and distribution of herbaceous communities, probably by setting up unfavorable understory conditions; since their thick canopy cover limits light availability, and their persistent acidic litter reduces soil pH and nutrient availability to the plants (Chávez & Macdonald, 2010; Macdonald & Fenniak, 2007; Radhamoni et al., 2023). Inter-specific competition in conifers may not significantly impact understorey diversification patterns. Instead, these hostile environments serve as a filter, hindering light- and nutrient-demanding species, while encouraging the dominance of species which can survive and flourish under such conditions, e.g., shade-tolerant species (Hart & Chen, 2006; Macdonald & Fenniak, 2007). Decreasing canopy cover enhances solar radiation on the forest floor, accelerating plant growth by supplying energy for photosynthesis (Hofmeister et al., 2009; Tinya et al., 2009). However, others observed no effect (Augusto et al., 2003; Lenière & Houle, 2006).

5. Conclusions

Our study demonstrated empirically how cutting-edge multivariate analytical tools could help analyze plant communities and identify a list of statistically derived combination of indicator species. These results can be used to monitor changes under global change and frame better conservation and management planning in plant communities distributed across varied environmental gradients. Further, the workflow adopted might be helpful as a pathway for similar vegetation classification studies. Additionally, many distinct questions were addressed, and the information gathered enabled the following conclusions:

1. The data demonstrate that understorey herbs in Gurez valley are relatively heterogeneous across the complex landscape. TTP harbors higher species diversity attributed to favorable ecological conditions and thus should be prioritized for species conservation without undermining the importance of ITG and TTP communities (H1).

2. Further, the results showed that altitude, canopy, trails, grazing, moisture index, and invasive aliens regulate understory species distribution in forest ecosystems and that the association of plant species changed in response to environmental gradients.

3. Variations in species diversity along altitudinal gradient indicate a strong influence of environmental filtering in structuring plant diversity patterns, with canopy and altitude being the key limiting drivers for understory vegetation as expected (H2), while disturbance and moisture status impose additional restraints.

4. The impact of anthropogenic factors on species distribution is alarming, as disturbances may change the community structure and disrupt ecological processes. Initiatives, viz., regulating increasing disturbances, monitoring forest diversity change, promoting reforestation activities, facilitating awareness campaigns for sustainable resource use, and effective rangeland management should draw considerable attention.

Acknowledgments

Authors thank two anonymous reviewers for their insightful comments and constructive feedback, which greatly improved the quality of this manuscript. AAD thanks Chief Conservator of Forests, Jammu and Kashmir Forest Department, for permission to conduct the fieldwork. First author is indebted to Centre for Biodiversity and Taxonomy, University of Kashmir, for renumerated identification of plant specimens. Thanks to my three field assistants – Mushtaq Ahmad Dar, Bilal Ahmad Dar, and Mudasir Ahmad Wani and forest guard – Aijaz Ahmad Khan, without them, the fieldwork would have been almost impossible.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Supplementary material

Supplemental data for this article can be accessed online at https://doi.org/10.1080/24749508.2023.2216532.

Data availability statement

The data is available with the corresponding author and is not available on any online web platform.

Additional information

Funding

AAD received financial support from University Grant Commission through Junior Research Fellowship (Ref. No.: 3796/[NET-JULY 2018]).