Genome-wide identification, characterization and expression pattern analysis of APYRASE family members in response to abiotic and biotic stresses in wheat

Screening of gene sequences

For identification of APY gene family members in the wheat genome, the amino acid sequence of seven Arabidopsis thaliana APYs directly obtained from TAIR (https://www.arabidopsis.org/) were used to blast against the wheat transcriptome and genome databases (Appels, Eversole & Feuillet, 2018) using the tBLASTN program with an e-value of 1 × e⁻⁵⁰ as the threshold. The APY candidates were accepted only if the protein contained the conserved Apyrase Conserved Region (ACR) (Steinebrunner et al., 2000).

Sequence alignment and phylogenetic analyses of APYs

The protein sequences of APYs from other species were obtained from the NCBI database. The sequence alignment and phylogenic analysis were both carried out using MEGA5 software, with maximum likelihood method with 1,000 bootstrap replicates and other parameters set as default, as shown previously (Jones, Taylor & Thornton, 1992; Tamura et al., 2011). The protein theoretical molecular weight (MW) and isoelectric point (IP) were predicted using online tools (http://au.expasy.org/tools). The exon/intron structure analysis was carried out by comparing the APY CDSs and their corresponding genomic sequences using the Gene Structure Display Server, as previously described (Hao & Qiao, 2018).

Structure analysis of TaAPYs gene and amino acid sequences

The DNAman software was used to analyze the conservation property of CDS and amino acid (AA) sequences. The secondary structures of these proteins were predicted using the online tool NPS@SOPMA (https://npsa-prabi.ibcp.fr/) with default settings. The structure analysis included the percentage of each amino acid, the position and the alpha helix number, Beta Bridge, and Random coli, shown in different colors. The motif analysis was carried out using the MEME motif analysis with motif number set as 16 and other parameters as default. The membrane spanning motif method was analyzed using the online tool TMHMM (Jianyi, Ambrish & Yang, 2013) with default settings. The 3D models of tertiary structures were simulated using the Swiss-Model (https://www.swissmodel.expasy.org/) which is based on the automatic ExPASy (Expert Protein Analysis System) web server (Bertoni et al., 2017; Guex, Peitsch & Schwede, 2010; Pascal, Marco & Torsten, 2011; Waterhouse et al., 2018; Waterhouse et al., 2016).

Abiotic and biotic stress treatment

Two-week-old wheat seedlings (Yangmai 158) were used for stress treatment. As previously described, 300 mM NaCl, 42 °C, 200 mM CdCl₂, and 300 mM mannitol were separately used as salt, heat, heavy metal and drought treatment (Ni et al., 2018). The wheat leaf and root were collected at 1, 3, 6, 12 and 24 h post treatment. The Bgt spores were inoculated on the wheat leaves as previously described (Liang et al., 2019). The leaf sample was taken at 24, 48, 72 and 96 h post inoculation. The samples were immediately frozen in liquid nitrogen and stored at −80 °C. Each sample was prepared with three biological replicates.

RNA isolation and quantitative real-time PCR

Total RNA was isolated using The Plant RNA kit (Omega, Shanghai, China). The RNA quantity and quality were determined using gel electrophoresis and a Scandrop spectrophotometer (Analytikjena, Jena, Germany), and cDNA was synthesized using the TransScript one-step gDNA Removal and cDNA synthesis SuperMix (Transgen, Beijing, China). qPCR was carried out using the SYBR Green PCR Kit (Qiagen, China) on the Lightcycler96 system (Roche, Swiss). The qPCR program was set as follows: preheating: 95 °C for 10 min, one cycle; Amplification: 95 °C for 10 s, 60 °C for 20 s, and 72 °C for 20 s, 45 cycles; Melting curve: 95 °C for 2 min, 60 °C for 30 s, then continuously increased to 95 °C as previously described (Ni et al., 2018). Raw data were calculated using the software given in the Lightcycler96 system. The primer information is listed in Table S3.

Vector construction, recombinant protein expression and protein purification

The CDS of TaAPY3-1, with a removed region of the membrane spanning domain, was cloned intro the pET22a vector. The expression vector was then transformed into Escherichia coli BL21. The transformed cells were cultured in Luria-Bertani (LB) broth at 4 °C until the OD600 reached 0.8, and then were induced with 0.5 mM isopropyl beta-D-1-thiogalactopyranoside (IPTG). The induced cells were cultured at 37 °C for 5 h. The cells were then harvested by centrifugation at 10,000 ×g and resuspended with 20 ml lysis buffer (20 mM Tris, 500 mM NaCl, pH 7.6). The cells were then lysed by sonication as previously described (Ye et al., 2015).

The recombinant protein was further isolated from the denatured inclusion body using the Inclusion Body Protein Extraction Kit (Shengong, China). Briefly, the denatured protein was renatured by urea gradient dialysis in 500 ml renaturing buffers (10 mM Tris–HCl, 100 mM NaCl, 2 mM reduced glutathione, 0.2 mM oxidized glutathione, pH 8.5). The concentrations of urea in the renaturing buffers were 6.0, 4.0, 2.0, 1.0, and 0 M, respectively. Every dialysis step was carried out at 4 °C for 12 h. After the refolding process, the insoluble protein was removed by centrifugation at 10,000 ×g for 30 min at 4 °C. The amount and purity of recombinant protein were assessed using the BCA Protein Assay Kit (Shengong, China).

Protein catalytic activity assay

The reaction system was set as: 50 mM Tris–HCl, 8 mM CaCl₂, 0.25 ng/ul BSA, 2.5 mM DTT, 150 mM NaCl, and 0.05% Tween 20 (100 µl) as previously described (Dong et al., 2012). Then, 10 mM ATP (or ADP, TTP, CTP and GTP) and 4 µg purified recombinant APY were added into this system. The catalytic activity was analyzed under different conditions (temperature, pH, substrates, and ions) using the Phosphate Assay Kit (Jiancheng, China).

Genome-wide identification and phylogenetic analysis of the APY family members in wheat

The protein sequences of the seven Arabidopsis APYs and the conserved ACR domains were used as the query sequences to BLAST against the recently published wheat genome and transcriptome database (Appels, Eversole & Feuillet, 2018). After careful validation of the candidates, a total of 27 APY members were identified with top hits for the AtAPY homologs (AtAPY1-7) in the wheat genome. These wheat APY candidates exhibited high sequence similarity and all contained five Apyrase Conserved Region (ACR) domains (Fig. S1). The 27 wheat APYs were further divided into nine groups, each group with three homologs located at different genome sets (A, B and D) (Table 1). The CDS information of these APYs is listed in Table S1. Based on the sequence similarity to the Arabidopsis APY homologs, the identified wheat APYs were separately named as TaAPY1, TaAPY2, TaAPY3-1, TaAPY3-2, TaAPY3-3, TaAPY3-4, TaAPY5, TaAPY6 and TaAPY7. As shown in Table 1, the APY genes were predicted to encode 430 to 706 amino acids in length, with putative molecular weights ranging from 46.446 to 77.557 kDa, and the protein isoelectric points (PIs) from 5.93 to 9.2 (Table 1). As the three copies of the APYs from different wheat genome sets (A, B and D) had very high CDS sequence similarity (Table 1), the APYs from genome set A were used in the following bioinformatic and biochemical analysis.

To investigate the evolutionary relationships among the APYs, a phylogenetic tree was constructed using the APY homologs from Arabidopsis thaliana, Zea mays, Oryza sativa, Aegilops tauschii, Phoenix dacylifera, and some other species (Table S2). The APYs can be mainly divided into three distinct groups: I, II, and II. Specifically, TaAPY1, TaAPY2, TaAPY3-1, TaAPY3-2, TaAPY3-4, and TaAPY3-3 were clustered in Group I, TaAPY5 and TaAPY6 in Group II, and TaAPY7 in Group III (Fig. 1).

Gene structure and conserved motif analysis of the APY genes in wheat

To investigate the structural diversity of the APYs, the online structural analysis tool NPS @ SOPMA (Deléage, 2017) was used for conserved motif analysis. As showed in Fig. S1, the amino acid sequences of those nine APYs were highly conserved. Additionally, a total of 16 motifs can be detected among the APYs (Fig. S2). Generally, the APYs can be mainly divided into two groups (I, II) by the motif analysis. Nine wheat APYs all contained conserved motif 1 and 8 (Fig. 2). APYs in the same group had specific motifs, such as motif 2, 3, 16, 19 to group I, motif 12 and 13 to group II, motif 15 to TaAPY3, motif 10 to TaAPY3-1 and TaAPY3-3, and motif 5 to TaAPY7 (Fig. 2). Further, the membrane spanning motif (MSM) analysis showed that Group I APYs (TaAPY1, 2, 3-1, 3-2 and 3-3) were predicted to contain only one MSM at N-terminal, whereas Group II APYs (TaAPY5 and 6) had two MSMs, separately located at the N- and C-terminals (Fig. 3). Interestingly, it was predicted that TaAPY3-4 contained no MSMs, while TaAPY7 contained three MSMs (Fig. 3). To compare, the MSM of the seven APY members from Arabidopsis were also analyzed. The results showed that except for AtAPY5 and 7, which separately had two and three MSMs, the others only exhibited one MSM at the N-terminal, which was similar to the wheat APYs (Fig. S3). The membrane-spanning was closely related to protein allocation and transportation. Thus, these results provided important information of their potential cellular roles. Moreover, the 3D structure analysis of the nine proteins showed that TaAPY5, 6 and 7 contains four subunits while other TaAPYs only have two. Two similar subunits of the APY were linked with an extended strand surrounded by an alpha helix (Fig. 4), which was a signature character of the GDA1-CD39 nucleoside phosphatase super family. Based on the results of 3D protein structure simulation (Fig. 4), the nine wheat APYs can also be divided into two groups: Group I (TaAPY1, 2, 3-1, 3-2, 3-3 and 3-4), and Group II (TaAPY5, 6 and 7), similar to the motif analysis. Although TaAPY5, TaAPY6, and TaAPY7 were categorized into different groups in the phylogenetic trees (Fig. 1), they shared high 3D structure similarity (Fig. 4).

Gene expression pattern of the APYs in wheat seedlings in response to abiotic and abiotic stresses

To investigate the expression pattern of the nine APYs in wheat, the gene expression in shoot and root, and their responses to abiotic stresses (salt, heat, heavy metal and osmotic stresses) were further analyzed using quantitative real-time PCR. The results showed that all nine APYs had similar expression levels at the seedling leaf and root, with TaAPY3-4, TaAPY3-1 and TaAPY3-3 having the highest expression level, and TaAPY6 the lowest both in the shoot and root (Fig. 5), suggesting TaAPY3-1, TaAPY3-3 and TaAPY3-4 could be the predominant APYs in wheat.

Previous studies have shown that several APYs could be involved in the regulation of abiotic stress adaptation (Clark et al., 2014). Thus, we further investigated the time-course expression profiles of TaAPYs in the leaf and root of wheat seedlings in response to different abiotic stresses. The results showed that in leaves, most of the TaAPYs could be upregulated after being subjected to cadmium treatment (Fig. 6A). Specifically, TaAPY3-4 exhibited high expression level both in the leaf and root under cadmium stress, while TaAPY3-1 reached a peak at 6 h in the leaves. The expression of the other TaAPYs in the root was not as sensitive as in the leaf, and only TaAPY6 exhibited a significant upregulation at 6 h (Fig. 6A). Further, mannitol treatment was used to produce an artificial drought stress condition in the wheat seedlings. The results showed that all the TaAPYs could be up-regulated within 24 h, among which, TaAPY1, TaAPY3-4, and TaAPY6 reached an extremely high expression level in the leaves (Fig. 6B). On the contrary, very few genes exhibited significant changes at different times post-mannitol-treatment in the root (Fig. 6B), suggesting that these TaAPYs also regulated the drought responses in the shoot, not only in the root. Under heat stress, the expression of TaAPY3-2, TaAPY3-4, TaAPY5, TaAPY6, and TaAPY7 began to increase both in the leaves and root at 12 h, with the highest increase fold of TaAPY3-2 (Fig. 6C). For salt stress, the expression of all the TaAPYs significantly increased in the leaves at 12 h, but only TaAPY1, TaAPY3-4, TaAPY7, and TaAPY5 were shortly up-regulated in the root, while others remained unaffected throughout (Fig. 6D).

To test for biotic stress, powdery mildew pathogen was used in the wheat seedlings. The results showed that seven TaAPY s (except TaAPY1) showed significant sensitivity to Blumeria graminis infection. Specifically, the significant up-regulation of TaAPY2, TaAPY3-1, TaAPY3-2, TaAPY5 and TaAPY6 could be detected at the pre-penetration stage (24 h), whereas TaAPY3-3, TaAPY3-4, and TaAPY7 were significantly up-regulated at the late infection stages (Fig. 7). These results suggested that the wheat APYs could be involved in the regulation of biotic stress responses. However, different APYs may have diverse roles at different infection stages.

Enzymatic analysis of the recombinant APY3-1 in wheat

To further validate the enzymatic activity of the wheat APYs, the recombinant APY3-1, without the cross-membrane domain, was cloned and purified using the Escherichia coli expression system (Figs. 8A and 8B). The production of the inorganic phosphate in the system was used to determine ATP hydrolyzation as previously described (Dong et al., 2012). The results showed that the recombinant TaAPY3-1 exhibited high enzymatic activity with a relatively wide temperature range from 25 to 52 °C, and had the highest catalytic activity at 37 °C (Fig. 8C). Further, the APY3-1 had relatively high enzyme activities at the acid conditions, with maximal activity detected at pH 4.5 to 5.5 (Fig. 8D). Moreover, the hydrolyzation efficiency of APY3-1 on different substrates was also evaluated. The results showed that the APY3-1 exhibited slightly lower enzymatic activity when degrading the ADP compared with ATP, but had very low activity during the degradation of TTP, GTP and CTP, suggesting that TaAPY3-1 has a high substrate preference (Fig. 8E). Moreover, it has been demonstrated that all NTPDase family members require divergent metal ions as cofactors for their enzymatic activity. Thus, we also evaluated the enzymatic activity of different ions in degrading the ATP. The results showed that Ca²⁺ was the most effective cofactor. The preference order was as follows: Ca²⁺ > Mg²⁺ > Zn²⁺ (Fig. 8F). Without the ions, the APY3-1 had no enzymatic activity, suggesting that apyrase activity is dependent on ions as cofactors, with a preference for Ca²⁺. The catalytic activity of TaAPY3-1 can reach a peak of Vmax = 61 (Pi µM/h/µg protein), Km = 8.7 mM under the most appropriate conditions (37 °C, pH 5.5, 8 mM Ca²⁺) (Fig. 8G). Conclusively, these results suggested that the APY3-1 had high and specific ATP and ADP degradation activities, which was consistent with the functions of the APY homologs reported in other species.