DNA barcoding and phylogenetic insights into the selected endemic flora of the Western Himalayas

The Himalayan “Him and Alaya” means “abodes of snow” mountain system, extending over approximately 2400 km across South and Central Asia. It constitutes one of the most significant global biodiversity hotspots¹. Its vast altitudinal gradient, ranging from subtropical foothills to nival zones above 7000 m, generates an extraordinary diversity of ecological niches and microhabitats². The western Himalayan region, encompassing northern Pakistan, northwestern India, and parts of western Nepal, is distinguished by its complex orogeny, varied edaphic conditions, and pronounced climatic heterogeneity³. These factors have facilitated high levels of speciation and endemism, rendering the region a repository of unique genetic resources and an invaluable center of plant diversity⁴.

Floral diversity in the western Himalaya area is the basis of ecosystem stability and is the source of abundant ecosystem services⁵. A large number of species are used as a source of medicine, fragrance, decoration, and food supplements, which are the main ingredients of traditional ethnobotanical knowledge and are the source of rural livelihoods^6,7. In fact, endemic taxa are mainly examples of unique evolutionary lineages which are irreplaceable⁸. Often, they have acquired specialized adaptations for very narrow ecological niches⁹. They have very limited distribution ranges and are thus very susceptible to environmental changes, including human-induced habitat degradation, overharvesting, invasion by alien species, and climate change-driven shrinking of their ranges¹⁰. Their extinction would not only erode global biodiversity but also diminish future opportunities for biotechnological, pharmaceutical, and ecological applications¹¹.

Despite the vast botanical importance of the region, there still exists a major taxonomic and molecular knowledge gap. Historically, conventional plant identification in the western Himalaya has been based on morphological taxonomy that, although being the basis, has limitations such as phenotypic plasticity, developmental variation, seasonal morphology, and the presence of cryptic or closely allied species^12,13. In addition, herbarium collections and floristic inventories for most of the endemic taxa are incomplete or have not been updated, whereas molecular reference sequences for these species are lacking or not at all available in public databases¹⁴. As a result, a great number of plant species in the area were misidentified, undocumented, or totally unrecognized in global biodiversity assessments, thereby hindering the conservation prioritization process¹⁵.

In this context, DNA barcoding is noted as a giant leap for plant systematics and biodiversity research¹⁶. This molecular diagnostic tool takes the use of short, standardized genetic loci typically chloroplast regions such as matK and rbcL, and nuclear regions such as ITS1 and ITS2 to produce repeatable sequence signatures for species-level identification¹⁷. DNA barcoding has several advantages over traditional approaches, such as morphological, physiological, and anatomical. It can provide the correct identification of a sample even if it is a minute or fragmentary material, it can distinguish morphologically similar or cryptic taxa and it also allows integration of regional biodiversity data into global repositories such as the Barcode of Life Data System (BOLD), and GenBank^18,19,20. Moreover, when combined with phylogenetic analysis, DNA barcoding can elucidate evolutionary relationships, assess genetic distinctiveness, and detect novel previously unreported taxa²¹.

Given the limited molecular characterization of endemic flora from the western Himalayan region, there is a critical need to integrate DNA based approaches for accurate species identification and phylogenetic assessment. Therefore, the main objective of this study was to authenticate selected endemic and rare plant species from the western Himalayas using a multi-locus DNA barcoding strategy and to evaluate their phylogenetic relationships. This study generated and analyzed sequence data from chloroplast and nuclear barcode markers to confirm taxonomic identity, assess genetic relatedness among species, and provide reliable molecular reference records for regionally important endemic taxa. The outcome of this research aims to strengthen molecular taxonomic knowledge of western Himalayan flora and support biodiversity conservation and future genetic resource management in this ecologically sensitive region.

Introduction to the study area

The current study was carried out in the Himalayan region of Pakistan. The Western Himalayan region in Pakistan, part of the broader Himalayan system, is a key biodiversity hotspot that extends across the northern areas and Azad Jammu and Kashmir, covering approximately 50,916 km²²². It stretches approximately 320 km into Pakistan and is bordered by Himachal Pradesh in India. This region includes three major ranges: Ladakh, Pir Panjal, and Zaskar²³. The area is known for its diverse habitats, primarily broadleaf and coniferous forests, which support a wide range of plant biodiversity, including many medicinal species²⁴. Protected areas such as Ayubia National Park, Gol National Park, Chitral, and Machiara National Park preserve these ecosystems²⁵. Despite its ecological importance, the region remains underexplored due to its remoteness, harsh terrain, and challenging climate. As a result, plant species distribution and diversity are poorly documented²⁶. Comprehensive research is needed to fill these gaps and support conservation efforts. The current study was conducted to explore and document the rich but understudied flora of this significant region (Fig. 1).

