BioSankey: Visualization of Microbial Communities Over Time

2.1 Import of Data into BioSankey and Generation of the Project-Specific Website

The absolute or relative abundance of an element type (e.g. normalized microbial species abundance or gene expression) is provided by the user as a comma separated file, where a row depicts the microbial species or gene and the second column contains the description of the gene or operational taxonomic unit (OTUs) followed by the normalized read counts per condition or optionally, the expression when genes should be analyzed. These identifiers can be used to query the dataset in the website. Optionally, time series data of the fluctuation of up- and down-regulated DEGs can be provided by the user (e.g. first versus second time point), which we will describe in use case 4. The whole functionality of BioSankey is summarized in Figure 1, starting from the specification of the possible input files, the integration of the information into a web page and finally the generation of the web page. We have created a graphical user interface (GUI) to guide the user throughout the configuration process and to automatically execute the python script for generating the project-specific web page. The python scripts are implemented under version 3.6 and require no additional packages and should work both on Linux and Windows.

Figure 1:

The workflow of the tool BioSankey: the green items on top, “microbial abundances” and “gene expression”, are alternative inputs, the yellow item on the left is optional information and the blue sheets at the bottom are the output, where the dashed arrow to “PDF” is an optional step.

2.2 Visualization of Elements in the HTML Site

The user can select different numbers of genes to be shown by the filter criteria. All genes are highlighted as Sankey diagrams, whereas if more than the defined amount of genes are searched, genes must be selected from a selection box. The search panel contains three visualization modes: In the “METAGENOME” mode, microbial communities are visualized from the highest taxonomic unit to the lowest taxonomic unit, a feature which is inspired by tools such as Krona and MEGAN to allow to analyze the abundances of each bacterium and to visualize the abundance over time. In addition, if the data is provided, a user can select “DEG categories”, where an overview of transitions of up- and down-regulated genes between consecutive time points is given, whereas in the mode “GENES” all genes are visualized and can be selected.

3.1 General use of BioSankey

BioSankey can visualize input data from two different data sources: (i) microbial data and (ii) gene expression data. In general, the tool can be used for any other dataset that contains quantitative data from different time points. When integrated into BioSankey, it is possible to infer the fluctuations in abundances in microbial species originating from metagenome projects. These different input data can be then used for generating a project-specific website by making use of JavaScript and the Google API to enable an intuitive and interactive selection and visualizations of all data elements. Furthermore, filter criteria are enabled by functional descriptions when genes instead of microbial species were integrated. The whole webpage can be exchanged with collaboration partners because no web server or database is required. Figure 1 provides the workflow of the tool where either the gene expression or microbial abundance can be integrated. When the webpage is generated, we offer an export function for the tool to generate Sankey diagrams in the PDF format.

3.2 Comparability to Other Tools

We have compared BioSankey to other tools, which allow to visualize the abundances of species or taxa, such as Krona [1] and iTOL [3]. An overview of the functionality that BioSankey provides in comparison to the two other tools is given in Table 1. While BioSankey allows to visualize taxonomical and time-series data, Krona and iTOL do not allow to visualize time-series data but offer a broad range of export functionalities. The advantage of BioSankey is to visualize the expression of single genes or selected species, which is so far not supported in the other tools.

3.3 Microbial Community Analysis

When metagenomic studies are considered, reads are often assembled with tools such as SSpace [26] and then enter a binning approach with help of tools such as Maxbin [27] or Concoct [28], which group contigs into Bins based on their tetranucleotide frequencies and additional intrinsic features. As an alternative, more often, the sequencing of the 16S rRNA gene is used to assess the abundance of the species in a metagenomics dataset. Thereby, sequences are clustered based on 97% sequence identity to form OTUs, which can be analyzed with resources such as the SILVA server [29] after being processed with Mothur [30] and QIIME [31]. To demonstrate BioSankey for analyzing microbial communities, we have used data generated with QIIME from [32], comprising time-series of microbial communities from different human tissues. The goal of this study was to get insights about variability in these tissues and to define potential core microbiomes. The authors of this paper use line- and pie-charts, scatter-plots, principal component analysis (PCA) and area charts, the latter two especially for showing time variation. PCAs over time were visualized in a video. In addition to this extensive collection of visualization, we added two alternative visualization options. First, we selected the tongue tissue and extracted all information on genus level and used BioSankey to generate a project-specific HTML site. Out of 373 genera, 250 of them had support by at least one read and were used for BioSankey. This is shown in Figure 2A. In comparison, we made also a Krona plot for the same data (Figure 2B). With the Sankey diagrams however, we can also visualize the changes of the microbiome over time in one diagram, while Krona plots have the advantage, that broader or very detailed taxonomic hierarchies can be depicted.

