Evolution and function of the extended miR-2 microRNA family.

MicroRNAs are essential post-transcriptional regulators. Many animal microRNAs are clustered in the genome, and it has been shown that clustered microRNAs may be transcribed as a single transcript. Polycistronic microRNAs are often members of the same family, suggesting a role of tandem duplication in the emergence of clusters. The mir-2 microRNA family is the largest in Drosophila melanogaster, with 8 members that are mostly clustered in the genome. Previous studies suggest that the copy number and genomic distribution of mir-2 family members has been subject to significant change during evolution. The effects of such changes on their function are still unknown. Here we study the evolution of function in the mir-2 family. Our analyses show that, in spite of the change in number and organization among invertebrates, most mir-2 loci produce very similar mature microRNA products. Multiple mature miR-2 sequences are predicted to target genes involved in neural development in Drosophila. These targeting properties are conserved in the distant species Caenorhabditis elegans. Duplication followed by functional diversification is frequent during protein-coding gene evolution. However, our results suggest that the production of microRNA clusters by gene duplication rarely involves functional changes. This pattern of functional redundancy among clustered paralogous microRNAs reflects birth-and-death evolutionary dynamics. However, we identified a small number of mir-2 sequences in Drosophila that may have undergone functional shifts associated with genomic rearrangements. Therefore, redundancy in microRNA families may facilitate the acquisition of novel functional features.

Introduction

MicroRNAs, crucial regulators of gene expression at the post-transcriptional level, are often clustered in the genome. According to miRBase, more than a quarter of both Drosophila and human microRNAs are less than 10 kb away from other microRNAs. These clustered microRNAs are often co-expressed, suggesting that they are produced from a single transcript.- The majority of microRNA clusters contain members of the same family, indicating a major role of tandem duplication in cluster formation.- In the case of protein-coding genes, duplication is acknowledged as the main source of functional innovation, since duplicates are free to diversify in their functions. Similarly, duplicated microRNAs may acquire new targets and therefore novel functions. However, microRNAs processed from the same transcript are linked by their expression pattern, imposing a functional constraint on their evolutionary diversification. Whether microRNA tandem duplications facilitate the emergence of new functions or generate redundant products remains to be explored.

Mir-2 is the largest microRNA family in Drosophila melanogaster and one of the first to be discovered.- The mir-2 family has 8 members in the D. melanogaster genome (mir-2a-1, mir-2a-2, mir-2b-1, mir-2b-2, mir-2c, mir-13a, mir-13b-1 and mir-13b-2), six of which are organized in two clusters. In most other studied insects, there are five mir-2 sequences encoded by a single transcript (see ref. 15 and references therein). Caenorhabditis elegans has only one mir-2 sequence.,

Here we study the mir-2 family to investigate the impact of microRNA family expansions on functional diversification. We combine comparative genomics with expression data analyses and functional annotation of predicted targets to compare the functional features of mir-2 sequences. Our results will help us to understand the role of tandem microRNA duplications in the evolution of gene regulation.

Results

Mir-2 is a conserved microRNA family in invertebrates

In order to characterize mir-2 family members, we performed comprehensive sequence similarity searches against multiple sequenced organisms (see Materials and Methods). We detected mir-2 hairpin precursor sequences in many invertebrates (Fig. 1A; ) but none in vertebrate species. The 3′ arm of the hairpin is highly conserved, although the many changes in the 5′ arm are fully consistent with the precursor hairpin structure (Fig. 1A). The gene copy number is highly variable among species (from one in C. elegans to eight in D. melanogaster) suggesting that the mir-2 content of each lineage is the product of multiple birth-and-death events.

Figure 1. Sequence conservation in the mir-2 family. (A) The alignment of mir-2 precursor sequences in representative genomes, shadowed by sequence conservation (visualized using Ralee), where darker tones reflect higher conservation. Structure of the consensus sequence is shown below the alignment in dot-bracket annotation. The open white box over the alignment indicates the canonical mature product, with the seed sequence highlighted (black). (B) Consensus structure of the mir-2 precursor in invertebrates, colored with VARNA according to sequence conservation. The canonical and non-canonical mature products produced by some mir-2 precursors are also indicated.

