mbQTL: an R/Bioconductor package for microbial quantitative trait loci (QTL) estimation.

MOTIVATION: In recent years, significant strides have been made in the field of genomics, with the commencement of large-scale studies aimed at collecting host mutational profiles and microbiome data. The amalgamation of host gene mutational profiles in both healthy and diseased subjects with microbial abundance data holds immense promise in providing insights into several crucial research questions, including the development and progression of diseases, as well as individual responses to therapeutic interventions. With the advent of sequencing methods such as 16s ribosomal RNA (rRNA) sequencing and whole genome sequencing, there is increasing evidence of interplay of human genetics and microbial communities. Quantitative trait loci associated with microbial abundance (mbQTLs), are genetic variants that influence the abundance of microbial populations within the host. RESULTS: Here, we introduce mbQTL, the first R package integrating 16S ribosomal RNA (rRNA) sequencing and single-nucleotide variation (SNV) and single-nucleotide polymorphism (SNP) data. We describe various statistical methods implemented for the identification of microbe-SNV pairs, relevant statistical measures, and plot functionality for interpretation. AVAILABILITY AND IMPLEMENTATION: mbQTL is available on bioconductor at https://bioconductor.org/packages/mbQTL/.

Neonatal Paenibacilliosis: Paenibacillus Infection as a Novel Cause of Sepsis in Term Neonates With High Risk of Sequelae in Uganda.

BACKGROUND: Paenibacillus thiaminolyticus may be an underdiagnosed cause of neonatal sepsis. METHODS: We prospectively enrolled a cohort of 800 full-term neonates presenting with a clinical diagnosis of sepsis at 2 Ugandan hospitals. Quantitative polymerase chain reaction specific to P. thiaminolyticus and to the Paenibacillus genus were performed on the blood and cerebrospinal fluid (CSF) of 631 neonates who had both specimen types available. Neonates with Paenibacillus genus or species detected in either specimen type were considered to potentially have paenibacilliosis, (37/631, 6%). We described antenatal, perinatal, and neonatal characteristics, presenting signs, and 12-month developmental outcomes for neonates with paenibacilliosis versus clinical sepsis due to other causes. RESULTS: Median age at presentation was 3 days (interquartile range 1, 7). Fever (92%), irritability (84%), and clinical signs of seizures (51%) were common. Eleven (30%) had an adverse outcome: 5 (14%) neonates died during the first year of life; 5 of 32 (16%) survivors developed postinfectious hydrocephalus (PIH) and 1 (3%) additional survivor had neurodevelopmental impairment without hydrocephalus. CONCLUSIONS: Paenibacillus species was identified in 6% of neonates with signs of sepsis who presented to 2 Ugandan referral hospitals; 70% were P. thiaminolyticus. Improved diagnostics for neonatal sepsis are urgently needed. Optimal antibiotic treatment for this infection is unknown but ampicillin and vancomycin will be ineffective in many cases. These results highlight the need to consider local pathogen prevalence and the possibility of unusual pathogens when determining antibiotic choice for neonatal sepsis.

GameRank: R package for feature selection and construction.

MOTIVATION: Building calibrated and discriminating predictive models can be developed through the direct optimization of model performance metrics with combinatorial search algorithms. Often, predictive algorithms are desired in clinical settings to identify patients that may be high and low risk. However, due to the large combinatorial search space, these algorithms are slow and do not guarantee the global optimality of their selection. RESULTS: Here, we present a novel and quick maximum likelihood-based feature selection algorithm, named GameRank. The method is implemented into an R package composed of additional functions to build calibrated and discriminative predictive models. AVAILABILITY AND IMPLEMENTATION: GameRank is available at https://github.com/Genentech/GameRank and released under the MIT License.

mirTarRnaSeq: An R/Bioconductor Statistical Package for miRNA-mRNA Target Identification and Interaction Analysis.

We introduce mirTarRnaSeq, an R/Bioconductor package for quantitative assessment of miRNA-mRNA relationships within sample cohorts. mirTarRnaSeq is a statistical package to explore predicted or pre-hypothesized miRNA-mRNA relationships following target prediction.We present two use cases applying mirTarRnaSeq. First, to identify miRNA targets, we examined EBV miRNAs for interaction with human and virus transcriptomes of stomach adenocarcinoma. This revealed enrichment of mRNA targets highly expressed in CD105+ endothelial cells, monocytes, CD4+ T cells, NK cells, CD19+ B cells, and CD34 cells. Next, to investigate miRNA-mRNA relationships in SARS-CoV-2 (COVID-19) infection across time, we used paired miRNA and RNA sequenced datasets of SARS-CoV-2 infected lung epithelial cells across three time points (4, 12, and 24 hours post-infection). mirTarRnaSeq identified evidence for human miRNAs targeting cytokine signaling and neutrophil regulation immune pathways from 4 to 24 hours after SARS-CoV-2 infection. Confirming the clinical relevance of these predictions, three of the immune specific mRNA-miRNA relationships identified in human lung epithelial cells after SARS-CoV-2 infection were also observed to be differentially expressed in blood from patients with COVID-19. Overall, mirTarRnaSeq is a robust tool that can address a wide-range of biological questions providing improved prediction of miRNA-mRNA interactions.

