Meng-Shi Yang , Xiao-Jian Xu , Bin Zhang Fei Niu Bai-Yun Liu

Abstract The heterogeneity of traumatic brain injury (TBI)-induced secondary injury has greatly hampered the development of effective treatments for TBI patients. Targeting common processes across species may be an innovative strategy to combat debilitating TBI. In the present study,a cross-species transcriptome comparison was performed for the first time to determine the fundamental processes of secondary brain injury in Sprague-Dawley rat and C57/BL6 mouse models of TBI, caused by acute controlled cortical impact. The RNA sequencing data from the mouse model of TBI were downloaded from the Gene Expression Omnibus (ID: GSE79441) at the National Center for Biotechnology Information. For the rat data, peri-injury cerebral cortex samples were collected for transcriptomic analysis 24 hours after TBI. Differentially expressed gene-based functional analysis revealed that common features between the two species were mainly involved in the regulation and activation of the innate immune response, including complement cascades as well as Toll-like and nucleotide oligomerization domainlike receptor pathways. These findings were further corroborated by gene set enrichment analysis. Moreover, transcription factor analysis revealed that the families of signal transducers and activators of transcription (STAT), basic leucine zipper (BZIP), Rel homology domain (RHD),and interferon regulatory factor (IRF) transcription factors play vital regulatory roles in the pathophysiological processes of TBI, and are also largely associated with inflammation. These findings suggest that targeting the common innate immune response might be a promising therapeutic approach for TBI. The animal experimental procedures were approved by the Beijing Neurosurgical Institute Animal Care and Use Committee (approval No. 201802001) on June 6, 2018.

Key Words: cognitive impairment; cross-species comparison; gene set enrichment analysis; inflammation; innate immune; neurodegenerative disease; secondary injury; transcription factor; transcriptome; traumatic brain injury

Introduction

Traumatic brain injury (TBI) affects more than 50 million people globally each year and imposes a large burden on individuals and society (Maas et al., 2017; Tang et al., 2020). The inherently heterogeneous nature of TBI and the multiplicity of its pathophysiological processes hamper the development of effective therapies (Lipponen et al., 2018). Until recently,improvements in outcome for TBI patients have largely been achieved by the management of elevated intracranial pressure and cerebral edema via treatments such as decompressive craniectomy and mannitol; however, treatments targeting specific perturbed molecular biological processes remain scarce (Marehbian et al., 2017). Given the preventive but not therapeutically influenced treatment options for primary injury in TBI, delayed secondary injury is the most promising therapeutic window for preserving and improving neuronal function after TBI (Zhong et al., 2017; Clervius et al., 2019). A better understanding of the molecular mechanisms underlying secondary injury is therefore important for identifying potential targets for TBI treatment.

Rodent models are powerful tools in the study of human biology and disease, and have provided invaluable knowledge,especially in neuroscience research (Ellenbroek and Youn, 2016). It is estimated that approximately 35–40%of all scientific publications concerning the mammalian nervous system use rats and mice (Keifer and Summers,2016). However, in light of the intrinsic characteristics of a specific model (such as species, age, and sex), as well as the heterogeneity of complex human diseases, no single animal model can completely mimic all of the pathological and clinical features of human disorders (Jassam et al., 2017). Fortunately,cross-species omics comparisons, particularly using the increasingly well-established transcriptomics, make up for the defects of single animal models to some extent (Breschi et al.,2017). Increasing evidence indicates that comparative gene expression analysis between species is an important strategy for identifying the crucial gene functional networks in both physiological and pathological conditions (Uosaki and Taguchi,2016). Moreover, the remarkably rapid advance of nextgeneration sequencing, coupled with its plummeting price,has accelerated the application of cross-species comparisons in biomedical studies and made them practical for use in everyday research (Lachmann et al., 2018).

