Genome-wide identification and expression analysis of the IQD gene family in moso bamboo (Phyllostachys edulis)

Members of the plant-specific IQ67-domain (IQD) protein family are involved in various aspects of normal plant growth and developmental processes as well as basal defence response. Although hundreds of IQD proteins have been identified, only a small number of IQDs have been functionally characterized. Moreover, no systematic study has been performed on moso bamboo. In this study, we performed for the first time a genome-wide identification and expression analysis of the IQD gene family in moso bamboo. We identified 29 non-redundant PeIQD encoding genes. Analysis of the evolutionary patterns and divergence revealed that the IQD genes underwent a large-scale event around 12 million years ago and the division times of IQD family genes between moso bamboo and rice, and, between moso bamboo and Brachypodium, were found to be 20–35 MYA and 25–40 MYA, respectively. We surveyed the putative promoter regions of the PeIQD genes, which showed that largely stress-related cis-elements existed in these genes. The expression profiles of the IQD genes shed light on their functional divergence. Additionally, a yeast two-hybrid assay proved that PeIQD8 can interact with PeCaM2 and that IQ or I in the IQ motif is required for PeIQD8 to combine with CaM2.


Results
Identification of IQD gene family in moso bamboo. In order to conduct a genome-wide identification of the IQD gene families in moso bamboo, the HHM profile of IQ domain (PF00612) was exploited as query to identify the IQD genes in the moso bamboo genome (http://www.ncgr.ac.cn/bamboo, accessed February 2016). We identified 35 putative IQD genes, then these genes were conducted to confirm the existence of the conserved IQD domain with InterproScan. After removing redundant sequences, 29 IQD genes were identified, and all of these predicted that IQD proteins had a typical "IQ calmodulin-binding motif " domain (PF00612) which is a major calcium (Ca 2+ ) regulator. We named these 29 IQD genes -IQD1 to IQD29, according to their physical locations (from top to bottom) on chromosomes ( Table 1). The identified IQD genes in moso bamboo encode proteins ranging from 190 to 940 amino acids (aa) in length with an average of 486 aa. The other characteristics of the IQD genes -including isoelectric point (pI), molecular weight (MW) chromosome location, and exons -are listed in Table 1. The open reading frame (ORF) length ranged from 573 bp (IQD6) to 2823 bp (IQD17), the MW from 20417.05 (IQD6) to 105765.41 (IQD17) Da, and the pI from 5.02 (IQD5) to 11.12 (IQD9).
Phylogenetic analysis of the IQD gene family. To examine the phylogenetic relationship of the IQD domain proteins among the three grass subfamilies -including Pooideae (Brachypodium distachyon), Ehrhartoideae (rice), Bambusoideae (Phyllostachys edulis) and four dicots (Arabidopsis thaliana, soybean, tomato and the woody plant poplar), a rooted tree was constructed from alignments of the full-length IQD protein sequences (Fig. 1A). The phylogenetic tree was constructed using MEGA 6.0 software by employing the Nneighbor-joining (NJ), minimal evolution and maximum parsimony methods. The tree topologies produced by the three algorithms were largely comparable with only minor modifications at interior branches (data not shown). Therefore, only the NJ phylogenetic tree was further analyzed in our study. The NJ phylogenetic tree comprised 226 full-length IQD protein sequences from Oryza sativa (28), Brachypodium distachyon (23), Arabidopsis thaliana (34), Populus trichocarpa (40), Glycine max (67) and Solanum lycopersicum (34) 6,12,13,16,17 . The characteristics of the IQD genes-including pI, MW, chromosome location, ORF length and amino acids -are showed in Table S1. As indicated in Fig. 1A, the phylogenetic tree highlighted that the IQD family of proteins could be divided into four well-conserved subfamilies, IQD I-VI, as described previously and with significant statistical support 6,12,13 . We further examined each of the subfamilies containing the IQD genes. The IQDI subfamily Scientific RepoRts | 6:24520 | DOI: 10.1038/srep24520 was divided into four clades designated as clade a, b, c and d, which was consistent with the nomenclature in previous studies of Arabidopsis, rice and soybean 6,12,13 . The other three subfamilies were separated into two clades (a and b) (Fig. 1A). The phylogenetic tree revealed that the plant IQD sequence distribution predominates with species bias (Fig. 1B). The largest subfamilies of the plant species were found to contain mostly IQD I genes . In contrast, the IQD II subfamily contained the fewest IQD genes -except in the case of Brachypodium and tomato, where IQD IV contained the fewest. It also appeared that more IQD IV genes were present in moso bamboo and soybean IQD IV than in other species. In fact, we found nine and ten IQD IV genes to be present in moso bamboo and soybean, respectively, while only three and two IQD IV genes were present in Brachypodium and tomato, respectively. The percentage distribution of the IQD proteins among each subfamily was calculated for all seven species (Fig. 1C). Remarkably, there were fewer IQD genes in each subfamily for monocotyledonous than for dicotyledonous. Gene structure and conserved motifs of moso bamboo IQD genes. To gain further insights into the structural diversity of moso bamboo IQD genes, we first constructed a separate phylogenetic tree exclusively using the full-length IQD protein sequences of moso bamboo. Moso bamboo proteins were also classified into four independent subfamilies, in good agreement with that described above for the seven plant species (Figs 1A and 2A). It is well known that genetic structural diversity is a possible mechanism for the evolution of multigene families. To gain further insights into the structural diversity of IQD genes, we compared the exon/intron organization in the coding sequences of individual IQD genes in moso bamboo (Fig. 2B). It is interesting to note that we found the intron/exon structure of most sister gene pairs to be conserved, but there are also some differences. For example, the six sister gene pairs (PeIQD24/-14, PeIQD19/-23, PeIQD10/-12, PeIQD11/-16, PeIQD20/-25 and PeIQD2/-29) were found to have the same exon/intron number and intron phase. However, their intron lengths showed great variability, ranging from a few tens of bp to − 2 kb. In addition, four sister gene pairs showed greater changes in their structural organizations (PeIQD4/-9, PeIQD18/-21, PeIQD26/-5, PeIQD22/-6) and varying numbers of exons and introns.
We further searched for the conserved motifs in PeIQD proteins by MEME program to obtain more insights into the diversity of motif compositions among PeIQDs. As shown in Fig. 3, a total of 10 conserved motifs designated as motif 1 to motif 10 were identified. The details of the 10 putative motifs were referred in Table S2. Motif 1 was found to encode the IQ domain, while motif 2 encoded proteins of unknown function (DUF4005). The specific motifs of the other subfamilies have not been functionally annotated. As expected, most of the closely related members had common motif compositions, suggesting functional similarities among IQD proteins within the same subfamily (Fig. 3). The most common motif is motif 1, found in all 29 moso bamboo IQD genes (Fig. 3).
To predict calmodulin-binding sites, we searched the Calmodulin Target Database, which provides various structural and biophysical parameters for the 29 PeIQD protein sequences. Consecutive strings of amino acid residues with high values of ≥ 7 indicated the locations of the putative binding sites. The search results showed that all of the predicted PeIQD proteins contained at least one CaM-binding site. These CaM-binding sites are described in Table 2.

Evolutionary patterns and divergence of the IQD gene family in moso bamboo, rice and
Brachypodium. In comparative genomics, the phylogeny-based and bidirectional best-hit methods are popular strategies for identifying possible paralogous or orthologous genes. Using these two methods, we found nine putative paralogous pairs in the moso bamboo genome, five orthologous pairs between OsIQD and PeIQD, and two orthologous pairs between BdIQD and PeIQD. All gene-pairs are listed in Table S3. To evaluate the divergence times among these three monocotyledon and gramineous plants, we used a relative Ks measure as a proxy for time. Figure 4(A,B) shows the frequency distributions of the relative Ks values obtained from duplicated and paralogous gene-pairs in the moso bamboo, and from thorthologous pairs between moso bamboo and rice, and, between moso bamboo and Brachypodium genomes. The relative Ks distribution peaks around 0.15 in moso bamboo suggested a large-scale event around 12 million years ago (MYA). Peng et al. have found that bamboo underwent whole-genome duplication 7-12 MYA 14 , according to analyses of clustered gene families, though our results, indicated that IQD family genes underwent a longer large-scale event. Similarly, the relative Ks distribution peaked at 0.25 for the seven duplicated orthologous gene-pairs, indicating division within the three groups of IQD genes at 20 MYA. A previous study estimated that the divergence time of rice and moso bamboo was 48.6 MYA, and for Brachypodium and moso bamboo was 46.9 MYA , and the results of the study by Peng, et al. 14 -combined with those of our own study -revealed that the IQD family have undergone gene evolution after separation of the three progenitors.
To determine the Ka/Ks ratios of different loci in the coding sequences, we performed a sliding-window analysis of the 16 gene pairs ( Figure S2). The results showed that the IQD domains have undergone strongly positive selection  (Ka/Ks 1), while most loci and regions have undergone moderately or strongly negative selection, as predicted by the overall Ka/Ks ratio.

Protein interaction analysis of PeIQDs with PeCaM2.
Previous studies have reported that IQD1 and IQD20 can interact with CaM in Arabidopsis 6,7 . Recently, it was proposed that ZmIQD15 can interact with CaM in vitro 18 . To attest whether this interaction also occurs in moso bamboo, yeast two-hybrid analysis was adopted to examine the interaction between PeIQDs and PeCaM2. It was verified that expression of BD-PeIQDs and AD-PeCaM2 was not toxic for AH109 yeast cells and that these fusion proteins have no transcriptional activation ability when expressed separately ( Figure S3). Our results showed that all the transformants tested were found to grow well on SD/-Leu/-Trp medium. When transferred onto SD/-Trp/-Leu/-Ade/-His/X-α -Gal plates for 5 days, both the positive control and the experimental group turned blue. In contrast, the negative control displayed no α -galactosidase activity (Fig. 5). The results suggested that PeIQDs can interact with PeCaM2 in vitro.

Mutations in the IQ motif of PeIQD8 result in the loss of affinity CaM binding.
When the recombinant plasmid pGBKT7-PeIQD8 Del112-113 and pGBKT7-PeIQD8 I112T were transformed into the yeast AH109 strain, the self-activation was lost ( Figure S3). From Fig. 6, we can clearly see that only the pGBKT7-PeIQD8 and pGADT7-PeCaM2 were co-transformed into the AH109 yeast strain, the selective agar plates (SD/-Trp/-Le u/-Ade/-His/X-α -Gal) turned blue. However, the yeasts that were co-transformed with pGBKT7-PeIQD8 Del112-113 Scientific RepoRts | 6:24520 | DOI: 10.1038/srep24520 and pGADT7-PeCaM2 or pGBKT7-PeIQD8 I112T and pGADT7-PeCaM2 were not able to grow on the selective medium. The results suggested that the IQ or I in the IQ motif is required for PeIQD8 to combine with CaM2.
Analysis of the putative promoter regions of the IQD gene family. Cis-regulatory elements largely determine the tissue-specific or stress responsive expression patterns of genes 19 , and multi-stimulus responsive genes are closely correlated with cis-regulatory elements in the promoter regions 20,21 . The function of the IQD gene family has been studied for many species (Populus, soybean, maize, and others); such genes have been found to be mainly regulated by drought and MeJA stresses. Two types of cis-elements, including the dehydration-responsive element and the wound element were detected in the current study. This observation prompted us to investigate possible stress-responsive cis-elements in the promoter regions of the moso bamboo   Table 2. Predicted calmodulin binding sites in moso bamboo IQD proteins. Predicted calmodulin binding sites obtained from the Calmodulin Target Database are shown for strings of amino acid residues with a score of at least 7. Residues with the highest score (9) are highlighted in bold. Numbers before strings indicate the location of the first amino acid residues of the strings in moso bamboo IQD protein sequences. IQD genes by searching the PLACE database. Cis-elements located upstream of genes comprising up to 2000 bp have marked effects on binding to target genes 22,23 , so the 2000-bp putative promoter regions were used to identify putative stress-responsive cis-regulatory elements. The results (Table S5) showed that numerous abiotic stress cis-elements -S000176 and S000415 for drought stress, and S000457 for wound stress -were found widely in the promoter regions of the IQD genes in moso bamboo. This clearly showed that IQD gene family respond to abiotic stress and have potential functions in enhancing abiotic stress resistance. For instance, PeIQD13 possessed up to 22 drought-stress elements (S000415) and PeIQD9 had up to15 wound-stress elements (S000457). Further analysis of the putative promoter regions of the PeIQD gene family will be important in advancing our understanding of the stress tolerance mechanism in moso bamboo.   From the heat map of the real-time quantitative PCR (qRT-PCR) analysis results for the PeIQD genes, we can see that all 29 genes were found to be MeJA responsive, but that some differences were observed among these genes (Fig. 8). Similarly, the expression patterns of 29 PeIQD genes under MeJA stress were revealed by qRT-PCR

Discussion
The plant-specific IQD gene family has been comprehensively analyzed in many plants, including Arabidopsis, Oryza sativa, Brachypodium, Populus, soybean and tomato, but no IQD genes have been found in algae, suggesting that IQD proteins belong to an ancient family of CaM/CML-binding proteins and developed during the early evolution of land plants. Predicted IQD-like gene in Physcomitrella patens showed that the IQD gene family evolved during the early evolution of land plants 450-700 Mya ago 24,25 . In our study, we identified and characterized 29 IQ67 domain-encoding genes in moso bamboo using genome wide analysis, and compared these with 34 AtIQDs, 28 OsIQDs, 23 BdIQDs, 40 PtIQDs, 67 GmIQDs and 34 SlIQDs, indicating that the number of IQD genes in moso bamboo (29) is higher than that in Oryza sativa and Brachypodium distachyon. but lower than that in dicotyledons, such as Arabidopsis thaliana, Populus trichocarpa, Glycine max and tomato 6,12,13,16 . In these plants, genome sizes vary enormously; for example, the genome size for moso bamboo is 2021 Mb 14 29 . There are some duplication events may help us to understand the higher number of IQD genes within the smaller genome size for Arabidopsis, soybean, Populus and so on. At least four different large-scale duplication events were shown to have occurred in the Arabidopsis genome 100-200 million years ago, and 17% of all genes were shaped in tandem arrays 30,31 . The genome duplications of soybean occurred at approximately 59 and 13 million years ago, resulting in a highly duplicated genome with nearly 75% of the genes present in multiple copies 28 . And the Populus genome has undergone at least three rounds of whole genome duplications in its evolutionary history 27 . The 226 IQD proteins were found to be divided into four subfamilies, with every subfamily containing both monocotyledons and dicotyledons members, indicating that IQD genes had diversified before the monocot-dicot split. Noticeably, IQD genes with the same functions showed a tendency to be grouped into one subfamily, similar to findings of the previous reports 6,12,13,22 . For instances, subfamily IQDIII encompassed the IQD proteins involved in drought and MeJA stresses (GmIQDIII and PtIQDIII) 12,13 . Moreover, we sought additional evidence -such as gene structure, motif compositions and expression patterns as described above -to support the reliability of the subfamily classification.
In order to research the occurrence of the whole-genome duplication event and the tetraploid origin of bamboo, we investigate the gene collinearity in moso bamboo, rice and Brachypodium, as the most recent whole-genome duplication was likely linked to polyploidy events 14 . Recent gene duplication events, which play an important role in the rapid expansion and evolution of gene families, result in many paralogous pairs in different species 32 . And large-scale duplication events are defined as simultaneous duplications of genes. To better explain the patterns of macroevolution, estimates of the evolutionary rates are extremely useful. We estimated the Ks and Ka models of paralogous genes (Pe-Pe) and orthologous genes (Pe-Os and Pe-Bd), than the Ks value was calculated for each gene-pair. Lastly, we use the formula (T = Ks/2λ ) to calculate the approximate date of the duplication event, assuming a λ value of synonymous substitution of 6.5 × 10 -9 for moso bamboo, rice and Brachypodium 14,33,34 . We estimated the large-scale event occurred approximately was 5-12 MYA in moso bamboo, and the divergence times for orthologous genes (Pe-Os and Pe-Bd) to be 20-35 MYA and 25-40 MYA . Peng et al. estimated that the divergence time for moso bamboo and rice was 48.6 MYA, and for moso bamboo and Brachypodium was 46.9 MYA 14 . According to the ratio of nonsynonymous to synonymous substitutions (Ka/Ks), the history of selection acting on coding sequences can be measured. In general, Ka/Ks ratio greater than 1, less than 1 and equal to 1 represents positive selection, negative or stabilizing selection and neutral selection, respectively 32,35 . In the current study, the Ka/Ks ratios were less than 1 in moso bamboo, indicating that the nine gene pairs experienced purifying selection during the process of evolution. The Ka/Ks ratios were less than 0.9 for Pe-Os and Pe-Bd, indicating that the seven gene pairs experienced negative or purifying selection during the process of evolution.
The characteristic of the IQ67 domain is remarkable: it is characterized by a unique and repetitive arrangement of three different calmodulin recruitment motifs, known as the IQ, 1-5-10, and 1-8-14 motifs. In addition, these motifs contain some basic and hydrophobic amino acid residues 6,36 . All predicted PeIQD proteins have a typical "IQ calmodulin-binding motif " domain (PF00612) which is a major calcium (Ca 2+ ) regulator (Fig. 3). Additional, PeIQD proteins are constructed by multiple sequence alignment of their IQ domains, which span Scientific RepoRts | 6:24520 | DOI: 10.1038/srep24520 approx 67 amino acids ( Figure S1). Previous studies have revealed that at least five Arabidopsis protein families play a role in the calcium signalling pathway because they contain common IQ motifs, such as the cyclic nucleotide gated channels (CNGC) family, the IQ-Motif (IQM) family, the CaM-binding transcriptional activator (CAMTA) family, the myosin family, and the IQD family [37][38][39] , which contain one, one, two, five and up to three IQ motifs, respectively. In our study, the PeIQD proteins contained three IQ motifs. Furthermore, the presence of possible CaM binding sites was predicted by using the Calmodulin Target Database for all the PeIQD proteins; all of the predicted PeIQD proteins were found to contain at least one CaM-binding site ( Table 2), suggesting that these proteins are a typical class of CaM targets. We also found that every PeIQD protein contains CaMBDs (CaM binding domains) of different primary structures, indicating that CaMBD is not a conserved domain. In the IQD family, some proteins have function in calcium signaling pathway were testified. AtIQD1 can control glucosinolate accumulation through binding to CaM in a Ca 2+ -dependent fashion under biotic stresses 8 . AtIQD20 also interacts with CaM in a Ca 2+ -independent manner 6 . AtIQD26 and AtCaM2 (AT2G27030) can be shown to interact by means of the yeast two-hybrid experiment 36 . It was noteworthy that the amino acid sequences for AtCaM2 (AT2G27030) and PH01000202G0880 were 97.32% similarity ( Figure S4). We used the InterproScan program to examine that PH01000202G0880 has the CaM domain and named it PeCaM2. Next, we used yeast two-hybrid analysis to detect the interaction between CaM and PeIQDs, and we randomly selected PeIQD8 to be tested. Our results showed that PeIQD8 can interact with PeCaM2, which provides further favourable evidence that PeIQDs may perform their functions by interacting with CaM. Zhou et al. proved that AtIQM1 (IQ motif-containing protein) can bind with CaM5 via its IQ-motif 40 . To investigate the important role of the IQ motif in the interaction between IQD proteins and CaMs, we designed two mutant proteins (PeIQD8 Del112-113 and PeIQD8 I112T ) to verify the binding characteristics of PeIQD8 with CaM2. Finally, we found that the two mutant proteins could not interact with PeCaM2. These results implied that the IQ motif is the key structural domain and that the I in the IQ motif is the key amino acid residue determining this binding activity. Furthermore, the interacting partners within IQD proteins promote the reconstruction of signalling pathways that involve IQD proteins.
In Arabidopsis, Jasmonic acid (JA) treatment leads to elevating levels of specific glucosinolate 41,42 , and glucosinolates are a small but diverse class of defense related secondary metabolites, which play an important roles in plant defence 43 . The overexpression of AtIQD1 causes the accumulation of glucosinolates and increases resistance against herbivory by augmenting and fine-tuning glucosinolate accumulation 7 . In addition, Feng et al. verify that 24 soybean IQD III genes are regulated by MeJA stress 12 . Similarly, the expression of 12 selected PeIQD members are regulated by MeJA. The IQD genes are also regulated by drought treatment. For example, the PtIQDIII members and ZmIQD genes have been testified 12,18 . Base on the phylogenetic analysis of IQD members, we used further qRT-PCR to investigate the expression patterns of PeIQD family genes under MeJA and drought treatments. The results demonstrated that the PeIQDs genes were either increased or repressed under the PEG and MeJA treatments, and we speculated that moso bamboo IQD genes might have the similar biological function in defence to insect herbivory and drought stress. By considering all of the results for abiotic stresses, we deduced that most genes within the same subclass of the phylogenetic tree showed different expression patterns. For example, three pairs of duplicated genes (PeQD5-PeIQD26, PeIQD18-PeIQD21 and PeIQD19-PeIQD23) showed the different expression patterns under drought (PEG) stress; similarly, six pairs of duplicated genes (PeIQD10-PeIQD12, PeIQD20-PeIQD25, PeIQD2-PeIQD29, PeIQD18-PeIQD21, PeIQD5-PeIQD26 and PeIQD6-PeIQD22) showed similar expression patterns under MeJA stress. This is because the regulatory sequences that responded to the stress conditions had diverged dramatically along with the evolution of each gene after duplication. In contrast, some pairs of duplicated genes expressed similar patterns, indicating that these duplicated genes in the same subclass share a high sequence similarity. Here, the expression patterns were similar, revealing that the regulatory sequences that responded to the stress conditions diverged to a much lesser extent the evolution of each gene after duplication. The new information obtained in this study may aid in the selection of appropriate candidate genes for further functional characterization.

Materials and Methods
Identification of IQD family genes in the moso bamboo genome. The conserved IQ domain (PF00612) was originally applied as a probe to search the National Center for Gene Research database (http:// www.ncgr.ac.cn/bamboo) 14 . All redundant sequences were discarded from further analysis based on cluster W alignment results 44 , sequence identification numbers and chromosome location. Furthermore, to verify the reliability of the initial results, all non-redundant candidate IQD sequences were analyzed to confirm the presence of the conserved IQ domain using the InterproScan program 45 Table S1. Moso bamboo IQD gene information, including the number of amino acids, ORF lengths and chromosome locations were obtained from the National Center for Gene Research (http://www.ncgr.ac.cn/bamboo). Physicochemical parameters including the molecular weight (MW) and isoelectric point (pI) of each gene product were calculated using compute the pI/Mw tool from ExPASy (http://www.expasy.org/tools/) and the parameter (resolution) was set to average 46 . We used the Gene Structures Display Server (GSDS) (http://gsds.cbi.pku.edu.cn/) to illustrate the exon/intron structure for individual IQD genes by comparing their cDNAs and the corresponding genomic DNA sequences 47 . Identification of the conserved motifs and putative calmodulin-binding sites. The conserved motifs were analyzed using MEME version 4.11.1 48,49 . The identified protein motifs were further annotated with ScanProsite 50 . All IQD protein sequences were examined against the Calmodulin Target Database (http://calcium. uhnres.utoronto.ca/ctdb/ctdb/home.html) to predict putative calmodulin-binding sites 51 .
Evolutionary patterns and divergence of the IQD gene family in moso bamboo, rice and Brachypodium analysis. Pairwise alignment of IQD gene encoding sequences of the orthologous and paralogous pairs was performed using ClustalX 2.11 software. Paralogous IQD gene pairs in moso bamboo were identified on the basis of alignment results. The criteria described in previous studies were adopted 52,53 : the shorter sequences covers over 70% of the longer sequence after alignment and the minimum identity of aligned regions is 70%. To identify putative orthologues between two species (A and B), each sequence from species A was searched against all sequences from species B using the BLASTN tool. Additionally, each sequence from species B was searched against all sequences from species A. The two sequences were defined as orthologues if each of them was the best hit of the other, and if more than 300 bp of the two sequences aligned. In addition, to further analyze gene duplication events, the synonymous substitution rate (Ks) and nonsynonymous substitution rate (Ka) were calculated using the software DnaSp 54,55 . The date of duplication events was subsequently estimated according to the equation T = Ks/2λ . The approximate value for clock-like rate was 6.5 synonymous substitutions per 10 −9 years for moso bamboo, rice and Brachypodium distachyon 14,33 . A sliding window analysis of Ka/Ks ratios was performed with the following parameters: window size, 150 bp; step size, 9 bp.
Yeast two-hybrid assay. The Matchmaker GAL4 two-hybrid system (Clontech, Palo Alto, CA) was used to test the interaction between PeIQDs and CaM. The full length of PeIQD8 and PeCaM2 (PH01000202G0880) cDNAs (Table S4) were cloned into the pGBKT7 bait vector and pGADT7 prey vector, respectively. Next, they were co-transformed into yeast strain AH109. The transformed yeast cells were selected on SD/-Trp/-Leu and SD/-Trp/-Leu/-Ade/-His/X-α -Gal plates at 30 °C for 3-5 days to determine the protein-protein interaction. We used pGBKT7-53 and pGADT7-T as a positive control, while the co-transformants with pGBKT7-Lam and pGADT7-T were used as negative controls in the yeast clone experiments.

Effects of mutations in the IQ motif of PeIQD8 on its CaM binding capacity.
To confirm the effects of mutations of key amino acid residues in the IQ motif of PeIQD8 on its CaM binding capacity, the IQ was deleted or the I was mutated to T using mutagenesis technology in vitro. Two pairs of primers were designed to amplify the IQD8 Del112-113 (deletion of IQ112-113 in the PeIQD8 ORF) and IQD8 I112T (where I112 was mutated to T112), using the sequencing-confirmed plasmid pEASY-T1 Simple-PeIQD8 as a template (Table S6). The PCR products were subsequently treated with Dpn I (TaKaRa, Japan) to eliminate methylation. After verification by sequencing, the IQD8 Del112-113 and IQD8 I112T fragments were inserted into KpnI and XbaI-digested pGBKT7 vector. Finally, the CaM binding capacity of the mutant proteins was analyzed via a yeast two-hybrid assay, as described above.
Analysis of the putative promoter regions of the IQD gene family in moso bamboo. The 2000-bp upstream sequences of the transcriptional start site of the PeIQDs were chosen to identify the cis-elements in the putative promoter regions. The PLACE website (http://www.dna.affrc.go.jp/PLACE/signalscan.html) was adopted to identify putative cis-elements among the promoter sequences 56 .
Plant material and growth conditions. Moso bamboo seeds were collected from Guilin, Guang Xi Province, China. The seeds from individual plants were germinated on sterile filter papers in culture dishes, while keeping the filter papers moist and in darkness at 25 °C. The seedlings were transferred to plastic pots containing vermiculite and grown in an illuminated incubator with 16/8 h of light/dark at 25/18 °C and 80% humidity, and watered with Hoagland nutrient solution every week. Plants were cultivated for three months. For the stress treatments, young leaves were sprayed with either 100 uM MeJA or 20% PEG-6000 solution and sampled at five time points (1, 3, 6, 12 and 24 h) after treatment. All samples were immediately frozen in liquid nitrogen and stored at − 80 °C for RNA extraction after collection. The untreated samples were used as the control (0 h).

RNA isolation and qRT-PCR.
The total RNA from young leaf was extracted using TRIzol reagent (Invitrogen, Ca, USA) according to the manufacturer's instructions. The total RNA was extracted from frozen samples using an RNAprep Pure Plant Kit (Tiangen) according to the manufacturer's instructions. The first-strand cDNA was then synthesized using a PrimeScriptTM RT Reagent Kit (TaKaRa). Gene-specific primers were designed using Primer Express 3.0 and Tonoplast intrinsic protein 41 gene (TIP41) was used as the reference gene 57 . Real-time PCR was performed on an ABI 7300 Real-Time system (Applied Biosystems). Each reaction was carried out in a final volume of 20 ul containing 12.5 μl of SYBR Green Master Mix reagent (Applied Biosystems), 1.5 μl of cDNA sample, and 1μl gene-specific primers. Each primer pair of PeIQD genes was listed in Table S6. The qPCR reaction conditions were as follows: 50 °C for 2 min, 95 °C for 10 min, 40 cycles of 95 °C for 15 s, and annealing at 55-60 °C for 30 s. A melting curve was generated to analyze the specificity of the reactions, and three biological replicates were made for each biological replicate. The relative expression level was calculated as 2 −ΔΔCT [ΔC T = C T , Target − C T , CYP2 . ΔΔCT = ΔC T , treatment − ΔC T , CK (0 h)].The relative expression level [2 −ΔΔCT, CK (0 h) ] in the control plants (without treatment) was normalized to 1 as described previously 58 .

Conclusions
In the current study, we systematically identified and characterized by bioinformatics the plant-specific IQD gene family of putative calmodulin target proteins in the only Bambusoideae plant, moso bamboo (Phyllostachys edulis). We explored 29 IQD genes in the moso bamboo genome and explored their expression profiles under drought (PEG) and MeJA conditions, which were investigated using quantitative real-time PCR (qRT-PCR). The qRT-PCR results elucidated the precise role of the individual PeIQD gene. Yeast two-hybrid analysis revealed that PeIQD8 can interact with PeCaM2. This study presents a thorough overview of the moso bamboo IQD gene family and provides a new perspective on the evolution of this gene family. These data will provide an insight into further understanding of the functions of IQD members and their roles bamboo in moso bamboo growth and development.