MicrobiomeExplorer: an R package for the analysis and visualization of microbial communities.

SUMMARY: We developed the MicrobiomeExplorer R package to facilitate the analysis and visualization of microbial communities. The MicrobiomeExplorer R package allows a user to perform typical microbiome analytic workflows and visualize their results, either through the command line or an interactive Shiny application included with the package. In addition to applying common analytical workflows, the application enables automated analysis report generation. AVAILABILITY AND IMPLEMENTATION: Available at https://github.com/zoecastillo/microbiomeExplorer. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Prognostic mutational subtyping in de novo diffuse large B-cell lymphoma.

BACKGROUND: Diffuse large B-cell lymphoma (DLBCL) is a heterogeneous disease defined using a number of well-established molecular subsets. Application of non-negative matrix factorization (NMF) to whole exome sequence data has previously been used to identify six distinct molecular clusters in DLBCL with potential clinical relevance. In this study, we applied NMF-clustering to targeted sequencing data utilizing the FoundationOne Heme® panel from the Phase III GOYA (NCT01287741) and Phase Ib/II CAVALLI studies (NCT02055820) in de novo DLBCL. Biopsy samples, survival outcomes, RNA-Seq and targeted exome-sequencing data were available for 423 patients in GOYA (obinutuzumab [G]-cyclophosphamide, doxorubicin, vincristine, and prednisone [CHOP] vs rituximab [R]-CHOP) and 86 patients in CAVALLI (venetoclax+[G/R]-CHOP). RESULTS: When the NMF algorithm was applied to samples from the GOYA study analyzed using a comprehensive genomic profiling platform, four of the six groups previously reported were observed: MYD88/CD79B, BCL2/EZH2, NOTCH2/TNFAIP3, and no mutations. Mutation profiles, cell-of-origin subset distributions and clinical associations of MYD88/CD79B and BCL2/EZH2 groups were similar to those described in previous NMF studies. In contrast, application of NMF to the CAVALLI study yielded only three; MYD88/CD79B-, BCL2/EZH2-like clusters, and a no mutations group, and there was a trend towards improved outcomes for BCL2/EZH2 over MYD88/CD79B. CONCLUSIONS: This analysis supports the utility of NMF used in conjunction with targeted sequencing platforms for identifying patients with different prognostic subsets. The observed trend for improved overall survival in the BCL2/EZH2 group is consistent with the mechanism of action of venetoclax, suggesting that targeting sequencing and NMF has potential for identifying patients who are more likely to gain benefit from venetoclax therapy.

Multivariable association discovery in population-scale meta-omics studies.

It is challenging to associate features such as human health outcomes, diet, environmental conditions, or other metadata to microbial community measurements, due in part to their quantitative properties. Microbiome multi-omics are typically noisy, sparse (zero-inflated), high-dimensional, extremely non-normal, and often in the form of count or compositional measurements. Here we introduce an optimized combination of novel and established methodology to assess multivariable association of microbial community features with complex metadata in population-scale observational studies. Our approach, MaAsLin 2 (Microbiome Multivariable Associations with Linear Models), uses generalized linear and mixed models to accommodate a wide variety of modern epidemiological studies, including cross-sectional and longitudinal designs, as well as a variety of data types (e.g., counts and relative abundances) with or without covariates and repeated measurements. To construct this method, we conducted a large-scale evaluation of a broad range of scenarios under which straightforward identification of meta-omics associations can be challenging. These simulation studies reveal that MaAsLin 2's linear model preserves statistical power in the presence of repeated measures and multiple covariates, while accounting for the nuances of meta-omics features and controlling false discovery. We also applied MaAsLin 2 to a microbial multi-omics dataset from the Integrative Human Microbiome (HMP2) project which, in addition to reproducing established results, revealed a unique, integrated landscape of inflammatory bowel diseases (IBD) across multiple time points and omics profiles.

metagenomeFeatures: an R package for working with 16S rRNA reference databases and marker-gene survey feature data.

