Identification and comparative analysis of drought-associated microRNAs in two cowpea genotypes

Drought is one of the main abiotic factors that cause reduction or total loss of crop production. Because water is becoming limited for agriculture in many areas of the world, the investigation of natural mechanisms of drought tolerance is an important strategy for understanding the biological basis of response to drought stress and for selection of plants with improved drought tolerance [1, 2]. Cowpea [Vigna unguiculata (L.) Walp.] is an economically important crop in semi-arid and arid tropical regions in Africa, Asia, and Central and South America, where cowpea is consumed as human food and nutritious fodder to livestock [3, 4]. As a leguminous species, cowpea belongs to the same tribe (Phaseoleae) as common bean and soybean. Compared to these close relatives and most other crops, cowpea is well adapted to these regions because of its ability to fix nitrogen in poor soil and greater drought tolerance [4, 5]. Therefore, cowpea is an excellent system for investigating the genetic basis of drought tolerance.

Efforts have been made to identify genetic elements that are involved in drought stress response in cowpea. For example, over a dozen genes have been shown to be associated with drought stress response through cloning and characterization of cDNAs [6–12]. In addition, ten drought tolerance quantitative trait loci (QTL) associated with tolerance in seedlings have been mapped in cowpea [13]. However, it is largely unknown how the expression of drought-associated cowpea genes or loci is regulated and how small RNAs are involved in the regulation.

MicroRNAs (miRNAs) are 20-24 nt single-stranded RNA molecules that are processed from RNA precursors that fold into stem-loop structures [14]. MiRNAs regulate gene expression of target mRNAs at the posttranscriptional level, which are recognized by nearly perfect base complementarity. Upon miRNA-target recognition, typically the target is negatively regulated via mRNA cleavage or translational repression [15]. Functional analyses have demonstrated that miRNAs are involved in a variety of developmental processes in plants [16]. In addition, miRNAs play critical roles in plant resistance to various abiotic and biotic stresses [17–19]. In particular, several approaches have been employed to study miRNAs that are involved in drought stress tolerance in plants. In one of the pioneering studies on stress-responsive miRNAs, Sunkar and Zhu [20] used small RNA cloning techniques to identify 26 novel miRNAs, among which miR393, miR397b, and miR402 were upregulated by dehydration and miR389a downregulated. Another miRNA family, miR169, was found to be downregulated by drought stress in an ABA-dependent pathway. The repression of miR169 leads to higher expression of its target gene NFYA5, which in turn enhances the drought resistance of the plant [21]. Many more miRNAs that are up- or down-regulated in drought condition were discovered by global miRNA expression profiling experiments with either microarray hybridization or small RNA deep sequencing [22–25].

Although numerous miRNAs have been identified in many plant species, including leguminous plants Medicago truncatula [26, 27], soybean [28], and common bean [29], only two sequences have been reported for cowpea in the miRBase registry. Recently, 47 potential miRNAs belonging to 13 miRNA families were predicted in cowpea [30]. In another study, 18 conserved miRNAs belonging to 16 families were identified [31]. Both studies used a homology search approach to identify cowpea miRNAs that are conserved in other plants. In this study, we used Illumina deep sequencing technology to generate small RNA reads and used these reads to identify miRNAs in cowpea, especially cowpea-specific miRNAs and those associated with drought tolerance. To our knowledge, this is the first report of miRNAs identified through direct small RNA cloning in cowpea.

Despite inherent drought tolerance, cowpea varieties display significantly different levels of drought tolerance [32–34]. The study and comparison of plant genotypes differing in sensitivity to drought is a promising approach to discover natural tolerance mechanisms [35]. In order to gain insight into the role of miRNAs in tolerance to drought, we used two representative cowpea genotypes: California Blackeye No. 46 (CB46) and IT93K503-1. The drought-sensitive CB46 is the most widely grown blackeye-type cultivar in the United States and was developed at the University of California, Davis [36]. IT93K503-1 is a drought-tolerant breeding line developed by the International Institute of Tropical Agriculture (IITA) in Ibadan, Nigeria. We grew these two genotypes in well-watered and drought stress conditions and used leaves from the vegetative stage to construct four small RNA libraries. Using small RNA reads from these libraries, we identified 157 candidate miRNAs. Comparison of the expression pattern of miRNAs among libraries indicates that some miRNAs display different levels of expression in different genotypes, and thus may be a key factor to their different levels of drought tolerance.

Regulation of gene expression through sequence-specific interaction between miRNAs and their target mRNAs offers an accurate and inheritable mechanism for plants to respond to environment stimuli [18]. Due to water limitations, drought is a major stress that limits the geographic distribution and yield of many crops. Therefore, extensive effort has been made for discovering genetic elements and mechanisms of drought tolerance, including the discovery of drought-associated miRNAs. As an important drought-tolerant crop in semi-arid and arid areas, cowpea offers a good system for the study of drought tolerance. Here we used deep sequencing of small RNA libraries from two cowpea genotypes and identified 157 miRNAs. By comparing the expression level of miRNAs in drought-stressed sample to control sample, we also identified 30 miRNAs that were upregulated in drought condition and 14 downregulated. This list of drought-associated miRNAs includes miRNA families that were known to be associated with drought in other plant species, indicating that they are involved in conserved drought response pathways. Some miRNA families, including some cowpea-specific miRNAs, were found to be associated with drought for the first time, suggesting that they may be involved in lineage- or species-specific stress response pathways and functions.

