A general modular framework for gene set enrichment analysis

Historic breeding practices contribute to germplasm divergence in leaf specialized metabolism and ecophysiology in cultivated sunflower (Helianthus annuus)

PREMISE

The use of hybrid breeding systems to increase crop yields has been the cornerstone of modern agriculture and is exemplified in the breeding and improvement of cultivated sunflower (Helianthus annuus). However, it is poorly understood what effect supporting separate breeding pools in such systems, combined with continued selection for yield, may have on leaf ecophysiology and specialized metabolite variation.

METHODS

We analyzed 288 lines of cultivated H. annuus to examine the genomic basis of several specialized metabolites and agronomically important traits across major heterotic groups.

RESULTS

Heterotic group identity supports phenotypic divergences between fertility restoring and cytoplasmic male-sterility maintainer lines in leaf ecophysiology and specialized metabolism. However, the divergence is not associated with physical linkage to nuclear genes that support current hybrid breeding practices in cultivated H. annuus. Additionally, we identified four genomic regions associated with leaf ecophysiology and specialized metabolism that colocalize with previously identified QTLs for quantitative self-compatibility traits and with S-protein homolog (SPH) proteins, a recently discovered family of proteins associated with self-incompatibility and self/nonself recognition in Papaver rhoeas (common poppy) with suggested conserved downstream mechanisms among eudicots.

CONCLUSIONS

Further work is necessary to confirm the self-incompatibility mechanisms in cultivated H. annuus and their relationship to the integrative and polygenic architecture of leaf ecophysiology and specialized metabolism in cultivated sunflower. However, because self-compatibility is a derived quantitative trait in cultivated H. annuus, trait linkage to divergent phenotypic traits may have partially arisen as a potential unintended consequence of historical breeding practices and selection for yield.

Regulatory Changes in the Fatty Acid Elongase eloF Underlie the Evolution of Sex-specific Pheromone Profiles in Drosophila prolongata

Pheromones play a key role in regulating sexual behavior throughout the animal kingdom. In Drosophila and other insects, many cuticular hydrocarbons (CHCs) are sexually dimorphic, and some are known to perform pheromonal functions. However, the genetic control of sex-specific CHC production is not understood outside of the model species D. melanogaster. A recent evolutionary change is found in D. prolongata, which, compared to its closest relatives, shows greatly increased sexual dimorphism in both CHCs and the chemosensory system responsible for their perception. A key transition involves a male-specific increase in the proportion of long-chain CHCs. Perfuming D. prolongata females with the male-biased CHCs reduces copulation success, suggesting that these compounds function as sex pheromones. The evolutionary change in CHC profiles correlates with a male-specific increase in the expression of multiple genes involved in CHC biosynthesis, including fatty acid elongases and reductases and other key enzymes. In particular, elongase F, which is responsible for producing female-specific pheromones in D. melanogaster, is strongly upregulated in D. prolongata males compared both to females and to males of the sibling species. Induced mutations in eloF reduce the amount of long-chain CHCs, resulting in a partial feminization of pheromone profiles in D. prolongata males while having minimal effect in females. Transgenic experiments show that sex-biased expression of eloF is caused in part by a putative transposable element insertion in its regulatory region. These results reveal one of the genetic mechanisms responsible for a recent evolutionary change in sexual communication.

ReactomeGSA: new features to simplify public data reuse

MOTIVATION

ReactomeGSA is part of the Reactome knowledgebase and one of the leading multi-omics pathway analysis platforms. ReactomeGSA provides access to quantitative pathway analysis methods supporting different 'omics data types. Additionally, ReactomeGSA can process different datasets simultaneously, leading to a comparative pathway analysis that can also be performed across different species.

RESULTS

We present a major update to the ReactomeGSA analysis platforms that greatly simplifies the reuse and direct integration of public data. In order to increase the number of available datasets, we developed the new grein_loader Python application that can directly fetch experiments from the GREIN resource. This enabled us to support both EMBL-EBI's Expression Atlas and GEO RNA-seq Experiments Interactive Navigator within ReactomeGSA. To further increase the visibility and simplify the reuse of public datasets, we integrated a novel search function into ReactomeGSA that enables users to search for public datasets across both supported resources. Finally, we completely re-developed ReactomeGSA's web-frontend and R/Bioconductor package to support the new search and loading features, and greatly simplify the use of ReactomeGSA.

AVAILABILITY AND IMPLEMENTATION

The new ReactomeGSA web frontend is available at https://www.reactome.org/gsa with an built-in, interactive tutorial. The ReactomeGSA R package (https://bioconductor.org/packages/release/bioc/html/ReactomeGSA.html) is available through Bioconductor and shipped with detailed documentation and vignettes. The grein_loader Python application is available through the Python Package Index (pypi). The complete source code for all applications is available on GitHub at https://github.com/grisslab/grein_loader and https://github.com/reactome.

Unraveling the molecular relevance of brain phenotypes: A comparative analysis of null models and test statistics

Correlating transcriptional profiles with imaging-derived phenotypes has the potential to reveal possible molecular architectures associated with cognitive functions, brain development and disorders. Competitive null models built by resampling genes and self-contained null models built by spinning brain regions, along with varying test statistics, have been used to determine the significance of transcriptional associations. However, there has been no systematic evaluation of their performance in imaging transcriptomics analyses. Here, we evaluated the performance of eight different test statistics (mean, mean absolute value, mean squared value, max mean, median, Kolmogorov-Smirnov (KS), Weighted KS and the number of significant correlations) in both competitive null models and self-contained null models. Simulated brain maps (n = 1,000) and gene sets (n = 500) were used to calculate the probability of significance (Psig) for each statistical test. Our results suggested that competitive null models may result in false positive results driven by co-expression within gene sets. Furthermore, we demonstrated that the self-contained null models may fail to account for distribution characteristics (e.g., bimodality) of correlations between all available genes and brain phenotypes, leading to false positives. These two confounding factors interacted differently with test statistics, resulting in varying outcomes. Specifically, the sign-sensitive test statistics (i.e., mean, median, KS, Weighted KS) were influenced by co-expression bias in the competitive null models, while median and sign-insensitive test statistics were sensitive to the bimodality bias in the self-contained null models. Additionally, KS-based statistics produced conservative results in the self-contained null models, which increased the risk of false negatives. Comprehensive supplementary analyses with various configurations, including realistic scenarios, supported the results. These findings suggest utilizing sign-insensitive test statistics such as mean absolute value, max mean in the competitive null models and the mean as the test statistic for the self-contained null models. Additionally, adopting the confounder-matched (e.g., coexpression-matched) null models as an alternative to standard null models can be a viable strategy. Overall, the present study offers insights into the selection of statistical tests for imaging transcriptomics studies, highlighting areas for further investigation and refinement in the evaluation of novel and commonly used tests.

bootGSEA: a bootstrap and rank aggregation pipeline for multi-study and multi-omics enrichment analyses

Introduction: Gene set enrichment analysis (GSEA) subsequent to differential expression analysis is a standard step in transcriptomics and proteomics data analysis. Although many tools for this step are available, the results are often difficult to reproduce because set annotations can change in the databases, that is, new features can be added or existing features can be removed. Finally, such changes in set compositions can have an impact on biological interpretation. Methods: We present bootGSEA, a novel computational pipeline, to study the robustness of GSEA. By repeating GSEA based on bootstrap samples, the variability and robustness of results can be studied. In our pipeline, not all genes or proteins are involved in the different bootstrap replicates of the analyses. Finally, we aggregate the ranks from the bootstrap replicates to obtain a score per gene set that shows whether it gains or loses evidence compared to the ranking of the standard GSEA. Rank aggregation is also used to combine GSEA results from different omics levels or from multiple independent studies at the same omics level. Results: By applying our approach to six independent cancer transcriptomics datasets, we showed that bootstrap GSEA can aid in the selection of more robust enriched gene sets. Additionally, we applied our approach to paired transcriptomics and proteomics data obtained from a mouse model of spinal muscular atrophy (SMA), a neurodegenerative and neurodevelopmental disease associated with multi-system involvement. After obtaining a robust ranking at both omics levels, both ranking lists were combined to aggregate the findings from the transcriptomics and proteomics results. Furthermore, we constructed the new R-package "bootGSEA," which implements the proposed methods and provides graphical views of the findings. Bootstrap-based GSEA was able in the example datasets to identify gene or protein sets that were less robust when the set composition changed during bootstrap analysis. Discussion: The rank aggregation step was useful for combining bootstrap results and making them comparable to the original findings on the single-omics level or for combining findings from multiple different omics levels.

Blood and urine multi-omics analysis of the impact of e-vaping, smoking, and cessation: from exposome to molecular responses