Field visit and sample collection

Multiple field visits were conducted across various sites in the Western Himalayan region to document and collect endemic plant species (Fig. 1). All plant material used in this study was collected under a permit issued by the Wildlife Department. Collection, handling, and experimental procedures complied with relevant institutional, national, and international guidelines and legislation for research involving plants. Species were identified based on literature and the Flora of Pakistan, focusing on those with restricted distributions and specific habitat requirements, and their reported locations were verified before fieldwork. During field visits, herbaceous plants were carefully uprooted, while branches with leaves were collected from trees, ensuring the collection of flowers, fruits, leaves, and, when necessary, roots for accurate morphological identification. Collected plant tissues were preserved in plastic bags to maintain specimen quality during transport. The plant samples were collected from the Himalayan region of Pakistan and deposited in the Herbarium of Hazara University for future reference. The voucher identification numbers listed in Table 1 were obtained with the permission of the relevant institutional ethical research cmmittee. Although species identifications are available in online floras with corresponding references, all specimens included in the present study were formally identified and authenticated by Dr. Abdul Majid, Department of Botany, Hazara University, Mansehra, who also provided the voucher specimens (Table 1).

DNA extraction and gel documentation

The CTAB (Cetyl trimethyl ammonium bromide) method was used to extract total genomic DNA with some modifications⁴⁴. A fine powder of the weighed samples was made by grinding in a mortar and pestle. For each 100 mg of ground tissue, 800 µL of pre-warmed CTAB Extraction Buffer was used at 65 °C. The mixture was vortexed thoroughly. The homogenate was then incubated for 2 h at 65 °C. After incubation, 600 µl of PCI was added, and the homogenate was centrifuged for about 20 min at 13,000 rpm. The supernatant, along with 500 µL of ice cold iso-propanol, was frozen overnight⁴⁵. The samples were again centrifuged for 20 min at 13,000 rpm, and then isopropanol was added. 500 µL of 70% ethanol was added, and the samples were centrifuged for 10 min at 13,000 rpm. The tubes were inverted for drying to remove ethanol completely⁴⁶. To dilute the dried sample, 60–80 µL of ddH₂O was added. The concentration and quality of extracted DNA were checked on 1% agarose gel electrophoresis. Added 1 gram of agarose powder in combination with 2 ml TAE in 98 ml dH2O. Preheated the solution in a microwave oven till solution became transparent. 20 uL ethidium bromide was added to the mixture, shaken well, poured into the gel tray, and combs were fitted for well production. The gel was transferred into the gel tank, and the DNA sample of 2 uL was loaded with the addition of 2ul loading dye, connected the power supply, and ran the sample for 20 min. The gel product was confirmed by UV Trans-illuminator⁴⁷.

Primer selection and PCR amplification

The most-studied and well-known primers for the nuclear and chloroplast genomes of plants were selected, as shown in Table 2. The selected primers were amplified according to the prescribed conditions. The volume of 25 µL of PCR reaction mixture was prepared in a 200 µL PCR tube, and each tube contained approximately 14 µl ddH2O, 2 µl of DNA template, 2.5 µl of 10× PCR buffer, 2 µl MgCl2, 2 µl of dNTPs 0.5 units of Taq Polymerase kit (Catalog no. K0171), and 2 µl of each forward primer and reverse primer. Amplification was performed on an Applied Biosystems 2720 Thermal Cycler. The initial step for 5 min at 94 °C was followed by 35 cycles of 30 s at 94 °C, 40 s at 52 °C, 35 s at 72 °C, and 1 cycle at 72 °C for 10 min. The amplified products were electrophoresed on 1% TAE agarose gel.

Nucleotide sequencing and analysis

