Genome-wide analysis of miR172-mediated response to heavy metal stress in chickpea (Cicer arietinum L.): physiological, biochemical, and molecular insights

Abstract

Background

Chickpea (Cicer arietinum L.), a critical diploid legume in the Fabaceae family, is a rich source of protein, vitamins, and minerals. However, heavy metal toxicity severely affects its growth, yield, and quality. MicroRNAs (miRNAs) play a crucial role in regulating plant responses to both abiotic and biotic stress, including heavy metal exposure, by suppressing the expression of target genes. Plants respond to heavy metal stress through miRNA-mediated regulatory mechanisms at multiple physiological, biochemical, and molecular levels. Although the Fabaceae family is well represented in miRNA studies, chickpeas have been notably underrepresented. This study aimed to investigate the effects of heavy metal-induced stress, particularly from 100 µM concentrations of cadmium (Cd), chromium (Cr), nickel (Ni), lead (Pb), and 30 µM arsenic (As), on two chickpea varieties: ILC 482 (sensitive) and Azkan (tolerant). The assessment focused on physiological, biochemical, and molecular parameters. Furthermore, a systematic characterization of the miR172 gene family in the chickpea genome was conducted to better understand the plant’s molecular response to heavy metal stress.

Results

Variance analysis indicated significant effects of genotype (G), treatment (T), and genotype-by-treatment (GxT) interactions on plant growth, physiological, and biochemical parameters. Heavy metal stress negatively impacted plant growth in chickpea genotypes ILC 482 and Azkan. A reduction in chlorophyll content and relative leaf water content was observed, along with increased cell membrane damage. In ILC 482, the highest hydrogen peroxide (H₂O₂) levels in shoot tissue were recorded under As, Cd, and Ni treatments, while in Azkan, peak levels were observed with Pb treatment. Malondialdehyde (MDA) levels in root tissue were highest in ILC 482 under Cd and Ni exposure and in Azkan under As, Cr, and Cd treatments. Antioxidant enzyme activities, including superoxide dismutase (SOD), catalase (CAT), peroxidase (POD), and ascorbate peroxidase (APX), were significantly elevated under heavy metal stress in both genotypes. Gene expression analysis revealed upregulation of essential antioxidant enzyme genes, such as SOD, CAT, and APX, with APX showing notable increases in both shoot and root tissues compared to the control. Additionally, seven miR172 genes (miR172a, miR172b, miR172c, miR172d, miR172e, miR172f, and miR172g) were identified in the chickpea genome, distributed across five chromosomes. All genes exhibited conserved hairpin structures essential for miRNA functionality. Phylogenetic analysis grouped these miR172 genes into three clades, suggesting strong evolutionary conservation with other plant species. The expression analysis of miR172 and its target genes under heavy metal stress showed varied expression patterns, indicating their role in enhancing heavy metal tolerance in chickpea.

Conclusions

Heavy metal stress significantly impaired plant growth and physiological and biochemical parameters in chickpea genotypes, except for cell membrane damage. The findings underscore the critical role of miR172 and its target genes in modulating chickpea’s response to heavy metal stress. These insights provide a foundational understanding for developing stress-tolerant chickpea varieties through miRNA-based genetic engineering approaches.

Background

The world population is expected to reach 9.7 billion by 2050, necessitating increased food production, particularly cereals and grain legumes [1]. Moreover, there is a progressive rise in the inclination toward plant-based proteins, as evidenced by the current surge in demand for pulses [2]. Common vegetable protein sources include oilseeds, grains, and legumes cultivated commercially [3]. Chickpea (Cicer arietinum L.) is a major leguminous crop with a high concentration of carbs, starch, protein, crude fiber, vitamins, and minerals [4]. Chickpea also plays a crucial role in agriculture by facilitating the binding of free nitrogen in the air to the soil through a symbiotic relationship with Rhizobium bacteria [5]. Based on 2020 FAO data, India, Australia, Turkey, and Ethiopia, are among the top five countries in chickpea production, producing the most chickpeas for global markets with about 15.7 million tons [6]. However, the collective efficiency of chickpea production has significantly decreased in recent years because of escalating levels of heavy metals in the soil [7, 8].