Cigarette smoking is a major preventable cause of morbidity and mortality. While quitting smoking is the best option, switching from cigarettes to non-combustible alternatives (NCAs) such as e-vapor products is a viable harm reduction approach for smokers who would otherwise continue to smoke. A key challenge for the clinical assessment of NCAs is that self-reported product use can be unreliable, compromising the proper evaluation of their risk reduction potential. In this cross-sectional study of 205 healthy volunteers, we combined comprehensive exposure characterization with in-depth multi-omics profiling to compare effects across four study groups: cigarette smokers (CS), e-vapor users (EV), former smokers (FS), and never smokers (NS). Multi-omics analyses included metabolomics, transcriptomics, DNA methylomics, proteomics, and lipidomics. Comparison of the molecular effects between CS and NS recapitulated several previous observations, such as increased inflammatory markers in CS. Generally, FS and EV demonstrated intermediate molecular effects between the NS and CS groups. Stratification of the FS and EV by combustion exposure markers suggested that this position on the spectrum between CS and NS was partially driven by non-compliance/dual use. Overall, this study highlights the importance of in-depth exposure characterization before biological effect characterization for any NCA assessment study.

Morphological connectivity differences in Alzheimer's disease correlate with gene transcription and cell-type

Alzheimer's disease (AD) is one of the most prevalent forms of dementia in older individuals. Convergent evidence suggests structural connectome abnormalities in specific brain regions are linked to AD progression. The biological basis underpinnings of these connectome changes, however, have remained elusive. We utilized an individual regional mean connectivity strength (RMCS) derived from a regional radiomics similarity network to capture altered morphological connectivity in 1654 participants (605 normal controls, 766 mild cognitive impairment [MCI], and 283 AD). Then, we also explored the biological basis behind these morphological changes through gene enrichment analysis and cell-specific analysis. We found that RMCS probes of the hippocampus and medial temporal lobe were significantly altered in AD and MCI, with these differences being spatially related to the expression of AD-risk genes. In addition, gene enrichment analysis revealed that the modulation of chemical synaptic transmission is the most relevant biological process associated with the altered RMCS in AD. Notably, neuronal cells were found to be the most pertinent cells in the altered RMCS. Our findings shed light on understanding the biological basis of structural connectome changes in AD, which may ultimately lead to more effective diagnostic and therapeutic strategies for this devastating disease.

Thymopentin ameliorates experimental colitis via inhibiting neutrophil extracellular traps

BACKGROUND

The long-term prognosis of Crohn's disease (CD) remains unsatisfactory. Therefore, we assessed the therapeutic effect of thymopentin (TP5) in a mouse model of 2,4,6-trinitrobenzene sulfonic acid (TNBS)-induced colitis, which mimics CD, and analyzed its impact on neutrophil extracellular traps (NETs).

METHODS

NET markers, including myeloperoxidase (MPO), neutrophil elastase (NE), citrullinated histone H3 (CitH3), peptidyl arginine deiminase IV (PAD4), and double-stranded DNA (dsDNA) were assessed by immunostaining and enzyme-linked immunosorbent assay. NET formation was evaluated in vitro. Neoseptin 3, a specific NET agonist, was used to reverse the effect of TP5 on TNBS-induced colitis. The action mechanism of TP5 was investigated using RNA-seq.

RESULTS

TP5 ameliorated weight loss (P < 0.001), disease activity index (DAI) (P = 0.05), colon shrinkage (P = 0.04), and elevated levels of tumor necrosis factor-alpha (TNF-α), interleukin (IL)-1β, IL-6, and neutrophils in the TNBS group. The TNBS group exhibited increased MPO, NE, CitH3, PAD4, dsDNA and MPO-DNA levels (all P < 0.001), which decreased after TP5 administration (P = 0.01, P < 0.001, P < 0.001, P < 0.001, P = 0.02, and P = 0.02 respectively). Tissue CitH3 levels were positively correlated with DAI and TNF-α levels (P < 0.05). Furthermore, phorbol 12-myristate 13-acetate-stimulated NET formation increased by 1.8-, 2.8-, and 2.3-fold in vitro in the control, TNBS + saline, and TNBS + TP5 groups, respectively. Neoseptin 3 significantly reversed the effect of TP5. RNA-seq revealed potential pathways underlying the effect of TP5.

CONCLUSION

TP5 effectively ameliorated colitis by suppressing NETs in the experimental CD model.

Detection of clinically important BRCA gene mutations in ovarian cancer patients using next generation sequencing analysis

Ovarian cancer, a complex and aggressive malignancy, remains a significant challenge in clinical oncology due to its heterogeneous nature and limited therapeutic options. In this study, across Pakistani ovarian cancer patients, we conducted a comprehensive analysis of mutations within the BRCA1 and BRCA2 genes to elucidate their potential implications in ovarian cancer susceptibility and progression. Employing Next-Generation Sequencing (NGS), we conducted a comprehensive mutational analysis of BRCA1/2 genes. Kaplan Meier analysis was used to analyze the effect of pathogenic mutations on the survival outcomes of ovarian cancer patients. Reverse transcription-quantitative polymerase chain reaction (RT-qPCR) and Immunohistochemistry (IHC) analyses were conducted to analyze the downstream effect of the pathogenic mutations. Targeted bisulfite sequencing (bisulfite-seq) analysis facilitated the investigation of epigenetic contributions to gene expression regulation. Enrichment analysis was conducted to uncover significant Gene Ontology (GO) terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways associated with BRCA1/2. Exploring DrugBank, we identified potential drugs capable of modulating BRCA1/2 expression regulation. NGS analysis identified three clinically significant pathogenic mutations within the BRCA1 gene and two within the BRCA2 gene, shedding light on their potential involvement in ovarian cancer susceptibility and progression. Kaplan Meier analysis unveiled poor overall survival (OS) associated with the identified pathogenic mutations, accentuating their prognostic value. Expression analysis using reverse transcription-quantitative polymerase chain reaction (RT-qPCR) and IHC demonstrated a significant up-regulation of BRCA1 and BRCA2 genes in ovarian cancer samples harboring pathogenic mutations. Bisulfite-seq revealed a significant hypomethylation within promoter regions of mutated BRCA1 and BRCA2 genes in ovarian cancer samples, compared to non-mutated cases with pathogenic mutations, indicating the role of epigenetics in expression dysregulation as well. By uncovering clinically significant pathogenic mutations in BRCA1/2 genes and establishing their link with up-regulated gene expression, this study significantly advances our understanding of ovarian cancer's underlying causes in the Pakistani population.

Bayesian Statistics Improves Biological Interpretability of Metabolomics Data from Human Cohorts

Univariate analyses of metabolomics data currently follow a frequentist approach, using p-values to reject a null hypothesis. We here propose the use of Bayesian statistics to quantify evidence supporting different hypotheses and discriminate between the null hypothesis versus the lack of statistical power. We used metabolomics data from three independent human cohorts that studied the plasma signatures of subjects with myalgic encephalomyelitis/chronic fatigue syndrome (ME/CFS). The data are publicly available, covering 84-197 subjects in each study with 562-888 identified metabolites of which 777 were common between the two studies and 93 were compounds reported in all three studies. We show how Bayesian statistics incorporates results from one study as "prior information" into the next study, thereby improving the overall assessment of the likelihood of finding specific differences between plasma metabolite levels. Using classic statistics and Benjamini-Hochberg FDR-corrections, Study 1 detected 18 metabolic differences and Study 2 detected no differences. Using Bayesian statistics on the same data, we found a high likelihood that 97 compounds were altered in concentration in Study 2, after using the results of Study 1 as the prior distributions. These findings included lower levels of peroxisome-produced ether-lipids, higher levels of long-chain unsaturated triacylglycerides, and the presence of exposome compounds that are explained by the difference in diet and medication between healthy subjects and ME/CFS patients. Although Study 3 reported only 92 compounds in common with the other two studies, these major differences were confirmed. We also found that prostaglandin F2alpha, a lipid mediator of physiological relevance, was reduced in ME/CFS patients across all three studies. The use of Bayesian statistics led to biological conclusions from metabolomic data that were not found through frequentist approaches. We propose that Bayesian statistics is highly useful for studies with similar research designs if similar metabolomic assays are used.

Transcriptome analysis reveals key regulatory genes for root growth related to potassium utilization efficiency in rapeseed (Brassica napus L.)

Root system architecture (RSA) is the primary predictor of nutrient intake and significantly influences potassium utilization efficiency (KUE). Uncertainty persists regarding the genetic factors governing root growth in rapeseed. The root transcriptome analysis reveals the genetic basis driving crop root growth. In this study, RNA-seq was used to profile the overall transcriptome in the root tissue of 20 Brassica napus accessions with high and low KUE. 71,437 genes in the roots displayed variable expression profiles between the two contrasting genotype groups. The 212 genes that had varied expression levels between the high and low KUE lines were found using a pairwise comparison approach. The Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) functional classification analysis revealed that the DEGs implicated in hormone and signaling pathways, as well as glucose, lipid, and amino acid metabolism, were all differently regulated in the rapeseed root system. Additionally, we discovered 33 transcription factors (TFs) that control root development were differentially expressed. By combining differential expression analysis, weighted gene co-expression network analysis (WGCNA), and recent genome-wide association study (GWAS) results, four candidate genes were identified as essential hub genes. These potential genes were located fewer than 100 kb from the peak SNPs of QTL clusters, and it was hypothesized that they regulated the formation of the root system. Three of the four hub genes' homologs-BnaC04G0560400ZS, BnaC04G0560400ZS, and BnaA03G0073500ZS-have been shown to control root development in earlier research. The information produced by our transcriptome profiling could be useful in revealing the molecular processes involved in the growth of rapeseed roots in response to KUE.