Since mir-2 sequences are short and very similar, their genomic contexts can improve our ability to annotate and explore their evolutionary origins. The genomic organization of mir-2 family members across phyla (Fig. 2) suggests that the ancestral mir-2 microRNA was clustered with mir-71, an evolutionarily unrelated microRNA. Mir-71 itself is found in protostomes, but also in cephalochordates, hemichordates and echinoderms., The origin of mir-71 therefore pre-dates the split of protostomes and deuterostomes, although it has been lost in chordates. Mir-2 arose later, most likely before the last common ancestor of protostomes. Although mir-71 and mir-2 are still linked in most species, mir-71 has been lost independently in two dipteran lineages. The expansion of the mir-2 family by tandem duplication and deletion has generated mir-2 clusters of different lengths in different species. The mir-13 subfamily has a conserved characteristic one-nucleotide deletion in its 3′ arm (Fig. 1A), indicating that these sequences originated from duplicated mir-2 locus in the common ancestor of insects. Combined analysis of sequence conservation and cluster structure (Figs. 1A and 2) suggests that the ancestral insect cluster split in two in the Drosophila lineage, with subsequent additional duplications. As a consequence, different mir-2 copies in Drosophila are under the transcriptional control of different regulatory sequences.

Figure 2. Copy distribution of mir-2 sequences. Phylogenetic tree of invertebrate species and genomic organization of mir-2 sequences. Divergence times were extracted from ref. 47. Black arrows depict mir-2 family members, and white arrow mir-71 sequences. Arrows linked by the same straight line indicate microRNAs linked in the genome by less than 10 kb.

Functional conservation and redundancy of mir-2 products

The pattern of sequence conservation in the mir-2 family sequences shown in Figure 1 suggests that the dominant mature microRNA is produced from the 3′ arm of mir-2 precursors. Our re-analysis of deep-sequencing data from D. melanogaster, Tribolium castaneum and C. elegans confirms that the 3′ arm is highly expressed compared with the 5′ arm in most mir-2 family members (). Deep-sequencing analyses from honeybee and silkworm also reveal the same expression pattern., We observe that most mir-2 sequences conserve the location of the Drosha and Dicer cleavage sites. This position determines the first nucleotides of the microRNA, and hence the seed sequence. The seed is defined as nucleotides 2 to 7 of a mature microRNA, and it is crucial for transcript targeting. Since sequence conservation is very high in the 3′ arms, seed sequences are the same for all mir-2 family products in which the Dicer cleavage site is conserved (Fig. 1).

Functional shifts in mir-2 products

Deep sequencing data from Drosophila suggest that the 3′ arm of mir-2a produces two alternative mature products, in contrast to the majority of mir-2 family members. Each accounts for a significant proportion of the reads produced by mir-2a loci (47% and 28%), and they are offset from one another by 2 nucleotides. The first of these products (the 5′-most) is processed identically to the conserved mature sequence produced from the majority of mir-2 family members, termed the ‘canonical' product here (Fig. 1B). The second is offset by 2 nucleotides in the 3′ direction, and is termed the 'non-canonical' product (Fig. 1B). Both of these products map exactly to two alternative hairpin precursors called mir-2a-1 and mir-2a-2, suggesting that both products could potentially be made from either locus. However, the 5′ arms of these two hairpins are not identical in sequence, and therefore reads mapping to the 5′ arms can be assigned to one or other hairpin. It has been previously reported that the characteristic pattern of two nucleotide overhang at the 3-prime ends of mature microRNA duplexes allows the assignment of reads from the 3′ arm to one or other hairpin, even though the 3′ arm sequences are identical.,, This approach predicts that the non-canonical mature sequence, offset by 2 nucleotides, is produced overwhelmingly from the mir-2a-2 locus, whereas mir-2a-1 is processed identically to the other mir-2 family members. Analysis of deep sequencing data from an RNA immunoprecipitation (RIP-seq) of Argonaute proteins shows that both canonical and non-canonical mature products are loaded into the RNA induced silencing complex (RISC), and are therefore likely to be functional. The seed sequences of canonical and non-canonical mature microRNAs are offset, and hence differ in sequence, suggesting that they regulate different targets. Drosophila mir-2c also produces an offset, non-canonical, mature product. However this microRNA is expressed at a very low level and is not found in the AGO RIP-seq data set., Our data show that a significant fraction of non-canonical mir-2 products are also expressed from mir-2 loci in T. castaneum () and in the honeybee Apis mellifera (data not shown). However, the strategy described above cannot be applied to assign reads to a single locus.

