The integrated analysis of RNA-seq and microRNA-seq depicts miRNA-mRNA networks involved in Japanese flounder (Paralichthys olivaceus) albinism

Abstract

Albinism, a phenomenon characterized by pigmentation deficiency on the ocular side of Japanese flounder (Paralichthys olivaceus), has caused significant damage. Limited mRNA and microRNA (miRNA) information is available on fish pigmentation deficiency. In this study, a high-throughput sequencing strategy was employed to identify the mRNA and miRNAs involved in P. olivaceus albinism. Based on P. olivaceus genome, RNA-seq identified 21,787 know genes and 711 new genes by transcripts assembly. Of those, 235 genes exhibited significantly different expression pattern (fold change ≥2 or ≤0.5 and q-value≤0.05), including 194 down-regulated genes and 41 up-regulated genes in albino versus normally pigmented individuals. These genes were enriched to 81 GO terms and 9 KEGG pathways (p≤0.05). Among those, the pigmentation related pathways-Melanogenesis and tyrosine metabolism were contained. High-throughput miRNA sequencing identified a total of 475 miRNAs, including 64 novel miRNAs. Furthermore, 33 differentially expressed miRNAs containing 13 up-regulated and 20 down-regulated miRNAs were identified in albino versus normally pigmented individuals (fold change ≥1.5 or ≤0.67 and p≤0.05). The next target prediction discovered a variety of putative target genes, of which, 134 genes including Tyrosinase (TYR), Tyrosinase-related protein 1 (TYRP1), Microphthalmia-associated transcription factor (MITF) were overlapped with differentially expressed genes derived from RNA-seq. These target genes were significantly enriched to 254 GO terms and 103 KEGG pathways (p<0.001). Of those, tyrosine metabolism, lysosomes, phototransduction pathways, etc., attracted considerable attention due to their involvement in regulating skin pigmentation. Expression patterns of differentially expressed mRNA and miRNAs were validated in 10 mRNA and 10 miRNAs by qRT-PCR. With high-throughput mRNA and miRNA sequencing and analysis, a series of interested mRNA and miRNAs involved in fish pigmentation are identified. And the miRNA-mRNA regulatory network also provides a solid starting point for further elucidation of fish pigmentation deficiency.

Fig 1. The volcanoplot of differentially expressed genes in RNA-seq.

The differentially expressed genes including 41 up-regulated genes and 194 down-regulated genes in PO_alb vs PO_con were illustrated in a volcanoplot (fold change≥2 and q-value≤0.05).

Fig 2. The expression levels of ten genes in qRT-PCR and RNA-seq.

The expression fold change (FC) 10 P. olivaceus genes in PO_alb versus PO_con detected by qRT-PCR and RNA-seq were calculated by 2-ddCt and readcount, respectively. And these genes’ log2FC values of qRT-PCR and RNA-seq are shown in A and B, respectively. In qRT-PCR analysis, β-actin expression levels were used for the internal control, and values are indicated as means ± standard error (SE) derived from triplicate experiments.

Fig 3. KEGG pathways enrichment in RNA-seq (p<0.05).

Gene number: number of target genes in each pathway. Rich factor: the ratio of the number of target genes divided by the number of all the genes in each pathway.

Fig 4. Length distribution of small RNAs.

The clean reads of small RNAs were distributed from 18 nt to 26 nt in six libraries including PO_con1, PO_con2, PO_con3, PO_alb1, PO_alb2, and PO_alb3. The reads of 20–24 nt accounted for 84–93% of small RNAs.

Fig 5. Conservation profile of the identified miRNA.

The frequency of identified miRNAs is calculated by mapping with 37 species (horizontal axis) in miRbase, and the miRNA counts are shown on the vertical axis.

Fig 6. The illustration of four novel miRNAs in P. olivaceus.

The secondary structures formed by four novel miRNA precursors are illustrated by RNAfold. The colour bar represents base-pair probabilities from 0 (blue) to 1 (red).

Fig 7. Clustering of expression patterns of 33 differentially expressed miRNAs.

The expression patterns of 33 differentially expressed miRNAs (p<0.05) in the PO_con2, PO_con3, PO_alb2, and PO_alb3 libraries are displayed in a pheatmap. Each column represents one library, and the colour bar indicates relative expression level from high (red) to low (green).

Fig 8. The expression levels of 10 randomly selected miRNAs in qRT-PCR and miRNA-seq.

The expression fold changes of 10 P. olivaceus miRNAs in PO_alb versus PO_con detected by poly (A)-tailed qRT-PCR and high-throughput sequencing are shown in A and B, respectively. In qRT-PCR analysis, 5s rRNA expression levels were used for the internal control, and values are indicated as means ± standard error (SE) derived from triplicate experiments.

Fig 9. The top 20 enriched GO terms (A) and KEGG pathways (B).

Gene number: number of target genes in each term or pathway. Rich factor: the ratio of the number of target genes divided by the number of all the genes in each term or pathway.

Fig 10. A proposed network of putative interactions between miRNAs and mRNAs in tyrosine metabolism.

The regulation network of miRNAs and mRNAs involved in tyrosine metabolism is illustrated by Cytoscape. Yellow ellipses represent miRNAs, and red rectangles indicate their target genes associated with tyrosine metabolism. Abbreviations of the genes are detailed in the text.