Assessing Kidney Injury Induced by Mercuric Chloride in Guinea Pigs with In Vivo and In Vitro Experiments

Acute kidney injury, which is associated with high levels of morbidity and mortality, affects a significant number of individuals, and can be triggered by multiple factors, such as medications, exposure to toxic chemicals or other substances, disease, and trauma. Because the kidney is a critical organ, understanding and identifying early cellular or gene-level changes can provide a foundation for designing medical interventions. In our earlier work, we identified gene modules anchored to histopathology phenotypes associated with toxicant-induced liver and kidney injuries. Here, using in vivo and in vitro experiments, we assessed and validated these kidney injury-associated modules by analyzing gene expression data from the kidneys of male Hartley guinea pigs exposed to mercuric chloride. Using plasma creatinine levels and cell-viability assays as measures of the extent of renal dysfunction under in vivo and in vitro conditions, we performed an initial range-finding study to identify the appropriate doses and exposure times associated with mild and severe kidney injuries. We then monitored changes in kidney gene expression at the selected doses and time points post-toxicant exposure to characterize the mechanisms of kidney injury. Our injury module-based analysis revealed a dose-dependent activation of several phenotypic cellular processes associated with dilatation, necrosis, and fibrogenesis that were common across the experimental platforms and indicative of processes that initiate kidney damage. Furthermore, a comparison of activated injury modules between guinea pigs and rats indicated a strong correlation between the modules, highlighting their potential for cross-species translational studies.

Cortical profiles of numerous psychiatric disorders and normal development share a common pattern

The neurobiological bases of the association between development and psychopathology remain poorly understood. Here, we identify a shared spatial pattern of cortical thickness (CT) in normative development and several psychiatric and neurological disorders. Principal component analysis (PCA) was applied to CT of 68 regions in the Desikan-Killiany atlas derived from three large-scale datasets comprising a total of 41,075 neurotypical participants. PCA produced a spatially broad first principal component (PC1) that was reproducible across datasets. Then PC1 derived from healthy adult participants was compared to the pattern of CT differences associated with psychiatric and neurological disorders comprising a total of 14,886 cases and 20,962 controls from seven ENIGMA disease-related working groups, normative maturation and aging comprising a total of 17,697 scans from the ABCD Study® and the IMAGEN developmental study, and 17,075 participants from the ENIGMA Lifespan working group, as well as gene expression maps from the Allen Human Brain Atlas. Results revealed substantial spatial correspondences between PC1 and widespread lower CT observed in numerous psychiatric disorders. Moreover, the PC1 pattern was also correlated with the spatial pattern of normative maturation and aging. The transcriptional analysis identified a set of genes including KCNA2, KCNS1 and KCNS2 with expression patterns closely related to the spatial pattern of PC1. The gene category enrichment analysis indicated that the transcriptional correlations of PC1 were enriched to multiple gene ontology categories and were specifically over-represented starting at late childhood, coinciding with the onset of significant cortical maturation and emergence of psychopathology during the prepubertal-to-pubertal transition. Collectively, the present study reports a reproducible latent pattern of CT that captures interregional profiles of cortical changes in both normative brain maturation and a spectrum of psychiatric disorders. The pubertal timing of the expression of PC1-related genes implicates disrupted neurodevelopment in the pathogenesis of the spectrum of psychiatric diseases emerging during adolescence.

MAPKAPK2-centric transcriptome profiling reveals its major role in governing molecular crosstalk of IGFBP2, MUC4, and PRKAR2B during HNSCC pathogenesis

Transcriptome analysis of head and neck squamous cell carcinoma (HNSCC) has been pivotal to comprehending the convoluted biology of HNSCC tumors. MAPKAPK2 or MK2 is a critical modulator of the mRNA turnover of crucial genes involved in HNSCC progression. However, MK2-centric transcriptome profiles of tumors are not well known. This study delves into HNSCC progression with MK2 at the nexus to delineate the biological relevance and intricate crosstalk of MK2 in the tumor milieu. We performed next-generation sequencing-based transcriptome profiling of HNSCC cells and xenograft tumors to ascertain mRNA expression profiles in MK2-wild type and MK2-knockdown conditions. The findings were validated using gene expression assays, immunohistochemistry, and transcript turnover studies. Here, we identified a pool of crucial MK2-regulated candidate genes by annotation and differential gene expression analyses. Regulatory network and pathway enrichment revealed their significance and involvement in the HNSCC pathogenesis. Additionally, 3'-UTR-based filtering recognized important MK2-regulated downstream target genes and validated them by nCounter gene expression assays. Finally, immunohistochemistry and transcript stability studies revealed the putative role of MK2 in regulating the transcript turnover of IGFBP2, MUC4, and PRKAR2B in HNSCC. Conclusively, MK2-regulated candidate genes were identified in this study, and their plausible involvement in HNSCC pathogenesis was elucidated. These genes possess investigative values as targets for diagnosis and therapeutic interventions for HNSCC.

Convergent genomics of longevity in rockfishes highlights the genetics of human life span variation

Longevity is a defining, heritable trait that varies dramatically between species. To resolve the genetic regulation of this trait, we have mined genomic variation in rockfishes, which range in longevity from 11 to over 205 years. Multiple shifts in rockfish longevity have occurred independently and in a short evolutionary time frame, thus empowering convergence analyses. Our analyses reveal a common network of genes under convergent evolution, encompassing established aging regulators such as insulin signaling, yet also identify flavonoid (aryl-hydrocarbon) metabolism as a pathway modulating longevity. The selective pressures on these pathways indicate the ancestral state of rockfishes was long lived and that the changes in short-lived lineages are adaptive. These pathways were also used to explore genome-wide association studies of human longevity, identifying the aryl-hydrocarbon metabolism pathway to be significantly associated with human survival to the 99th percentile. This evolutionary intersection defines and cross-validates a previously unappreciated genetic architecture that associates with the evolution of longevity across vertebrates.

Exploration of the Mechanism of Linoleic Acid Metabolism Dysregulation in Metabolic Syndrome

We aimed to explore the mechanism of the linoleic acid metabolism in metabolic syndrome (MetS). RNA-seq data for 16 samples with or without MetS from the GSE145412 dataset were collected. Gene set variation analysis (GSVA), gene set enrichment analysis (GSEA), and gene differential expression analysis were performed. Expression data of differentially expressed genes (DEGs) involved in the linoleic acid metabolism pathway were mapped to the pathway by using Pathview for visualization. There were 19 and 10 differentially expressed biological processes in the disease group and healthy group, respectively. 9 KEGG pathways were differentially expressed in the disease group. Linoleic acid metabolism was the only differentially expressed pathway in the healthy group. The GSVA enrichment score of the linoleic acid metabolism pathway in the disease group was markedly lower than that in the healthy group. The GSEA result showed that the linoleic acid metabolism pathway was significantly downregulated in the disease group. JMJD7-PLA2G4B, PLA2G1B, PLA2G2D, CYP2C8, and CYP2J2 involved in the pathway were significantly downregulated in the disease group. This study may provide novel insight into MetS from the point of linoleic acid metabolism dysregulation.

Revealing key regulators of neutrophil function during inflammation by re-analysing single-cell RNA-seq

Excessive neutrophil infiltration and dysfunction contribute to the progression and severity of hyper-inflammatory syndrome, such as in severe COVID19. In the current study, we re-analysed published scRNA-seq datasets of mouse and human neutrophils to classify and compare the transcriptional regulatory networks underlying neutrophil differentiation and inflammatory responses. Distinct sets of TF modules regulate neutrophil maturation, function, and inflammatory responses under the steady state and inflammatory conditions. In COVID19 patients, neutrophil activation was associated with the selective activation of inflammation-specific TF modules. SARS-CoV-2 RNA-positive neutrophils showed a higher expression of type I interferon response TF IRF7. Furthermore, IRF7 expression was abundant in neutrophils from severe patients in progression stage. Neutrophil-mediated inflammatory responses positively correlate with the expressional level of IRF7. Based on these results, we suggest that differential activation of activation-related TFs, such as IRF7 mediate neutrophil inflammatory responses during inflammation.

Identification of Schistosoma mansoni miracidia attractant candidates in infected Biomphalaria glabrata using behaviour-guided comparative proteomics