Unlike other mir-2 members, the mir-2a-2 precursor produces approximately equal amounts of mature sequences from each arm of the hairpin., Nevertheless, mature sequences derived from the 5′ arm are not observed in AGO RIP-seq experiments, and are not, therefore, predicted to be loaded into the RISC complex. This further supports a dominant role of the mature sequence from the 3′ arm across the mir-2 family.

Mir-2 products are likely to target neural genes

We have shown that mature products from mir-2 loci are highly conserved and are likely to have the same targeting properties. Do mir-2 sequences therefore conserve their targets throughout evolution? We address this question by comparing the targets of D. melanogaster and C. elegans miR-2 mature sequences. We used the canonical seed method to predict transcripts whose 3′UTR are targeted by all miR-2 family members in D. melanogaster and the only miR-2 sequence in C. elegans (see Methods). All but two miR-2 sequences in Drosophila have identical seeds and therefore identical predicted target sets (Fig. 1A).The two microRNAs with different targets were miR-2a-2 and miR-2c, which are offset with respect to the canonical mir-2 products (Fig. 1B).

We mapped Gene Ontology terms to the predicted targets of miR-2 family members, and analyzed the set of terms that were statistically enriched in the targeted gene set (see Materials and Methods). We focused on terms within the ‘Developmental process’ category, which is particularly informative for development and tissue specificity. We detected 675 genes targeted by Drosophila the miR-2 canonical sequence, and 979 for the functional Caenorhabditis miR-2 product. For both Drosophila and Caenorhabditis, we observed an enrichment in genes involved in neural development (Table 1). We therefore predict a role for mir-2 in neural function. Indeed, expression data from deep-sequencing analyses in Drosophila indicate that mir-2 products are highly expressed in adult heads. We also studied the targets of the non-canonical products from miR-2a-2 and miR-2c in Drosophila. Both miR-2a-2 and miR-2c are predicted to target 286 genes. In these cases, we did not find any significantly enriched functional classes (not shown).

Table 1. Top 20 enriched GO terms in the developmental process category

Species	Enriched GO term	# genes1	q-value2
Drosophila	multicellular organismal development	121	0.0000
	nervous system development	67	0.0000
	central nervous system development	25	0.0000
	sensory organ development	37	0.0000
	anatomical structure morphogenesis	91	0.0000
	organ morphogenesis	48	0.0000
	neurogenesis	54	0.0000
	cell differentiation	84	0.0000
	neuron differentiation	46	0.0000
	developmental process	129	0.0000
	cell fate commitment	31	0.0000
	organ development	79	0.0000
	generation of neurons	53	0.0000
	system development	106	0.0000
	anatomical structure development	123	0.0000
	cellular developmental process	85	0.0000
	brain development	16	0.0011
	eye development	30	0.0022
	neuron development	38	0.0023
	regionalization	37	0.0023
Caenorhabditis	cellular component morphogenesis	38	0.0047
	anatomical structure morphogenesis	112	0.0062
	neurogenesis	22	0.0164
	generation of neurons	22	0.0164
	neuron development	20	0.0165
	cell morphogenesis	22	0.0167
	neuron differentiation	21	0.0167
	muscle structure development	23	0.0187
	muscle organ development	6	0.0191
	nervous system development	22	0.0204
	neuron projection morphogenesis	18	0.0209
	organ morphogenesis	10	0.0226
	axonal fasciculation	11	0.0232
	neuron projection development	18	0.0233
	anatomical structure formation involved in morphogenesis	25	0.0238
	cell projection morphogenesis	19	0.0245
	syncytium formation by plasma membrane fusion	3	0.0334
	syncytium formation	3	0.0334
	cell part morphogenesis	19	0.0347
	neuron recognition	11	0.0356

1Number of genes with predicted canonical seed targets (see Methods) annotated to a GO term; 2q-value is the p-value corrected for a false discovery rate of 0.05 (ref. 46)