The PCR-amplified products were sequenced using the Sanger method at Macrogn Inc., South Korea. And after successful sequencing, to ensure the resulting sequences accuracy and to recover numerous defined sequences, the BLAST analysis in GenBank was performed to retrieve the most similar sequences based on Max score, total aligned sequence length, query cover, E. Value, and species identity success (https://blast.ncbi.nlm.nih.gov/Blast.cgi)⁶⁵. Sequences of composed data were arranged via Bio Edit (https://bioedit.software.informer.com/7.2/) and multiple sequence alignment (MUSCLE) software (https://www.ebi.ac.uk/jdispatcher/msa/muscle) aligned by ClustalW (https://www.genome.jp/tools-bin/clustalw)^66,67,68. The resultant data was processed to compute the Kimura-2-parameter (K2P) distances for each region by MEGA software (https://www.megasoftware.net)⁶⁹.

PCR amplification and sequencing success of DNA barcode regions

Initially, 9 different primers targeting the nuclear and chloroplast genomes were used to amplify and sequence 32 species from 21 families to assess amplification, sequencing, and species identification success rates. Among them, the rbcL region showed 100% PCR amplification and 96% sequencing success rates. The matK region showed a 92% amplification rate, with a lower sequencing success rate of 56.5%. The trnH-psbA primer demonstrated 100% amplification and 80% sequencing success, whereas the ITS region showed lower rates: 80.8% amplification and 52.4% sequencing success. Additionally, five barcode regions (trnV, trnA, rbcLa, rpoB, and ycf3) were tested on six endemic species. For rbcLa, amplification success was 100% with 66.6% sequencing success. The ycf3 region showed 66.6% amplification and 75% sequencing rates, while rpoB achieved 83.3% amplification and 60% sequencing success. The trnA and trnV regions showed 83.3% amplification success, with sequencing rates of 66.6% and 83%, respectively. Overall, the highest sequence recovery was observed for rbcL, followed by trnH-psbA and matK, while ITS had the lowest amplification and sequencing success among the tested regions (Table 3).

Sequencing success at the family level showed a 100% success rate in Oleaceae, Rhamnaceae, and Balsaminaceae, while several families exhibited moderate to low success rates (Table 3). Among the tested markers, rbcL demonstrated the highest sequencing success rate, with 21 families (91%) successfully sequenced, followed by trnH-psbA, with 13 families. The matK marker showed a 50% success rate in 12 families, while ITS2 had the lowest efficiency, with a 26% success rate across six families. Additionally, trnA, trnV, and rbcLa showed an 80% success rate at the identification level, whereas ITS2 consistently exhibited the lowest performance among the tested barcode regions (Table 4). The obtained nucleotide sequences were further submitted to the online GenBank NCBI database https://www.ncbi.nlm.nih.gov/genbank/ (Table 5).

Authentication based on BLAST method

BLAST, a similarity-based algorithm, was used to authenticate species by matching query sequences to reference databases, with identification considered successful at ≥ 97% sequence identity and query coverage. The BLAST results helped refine morphological identifications, with our samples showing at least 98% similarity to reference sequences. In this study, newly generated ITS2, rbcL, trnH-psbA, trnA, and rpoB sequences demonstrated 100% query coverage, confirming high accuracy in species identification (Table 6).

Species discrimination based on Taxon DNA

Species discrimination and identification success were further assessed using the Best Match (BM) and Best Close Match (BCM) functions in Taxon DNA (Species Identifier 1.7.7), based on sequence similarity across all recovered barcode regions. Results showed that rbcL and trnH-psbA provided the highest correct recognition rates, with BM values of 84.6% and 75.0%, respectively (Table 7).

Genetic divergence of candidate barcode regions

The interspecific and intraspecific distances were calculated for distance-based analysis by using ExcaliBar software. DNA barcoding is considered efficient when there is a significant difference between interspecific and intraspecific divergence, and it is also known as the barcoding gap. The barcoding gap occurs when the maximum intraspecific distance is higher than the minimum interspecific distance. In the present study, a clear barcoding gap was found in trnH-psbA, followed by matK and rbcLa (Table 8).

Phylogenetic analysis

For phylogenetic analysis, some reported reference sequences were retrieved from the GenBank database (http://www.ncbi.nlm.nih.gov/genbank/). The phylogenetic analysis was performed using the Neighbour-joining method in MEGA 6. Overall, the phylogenetic tree based on various tested barcode regions showed strong node support and had a powerful resolution between species to differentiate the relationship. In the Neighbour joining tree method, all the candidate barcodes had different degrees of resolution power, but some had similar topology. Among all tested barcode regions, the ITS2 showed slightly low level of resolution power. In some cases, the tested barcode regions showed strong discrimination power but failed to delineate true species boundaries. Some species within the same genus were found in clusters, while others fell outside and occurred in a dispersed form, representing the studied species as an outgroup. The Pairwise Genetic Distance Matrix represented the Evolutionary divergence of the targeted species and the NCBI retrieved species.

Phylogeny of query species based on ITS2 sequence

The ITS2 sequence based phylogeny showed the relationship of query species with reference database sequences, among (LC511740.1-Anemone obtusiloba) clustered with (LC554188.1-Anemone obtusiloba) with bootstrap of 89, (LC510570.1-Epimedium elatum) with (JN010975.1-Epimedium elatum), (LC527449.1-Impatiens edgeworthii) with (JX524787.1-Impatiens glandulifera), (LC527456.1-Phytolacca latbenia) with (OL824848.1-Phytolacca latbenia Walter), (LC528221.1- Rhamnus parvifolius) with (AY626439.1-Rhamnus purpurea), and (LC532160.1 zanthoxylum armatum) with (MH016484.1-zanthoxylum armatum) with a strong bootstrap support of 69, 75, 82, 86, and 87 showed in Fig. 2.

Phylogeny of query species based on matK sequence

The matK based phylogeny showed the relationship of (LC521897.1-Bistorta amplexicaulis) with (KP172064.1-Rydingia limbata), (LC521900.1- Bupleurum lanceolatum) with (MF694831.1-Bupleurum falcatum), and (MK926096.1-Bupleurum rotundifolium), (LC527407.1-Deutzia staminea) with (MH659280.1-Deutzia parviflora), (LC527409.1-Dioscorea balcanica) with (LC327678.1-Dioscorea bulbifera), (LC527451.1-Impatiens edgeworthii) with (KX677366.1-Impatiens glandulifera), (LC527454.1-Phytolacca latbenia) with (MH659547.1-Phytolacca acinosa), and (LC528225.1-Rhododendron lepidotum) with (JF956003.1-Rhododendron lepidotum), while (LC527447.1-Galium tetraphyllum), and (LC528137.1-Potentilla erecta) with the member of genus Rhamnus with bootstraps of 44, 57, 96, 98, 99, and 100 showed in Fig. 3.

Phylogeny of query species based on trnH-psbA sequence

Similarly, phylogeny based on trnH-psbA sequences, the query sequences of (LC511742.1-Anemone obtusiloba) share clade with (LC554187.1-Anemone obtusiloba), (LC516690.1-Aesculus indica) with the member of genus Aesculus, (LC521898.1-Bistorta amplexicaulis) with (EF633739.1-Bistorta amplexicaulis), (LC527403.1-Delphinium cashmerianum) with (OK148444.1-Delphinium montanum), (LC511739.1-Epimedium elatum) with Epimedium spp, (LC528220.1-Rhamnus parvifolius) with (KP299593.1-Rhamnus sphaerosperma), (LC528224.1-Rhododendron lepidotum) with (JN046857.1-Rhododendron lepidotum), (LC528387.1-Swertia ciliata) with Swertia spp, and (LC528389.1-Thymus linearis) with bootstrap of 84 to 100 showed in Fig. 4.

Phylogeny of query species based on the rbcL sequence

In phylogeny based on rbcL sequence of query sequences showed relationship of (LC516689.1-Asculus indica) with (PP809764.1-Asculus polyneura) with bootstrap of 84, (LC521896.1-Arisaema utile) with (NC064687.1-Arisaema decipiens) and (KJ773277.1-Arisaema dracobtium) with bootstrap of 60, (LC521899.1-Bistorta amplexicaulis) with (EF653765.1-Bistorta amplexicaulis), and (FM883606.1-Bistorta amplexicaulis) with bootstrap of 64, (LC521901.1-Bupleurum lanceolatum) with (OR508826.1-Bupleurum marginatum) with bootstrap of 83, (LC527402.1-Delphinium cashmerianum) with (OK148444.1-Delphinium montanum) with bootstrap of 11, (LC527404.1-Dephne papyracea) with (MG833726.1-Dephne bholua) with bootstrap of 72, (LC527408.1-Deutzia staminea) with (NC057286.1-Deutzia glabrata) with bootstrap of 09, (LC511738.1-Epimedium elatum) with (MW483094.1-Epimedium elatum) with bootstrap of 11, (LC527448.1-Heracleum polyadenum) with (HK518827.1-Heracleum maximum) with bootstrap of 53, (LC527450.1-Impatiens edgeworthii) with (AB043532.1-Impatiens amplexicaulis) with bootstrap of 34, (LC527452.1-Incarvillea emodi) with (JQ933368.1-Incarvillea emodi) with bootstrap of 53, (LC527453.1-Lonicera japonica) with (OP388439.1-Lonicera angustifolia) with bootstrap of 45, (LC527455.1-Phytolacca latbenia) with (NC041113.1-Phytolacca insularis) with bootstrap of 44, (LC528323.1-Rhamnus parvifolius) with (AM235104.1-Rhamnus prinoides) and (ON009012.1- Rhamnus globose) with bootstrap of 57, (LC528386.1-Rhododendron lepidotum) with Rhododendron spp, (LC527457.1-Swertia cilita) with Swertia spp, (LC528388.1-Thalictrum rochebruneanum) with (JX258372.1-Thalictrum macrocarpum) with bootstrap of 56, (LC528390.1-Thymus linearis) with Thymus spp, and (LC527406.1-Ototropis elegans) with (MN267864.1-Semenovia transiliensis) and (NC064353.1-Tordyliopsis brunonis) with bootstraps of 32 showed in Fig. 5.

The present study represents the first comprehensive multilocus DNA barcoding effort targeting endemic plant species of the Western Himalayan biodiversity hotspot in Pakistan. This region, although floristically rich, is highly vulnerable to climate change, glacier melting, soil erosion, landslides, habitat fragmentation, and anthropogenic pressure^70,71. Endemic flora are particularly susceptible due to their restricted geographic distribution, medicinal exploitation, and limited conservation management⁷². In this context, accurate molecular identification is essential for documenting and safeguarding regional plant diversity⁷³.

In the current study, a total of 32 endemic species of 31 genera belonging to 24 families were molecularly characterized to establish a reference DNA barcode library. DNA barcoding provides a rapid and standardized approach for species authentication using short genomic rather than whole genomes⁷⁴. This technique has proved to be a rapid, cost effective, and accurate method for species identification⁷⁵. This technique was initially applied in biodiversity records, molecular systematics, and wildlife forensic^76,77.

Barcode performance

Nine loci (trnH-psbA, rbcL, matK, ITS2, rbcLa, trnA, trnV, ycf3, and rpoB) were evaluated amplification efficiency, sequencing success, and discriminatory power. Among these, rbcL and trnH-psbA showed the highest amplification and sequencing success rates, followed by matK. ITS2 showed comparatively lower sequencing efficiency, consistent with challenges commonly reported for nuclear ribosomal in plants.

Potential DNA barcode selection based on amplification, sequencing, and resolving the taxonomic uncertainties is a crucial step in accurate species identification⁷⁸. Quantitative species discrimination analysis demonstrated that rbcL and trnH-psbA achieved the highest identification success under Best Match (BM) and Best Close Match (BCM) criteria. The combined loci e.g., rbcL + matK improved species resolution compared to single markers, supporting the CBOL recommendation of rbcL + matK as core barcode^79,80.

The non-coding intergenic sequence trnH-psbA has widely used in phylogenetic analysis because of its high substitution rate⁸¹. In the present study, trnH-psbA was observed to be an effective marker after rbcL for amplification and sequencing, as well as for sequence recovery. Similarly, trnH-psbA showed 100% success in amplification and sequencing compared with the ITS and matK gene regions⁸². Similarly, in previous studies, the trnH-psbA revealed 97.05% amplification success rate and 92.02% the sequencing success rate⁸³.

Analysis of Kimura 2-parameter (K2P) distances revealed clear interspecific divergence exceeding intraspecific variation in most taxa, indicating the presence of barcode gap for the majority of species. However, a few closely related congeners exhibited overlapping intra-and interspecific distances, explaining cases of reduced resolution⁸⁴.

Taxonomic resolution and clarification of problematic taxa

The molecular data helped clarify identification in several morphologically similar taxa within families such as Ranunculaceae, Apiaceae, Rosaceae, and Caprifoliaceae. In certain cases, preliminary morphology-based identification required re-evaluation after BLAST and phylogenetic analysis indicated alternative placements consistent with reference sequences.

For example, closely related taxa within Apiaceae and Rosaceae that exhibit overlapping vegetative characters were more reliably separated using multilocus data. Instances of low interspecific divergence suggest possible recent speciation events or the need for additional markers in specific genera. No strong evidence of deep cryptic divergence was detected, although a few taxa warrant further population-level study.

These finding demonstrate how DNA barcoding reduces misapplied names, strengthens herbarium voucher authentication, and refines species boundaries in morphologically complex Himalayan taxa.

Phylogenetic reconstruction

Phylogenetic analysis conducted using Neighbour-Joining method through Mega software version 11⁸⁵, with appropriate substitution models and bootstrap support values, confirmed clustering patterns largely consistent with current taxonomic classification at family and genus level. Nodes with bootstrap support > 50% were retained for interpretation, and several clades showed strong statistical support. However, the plastid loci rbcL and trnH-psbA provided stable backbone resolution, while ITS2 contributed additional variation at lower taxonomic levels when successfully amplified. A concatenated multilocus dataset improved overall tree resolution compared to single-locus trees⁸⁶.

DNA barcoding is a widely accepted method to evaluate a suitable barcode region⁸⁷. However, in Best Match, Best Close Match, and species barcode, the highest rate was observed for rbcL and trnH-psbA, followed by matK and rbcLa. The rbcL was the most suitable DNA barcode region at genus level identification, followed by trnH-psbA in species identification, while other barcode regions showed moderate results in identification, consistent with the previous study of⁸⁸.

Species discrimination

A total of 32 endemic plants of 31 genera belonging to 24 families were molecularly characterized and identified at genus, species, and family-level through DNA barcoding, such as Thalictrum secundum, Aquilegia pubiflora, Anemone obtusiloba, and Delphinium cashmerianum of family (Ranunculaceae), Pimpinella stewartii (Umbeliferaceae⁸⁹, Bistorta amplexicualis (Polygonaceae)⁹⁰. Abelia triflora and Lonicera japonica (Caprifoliaceae)⁹¹, Aesculus indica of family (Spindaceae), Arisaema utile of (Arecaceae)⁹², Coriaria and Coriari nepalensis (Coriaracea), Desmodium elegans (Papilionacea), Dephni papyraceae of family (Thymelaceae)⁹³, , Epimedium elatum (Berberidaceae)⁹⁴, Galium tetraphyllum (Rubiaceae)⁹⁵. Bupleurum lanceolatum and Heracleum polyadenum (Apiaceae)⁹⁶, Impatiens edgeworthii, (Balsiminaceae), Incarvelea emodi (Bigononiacea). Potentilla erecta, Cotoneaster Cotoneaster frigidus var. affinis and Spirae hazarica from family (Rosaceae)⁹⁷, Phytolaca latbenia, (Phytolacacea), Rhamnus parvifolius (Rhamnaceae)⁹⁸, Thymus linearis (Labiateae)⁹⁹, Jasminum leptophyllum (Oleaceae)¹⁰⁰, Rhododendron lepidotum (Ericacea)¹⁰¹, Swertia ciliata (Gentianaceae)¹⁰², Zanthoxylum armatum (Rutaceae)¹⁰³, Dioscorea balcanica (Dioscoreaceae)¹⁰⁴, and Deutzia staminea (Hydrangeaceae)¹⁰⁵.

Conservation implications

Accurate species identification is fundamental for conservation planning. Misidentification of endemic or threatened plants can lead to flawed red list assessments, inappropriate conservation priorities, and ineffective protected-area management. By establishing a validated barcode reference library for 32 endemic species, this contributes molecular tools for (i) verification of threatened taxa in Red List evaluations, (ii) Monitoring of medicinal plant trade and prevention of adulteration, (iii) supporting protected-area biodiversity inventories, and (iv) detecting taxonomic ambiguities before conservation decisions. In a biodiversity hotspot such as the western Himalaya, integrating molecular authentication with traditional taxonomy enhances long-term conservation strategies and strengthens evidence-based biodiversity management.

In conclusion, the current study is the first molecular systematic investigation of the endemic flora of the western Himalayan region of Pakistan. Thirty two species were investigated, and different multi-locus markers such as ITS2, matK, rbcL, rpoB, trnA, trnV, trnH-psbA, and ycf3 were applied to check the basic criterion of amplification, sequencing, and species discrimination. Among the tested candidate barcode rbcL showed the highest success rate in amplification and sequencing, followed by trnH-psbA, and matK. At the same time, the remaining primers showed amplification rates greater than 60% in amplification, sequencing, and identification. Overall, DNA barcoding has demonstrated the taxonomic uncertainties in the available flora at the family, genus, and species levels.