Schistosomiasis, caused by infection with Schistosoma digenetic trematodes, is one of the deadliest neglected tropical diseases in the world. The Schistosoma lifecycle involves the miracidial infection of an intermediate freshwater snail host, such as Biomphalaria glabrata. Dispersing snail host-derived Schistosoma miracidia attractants has been considered a method of minimising intermediate host infections and, by extension, human schistosomiasis. The attractiveness of B. glabrata to miracidia is known to be reduced following infection; however, the relationship between duration of infection and attractiveness is unclear. Excretory-secretory proteins (ESPs) most abundant in attractive snail conditioned water (SCW) are key candidates to function as miracidia attractants. This study analysed SCW from B. glabrata that were naïve (uninfected) and at different time-points post-miracidia exposure (PME; 16h, 1-week, 2-weeks and 3-weeks PME) to identify candidate ESPs mediating Schistosoma mansoni miracidia behaviour change, including aggregation and chemoklinokinesis behaviour (random motion, including slowdown and increased turning rate and magnitude). Miracidia behaviour change was only observed post-addition of naïve and 3W-PME SCW, with other treatments inducing significantly weaker behaviour changes. Therefore, ESPs were considered attractant candidates if they were shared between naïve and 3W-PME SCW (or exclusive to the former), contained a predicted N-terminal signal peptide and displayed low identity (<50%) to known proteins outside of the Biomphalaria genus. Using these criteria, a total of 6 ESP attractant candidates were identified, including acetylcholine binding protein-like proteins and uncharacterised proteins. Tissue-specific RNA-seq analysis of the genes encoding these 6 ESPs indicated relatively high gene expression within various B. glabrata tissues, including the foot, mantle and kidney. Acetylcholine binding protein-like proteins were highly promising due to their high abundance in naïve and 3W-PME SCW, high specificity to B. glabrata and high expression in the ovotestis, from which attractants have been previously identified. In summary, this study used proteomics, guided by behavioural assays, to identify miracidia attractant candidates that should be further investigated as potential biocontrols to disrupt miracidia infection and minimise schistosomiasis.

SEMgsa: topology-based pathway enrichment analysis with structural equation models

BACKGROUND

Pathway enrichment analysis is extensively used in high-throughput experimental studies to gain insight into the functional roles of pre-defined subsets of genes, proteins and metabolites. Methods that leverages information on the topology of the underlying pathways outperform simpler methods that only consider pathway membership, leading to improved performance. Among all the proposed software tools, there's the need to combine high statistical power together with a user-friendly framework, making it difficult to choose the best method for a particular experimental environment.

RESULTS

We propose SEMgsa, a topology-based algorithm developed into the framework of structural equation models. SEMgsa combine the SEM p values regarding node-specific group effect estimates in terms of activation or inhibition, after statistically controlling biological relations among genes within pathways. We used SEMgsa to identify biologically relevant results in a Coronavirus disease (COVID-19) RNA-seq dataset (GEO accession: GSE172114) together with a frontotemporal dementia (FTD) DNA methylation dataset (GEO accession: GSE53740) and compared its performance with some existing methods. SEMgsa is highly sensitive to the pathways designed for the specific disease, showing low p values ([Formula: see text]) and ranking in high positions, outperforming existing software tools. Three pathway dysregulation mechanisms were used to generate simulated expression data and evaluate the performance of methods in terms of type I error followed by their statistical power. Simulation results confirm best overall performance of SEMgsa.

CONCLUSIONS

SEMgsa is a novel yet powerful method for identifying enrichment with regard to gene expression data. It takes into account topological information and exploits pathway perturbation statistics to reveal biological information. SEMgsa is implemented in the R package SEMgraph, easily available at https://CRAN.R-project.org/package=SEMgraph .

Chimeric GPCRs mimic distinct signaling pathways and modulate microglia responses

G protein-coupled receptors (GPCRs) regulate processes ranging from immune responses to neuronal signaling. However, ligands for many GPCRs remain unknown, suffer from off-target effects or have poor bioavailability. Additionally, dissecting cell type-specific responses is challenging when the same GPCR is expressed on different cells within a tissue. Here, we overcome these limitations by engineering DREADD-based GPCR chimeras that bind clozapine-N-oxide and mimic a GPCR-of-interest. We show that chimeric DREADD-β2AR triggers responses comparable to β2AR on second messenger and kinase activity, post-translational modifications, and protein-protein interactions. Moreover, we successfully recapitulate β2AR-mediated filopodia formation in microglia, an immune cell capable of driving central nervous system inflammation. When dissecting microglial inflammation, we included two additional DREADD-based chimeras mimicking microglia-enriched GPR65 and GPR109A. DREADD-β2AR and DREADD-GPR65 modulate the inflammatory response with high similarity to endogenous β2AR, while DREADD-GPR109A shows no impact. Our DREADD-based approach allows investigation of cell type-dependent pathways without known endogenous ligands.

Understanding the function of GPCRs requires stimulation with their specific ligands. Here, the authors design chemogenetic G-protein coupled receptors that allows for the study of receptors without knowing the immediate ligand, and demonstrate its use for the β2-adrenergic receptor in microglia.

Prokineticins as a Prognostic Biomarker for Low-Grade Gliomas: A Study Based on The Cancer Genome Atlas Data

Lower-grade glioma (LGG) is a crucial pathological type of glioma. Prokineticins have not been reported in LGG. Prokineticins as a member of the multifunctional chemokine-like peptide family are divided into two ligands: PROK1 and PROK2. We evaluated the role of PROK1 and PROK2 in LGG using TCGA database. We downloaded the datasets of LGG from TCGA and evaluated the influence of prokineticins on LGG survival by survival module. Correlations between clinical information and prokineticins expression were analyzed using logistic regression. Univariable survival and multivariate Cox analysis was used to compare several clinical characteristics with survival. Correlation between prokineticins and cancer immune infiltrates was explored using CIBERSORT and correlation module of GEPIA. We analyzed genes of PROK1 and PROK2 affecting LGG, screened differentially expressed genes (DEGs), interacted protein-protein with DEGs through the STRING website, then imported the results into the Cytospace software, and calculated the hub genes. To analyze whether hub genes and prokineticins are related, the relationship between PROK1 and PROK2 and hub genes was assessed and shown by heat map. In addition, gene set enrichment analysis (GSEA) was performed using the TCGA dataset. The univariate analysis using logistic regression and PROK1 and PROK2 showed opposite expression differences between tumor and normal tissues ( $p<0.05$ ). PRO1 and PROK2 expressions showed significant differences in tumor grade, age, Iiscitrate DeHydrogenase (IDH) status, histological type, and 1P/19q codeletion. Multivariate analysis revealed that the up-regulated PROK1 and PROK2 expression is an independent prognostic factor for bad prognosis. Specifically, prokineticin expression level has significant correlations with infiltrating levels of Th1 cells, NK CD 56bright cells, and Mast cells in LGG. We screened 21 DEGs and obtained 5 hub genes (HOXC10, HOXD13, SOX4, GATA4, HOXA9). GSEA-identified FCMR activation, creation of C4 and C2 activators, and CD22-mediated BCR regulation in gene ontology (GO) were differentially enriched in high PROK1 and PROK2 expression phenotype pathway, cytoplasmic ribosomal proteins, and ribosome and were differentially enriched in the low PROK1 and PROK2 expression phenotype pathway. Prokineticins are a prognostic biomarker and the correlation between hub genes and LGG requires further attention.

WGCNA Analysis Identifies the Hub Genes Related to Heat Stress in Seedling of Rice (Oryza sativa L.)

Frequent high temperature weather affects the growth and development of rice, resulting in the decline of seed–setting rate, deterioration of rice quality and reduction of yield. Although some high temperature tolerance genes have been cloned, there is still little success in solving the effects of high temperature stress in rice (Oryza sativa L.). Based on the transcriptional data of seven time points, the weighted correlation network analysis (WGCNA) method was used to construct a co–expression network of differentially expressed genes (DEGs) between the rice genotypes IR64 (tolerant to heat stress) and Koshihikari (susceptible to heat stress). There were four modules in both genotypes that were highly correlated with the time points after heat stress in the seedling. We further identified candidate hub genes through clustering and analysis of protein interaction network with known–core genes. The results showed that the ribosome and protein processing in the endoplasmic reticulum were the common pathways in response to heat stress between the two genotypes. The changes of starch and sucrose metabolism and the biosynthesis of secondary metabolites pathways are possible reasons for the sensitivity to heat stress for Koshihikari. Our findings provide an important reference for the understanding of high temperature response mechanisms and the cultivation of high temperature resistant materials.

An Interaction-Based Method for Refining Results From Gene Set Enrichment Analysis

Purpose: To demonstrate an interaction-based method for the refinement of Gene Set Enrichment Analysis (GSEA) results.