SUMMARY: We developed the metagenomeFeatures R Bioconductor package along with annotation packages for three 16S rRNA databases (Greengenes, RDP and SILVA) to facilitate working with 16S rRNA databases and marker-gene survey feature data. The metagenomeFeatures package defines two classes, MgDb for working with 16S rRNA sequence databases, and mgFeatures for marker-gene survey feature data. The associated annotation packages provide a consistent interface to the different databases facilitating database comparison and exploration. The mgFeatures-class represents a crucial step in the development of a common data structure for working with 16S marker-gene survey data in R. AVAILABILITY AND IMPLEMENTATION: https://bioconductor.org/packages/release/bioc/html/metagenomeFeatures.html. SUPPLEMENTARY INFORMATION: Supplementary material is available at Bioinformatics online.

Prognostic impact of somatic mutations in diffuse large B-cell lymphoma and relationship to cell-of-origin: data from the phase III GOYA study.

Diffuse large B-cell lymphoma represents a biologically and clinically heterogeneous diagnostic category with well-defined cell-of-origin subtypes. Using data from the GOYA study (NCT01287741), we characterized the mutational profile of diffuse large B-cell lymphoma and evaluated the prognostic impact of somatic mutations in relation to cell-of-origin. Targeted DNA next-generation sequencing was performed in 499 formalin-fixed paraffin-embedded tissue biopsies from previously untreated patients. Prevalence of genetic alterations/mutations was examined. Multivariate Cox regression was used to evaluate the prognostic effect of individual genomic alterations. Of 465 genes analyzed, 59 were identified with mutations occurring in at least 10 of 499 patients (≥2% prevalence); 334 additional genes had mutations occurring in ≥1 patient. Single nucleotide variants were the most common mutation type. On multivariate analysis, BCL2 alterations were most strongly associated with shorter progression-free survival (multivariate hazard ratio: 2.6; 95% confidence interval: 1.6 to 4.2). BCL2 alterations were detected in 102 of 499 patients; 92 had BCL2 translocations, 90% of whom had germinal center B-cell-like diffuse large B-cell lymphoma. BCL2 alterations were also significantly correlated with BCL2 gene and protein expression levels. Validation of published mutational subsets revealed consistent patterns of co-occurrence, but no consistent prognostic differences between subsets. Our data confirm the molecular heterogeneity of diffuse large B-cell lymphoma, with potential treatment targets occurring in distinct cell-of-origin subtypes. clinicaltrials.gov identifier: NCT01287741.

Metaviz: interactive statistical and visual analysis of metagenomic data.

Large studies profiling microbial communities and their association with healthy or disease phenotypes are now commonplace. Processed data from many of these studies are publicly available but significant effort is required for users to effectively organize, explore and integrate it, limiting the utility of these rich data resources. Effective integrative and interactive visual and statistical tools to analyze many metagenomic samples can greatly increase the value of these data for researchers. We present Metaviz, a tool for interactive exploratory data analysis of annotated microbiome taxonomic community profiles derived from marker gene or whole metagenome shotgun sequencing. Metaviz is uniquely designed to address the challenge of browsing the hierarchical structure of metagenomic data features while rendering visualizations of data values that are dynamically updated in response to user navigation. We use Metaviz to provide the UMD Metagenome Browser web service, allowing users to browse and explore data for more than 7000 microbiomes from published studies. Users can also deploy Metaviz as a web service, or use it to analyze data through the metavizr package to interoperate with state-of-the-art analysis tools available through Bioconductor. Metaviz is free and open source with the code, documentation and tutorials publicly accessible.

Exploring regulation in tissues with eQTL networks.

Characterizing the collective regulatory impact of genetic variants on complex phenotypes is a major challenge in developing a genotype to phenotype map. Using expression quantitative trait locus (eQTL) analyses, we constructed bipartite networks in which edges represent significant associations between genetic variants and gene expression levels and found that the network structure informs regulatory function. We show, in 13 tissues, that these eQTL networks are organized into dense, highly modular communities grouping genes often involved in coherent biological processes. We find communities representing shared processes across tissues, as well as communities associated with tissue-specific processes that coalesce around variants in tissue-specific active chromatin regions. Node centrality is also highly informative, with the global and community hubs differing in regulatory potential and likelihood of being disease associated.

Smooth quantile normalization.