The complex diversity of TBI-induced secondary injury greatly hinders the development of medical treatment for TBI patients (Chen et al., 2019). However, conserved biological processes across species often play a critical role in homeostatic regulation. Therefore, the targeting of conserved or common processes may be an innovative strategy to combat debilitating TBI. To address the question of which processes are shared between rats and mice after TBI, a comparative transcriptomic analysis was conducted to identify the common molecular and cellular signatures.

Materials and Methods

Animals

Adult male Sprague-Dawley rats (aged 6–7 weeks, weighing 280–300 g,n= 30) were obtained from Beijing Vital River Experimental Animals Technology Co., Ltd. (Beijing, China,SCXK (Jing) 2016-0006). All animals were housed individually in temperature- and humidity-controlled animal quarters with free access to food and water. Rats were maintained on a 12-hour light/dark cycle to minimize discomfort. All animal procedures were approved by the Beijing Neurosurgical Institute Animal Care and Use Committee (approval No.201802001) on June 6, 2018.

Controlled cortical impact models

The rats were randomly divided into two groups (n= 15 per group): the control group and the TBI group. As in our previous study (Zhang et al., 2020a), rats were anesthetized and fixed in a stereotaxic frame (RWD Life Science Co., Shenzhen, China). In the surgical procedure, each rat’s temperature was maintained at 37.0 ± 0.5°C using a thermal plate. A midline scalp incision was made to expose the skull bone. Next, an electric drill(RWD Life Science Co.) was used to perform a craniotomy(approximately 5 mm in diameter) on the right parietal cortex without destroying the dura. Rats were subjected to controlled cortical impact (CCI) with a 5-mm diameter flat tip using a computer CCI device (PinPoint Precision Cortical Impactor, PCI 3000; Hatteras Instruments, Cary, NC, USA) with the following parameters: velocity 3 m/s, dwell time 300 ms, and depth 2 mm. The CCI model was considered successful when a cerebral cortical contusion was noticeable (Song et al., 2019).

The RNA sequencing (RNA-seq) data from the mouse model of CCI were downloaded from the Gene Expression Omnibus(GEO) (GSE79441) at the National Center for Biotechnology Information (NCBI). The impact parameters were as follows:velocity 5 m/s, dwell time 100 ms, and depth 2 mm. After impact, the skull was replaced and sealed with bone wax.The control group rats underwent the same process without impact.

RNA extraction, cDNA library construction, and sequencing

At 24 hours after TBI, rats were anesthetized using isoflurane inhalation (RWD Life Science Co.) and perfused with cold 0.9%saline solution via the left ventricle. The brain tissue was then isolated and the cerebral perilesional cortices were dissected(Zhang et al., 2019). After residual blood was removed, cortical RNA was obtained using Trizol reagent (Invitrogen, Carlsbad,CA, USA) according to the manufacturer’s instructions. Pooled samples can reduce variability and compensate for the loss of the number of replicates (Takele Assefa et al., 2020).Additionally, pooled samples can maximize biological diversity(Zhong et al., 2016; Patel et al., 2018). Therefore, following RNA isolation, the lysates of five rats were pooled into one sample, and each group contained three replicates. For the mouse samples (GSE79441), three mice in each group were pooled into one sample, producing three biological replicates in each group. The RNA purity and integrity were evaluated using a Bioanalyzer 2100 and an RNA 1000 Nano LabChip Kit (Agilent Technologies, Santa Clara, CA, USA) with RNA integrity number ＞ 7.0. The mRNA was then fragmented into small pieces and reverse-transcribed to construct a cDNA library using the manufacturer’s protocol for the mRNA-seq Sample Preparation Kit (Illumina, San Diego, CA, USA). The average insert size for the paired-end libraries was 300 bp (±50 bp). Next, paired-end sequencing was performed using a Hiseq4000 (Illumina).

RNA-seq data processing