Method: Intravitreal injection of miR-124-3p antagomir was used to knockdown the expression of miR-124-3p in mouse retina at postnatal day 3 (P3). Whole retinal RNA was extracted for mRNA transcriptome sequencing at P9. After preprocessing the dataset, GSEA was performed, and the leading-edge subsets were obtained. The Apriori algorithm was used to identify the frequent genes or gene sets from the union of the leading-edge subsets. A new statistic d was introduced to evaluate the frequent genes or gene sets. Reverse transcription quantitative PCR (RT-qPCR) was performed to validate the expression trend of candidate genes after the knockdown of miR-124-3p.

Results: A total of 115,140 assembled transcript sequences were obtained from the clean data. With GSEA, the NOD-like receptor signaling pathway, C-type-like lectin receptor signaling pathway, phagosome, necroptosis, JAK-STAT signaling pathway, Toll-like receptor signaling pathway, leukocyte transendothelial migration, chemokine signaling pathway, NF-kappa B signaling pathway and RIG-I-like signaling pathway were identified as the top 10 enriched pathways, and their leading-edge subsets were obtained. After being refined by the Apriori algorithm and sorted by the value of the modulus of d, Prkcd, Irf9, Stat3, Cxcl12, Stat1, Stat2, Isg15, Eif2ak2, Il6st, Pdgfra, Socs4 and Csf2ra had the significant number of interactions and the greatest value of d to downstream genes among all frequent transactions. Results of RT-qPCR validation for the expression of candidate genes after the knockdown of miR-124-3p showed a similar trend to the RNA-Seq results.

Conclusion: This study indicated that using the Apriori algorithm and defining the statistic d was a novel way to refine the GSEA results. We hope to convey the intricacies from the computational results to the low-throughput experiments, and to plan experimental investigations specifically.

On the influence of several factors on pathway enrichment analysis

Pathway enrichment analysis has become a widely used knowledge-based approach for the interpretation of biomedical data. Its popularity has led to an explosion of both enrichment methods and pathway databases. While the elegance of pathway enrichment lies in its simplicity, multiple factors can impact the results of such an analysis, which may not be accounted for. Researchers may fail to give influential aspects their due, resorting instead to popular methods and gene set collections, or default settings. Despite ongoing efforts to establish set guidelines, meaningful results are still hampered by a lack of consensus or gold standards around how enrichment analysis should be conducted. Nonetheless, such concerns have prompted a series of benchmark studies specifically focused on evaluating the influence of various factors on pathway enrichment results. In this review, we organize and summarize the findings of these benchmarks to provide a comprehensive overview on the influence of these factors. Our work covers a broad spectrum of factors, spanning from methodological assumptions to those related to prior biological knowledge, such as pathway definitions and database choice. In doing so, we aim to shed light on how these aspects can lead to insignificant, uninteresting or even contradictory results. Finally, we conclude the review by proposing future benchmarks as well as solutions to overcome some of the challenges, which originate from the outlined factors.

CBEA: Competitive balances for taxonomic enrichment analysis

Research in human-associated microbiomes often involves the analysis of taxonomic count tables generated via high-throughput sequencing. It is difficult to apply statistical tools as the data is high-dimensional, sparse, and compositional. An approachable way to alleviate high-dimensionality and sparsity is to aggregate variables into pre-defined sets. Set-based analysis is ubiquitous in the genomics literature and has demonstrable impact on improving interpretability and power of downstream analysis. Unfortunately, there is a lack of sophisticated set-based analysis methods specific to microbiome taxonomic data, where current practice often employs abundance summation as a technique for aggregation. This approach prevents comparison across sets of different sizes, does not preserve inter-sample distances, and amplifies protocol bias. Here, we attempt to fill this gap with a new single-sample taxon enrichment method that uses a novel log-ratio formulation based on the competitive null hypothesis commonly used in the enrichment analysis literature. Our approach, titled competitive balances for taxonomic enrichment analysis (CBEA), generates sample-specific enrichment scores as the scaled log-ratio of the subcomposition defined by taxa within a set and the subcomposition defined by its complement. We provide sample-level significance testing by estimating an empirical null distribution of our test statistic with valid p-values. Herein, we demonstrate, using both real data applications and simulations, that CBEA controls for type I error, even under high sparsity and high inter-taxa correlation scenarios. Additionally, CBEA provides informative scores that can be inputs to downstream analyses such as prediction tasks.

Gene-Interaction-Sensitive enrichment analysis in congenital heart disease

Background

Gene set enrichment analysis (GSEA) uses gene-level univariate associations to identify gene set-phenotype associations for hypothesis generation and interpretation. We propose that GSEA can be adapted to incorporate SNP and gene-level interactions. To this end, gene scores are derived by Relief-based feature importance algorithms that efficiently detect both univariate and interaction effects (MultiSURF) or exclusively interaction effects (MultiSURF*). We compare these interaction-sensitive GSEA approaches to traditional χ² rankings in simulated genome-wide array data, and in a target and replication cohort of congenital heart disease patients with conotruncal defects (CTDs).

Results

In the simulation study and for both CTD datasets, both Relief-based approaches to GSEA captured more relevant and significant gene ontology terms compared to the univariate GSEA. Key terms and themes of interest include cell adhesion, migration, and signaling. A leading edge analysis highlighted semaphorins and their receptors, the Slit-Robo pathway, and other genes with roles in the secondary heart field and outflow tract development.

Conclusions

Our results indicate that interaction-sensitive approaches to enrichment analysis can improve upon traditional univariate GSEA. This approach replicated univariate findings and identified additional and more robust support for the role of the secondary heart field and cardiac neural crest cell migration in the development of CTDs.

Supplementary Information

The online version contains supplementary material available at (10.1186/s13040-022-00287-w).

Immunity-associated volatile emissions of β-ionone and nonanal propagate defence responses in neighbouring barley plants

Plants activate biochemical responses to combat stress. (Hemi-)biotrophic pathogens are fended off by systemic acquired resistance (SAR), a primed state allowing plants to respond faster and more strongly upon subsequent infection. Here, we show that SAR-like defences in barley (Hordeum vulgare) are propagated between neighbouring plants, which respond with enhanced resistance to the volatile cues from infected senders. The emissions of the sender plants contained 15 volatile organic compounds (VOCs) associated with infection. Two of these, β-ionone and nonanal, elicited resistance upon plant exposure. Whole-genome transcriptomics analysis confirmed that interplant propagation of defence in barley is established as a form of priming. Although gene expression changes were more pronounced after challenge infection of the receiver plants with Blumeria graminis f. sp. hordei, differential gene expression in response to the volatile cues of the sender plants included an induction of HISTONE DEACETYLASE 2 (HvHDA2) and priming of TETRATRICOPEPTIDE REPEAT-LIKE superfamily protein (HvTPL). Because HvHDA2 and HvTPL transcript accumulation was also enhanced by exposure of barley to β-ionone and nonanal, our data identify both genes as possible defence/priming markers in barley. Our results suggest that VOCs and plant-plant interactions are relevant for possible crop protection strategies priming defence responses in barley.

The misuse of distributional assumptions in functional class scoring gene-set and pathway analysis

Gene-set analysis (GSA) is a standard procedure for exploring potential biological functions of a group of genes. The development of its methodology has been an active research topic in recent decades. Many GSA methods, when newly proposed, rely on simulation studies to evaluate their performance with an implicit assumption that the multivariate expression values are normally distributed. This assumption is commonly adopted in GSAs, particularly those in the group of functional class scoring (FCS) methods. The validity of the normality assumption, however, has been disputed in several studies, yet no systematic analysis has been carried out to assess the effect of this distributional assumption. Our goal in this study is not to propose a new GSA method but to first examine if the multi-dimensional gene expression data in gene sets follow a multivariate normal (MVN) distribution. Six statistical methods in three categories of MVN tests were considered and applied to a total of 24 RNA data sets. These RNA values were collected from cancer patients as well as normal subjects, and the values were derived from microarray experiments, RNA sequencing, and single-cell RNA sequencing. Our first finding suggests that the MVN assumption is not always satisfied. This assumption does not hold true in many applications tested here. In the second part of this research, we evaluated the influence of non-normality on the statistical power of current FCS methods, both parametric and nonparametric ones. Specifically, the scenario of mixture distributions representing more than one population for the RNA values was considered. This second investigation demonstrates that the non-normality distribution of the RNA values causes a loss in the statistical power of these GSA tests, especially when subtypes exist. Among the FCS GSA tools examined here and among the scenarios studied in this research, the N-statistics outperform the others. Based on the results from these two investigations, we conclude that the assumption of MVN should be used with caution when evaluating new GSA tools, since this assumption cannot be guaranteed and violation may lead to spurious results, loss of power, and incorrect comparison between methods. If a newly proposed GSA tool is to be evaluated, we recommend the incorporation of a wide range of multivariate non-normal distributions or sampling from large databases if available.

