De novo assembly and sex-specific transcriptome profiling in the sand fly Phlebotomus perniciosus (Diptera, Phlebotominae), a major Old World vector of Leishmania infantum

Abstract

Background: The phlebotomine sand fly Phlebotomus perniciosus (Diptera: Psychodidae, Phlebotominae) is a major Old World vector of the protozoan Leishmania infantum, the etiological agent of visceral and cutaneous leishmaniases in humans and dogs, a worldwide re-emerging diseases of great public health concern, affecting 101 countries. Despite the growing interest in the study of this sand fly species in the last years, the development of genomic resources has been limited so far. To increase the available sequence data for P. perniciosus and to start studying the molecular basis of the sexual differentiation in sand flies, we performed whole transcriptome Illumina RNA sequencing (RNA-seq) of adult males and females and de novo transcriptome assembly.

Results: We assembled 55,393 high quality transcripts, of which 29,292 were unique, starting from adult whole body male and female pools. 11,736 transcripts had at least one functional annotation, including full-length low abundance salivary transcripts, 981 transcripts were classified as putative long non-coding RNAs and 244 transcripts encoded for putative novel proteins specific of the Phlebotominae sub-family. Differential expression analysis identified 8590 transcripts significantly biased between sexes. Among them, some show relaxation of selective constraints when compared to their orthologs of the New World sand fly species Lutzomyia longipalpis.

Conclusions: In this paper, we present a comprehensive transcriptome resource for the sand fly species P. perniciosus built from short-read RNA-seq and we provide insights into sex-specific gene expression at adult stage. Our analysis represents a first step towards the identification of sex-specific genes and pathways and a foundation for forthcoming investigations into this important vector species, including the study of the evolution of sex-biased genes and of the sexual differentiation in phlebotomine sand flies.

Fig. 1

Analysis workflow. Flowchart illustrating the main steps of the bioinformatics pipeline: a assembly, (b) annotation and (c) differential expression analysis. FDR = false discovery rate, FC = fold change

Fig. 2

Summary of the functional annotation results. The bar plot represents the number of transcripts of the PERNI data set for each annotation category

Fig. 3

GO term distribution of the annotated transcripts of the PERNI data set compared with those of P. papatasi and L. longipalpis. GO categories are shown within the divisions of Biological Process (BP), Molecular Function (MF), and Cellular Component (CC). Column heights reflect the percentage of annotated transcripts in each assembly that mapped to a given GO term. Asterisks indicate the statistical significance of the differences observed in the pair-wise comparison between P. perniciosus and P. papatasi and between P. perniciosus and L. longipalpis (* = P < 0.05. ** = P < 0.01)

Fig. 4

Differential expression analysis. Differential Gene Expression Analysis of adult male vs. female whole bodies of P. perniciosus. a Volcano plot showing the relative expression levels of the transcripts in adult male versus adult female whole body. The x-axis represents the log2 of the expression fold change (males vs. females) for each transcript of the PERNI data set. The y-axis represents the negative log10 of the adjusted p-value (FDR) of the Fisher’s Exact test. Blue data points (n = 8590) represent the transcripts that were statistically significant sex-biased (FDR < 0.01; dotted horizontal line). Dotted vertical lines represent the biological cut-off of 2-fold change in expression between males and females. b MA-plot showing the relative expression levels of the transcripts in adult male versus adult female whole body. The x-axis represents the log2 of the counts per million of mapped reads (CPM) for each transcript of the PERNI data set. The y-axis represents the log2 of the expression fold change (FC) for each transcript of the PERNI data set. Blue data points (n = 8590) represent the transcripts that were statistically significant sex-biased (FDR < 0.01; dotted horizontal line). c Histogram of the log2 of the expression fold change (FC) values between the 8590 sex-biased PERNI transcripts. Dotted orange vertical lines represent the biological cut-off of 2-fold change in expression between males and females. d Stack plot of the distribution of the annotation types of the sex-biased PERNI transcripts

Fig. 5

Semi-quantitative RT-PCR analysis of the top sex-biased PERNI transcripts. Semi-quantitative RT-PCR analysis of selected transcripts of male and female adults of P. perniciosus. a Normalization control using primer pairs for the sod gene of P. perniciosus at 20, 25 and 30 PCR amplification cycles. b Amplification of female-biased PERNI transcripts. For the transcript F4 we observed an extra male-specific amplification signal of larger molecular size. c Amplification of male-biased PERNI transcripts. For the transcript M1 we observed an extra male-specific amplification signal of larger molecular size and for the transcripts M5a, M6, M8, M15, M17 and M19 we observed extra female-specific amplification signals of larger molecular size. For the transcript M9 we observed an extra non-sex-specific amplification signals of larger molecular size. λ indicates the molecular weight marker (λ genomic DNA EcoRI-HindIII digested)

Fig. 6

Ka/Ks plot. Plot representing the distribution of the Ka/Ks values between P. perniciosus/L. longipalpis and P. perniciosus/P. papatasi. Blue dots represent the ortholog pairs showing relaxation of selective constrains in the comparison P. perniciosus/P. papatasi