Between-sample normalization is a critical step in genomic data analysis to remove systematic bias and unwanted technical variation in high-throughput data. Global normalization methods are based on the assumption that observed variability in global properties is due to technical reasons and are unrelated to the biology of interest. For example, some methods correct for differences in sequencing read counts by scaling features to have similar median values across samples, but these fail to reduce other forms of unwanted technical variation. Methods such as quantile normalization transform the statistical distributions across samples to be the same and assume global differences in the distribution are induced by only technical variation. However, it remains unclear how to proceed with normalization if these assumptions are violated, for example, if there are global differences in the statistical distributions between biological conditions or groups, and external information, such as negative or control features, is not available. Here, we introduce a generalization of quantile normalization, referred to as smooth quantile normalization (qsmooth), which is based on the assumption that the statistical distribution of each sample should be the same (or have the same distributional shape) within biological groups or conditions, but allowing that they may differ between groups. We illustrate the advantages of our method on several high-throughput datasets with global differences in distributions corresponding to different biological conditions. We also perform a Monte Carlo simulation study to illustrate the bias-variance tradeoff and root mean squared error of qsmooth compared to other global normalization methods. A software implementation is available from https://github.com/stephaniehicks/qsmooth.

Estimating gene regulatory networks with pandaR.

CONTACT: johnq@jimmy.harvard.edu or dschlauch@fas.harvard.edu. AVAILABILITY AND IMPLEMENTATION: PandaR is provided as a Bioconductor R Package and is available at bioconductor.org/packages/pandaR.

Tissue-aware RNA-Seq processing and normalization for heterogeneous and sparse data.

BACKGROUND: Although ultrahigh-throughput RNA-Sequencing has become the dominant technology for genome-wide transcriptional profiling, the vast majority of RNA-Seq studies typically profile only tens of samples, and most analytical pipelines are optimized for these smaller studies. However, projects are generating ever-larger data sets comprising RNA-Seq data from hundreds or thousands of samples, often collected at multiple centers and from diverse tissues. These complex data sets present significant analytical challenges due to batch and tissue effects, but provide the opportunity to revisit the assumptions and methods that we use to preprocess, normalize, and filter RNA-Seq data - critical first steps for any subsequent analysis. RESULTS: We find that analysis of large RNA-Seq data sets requires both careful quality control and the need to account for sparsity due to the heterogeneity intrinsic in multi-group studies. We developed Yet Another RNA Normalization software pipeline (YARN), that includes quality control and preprocessing, gene filtering, and normalization steps designed to facilitate downstream analysis of large, heterogeneous RNA-Seq data sets and we demonstrate its use with data from the Genotype-Tissue Expression (GTEx) project. CONCLUSIONS: An R package instantiating YARN is available at http://bioconductor.org/packages/yarn .

Regulatory network changes between cell lines and their tissues of origin.

BACKGROUND: Cell lines are an indispensable tool in biomedical research and often used as surrogates for tissues. Although there are recognized important cellular and transcriptomic differences between cell lines and tissues, a systematic overview of the differences between the regulatory processes of a cell line and those of its tissue of origin has not been conducted. The RNA-Seq data generated by the GTEx project is the first available data resource in which it is possible to perform a large-scale transcriptional and regulatory network analysis comparing cell lines with their tissues of origin. RESULTS: We compared 127 paired Epstein-Barr virus transformed lymphoblastoid cell lines (LCLs) and whole blood samples, and 244 paired primary fibroblast cell lines and skin samples. While gene expression analysis confirms that these cell lines carry the expression signatures of their primary tissues, albeit at reduced levels, network analysis indicates that expression changes are the cumulative result of many previously unreported alterations in transcription factor (TF) regulation. More specifically, cell cycle genes are over-expressed in cell lines compared to primary tissues, and this alteration in expression is a result of less repressive TF targeting. We confirmed these regulatory changes for four TFs, including SMAD5, using independent ChIP-seq data from ENCODE. CONCLUSIONS: Our results provide novel insights into the regulatory mechanisms controlling the expression differences between cell lines and tissues. The strong changes in TF regulation that we observe suggest that network changes, in addition to transcriptional levels, should be considered when using cell lines as models for tissues.

Privacy-preserving microbiome analysis using secure computation.

MOTIVATION: Developing targeted therapeutics and identifying biomarkers relies on large amounts of research participant data. Beyond human DNA, scientists now investigate the DNA of micro-organisms inhabiting the human body. Recent work shows that an individual's collection of microbial DNA consistently identifies that person and could be used to link a real-world identity to a sensitive attribute in a research dataset. Unfortunately, the current suite of DNA-specific privacy-preserving analysis tools does not meet the requirements for microbiome sequencing studies. RESULTS: To address privacy concerns around microbiome sequencing, we implement metagenomic analyses using secure computation. Our implementation allows comparative analysis over combined data without revealing the feature counts for any individual sample. We focus on three analyses and perform an evaluation on datasets currently used by the microbiome research community. We use our implementation to simulate sharing data between four policy-domains. Additionally, we describe an application of our implementation for patients to combine data that allows drug developers to query against and compensate patients for the analysis. AVAILABILITY AND IMPLEMENTATION: The software is freely available for download at: http://cbcb.umd.edu/∼hcorrada/projects/secureseq.html SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online. CONTACT: hcorrada@umiacs.umd.edu.