Systems toxicology study reveals reduced impact of heated tobacco product aerosol extract relative to cigarette smoke on premature aging and exacerbation effects in aged aortic cells in vitro

Aging and smoking are major risk factors for cardiovascular diseases (CVD). Our in vitro study compared, in the context of aging, the effects of the aerosol of Tobacco Heating System 2.2 (THS; an electrically heated tobacco product) and 3R4F reference cigarette smoke (CS) on processes that contribute to vascular pathomechanisms leading to CVD. Young and old human aortic smooth muscle cells (HAoSMC) were exposed to various concentrations of aqueous extracts (AE) from 3R4F CS [0.014-0.22 puffs/mL] or THS aerosol [0.11-1.76 puffs/mL] for 24 h. Key markers were measured by high-content imaging, transcriptomics profiling and multianalyte profiling. In our study, in vitro aging increased senescence, DNA damage, and inflammation and decreased proliferation in the HAoSMCs. At higher concentrations of 3R4F AE, young HAoSMCs behaved similarly to aged cells, while old HAoSMCs showed additional DNA damage and apoptosis effects. At 3R4F AE concentrations with the maximum effect, the THS AE showed no significant effect in young or old HAoSMCs. It required an approximately ten-fold higher concentration of THS AE to induce effects similar to those observed with 3R4F. These effects were independent of nicotine, which did not show a significant effect on HAoSMCs at any tested concentration. Our results show that 3R4F AE accelerates aging in young HAoSMCs and exacerbates the aging effect in old HAoSMCs in vitro, consistent with CS-related contributions to the risk of CVD. Relative to 3R4F AE, the THS AE showed a significantly reduced impact on HAoSMCs, suggesting its lower risk for vascular SMC-associated pathomechanisms leading to CVD.

Silver: Forging almost Gold Standard Datasets

Gene set analysis has been widely used to gain insight from high-throughput expression studies. Although various tools and methods have been developed for gene set analysis, there is no consensus among researchers regarding best practice(s). Most often, evaluation studies have reported contradictory recommendations of which methods are superior. Therefore, an unbiased quantitative framework for evaluations of gene set analysis methods will be valuable. Such a framework requires gene expression datasets where enrichment status of gene sets is known a priori. In the absence of such gold standard datasets, artificial datasets are commonly used for evaluations of gene set analysis methods; however, they often rely on oversimplifying assumptions that make them biased in favor of or against a given method. In this paper, we propose a quantitative framework for evaluation of gene set analysis methods by synthesizing expression datasets using real data, without relying on oversimplifying or unrealistic assumptions, while preserving complex gene-gene correlations and retaining the distribution of expression values. The utility of the quantitative approach is shown by evaluating ten widely used gene set analysis methods. An implementation of the proposed method is publicly available. We suggest using Silver to evaluate existing and new gene set analysis methods. Evaluation using Silver provides a better understanding of current methods and can aid in the development of gene set analysis methods to achieve higher specificity without sacrificing sensitivity.

Comparing enrichment analysis and machine learning for identifying gene properties that discriminate between gene classes

Biologists very often use enrichment methods based on statistical hypothesis tests to identify gene properties that are significantly over-represented in a given set of genes of interest, by comparison with a 'background' set of genes. These enrichment methods, although based on rigorous statistical foundations, are not always the best single option to identify patterns in biological data. In many cases, one can also use classification algorithms from the machine-learning field. Unlike enrichment methods, classification algorithms are designed to maximize measures of predictive performance and are capable of analysing combinations of gene properties, instead of one property at a time. In practice, however, the majority of studies use either enrichment or classification methods (rather than both), and there is a lack of literature discussing the pros and cons of both types of method. The goal of this paper is to compare and contrast enrichment and classification methods, offering two contributions. First, we discuss the (to some extent complementary) advantages and disadvantages of both types of methods for identifying gene properties that discriminate between gene classes. Second, we provide a set of high-level recommendations for using enrichment and classification methods. Overall, by highlighting the strengths and the weaknesses of both types of methods we argue that both should be used in bioinformatics analyses.

GeneTrail: A Framework for the Analysis of High-Throughput Profiles

Experimental high-throughput techniques, like next-generation sequencing or microarrays, are nowadays routinely applied to create detailed molecular profiles of cells. In general, these platforms generate high-dimensional and noisy data sets. For their analysis, powerful bioinformatics tools are required to gain novel insights into the biological processes under investigation. Here, we present an overview of the GeneTrail tool suite that offers rich functionality for the analysis and visualization of (epi-)genomic, transcriptomic, miRNomic, and proteomic profiles. Our framework enables the analysis of standard bulk, time-series, and single-cell measurements and includes various state-of-the-art methods to identify potentially deregulated biological processes and to detect driving factors within those deregulated processes. We highlight the capabilities of our web service with an analysis of a single-cell COVID-19 data set that demonstrates its potential for uncovering complex molecular mechanisms. GeneTrail can be accessed freely and without login requirements at http://genetrail.bioinf.uni-sb.de.

Refining Convergent Rate Analysis with Topology in Mammalian Longevity and Marine Transitions

The quest to map the genetic foundations of phenotypes has been empowered by the modern diversity, quality, and availability of genomic resources. Despite these expanding resources, the abundance of variation within lineages makes it challenging to associate genetic change to specific phenotypes, without an a priori means of isolating the changes from background genomic variation. Evolution provides this means through convergence-i.e., the shared variation that may result from replicate evolutionary experiments across independent trait occurrences. To leverage these opportunities, we developed TRACCER: Topologically Ranked Analysis of Convergence via Comparative Evolutionary Rates. Compared to current methods, this software empowers rate convergence analysis by factoring in topological relationships, because genetic variation between phylogenetically proximate trait changes is more likely to be facilitating the trait. Comparisons are performed not with singular branches, but with the complete paths to the most recent common ancestor for each pair of lineages. This ensures that comparisons represent a single context diverging over the same timeframe while obviating the problematic requirement of assigning ancestral states. We applied TRACCER to two case studies: mammalian transitions to marine environments, an unambiguous collection of traits which have independently evolved three times; and the evolution of mammalian longevity, a less delineated trait but with more instances to compare. By factoring in topology, TRACCER identifies highly significant, convergent genetic signals, with important incongruities and statistical resolution when compared to existing approaches. These improvements in sensitivity and specificity of convergence analysis generates refined targets for downstream validation and identification of genotype-phenotype relationships.

The RNA structurome in the asexual blood stages of malaria pathogen plasmodium falciparum

Plasmodium falciparum is a deadly human pathogen responsible for the devastating disease called malaria. In this study, we measured the differential accumulation of RNA secondary structures in coding and non-coding transcripts from the asexual developmental cycle in P. falciparum in human red blood cells. Our comprehensive analysis that combined high-throughput nuclease mapping of RNA structures by duplex RNA-seq, SHAPE-directed RNA structure validation, immunoaffinity purification and characterization of antisense RNAs collectively measured differentially base-paired RNA regions throughout the parasite's asexual RBC cycle. Our mapping data not only aligned to a diverse pool of RNAs with known structures but also enabled us to identify new structural RNA regions in the malaria genome. On average, approximately 71% of the genes with secondary structures are found to be protein coding mRNAs. The mapping pattern of these base-paired RNAs corresponded to all regions of mRNAs, including the 5' UTR, CDS and 3' UTR as well as the start and stop codons. Histone family genes which are known to form secondary structures in their mRNAs and transcripts from genes which are important for transcriptional and post-transcriptional control, such as the unique plant-like transcription factor family, ApiAP2, DNA-/RNA-binding protein, Alba3 and proteins important for RBC invasion and malaria cytoadherence also showed strong accumulation of duplex RNA reads in various asexual stages in P. falciparum. Intriguingly, our study determined stage-specific, dynamic relationships between mRNA structural contents and translation efficiency in P. falciparum asexual blood stages, suggesting an essential role of RNA structural changes in malaria gene expression programs.

Over-Generalizing About GC (Hypoxia): Pitfalls of Limiting Breadth of Experimental Systems and Analyses in Framing Informatics Conclusions

Accumulating evidence suggests that many immune responses are influenced by local nutrient concentrations in addition to the programming of intermediary metabolism within immune cells. Humoral immunity and germinal centers (GC) are settings in which these factors are under active investigation. Hypoxia is an example of how a particular nutrient is distributed in lymphoid follicles during an antibody response, and how oxygen sensors may impact the qualities of antibody output after immunization. Using exclusively a bio-informatic analysis of mRNA levels in GC and other B cells, recent work challenged the concept that there is any hypoxia or that it has any influence. To explore this proposition, we performed new analyses of published genomics data, explored potential sources of disparity, and elucidated aspects of the apparently conflicting conclusions. Specifically, replicability and variance among data sets derived from different naïve as well as GC B cells were considered. The results highlight broader issues that merit consideration, especially at a time of heightened focus on scientific reports in the realm of immunity and antibody responses. Based on these analyses, a standard is proposed under which the relationship of new data sets should be compared to prior “fingerprints” of cell types and reported transparently to referees and readers. In light of independent evidence of diversity within and among GC elicited by protein immunization, avoidance of overly broad conclusions about germinal centers in general when experimental systems are subject to substantial constraints imposed by technical features also is warranted.