After sequencing, the raw sequencing reads were filtered according to the following criteria: reads containing sequencing adaptors and primer were removed, and nucleotides with quality (Q) scores lower than 20 were also removed. The clean reads were then aligned to the rat genome using the HISAT (2.0) package (http://www.ccb.jhu.edu/software/hisat/index.shtml). The mapped reads were assembled using StringTie (1.3.0)(http://ccb.jhu.edu/software/stringtie/). The RNA-seq data from mouse peri-injury cortices after TBI were retrieved from the NCBI GEO database (Zhong et al., 2016). Rat and mouse protein-coding genes were used for the subsequent analysis. Cross-species expression correlations were estimated based on Pearson’s correlation analysis using R software (version: 3.6.2). Differentially expressed genes (DEGs) were determined using the DESeq2 package (version: 1.26.0). To correct for multiple testing, the false discovery rate (FDR) was calculated using the Benjamini–Hochberg (BH) procedure. DEGs were defined as those that had |log2 (fold change)| ＞1 and FDR ＜ 0.05. The DEGs that met this cutoff were selected for downstream analysis.

Gene Ontology and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analysis

Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis was used to identify significantly over-represented biological pathways for DEGs (Chen et al., 2018). Gene Ontology (GO) provides structured vocabularies for annotating genes. Generally, GO describes gene products within three independent categories: biological process (BP), molecular function (MF), and cellular component (CC). For the GO and KEGG analyses, a list of DEGs were imported to the web-based tool Metascape (https://metascape.org/gp/index.html#/main/step1), which is a biologist-oriented resource for the analysis of systems-level databases (Zhou et al., 2019). FDRs were calculated using the Benjamini-Hochberg procedure to account for multiple testing. Enriched terms with FDR ＜ 0.05 were considered to be significantly enriched.

Protein-protein interaction network and key gene analysis

Gene-encoded proteins and the protein-protein interaction(PPI) network were obtained from the Search Tool for Recurring Instances of Neighboring Genes (STRING) database,which contains both physical and functional interactions (Li et al., 2020a). The interaction network was then imported into Cytoscape, which is a public source software for visualizing and analyzing molecular interaction networks. The application cytoHubba in Cytoscape was used to identify key genes based on the node degree method (Zhang et al., 2020b). This method is a locally based method that only considers the direct neighborhood of a node.

Gene set enrichment analysis

Gene set enrichment analysis (GSEA) is a robust threshold-free computational method that determines whether an a priori defined set of genes are differentially expressed in different phenotypes (Subramanian et al., 2005). The Molecular Signatures Database (MSigDB) is a collection of annotated gene sets for use with GSEA software (https://www.gseamsigdb.org/gsea/index.jsp, version 4.0.3). The MSigDB gene sets are divided into eight major collections: H and C1–7. The C3 collection contains the regulatory target gene sets, which include subsets of microRNA targets and transcription factor targets (TFT). The C2 collection is divided into the following two sub-collections: chemical and genetic perturbations, and canonical pathways. This software was used to determine canonical pathway enrichment (the C2 KEGG subset of canonical pathways) (Kar et al., 2017). Rat and mouse genes were assigned to corresponding human orthologs using the GSEA “Collapse dataset to gene symbols” feature with the respective symbol remapping chip (Rat(mouse)_Gene_Symbol_Remapping_MSigDB.v7.0.chip). The permutation parameter was 1000 (Reisinger et al., 2020). The normalized enrichment score (NES) reflects the enrichment degree in gene expression data (Li et al., 2020b). The FDR is the estimated probability that an enriched gene set with a given NES represents a false-positive finding. Gene sets with FDR ＜0.25 were considered to be significantly enriched.

Transcription factor regulatory analysis

The Animal Transcription Factor DataBase (AnimalTFDB)can be used to annotate and predict animal transcription factor (TFs) (Hu et al., 2019). Rat and mouse TF lists were downloaded from AnimalTFDB. TFs with reads/fragments per kilobase of transcript per million RNA-seq mapped reads(RPKM/FPKM) ＞ 0 in at least three samples were set as expressed. For the PPI network, common significantly changed TFs between rats and mice were imported into STRING. The PPI network was visualized using Cytoscape. Hub TFs with the highest node degrees were determined using cytoHubba in Cytoscape. To further investigate the TF regulatory network,TFT analysis was performed using GSEA software. The FDR is the estimated probability that an enriched gene set with a given NES represents a false-positive finding. The number of permutations was 1000, and FDR ＜ 0.25 was considered to be statistically significant.

Results

Transcriptomic profiling in the rat cortex after TBI and its correlation with a mouse model of TBI

To explore the transcriptome landscape of the peri-injured cortex in a rat model of TBI, six cDNA libraries were generated and sequenced. On average, 51.7 million raw reads per cerebral cortex sample were obtained. Approximately 48.6 million (95.2%) reads per sample were aligned to the Rattus norvegicus reference genome (Genome build: rn6), where an average of 88.8% were mapped within exons, 4.2% mapped to introns, and 7.0% mapped to intergenic regions. An average of 45.7 million reads, representing about 89.7% of the total filtered reads, were mapped as paired-reads and used in the subsequent analysis.

Following the expression quantification of protein-coding genes using FPKM, the expression patterns of protein-coding genes in different samples were compared. There was a similar distribution of gene expression for protein-coding genes in all samples in both the control and TBI groups (Figure 1A). To compare transcriptomic responses to TBI between rats and mice, a transcriptome profile of mouse peri-injury cortices was downloaded from the GEO database (GSE79441). In total,14 346 expressed genes were shared by rats and mice. The expression of common genes correlated well across the two species (Pearson’s correlation analysis,r= 0.88,P＜ 2.2 × 10–16)(Figure 1B).

Identification of common differentially expressed genes in rat and mouse models of TBI

To investigate the conserved biological processes that occur in response to TBI injury in both rats and mice, common DEGs were first identified. A total of 2016 genes were differentially expressed in rats at 24 hours after TBI, of which 1413 were upregulated and 603 were downregulated. In mice, 1263 genes were differentially expressed after TBI (Figure 2A). Likewise,the majority of DEGs (89.9% of all DEGs) were upregulated.A cross-species comparison revealed that 639 differentially expressed genes were shared by both rats and mice, with 591 upregulated genes and 48 downregulated genes, respectively(Figure 2B). A scatter plot revealed a strong correlation of common DEGs between the two species when comparing the log2 fold change (Pearson’s correlation analysis,r= 0.7,P＜2.2 × 10–16) (Figure 2C). Furthermore, common DEGs were not only able to distinguish the TBI group from the control group in all rat and mouse samples, but were also able to distinguish between samples from the two species in the same group (i.e.,either the TBI or control group) (Figure 2D).

Functional classification of commonly expressed DEGs by KEGG and GO

To characterize the biological functions of the 639 common DEGs between rats and mice after TBI, KEGG analysis was performed using Metascape, a web-based tool designed for the analysis and interpretation of omics data. KEGG enrichment analysis revealed that the significantly enriched pathways were primarily involved in inflammatory processes,including proinflammatory cytokines, innate pattern recognition receptor (PRR) signaling, and inflammationrelated TF cascades (Figure 3A). Among these pathways, the tumor necrosis factor α (TNFα; the major cytokine during acute inflammation) pathway was most enriched, suggesting its critical role in secondary injury. Upon sensing pathogenassociated molecular patterns or damage-associated molecular patterns, PRRs initiate the upregulation of genes involved in inflammatory responses. Notably, both plasma membrane Toll-like receptor (TLR) signaling pathways and cytoplasmic nucleotide oligomerization domain (NOD)-like receptor signaling pathways were significantly enriched in common DEGs related to biological events. Additionally, the nuclear factor κB (NFκB), Janus kinase (JAK)–signal transducer and activator of transcription (STAT), and phosphoinositide 3-kinases (PI3K)–protein kinase B (AKT) signaling pathways were also observed in the enriched processes.

To obtain further global insights into the physiological roles of common DEGs after TBI, GO analysis was used. In this analysis,the characteristics of genes can be investigated according to their biological process, molecular function, and cellular component. GO annotation results revealed that inflammatory response, cytokine production, and leukocyte migration were the three most enriched terms in the biological process category, while cytokine binding and adhesion were foremost in that of molecular function. Within the cellular component category, the enriched terms were mainly associated with the extracellular matrix and membrane microdomain (Figure 3B).All of these GO terms further implicate neuroinflammation in the pathophysiology of TBI, as was suggested by the KEGG analysis.

Figure 1｜The distribution of rat protein-coding genes after traumatic brain injury (TBI) and their correlation with those of a mouse model of TBI.

Figure 2｜ Identification of common differentially expressed genes (DEGs) between rats and mice.

PPI network analysis and the characterization of core genes

PPI network analysis is a powerful tool for investigating biomolecular interaction networks and identifying key proteins among a cluster of molecules. To determine which of the common 639 DEGs play crucial roles in secondary injury after TBI, a PPI network was constructed using STRING and the cytoHubba application of Cytoscape. A total of 638 nodes and 6919 edges were obtained from the PPI network. The top 10 hub genes (with the highest degrees) were TNFα; interleukin(IL)-6; protein tyrosine phosphatase, receptor type, C (PTPRC);fibronectin 1 (FN1); STAT3; IL-1b; TLR4; chemokine (C-C motif)ligand 2 (CCL2); matrix metallopeptidase 9 (MMP9); and TLR2.The high degrees of these hub genes indicated that these proteins, such as TNFα, chemokines, and Toll-like genes, are principally involved in inflammatory processes. This finding is in accordance with the results of the functional analysis using KEGG and GO (Figure 3C).

Identification of enriched biological processes by GSEA

DEG-based functional analysis depends heavily on an arbitrary cutoff to identify significant genes, which might underestimate the biological importance of genes outside of the DEG selection criteria. In contrast, GSEA uses all of the information about the genes in an experiment to investigate the underlying biological mechanisms. To validate the findings from the DEG-based functional analysis, GSEA was performed for the enrichment of KEGG pathways. After converting the genes of rats and mice to corresponding human orthologs,pathway enrichment was carried out against the C2 KEGG subset of canonical pathways gene sets. The results revealed 32 and 22 enriched pathways in rats and mice, respectively.Of these, 21 canonical pathways were shared by both rats and mice (Figure 4A). Similar to the KEGG enrichment results, the common enriched pathways by GSEA were also largely implicated in inflammation, and included the JAKSTAT signaling pathway, TLR signaling pathway, and NODlike receptor signaling pathway, which had been previously observed in the DEG-based KEGG analysis (Figure 4B).Furthermore, 11 pathways were exclusively enriched in rat models of TBI, of which four pathways were associated with the DNA repair process (Figure 4C). Only RIG-I-like receptor signaling pathway was significantly enriched in mouse models compared with rat models (Figure 4D).

Cross-species expression of TFs in the cortex 24 hours after TBI

TFs are the master regulators for decoding the genome,thereby underlying multiple aspects of physiological and pathological processes. To get a better insight into TBIassociated mechanisms, abnormalities in TFs—especially common changes in the cross-species comparisons—were examined. Using the criteria that the expression of TFs was ＞0 in at least three samples, 1204 and 1323 TFs were obtained for rats and mice, respectively. Cross-species comparisons of TFs indicated that 1012 TFs were shared by both rats and mice (Figure 5A). Of these commonly expressed TFs, 111 TFs were significantly altered in the rat, of which 87 were upregulated and 24 were downregulated. For the mouse TFs, 79 were significantly altered, with 69 upregulated TFs and 10 downregulated TFs (Figure 5B). Moreover, 42 TFs were common DEGs between rats and mice, of which 39 were upregulated and 3 were downregulated. An interaction network of the 42 common TFs was constructed using STRING(Figure 5C). Based on the node degree, the top 10 hub genes were STAT3, STAT5a, interferon regulatory factor (IRF) 8, IRF1,activating transcription factor 3 (ATF3), JUN, RELB, runt-related transcription factor 3 (RUNX3), MYC, and CCAAT/enhancerbinding protein beta (CEBPB) (Figure 5D).

TF enrichment analysis by GSEA

To further understand the role of TFs in pathway perturbations after TBI, TF enrichment analysis was performed with GSEA using the TFT subset of C3: regulatory target gene sets. There were 73 and 31 enriched TFs in rats and mice, respectively.Furthermore, there were 25 common TFs between rats and mice, which were classified into seven families: STAT, basic leucine zipper (BZIP), E26 transformation-specific (ETS),Rel homology domain (RHD), IRF, paired box (PAX), and runt-related transcription factor (RUNT) (Figure 6AandB).Additionally, PPI network analysis suggested a core role of STAT3, IRF8, and IRF3 in the seven enriched families (Figure 6C).

Discussion

Cross-species transcriptome comparison not only helps us to understand the evolution of biological systems and infer gene functions, but also reveals fundamental events underlying diseases (Zhou and Gibson, 2004; Gerstein et al., 2014).Moreover, the rapid accumulation of next-generation highthroughput sequencing data from multiple species accelerates the application of cross-species comparisons in the exploration of disease-associated processes (Uosaki and Taguchi, 2016).We therefore hypothesized that the investigation of TBI might benefit from a cross-species transcriptome comparison. To this end, a cross-species transcriptome comparison was performed for the first time to explore the fundamental processes that occur after acute TBI in rats and mice. Our findings indicate that the common features between the two species are primarily implicated in the regulation and activation of the innate immune response, including complement cascades and Toll- and NOD-like receptor pathways, in acute TBI.Furthermore, our results suggest that STAT- and IRF-related TFs may play a critical regulatory role in acute TBI-induced cellular and molecular perturbations.

The innate immune system is evolutionarily conserved,and is the first line of host defense against exogenous and endogenous danger signals (Dunkelberger and Song, 2010).Innate immune defenses comprise cellular effectors and complement cascades, as well as anatomical barriers such as skin and mucous membranes (Hato and Dagher, 2015). A growing body of evidence indicates that the innate immune system greatly influences the progress and clinical outcomes of TBI (Needham et al., 2019). Abnormal complement activation can lead to the aberrant generation of opsonins and anaphylatoxins, and the consequent assembly of the terminal membrane attack complex (MAC), which has been reported to contribute to the development of secondary injury in both experimental and human TBI (Hammad et al.,2018). It has been demonstrated that the levels of C3, factor B, and sC5b-9 are significantly elevated in the ventricular cerebrospinal fluid of patients with TBI (Kossmann et al.,1997; Stahel et al., 2001). Furthermore, the complement components C1q, C3 fragments, and C5b-9 have been observed in the penumbra of cortical contusions in the human brain (Bellander et al., 2001). Interestingly, mRNA of C3 exists in the penumbra, as determined by in situ hybridization,suggesting that that complement components may be locally synthesized. Moreover, components of the lectin pathway,such as mannose-binding lectin (MBL), ficolin-2, and mannanbinding lectin serine protease 2 (MASP-2), were found to be present inside and outside brain vessels in both human and experimental cerebral contusions (Longhi et al., 2014; De Blasio et al., 2019). By taking advantage of preclinical animal models, the roles of complements in secondary injury after TBI have been validated (Leinhase et al., 2007; Ruseva et al., 2015;Rich et al., 2016). Further investigations have demonstrated that MAC is the major contributor to complement-mediated secondary injury in acute TBI (Stahel et al., 2009; Fluiter et al., 2014; Ruseva et al., 2015), whereas upstream products of complements, primarily from alternative pathways, trigger chronic neuroinflammation and cause neurodegeneration by microglial activation (Alawieh et al., 2018).

Damage-associated molecular patterns, such as ATP, high mobility group box 1 protein (HMGB1), monosodium urate,and heat shock proteins (HSPs), are released by damaged cells(Jassam et al., 2017; Krieg et al., 2017). Upon being released into the extracellular space, these molecules can trigger neuroinflammation by interacting with pattern recognition receptors, including TLRs and NOD-like receptors (Liesz et al., 2015; O’Brien et al., 2020). After brain injury, excessive activation of TLR4 (the most extensively investigated TLR in TBI) can evoke an immune response and promote cerebral edema, cumulatively resulting in cerebral secondary damage and cognitive impairment (Laird et al., 2014; Shi et al.,2019). TLR4 is mainly expressed in microglial cells. Genetic ablation of TLR4 in mice promotes microglia to switch to the M2 (anti-inflammatory) phenotype and ameliorates neurological impairment after TBI (Yao et al., 2017). Moreover,the absence or knockdown of TLR4 in experimental models weakens TBI-induced astroglial activation and astrogliosis,and mitigates the inflammatory response and subsequent brain injury, thus demonstrating that TLR4 is crucial for TBImediated inflammatory cascades in astrocytes (Ahmad et al.,2013; Jiang et al., 2018). Notably, silencing TLR4 suppresses autophagy, thereby ameliorating the neuroinflammatory response and improving neurological deficits (Jiang et al.,2018). In addition, Vaibhav et al. (2020) found that TLR4 and the downstream kinase peptidylarginine deiminase 4 (PAD4)induced neutrophil extracellular trap (NET) formation and cerebrovascular dysfunction after TBI, which may provide a promising target for improving TBI outcomes.

Figure 3｜Functional analysis using the Kyoto Encyclopedia of Genes and Genomes (KEGG) and Gene Ontology (GO) and the interaction network of common DEGs.

Figure 4｜GSEA based on the KEGG gene set.

Figure 5｜Transcription factor (TF)expression between rats and mice after traumatic brain injury.

Figure 6｜Gene set enrichment analysis based on transcription factor targets.

The NOD-like receptor protein-3 (NLRP3) inflammasome can trigger pyroptosis, an alternative form of programmed necrotic death, which is a pro-inflammatory and lytic mode of cell death that is usually associated with the secretion of IL-1β and IL-18 (Xu et al., 2017; Yang et al., 2019) . Recently,the NLRP3 inflammasome has received increasing attention in the pathophysiological processes of TBI (O’Brien et al.,2020). TBI upregulates NLRP3, apoptosis-associated specklike adapter protein (ASC), cleaved caspase-1, and cleaved IL-1β, and causes a cerebral inflammatory response whereby the NLRP3 inflammasome exacerbates neurological deficits after TBI (Ismael et al., 2018). Moreover, Ge et al. (2018) reported that NLR and AIM2 inflammasome-induced pyroptosis worsens blood–brain barrier damage in a CCI mouse model.In contrast, inhibition of the NLRP3 inflammasome by the selective inhibitor MCC950 improves acute and long-term neurological outcomes in experimental models of TBI by mitigating neuroinflammation, reducing brain edema, and alleviating TBI-induced blood-brain barrier leakage (Ismael et al., 2018; Xu et al., 2018).

Homeostatic imbalances can evoke TF-mediated transcription of specific genes, which may cause or offset biological dysfunctions (Lee and Young, 2013). The regulatory effects of TFs can be evaluated at two levels: the expression level and the activity of TFs. To understand the common characteristics of TFs between rats and mice after TBI, we first performed differential expression-based analysis. Cross-species analysis demonstrated that inflammation-associated TFs, such as STAT3, IRF8, IRF1, ATF3, JUN, RELB, and CEBPB, play a vital role in the pathophysiological processes of TBI. These findings were then corroborated by targeted analysis using GSEA.

Despite the high degree of homology between rats and mice,interspecies differences exist in gene expression patterns,social behaviors, and cognitive processes (?akowski, 2020).Notably, species-specific genes and related pathways are involved in neurogenic processes and neurological disorders(Prodromidou and Matsas, 2019). The investigation of speciesspecific biological events will deepen our understanding of the pathophysiology of secondary brain injury. However, in the present study, the rat and mouse models of TBI were induced by different treatments, such as impact parameters,which limits our interpretation of species-specific differences.Furthermore, although human and rodent models share many indispensable biological processes, the divergence of genomic and protein levels between humans and rodent models hinders the translation of findings from animal models to humans (Breschi et al., 2017). Additionally, the cerebral cortex possesses a complex cytoarchitecture (Lodato and Arlotta,2015), which is a challenge for bulk RNA-seq to dissect.Conventional bulk RNA-seq can only obtain the average expression of genes, which masks critical signals from specific cell types (Hwang et al., 2018). The responses of each cell type to TBI and the signaling pathway perturbations in each cell type remain unknown.

RNA-seq on cortical tissue as a whole may miss changes in specific cell types; in future studies we should consider the use of single-cell RNA-seq to investigate such changes

In conclusion, a lack of effective therapies for TBI has prompted intense investigation into innovative therapeutics for TBI treatment. Cross-species comparisons of the transcriptome may provide an alternative option for understanding the underlying mechanisms of TBI. The targeting of conserved or common processes may be an innovative strategy to combat debilitating TBI.

Acknowledgments:The authors acknowledge Pro. Shu-Min Wang and Qi-Qi Li from China Agricultural University and Qi-Lin Wang from Beijing Institute of Genomics, Chinese Academy of Sciences, China for experiment technical support as well as Xiang-Ning Feng from University of Chinese Academy of Sciences for English polishing.

Author contributions:Study design: BYL, XJX and MSY; animal experiment procedures: MSY, XJX and BZ; data collection: BZ and FN; statistical analysis: MSY, XJX and FN; manuscript writing: MSY; manuscript revising:XJX. All authors approved the final version of this paper.

Conflicts of interest:The authors declare that there are no conflicts of interest associated with this manuscript.

Financial support:This work was supported by the National Natural Science Foundation of China, No. 81471238, 81771327 (both to BYL),Construction of Central Nervous System Injury Basic Science and Clinical Translational Research Platform, Budget of Beijing Municipal Health Commission 2020, No. PXM2020_026280_000002 (to BYL). The funding bodies played no role in the study design, in the collection, analysis and interpretation of data, in the writing of the paper, and in the decision to submit the paper for publication.

Institutional review board statement:The animal experiment procedures were approved by the Beijing Neurosurgical Institute Animal Care and Use Committee (approval No. 201802001) on June 6, 2018. The experimental procedure followed the United States National Institutes of Health Guide for the Care and Use of Laboratory Animals (NIH Publication No. 85-23, revised 1996).

Copyright license agreement:The Copyright License Agreement has been signed by all authors before publication.

Data sharing statement:Datasets analyzed during the current study are available from the corresponding author on reasonable request.

Plagiarism check:Checked twice by iThenticate.

Peer review:Externally peer reviewed.

Open access statement:This is an open access journal, and articles are distributed under the terms of the Creative Commons Attribution-NonCommercial-ShareAlike 4.0 License, which allows others to remix,tweak, and build upon the work non-commercially, as long as appropriate credit is given and the new creations are licensed under the identical terms.

Open peer reviewer:Kathryn Leigh Wofford, University of Pennsylvania,USA; Clarissa Cavarsan, Scott and White Hospital Temple, USA.

Additional file:Open peer review report 1 and 2.