Investigation of the AQP Family in Soybean and the Promoter Activity of TIP2;6 in Heat Stress and Hormone Responses

Aquaporins (AQPs) are one diverse family of membrane channel proteins that play crucial regulatory roles in plant stress physiology. However, the heat stress responsiveness of AQP genes in soybean remains poorly understood. In this study, 75 non-redundant AQP encoding genes were identified in soybean. Multiple sequence alignments showed that all GmAQP proteins possessed the conserved regions, which contained 6 trans-membrane domains (TM1 to TM6). Different GmAQP members consisted of distinct Asn-Pro-Ala (NPA) motifs, aromatic/arginine (ar/R) selectivity filters and Froger’s positions (FPs). Phylogenetic analyses distinguished five sub-families within these GmAQPs: 24 GmPIPs, 24 GmTIPs, 17 GmNIPs, 8 GmSIPs, and 2 GmXIPs. Promoter cis-acting elements analyses revealed that distinct number and composition of heat stress and hormone responsive elements existed in different promoter regions of GmAQPs. QRT-PCR assays demonstrated that 12 candidate GmAQPs with relatively extensive expression in various tissues or high expression levels in root or leaf exhibited different expression changes under heat stress and hormone cues (abscisic acid (ABA), l-aminocyclopropane-l-carboxylic acid (ACC), salicylic acid (SA) and methyl jasmonate (MeJA)). Furthermore, the promoter activity of one previously functionally unknown AQP gene-GmTIP2;6 was investigated in transgenic Arabidopsis plants. The beta-glucuronidase (GUS) activity driven by the promoter of GmTIP2;6 was strongly induced in the heat- and ACC-treated transgenic plants and tended to be accumulated in the hypocotyls, vascular bundles, and leaf trichomes. These results will contribute to uncovering the potential functions and molecular mechanisms of soybean GmAQPs in mediating heat stress and hormone signal responses.