Gene-set distance analysis (GSDA): a powerful tool for gene-set association analysis

Background

Identifying sets of related genes (gene sets) that are empirically associated with a treatment or phenotype often yields valuable biological insights. Several methods effectively identify gene sets in which individual genes have simple monotonic relationships with categorical, quantitative, or censored event-time variables. Some distance-based methods, such as distance correlations, may detect complex non-monotone associations of a gene-set with a quantitative variable that elude other methods. However, the distance correlations have yet to be generalized to associate gene-sets with categorical and censored event-time endpoints. Also, there is a need to determine which genes empirically drive the significance of an association of a gene set with an endpoint.

Results

We develop gene-set distance analysis (GSDA) by generalizing distance correlations to evaluate the association of a gene set with categorical and censored event-time variables. We also develop a backward elimination procedure to identify a subset of genes that empirically drive significant associations. In simulation studies, GSDA more effectively identified complex non-monotone gene-set associations than did six other published methods. In the analysis of a pediatric acute myeloid leukemia (AML) data set, GSDA was the only method to discover that event-free survival (EFS) was associated with the 56-gene AML pathway gene-set, narrow that result down to 5 genes, and confirm the association of those 5 genes with EFS in a separate validation cohort. These results indicate that GSDA effectively identifies and characterizes complex non-monotonic gene-set associations that are missed by other methods.

Conclusion

GSDA is a powerful and flexible method to detect gene-set association with categorical, quantitative, or censored event-time variables, especially to detect complex non-monotonic gene-set associations. Available at https://CRAN.R-project.org/package=GSDA.

Supplementary information

The online version contains supplementary material available at 10.1186/s12859-021-04110-x.

A 7-month inhalation toxicology study in C57BL/6 mice demonstrates reduced pulmonary inflammation and emphysematous changes following smoking cessation or switching to e-vapor products

Cigarette smoking causes serious diseases, including lung cancer, atherosclerotic coronary artery disease, peripheral vascular disease, chronic bronchitis, and emphysema. While cessation remains the most effective approach to minimize smoking-related disease, alternative non-combustible tobacco-derived nicotine-containing products may reduce disease risks among those unable or unwilling to quit. E-vapor aerosols typically contain significantly lower levels of smoke-related harmful and potentially harmful constituents; however, health risks of long-term inhalation exposures are unknown. We designed a 7-month inhalation study in C57BL/6 mice to evaluate long-term respiratory toxicity of e-vapor aerosols compared to cigarette smoke and to assess the impact of smoking cessation (Cessation group) or switching to an e-vapor product (Switching group) after 3 months of exposure to 3R4F cigarette smoke (CS). There were no significant changes in in-life observations (body weights, clinical signs) in e-vapor groups compared to the Sham Control. The 3R4F CS group showed reduced respiratory function during exposure and had lower body weight and showed transient signs of distress post-exposure. Following 7 months of exposure, e-vapor aerosols resulted in no or minimal increase in pulmonary inflammation, while exposure to 3R4F CS led to impairment of lung function and caused marked lung inflammation and emphysematous changes. Biological changes observed in the Switching group were similar to the Cessation group. 3R4F CS exposure dysregulated the lung and nasal tissue transcriptome, while these molecular effects were substantially lower in the e-vapor group. Results from this study demonstrate that in comparison with 3R4F CS, e-vapor aerosols induce substantially lower biological responses including pulmonary inflammation and emphysematous changes, and that complete switching from CS to e-vapor products significantly reduces biological changes associated with CS in C57BL/6 mice.

TOXPANEL: A Gene-Set Analysis Tool to Assess Liver and Kidney Injuries

Gene-set analysis is commonly used to identify trends in gene expression when cells, tissues, organs, or organisms are subjected to conditions that differ from those within the normal physiological range. However, tools for gene-set analysis to assess liver and kidney injury responses are less common. Furthermore, most websites for gene-set analysis lack the option for users to customize their gene-set database. Here, we present the ToxPanel website, which allows users to perform gene-set analysis to assess liver and kidney injuries using activation scores based on gene-expression fold-change values. The results are graphically presented to assess constituent injury phenotypes (histopathology), with interactive result tables that identify the main contributing genes to a given signal. In addition, ToxPanel offers the flexibility to analyze any set of custom genes based on gene fold-change values. ToxPanel is publically available online at https://toxpanel.bhsai.org. ToxPanel allows users to access our previously developed liver and kidney injury gene sets, which we have shown in previous work to yield robust results that correlate with the degree of injury. Users can also test and validate their customized gene sets using the ToxPanel website.

Integrating biological knowledge and gene expression data using pathway-guided random forests: a benchmarking study

MOTIVATION

High-throughput technologies allow comprehensive characterization of individuals on many molecular levels. However, training computational models to predict disease status based on omics data is challenging. A promising solution is the integration of external knowledge about structural and functional relationships into the modeling process. We compared four published random forest-based approaches using two simulation studies and nine experimental datasets.

RESULTS

The self-sufficient prediction error approach should be applied when large numbers of relevant pathways are expected. The competing methods hunting and learner of functional enrichment should be used when low numbers of relevant pathways are expected or the most strongly associated pathways are of interest. The hybrid approach synthetic features is not recommended because of its high false discovery rate.

AVAILABILITY AND IMPLEMENTATION