We predicted target genes for 32 out of 44 drought-associated miRNAs. The predicted target mRNAs encode proteins of diverse function, most of them being transcription factors (Additional file 4). For most of the conserved miRNAs, it is expected that their targets are also conserved. For example, our results showed that miR156 was upregulated in response to drought in cowpea. MiR156 has been known to be responsive to abiotic stresses and targets SPB transcription factors in Arabidopsis, maize, rice and wheat [24, 39, 41–43]. This miRNA is also involved in the regulation of development during vegetative phase change [42], indicating that reprogramming of development is a crucial step in plants to cope with drought stress. Another miRNA, miR169, was downregulated in both cowpea genotypes. In Arabidopsis, miR169 was downregulated and its target, a Nuclear Factor Y transcription factor NFYA5, was induced by drought stress [21]. MiR169 most likely functions in a similar way in cowpea to enhance drought tolerance by inducing the expression of NFYA5 orthologs.

The cowpea genotypes studied in this work have different abilities of drought tolerance. Because the two genotypes are highly similar to each other in their genetic composition, their phenotypic variations such as drought tolerance are most likely caused by changes in regulatory processes, rather than changes in proteins [44]. Due to their different geographical origins, the two genotypes are adapted to the particular environmental conditions in their natural habitats. It is thus expected to find constitutive differences, which could be related to metabolism, use of energetic resources, mobilization of biomass, structure of radical system, wax deposition in leaves, membrane stability or density of stomata, among other characteristics. We found that nine miRNAs were predominantly or exclusively expressed in only one genotype, regardless of the treatments. On the other hand, 11 miRNAs were found to be differentially expressed under drought stress in one genotype, but not the other. Changes in miRNA expression are expected to cause changes in the expression of target genes between the two genotypes.

Among miRNAs that had genotype-specific regulation, miR160a and miR160b were upregulated in response to drought in the tolerant, but not in the sensitive cultivar (Table 2). Their putative targets are members of the family of Auxin Response Factors (ARFs). ARFs are key elements in regulation of physiological and morphological mechanisms mediated by auxins that may contribute to stress adaptation [45]. Moreover, negative regulation of ARF10 by miR160 was demonstrated to be critical during seed germination in Arabidopsis thaliana through the crosstalk between auxin and ABA-dependent pathways [46]. On the other hand, two members of the miR2111 family were upregulated by drought in the sensitive, but not in the tolerant cultivar (Table 3). Their putative targets are Kelch repeat-containing F-box proteins that belong to a large family with members known to be involved in response to biotic and abiotic stresses [47]. Furthermore, F-box proteins containing Kelch repeats have been found to be responsive to drought in chickpea, a close relative of cowpea [48]. This suggests that genotype-specific regulation of miRNAs might be part of the reason why some cowpea genotypes have stronger drought tolerance than others.

Among the new miRNA candidates that were identified in this study, ten were regulated by drought stress and target genes were predicted for five of them. For instance, vun_cand030 was downregulated by drought and putatively targets a zinc finger protein. Zinc finger proteins are known to be involved in a variety of functions in development and stress response [49]. Moreover, vun_cand015 was upregulated by drought in the tolerant cultivar and putatively targets a basic-helix-loop-helix (bHLH) transcription factor. These proteins have roles in response to abiotic stresses, such as iron deficiency [50], freezing, and salt stress [51]. This suggests these new miRNAs may be indeed an integral component of drought response in cowpea.

For many miRNAs there were more than one target predicted. The possibility of a miRNA to have multiple targets is commonly observed. To confirm these predicted targets, we need to perform detailed analysis of cleavage of mRNA targets at the miRNA recognition site by experimental approaches, such as RACE and degradome analysis [52–54]. Once we validate the targets of drought-associated miRNAs, we will be in a better position to link the expression changes of miRNAs and their targets to differences of drought tolerance in cowpea.

Because we do not have the complete cowpea genome sequence, some miRNA genes were not identified, even though they had significant expression in our small RNA libraries. To find out how many miRNA families have been missed, we mapped unique small RNA reads to plant miRNA precursors in the miRBase, allowing up to 2 mismatches. Although we did not miss a large number of miRNAs, we did find that miR2118, miR2911, and miR529 had significant expression in our libraries (Additional file 5). The latter two were also induced by drought stress. MiR529 was identified as drought-associated miRNA in rice [25]. However, contrary to the pattern that we found in cowpea, it was downregulated under drought stress in rice. It is not clear whether it was caused by different sampling time or tissue, or species-specific stress response mechanisms.

Like protein coding genes, many miRNA families possess more than one miRNA gene and miRNA genes from the same family may have either identical or similar but different mature miRNA sequences. During evolutionary process, homologous miRNA genes may functionally diverge from each other. In the set of miRNAs that we identified in cowpea, members from miR166 and miR167 families showed clear evidence for functional diversification. While one member miRNA gene (miR166a, miR167b) was induced by drought stress, another miRNA from the same family (miR166b, miR167a) was significantly downregulated (Additional file 4).

Using deep sequencing technology, we identified 157 miRNAs in cowpea, including 44 miRNAs that are drought-associated. By comparing mature miRNA counts in different genotypes and growth conditions, we found 9 miRNAs that were almost exclusively expressed in only one genotype and 11 miRNAs that were regulated by drought stress in one genotype, but not the other. Our study demonstrated that deep sequencing of small RNAs is a cost-effective way for miRNA discovery and expression analysis. Compared to the homology search method, deep sequencing allowed the detection of species-specific miRNAs and digital expression analysis. Our findings demonstrate that expression patterns of some miRNAs may be very different even between two genotypes of the same species. Further characterization of the targets of drought-associated miRNAs will help understand the details of response and tolerance to drought in cowpea.