Characterisation of microRNA promoters in the murine genome

<p>This poster was originally presented at the FEBS EMBO meeting held in Paris from 30th August to 4th September 2014.</p> <p><strong>Summary</strong></p> <p>The promoter of protein coding genes is conventionally considered to be located a few kilobases around their transcriptional start site (TSS). However – given the short half-life of miRNA primary transcripts (pri-miRNAs) – the TSS of miRNA is in most cases not annotated, therefore the position of the promoter is usually not know and potentially located tens or hundreds of kilobases upstream of the annotated miRNA gene (pre-miRNA). Several attempts have been made to identify miRNA TSSs and promoters in human, mouse, worm and plants. Nevertheless, to our knowledge, most of the studies in mouse only made use of a limited number of experimental techniques and/or sources of information available at the time of publication. To overcome these limitations we implemented an integrative approach based on computational predictions, sequence annotations and high-throughput chromatin modification data generated by the ENCODE project to produce an accurate and comprehensive annotation of miRNA promoters in the mouse genome. 
For each murine miRNA annotated in miRbase version 18 we scanned a genomic region of 100 kb upstream of it looking for features indicative of promoters. We found on average 15.7 candidate promoters for each miRNA with a mean size of ~1.4 kb in accordance with the typical size of known promoters. We also found that they were typically located 5 kb upstream of the miRNA, supporting the validity of our approach and in accordance with previous reports. When we inspected what chromatin modifications supported each of the candidate promoters we found that nearly all of them have a DNase Hypersensitivity Site (DHS), more than 50% have both a DHS and H3K4 trimethylation (H3K4me3), and more than 30% have a DHS mark, an H3K4me3 mark and a CpG island. We also found that the highest scoring promoters are on average more conserved than all other promoters and more conserved than random genomic locations of the same size, strongly supporting the validity of our approach. Finally, by chromatin immunoprecipitation assay (ChIP) and methyl-CpG immunoprecipitation (mCIP) we validated H3K4me3 and CpG methylation of the predicted promoters of selected miRNAs, confirming the validity of our approach. Taken together, these results provide a valuable resource that will help to shed light on the transcriptional regulation of miRNAs in the murine genome.</p> <p> </p> <p><strong>Methods</strong></p> <p>- CpG islands were downloaded from the UCSC genome browser in BED format and intersected with the miRbase annotation of murine miR to retain only those within 100kb upstream of any miR.</p> <p>- DNAse Hypersensitivity Sites (DHS): data were downloaded from the DNaseI Hypersensitivity by Digital DNaseI from ENCODE/University of Washington (UW DNaseI HS) track of the UCSC genome Browser. We used Hotspots for replicate 1 and 2 for all available samples. Data were downloaded in BED format using the Subtrack Merge function. All datasets in the track were merged (48 cell lines) and peaks common to at least 12 samples were retained.</p> <p>- H3K4me3 data was downloaded from the UCSC genome browser; track ENCODE-PSU/LICR/Caltech (all datasets, 32 different cell lines). Additionally, we also used H3K4me3 ChIP-Seq data published in (Mikkelsen et al., 2007). All datasets were merged and peaks common to at least 10 samples were retained.</p> <p>- Eponine predictions of TSSs. For each miR annotated in miRbase 18 we obtained the DNA sequence of the region 100kb upstream of the 5’ end of each miR using bedtools v2.16.2 (Quinlan and Hall, 2010) and we scanned it with Eponine (threshold 0.990) to identify candidate TSS.</p> <p>- Annotated TSS of host gene for intronic miRs. For each miR annotated as intronic according to miRbase 18 we extracted the coordinates of the 5’ end of the host transcript from ENSEMBL v64. The 5’ coordinates were expanded 50nt each side to obtain a TSS region of 100nt.</p> <p>- Genome-wide chromatin segmentation using chromHMM. The chromHMM track was downloaded from the UCSC genome browser. From the track we extracted all regions marked as promoters, merged them with bedtools and retained only peaks within 100kb of any miRbase annotated human miR. The coordinates of each human peak were then converted to the corresponding syntenic regions in the murine genome using liftover (minMatch=0.6) and intersected with bedtools to retain only those within 100kb of any murine miR.</p> <p>Of the above feature, all those bigger than the 3rd quartile of the size distribution plus 1.5 time the inter-quantile range were fragmented into n fragments of the mean feature size. After this step, all features were clustered using bedtools to produce an annotation of candidate promoters. Each resulting cluster within 100kb upstream of any miR was then intersected with the feature annotation using bedtools annotate to calculate how many features supported it. We then assigned a score to each cluster by summing the number of features that mapped to it (all features had the same weight, with the exception of DHS that weighed 20% more than the others) and weighed the resulting score according to the distance of the cluster from the corresponding miR using the following formula: Score=nFeature * (1 - ( sqrt(distance)/ sqrt(10^6) ) )</p> <p>For each miR the cluster with the highest score was the denominated the “best” cluster for that miR.</p> <p>To measure the promoter conservation the coordinates of each cluster were intersected with the PhasCons30wayPlacental (Siepel et al.) track of the UCSC genome browser using the R package rtracklayer (Lawrence et al., 2009).</p> <p>The cluster annotation was then imported in R where we calculated a z-score for each cluster. All plots were produced in R using ggplot2.</p