An R package providing functions for data analysis and simulation is available at GitHub (https://github.com/szymczak-lab/PathwayGuidedRF). An accompanying R data package (https://github.com/szymczak-lab/DataPathwayGuidedRF) stores the processed and quality controlled experimental datasets downloaded from Gene Expression Omnibus (GEO).

SUPPLEMENTARY INFORMATION

Supplementary data are available at Bioinformatics online.

CTD: An information-theoretic algorithm to interpret sets of metabolomic and transcriptomic perturbations in the context of graphical models

We consider the following general family of algorithmic problems that arises in transcriptomics, metabolomics and other fields: given a weighted graph G and a subset of its nodes S, find subsets of S that show significant connectedness within G. A specific solution to this problem may be defined by devising a scoring function, the Maximum Clique problem being a classic example, where S includes all nodes in G and where the score is defined by the size of the largest subset of S fully connected within G. Major practical obstacles for the plethora of algorithms addressing this type of problem include computational efficiency and, particularly for more complex scores which take edge weights into account, the computational cost of permutation testing, a statistical procedure required to obtain a bound on the p-value for a connectedness score. To address these problems, we developed CTD, "Connect the Dots", a fast algorithm based on data compression that detects highly connected subsets within S. CTD provides information-theoretic upper bounds on p-values when S contains a small fraction of nodes in G without requiring computationally costly permutation testing. We apply the CTD algorithm to interpret multi-metabolite perturbations due to inborn errors of metabolism and multi-transcript perturbations associated with breast cancer in the context of disease-specific Gaussian Markov Random Field networks learned directly from respective molecular profiling data.

Gene-set Enrichment with Mathematical Biology (GEMB)

BACKGROUND

Gene-set analyses measure the association between a disease of interest and a "set" of genes related to a biological pathway. These analyses often incorporate gene network properties to account for differential contributions of each gene. We extend this concept further-defining gene contributions based on biophysical properties-by leveraging mathematical models of biology to predict the effects of genetic perturbations on a particular downstream function.

RESULTS

We present a method that combines gene weights from model predictions and gene ranks from genome-wide association studies into a weighted gene-set test. We demonstrate in simulation how such a method can improve statistical power. To this effect, we identify a gene set, weighted by model-predicted contributions to intracellular calcium ion concentration, that is significantly related to bipolar disorder in a small dataset (P = 0.04; n = 544). We reproduce this finding using publicly available summary data from the Psychiatric Genomics Consortium (P = 1.7 × 10-4; n = 41,653). By contrast, an approach using a general calcium signaling pathway did not detect a significant association with bipolar disorder (P = 0.08). The weighted gene-set approach based on intracellular calcium ion concentration did not detect a significant relationship with schizophrenia (P = 0.09; n = 65,967) or major depression disorder (P = 0.30; n = 500,199).

CONCLUSIONS

Together, these findings show how incorporating math biology into gene-set analyses might help to identify biological functions that underlie certain polygenic disorders.

Using hidden Markov model to predict recurrence of breast cancer based on sequential patterns in gene expression profiles

A new approach is presented to predict breast cancer recurrence through gene expression profiles using hidden Markov models (HMM). In this regard, 322 genes were selected from 44 published gene lists related to breast cancer prognosis. Afterwards, using gene set enrichment analysis, 922 gene sets were found from subsets of genes with the same biological meaning. In order to extract the sequential patterns from gene expression data, we ranked the gene sets using appropriate criteria and used HMM in which the ranked gene sets considered as observation sequences and hidden states represented priority of gene sets for discriminating between expression profiles. In this experiment, seven publicly available microarray datasets, including 1271 breast tumor samples, were used to classify cancer patients into two groups according to risk of recurrence. Our experiments indicated the greater performance and more robustness of the proposed model compared with other widely used classification methods.

A proposal for score assignment to characterize biological processes from mass spectral analysis of serum

•

Biological process-associated scores generated from mass spectrometry of serum.

•

Scores demonstrated acceptable levels of reproducibility.

•

Scores associated with biological processes and clinical outcome in cancer patients.

•

Possible application to biomarker studies for treatment or monitoring of disease.

•

Multiple biological processes assessed from 3 µL of serum.

Introduction

Most diseases involve a complex interplay between multiple biological processes at the cellular, tissue, organ, and systemic levels. Clinical tests and biomarkers based on the measurement of a single or few analytes may not be able to capture the complexity of a patient’s disease. Novel approaches for comprehensively assessing biological processes from easily obtained samples could help in the monitoring, treatment, and understanding of many conditions.

Objectives

We propose a method of creating scores associated with specific biological processes from mass spectral analysis of serum samples.

Methods

A score for a process of interest is created by: (i) identifying mass spectral features associated with the process using set enrichment analysis methods, and (ii) combining these features into a score using a principal component analysis-based approach. We investigate the creation of scores using cohorts of patients with non-small cell lung cancer, melanoma, and ovarian cancer. Since the circulating proteome is amenable to the study of immune responses, which play a critical role in cancer development and progression, we focus on functions related to the host response to disease.

Results

We demonstrate the feasibility of generating scores, their reproducibility, and their associations with clinical outcomes. Once the scores are constructed, only 3 µL of serum is required for the assessment of multiple biological functions from the circulating proteome.

Conclusion

These mass spectrometry-based scores could be useful for future multivariate biomarker or test development studies for informing treatment, disease monitoring and improving understanding of the roles of various biological functions in multiple disease settings.

Gene Set Analysis: Challenges, Opportunities, and Future Research

Gene set analysis methods are widely used to provide insight into high-throughput gene expression data. There are many gene set analysis methods available. These methods rely on various assumptions and have different requirements, strengths and weaknesses. In this paper, we classify gene set analysis methods based on their components, describe the underlying requirements and assumptions for each class, and provide directions for future research in developing and evaluating gene set analysis methods.

A toxicogenomic approach to assess kidney injury induced by mercuric chloride in rats

Kidney injury caused by disease, trauma, environmental exposures, or drugs may result in decreased renal function, chronic kidney disease, or acute kidney failure. Diagnosis of kidney injury using serum creatinine levels, a common clinical test, only identifies renal dysfunction after the kidneys have undergone severe damage. Other indicators sensitive to kidney injury, such as the level of urine kidney injury molecule-1 (KIM-1), lack the ability to differentiate between injury phenotypes. To address early detection as well as detailed categorization of kidney-injury phenotypes in preclinical animal or cellular studies, we previously identified eight sets (modules) of co-expressed genes uniquely associated with kidney histopathology. Here, we used mercuric chloride (HgCl₂)-a model nephrotoxicant-to chemically induce kidney injuries as monitored by KIM-1 levels in Sprague Dawley rats at two doses (0.25 or 0.50 mg/kg) and two exposure lengths (10 or 34 h). We collected whole transcriptome RNA-seq data derived from five animals at each dose and time point to perform a toxicogenomics analysis. Consistent with documented injury phenotypes for HgCl₂ toxicity, our kidney-injury-module approach identified the onset of necrosis and dilation as early as 10 h after a dose of 0.50 mg/kg that produced only mild injury as judged by urinary KIM-1 excretion. The results of these animal studies highlight the potential of the kidney-injury-module approach to provide a sensitive and histopathology-specific readout of renal toxicity.

Enhancing reproducibility of gene expression analysis with known protein functional relationships: The concept of well-associated protein

Identification of differentially expressed genes (DEGs) is well recognized to be variable across independent replications of genome-wide transcriptional studies. These are often employed to characterize disease state early in the process of discovery and prioritize novel targets aimed at addressing unmet medical need. Increasing reproducibility of biological findings from these studies could potentially positively impact the success rate of new clinical interventions. This work demonstrates that statistically sound combination of gene expression data with prior knowledge about biology in the form of large protein interaction networks can yield quantitatively more reproducible observations from studies characterizing human disease. The novel concept of Well-Associated Proteins (WAPs) introduced herein-gene products significantly associated on protein interaction networks with the differences in transcript levels between control and disease-does not require choosing a differential expression threshold and can be computed efficiently enough to enable false discovery rate estimation via permutation. Reproducibility of WAPs is shown to be on average superior to that of DEGs under easily-quantifiable conditions suggesting that they can yield a significantly more robust description of disease. Enhanced reproducibility of WAPs versus DEGs is first demonstrated with four independent data sets focused on systemic sclerosis. This finding is then validated over thousands of pairs of data sets obtained by random partitions of large studies in several other diseases. Conditions that individual data sets must satisfy to yield robust WAP scores are examined. Reproducible identification of WAPs can potentially benefit drug target selection and precision medicine studies.

Assessing Chemical-Induced Liver Injury In Vivo From In Vitro Gene Expression Data in the Rat: The Case of Thioacetamide Toxicity

Consumers are exposed to thousands of chemicals with potentially adverse health effects. However, these chemicals will never be tested for toxicity because of the immense resources needed for animal-based (in vivo) toxicological studies. Today, there are no viable in vitro alternatives to these types of animal studies. To develop an in vitro approach, we investigated whether we could predict in vivo organ injuries in rats with the use of RNA-seq data acquired from tissues early in the development of toxicant-induced injury, by comparing gene expression data from RNA isolated from these rat tissues with those obtained from in vitro exposure of primary liver and kidney cells. We collected RNA-seq data from the liver and kidney tissues of Sprague-Dawley rats 8 or 24 h after exposing them to vehicle (control), low (25 mg/kg), or high (100 mg/kg) doses of thioacetamide, a known liver toxicant that promotes fibrosis; we used these doses and exposure times to cause only mild toxicant-induced injury. For the in vitro study, we treated two cell types from Sprague-Dawley rats, primary hepatocytes (vehicle; low, 0.025 mM; or high, 0.125 mM dose), and renal tube epithelial cells (vehicle; low, 0.125 mM; or high, 0.500 mM) dose) with the thioacetamide metabolite, thioacetamide-S-oxide, selecting in vitro doses and exposure times to recreate the early-stage toxicant-induced injury model that we achieved in vivo. RNA-seq data were collected 9 or 24 h after application of vehicle or thioacetamide-S-oxide. We found that our modular approach for the analysis of gene expression data derived from in vivo RNA-seq strongly correlated (R² > 0.6) with the in vitro results at two different dose levels of thioacetamide/thioacetamide-S-oxide after 24 h of exposure. The top-ranked liver injury modules in vitro correctly identified the ensuing development of liver fibrosis.

Heterogeneity in old fibroblasts is linked to variability in reprogramming and wound healing

Age-associated chronic inflammation (inflammageing) is a central hallmark of ageing¹, but its influence on specific cells remains largely unknown. Fibroblasts are present in most tissues and contribute to wound healing^2,3. They are also the most widely used cell type for reprogramming to induced pluripotent stem (iPS) cells, a process that has implications for regenerative medicine and rejuvenation strategies⁴. Here we show that fibroblast cultures from old mice secrete inflammatory cytokines and exhibit increased variability in the efficiency of iPS cell reprogramming between mice. Variability between individuals is emerging as a feature of old age^5-8, but the underlying mechanisms remain unknown. To identify drivers of this variability, we performed multi-omics profiling of fibroblast cultures from young and old mice that have different reprogramming efficiencies. This approach revealed that fibroblast cultures from old mice contain 'activated fibroblasts' that secrete inflammatory cytokines, and that the proportion of activated fibroblasts in a culture correlates with the reprogramming efficiency of that culture. Experiments in which conditioned medium was swapped between cultures showed that extrinsic factors secreted by activated fibroblasts underlie part of the variability between mice in reprogramming efficiency, and we have identified inflammatory cytokines, including TNF, as key contributors. Notably, old mice also exhibited variability in wound healing rate in vivo. Single-cell RNA-sequencing analysis identified distinct subpopulations of fibroblasts with different cytokine expression and signalling in the wounds of old mice with slow versus fast healing rates. Hence, a shift in fibroblast composition, and the ratio of inflammatory cytokines that they secrete, may drive the variability between mice in reprogramming in vitro and influence wound healing rate in vivo. This variability may reflect distinct stochastic ageing trajectories between individuals, and could help in developing personalized strategies to improve iPS cell generation and wound healing in elderly individuals.