Individual-specific changes in the human gut microbiota after challenge with enterotoxigenic Escherichia coli and subsequent ciprofloxacin treatment.

BACKGROUND: Enterotoxigenic Escherichia coli (ETEC) is a major cause of diarrhea in inhabitants from low-income countries and in visitors to these countries. The impact of the human intestinal microbiota on the initiation and progression of ETEC diarrhea is not yet well understood. RESULTS: We used 16S rRNA (ribosomal RNA) gene sequencing to study changes in the fecal microbiota of 12 volunteers during a human challenge study with ETEC (H10407) and subsequent treatment with ciprofloxacin. Five subjects developed severe diarrhea and seven experienced few or no symptoms. Diarrheal symptoms were associated with high concentrations of fecal E. coli as measured by quantitative culture, quantitative PCR, and normalized number of 16S rRNA gene sequences. Large changes in other members of the microbiota varied greatly from individual to individual, whether or not diarrhea occurred. Nonetheless the variation within an individual was small compared to variation between individuals. Ciprofloxacin treatment reorganized microbiota populations; however, the original structure was largely restored at one and three month follow-up visits. CONCLUSION: Symptomatic ETEC infections, but not asymptomatic infections, were associated with high fecal concentrations of E. coli. Both infection and ciprofloxacin treatment caused variable changes in other bacteria that generally reverted to baseline levels after three months.

Differential abundance analysis for microbial marker-gene surveys.

We introduce a methodology to assess differential abundance in sparse high-throughput microbial marker-gene survey data. Our approach, implemented in the metagenomeSeq Bioconductor package, relies on a novel normalization technique and a statistical model that accounts for undersampling-a common feature of large-scale marker-gene studies. Using simulated data and several published microbiota data sets, we show that metagenomeSeq outperforms the tools currently used in this field.

Microbiota that affect risk for shigellosis in children in low-income countries.

Pathogens in the gastrointestinal tract exist within a vast population of microbes. We examined associations between pathogens and composition of gut microbiota as they relate to Shigella spp./enteroinvasive Escherichia coli infection. We analyzed 3,035 stool specimens (1,735 nondiarrheal and 1,300 moderate-to-severe diarrheal) from the Global Enteric Multicenter Study for 9 enteropathogens. Diarrheal specimens had a higher number of enteropathogens (diarrheal mean 1.4, nondiarrheal mean 0.95; p<0.0001). Rotavirus showed a negative association with Shigella spp. in cases of diarrhea (odds ratio 0.31, 95% CI 0.17-0.55) and had a large combined effect on moderate-to-severe diarrhea (odds ratio 29, 95% CI 3.8-220). In 4 Lactobacillus taxa identified by 16S rRNA gene sequencing, the association between pathogen and disease was decreased, which is consistent with the possibility that Lactobacillus spp. are protective against Shigella spp.-induced diarrhea. Bacterial diversity of gut microbiota was associated with diarrhea status, not high levels of the Shigella spp. ipaH gene.

bacLIFE: a user-friendly computational workflow for genome analysis and prediction of lifestyle-associated genes in bacteria.

Bacteria have an extensive adaptive ability to live in close association with eukaryotic hosts, exhibiting detrimental, neutral or beneficial effects on host growth and health. However, the genes involved in niche adaptation are mostly unknown and their functions poorly characterized. Here, we present bacLIFE ( https://github.com/Carrion-lab/bacLIFE ) a streamlined computational workflow for genome annotation, large-scale comparative genomics, and prediction of lifestyle-associated genes (LAGs). As a proof of concept, we analyzed 16,846 genomes from the Burkholderia/Paraburkholderia and Pseudomonas genera, which led to the identification of hundreds of genes potentially associated with a plant pathogenic lifestyle. Site-directed mutagenesis of 14 of these predicted LAGs of unknown function, followed by plant bioassays, showed that 6 predicted LAGs are indeed involved in the phytopathogenic lifestyle of Burkholderia plantarii and Pseudomonas syringae pv. phaseolicola. These 6 LAGs encompassed a glycosyltransferase, extracellular binding proteins, homoserine dehydrogenases and hypothetical proteins. Collectively, our results highlight bacLIFE as an effective computational tool for prediction of LAGs and the generation of hypotheses for a better understanding of bacteria-host interactions.