In recent years, anthropogenic activities, rapid urbanization, and increased industrial activities to meet the growing population’s needs have led to heavy metal accumulation because of environmental pollution [9, 10]. Furthermore, by promoting their mobility and intensifying their hazardous effects, climate change’s direct influence may aid in releasing heavy metals from the structure of organic substances [11]. The toxicity of heavy metals is determined by their bioavailability, dosage, density, and mobility [12,13,14]. Density limits range from 3.5 g/cm3 to 7 g/cm³ [15]. Elements with a density of more than 5 g/cm³ can be categorized as non-essential and less/extremely toxic to plant growth and development [16,17,18]. The presence of toxic heavy metals such as arsenic [(As), density: 5.72 g/cm³], cadmium [(Cd), density: 8.65 g/cm3], chromium [(Cr), density: 7.2 g/cm³, nickel [(Ni), density: 8.90 g/cm³] and lead [(Pb), density: 11.35 g/cm³] in soil has become a significant environmental issue, posing high toxicity to living organisms [19,20,21,22,23]. Even low concentrations of heavy metals such as As, Cd, Pb, and mercury (Hg) disrupt biochemical and physiological processes in the plant and reduce productivity [24, 25]. Heavy metal stress negatively affects many essential enzymes for plant metabolism, resulting in the denaturation of proteins, loss of membrane integrity, decreased photosynthesis and respiration, and an imbalance of reactive oxygen species (ROS) such as superoxide radicals (O₂^•−), hydrogen peroxide (H₂O₂), hydroxyl radical (HO^•) singlet oxygen (¹O₂), and alkoxyl (RO^•) and peroxyl (ROO^•) radical [26, 27]. These extreme accumulations of ROSs have a role in causing oxidative damage by harming lipids in cell membranes and essential biomolecules such as nucleic acids, carbohydrates, proteins, and lipids, eventually triggering cell death [28,29,30,31].

Several studies have shown that metal toxicity changes the amount of chlorophyll in leaves, decouple photosynthetic enzymes, slows down growth and biomass, changes protein structure, changes cell redox state, and lowers antioxidant enzyme activity [32,33,34,35,36,37]. Plant antioxidative defense uses enzymatic and non-enzymatic components to scavenge or suppress ROS to prevent/delay cellular damage [38]. The antioxidant system of defense is composed of two groups: (i) enzymatic antioxidants, including superoxide dismutase (SOD), catalase (CAT), ascorbate peroxidase (APX), and peroxiredoxin (Prx), etc. (ii) non-enzymatic antioxidants such as vitamins [ascorbic acid (AA) and α-tocopherol], reduced glutathione (GSH), carotenoids, flavonoids and proline [39,40,41].

Plants perceive signals from their severe climate environment, regulate many signal-transduction processes, and direct and initiate their stress regulatory pathways to respond to stress [42]. These stress response regulatory mechanisms include physiological, biochemical (i.e. enzymatic antioxidants and amino acids), and molecular mechanisms. All signaling pathways have the potential to activate specific transcription factors (TFs), which in turn facilitate the transcription of essential genes responsible for maintaining plant homeostasis during periods of stress [43]. TFs are key regulators of gene expression that can simultaneously regulate many genes and play a role in the response to abiotic stress [44]. Plant-conserved mature microRNAs (miRNAs) are 20–25 nt long and generated from stem-loop sections of about 70 nt RNA precursors by RNase III Dicer-like enzymes [45, 46]. Various plant species have documented over 7,000 mature miRNAs and over 6,000 precursor miRNAs up to this present [47]. These molecules facilitate post-transcriptional silencing by binding to specific locations in the 3′ UTR region of the target messenger RNA (mRNA), leading to degradation or inhibition of mRNA translation [48]. Recent research on miRNAs and their targets has shown that plants have significant regulatory functions, particularly in adapting to different stressors, promoting plant growth and development, regulating gene interactions, and responding to hormonal changes [49,50,51].

miR172 is a conserved short non-protein-coding RNA molecule with 21 nt responsible for regulating diverse plant developmental processes and gene regulation [52]. miR172 is one of the earliest miRNA families in the kingdom of plants, having been found first in Arabidopsis by the cloning and sequencing of short RNAs [53]. There are three isoforms of miR172 in Arabidopsis encoded by the genes mir172a/b, c/d, and e [54]. The miR172 targets the mRNAs of the APETALA2-LIKE (AP2-LIKE) TF family, and this interactive relation is conserved in all spermatophytes [55]. The miR172 gene family regulates meristem length, root elongation, trichome initiation, shoot branching, and flower competence [56]. Many studies show that miR172 is essential for plant growth and response to environmental stress such as salinity, drought, heat, and heavy metals [57,58,59].

No study has explored the effects of As, Cd, Cr, Ni, and Pb on miR172 expression in chickpeas. This study examined chickpeas’ physiological, morphological, and molecular responses to five heavy metal stresses. Moreover, the miR172 gene family was identified, and chromosomal locations, 2D structures, motifs, and phylogenetic relationships were investigated.

Results

Growth parameters

SL and RL varied significantly by G and treatment T (Table 1). The observed GxT interaction effect for both SL and RL indicates a significant variation in the impact of heavy metals across different genotypes. Both genotypes (ILC 482 and Azkan) exhibited notable treatment-specific disparities in SL and RL compared to the control group in this study. The Azkan genotype showed the most significant reduction in shoot length in the As treatment group, while the ILC 482 genotype showed a significant shortening in the Cd treatment. The SL of the Ni-treated groups was more significant than that of the control group in both genotypes (Fig. S1A). Both genotypes exhibited a substantially decrease in RL when exposed to heavy metals compared to the control group. Moreover, both genotypes exhibited the lowest average RL when subjected to Cr and Cd treatments, as compared to other metal treatments (Fig. S1B).

Table 1 The variance analysis results of the parameters¹

Cell membrane damage (CMD)

According to the variance analysis results, the effect of genotype (G) on CMD was insignificant (p > 0.05), while the effect of treatment (T) and GxT interaction was very significant (p < 0.01) (Table 1). The results indicated that all heavy metal treatments significantly increased CMD levels compared to the control group in both genotypes. Notably, As treatment resulted in the highest CMD levels in both genotypes, while Cd treatment had relatively lower impacts, however, all metals showed statistically significant detrimental effects compared to the control (p < 0.001) (Fig. S2).

Relative water content (RWC)

Genotypes and treatments considerably impacted RWC. The observed interaction between genotype and treatment can be attributed to the variability in the RWC of different genotypes (Table 1). Heavy metal application resulted in a reduction in RWC in both genotypes compared to the control. Among all the metal treatments, the Ni and Pb treatments exhibited the most significant decrease in RWC in the ILC 482 genotype. On the other hand, in the Azkan genotype, Cr treatment was shown to be the heavy metal that caused the most significant reduction in RWC (Fig. S3).

Chlorophyll content (CC) (SPAD)

On chlorophyll content, the effect of genotype (G), treatment (T), and their interaction (GxT) was very significant (p < 0.01) (Table 1). The ILC 482 genotype exhibited a reduction in CC across all heavy metal treatments compared to the control. Among these treatments, Cd and Cr caused the most minor reductions, while Ni, Pb, and As treatments exhibited the lowest CC levels (Fig. S4). The chlorophyll content of the Azkan genotype varied according to the treatment. This genotype showed a significant increase in CC content in Ni, Cd, and As treatments compared to the control. At the same time, a decrease was determined in Cr and Pb treatments (Fig. S4).

Oxidative stress indicators

H₂O₂ and MDA levels were significantly affected by the treatments (Table 1). In both tissue types of both genotypes, the amount of H₂O₂ increased in heavy metal treatments compared to the control (Fig. S5). The Azkan genotype showed the maximum quantity of H₂O₂ in shoot tissue when treated with Pb, while the ILC 482 genotype showed the highest amounts in treatments with As, Cd, and Ni, respectively (Fig. S5A). For the ILC 482 genotype, the As treatment had the most significant quantity of H₂O₂ in the root tissue; for the Azkan genotype, the As, Cd, and Cr treatments had the highest amount (Fig. S5B). Moreover, both genotypes demonstrated a rise in MDA content in heavy metal treatments relative to the control, which is consistent with the results shown in H₂O₂. The ILC 482 genotype showed the greatest MDA level in shoot tissue regarding Cd and Cr treatments, whereas the Azkan genotype showed the highest MDA content in As and Cd treatments (Fig. S5C). When subjected to Cd and Ni treatments, the ILC 482 genotype exhibited the highest MDA concentration in root tissue. In contrast, the Azkan genotype showed the highest MDA concentration when treated with As, Cr, and Cd (Fig. S5D).

Antioxidant enzymes

The activities of SOD, CAT, POD, and APX enzymes were significantly affected by both genotype and treatment. Moreover, the fact that the activities of these enzymes in genotypes differed according to the treatments caused the GxT interaction to be significant (Table 1). Both genotypes showed an increase in SOD, POD, and APX in both shoot and root tissues under heavy metal treatments compared to the control. The most significant levels of SOD in shoot tissue of the genotypes were noticed in the Pb treatment for the ILC 482 genotype and the As and Cd treatments for the Azkan genotype (Fig. S6A). It was most remarkable in root tissue in the Ni and As treatments for the ILC 482 genotype and the As treatment for the Azkan genotype (Fig. S6B). The CAT enzyme activities for both genotypes (ILC-482 and Askan) show a significant increase in response to heavy metal treatments compared to control conditions in both tissues. Notably, the shoot tissue of the ILC 482 genotype exhibited the highest CAT activity following Cd and Cr treatments (Fig. 6SC). Similarly, in the Azkan genotype, Cd and Ni treatments resulted in a remarkable rise in CAT activity. In the root tissue, the ILC 482 genotype displayed the highest CAT activity in response to Cd, Pb, and As treatments, while in the Azkan genotype, the highest activity was observed following As treatments (Fig. S6D). The highest levels of POD activity were observed in shoot tissue treated with Pb in the ILC 482 genotype and As in the Azkan genotype (Fig. S6E). Root tissue exhibited the most significant levels of POD activity in all metal treatments, excluding Cd in the ILC 482 genotype and As in the Azkan genotype (Fig. S6F). Furthermore, after assessing the APX levels in the shoot and root tissues of the genotypes, it was found that the ILC 482 genotype exhibited the greatest APX concentration in both tissues under the As treatment. The Azkan genotype exhibited the greatest APX value in the shoot tissue when exposed to Ni and in the root tissue when treated with CR application (Fig. S6G and H).

Proline

The statistical analysis revealed a considerable change in proline content, which serves as a stress indicator, among different genotypes and treatments. In addition, the differences in proline content among different genotypes in the treatments resulted in significant GxT interaction (Table 1). Both genotypes exhibited a rise in proline content in reaction to heavy metal stress, with the same treatments resulting in a more notable increase than others. The genotype ILC 482 exhibited the maximum proline content in both tissues with Ni treatment and in the Azkan genotype under Cr treatment (Fig. S7A and B). Both genotypes had the lowest proline content in the control plants, suggesting that the metal treatments had a stress-inducing impact.

Correlation and PCA analysis

The length of the roots and shoots of plants and the RWC and CMD showed a negative correlation with the levels of MDA and H₂O₂. CC showed a positive relationship with MDA. In addition, H₂O₂ content exhibited a positive correlation with SOD, CAT, POD, and APX contents. Moreover, the proline content had a positive correlation with MDA, H₂O₂, SOD, CAT, POD, and APX contents (Fig. 1).

Fig. 1

Heatmap of Pearson’s correlation for different parameters of chickpea under heavy metal stress. Red indicates a positive correlation, with darker shades representing stronger correlations. Blue indicates a negative correlation, with darker shades representing stronger negative correlations. Asterisks (*) denote statistically significant correlations, with one asterisk indicating a p-value ≤ 0.05 and two indicating a p-value ≤ 0.01

PCA was conducted in the study to comprehend the fundamental structure of the dataset, which encompasses a range of biochemical and physiological indicators. Four principal components with Eigen values greater than one were identified, cumulatively explaining 97.24% of the total variance observed in the data (Table 2). PC1 accounted for 65.05% of the total variance, making it the most important contributor. The high negative loadings for several parameters, including CMD, MDA, H₂O₂, SOD, CAT, POD, and APX in both shoots and roots, suggest that PC1 primarily captures oxidative stress responses and antioxidant enzyme activities. This component likely reflects the intensity of oxidative damage and the plant’s antioxidative defense mechanisms. On the other hand, RWC showed a significant positive loading (0.945) on PC1, suggesting a contrasting pattern with oxidative stress indicators. CMD, one of the physiological characters with the highest negative charge, may be used as a heavy metal stress indicator. PC2 explained 13.85% of the variance. It exhibited positive loadings for CAT activity in shoots and root proline content, indicating a correlation between antioxidant activity and osmotic stress adaptation. The significant loadings suggest that PC2 represents the plant’s protective and adaptive responses under stress conditions, particularly involving osmolyte accumulation and specific enzymatic activity.

Table 2 Factor loadings of parameters based on principal components¹

According to biplot analysis, furthermore, the PCA biplot showed clear distinctions between the treatments, with heavy metal treatments clustering away from the control group (Fig. 2).

Fig. 2

Principal component analysis for parameters and treatments. #1 As, #2 Cd, #3 Control, #4 Cr, #5 Ni and #6 Pb. SL; Shoot length, RL; Root length, RWC; Relative water content, CMD; Cell membrane damage, CC; chlorophyll content, S; Shoot, R; Root, H₂O₂; Hydrogen peroxide, MDA; Malondialdehyde, SOD; Superoxide dismutase, POD; Peroxidase, and APX; Ascorbate peroxidase

PC1, which is mainly associated with growth-related parameters such as SL, RL, and relative RWC. These parameters are closely linked to the control group (#3), suggesting that plants in the control group exhibited better growth and water retention than those treated with heavy metals. Cd (#2) and Cr (#4) shift towards biochemical parameters like CAT and Proline, suggesting that these heavy metals are associated with oxidative stress responses, as indicated by the increased activity of these stress-related enzymes. Ni (#%5), Pb (#6), and As (#1) are positioned negatively along the PC2 axis, showing a negative correlation with root length and certain antioxidant enzyme activities (Fig. 2).

Expression analysis of SOD, CAT, and APX genes under heavy metal stress

When evaluating the expression profile of the SOD gene in ILC 482 and Azkan genotypes, we found that the mRNA expression of SOD in the shoot tissue of ILC-482 was slightly higher in all heavy metal treatments compared to the control, with significant differences between the control and treatments. The expression level of SOD for Azkan was much higher compared to ILC-482. Significant upregulation is observed under Cd and Ni treatments, particularly under Cd stress, where SOD mRNA expression is maximally induced (Fig. S8A). Furthermore, expression levels of SOD in roots of ILC 482 show a significant increase under all heavy metal treatments, especially in response to Cd, Pb, and As, with moderate increases in Ni and Pb conditions. Azkan exhibits similar trends, with AS and Pb treatments inducing the highest SOD expression in roots. Ni and Cr also had a notable induction (Fig. S8B).

All treatments significantly increased CAT mRNA expression in the shoot tissues of ILC482, with Cd and Cr showing the highest induction. The CAT expression in the shoots of Azkan follows a similar pattern, with highly significant increases in response to all treatments compared to control (Fig. S8C). Moreover, Cd, Cr, Pb, and As treatments strongly induce CAT expression in the roots, with the highest mRNA levels observed under Cd stress. The Ni treatment also showed significant upregulation. CAT mRNA expression in the root of Azkan showed a substantial increase in all the heavy metal treatments except for Cd (Fig. S8D).

APX expression in shoots of ILC 482 was significantly upregulated in all treatments, with Cr and Pb showing the highest induction. Azkan also indicated significant APX induction across all metal treatments, with Ni and Cr causing the most pronounced increases (Fig. S8E). All metal treatments highly upregulated APX in the root tissue of both ILC 482 and Azkan genotypes. The ILC 482 genotype showed the highest expression in the Pb application, while the Azkan genotype showed it in the As application (Fig. S8F).

Identification, chromosomal localization, secondary structure characterization of miR172 family genes in chickpea

Seven car-miR172 genes, namely as car-miR172a, car-miR172b, car-miR172c, car-miR172d, car-miR172e, car-miR172f, and car-miR172g, were identified in C. aritenium genome. Table 3 provides the mature sequence, chromosome information, and the start and end positions of the miR172 gene family members. The lengths of the mature sequences of the miR172 gene family were 21 nt. The car-miR172 family members were identified on five different chromosomes based on chromosomal localization data (Tables 3 and Fig. 3). The RNAfold web server was used to study and show the secondary structures of the members of the car-miR172 gene family. Upon examination, it was revealed that all genes had the hairpin structure necessary for miRNA production. All car-miRNAs were in the 3’ arm of the secondary stem-loop hairpin structure (Fig. 4).

Table 3 Mature sequence, sequence length, and chromosomal localization of Car-miR172 family members

Fig. 3

Chromosomal locations of car-miR172 gene family members. The chromosomes (CaChr1, CaChr2, CaChr4, and CaChr6) are represented as vertical bars, with their length scaled according to Megabase (Mb) pairs indicated on the left

Fig. 4

Stem-loop structures of the miR172 gene family in Chickpea. Each structure represents the precursor microRNA (pre-miRNA) sequences, with the mature miRNA sequence highlighted in red within the stem region of the hairpin

Identification of conserved motifs and phylogenetic analysis of the car-miR172 family members

As a result of the conserved region analysis of car-miR172 gene family members, it was determined that 18 nt were conserved. In Fig. 5, each letter represents a nucleotide (A, G, U, C). The height of each letter at a given position indicates the frequency and conservation of that nucleotide at that position across the miR172 family. The overall height of the stack at each position reflects the level of sequence conservation at that position, with higher stacks indicating more significant conservation. This phylogenetic tree shows the evolutionary relationships of the miR172 family among Cicer arietinum, Arabidopsis thaliana, Glycine max, Medicago truncatula, and Oryza sativa. The results showed the clustering of three major clades, A, B, and C. When the phylogenetic tree was analyzed, it was found that all Car-miR172 genes were located in clade C. Moreover, Osa-miR172a, Osa-miR172b, and Ath-miR172d genes were located in clade B, while Gma-miR172b and Gma-miR172k genes were located in clade A. (Fig. 6). These results supported the evolutionary conservation of miR172 family members in chickpea.

Fig. 5

The conserved domain analysis results of the mature sequence of the car-miR172 family

Fig. 6

Evolutionary relationship analysis of the car-miR172 gene family members. The colored symbols on the branches represent the plant species to which each microRNA belongs: Green square: Arabidopsis thaliana, Blue star: Cicer arietinum, Orange circle: Glycine max, Red triangle: Medicago truncatula, Brown triangle: Oryza sativa

Expression analysis of car-miR172 family members and their target genes under heavy metal stress

The expression status of the car-miR172 gene family members and its target genes were investigated under five heavy metal stresses in the sensitive genotype (ILC 482) and the resistant genotype (Azkan). The expression levels of car-miR172 gene family members and their target genes under heavy metal stress for each genotype are described below.

In the ILC 482 genotype, the expression level of car-miR172a was increased in all stresses except for Cd and Ni treatments in the root tissue, whereas an increase was determined only in Ni stress in shoot tissue. AP2 and LMBR1 genes in the ILC 482 genotype exhibited a high expression profile under all stresses in both tissues, unlike car-miR172a. In the Azkan genotype, the expression of car-miR172a increased in shoot tissue in all treatments, while in root tissue, it increased only in the Ni treatment. Again, target genes (AP2 and LMBR1) were increased in both tissue types of the Azkan genotype compared to car-miR172a (Fig. 7). car-miR172b exhibited an increased expression profile in all stresses in both tissue of both genotypes (Fig. 8). At the same time, the expression levels of Filamentation temperature-sensitive H (FtsH), AP2, and Arabidopsis Homolog of Trithorax1 (ATX1) genes, which are the target genes of car-miR172b, increased (Fig. 8).

Fig. 7

Expression patterns of car-miR172a and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

Fig. 8

Expression patterns of car-miR172b and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

The expression level of car-miR172c in shoot tissue of genotype ILC 482 increased under Cr and As treatments compared to control, while there was an increase in root tissue, especially Cd and As. While the expression of the car-miR172c gene was high in shoots of the Azkan genotype, it was considerably higher in Cr treatment in root tissue. Both genotypes for both tissues show increased gene expression under metal stress conditions AP2, P-type ATPases (P4-ATPase), and V-type ATPases (V-ATPase). Moreover, the expression levels vary by gene and treatment, with notable increases in Pb and As across both genotypes and tissues. In terms of these genes, Azkan has higher responses to metal treatments in both shoots and roots, while ILC 482 shows more pronounced expression in Pb and Ni treatments (Fig. 9).

Fig. 9

Expression patterns of car-miR172c and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

Upon analyzing the expression levels of the car-miR172d gene, it was observed that its expression increased in both shoot and root tissues of the ILC 482 genotype across all treatments, except for the Pb treatment in shoot tissue. The expression analysis of the car-miR172d gene revealed a significant upregulation in the shoot tissues of the Azkan genotype across all treatments, while one increase was observed in the Cr treatment. The genes AP2, LMBR1, and V-ATPase also exhibited increased expression in both genotypes and tissue types. Notably, the LMBR1 gene demonstrated the highest level of upregulation among all the genes analyzed (Fig. 10).

Fig. 10

Expression patterns of car-miR172d and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

In the ILC 482 genotype, the expression of the car-miR172e gene increased across all treatments in both shoot and root tissues, except the As treatment in shoot tissue. In contrast, the Azkan genotype exhibited an elevated expression level of the miR172e gene in shoot tissue under all heavy metal treatments, while in root tissue, a significant increase was observed only under Cr treatment. Under metal stress conditions, both genotypes showed enhanced gene expression for Ring Finger, AP2, and human UDP-glucose: glycoprotein glucosyltransferase (HUGT) (Fig. 11).

Fig. 11

Expression patterns of car-miR172e and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

In the ILC 482 genotype, the expression of car-miR172f was significantly elevated in shoot tissues under all treatments. In contrast, the increase in miR172f expression in root tissues was particularly prominent under Cd treatment. For the Azkan genotype, significant increases in miR172f expression were observed under Cd, Ni, and Pb treatments in shoot tissues, whereas root tissues exhibited notable increases specifically under Cr treatment. Furthermore, both the ILC 482 and Azkan genotypes in both tissue types showed increased expression of the AP2 and LMBR1 genes under metal stress conditions, with LMBR1 showing the most substantial upregulation (Fig. 12).

Fig. 12

Expression patterns of car-miR172f and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

The car-miR172g gene expression in shoot tissue of the ILC 482 genotype was markedly higher in treatments of Cr, Ni, and Pb than in control. Under Cr and Ni treatments, a significant increase was noted in root tissue. car-miR172g expression in shoot tissue of the Azkan genotype increased significantly only in response to As treatment. Under Pb and As treatments, higher expression levels were observed in root tissue. In addition, both genotype tissue types showed a significant increase in the expression level of AP2 and P4-ATPases, the target genes of mir172g. Its expression in shoot tissue of the ILC 482 genotype showed variability depending on the treatments (Fig. 13).

Fig. 13

Expression patterns of car-miR172fg and its target genes under heavy metal stress. The bars are presented as mean ± SE (n = 4)

Discussion

The impact of environmental stress on global agricultural production losses due to abiotic stress is becoming increasingly significant, resulting in a 70% decrease in output [60]. Heavy metal toxicity is a major abiotic stress that reduces plant biomass yield [61]. Non-essential hazardous metalloid chemicals such as As, Cd, Cr, Ni, and Pb have detrimental effects on crop productivity, quality, and human health [62,63,64]. Heavy metal toxicity has many effects on plants at the cellular and molecular levels. Heavy metal pollution harms plant root morphology, reducing plant development, chlorophyll content, levels of enzyme antioxidants, plant quality, and yield [65, 66].

Roots are the initial organs that encounter soils polluted with toxic metals [67]. The plants’ expanded root systems could assimilate the metals due to the roots’ water-absorbing capacity, which ultimately had a more significant negative impact on the roots [68]. The Azkan genotype had the highest shoot length loss in the As application group, while the ILC 482 genotype had a substantial shortening in the Cd application group. The shoot length of the Ni-treated groups was larger than that of the control group in both genotypes (Fig. S1A). In addition, the results revealed that heavy metals induced toxicity and significantly reduced the root length of the ILC 482 and Azkan cultivars. Both genotypes had the lowest mean root length after Cr and Cd treatments (Fig. S1B). The results that were obtained agree with the findings of other authors. Naz et al. [69] also reported a significant reduction in root and shoot length in chickpeas grown under increased doses of Cd, Cr, Zn, Ni, and Pb stress and all heavy metal combinations [68]. indicated that the toxicity of Cd, Cr, and Ni heavy metals in corn increased with dose, and while all heavy metals reduced root length, Cd was the most effective on shoot growth, and Cr and Ni the least. As the As dose increased, shoot and root length gradually decreased [70]. In research on Cd toxicity in three chickpea types, Ullah et al. [71] found that 50 µM Cd stress reduced root and shoot length and fresh and dry weights.

One of the main effects of heavy metal toxicity on plants is damage to cell membranes, particularly the plasma membrane [72]. Heavy metals can damage membranes by oxidizing and cross-linking protein thiols, easily by binding the hydroxyl part of phospholipids, by replacing calcium ions at key locations on cell membranes, by stopping essential proteins like H⁺-ATPase from working, or by changing the structure and flexibility of membrane lipids [73,74,75,76]. All those things lead to a drop in particular transporting activities and an increase in non-specific membrane permeability, which disturb ionic homeostasis and, in turn, the activities of several enzymes essential for fundamental cell metabolism [77]. All heavy metal applications increased CMD in both genotypes compared to the control group (Fig. S2).

Heavy metal toxicity can influence plasma membrane permeability, leading to a decrease in water content. In parallel with the results of CMD, it was determined that heavy metal applications caused a decline in RWC in both genotypes. Different research has declared that Cd, Cr, Ni, Pb, and As stress induces changes in plant water status in various plant species [78,79,80,81,82]. Furthermore, heavy metals can lead to severe dehydration in shots by limiting water movement from root parts to the shoot sections of plants [83]. In the Azkan genotype, the Cr application caused the highest decrease in RWC, while in the ILC 482 genotype, the Ni and Pb application caused the highest reduction in RWC (Fig. S3). Furthermore, RWC was a growth-related factor that affected PC1 (Fig. 2). This suggests that plants in control group (#3) grew and retained water better than those treated with heavy metals. The stability of water relations depends on the balance between transpiration and water uptake [84]. Ni toxicity decreased the area of plants’ transpiring surfaces [85]. Excess Ni in plant tissues can cause chlorosis, necrosis, water potential and transpiration decreases, growth inhibition, and mortality [86]. Jagetiya et al. [87] indicated that leaves of green gram supplied with excess Ni (10, 100, and 1,000 µM) displayed a reduction in water potential (ψw) and RWC. Vezza et al. [88] reported that the treatment of soybean plants with As did decrease their leaf RWC and root WC compared to control plants. Zhang et al. [89] show that water deficiency reduces CC and photosynthetic rate in soybeans, decreasing plant height, biomass, and seed yield.

A high CC in plants is the basis for preserving the normal photosynthesis of mesophyll cells under stress [90]. The study observed a decrease in CC in all heavy metal applications in the sensitive genotype variety. In contrast, a decrease was observed in only Cr and Pb applications in the tolerant variety. While Cd, Ni, and As reduced the CC content in the sensitive genotype, they increased the CC in the tolerant variety (Fig. S4). Piotto et al. [91] used the CC (SPAD index) as one of the metrics to assess the tolerance or sensitivity of tomatoes to Cd toxicity. Toxic Cd causes a drop in the number of chloroplasts or damage to the structure of the chloroplast, which speeds up the breakdown of chlorophyll and changes the plant’s photosynthetic capacity [35]. Also, plants under Pb toxicity have their photosynthetic pathways negatively affected because it disrupts the ultrastructure of chloroplasts and blocks the synthesis of essential pigments, including chlorophyll and carotenoids, in addition to plastoquinone [73]. Cu and Pb at 800 µM inhibit chlorophyll and carotenoid biosynthesis and their incorporation into photosynthetic machinery [92]. Navabpour et al. [93] indicated that the leaf CC was reduced significantly in wheat cultivars under Pb stress. Reduced chlorophyll (a + b) levels may account for the decline in photosynthetic rates. Heavy metals capture electrons from the photosynthetic electron transport chain, generating harmful reactive oxygen species (ROS) [92].

Heavy metal toxicity leads to the excessive buildup of ROS within the cell [94]. In this context, heavy metals are divided into two: redox-active and redox-active [95]. Redox-active heavy metals, such as Fe, Cu, Cr, and Co, play a direct role in the redox processes that occur in cells. These metals are directly involved in the redox reaction that happens inside cells and leads to the production of O₂^•− and then H₂O₂ and OH through the Haber-Weiss and Fenton reactions [74]. Redox-inactive heavy metals such as Cd, Zn, Ni, and Al can indirectly produce ROS by inducing the expression of lipoxygenase (LOX) in plant tissues and causing the oxidation of indirect polyunsaturated fatty acids [75].

This study demonstrated that Cd, Cr, Ni, Pb, and As induce oxidative stresses by overproducing H₂O₂ (Fig. S5A and B) in the roots and leaves of chickpeas. This is in agreement with the results of the presence of Cd, Cr, Ni, Pb, and As, which induces oxidative stress in the developing parts of Cicer arietinum by increasing the production of ROSs [8, 96,97,98,99]. It is generally known that Pb stress binds thiol groups and alters cell redox state, increasing ROS generation [100]. Increased ROS accumulation through membrane lipid and protein disintegration leads to lipid peroxidation, one of the main causes of abiotic stress-induced cell damage [101, 102]. MDA is a lipid peroxidation product that indicates oxidative damage [103]. The increase in MDA content is directly related to increased H₂O₂ production, suggesting oxidative damage caused by excessive ROS production [104]. As a matter of fact, in our study, a positive correlation was determined between H₂O₂ and MDA amounts in both tissue types (r = 0.334 at p < 0.05 and r = 0.624 at p < 0.01 for shoot and root, respectively) (Fig. 1).

The process of priming or buildup of H₂O₂ in plant tissues enhances their tolerance to heavy metals by stimulating the antioxidant mechanism in their tissues [67]. Plants have a complex antioxidant system to maintain cellular redox balance; this protection’s strength and potential depend primarily on the species and genotype [105]. An increase in SOD, CAT, and GR activities has been reported in Cd, Pb, Cu, and Zn-treated Pisum sativum root [106], SOD CAT, APX, GR, POD, GST, GPX, DHAR, and MDHAR in V, Cr, Ni, and Cd-treated pepper leaves [107], SOD, CAT, POD, and APX in Cd, Pb, and Ni-treated common bean [108].

All antioxidative enzymes showed parallel activity patterns in shoots and roots globally. The study found that the SOD enzyme activity was typically higher in root and shoot organs of two different genotypes (Fig. S6A and B). This outcome was expected, as other studies have already shown similar findings [101, 109,110,111,112]. Results indicated a positive correlation between the increase in SOD activity and H₂O₂ levels (Fig. 1). This may be explained by the fact that an increase in O₂^•− levels increases SOD activity, leading to the acceleration of O₂^•− degradation and consequently elevating the H₂O₂ content [113]. The application of Ni led to the highest increase in the root tissue of the sensitive genotype (Fig. S6A and B). Ni is involved in the composition of enzymes such as superoxide dismutases, glyoxalase, peptide deformylase, methyl-CoM reductases, urease, and some hydrogenases [73].

Exposure to Cd, Cr, Ni, Pb, and As increased CAT activity in both tissues of ILC 482 and Azkan genotype compared to control plants (Fig. S6C and D). A positive correlation was determined between CAT and H₂0₂, especially root tissue exposed to heavy metal (Fig. 1). CAT is a crucial enzyme that defends biomolecules from oxidative damage by decomposing H₂O₂ into water and oxygen [114]. Furthermore, the treatments with Cd and Cr showed a strong correlation with PC1 according to the results of the PC analysis (Fig. 2). These treatments were explicitly associated with CAT and Proline, which are oxidative stress responses. CAT is closely related to PODs, both structurally and functionally. It has two functions: it reacts “peroxidatically” at lower concentrations of peroxide and “catalytically” at higher concentrations [115]. Figure 1 has been examined, revealing a positive correlation between CAT and POD. Other studies supported the increase in CAT activities with Cd, Cr, Ni Pb, and As applications compared to the control [116,117,118,119].

POD is found in many plant tissues and organs, removes peroxides, inhibits membrane lipid peroxidation, and prevents membrane damage [120]. The levels of POD were increased after all heavy metal treatment (Fig. S6C and D). Upon metal ions’ entry into the plant body, a rapid accumulation of toxic compounds ensues, leading to an increase in the substrate of the POD enzyme [121]. Kang et al. [122] observed increased POD levels in potatoes under Cd stress, which supports our results. Treatments with concentrations of As 5, 25, and 125 µM in broccoli considerably enhanced the activity of CAT and POD enzymes [123].

APXs catalyze the reduction of H₂O₂ to water and serve a crucial function in plants’ antioxidant systems [124]. In our findings, all treated increased the activity of APX in chickpeas. Shortly, the highest APX enzyme activity was seen in the shoots/roots of chickpeas treated with As in the sensitive genotype (ILC 482) compared to the control (Fig. S6E and F). Ibrahim et al. [125] indicated higher APX levels in wheat seedlings under As stress, which supports our findings. Moreover, there was a positive and significant correlation (r = 0.599 and r = 0.594 at p < 0.01 for shoot and root tissues, respectively) between H₂O₂ and APX levels in both tissues (Fig. 1). In response to abiotic stress, plants release the intracellular signal H₂O₂, which effectively activates APX, contributing to the plant’s ability to cope with and adapt to stress [126].

Plants not only produce antioxidants to protect themselves, but they also create a variety of cellular biomolecules that serve as compatible and metabolic osmolytes. These biomolecules, like proline, help maintain osmotic equilibrium and counteract the harmful effects of metal toxicity [4]. Excessive cell proline accumulation is a sign of stress, and plants exposed to biotic and abiotic challenges frequently accumulate proline as a protective mechanism [127]. Proline is an osmoregulatory component and acts as a chaperone molecule to stabilize proteins and scavenge ROS [128, 129]. Plants’ ability to accumulate proline under heavy metal stress might be influenced by the direct effects of metal ions and water shortage [130]. This study showed proline accumulation in chickpea shoots and roots under the influence of heavy metals, and this accumulation was found to depend on the applied heavy metal dose and organs. Bushra et al. [131] show that the proline concentration of chickpea plants significantly increased in the highest Cd (100 mg kg^− 1 soil) level as compared with untreated plants, as that of [132] in chickpeas under Cd and Pb, of [133] in finger millet under As and Ni. The tolerant Azkan genotype showed higher proline accumulation in other heavy metals, except for Ni and Pb heavy metals. Plants that are tolerant to stress have a higher concentration of proline than those that are sensitive to stress [134]. Under Cr (120 µM) stress, two varieties of chickpeas showed increases in proline accumulation in their roots and leaves compared to control plants, with the tolerant genotype showing these more pronouncedly than the sensitive genotype (Fig. S7A and B).

Plants have evolved complex molecular defense mechanisms that include physiological mechanisms and gene regulation to cope with abiotic stress [135]. This defense mechanism can be effectively achieved through the interaction of multiple genes, cofactors, and miRNAs [136]. The mRNA profiles of the antioxidant enzymes encoding POD, CAT, and APX were studied to evaluate the antioxidative responses to heavy metal stresses applied in the study. Compared to control plants, both genotypes’ mRNA expression levels of antioxidant enzyme genes rose in response to Cd, Cr, Ni, Pb, and As stress (for the CAT gene, except for the Azkan genotype’s Cd application in root tissue) (Fig. S8). However, the expression of antioxidant enzymes such as SOD, CAT, and APX in heavy metals stressed environments has been shown to play an essential role in reducing heavy metal toxicity, as reported in studies under various plant species. In this connection, Sharma et al. [137] and Alwutayd et al. [138] reported the increased activities of SOD, CAT, and APX at Cr treatment in rice [93]. observed increased antioxidant enzyme activity (SOD, CAT, GPX, and APX) and gene overexpression in wheat leaf and root tissues under Pb stress. After 14 days of Ni stress in Leersia hexandra the expression of APX7, SODCP, SOD1, and CAT genes were significantly upregulated. Still, the expression of the glutathione metabolic cycle was reported to be down-regulated [139]. The treatments of 160, 320, and 1280 µM of Cd^+ 2 in tomato increased the expression level of the Cat2 gene compared to the control [140]. Transcript levels of LOX and SOD were significantly (p < 0.05) up-regulated during 25 µM As exposure in rice [141].

Understanding these mechanisms is crucial for developing strategies to improve plants’ stress tolerance. Despite gene function being regulated at various steps, miRNAs’ post-transcriptional transcript level regulation adds complexity to the mechanism [142]. miR172 was first comprehensively revealed in chickpeas under abiotic stress, shedding light on previously undiscovered genetic factors in this important plant species. The miR172 family is found in several plants; however, the number of members differs within plant kingdoms. In this study, 7 miR172 genes were identified in chickpeas (Table 3). Previous studies have identified 5 miR172 gene family members in A. thaliana [56, 143, 144], 4 in O. sativa [145, 146], 3 in H. vulgare [147], and 12 in G. max [148]. The car-miR172 gene family members are on five chickpea chromosomes (Fig. 3). Bansal et al. [57] demonstrated that the miR172 genes are unevenly distributed among only five (chromosomes 4, 5, 6, 9, and 11) out of the twelve chromosomes in tomato. In A. thaliana, miR172 genes have been identified on chromosomes 2, 3, and 5 [56]. Secondary structural determinants of the precursors, such as the sequences surrounding the miRNA region and stem-loops bulges, are crucial for miRNA processing [149,150,151]. This study revealed that every gene possesses the essential hairpin structure for the formation of miRNA (Fig. 4). miR172 gene family members in Elaeis guineensis also have a hairpin structure, and similarly, mature sequences are located close to the ends [152]. The four precursor sequences’ bna-mirl72 gene transcript secondary structures created varied stem-loop structures, which may affect member function [153]. The alignment of mature miR172 sequences demonstrates that seven mature miR172 sequences are highly conserved, with variations occurring only in the first, last, and medium bases (Fig. 5). Differences or similarities in mature sequences may lead to variations in functional abilities [154]. The results of the motif composition of car-miR172 strongly supported its classification through phylogenetic analysis. Analyzing the amino acid sequences of the conserved miR172 domain from C. arietinum, A. thaliana, G. max, O. Sativa, and M. truncatula revealed distinct classifications into three different groups, providing insights into their evolutionary relationships (Fig. 6). The car-miR172 genes belong to the same group and are closely related to other species and each other. Similar results were found in a study where a phylogenetic tree was constructed for all leader sequences of miR172 in the AA, BB, and CC genomes of Brassica and divided into three classes [155]. The phylogenetic analysis of the pioneer regions of miR172 between Oryza sativa and six wild cousins showed that three evolutionary branches emerged [156].

Numerous conserved and species-specific miRNAs linked to stress responses have been discovered in various crops [157,158,159]. When plant organs are exposed to heavy metal stress, many dysregulated miRNAs and the genes they target play crucial roles [160]. In plants, miRNAs responsive to metal stress are identified in multiple sets along with their targets [161]. Nevertheless, the miRNA and their target gene expression profiles in chickpea organs still need investigation in response to heavy metal stresses. During heavy metal stress, miRNAs selectively target different genes engaged in metabolic pathways related to metal stress response. These processes include the signaling of phytohormones, the absorption and distribution of sulfates, the activation of antioxidant scavengers, and the synthesis of miRNAs [162, 163]. In accordance with previous studies, the majority of the predicted target genes of various conserved and novel miRNAs were found to encode TFs such as NAC, AP2, NFYA5, ARF, MADS, bHLH, MYB, WRKY, NAC, DREB, and SPL [154, 164,165,166].

The results grouped the expression of miR172 members and their target genes (FtsH, Ring Finger, AP2, P4-ATPase, HUGT, ATX1, V-ATPase, and LMBR1) into up-regulated, down-regulated, and genes with no significant change in gene expression (Figs. 7, 8, 9, 10, 11, 12 and 13). Seven members of the car-miR172 family, namely car-miR172a, b, c, d, e, f, and g, potentially target AP2 genes (Tables S2). Moreover, these genes also showed different expression levels in the root and shoot tissues of the two genotypes (sensitive and tolerant) when under heavy metal stress. miRNAs were expressed differently between 2 genotypes against Cr stress in tobacco, which indicates the diversity in their genetic structures and different responses [167]. Additionally, the expression of target genes was determined to be higher than car-miR172a (Fig. 7). It is known that upregulation of miRNAs leads to downregulation of their target mRNAs, and the opposite can occur in a stressful environment [168]. Yu et al. [169] identified different As-responsive miRNAs in rice, revealing that these genes are involved in transport, signaling, and metabolic processes.

Lower expression of miR172a, miR172c, and miR172d in Azkan than in ILC 482 may partly explain the high As tolerance of Azkan. Under As (III) stress, decreased expression of miR172 with concomitant increased expression of AP2-like transcripts has been reported [170]. In Aquilaria sinensis, the expression level of AsERF1 increased under Cd stress, reaching 3.5-fold at the highest point compared to the control group [171]. Interestingly, miR172b and FtsH, AP2, and ATX1 target genes were increased in both tissues of both genotypes in all treatments (Fig. 8). Inal et al. [172] investigated the expression of the PvFtsH gene in response to drought and salt stress in beans. They revealed that it can respond to environmental and physiological stress factors. The expression of miR172 in cotton was decreased by 50 µM Pb but increased by 100 and 200 µM Pb treatments. It has also been reported that specific miRNAs and target genes (miR414–RAX1, miR833a–RAX2 and miR5658–BIT1) showed a positive correlation in root tissue except at 100 µM Pb concentration [173]. Both miR319 and miR159 showed upregulation under Ni stress conditions in castor beans.

In all experimental conditions involving heavy metal treatments, the expression of V-ATPase, a gene targeted by car-miR172c and car-miR172d, exhibited increased activity across all genotypes and tissue types (Figs. 9 and 10). V-ATPases enhance gene expression and transport activity in cucumber and barley when exposed to heavy metals [174]. The expression levels of Ring Finger, AP2, and HUGT target genes, which car-miR172e regulates, were elevated, regardless of car-miR172e expression in all treatments (Fig. 11). In transgenic lines overexpressing the Ring Finger gene in S. lycopersicum exposed to Cd stress, Cd accumulation in shoots and roots decreased, and antioxidant activity genes like CAT, DHAR, and MDHAR were increased, increasing the plant’s tolerance to heavy metal stress [175]. In the ILC 482 genotype, car-miR172f showed an increased expression profile in shoot tissue under all stresses, while it showed a decreased expression profile in root tissue under all stresses except Cd stress. In addition, the LMBR1 gene exhibited the highest expression profile in Cd application in root tissue (Fig. 12). miR172f, miR398, and miR857 were reported to be up-regulated under 80 µM CdCl₂ stress in Brassica napus [176]. car-miR172g expression was down-regulated in both genotypes and tissues (Fig. 13). miR172, which targets TFs such as AP2, has also been reported to be down-regulated by Cd in Solanum torvum roots [177]. This strong evidence suggests that miRNA regulates the response to environmental factors, contributing to the plant’s overall resilience. With its unique perspective, this research will be a starting point for understanding the regulatory mechanisms and biological functions of the identified miRNAs and their target genes in the future.

Conclusion

This study describes in great detail how heavy metal stress affects two types of chickpeas, ILC 482 and Azkan. It shows that growth, physiological, and biochemical parameters significantly differ depending on the genotype and treatment. When root length was evaluated, Cr and Cd applications exhibited the lowest root length averages compared to other metal applications in both genotypes. Metal stressors caused cell membrane damage and decreased relative water and chlorophyll content, indicating significant physiological disturbances. Biochemically, H₂O₂ and MDA levels were elevated, along with antioxidant enzymes like SOD, CAT, POD, and APX and non-enzymatic molecules like proline, indicating a robust antioxidative response to oxidative stress. The study found seven miR172 genes in the chickpea genome: car-miRNA172a, b, c, d, e, f, and g. The genes are spread across five chromosomes and share hairpin structures necessary for miRNA processing. These miR172s are involved in heavy metal stress and identified their possible target genes, such as AP2, FtsH, ATX1, P4-ATPase, V-ATPase Ring Finger, HUGT, and LMBR1, which may serve as TFs in various stress response pathways.

miR172 plays a significant role in the molecular response to heavy metal stress in plants, but further functional validation using gene knockdown or overexpression experiments is necessary to illuminate its precise regulatory mechanisms. Advanced techniques like CRISPR/Cas9 for genome editing can provide deeper insights into miR172’s role and uncover additional pathways involved in stress tolerance, ultimately helping produce heavy metal-resistant chickpeas. The results can be used in breeding programs to develop stress-tolerant chickpea varieties, especially for agricultural areas contaminated with heavy metals. miRNA-based genetic engineering methods are an essential tool in such breeding studies. Moreover, extending research through field experiments and cross-species analysis can help translate laboratory findings into practical agronomic solutions. Exploring the interactions between miR172-regulated pathways and phytohormones like ABA, ethylene, and auxin is essential for a deeper understanding of stress modulation in chickpeas. To better understand the response of miR172 and its pathways to multiple stress factors, it is important to understand its role in other abiotic stresses besides heavy metal stress.

Methods

Plant material and design of experiment

This study used seeds of Azkan (tolerant) and ILC 482 (sensitive) genotypes as plant materials. Azkan and ILC 482 genotypes were obtained from the Transitional Zone Agricultural Research Institute in Eskisehir, Turkey, and the GAP International Agricultural Research and Training Center in Diyarbakir, Turkey, respectively.

The seeds of these two genotypes were sterilized with sodium hypochlorite (10%) for 5 min and washed three times with pure water. Sterilized seeds were grown in 70 × 50 × 15 cm sized perlite-containing germination boxes for 12 days under controlled conditions in a plant growth chamber with 250 mmol at 25 °C and a 16-hour photoperiod. The plant seedlings were raised for three weeks in a hydroponic system (pH = 5.5–5.7) using a half-strength Hoagland solution [71, 178]. Then, treatments [control (no heavy metals), 100 µM Cd(SO₄)₂, 100 µM K₂Cr₂O₇, 100 µM NiSO₄.6H₂O, 100 µM Pb(NO₃)₂, and 30 µM NaAsO₂] were applied with half-strength Hoagland solution. One preliminary experiment was performed to determine the heavy metal concentrations of chickpeas. For this, the standard LC₅₀ method was employed, which represents the concentration at which half of the maximal growth of the chickpea root is inhibited. The research was established according to a randomized trial plan with three replications and ten plants in each replication. Plants were harvested seven days [128] after treatment and shoot length (SL) and root length (RL) were measured. Plant samples (shoots and roots) were frozen with liquid nitrogen and stored at -80 ◦C until biochemical and molecular analysis.

Cell membrane damage (CMD %)

CMD was performed according to Lutts et al. [179]. The calculation of cell membrane damage in the treatment groups was determined using the following formula.

$$\text{CMD}(\%)=(\text{EC}1/\text{EC}2)\times 100.$$

Relative water content (RWC %)

To determine the RWC ratio in leaf tissue, four leaf samples were selected from each replicate group of treatments, and their fresh weights (FW) were measured. Then, the leaf samples were soaked in pure water for 4 h, and their turgor weights (TW) were determined. They were kept in an oven at 80 °C for 48 h, and their dry weights (DW) were measured. The RWC (%) was determined using the following formula [180].

RWC (%): [(FW – DW) / (TW – DW)] x 100.

Chlorophyll content (CC) (SPAD)

SPAD-502 (Chlorophyll Meter; Konica Minolta, Tokyo, Japan) was used to determine CC. This device measured leaf transmittance at 650/940 nm electromagnetic spectrum wavelengths. Readings were taken by chance in the middle parts of 5 plants’ upper leaves using the SPAD-502 device [181, 182].

Hydrogen peroxide (H2O2), malondialdehyde (MDA), and antioxidant activities

H₂O₂ and MDA levels were determined in shoot and root tissues using the methodology described by Shams et al. [183]. The H₂O₂ concentration was measured by employing a standard calibration curve, which involved utilizing various concentrations of H₂O₂. The MDA concentration was determined by measuring the levels of thiobarbituric acid-reactive chemicals and calculating it from the absorbance curve using an extinction coefficient of 155 mmol L^− 1 cm^− 1. The antioxidant enzyme activities (SOD, CAT, POD, and APX) were measured on fresh leaf samples using the techniques of Shams et al. [183].

Proline content

Proline content was determined according to the method used by [184, 185]. Fresh shoot and root tissues (0.1 g) were homogenized with sulphosalicylic acid (3%) (2 mL) and heated at 100 °C by adding glacial acetic acid (1 mL) and Ninhydrin reagent (1 mL). Then, toluene was added to the cooled samples, and readings were taken at 520 nm wavelength. The proline level in shoot and root tissues was evaluated with a standard curve constructed from pure proline.

Statistical analysis of morphological, physiological, and biochemical characters

Data was subjected to a two-way analysis of variance using the Minitab 19.2.0 program according to a completely randomized experimental arrangement in 2 (Genotype; G) x 6 (Treatment; T) factorial design. Tukey’s HSD (honestly significant difference) test was used to identify differences between means at 0.05 level. A simple main effect explained the GxT interaction effect. Correlation analysis used the Pearson correlation test (two-tailed) to reveal the relationships among examined variables. Graphical analysis was performed using Origin 2021b software (OriginLab Cooperation, Northampton, MA, USA). Using a correlation matrix, a principal component analysis (PCA) was used to detect loading values using GraphPad Prism software, version 9.5.1. The PCA results were visually represented by constructing biplots generated using GraphPad Prism software, version 9.5.1.

Identification and classification of miR172 genes family in chickpea

miR172 precursor and mature sequences in all plants (C. arietinum, A. thaliana, G. max, O. sativa, and M. truncatula) were searched for using sRNAanno (http://www.plantsrnas.org) and PmiREN 2.0 database (Table S1) [186, 187]. The positions of miR172 gene family members on chickpea chromosomes were scanned on the genome of C. arietinum (https://phytozome-next.jgi.doe.gov/info/Carietinum_v1_0) using the BLAST tool in the Phytozome v13 database [188]. The obtained chromosome information was visualized using the MG2C v2.1 (http://mg2c.iask.in/mg2c_v2.1/) web tool [189]. The RNAfold Web Server (http://rna.tbi.univie.ac.at/cgi-bin/RNAWebSuite/RNAfold.cgi) web tool was used to predict the secondary structures of car-miR172 gene family members [190].

Identification of conserved motifs and phylogenetic analysis of the car-miR172 family members

Using the mature sequence information of car-miR172 family members, conserved regions were visualized with the SeqLogo tool in TBtools [191]. To construct a phylogenetic tree among C. arietinum, A. thaliana, G. max, O. Sativa, and M. truncatula plants, the leader sequences of miR172 genes obtained from the sRNAanno database were aligned using the ClustalW tool [192]. Then, the MEGA 11 program was used [193]. The created phylogenetic tree was visualized with the help of the Evolview v2 web tool [194].

Prediction of genes that may be targets of car-miR172

The target genes of car-miR172 were predicted with the help of the psRNATarget [195] database using default settings (Table S2). The settings in the psRNATarget database were: Maximum expectation 4.0 (Range: 0–5), length for complementarity scoring 20 (Range: 15–30 bp), Target accessibility 25.0 (Range: 0-100), Flanking length around target site for target accessibility analysis 17, Range of central mismatch leading to translational inhibition 9–11 nt. The matching gene car-miR172 gene family members were targeted, and the miRNA-target gene relationship was determined.

Expression analysis of the car-miR172 family, target, APX, CAT, and SOD genes in chickpea

Total RNA isolation in shoot and root tissues was performed with Trizol^® Reagent. According to the manufacturer’s procedure, the SuScript cDNA Synthesis Kit (CAT: RT01A025) was used to synthesize cDNA (complementary DNA). Furthermore, unlike cDNA synthesis of other genes, miRNA-cDNA synthesis was performed using miRNA-specific stem-loop RT primer (Tables S1) [143]. qRT-PCR analysis was conducted with the RotorGene^® Q Real-Time PCR System. 2X SuYBRGreen qPCR Mastermix (CAT: PCR01C0252, Ankara, Turkey) mix was 20 µl in total, including cDNA (100 ng) 4 µL, primer (forward + reverse) 1 µl, SYBR green master mix 10 µL, and dH₂O 5 µl. The reaction was carried out as follows: 95 °C for 5 min, followed by 40 cycles of 95 °C for 10 s, 55 °C for 15 s, and 72 °C for 20 s. Genes were amplified from four biological replicates, with three technical replicates for each one.

Sequence information of the primers used in the study is shown in Tables S1 and S2. Primers were designed using the web-based Primer3 (http://bioinfo.ut.ee/primer3-0.4.0/) program. The specificity of the primers was tested using PCR, which was visible on a 1.2% agarose gel, as well as their product qRT-PCR melting curves. To evaluate the normalization of expression, the actin gene was used for APX, CAT, and SOD genes, while the U6 gene [196] was used for target genes. At the same time, expression measurements of the miR172 gene family were normalized according to the weighted average normalization method [197]. The 2^−ΔΔCT normalization formula developed by [198] was used to analyze fold changes in gene expression. The data’s averages and standard deviation values were obtained in the study using four replicates, and Data were visualized with GraphPad Prizm 8. Means were compared to the control group using Dunnett’s test.