Arbuscular mycorrhizal fungi enhance drought resistance in Bombax ceiba by regulating SOD family genes

Experimental design

We employed a fully randomized block design with two variables: (1) AMF treatment: presence (+AMF) or absence (−AMF) of Rhizophagus irregularis; and (2) soil water regime: well-watered (WW) at 60% water holding capacity or drought stress (DS) at 30% water holding capacity. Each replicate unit consisted of three pots, and each treatment was replicated four times.

Mycorrhizal inoculum

The AMF treatment group was inoculated with R. irregularis fungus (Glomeromycetes; Glomerales; Glomeraceae; Rhizophagus), which was obtained from the Beijing Academy of Agriculture and Forestry Sciences. It was propagated with corn (Zea mays L.) and Trifolium repens L. at the Center for Yunnan Plateau Biological Resources Protection and Utilization, Qujing Normal Uni-versity (Qujing, China). Propagule density was estimated at 226 propagules/mL using the most likely number method (Feldmann & Idczak, 1992). During transplanting, seedling roots were treated with 10 mL mycorrhizal inoculum. Control seedings received 10 mL of sterilized inoculum and 10 mL of inoculum washing solution obtained from live inoculum, which had previously been filtered through a 1 μm nylon mesh.

Plant culture

After collection from a hot-dry valley (25°40′50.06″N and 101°53′27.76″E), the surfaces of B. ceiba seeds were sterilized by soaking in a 30% hydrogen peroxide solution for 30 min and then rinsed five times with sterile water. To initiate pre-germination, sterilized seeds were incubated on sterile gauze in Petri dishes (15 cm) at ambient temperature. Seeds were washed twice daily with sterilized water. Following germination, seeds were transferred to incubation plates containing autoclaved vermiculite and incubated at 25 °C, with a photoperiod of 14 h of light and 10 h of darkness. Similarly sized seedlings (~6–7 cm) were selected and transplanted into plastic pots (26 cm diameter, 19.5 cm height). Each pot contained 8 kg of sand, which had been washed five times with tap water and autoclaved at 121 °C for 2 h for use as a substrate.

Soil water regimes

Soil moisture in WW and DS treatments was maintained at 60% and 30% of substrate water holding capacity, respectively, after Li et al. (2022). Prior to drought simulation, well-watered conditions were maintained for all pots. We began the drought simulation 4 weeks after seedlings were transplanted. Pots were weighed daily and watered as needed to maintain stable soil moisture levels for 45 days. Plants were grown in a controlled environment with temperatures ranging between 24 °C and 32 °C. They were exposed to 14 h of daylight, and humidity was maintained between 40% and 60%. Every 10 days, all pots received 100 mL of Hoagland solution for irrigation (Hoagland & Arnon, 1950).

Plant sampling and biomass measurement

At the end of the 45-day drought simulation, seedling height and basal diameter were measured using a ruler and a vernier caliper, respectively, and stems and roots were collected from each plant. To evaluate AMF colonization, roots were washed with tap water, and a little of them were conserved in FAA fixative (37% formaldehyde, glacial acetic acid, and 95% ethanol, 9:0.5:0.5 (v:v:v)). The remaining portion was flash frozen using liquid nitrogen and stored at −80 °C for future analysis.

Root staining and evaluation of mycorrhizal colonization

Root samples were stained using trypan blue, after Phillips & Hayman (1970). Prior to staining, roots were acid treated. They were first submerged in 10% KOH solution at 90 °C for 30 min before treatment with 10% H₂O₂ at ambient temperature for 10 min, followed by immersion in 2% hydrochloric acid for 5 min. Roots were then stained using trypan blue for 30 min at 90 °C. To evaluate rates of mycorrhizal colonization of hyphal, arbuscular, vesicle, and spore structures, the magnified cross-sections were utilized (McGonigle et al., 1990) using a Leica DM2500 microscope (CMS, GmbH, Wetzlar, Germany).

ROS and lipid peroxidation measurement

