dataesr/scientific-paragraphs-categorization huggingface dataset:
Literature-based TF-target relationships facilitate TF inference To establish a comprehensive gene set database containing TF-target relationships, we considered data from different sources (Additional file 3: Table  S1 , also see Additional file 1: Supplementary Note 1). They are (D1) a literature-based database, consisting of data from TRRUST  [40] , RegNetwork  [41] , TFactS  [42] , and TRED  [43] ; (D2) a gene regulatory network database FANTOM5  [44] , whose interactions are extracted from networks constructed using RNA expression data from 394 individual tissues; (D3) a database derived from resources of putative TF binding targets, including ChEA  [45] , TRANSFAC  [46] , JASPAR  [47] , and ENCODE  [48] ; and (D4) a database derived from motif-enrichment analysis, RcisTarget  [35]. These databases have been frequently used to study the transcriptional regulations and have already been utilized for network construction  [29, 49].
We evaluated the performance of these databases by GSEA on a benchmark gene expression dataset. GSEA is a popular statistical method that can be used to evaluate the significant overlapping between a set of genes and differentially expressed genes between two experimental conditions. Using various types of TF-target databases, our goal is to find the best version of the database, so that GSEA can detect the target gene sets of the relevant TFs to be statistically significant. The benchmark dataset, denoted as set B, consists of a compilation of 12 microarray and 32 RNA-seq gene expression data (Additional file 3: Table  S2 ). Each of these datasets contains at least three samples under the normal condition (control) and three samples under the treatment condition in which a specific TF is treated by knockdown (KD).
We applied GSEA (with slight modifications, details in the "Methods" section) on set B to evaluate whether the enrichment analysis can detect the perturbed TFs. The underlying assumption is that, with a better TF-target gene set database, GSEA will be more likely to detect the corresponding perturbed TFs. For each TF-target database and each gene expression data in set B, we calculated the q-values of all the TFs in the database by GSEA to determine whether the target genes of the perturbed TF are enriched in the differentially expressed genes. We found that more significant q-values are usually associated with relatively larger number of targets for each TF; however, too many (e. g. , greater than 2000) targets will result in non-significant q-values. The summary statistics, such as the total number of TFs and the average number of target genes per TF, are summarized in Additional file 3: Table  S1.
Furthermore, these corresponding q-values from all the gene expression data are converted to specificity and sensitivity values (see the "Methods" section), and different databases are compared based on the area under the sensitivity-specificity curves (Fig. 1d ). We found that the literature-based database has the best overall performance; thus, we used this database for further analyses. Our results are in line with a previous benchmark study  [50]  that literature-based TF-target database outperforms others in capturing transcriptional regulation.
Inferring TF activity without using TF expression NetAct can accurately infer TF activity for an individual sample directly from the expression of genes targeted by the TF (see the "Methods" section). Here, we will illustrate how NetAct infers TF activity on two cases of microarray KD experiments-one case for shRNA KD of FOXM1 and shRNA KD of MYB in lymphoma cells (GEO: GSE17172  [51] ), and another case for KD of BCL6 on both OCI-Ly7 and Pfeiffer GCB-DLBCL cell lines (GEO: GSE45838  [2] ). NetAct first successfully identified the TFs that undergo knockdown in each case, i. e. , FOXM1, MYB, and BCL6, by applying GSEA on the optimized TF-target database (q-value < 0. 15). Next, for each identified TF, NetAct calculates its activity using the mRNA expression of the direct targets of the TF. We first constructed a Spearman correlation matrix from the expression of the targeted genes.
Evaluating activity inference and network construction in a simulation benchmark To evaluate the accuracy and robustness of inferred TF activity, we performed extensive benchmark tests to compare NetAct with other existing methods. We first performed the benchmark tests on simulated data because TF activity is usually not directly measurable. The activity of a TF can be related to its protein level or the level of a particular posttranslational modification, such as phosphorylation. Therefore, it is very difficult to obtain the ground truth of TF activity from an experimental dataset. Thus, in this benchmark test, we rely on mathematical modeling to simulate both the expression and activity of each TF from a synthetic TF-target network. With this simulated data, we benchmark NetAct against other methods.
To establish the simulated benchmark dataset, we first constructed a synthetic TF-target network with a total of 30 TFs. Each TF has 20 target genes randomly selected with replacement from a pool of 1000 genes. In addition, each TF also regulates two (randomly selected) of the 30 TFs. This synthetic network has a hierarchical structure, where a target gene may be co-regulated by multiple TFs. The type of each TF-to-TF regulation is either excitatory, inhibitory, or signaling, with a chance of 25%, 25%, and 50%, respectively; the type of each TF-to-target regulation is either excitatory or inhibitory with a 50% chance for each. Here, the signaling regulation changes the activity of a TF without changing its expression, whereas the excitatory or inhibitory interactions change both the activity and expression.
4a  (also Additional file 2: Figs. S3-S6 ), NetAct significantly outperforms each of the other methods in reproducing the simulated activity profiles at each perturbation level. As expected, the performance of NetAct is decreased by increasing the perturbation levels of the regulon data; however, NetAct still performs reasonably well even when only 25% of the actual target genes are kept in the regulon data. The results indicate that NetAct can robustly and accurately infer TF activity even with a noisy TF-target database. Furthermore, we tested another scenario where the test data contains simulated data from two experimental conditions, e. g. , one representing an unperturbed condition and the other representing a perturbed condition. Here, we used the same synthetic network but compiled 40 expression and activity data from the abovementioned simulation (unperturbed condition), together with 43 expression and activity data from the simulations in which a specific TF (TF9) is knocked down (perturbed condition).
We then performed a similar test as above and found that NetAct outperformed each of the other methods (Additional file 2: Fig. S2 , Additional file 2: Fig. S7a ). The notable performance gain of NetAct mainly emanates from the removal of incoherent (or noisy) targets of a TF before the activity calculation in NetAct (see the "Methods" section). In addition, we performed a network construction benchmark of NetAct and a few other network construction algorithms using the in silico simulation dataset, as shown in Fig. 4b-d. NetAct, using the TF activity inferred from the original regulon database, outperforms not only network construction methods using gene expression, such as GENIE3  [56] , GRNBoost2  [57] , and ppcor  [58, 59] , but also GENIE3 using the TF activity inferred by AUCell (Fig.
4b ). The last approach was presented to mimic a popular method SCENIC. Moreover, we evaluated the performance of NetAct when using a perturbed regulon database. We found that NetAct remains performing well when the perturbation level is as large as 50%, when evaluated by all the ground truth interactions (Fig. 4c ) and by those not presented in the regulon database (Fig. 4d ). The latter case was designed to evaluate the capability of NetAct in predicting novel interactions. We observed similar outcomes for the case of the second scenario of the simulation data from two conditions (Additional file 2: Fig. S7b-d , see Additional file 1: Supplementary Note 6 for details of the benchmark Fig. 4  The performance of activity and network inference from a simulation benchmark. a TF activity inference.
Characterizing cellular state transitions by GRN construction and modeling In the previous sections, we demonstrated the capability of NetAct in identifying the key TFs and predicting TF activity. With these data, NetAct further constructs a TF-based GRN using the mutual information (MI) of the activity from the identified TFs (details in the "Methods" section). We then applied RACIPE to the constructed network to check whether the simulated network dynamics are consistent with the experimental observations. Below, we show the utility of NetAct with two biological examples: epithelialmechanical transition (EMT) and macrophage polarization. In the first case (EMT), we analyzed a set of time-series microarray data on A549 epithelial cells undergoing TGF-β-induced epithelial-mesenchymal transition (EMT) (GEO: GSE17708)  [60].
According to the overall structure of the transcriptomics profiles, we arranged the samples from different time points into three groups-early stage (time points 0 h, 0. 5 h, and 1 h), middle stage (time points 2 h, 4 h, and 8 h), and late stage (time points 16 h, 24 h, and 72 h). We then performed three-way GSEA with our human literature-based TF-target database to identify the enriched TFs that are active between early-middle, early-late, and middle-late time points. Forty-one TFs (q-value cutoff 0. 01) were identified including many major transcriptional master regulators, such as BRCA1, CTNNB1, MYC, TWIST1, TWIST2, and ZEB1, and factors that are directly associated with TGF-β signaling pathways, such as SMAD3  [61] , FOS, and JUN  [62]. The hierarchical clustering analysis (HCA) of the expression and activity profiles for these TFs is shown in Fig.
These predictions are all consistent with existing experimental evidence (Additional file 3: Table  S4 ). Compared to a previous model of the EMT network based on an extensive literature survey  [19] , the GRN constructed by NetAct identified some of the same regulators induced by the TGF-β pathway, such as SMAD3/4, TWIST2, ZEB1, CTNNB1, NFKB1, RELA, FOS, and EGR1. Because of the lack of microRNAs and protein-protein interactions in the database, NetAct did not identify factors like miR200 and signaling molecules like PI3K. Interestingly, the NetAct model identifies STAT1/3, which was connected to other signaling pathways, such as HGF, PDGF, IGF1, and FGR, but not TGF-β in the previous network model. In addition, the NetAct model identified regulators in other important pathways in TGF-β-induced EMT in cancer cells, e.
g. , cell cycle pathway (RB1 and E2F1) and DNA damage pathway (P53). In the second case, we studied the macrophage polarization program in mouse bone marrow-derived macrophage cells using time series RNA-seq data (GEO: GSE84517)  [64]. In this experiment, macrophage progenitor cells (denoted as UT condition) were treated with (1) IFNγ to induce a transition to the M1 state, (2) IL4 to induce a transition to the M2 state, and (3) both IFNγ and IL4 to induce a transition to a hybrid M state. Here, we reprocessed the raw counts of RNA-seq with a standard protocol (details in Additional file 1: Supplementary Note 2). From principal component analysis (PCA) on the whole transcriptomics (Fig. 6b ), we found that the gene expression undergoes distinct trajectories when macrophage cells were treated with either IFNγ (M1 state) or IL4 (M2 state).
When both IFNγ and IL4 were administered, the gene expression trajectories are in the middle of the previous two trajectories, suggesting that cells are in a hybrid state (hybrid M state). We aim to use NetAct to elucidate the crosstalk in transcriptional regulation downstream of cytokine-induced signaling pathways during macrophage polarization. Here, we applied GSEA on six comparisons-untreated versus IFNγ-treated samples (one comparison between the untreated and the treated after 2 h, another between the untreated and the treated after 4 h, same for the other comparisons), untreated versus IL4-treated samples, and untreated versus IFNγ + IL4-treated samples. Using our mouse literature-based TF-target database, we identified 79 TFs (q-value cutoff 0. 05 for UT vs IL4-2 h and 0. 01 for all others).
S10 ). Taken together, we showed that the NetAct-constructed GRN model captures the multiple cellular state transitions during macrophage polarization. In conclusion, we show that NetAct can identify the core TF-based GRN using both the literature-based TF-target database and the gene expression data. We also demonstrate how RACIPE-based mathematical modeling complements NetAct-based GRN inference in elucidating the dynamical behaviors of the inferred GRNs. Together, these two methods can be applied to infer biologically relevant regulatory interactions and the dynamical behavior of biological processes.
Discussion In this study, we have developed NetAct-a computational platform for constructing and modeling core transcription factor (TF)-based regulatory networks. NetAct takes a data-driven approach to establish gene regulatory network (GRN) models directly from transcriptomics data and takes a mathematical modeling approach to characterize cellular state transitions driven by the inferred GRN. The method specifically integrates both literature-based TF-target databases and transcriptomics data of multiple experimental conditions to accurately infer TF transcriptional activity based on the expression of their target genes. Using the inferred TF activity, NetAct further constructs a TF-based GRN, whose dynamics can then be evaluated and explored by mathematical modeling. Our approach in combining top-down and bottom-up systems biology approaches will contribute to a better understanding of the gene regulatory mechanism of cellular decision-making.
One of the key components of NetAct is a pre-compiled TF-target gene set database. Here, we have evaluated different types of TF-target databases in identifying knockeddown TFs using publicly available transcriptomics datasets. In this test, we have considered databases derived from the literature, gene co-expression, cis-motif prediction, and TF-binding motif data. Our benchmark tests suggest that the literature-based database clearly outperformed the other databases. The literature-based database usually contains a small (~ 30) number of target genes for each TF, but these data have direct experimental evidence, therefore being more reliable than those from the other sources. However, the literature-based database for sure has missing regulatory interactions, therefore maybe limiting the overall performance of NetAct.
One way to address this issue is to further update the literature-based database, once new information is available. Another potential approach is to compile a database by combining different types of databases together. However, this might be quite challenging as different databases have data of very different sizes (the number of target genes) and quality. Future investigations on this direction can help to expand our knowledge of transcriptional regulation and meanwhile improve the performance of the algorithm. NetAct also has a unique approach to infer the TF activity from the gene expression of the target genes with the consideration of activation/inhibition nature. From our in silico benchmark tests, we found that NetAct outperforms major activity inference methods, owing to the design of the filtering step and the use of a high-quality TF-target database.
NetAct is also robust against some inaccuracy in the TF-target database and noises in gene expression data, because of its capability of filtering out irrelevant targets as well as remaining key targets. One potential issue is the assignment of the sign of TF activity, as it is algorithmically assigned according to the correlation with TF expression. In the case where the TF expression is very noisy or the expression is completely unrelated to TF activity, the sign assignment might be inaccurate. To deal with this issue, we have devised a semi-manual approach that identifies the sign of TF activity according to the sign of other interacting TFs. Another potential issue is that some TFs from the same family may have very similar target genes; therefore, NetAct will have difficulty in identifying exactly which TF from the family is most relevant.
Additional data resources, such as epigenomics  [68] , TF-binding data  [35] , and Hi-C data  [69] , will be helpful to address this problem. One of the future directions is to design methods to integrate these data resources. Lastly, instead of constructing a global transcriptional regulatory network, NetAct focuses on modeling a core regulatory network with only interactions between key TFs. The underlying hypothesis is that these TFs and the associated regulatory interactions play major roles in controlling the gene expression of different cellular states and the patterns of state transitions. With the core network identified using NetAct, we can further perform simulations with mathematical modeling algorithms, such as RACIPE, to analyze the control mechanism of the core network. These simulations allow us to generate new hypotheses, which can be further tested experimentally.
Conclusions We developed NetAct, a computational platform for constructing and modeling core transcription factor regulatory networks using both transcriptomics data and literaturebased transcription factor-target gene databases. Utilizing both types of resources allows us to identify regulatory genes and links specific to the data and fully take advantage of the existing knowledgebase of transcriptional regulation. Our method of combining topdown and bottom-up systems biology approaches contributes to a better understanding of the mechanism of gene regulation driving cellular state transitions.
Selecting enriched TFs For a comparison between two experimental conditions, we obtained a ranked gene list quantified by the absolute value of the test statistics (t statistics in microarray and Wald test statistics in RNA-seq) from differential expression (DE) analysis  [70] , followed by gene set enrichment analysis (GSEA)  [39]  using our optimized transcription factor (TF)target gene set database. Here, for each TF, the corresponding gene set consists of all its target genes. GSEA identifies important TFs whose targets are enriched in DE genes between the two conditions. The significance test is achieved through 10,000 permutations of the gene list names and TFs are kept for further analysis when the q-value is below a certain threshold cutoff (0. 05 by default). A C++ implementation of this version of GSEA, specifically for gene name permutations, has been provided in NetAct for fast computation.
For multiple comparisons, a set of enriched TFs are first identified from each pairwise comparison and then a union of the multiple sets of TFs is considered. In the database benchmark test, for each database, we computed the sensitivity and specificity values for different q-value cutoffs. Here, for each cutoff value, we defined the sensitivity as the proportion of datasets where the gene sets for the KD TFs were enriched with q-values below the cutoff value. We also defined specificity as the fraction of cases where the gene sets for the other TFs (non-KD TFs in the benchmark) were not enriched with q-values above the cutoff value. We then computed the area under the ROC curve (AUC) using the DescTools R package  [71].
Inferring TF activity TF activity is inferred from the expression of target genes retrieved from the TF-target database. NetAct defines the activity of the selected TFs using two different schemesone using only the expression of target genes and the other using the expression of both the TF and its target genes. The second scheme is only used for the situation of noisy target gene expression. For each TF, the algorithm selects the better scheme according to its performance, as described below.
Sign assignment for non-DE TF If the expression patterns of the identified TFs fail to show the significant differences between cell type conditions, a semi-manual method to assign the sign of activity can be adopted. Putative interaction partners between DE and non-DE TFs in the inferred network are identified using Fisher's exact test between TF targets in the NetAct TF-target database. The most significant pairs are then cross-referenced with the STRING database (https:// string-db. org) to identify instances of protein-protein interactions (PPIs). A literature search is then performed to identify the nature of the PPI, and the sign of the non-DE TF is adjusted based on the DE TF and the type of PPI. Note that the last step needs to be done manually for each modeling application. Additional file 3: Table  S3  shows the details of TF sign flipping and supported experimental evidence for the two network modeling applications.
Network construction and mathematical modeling NetAct constructs a TF regulatory network using both the TF-TF regulatory interactions from the TF-target database and the activity values. (1) The network is constructed using mutual information between the activity values of two TFs. (2)  Interactions are filtered out if they cannot be found in the TF-target regulatory database (i. e. , D1). (3) The sign of each link is determined by the sign of the Spearman correlation between the activity of two TFs. (4)  We keep the interaction between two TFs if their mutual information is higher than a threshold cutoff. With different cutoff values for mutual information, NetAct establishes networks of different sizes. To identify the best network model capturing gene expression profiles, we apply mathematical modeling to each of the TF networks using RACIPE  [29].
Availability of data and materials The information of the TF-target gene set databases is listed in Additional file 3: Table  S1. The public gene expression datasets for algorithm optimization and benchmark are listed in Additional file 3: Table  S2. The datasets and computational scripts for in silico benchmark; the network modeling scripts, including those for data processing, network construction, and network simulations; and the inferred network topology files are available in GitHub  [76]  and in Zenodo  [77]. The NetAct software is available in GitHub  [78]  and in Zenodo  [79]. NetAct is platform-independent, written in R with a partial of codes in C++ for improved performance. NetAct is licensed under the MIT License.
7%, mean relative abundance = 1. 1%) corresponding to P. viridiflava for which 74 isolated strains from the natural populations of A. thaliana were confirmed to be pathogenic  (Bartoli et al. , 2018). We therefore focused on the remaining 12 OTUs (Figure  1A ) for which taxonomic affiliation based on gyrB amplicon was identified according to GenBank BLASTN (Supplementary Data Set 1).
4 µl of 10 mM dNTP, 0. 4 µl of each 10 µM primer, 0. 1 µl of Tween 20, 0. 2 µl of 5 u/µl GoTaq G2 Polymerase (Promega), and 14. 5 µl of PCR-grade water] and mixed by pipetting up and down 3 times. Afterward, the same tip was dipped and mixed by pipetting up and down 3 times in a well of a 96-well 0. 8 ml deep well plate filled with 300 µl liquid 3 g/L TSB medium. After incubation at 28 • C for 48 h, 150 µl of 60% glycerol was added and plates were stored at -80 • C. An in-house standardized amplification was then performed as follows: 94 • C (2 min) followed by 40 cycles of amplification at 94 • C (30 s), 50 • C (60 s), and 68 • C (90 s), with a final extension step of 10 min at 68 • C. PCR products were checked on a 2% agarose gel. In order to multiplex the 7,259 amplicons, corresponding to 7,259 picked colonies, on a single Miseq run, seven forward and seven reverse gyrB primers containing different internal tags were combined for each PCR plate (Supplementary Data Set 2).