The seed model for microRNA targets predicts that offset mature products from the mir-2a-1 and mir-2a-2 loci will target different sites. However, it is well established that sequence complementarity outside the seed motif is important (and perhaps even sometimes sufficient) for target recognition (reviewed in ref. 25). To explore whether offset microRNAs with the same nucleotide sequence may have different targeting properties, we predicted targets with a different tool, miRanda, which places less weight on the microRNA seed and accounts more fully for the hybridization energy between the microRNA and the target. In this particular case, Drosophila miR-2a-1 is predicted to target 553 transcripts, and miR-2a-2 putatively binds to 788, with 368 targeted genes common to both. The overlap of target genes is greater than expected by chance (p = 0.008, see Materials and Methods). This suggests that, although the seed shifting between miR-2a-1 and miR-2a-2 may induce functional changes, the two microRNAs likely conserve partially redundant targeting properties.

Discussion

The evolutionary history of microRNA families is characterized by frequent duplications, losses and rearrangements.,,, Here we describe the evolution of the largest conserved insect microRNA family: mir-2. We showed that this family is widely represented in invertebrates, and the copy number and genomic distribution varies greatly between species. Deep-sequencing data reveal that all mir-2 family members produce their dominant mature microRNAs from the 3′ arm, whose sequence is highly conserved (Fig. 1). Moreover, most mir-2 precursors have the same Dicer cleavage site, thus producing functional mature miR-2 sequences with the same seed region and predicted targets. According to the available deep-sequencing data, most mir-2 loci within the same species produce redundant products. In Drosophila, antisense-mediated inactivation of mir-2 sequences shows that multiple mir-2 loci have similar (if not identical) functions.,

It is well-established that pairs of protein-coding loci resulting from gene duplication rapidly diverge in their sequence and/or expression pattern, since functional redundancy is generally a transient situation. Duplication has been also proposed as a mechanism of microRNA functional diversification, although there is no direct evidence of this pattern so far. The mir-2 family suggests that microRNA families may tolerate a situation of functional redundancy in the longer term, as multiple almost identical copies are present in each invertebrate genome. One possible explanation is that mir-2 products are required at high levels and local tandem duplications produce a net increase in the expression level. This is supported by a previous observation that increased expression levels are associated with an increase in microRNA copy number. On the other hand, the presence of redundant mir-2 paralogs might reflect essentiality (see discussion in ref. 29). Functional redundancy in clustered paralogous microRNAs has been previously reported, and may simply reflect high turnover and birth-and-death evolutionary dynamics.,, These processes will generate clusters of very similar sequences, and account for the copy number differences between different species. The data strongly suggest that mir-2 family evolution is dominated by high turnover and birth-and-death dynamics mostly driven by random drift.

Clustered paralogous microRNAs are evolutionarily constrained since their expression pattern is linked. However, mir-2 family members in the Drosophila genus are located in two clusters and two single loci. This decoupling of their regulatory sequences may have facilitated functional changes. Indeed, we observe that the identical 3′ arms of the mir-2a-1 and mir-2a-2 hairpin precursors produce different offset mature sequences, which we call here canonical and non-canonical miR-2 products (Fig. 1B). This phenomenon is called “seed shifting,” and has been described to induce functional changes between orthologous microRNAs., Experiments in Drosophila suggest that mir-2 products are expressed in brain and have (at least partially) redundant functions., However, in situ hybridizations show that the three clusters mir-2b-2~mir-2a-2~mir-2a-1, mir-2c~mir-13a~mir-13b-1 and mir-13b-2 have different spatial expression patterns during early development. We suggest that genomic reorganizations breaking the linkage between mir-2 loci in Drosophila triggered a subfunctionalization event. Interestingly, in the flatworm Schistosoma mansoni we observe a duplication of the entire ancestral mir-2 cluster (Fig. 2 and ref. 35). The functional analysis of the mir-2 family in this parasitic species might shed light on the evolutionary dynamics of clustered microRNAs.