Powdered samples were mixed with extraction buffer (1% PVP, 1 mM EDTA, and 50 mM potassium phosphate buffer to achieve uniformity) and centrifuged at 14,000 g for 30 min at 4 °C. Soluble protein and SOD activity in the resulting supernatant was quantified using the nitroblue tetrazolium reduction approach with photochemical inhibition using the procedure described by Beyer & Fridovich (1987).

Identification of SOD family genes across the B. ceiba genome

The entire genome of B ceiba and associated protein sequences were retrieved from the GigaDB Dataset (http://gigadb.org/dataset/100445) (Gao et al., 2018). To identify SOD-encoding genes, nine reported SOD protein sequences from Arabidopsis thaliana were employed as query sequences in a local BLASTP search against the B. ceiba genome, with an E-value cutoff of 1.0 × e⁻⁵. Candidate SODs were selected based on E-values (<1.0 × e⁻¹⁰) and sequence homology values (>60%). Corresponding cDNA and protein sequences for the selected candidates were extracted. All candidate sequences that meeting the criteria were further validated using the InterPro (https://www.ebi.ac.uk/interpro/) database, the Hammer SMART (http://smart.embl-heidelberg.de/#) database, and NCBI Conserved Domain Search (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi) for conserved domain verification. Databases were queried for SOD domains, such as SOD_Fe_C-terminal (PF02777), iron/manganese SOD, SOD_Fe_N-alpha-hairpin (PF00081), and Copper/zinc SOD, Cu/Zn SOD (PF00080.23). Sequences containing the relevant domains were classified as BcSODs.

Redundant sequences were eliminated via alignment using ClustalW. Subsequent analyses employed unique sequences. Phylogenetic analysis was conducted for candidate SOD genes in B. ceiba and compared with SOD genes in Arabidopsis thaliana and various other plant species. We employed the naming conventions used for A. thaliana SOD genes (Table 1).

SOD protein properties, subcellular location, multiple sequence alignment, and phylogenetic analysis

Physicochemical properties, such as amino acid molecular weights, lengths, and isoelectric points of SOD family members were predicted using the ProtParam function in the Expasy online platform (http://web.expasy.org/protparam/). Subcellular localization of SOD protein family members was predicted using WoLFPSORT (https://wolfpsort.hgc.jp/).

Multiple amino acid sequences of BcSODs were aligned and colored using DNAMAN7.0 (http://www.lynnon.com/) (see Fig. S1). To analyze phylogenetic relationships between SOD proteins in different plants, full-length amino acid sequences of 35 SOD proteins from Arabidopsis, Oryza sativa, Populus trichocarpa, and Dimocarpus longan were used to construct a phylogenetic tree. Trees were constructed using ClustalW in MEGA 7 (http://www.megasoftware.net/) with default parameters, employing the maximum likelihood method. We combined the JTT matrix-based model with pairwise gap deletion. A total of 1,000 bootstrap replicates were carried out.

Protein-protein interaction analysis of BcSODs

The apparent physical protein-protein interactions of BcCSD, BcMSD and BcFSD were estimated using STRING v9.0 (http://string-db.org/), and the associated graphic was generated using cytoscape (version 3.10.1, https://Cytoscape.org/) (Otasek et al., 2019).

Conserved motif, gene structure, and promoter analysis of BcSODs

The discovery of conserved motifs among the 13 BcSOD proteins was performed using the Multiple Em for Motif Elicitation (MEME) (http://meme-suite.org/tools/meme) tool, with motif number set to 12 and motif width ranging from 6 to 50. Structure conservativeness was predicted for BcSODs using the NCBI website. Downloaded files were visualized using the “Visualize NCBI CDD Domain Pattern” function in TBtools 2.019 (Chen et al., 2020). To predict potential stress-responsive cis-elements in the promoter regions of BcSOD genes, 3,000 bp long sequences were extracted upstream from the transcription start sites for each gene using GFF3/GTF Manipulate in the Sequence Toolkit and Fasta Extract in TBtools 2.019. These sequences were considered to represent the promoter sequences. The distribution and types of cis-acting elements in the promoters were subsequently analyzed using the PlantCARE website. Results were visualized using the Biosequence Structure Illustrator function of TBtools 2.019.

Gene expression profile of BcSODs in B. ceiba

Total RNA was extracted from shoots and roots using a Plant RNA extraction kit (Omega Bio-Tek, Norcross, GA, USA) according to the manufacturer’s instructions. RNA concentrations were quantified using a NanoDrop 2000 (Thermo Scientific, Waltham, MA, USA), and RNA quality was confirmed using agarose gel electrophoresis. First-strand cDNA synthesis was performed using a PrimerScript^® RT Reagent Kit with gDNA Eraser (TaKaRa Bio, Dalian, China) according to the supplier’s protocol.

Root and shoot BcSOD expression in response to drought stress was profiled using quantitative real-time PCR (qRT-PCR) using a Roche LightCycle 96 machine (Roche, Basel, Germany) with SYBR Green qPCR kits (TaKaRa, Dalian, China) according to the manufacturer’s instructions. The Actin gene was used as the internal standard (Li et al., 2022). Gene-specific primers were designed for these amplifications (Table S1). Each reaction included 10.0 µL of SYBR R Premix Ex Taq TM (TaKaRa, Dalian, China), 4.0 µL of tenfold diluted cDNA as a template, 1 µL of each specific primer, and 4 µL of ddH2O, for a total volume of 20 µL. Thermal cycles were as follows: an initial step at 95 °C for 3 min, followed by 40 cycles of 95 °C for 20 s, 56 °C for 20 s, and 72 °C for 20 s. Each qRT-PCR was performed three times using four separate RNA extracts from four biological replicates to minimize inherent errors. The relative expression levels of all BcSOD genes were calculated using the 2^−ΔΔCT method (Livak & Schmittgen, 2001).

Root colonization and growth conditions

Prominent mycorrhizal formations, including arbuscules, vesicles, and inter-radical spores, were observable in the roots of B. ceiba plants from both water regime treatments (Figs. 1A–1C and Table 1), indicating that they had been colonized by AMF (R. irregularis). As anticipated, plants that were not inoculated (NMF) showed no indication of mycorrhizal colonization in either water regime treatment (Fig. 1A). Drought stress caused a significant reduction in hypha and arbuscule colonization (38.58% and 22.01%, respectively), while simultaneously increasing spore and vesicle colonization (44.25% and 28.67%, respectively) compared to well-watered (Table 1).

Seedling growth was significantly affected by the combined influence of drought stress and AMF colonization. Drought stress reduced seedling height, but AMF colonization significantly increased height and basal diameter in both the DS (49.39% and 45.59%, respectively) and the WW (45.64% and 59.89%, respectively) treatment groups compared to NMF (Figs. 2A–2C). In general, AMF plants had better growth than NMF plants in both DS and WW (Fig. 2A). Furthermore, SOD activity was higher in AMF roots and shoots in both water regime treatments. Under drought conditions, AMF colonization increased SOD notably in both roots and shoots. However, the influence of AMF colonization was minimal in the WW treatment group. Moreover, SOD activities were higher under DS relative to WW (Fig. 3).

Identification and physicochemical characterization of SOD proteins

Thirteen non-redundant putative BcSOD genes were successfully identified from the B. ceiba genome using local blastp searches. Potential matches, acquired via protein sequence analysis using documented SODs, were assessed for SOD domains using bioinformatics tools, including CDD, InterProScan, and Pfam. Our analysis identified a total of 13 SOD genes in B. ceiba (Table 2).

Amino acid identity in the Zn/Cu-SODs subgroup ranged between 21.20% and 94.10% and between 23.00% and 95.40% in the Fe/Mn-SODs subgroup (Table S2). BcSOD proteins ranged in size from 140 to 306 amino acids, with calculated molecular weights between 14.82 and 34.91 kDa. Protein isoelectric points ranged from 4.87 to 9.39. Except for BcCSD2a, all BcSODs had negative Grand Average of Hydropathicity (GRAVY) values, ranging from −0.57 to 0.064. This suggests that they were primarily hydrophilic (Gasteiger et al., 2005). Instability index values suggested that all BcSODs except for BcFSD1 and BcFSD1 were unstable (Table 2). According to WoLFPSORT predictions, three members (BcCSD1a, BcCSD1b, and BcCSD1c) were located in the cytoplasm, while five members (BcFSD1, BcFSD2, BcFSD3, BcCSD2a, and BcCSD2b) were present in chloroplasts. All three MnSODs were found in mitochondria, BcCSDc was found in the peroxisome, and BcCSD3 was found in the cytoskeleton (Table 2).

InterPro and CDD were used to identify protein domains. Conserved domains, namely Mn/Fe_SOD_N and Mn/Fe_SOD_C, were detected in sequences containing BcMSD1a, BcMSD1b, BcMSD 2, BcFSD1, BcFSD2, and BcFSD3 (Figs. 4A–4C and S1A ). Furthermore, the SOD_Cu_Zn_dom domain was detected in various sequences, including BcCSD1a, BcCSD1b, BcCSD1c, BcCSD2a, BcCSD2b, BcCSD3, and BcCSDc (Figs. 4C and S1B). BcSODs were categorized into three separate clusters using phylogenetic connections and domain recognition: MnSOD, FeSOD, and Zn/CuSOD (Fig. 4A). This categorization demonstrates the structural and functional diversity present in the BcSOD gene family in B. ceiba.

Phylogenetic analysis of BcSODs from B. ceiba and other plants

To analyze the evolutionary relationships between SODs in B. ceiba and other plants, complete BcSOD protein sequences were compared with SOD protein sequences from A. thaliana, Populus trichocarpa, Oryza sativa, and Dimocarpus longan. This comparison was used to construct an unrooted phylogenetic tree (Fig. 5 and Table S3), which indicated that the 46 SODs could be categorized into three separate clusters: Fe-SODs, Mn-SODs, and Cu/Zn-SODs. This division is based on BcSOD domains and is supported by strong bootstrap values.

Phylogenetic analysis offered a distinct understanding regarding BcSOD categorization, which was based on association with different metallic coenzymes. The Cu/Zn-SOD group included BcCSD1a, BcCSD1b, BcCSD1c, BcCSD2a, BcCSD2b, BcCSD3, and BcCSDc. BcFSD1, BcFSD2, and BcFSD3 were classified as part of the Fe-SOD group, whereas BcMSD1a, BcMSD1b, and BcMSD2 were placed in the Mn-SOD group (Fig. 5). Importantly, Cu/ZnSODs and Fe/MnSODs were clearly segregated within monocot and dicot clades. Additionally, subcellular forecasts for SODs are likely to help reinforce protein clustering. For example, AtCSD2, BcCSD2a, BcCSD2b, DLCSD2, and PtCSD2 constituted a closely associated group with a significant bootstrap value (100%), and all of these proteins were found in chloroplasts (Table 2).

Protein-protein interaction predictions within the SOD gene family

Our predictions of protein-protein interactions assume independence from gene context methods used by STRING to identify potential interacting proteins. The highest identity homologous AtSOD proteins were designate as STRING proteins based on orthologs in A. thaliana. Figure 6 and Table S5 show that all 13 BcSOD proteins were associated with nine known AtSOD proteins in the interaction network. BcCSD1a, BcCSD1b, and BcCSD1c were linked with AtCSD1; BcCSD2a and BcCSD2b were associated with AtCSD2; and AtCSD3 and AtCSDc were the highest homologous proteins of BcCSD3 and BcCSDc, respectively. For FSDs, BcFSD1, BcFSD2, and BcFSD3 corresponded to AtFSD1, AtFSD2, and AtFSD3, respectively. Additionally, BcMSD1a and BcMSD1b were associated with AtMSD1, and BcMSD2 was associated with AtMSD2. This is consistent with our phylogenetic classification, in which BcSOD proteins corresponding to three types of AtSOD proteins were grouped together in the phylogenetic tree (Fig. 5). The three types of SOD proteins form a robust interaction network, with each BcSOD putatively interacting directly or indirectly with all other BcSOD proteins, and vice versa (Fig. 6). This suggests that protein complexation has a potential regulatory role. Furthermore, the high level of interaction between SODs and catalase indicates a coordinated response of various enzymes within the plant antioxidant enzyme system in response to abiotic and biological stresses. The widespread sub-cellular localization of BcCSD, BcFSD, and BcMSD proteins indicates their reciprocal regulation and co-expression.

Analysis of cis-regulatory elements in the promoter of BcSODs

The examination of putative cis-acting elements through statistical analysis indicated a notable presence of environmental and hormone-responsive elements dispersed throughout the promoter region of BcSODs (Fig. 7 and Table S6). As depicted in Fig. 7, five phytohormone (salicylic acid (SA), gibberellin (GA), methyl jasmonate (MeJA), auxin, and abscisic acid (ABA)) correlated responsive elements including TGA-element, TGACG-motif, TCA-element, TATC-box, P-box, GARE-motif, CGTCA-motif, ABRE, etc., were identified (Fig. 7 and Table S6). The distribution of these elements varied among genes, emphasizing their pivotal role in phytohormone mediation. Moreover, environmental-responsive elements (drought, anaerobic, low-temperature, light, and defense and stress), including TCT-motif, TC-rich repeats, TCCC-motif, MBS, LTR, LAMP-element, I-box, GT1-motif, G-box, GATA-motif, ARE, AT1-motif, etc., were identified (Fig. 7 and Table S6). Notably, a multitude of light-responsive elements was widely distributed among all promoter region of BcSODs, with BcCSD1b featuring the lowest (11) and BcCSD1c the highest (26) number of light-responsive elements. Additionally, plant hormone-responsive elements (ABA, MeJA, GA, SA, and auxin) exhibited diverse distribution in the promoter region of BcSODs. For instance, BcCSDc featured only two hormone-responsive elements (CGTCA-motif and TGACG-motif) associated with MeJA-responsiveness. BcMSD1a, BcCSD1b and BcCSD1c contained all five types of hormone-responsive elements, while BcFSD1, BcMSD1b, BcMSD2 and BcCSD2b contained four types. Cluster analysis grouped specific BcSODs together, implying potential shared roles in response to similar abiotic stressors or conditions. Specifically, BcCSD1a, BcCSD2a, BcCSD1c, and BcMSD2 formed one cluster, while BcFSD1, BcFSD2, BcFSD3, BcCSDc, and BcCSD3 comprised another. BcMSD1a, BcMSD1b, BcCSD1b, and BcCSD2b fell into the same subgroup. These clusters suggest potential functional similarities among the grouped BcSODs. Overall, the results suggest that the expression levels of BcSODs may diverge under phytohormone and abiotic stress conditions.

Responses of BcSOD expression to AMF inoculation and drought stress

We measured the impact of drought stress and AMF on the mRNA levels of BcSOD. All seven Cu/ZnSODs in the root were significantly upregulated by AMF colonization under drought conditions, except for BcCSD3 (Fig. 8A). Furthermore, AMF colonization demonstrated a positive impact on the expression of BcMSD1a and BcMSD1b, while it led to a decrease in BcFSD1 transcript levels in the root under drought stress conditions (Figs. 8A and 8B). In well-watered conditions, AMF colonization exhibited a dual regulatory effect in root tissues: upregulation of BcCSDs and downregulation of BcFSD1 expression (Figs. 8C and 8D). In shoots tissues under drought stress, all seven BcCSDs were upregulated by AMF inoculation, but there was no significant effect on the BcFSDs and BcMSDs group (Fig. 8C). Under well-watered conditions, the impact of AMF colonization was relatively modest, except for the notable upregulation of BcCSD1a, BcCSD2a, and BcFSD2 in shoots tissues (Fig. 8D). In summary, AMF colonization emerges as a key modulator, particularly enhancing the expression of Cu/ZnSODs during drought stress in both root and shoot tissues.