Figure 2:

Illustration of the taxonomic distribution of the most abundant genera from the tongue tissue as described in Caporaso et al. [32]. On the left side (subfigure A) is the diagram generated by BioSankey, on the right side (subfigure B) is the diagram generated by the tool Krona (Ondov et al. [1]. The information in both diagrams is the same, the difference – except style, is that the Sankey diagram is read left to right from general to more specific (= to lower taxonomic ranks) and the Krona diagram from inside to outside.

Further, we used the first six time points on genus level of the tongue microbial inhabitants for BioSankey (Figure 3). This is a partly replacement of the figure on genus level in Additional File 10 in the paper of Caporaso et al. [32]. The Sankey plot visualization contains more information (the labeling) and is more appealing to the reader, with the drawback of being infeasible to show all 396 time points at once. In this case, the BioSankey plots might be the right choice to zoom in or show a selection of up to 20 time-points or higher, depending on the used computer screen or resolution.

Figure 3:

Illustration of four timepoints from the microbial community of the tongue as described in Caporaso et al. [32]. The taxonomy level is genus. The largest absolute genera and the biggest absolute and relative changes over time are visible: Prevotella as largest genus and Rothia, and Granulicatella with large changes (large bars, with small lines from left to right, e.g. “Rothia_2” and “Rothia_3”).

3.4 Differential Gene Expression Visualizations Embedded into BioSankey

To demonstrate BioSankey, we used published gene expression time-series data to describe the effects of Camptothecin in U87-MG cell lines by integrating up- and down-regulated genes over time [33]. Camptothecin is a drug, that specifically targets topoisomerase I. In this study, the authors showed the effects of Camptothecin to two glioblastoma cell lines (U87-MG and DBTRG-05). This is important to assess the use of this drug for malignant gliomas. With the gene expression in hand the authors could infer the affected pathways and assess the changes over time by considering expression data of six time points (2, 6, 16, 24, 48 and 72 h). By using BioSankey, we show, that Sankey plots can be used to visualize the reported ∼80% of down-regulated genes in U87-MG even in a time-resolved manner (Figure 4).

Figure 4:

Groups of differentially expressed genes between pairs of consecutive timepoints. The number of genes is the same for each timepoint. Genes, which are differentially expressed in only one time point are prefixed with “_LC” in all others (LC; low-change). This example is the response of U87-MG cell lines to Camptothecin which shows, that the effect is highest after 16 h and 72 h (the smallest amount of LC), the biggest change in the series is also from/to 16 h and 72 h and there are virtually no genes changing between up- and down-regulated (only the little line from “2 h_up” and “6 h_down”).

As a general overview of (mode: “DEG”) we provide a feature in BioSankey to highlight the amount of differential expressed genes as shown in Figure 4. The genes at each time point are filtered for up- and down-regulated expression (in the example a minimal fold-change of 2). In Figure 4, all transitions of genes between these states (up-, down- or low-regulated) are shown in each time point.

For a detailed analysis, the genes of a transition (e.g. all genes up-regulated at 2 h and at 6 h) can be selected and visualized separately. To obtain an overview visualization, a user must provide lists of DEG for each time point in a directory. From the visualization we can observe, that only a certain fraction of the genes (64 of 408) is up-regulated already at the first hours and only 41 of them are differentially expressed at the later time points. Furthermore, a user can select a time-point of genes that are up-, down- or not differentially expressed and can then extract the respective genes and visualize their expression to find a candidate gene or to obtain a general overview of the functionality of these genes.