Mir-2 loci are highly expressed in adult heads in Drosophila and in neurons in Caenorhabditis. We show that the predicted targets of mir-2 microRNAs in both Drosophila and Caenorhabditis are significantly enriched for transcripts with neural development functions (Table 1). Mir-2 has also been found to be highly expressed in heads of Bombyx mori. Antisense-mediated inactivation of mir-2 in Drosophila produces embryos with defects in head and posterior abdominal segments. Mir-2 has been shown to specifically target the pro-apoptotic genes rpr, grim and skl. Strikingly, these three genes are involved in the selective death by apoptosis of neuroblasts during the normal development of the nervous system. By targeting these pro-apoptotic genes, mir-2 can act as an anti-apoptotic factor in neurons. Indeed, the repression of rpr and grim by ABD-B prevent apoptosis in neural cells. In the light of these data, we speculate that mir-2 microRNAs have a fundamental role in neuron survival during development and adulthood.

Finally, we note that early works associate mir-6 and mir-11 sequences with the mir-2 family because they have identical (or very similar) seed sequences (e.g., ref. 29). However, there is no evidence of an evolutionary relationship between these three families. Moreover, mir-6 and mir-11 have a distinct expression pattern from mir-2, so functional overlap among these families is unlikely.,, We strongly encourage the use of the family name mir-2 to represent only mir-2/mir-13 sequences.

In summary, the mir-2 family is an invertebrate-specific family of microRNAs probably involved in neural development and maintenance. The number and genomic organization of mir-2 loci varies greatly between species, although the function of paralogous microRNAs is most often redundant. The retention of redundant sequences may be facilitated by the co-transcription of clustered microRNAs. In Drosophila, the ancestral mir-2 cluster has split into multiple independent transcripts, decoupling the transcriptional regulation among mir-2 loci. In this species we find evidence of potential functional shifts of some mir-2 family members.

Materials and Methods

We retrieved all mir-2 precursor sequences from miRBase (version 17) and used BLAST (w = 4, r = 2, q = -3) to search for homologous sequences in multiple genomes from NCBI (www.ncbi.nlm.nih.gov/genome): Drosophila melanogaster, D. virilis, D. willistoni. D. pseudoobscura, Aedes aegypti, Anopheles gambiae, Acyrthosiphon pisum, Bombyx mori, Apis mellifera, Tribolium castaneum, Capitella teleta, Daphnia pulex, Caenorhabditis elegans, Gallus gallus, Mus musculus and Homo sapiens. We aligned the putative microRNA hairpin sequences with CMfinder (n = 5, m = 30, M = 100), chose the output alignment that best reflects the microRNA hairpin pairing, and manually refined the alignment using RALEE. The consensus sequence of the alignment was built by taking the most abundant base for each column. All columns with more than 60% gaps were excluded. VARNA 3.7 was used to visualize the consensus microRNA structure.

Small RNA libraries, with accession numbers GSE7448 (D. melanogaster), GSE15169 (C. elegans) and GSE26036 (T. castaneum), were retrieved from the GEO database (www.ncbi.nlm.nih.gov/geo/). Reads were mapped to the reference genomes using a sequential trimming approach with the SeqTrimMap tool using default parameters, and microRNAs were detected as described previously. Briefly, predicted hairpin structures within the genome with reads mapped to both arms were first extracted. We only further considered high-quality predictions after careful visual inspection.

To identify potential targets of mir-2 sequences, we first extracted 3′UTR sequences from ENSEMBL (v. 62) via Biomart for D. melanogaster and C. elegans. Only the longest isoforms were considered. We then detected canonical seed targets in these 3′UTRs. Additional target predictions were performed with miRanda, using default parameters. The significance of the observed overlap of targets between dme-miR-2a-1 and dme-miR-2a-2 was assessed as follows: first, we calculate the target list overlap between both products (number of genes with common targets divided by the total number of targeted genes); then, we calculated the target list overlap between 1000 pairs of randomly selected microRNA products in D. melanogaster; finally, the associated p-value is estimated as the proportion of random overlap measures equal or greater than the actual overlap value. In order to detect functional categories enriched in genes targeted by mir-2 products, we analyzed the annotation of transcripts with putative mir-2 target sites in Gene Ontology with HT-GOMiner. We used a false discovery rate of 0.05, and focused specifically on enriched terms inside the “developmental process” GO class.