Genome-wide characterization and expression analysis of PP2CA family members in response to ABA and osmotic stress in Gossypium

Analysis of the PP2C family in four Gossypium species

To explore all the members of the PP2C family in Gossypium, the protein sequences of 80 AtPP2Cs were initially applied as queries to search against the databases of G. arboreum (BGI-CGB v2.0 assembly genome), G. raimondii (JGI assembly v2.0 data), G. hirsutum (NAU-NBI v1.1 assembly genome) (www.cottongen.org), and G. barbadense (http://database.chgc.sh.cn/cotton/index.html), respectively, using the BLAST program with default setting (E-value < e⁻¹⁰) (Camacho et al., 2009). After removing the redundant sequences from the data set, the putative Gossypium PP2Cs were then characterized using the PP2C model (PF00481) (http://pfam.xfam.org/) by the Hmmer software (http://hmmer.org/), and the proteins without a PP2C catalytic domain were deleted. The molecular weight (MW) and isoelectric point of the Gossypium PP2Cs were predicted by the online tool ExPaSy (http://web.expasy.org/protparam/), which can give various physico-chemical properties of a protein based on its amino acid sequence (the extinction coefficient and the absorbance of a native protein in water at 280 nm were used). The composition and position of exons and introns of the PP2CAs were obtained from CottonGen (https://www.cottongen.org/) and characterized by the Gene Structure Display Server tools (http://gsds.cbi.pku.edu.cn/) (Hu et al., 2015). The conserved domains of PP2CAs were validated in NCBI (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi) using the automatic mode (Marchler-Bauer et al., 2017). The MEME program (meme-suite.org/tools/meme) was applied to determine the motifs of PP2CAs in Gossypium (“any number of repetitions” to be distributed in sequences was set). The locations of Gossypium PP2CAs in chromosomes were assessed using MapInspect software (http://www.mybiosoftware.com/mapinspect-compare-display-linkage-maps.html).

Analysis of synteny and Ka/Ks of PP2CAs

The homologous regions of PP2CAs in Gossypium were identified by MCScanx software (http://chibba.pgml.uga.edu/mcscan2/) and syntenic blocks were determined by the CIRCOS program (http://www.circos.ca/). The syntenic maps of the PP2CAs were obtained using the circos-0.69±3 software with default parameters (http://www.circos.ca/). Some genes located within the same or adjacent intergenic region were regarded as tandem duplications. The nucleotide substitution parameters Ka (nonsynonymous) and Ks (synonymous) were assessed by the PAML program (http://abacus.gene.ucl.ac.uk/software/paml.html). Then, the ratio of Ka/Ks was calculated. Ka/Ks <1 means purifying selection; Ka/Ks = 1 indicates neutral selection, while Ka/Ks >1 represents positive selection (Hurst, 2002).

Phylogenetic analysis of PP2CAs

The PP2CA databases were downloaded for Arabidopsis thaliana (http://www.arabidopsis.org/), Theobroma cacao (http://cocoagendb.cirad.fr), Ricinus communis (http://castorbean.jcvi.org), Populus trichocarpa (http://www.phytozome.net/poplar), Glycine max (http://www.phytozome.net/soybean), B. distachyon (http://plants.ensembl.org/Brachypodium_distachyon/Info/Index), Oryza sativa (http://rapdb.dna.affrc.go.jp), and the four Gossypium species mentioned above. The amino acid sequences of PP2CAs were aligned by the MUSCLE software (Edgar, 2004), and a phylogenetic tree of the PP2CAs was generated using the IQ-TREE server (http://www.iqtree.org/) following the maximum likelihood (ML) method (Nguyen et al., 2015; Trifinopoulos et al., 2016). The best-fitting model was chosen using the ModelFinder (Kalyaanamoorthy et al., 2017) and the support was assessed using the ultrafast bootstrap (Minh, Nguyen & Von Haeseler, 2013). Evolutionary tree was visualized using the FigTree v1.4.4 software (available from https://github.com/rambaut/figtree/releases).

Measurements of GhPP2CAs expression in tissues and in response to ABA or osmotic stress

For measuring the expression of GhPP2CAs in tissues in each experiment, seeds of G. hirsutum L. acc. Texas Marker-1 (TM-1) were sown in pots containing the mixed nutrient soil (rich soil:vermiculite = 2:1, v/v) in a growth chamber. After 21 days, about two gram samples of roots, stems, or leaves were collected from 10 plants. About 20 flowers were collected 1 day post anthesis, and about five gram fibers were obtained from ovules 23 days post anthesis of cotton plants grown in the fields. For monitoring the expression of GhPP2CAs after ABA treatment or under osmotic stress, TM-1 seeds were germinated and planted in liquid 1/2 MS medium (Murashige & Skoog, 1962) in a growth chamber (the medium was aerated. day/night temperature cycle of 28/26 °C, 14 h light/10 h dark, and about 50% relative humidity). Three weeks later, the plants were sprayed with 100 μM ABA or treated with 10% PEG6000 (dissolved in medium) for 0, 3, 6, 12, and 24 h, respectively. Then, about two gram roots were sampled, frozen in liquid nitrogen and stored at −70 °C. Total RNA was extracted from some of the samples and cDNA was generated according to the method described previously (Ma et al., 2012; Zhang et al., 2017b).

Quantitative real-time RT-PCR (qRT-PCR) experiments were constructed in an ABI 7,500 real-time PCR amplifier using the cDNA, SYBR Green Master mix, the specific primers of GhPP2CA genes (Table S1). GhUBQ7 was used as an internal control (Lu et al., 2017). Experiments were independently repeated three times. The interval between two repeated experiments was 7–10 days.

Monitoring protein interaction by yeast-two hybrid method

The CDS sequences of GhPYLs (GhPYL2-2D, GhPYL6-2A, and GhPYL9-2A) and GhPP2CAs were amplified, and cloned into pGADT7 and pGBKT7 vectors, respectively, using gene specific primers (Table S2). After sequencing, the fused vectors were transformed into AH109. The cotransformants were plated on nonselective SD/-Leu/-Trp solid medium and selective SD/-Leu/-Trp/-His/-Ade solid medium as described previously (Lu et al., 2017; Zhang et al., 2017b).

Genome-wide analysis of PP2CAs in four Gossypium species

To identify the putative PP2CA family members in Gossypium, the amino acid sequences of 80 Arabidopsis PP2Cs (Xue et al., 2008) were used to survey the Gossypium databases. Putative PP2Cs were assigned to a total of 114, 116, 239, and 232 genomic sequences that were retrieved from G. arboreum, G. raimondii, G. hirsutum, and G. barbadense, respectively. They were individually denominated as GaPP2Cs, GrPP2Cs, GhPP2Cs, and GbPP2Cs (Table S3). According to the phylogenetic relationships of PP2Cs between Gossypium and Arabidopsis, the Gossypium PP2Cs could be clustered into 12 clades (A–L) (Figs. S1–S4). The PP2CAs possessed 14 GaPP2CAs, 13 GrPP2CAs, 27 GhPP2CAs, and 23 GbPP2CAs, respectively. They were named individually according to their gene identifiers (Table 1). In this report, we focused on the PP2CA family members in the four Gossypium species.

It was found that the predicted coded amino acid lengths of Gossypium PP2CAs ranged from 118 to 593, with an average of 420. These PP2CAs had MWs of 12.8 kDa (GhPP2CA27) to 66 kDa (GaPP2CA6). The mean theoretical pIs of PP2CAs were approximately 5.9 with a minimum of 4.65 (GrPP2CA5) and a maximum of 8.74 (GhPP2CA22) (Table 1).

Phylogenetic and structural analysis of PP2CAs in Gossypium

In order to understand the evolutionary relationship among GaPP2CAs, GrPP2CAs, GhPP2CAs, and GbPP2CAs, we conducted a phylogenetic tree using the protein sequences of the Gossypium PP2CAs (Fig. 1A). As expected, most of GaPP2CAs were individually clustered closely with their corresponding orthologs of GhPP2CAs (GhPP2CA1–13) and GbPP2CAs (GbPP2CA1–10) in A genomes, and a majority of GrPP2CAs individually clustered closely with their homologs of GhPP2CAs (GhPP2CA14–27) and GbPP2CAs (GbPP2CA11–23) in D genomes. Noteworthily, GaPP2CA10 clustered together with GaPP2CA11, and a similar case occurred between GbPP2CA21/GbPP2CA22 and GbPP2CA16/GbPP2CA23. Moreover, homologues of 14 GaPP2CAs and of 13 GrPP2CAs were found in the G. hirsutum At and Dt subgenomes, respectively; and homologs of 10 GaPP2CAs (except GaPP2CA1, GaPP2CA3, GaPP2CA4, and GaPP2CA14) and 11 GrPP2CAs (except GrPP2CA6 and GrPP2CA9) were detected in the G. barbadense At′ or Dt′ subgenomes, respectively. Additionally, three pairs of paralogues with high sequence similarity in GbPP2CAs including GbPP2CA8/GbPP2CA10, GbPP2CA16/GbPP2CA23, GbPP2CA21/GbPP2CA22 were clustered together. They were seemingly derived from GaPP2CA8, GrPP2CA7, and GrPP2CA1, respectively.

Most PP2CAs had three to four exons except that GaPP2CA6, GrPP2CA9, GhPP2CA22, GbPP2CA13 possessed five exons, and GhPP2CA27 had two exons. Among the PP2CA genes, GaPP2CA6, GrPP2CA9, GhPP2CA9, and GhPP2CA22 individually had a longer intron sequence than other genes did (Fig. 1B). These results indicate that the exon/intron structures of the Gossypium PP2CA genes were highly conserved.

The motif compositions of the PP2CA proteins were analyzed in the four Gossypium species. Twenty putative motifs named motif 1 to motif 20 were identified. Among those, motif 1, 2, 3, 4, 5, 6, and 7 existed in every cluster and the majority of the PP2CA members. Moreover, most orthologous PP2CA proteins in the four Gossypium plants had the same or very similar compositions and distributions of motifs, suggesting that the PP2CA members in the same cluster likely share similar functions (Fig. 1C).

Chromosomal distributions of PP2CAs in Gossypium

To determine the putative evolutionary relationships of the Gossypium PP2CA genes, the positions of the genes on chromosomes were analyzed. We found that the distributions of these PP2CAs were uneven. The 14 GaPP2CAs, 13 GrPP2CAs, 27 GhPP2CAs, and 23 GbPP2CAs were distributed on 8, 9, 16, and 14 chromosomes, respectively. Most of the chromosomes contained one PP2CA gene. By contrast, some chromosomes individually had two PP2CA genes. They were D11 and D13 in G. raimondii, At10, At13, Dt05, Dt10, and Dt13 in G. hirsutum, and Dt′07 and Dt′13 in G. barbadense. Besides, each of D09 and Dt′05 owned three PP2CAs. A13, At05, and At′05 separately possessed four PP2CAs. In contrast, GaPP2CA5, GaPP2CA8, GaPP2CA14, GhPP2CA18, GhPP2CA20, GhPP2CA27, GbPP2CA2, and GbPP2CA12 were located on scaffolds, in which contigs were not spliced into any chromosome in genomic mapping.

We compared the positions of the orthologs of GaPP2CAs, GrPP2CAs, and GhPP2CAs or GbPP2CAs in chromosomes. As expected, most homologs of GaPP2CAs and GrPP2CAs in G. hirsutum were located in their corresponding At subgenomes and Dt subgenomes. A similar situation also occurred in G. barbadense. However, homologous genes of GaPP2CAs and GrPP2CAs were barely located in their corresponding homoeologous chromosomes and collinear loci in G. hirsutum and G. barbadense (Fig. 2). For example, GaPP2CA7 was found in chromosome A09, however, its ortholog GhPP2CA3 and GbPP2CA9 were present in chromosome At05 and At′05, respectively. These results imply that specific, unique, and complex variation events in PP2CA-contained homoeologous chromosomes may happen within each of the two diploid and tetraploid species during genetic evolution.

Synteny analysis of PP2CA genes

During evolutionary processes, tandem and segmental duplications contribute to expanding gene family in plants (Cannon et al., 2004). We examined the duplication relationship of the PP2CAs among G. arboreum, G. raimondii, and G. hirsutum (the related database for G. barbadense was lacking). It was found that GaPP2CA10 and GaPP2CA11 joined together, and GbPP2CA16 and GbPP2CA23 clustered together in the chromosome. There are less than five genes between each pair of the genes, suggesting that the two pairs of genes are tandemly duplicated.

The synteny relationship of gene pairs was also explored among GaPP2CAs, GrPP2CAs, and GhPP2CAs. A total of 136 homologous gene pairs were observed in 133 collinearity blocks. Most of the blocks had one gene pair. Some blocks owned two gene pairs (GrPP2CA9/GhPP2CA22, GrPP2CA10/GhPP2CA23) between chromosomal D11 and Dt10. Another block harbored three gene pairs (GrPP2CA5/GhPP2CA2, GrPP2CA6/GhPP2CA3, GrPP2CA7/GhPP2CA4) between chromosomal D09 and At05 (Fig. 3). These findings imply that segmental duplication plays major roles in generating PP2CAs during evolution in Gossypium.

Analysis of Ka/Ks values of PP2CAs

To further understand the evolution processes among Gossypium PP2CAs, the effects of selection on duplication of PP2CA genes were determined. The Ka and Ks substitutions, and Ka/Ks values were calculated for the homologous gene pairs among GaPP2CAs, GrPP2CAs, and GhPP2CAs. The mean values of Ka/Ks for these gene pairs between species Ga/Gh, Gr/Ga, Gr/Gh, Ga/Ga, Gr/Gr, Gh/Gh were 0.22, 0.21, 0.22, 0.21, 0.19, and 0.21, respectively. All of them were less than 1, indicating that the formation of these genes were mainly under purifying selection during evolution. The Ka/Ks ratios for the two gene pairs GrPP2CA11/GhPP2CA15 and GrPP2CA3/GhPP2CA21 were higher than 1, suggesting that the two gene pairs were generated under positive selection and the selection likely has effects to change these genes during evolution (Fig. 4).

Phylogenetic analysis of PP2CAs in Gossypium and other plants

We constructed a phylogenetic tree of PP2CA proteins in G. arboreum, G. raimondii, G. hirsutum, G. barbadense, A. thaliana, T. cacao, R. communis, P. trichocarpa, G. max, B. distachyon, and O. sativa using the ML method, and analyzed the evolutionary relationships of these PP2CAs. It was found that the PP2CAs included both dicotyledonous and monocotyledonous members (Fig. 5). This suggests that these PP2CAs formed before the divergence of eudicots and monocots and are in general highly conserved. Indeed, the PP2CAs from the eudicots Gossypium, cacao, poplar, castor, soybean and Arabidopsis clustered more closely, and those of the monocots rice and distachyon clustered together. Moreover, many PP2CAs from Gossypium clustered more closely with those from cacao than from poplar, castor, soybean and Arabidopsis (Fig. 5), indicating that PP2CAs of Gossypium had closer relationship with those of cacao than those of other plants. As expected, PP2CAs in the four Gossypium species always clustered together, in line with their homologous evolutionary relationships (Fig. 5).

Expression patterns of GhPP2CA genes in different tissues

The transcript abundances of 27 GhPP2CAs in various tissues were measured by qRT-PCR to determine the putative functions of the PP2CAs in cotton. The results showed that all of the GhPP2CAs were highly expressed in flowers. GhPP2CA11 and GhPP2CA27 were also preferentially expressed in roots. Moreover, the transcriptional levels of GhPP2CA3, 11, 13, 27 were high in fibers. The transcripts of GhPP2CA4, 16, 22 were abundant in stems. These results imply that most cotton PP2CA members may function in reproductive development, and some PP2CAs also play roles in some specific tissues like roots, fibers, and stems (Fig. 6).

Transcriptional changes of GhPP2CAs in responses to ABA and osmotic stress

To gain insight into the roles of GhPP2CAs in ABA signaling, transcriptional abundances of GhPP2CAs in roots were detected after treatments with 100 μM ABA or 10% PEG6000 for indicated periods of time. We observed that the transcriptional levels of some GhPP2CA genes such as GhPP2CA5, 11, 18, 20, 25, 27 continually increased with the extension of ABA treatment time. In contrast, the expression levels of some members including GhPP2CA2–4, 8, 10, 24 had decreasing trends. The expression levels of some genes were decreased at 3 or 6 h but increased at 12 or 24 h. These genes included GhPP2CA1, 16, 17, 19, 26. The expression levels of several genes increased at 3 or 6 h but decreased at 12 or 24 h. These genes were GhPP2CA6, 7, 9, 12–15, 21–23 (Fig. 7). Treatment of cotton seedlings with PEG6000 also altered the expression of most GhPP2CA genes (Fig. 8). The majority of GhPP2CAs were upregulated after treatments with PEG for a short time period and downregulated afterward. For example, the transcriptional levels of GhPP2CA8 and GhPP2CA21 were prominently enhanced at 3 h, and then reduced at 6, 12, and 24 h, while those of GhPP2CA5, 6, 17, 23, 25, 26 were pronouncedly increased at 6 h and decreased at 12 and 24 h. By contrast, the expression of some genes was significantly elevated at 12 or 24 h post PEG treatment. These genes included GhPP2CA1, 3, 4, 7, 9, 12–14, 18, 20, 24, 27. The overall trend of expression of GhPP2CA10 and GhPP2CA15 was increased while that of GhPP2CA2, 11, 16 was decreased under osmotic stress (Fig. 8). Together, these data suggest that GhPP2CAs exhibit diverse expression patterns in responses to ABA and osmotic stress.

Many GhPP2CAs interact with GhPYL2-2D, GhPYL6-2A, and GhPYL9-2A

Clade A type 2C protein phosphatases have been documented to interact with ABA receptor PYLs in ABA signal pathway (Ma et al., 2009; Park et al., 2009). Accordingly, we investigated the interactions between GhPP2CAs and GhYPLs in the absence or presence of ABA by yeast-two hybrid method. A total of 11 GhPP2CAs were cloned, and three GhPYLs GhPYL2-2D (Gh_D08G2587), GhPYL6-2A (Gh_A06G1418), and GhPYL9-2A (Gh_A11G0870) were randomly selected and cloned. These genes were fused into yeast vectors, and yeast-two hybrid experiments were performed. In the absence of ABA, GhPP2CA2, and GhPP2CA25, respectively, interacted with GhPYL2-2D while multiple GhPP2CAs like GhPP2CA2, 3, 6, 10, 13, 15, 18, 19, 25 individually interplayed with GhPYL2-2D in the presence of ABA (Fig. 9). In contrast, several GhPP2CAs could, respectively, interact with GhPYL6-2A or GhPYL9-2A either with or without ABA. These GhPP2CAs included GhPP2CA2, 6, 10, 13, 15, 18, 19, 24, 25. Besides, GhPP2CA3 interact with GhPYL9-2A but not with GhPYL6-2A either in the presence or absence of ABA (Figs.10 and 11). These results imply that GhPP2CAs differentially interact with GhPYLs in responding to ABA in cotton.