Sma-miR-277 family predicted targets downregulated in developing female parasites.

<p>Fold change expression (Log₂) of high confidence targets of sma-miR-277 family during the development of male and female worms in two conditions: paired (solid line red) and unpaired (dashed green). Black lines represent the mean expression of genes in paired (solid black line) and unpaired (dashed black line) worms.</p

Sma-miR-277 and sma-miR-4989 belong to a gene cluster.

<p>(A) The genomic locus in Chromosome 4 of sma-miR-277 and sma-miR-4989 suggests they belong to a gene cluster. The average distance between genes (represented by coloured boxes) is 109 bases. Here the mature miRNAs (sma-miR2-277 and sma-miR-4989) and passenger miRNAs are represented with coverage plot and aligned reads from one of the small RNA libraries. (B) Predicted stem-loop structures for sma-miR-277 and sma-miR-4989 –individual cases. Mature miRNAs are located in the 3’-end of the hairpin. (C) Due to the cluster organisation of sma-miR-277/sma-miR-4989, it is likely that they are transcribed as one precursor RNA molecule. This figure represents the predicted stem-loop structure for sma-miR-277/sma-miR-4989 when arising from a larger transcript.</p

Sma-miR-4989 is expressed in the cells surrounding the oesophagus and cells of the tegument in adult worms.

<p>Whole mount <i>in situ</i> hybridisation for (A) <i>cathepsin B</i>, (B) sma-miR-124a-3p (<i>124a</i>), and (C) sma-miR-4989. (D) Fluorescence <i>in situ</i> hybridisation showing the colocalization of sma-miR-4989 with four co-expressed tegument-specific mRNAs (<i>calpain</i>, <i>npp-5</i>, <i>annexin</i> and <i>gtp-4</i>). Nuclei are stained with DAPI and shown in blue. Anterior of worms is to the left in A-C. Scale Bars: A-C 100 μm; D 10 μm.</p

Sma-miR-4989 is significantly up-regulated during male and female maturation.

<p>Fold change expression of sma-miR-4989 during development of juvenile to adult worms in male (blue bars) and females (red bars) as measured by RT-qPCR. Samples were collected at the time points (days post infection) indicated in the x-axis from murine hosts infected with pooled (mixed sex) cercariae. Each barplot represents the mean of three biological replicates. T-tests were performed between 49 d.p.i. and 28 d.p.i. and were both significant with p-value ≤ 0.01. Error bars show the standard error of the mean, based on three biological replicates.</p

MiRNA target prediction based on both miRNA-unaware and miRNA-guided approaches.