Key Structural Features of the AQP Proteins
To understand the possible physiological role and substrate specificity of soybean AQP proteins, the TM domains, NPA motifs, ar/R selectivity filters, and FPs were investigated (Table 2; Figures S1-S6). Protein structure analyses supported that all GmAQP proteins possessed the typically conserved regions, which contained 6 TM domains (TM1 to TM6) ( Figure S1). All GmPIPs and GmTIPs contained two conserved NPA motifs in LB and LE. In GmNIPs, the first NPA showed the same sequence as in PIPs and TIPs, except for GmNIP5;2, where A was replaced by S. The second NPA motif showed an A to V substitution in four NIPs (GmNIP1;6, GmNIP5;2, GmNIP6;2, and GmNIP6;3). GmSIPs showed a second NPA motif completely conserved with the other members. Instead, all of the first NPA motifs showed the replacement of A by T (GmSIP1;1, GmSIP1;2, GmSIP1;3, and GmSIP1;4), A by S (GmSIP1;5 and GmSIP1;6), and A by L (GmSIP2;1 and GmSIP1;2). In GmXIPs, A in the first NPA of GmXIP1;1 was changed to I, and N in the second NPA was changed to S. The second NPA of GmXIP1;2 was completely conserved, while the N and A residues in the first NPA were replaced by S and V. Table 2. Conserved amino acid residues (Asn-Pro-Ala (NPA) motifs, aromatic/arginine (ar/R) filters and Froger's positions (FPs)) and trans-membrane (TM) domains of AQP proteins in soybean.

Gene Name
Gene Symbol TM Number

NPA Motifs ar/R Selectivity Filters FPs
The ar/R positions (H2, H5, LE1, and LE2) of GmAQPs showed increased sub-family specificity compared to the two NPA motifs. In GmPIPs, all selectivity filters were F-H-T-R. In GmTIPs, the ar/

Evolutionary Characterization of the AQP Genes
To investigate the classification and evolutionary relationship of soybean AQP proteins, the phylogenetic tree was constructed with the full-length GmAQP protein sequences from Arabidopsis, Phaseolus vulgaris, Populus trichocarpa, and Lotus japonicus (Dataset S1; Dataset S2). Soybean AQPs grouped into five sub-families (PIPs, 24; TIPs, 24; NIPs, 17; SIPs, 8; XIPs, 2) ( Figure 2). Among the PIP sub-family, 24 members were divided into two groups: PIP1 with 10 members and PIP2 with 14 members. Five groups were found for the TIP sub-family (TIP1 to TIP5), with 10 members in the TIP1 group, 7 members in the TIP2 group, 4 members in the TIP3 group, 2 members in the TIP4 group, and 1 member in the TIP5 group. Seven groups belonged to the NIP sub-family (NIP1 to NIP7), with 6 members in the NIP1 group, 1 member in the NIP2 group, 1 member in the NIP3 group, 2 members in the NIP4 group, 2 members in the NIP5 group, 3 members in the NIP6 group, and 2 members in the NIP7 group. The SIP sub-family were composed of two groups (SIP1 and SIP2), with 6 members in the SIP1 group and 2 members in the SIP2 group. Only one group was identified for the XIP sub-family (XIP1), with 2 members.
tissues, heat stress obviously up-regulated the expression of GmAQPs during the early durations (1.5 hour) in roots, whereas heat stress slightly up-regulated expression of GmAQPs during the late durations (6.0 hour) in leaves (Figure 4; Table S1). In roots, most of the analyzed members (GmPIP1;7, GmPIP1; 8,GmPIP2;4,GmPIP2;5,GmPIP2;13,GmPIP2;14,GmTIP1;7,GmTIP2;2,and GmTIP2;6) were transcriptionally up-regulated, whereas GmTIP4;1 and GmSIP1;3 were extremely down-regulated after heat treatment. In leaves, the transcripts of GmPIP2;5, GmPIP2;6, GmTIP1;7, and GmTIP4;1 were inhibited by heat stress, while GmPIP1;7, GmPIP1;8, GmPIP2;4, GmPIP2;13, GmPIP2;14, GmTIP2;2, and GmSIP1;3 were firstly inhibited and then promoted. Amongst them, GmPIP1;8, GmPIP2;13, GmPIP2;14, and GmTIP2;6 were relatively dramatically induced in roots, but GmTIP4;1 was sharply repressed in both the roots and leaves under heat stress. Further, these candidate soybean AQP genes were subjected to qRT-PCR analyses to evaluate their roles in response to hormone signals ( Figure 5; Table S1). For ABA treatment, most transcripts of PIP genes underwent down-regulation after ABA treatment, except GmPIP1;8 and GmPIP2;13, whereas TIP genes shared up-regulation in both the roots and leaves. GmPIP1;8 was strongly up-regulated after 0.5 hour of ABA stress in roots and down-regulated in leaves. GmPIP2;13 was abundantly down-regulated after 0.5 hour of ABA stress in roots and up-regulated in leaves. For ACC treatment, most transcripts of PIP genes showed an increase in roots after 0.5 hour, 1.5 hour, and 6.0 hour of ACC stress, except GmPIP2;13. However, in leaves, a subset of PIP genes (GmPIP1;7, GmPIP1;8, GmPIP2;4, GmPIP2;5, GmPIP2;6, and GmPIP2;14) initially displayed up-regulation expression following gradual down-regulation expression. GmPIP2;13 was up-regulated only during early duration (0.5 hour) and then down-regulated in roots and continuously up-regulated in leaves. For TIP genes, all of them exhibited up-regulated expression after ACC treatment. GmSIP1;3 was up-regulated in roots and down-regulated in leaves. For SA treatment, a cluster of PIP genes (GmPIP1;7, GmPIP1;8, GmPIP2;4, and GmPIP2;5) presented down-regulation in both the roots and leaves. GmPIP2;6 was up-regulated in both the roots and leaves during early durations (0.5 and 1.5 hour) following gradual down-regulation. GmPIP2;13 was up-regulated after 0.5 hour of SA stress in roots and down-regulated in leaves. GmPIP2;14 displayed obvious down-regulation in roots and up-regulation in leaves after 0.5 hour of SA stress following down-regulation. For TIP genes, the transcriptional level of GmTIP1;7 was detected with up-regulation in both the roots and leaves. GmTIP2;2 showed sharp down-regulation after 1.5 hour of SA stress in the roots. GmTIP4;1 was up-regulated at 0.5 hour of SA stress and then down-regulated in both the roots and leaves. GmSIP1;3 displayed high transcript abundance in roots and low transcript abundance in leaves. For MeJA treatment, the expression levels of PIP genes (GmPIP1;7, GmPIP1;8, GmPIP2;4, GmPIP2;5, and GmPIP2;14) were promoted in roots and inhibited in leaves, while GmPIP2;13 was inhibited in roots and promoted in leaves. For TIP genes, GmTIP1;7 maintained up-regulation in both the roots and leaves. GmTIP2;2 showed up-regulation in roots and down-regulation in leaves. GmTIP4;1 presented up-regulation during early duration (0.5 hour) following gradual down-regulation in both the roots and leaves. The transcript of GmSIP1;3 maintained a high level in both the roots and leaves.

GUS Activity of the GmTIP2;6 Promoter
To characterize the function of AQP promoter in response to heat stress and hormone signals, the GmTIP2;6 promoter, with a relatively high number of HSE elements, was fused to the GUS reporter gene and transferred into Arabidopsis (Dataset S3). Under normal growth conditions, the expression pattern of the GUS gene driven by the GmTIP2;6 promoter was weakly detected in the hypocotyls ( Figure 7A). After heat treatment, GUS activity was remarkably induced and increased in hypocotyls, roots, leaf vascular bundles, and young leaf trichomes. Similarly, after ACC treatment, the responsiveness of the GmTIP2;6 promoter was enhanced in hypocotyls, roots, leaf vascular bundles, and young leaf trichomes ( Figure 7A). Further quantitative GUS activity analyses verified that the translational levels of GUS protein in the transgenic plants with heat and ACC treatments were evidently stronger than those without treatments ( Figure 7B). This result confirmed that GmTIP2;6 was one heat stress and ACC hormone inducible promoter.

GUS Activity of the GmTIP2;6 Promoter
To characterize the function of AQP promoter in response to heat stress and hormone signals, the GmTIP2;6 promoter, with a relatively high number of HSE elements, was fused to the GUS reporter gene and transferred into Arabidopsis (Dataset S3). Under normal growth conditions, the expression pattern of the GUS gene driven by the GmTIP2;6 promoter was weakly detected in the hypocotyls ( Figure 7A). After heat treatment, GUS activity was remarkably induced and increased in hypocotyls, roots, leaf vascular bundles, and young leaf trichomes. Similarly, after ACC treatment, the responsiveness of the GmTIP2;6 promoter was enhanced in hypocotyls, roots, leaf vascular bundles, and young leaf trichomes ( Figure 7A). Further quantitative GUS activity analyses verified that the translational levels of GUS protein in the transgenic plants with heat and ACC treatments were evidently stronger than those without treatments ( Figure 7B). This result confirmed that GmTIP2;6 was one heat stress and ACC hormone inducible promoter. hypocotyls, roots, leaf vascular bundles, and young leaf trichomes. Similarly, after ACC treatment, the responsiveness of the GmTIP2;6 promoter was enhanced in hypocotyls, roots, leaf vascular bundles, and young leaf trichomes ( Figure 7A). Further quantitative GUS activity analyses verified that the translational levels of GUS protein in the transgenic plants with heat and ACC treatments were evidently stronger than those without treatments ( Figure 7B). This result confirmed that GmTIP2;6 was one heat stress and ACC hormone inducible promoter.
Plant hormones were important signal molecules that controlled plant growth and development in response to heat stimulus, including ABA, SA, and MeJA [44]. However, evidence for the molecular mechanism of AQPs' involvement in the hormone response process remained scanty. The preliminary promoter element analyses indicated that different members of soybean GmAQPs possessed distinct hormone-related elements ( Figure S7). For instance, two GmAQPs (GmPIP1;1 and GmNIP4;2), three GmAQPs (GmPIP2;9, GmTIP2;6, and GmNIP1;2), eight GmAQPs (GmPIP1;8, GmTIP2;7, GmTIP4;1, GmNIP1;1, GmNIP1;6, GmNIP6;3, GmSIP1;1, and GmSIP2;1) and three GmAQPs (GmTIP5;1, GmNIP7;1, and GmSIP1;4) contained ABA, ET, SA, or MeJA-special element, respectively. In contrast, a combination of four hormone-related elements were observed in GmPIP1;10, GmPIP2;1, and GmTIP4;2. Furthermore, gene expression analyses showed that different GmAQPs displayed distinct transcriptional changes, up-regulation or down-regulation under different hormone treatments, as evidenced by the qRT-PCR assay ( Figure 5). For example, five highly up-regulated GmPIPs transcripts during ACC and MeJA hormone treatments in roots included GmPIP1;7, GmPIP1;8, GmPIP2;4, GmPIP2;5, and GmPIP2;14, and four abundantly accumulated GmTIPs during ABA and ACC hormone treatments in both the roots and leaves included GmTIP1;7, GmTIP2;2, GmTIP2;6, and GmTIP4;1. Moreover, the expression of GmSIP1;3 got enhanced and underwent significant changes during ACC, SA, and MeJA hormone treatments in roots. However, the transcript changes of four GmPIPs (GmPIP1;7, GmPIP1;8, GmPIP2;4, and GmPIP2;5) were observed with down-regulation under ABA, SA, and MeJA hormone treatments in leaves. In rice and oilseed rape, ABA highly enhanced the expression of OsTIP1;1 in the shoots and roots and BnTIP2 in the seeds [45,46]. In resurrection plant, ABA greatly decreased the expression of CpTIP in the callus [47]. In wheat, ethylene (ET) up-regulated the gene expression of wheat TaAQP8 (one PIP sub-family gene) under salt stress [48]. In rose, ET decreased the expression of RhTIP1;1 in the flower [49]. Systematical microarray analysis of AQPs under ABA, ACC, and MeJA based on the public database of Arabidopsis eFP Browser showed that most Arabidopsis AQPs were involved in the process ( Figure S9; Table S1) [16,43]. Other stress-related signals such as brassinolide (BR), gibberellin (GA), and auxin (IAA) also regulated the expression of AQPs [50][51][52][53] but by, as yet, unclear mechanisms. All this evidence gives clues to the role of AQPs in complicated hormone signal transduction systems in different tissues. In this regard, it will be of interest to decipher further how AQP genes work in the interactive signaling regulation network, especially under heat stress.
Promoters, the direct indication of gene expression patterns, were extensively involved in the responses of signal molecules and environmental elicitors. It was noteworthy that we validated that soybean GmTIP2;6 promoter was a typical inducible promoter, which responded to ACC and heat stresses in hypocotyls, vascular bundles, and leaf trichomes (Figure 7). This was consistent with the qRT-PCR result that GmTIP2;6 was up-regulated by heat stress and ACC hormone (Figures 4  and 5). In cotton, strong expression of the GUS gene driven by GhPIP2;7 promoter was detected in leaves of 5 to 10-day-old transgenic Arabidopsis seedlings, but GUS activity gradually became weak as the seedlings further developed. GUS activity driven by cotton GhPIP2;7 promoter was remarkably increased after mannitol treatment [54]. In Arabidopsis, AtNIP3;1 promoter-mediated GUS activity was specifically expressed in the roots [55]. The expression of the GUS gene driven by AtPIP2;7 promoter was strongly detected in cotyledons, emerging leaf primordia, and root elongation zones, and salt stress induced strong repression of AtPIP2;7 promoter activity [35]. In soybean, the GUS activity of GmTIP2;3 promoter was expressed in the root, stem, and leaf and preferentially expressed in the stele of root and stem [56]. These data indicated different promoters of AQP members played different roles in different tissues or development stages. In the continued study, it will be very meaningful to investigate the core elements of AQP promoters for quantitative and qualitative gene expression regulation.

Categorization of Soybean AQP Genes
The soybean genome sequences were retrieved from Phytozome V12.1 (http://phytozome.jgi.doe. gov/pz/portal.html). The keyword searches of aquaporin, Hidden Markov model profile (PF00230), and KEGG Orthology terms (PIPs, K09872; TIPs, K09873; NIPs, K09874; SIPs, K09875) were applied to identify the soybean candidate AQP members [57]. Known Arabidopsis, Zea mays, Oryza sativa, Populus trichocarapa, Phaseolus vulgaris and Lotus japonicas AQPs were also subjected to BlastP and BlastN against the soybean database with cut-off E-value of e −5 [58][59][60][61]. GmAQPs were named based on their sequence homology with known AQPs and soybean genome annotation. The decrease redundancy tool (http://web.expasy.org/decrease_redundancy/) was utilized to discard the redundant AQP sequences. Further, the resulting candidate sequences were checked for the presence of six TM domains and two NPA motifs by SMART (http://smart.embl-heidelberg.de/smart/batch.pl) and NCBI-CDD (http: //www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi) web servers. The numbers of phosphorylation sites of AQP proteins were predicted with NetPhos 2.0 (http://www.cbs.dtu.dk/services/NetPhos/). The molecular weight (MW) and isoelectric point (pI) of AQP proteins were calculated using ExPASy (http://web.expasy.org). The chromosomal positions of GmAQPs were mapped using MapInspect software based on the starting position of all genes on each chromosome. Tandem duplications were identified manually. Adjacent genes of the same sub-group tightly linked within 20 kb of each other and the identity of the genes ≥80% are considered as tandem duplicated genes [62].

Phylogenetic Tree
The GmAQP full length protein sequences were aligned using ClustalX2 software. As soon as the ALN file was generated, MEGA6 was carried out to construct the neighbor-joining (NJ) phylogenetic tree [63]. The criteria were adopted with pairwise deletion option and Poisson correction model. Bootstrap test was performed with 1000 replicates.

Heat Stress and Hormone Treatments
Soybean cultivar GMLN012012017, with the characteristic of heat tolerance, was used in this study. Soybean seeds were cultivated in pots in an illuminated incubator (PTC-300, Shanghai, China) adjusted to 22 • C temperature, 60% relative humidity, 16/8 h photoperiod, and25000 Lux light intensity. For high-temperature treatment, 21-day-old seedlings in pots were transferred to the illuminated incubator adjusted to 42 • C. The un-treated samples were used as the control (0 hour). For hormone treatments, the root systems of 21-day-old seedlings were washed gently with water to remove soil and then the plants were soaked into 200 mL solutions with 100 µM abscisic acid (ABA), 100 µM l-aminocyclopropane-l-carboxylic acid (ACC), 100 µM salicylic acid (SA) or 100 µM methyl jasmonate (MeJA). The samples soaked with water were used as the control (0 hour). Each single seedling was sampled at one time point (0, 0.5, 1.5, 6.0 or 12 hour). Then, whole leaves and roots of each control or treated seedling were collected and frozen in liquid nitrogen and stored at −80 • C for further analyses.

qRT-PCR
Total RNA was separately extracted from the frozen samples using a RNA Simple Total RNA Extraction Kit (Tiangen, Beijing, China) according to the manufacturer's protocol. Then, the cDNA was synthesized using a FastQuant RT Kit (Tiangen, Beijing, China). Gene-specific primers were designed using PrimerQuest Tool (http://sg.idtdna.com/PrimerQuest/Home/) (Table S2). GmActin11 (Glyma.18G290800) was selected as the internal reference gene [65][66][67][68]. The amplification reactions were performed on Applied Biosystems StepOnePlus TM Real-Time System using KAPA SYBR ® Fast qPCR Kit (Tiangen, Beijing, China) with the following parameters: initializing denaturation at 95 • C for 5 min, followed by 45 cycles of denaturation at 95 • C for 5 seconds, annealing at 58 • C for 5 seconds, and extension at 72 • C for 30 seconds. Three technical replicates were maintained for each sample. The relative expression levels were calculated as 2 −∆∆Ct [69]. The heatmaps for the expression profiles of GmAQP genes were generated with BAR HeatMapper Tool [16].

Promoter Element Prediction
1.5 kb promoter regions, upstream of the AQP gene coding sequences, were extracted from Phytozome V12.1. Promoter regions were subsequently analyzed using PlantCARE database (http: //bioinformatics.psb.ugent.be/webtools/plantcare/html/) to illustrate the number and composition of hormone and heat stress responsive elements.

Promoter Cloning and Arabidopsis Transformation
The promoter of GmTIP2;6 was isolated using the specific primer pairs (Table S2). To generate the proGmTIP2;6::GUS construct, the CaMV 35S promoter was replaced by the promoter of GmTIP2;6. The promoter of GmTIP2;6 was inserted into Sal I/Sma I sites and sub-cloned into the pBI121 vector whose Hind III site was replaced by three continuous sites (Hind III, Pst I, and Sal I) upstream of the CaMV 35S promoter ( Figure S10). The GUS fusion construct was then introduced into Arabidopsis (Col-0) by Agrobacterium-mediated floral-dip method [70,71]. Transformed seeds were selected on MS medium with 50 mg/L kanamycin (Kan). Homozygous lines of T3 were used for the following GUS activity assays.

GUS Activity Detection
To evaluate GUS activity, the proGmTIP2;6::GUS transgenic seeds were sowed on MS medium for 5 days, and then transferred to MS medium supplemented with 100 µM ACC for 7 days. For the high-temperature treatment, 5-day-old seedlings on MS medium were transferred to the temperature-controlled chamber adjusted to 37 • C for 7 days. Seedlings on MS medium without any additions were used as controls [72]. GUS staining was conducted using standard protocols. In brief, seedlings were incubated in 10 mL tubes with 1 mg/mL 5-bromo-4-chloro-3-indolylglucuronic acid, 5 mM potassium ferricyanide, 5 mM potassium ferrocyanide, 0.03% Triton X-100, and 0.1 M sodium phosphate buffer, pH 7.0 overnight at 37 • C. Then, seedlings were immersed in 70% ethanol to remove the chlorophyll and visualized on the Leica microscope. GUS activities were measured by monitoring the cleavage of GUS substrate 4-methylumbelliferyl glucuronide as reported previously [73]. Data analyses of variance were used to compare the statistical difference based on Student's t-test, at a significant level of p < 0.01.

Supplementary Materials:
The following are available online at http://www.mdpi.com/1422-0067/20/2/262/s1. Dataset S1. Fasta files of 75 AQP proteins in soybean. Dataset S2. Fasta files of AQP proteins in Arabidopsis, Phaseolus vulgaris, Populus trichocarpa and Lotus japonicus. Dataset S3. Fasta file of the promoter for GmTIP2;6 in soybean. Table S1. Gene expression numbers for GmAQPs and AtAQPs in response to heat stress and hormone signals. Table S2. List of primers used in this study. Figure S1. Predicted structures of 75 AQP proteins in soybean. Figure S2. Conserved amino acid residues (NPA motifs, ar/R filter, FPs) and TM domains of PIPs in soybean. Figure S3. Conserved amino acid residues (NPA motifs, ar/R filter, FPs) and TM domains of TIPs in soybean. Figure S4. Conserved amino acid residues (NPA motifs, ar/R filter, FPs) and TM domains of NIPs in soybean. Figure S5. Conserved amino acid residues (NPA motifs, ar/R filter, FPs) and TM domains of SIPs in soybean. Figure S6. Conserved amino acid residues (NPA motifs, ar/R filter, FPs) and TM domains of XIPs in soybean. Figure S7. Distribution of cis-acting elements in 75 soybean AQP gene promoters. 1500 bp adjacent to the AQP coding sequence. The elements are represented by different colors. The scale bar represents 500 bp. Figure S8. Expression profiles of Arabidopsis AQP genes in response to heat stress. Differential sub-family AQP gene expression in response to heat stress across eight time points (0.25, 0.5, 1, 3, 4, 6, 12 and 24 hour) in two tissues (roots and shoots). Differences in gene expression changes are shown in color in the green-red scale. Figure S9. Expression profiles of Arabidopsis AQP genes in response to ABA, ACC, and MeJA signals. Differential sub-family AQP gene expression in response to ABA, ACC, and MeJA hormones across three time points (0.5, 1, and 3 hour). Differences in gene expression changes were shown in color in the green-red scale. Figure S10.