<p>(A) Sylamer enrichment landscape plots for 7mers in male (top) and female (bottom) expression data. The x-axis represents a list of transcripts, ranked from more expressed in juveniles to more expressed in adults. The y-axis represents the significance values acquired for each 7mer at each position in the ranked list of transcripts. Coloured boxes represent the fraction of transcripts significantly (adjusted p-value < 0.01) differentially expressed between juvenile and adult worm as found using DESeq2. These transcripts were subsequently filtered based on the presence of the 7mers TGCATTT or GCATTTA as found by Sylamer. The resulting sets are referred to as Male and Female Sylamer genes. (B) Venn Diagram showing the intersection of Male and Female Sylamer genes with schistosome-conserved miRNA targets as found using TargetScan with conservation + miRanda. The overlap represents transcripts with highly conserved sma-miR-277 target sites across the three <i>Schistosome</i> spp (S. <i>mansoni</i>, <i>S</i>. <i>haematobium</i> and S. <i>japonicum</i>) that are also significantly down regulated during worm development.</p

High confidence targets of the sma-miR-277 family identified with a combined approach (Sylamer, miRanda and TargetScan with conservation).

<p>F = female; M = male.</p

Candidate target genes of miR-25.

a<p>number of sites in the 3′UTR that are complementary to a 7 nt miR-25 seed region (GTGCAAT or TGCAATA).</p>b<p>miR-25 transfection in Dgcr8 KO cells compared to a control transfection.</p>c<p>Dgcr8 KO cells compared to Dgcr8 Het cells.</p>d<p>25-iPS cells compared to 4F-iPS cells.</p>e<p>expression fold change.</p>f<p><i>p</i> values adjusted for multiple testing (Benjamini & Hochberg).</p

Characterization of 25-iPSCs.

<p><b>a.</b> Immunofluorescence-staining of iPSCs to detect SSEA1 and Nanog expression. DNA was stained with DAPI (blue). <b>b.</b> RT-PCR analysis of pluripotency gene expression in 25-iPSCs and 4F-iPSCs. MEFs and ES cells serve as the negative and positive controls, respectively. Beta-actin was used as the PCR control. dH2O: no DNA template. <b>c.</b> Bisulfide genomic sequencing of the promoter regions of <i>Oct4</i> and <i>Nanog</i> loci to detect CpG methylation. Open and filled circles represent unmethylated or methylated CpG dinucleotides, respectively. <b>d.</b> Normal karyotypes in a male 4F-iPS cell and a male 25-iPS cell. <b>e.</b> No exogenous Yamanaka factor expression in 25-iPSCs. RT-PCR was performed using primers to amply OCKS junctions (Primer sequences are in <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0040938#pone.0040938.s002" target="_blank">Table S2</a>). Robust expression from PB-TRE-OCKS was detected in the presence of Dox. No expression was found if Dox was removed from the media. +RT: reverse transcriptase; -RT: no reverse transcriptase in reverse transcription. Beta-actin was used as the PCR control. <b>f.</b> Haematoxylin eosin stained sections of teratomas derived from 25-iPSCs. Cells representing all three germ layers were readily identifiable. <b>g.</b> Germline-transmission pups (agouti) from two chimeras derived from 25-iPSCs that were crossed to wild-type C57BL6 (albino) females.</p

Identification of miR-25 targets.

<p><b>a.</b> Timeline of miR-25 mimic transfection into the <i>Dgcr8</i>-deficient ES cells for discovering transcripts that were repressed by miR-25. The <i>C. elegans</i> miRNA, cel-miR-239b, was used as the negative control. <b>b.</b> Sylamer plots following the analysis of microarray data. The x-axis represents the 3′UTR sequences of all transcripts, sorted from the most down-regulated after miR-25 transfection compared to the cel-miR-239b transfection. Sylamer was used to analyze the enrichment or depletion of all possible heptamers matching mouse miRNA seeds, at every 500 sequences of the gene list. The horizontal dashed lines denote Bonferroni corrected p-value thresholds of 0.05. Words matching the mmu-miR-25 and cel-miR-293b seeds are highlighted. <b>c.</b> Cumulative fraction plot of the fold-changes obtained by comparing 25-iPS <i>vs</i> 4F-iPS microarrays. The fold-changes for all genes are indicated by the black curve; the green curve represents the 54 miR-25 targets from <a href="http://www.plosone.org/article/info:doi/10.1371/journal.pone.0040938#pone.0040938.s004" target="_blank">Table S4</a>.</p