Identification of the Potential Molecular Mechanisms Linking RUNX1 Activity with Nonalcoholic Fatty Liver Disease, by Means of Systems Biology

Nonalcoholic fatty liver disease (NAFLD) is the most prevalent chronic hepatic disease; nevertheless, no definitive diagnostic method exists yet, apart from invasive liver biopsy, and nor is there a specific approved treatment. Runt-related transcription factor 1 (RUNX1) plays a major role in angiogenesis and inflammation; however, its link with NAFLD is unclear as controversial results have been reported. Thus, the objective of this work was to determine the proteins involved in the molecular mechanisms between RUNX1 and NAFLD, by means of systems biology. First, a mathematical model that simulates NAFLD pathophysiology was generated by analyzing Anaxomics databases and reviewing available scientific literature. Artificial neural networks established NAFLD pathophysiological processes functionally related to RUNX1: hepatic insulin resistance, lipotoxicity, and hepatic injury-liver fibrosis. Our study indicated that RUNX1 might have a high relationship with hepatic injury-liver fibrosis, and a medium relationship with lipotoxicity and insulin resistance motives. Additionally, we found five RUNX1-regulated proteins with a direct involvement in NAFLD motives, which were NFκB1, NFκB2, TNF, ADIPOQ, and IL-6. In conclusion, we suggested a relationship between RUNX1 and NAFLD since RUNX1 seems to regulate NAFLD molecular pathways, posing it as a potential therapeutic target of NAFLD, although more studies in this field are needed.


Introduction
Nonalcoholic fatty liver disease (NAFLD) is a condition characterized by excess fat in the liver, without alcohol implication in the onset of the disease. The term NAFLD comprehends a substantial number of liver conditions, ranging from simple steatosis (SS) to the more aggressive form of nonalcoholic steatohepatitis (NASH), which may lead to cirrhosis and hepatocellular carcinoma [1]. SS is defined as the presence of ≥5% hepatic steatosis without evidence of hepatocellular injury in the form of hepatocyte ballooning, inflammation [2], and fibrosis, three remarkable events in NASH pathology. The progress of the disease may vary from individuals, depending on the accumulated fat to the immunological and the oxidant stress responses [3,4].
NAFLD is the most prevalent chronic liver disease, with a global prevalence in adults between 23-25% [5,6]. Nevertheless, nowadays there is no definitive diagnostic test apart from invasive liver biopsy, and no specific approved treatment besides exercise and dietary interventions. Pharmacologic-based therapies for NAFLD are limited, but many clinical trials are in process [7]. For this reason, knowledge about NAFLD pathophysiology is continuously growing.
RUNX1 belongs to the runt-related transcription factor (RUNX) family of genes, and is also known as acute myeloid leukemia 1 [8]. RUNX1 regulates the differentiation of hematopoietic stem cells into mature blood cells [9,10]. It also plays a major role in the development of the neurons that transmit pain [11], and in angiogenesis and inflammation [12]. In addition, RUNX1 involvement in apoptotic processes has been reported on one hand to induce apoptosis and inhibit tumor progression in neuroblastoma [13] and leukaemia [14], while it contrarily seems to present an antiapoptotic effect in pancreatic and ovarian cancer [15,16].
Diseases associated with RUNX1 include platelet disorders with associated myeloid malignancy and blood platelet disease [17]. Related pathways include transport of glucose and other sugars, bile salts, organic acids, metal ions, and amine compounds, as well as transforming growth factor-beta (TGF-β) signaling pathways [18,19]. Recently, Kaur et al. reported a relationship between RUNX1 and NAFLD. Authors related its activity with the progression to NASH, since the interaction of RUNX1 and C-C motif chemokine 2 (CCL2), an important adhesion molecule, mediates the infiltration of pro-inflammatory and proangiogenic factors in NASH [20]. Thus, we previously wanted to study the role of RUNX1 mRNA and protein expression in NAFLD in a cohort of women with morbid obesity. We hypothesized that RUNX1 may play a protective role in NAFLD since its expression was enhanced in early stages of the disease and decrease along with the progression to NASH [21]. Given these controversies among our previous results and what was already known, the objective of the present work is to determine the proteins and the potential molecular mechanisms that could establish a link between the activity of RUNX1 and NAFLD pathogenesis by means of systems biology.

Bibliographic and Metadata Analysis in Databases
First, we built the molecular description of NAFLD pathophysiology through systematic searches and reviewing the most up-to-date scientific knowledge regarding this pathology (Supplementary Table S1). Accordingly, NAFLD was divided in specific pathophysiological processes-called motives-involved in SS, in NASH, or in both forms (Supplementary Table S1), and the corresponding molecular effectors (or key proteins) playing biological roles in these mechanisms were identified (Supplementary Table S2).
The interactome around RUNX1 was manually curated in order to better fit the mathematical models. Protein relationship databases including TRRUST database (Transcriptional Regulatory Relationships Unravelled by Sentence-based Text-mining) [22], BioGRID (The Biological General Repository for Interaction Datasets) [23], HPRD (Human Protein Reference Database) [24,25], INTACT (IntAct Molecular Interaction Database) [26], KEGG (Kyoto Encyclopedia of Genes and Genomes) [27], REACTOME (Reactome Pathway Database) [28], and available scientific literature were the sources used to identify and curate new direct interactors of RUNX1.

Mechanistic Model Generation
The compiled information was used to generate a mathematical model that simulate NAFLD pathophysiology by applying Therapeutic Performance Mapping System (TPMS) technology [29], which integrates all available biological, pharmacological, and medical knowledge to simulate human physiology in silico (Supplementary Table S3). Then, we used an artificial neural networks (ANNs) strategy [30,31] to analyze these models in order to establish the functional relationships between RUNX1 and NAFLD, considering the motives both together and individually. ANNs evaluate the relationship among protein sets or regions inside the Anaxomics network, providing a predictive score that quantifies the probability of the existence of a functional relationship between the evaluated regions. Each score is associated with a p-value that describes the probability of the result being a true positive. The ranking score has been divided into five categories: very high (ANN score > 92; p < 0.01), high (ANN score = 78-92; p = 0.01-0.05), medium-high (ANN score 71-78; p = 0.05-0.1), medium (ANN score 37-71; p = 0.1-0.25), low (ANN score < 37; >0.25).
Sampling methods-based mathematical models were then generated to determine the potential molecular mechanisms that could justify our hypothesis:

2.
Activation of RUNX1 promoting lipotoxicity and hepatic injury and liver fibrosis.
TPMS sampling-based methods trace the most probable mechanisms of action (MoA) or paths, both in biological and mathematical terms, which lead from a stimulus (e.g., activation of RUNX1) to a response (e.g., activation of IR) through the biological human protein network. In this way, it identifies the set of possible MoA that achieve a response when the system is stimulated with the specific stimulus. A population of possible solutions was obtained, and this variability was exploited and analysed to obtain a representation with the most represented paths among the set of possible solutions. A detailed description of the applied methodology was described elsewhere [29,32] and in Appendix A.

Functional Relationship between RUNX1 and NAFLD: ANNs Analysis
The possible functional relationship between RUNX1 and NAFLD, defined as the set of proteins included in its molecular characterization, has been evaluated by means of ANNs analysis. To deepen our insights, the analysis has also been performed individually for each pathophysiological motive included in NAFLD characterization: (1) increased body fat, (2) hepatic IR, (3) altered fatty acid metabolism, (4) lipotoxicity, and (5) hepatic injury and liver fibrosis. The first three pathophysiological processes occur in both SS and NASH, while the last two only happen in NASH pathophysiology or participate in the progression of NAFLD to NASH.
In this study, the relationship between RUNX1 and NAFLD or individual NAFLD motives has been evaluated, assuming that a possible functional relationship could indicate a participation of RUNX1 in NAFLD pathophysiology, either in promoting or reverting the process, since ANNs only indicate the existence of a possible relationship but not its direction. As shown in Table 1, the results obtained suggest a medium relationship of RUNX1 with the global NAFLD, considering all motives simultaneously. When considering the motives separately, however, RUNX1 seems to show a high relationship with hepatic injury and liver fibrosis, and a medium relationship with both lipotoxicity and hepatic IR.
The different columns show the ANNs score obtained for NAFLD globally and for each individual pathophysiological motive, some involved in SS and NASH stages, while others are only implicated in NASH. Category splitting was based on p-value breaks. RUNX1, runt-related transcription factor 1; NAFLD, nonalcoholic fatty liver disease; SS, simple steatosis; NASH, nonalcoholic steatohepatitis. A darker color indicates a higher ANN score.

Mechanisms of Action of RUNX1
Then, TPMS sampling methods-based mathematical models were generated simulating NAFLD pathophysiology to identify the key proteins and the most probable paths that link the activation of RUNX1 with the most strongly-related motives according to ANNs analysis (hepatic IR, and lipotoxicity, hepatic injury, and liver fibrosis). To provide new insights on the different disease stages (SS or NASH), we studied two independent MoA, considering whether the motive occurs in early or later stages of the disease: (1) RUNX1 promoting IR and (2) RUNX1 promoting lipotoxicity and hepatic injury and fibrosis, respectively. Figure 2 summarizes some of the most interesting pathways that could be regulated by RUNX1 in the context of promotion of IR in NAFLD, including the modulation of genes such as CCAAT/enhancer-binding protein alpha (CEBPA), histone deacetylase 1 (HDAC1), the transcription factor c-JUN, nuclear factor kappa B (NFκB), and some types of protein kinase C (PKCβ and PKCε). TNF, tumour necrosis factor; NFκB, nuclear factor kappa B; IL6, interleukin 6; LCN2, neutrophil gelatinase-associated lipocalin 2; CEBPA, CCAAT/enhancer-binding protein beta; SOCS3, suppressor of cytokine signaling 3; PKC, protein kinase C; JNK1, c-Jun N-terminal kinase 1; MTOR, mammalian target of rapamycin serine/threonine-protein kinase; IRS, insulin receptor substrate. This picture was generated using Graphviz software. Table 2 shows the IR effector proteins that are regulated by the activation of RUNX1 to promote this motive (considering that the activity values of the proteins in our models range from 1 to -1, only proteins with activation state > 0.1 are shown); the table contains all modulated proteins, not only those highlighted by the most represented paths. RUNX1 could be promoting IR through the regulation of 64.10% of the effector proteins involved in this motive. The IR effector proteins most activated by RUNX1-dependent downstream pathways are, in decreasing order: NFκB, JNK, PKCε, tumour necrosis factor (TNF), inhibitor of nuclear factor kappa B kinase subunit beta (IKBKB), and prostaglandin G/H synthase 2 (PTGS2); while the most inhibited ones are insulin receptor substrate (IRS)-1, phosphatase and tensin homolog (PTEN), IRS2 and sirtuin 1 (SIRT1). NAFLD, nonalcoholic fatty liver disease; NAFLD, nonalcoholic fatty liver disease; RUNX1, runt-related transcription factor 1; MoA, mechanism of action. Green color indicates a positive interaction between the effector protein and RUNX1, while red color indicates a negative one. A more intense color indicates a higher intensity of activation/inhibition.

Mechanism of Action of RUNX1 Promoting Lipotoxicity and Hepatic Injury-Liver Fibrosis
As shown in Figure 3, most of the molecular pathways that may justify the potential role of RUNX1 promoting lipotoxicity and fibrosis-related processes are shared with those involved in the motive IR. PKCε, protein kinase C epsilon; CEBPA, CCAAT/enhancer-binding protein alpha; SPP1, osteopontin; JNK1, c-Jun N-terminal kinase 1; BAX, BCL2 Associated X; CCL2, C-C motif chemokine 2; IL, interleukin; MMP2, matrix metalloproteinase-2; NLRP3, NACHT, LRR and PYD domains-containing protein 3; TLR, toll-like receptor; NOS2, inducible nitric oxide synthase; LCN2, neutrophil gelatinaseassociated lipocalin; PLIN1, perlipin; NOX, NADPH oxidase; IR, insulin resistance; LIPO, lipotoxicity; HILF, hepatic injury and liver fibrosis. This picture was generated using Graphviz software. Table 3 describes the lipotoxicity and fibrosis effector proteins that are regulated by the activation of RUNX1 (considering that the activity values of the proteins in our models range from 1 to −1, only proteins with activation state >0.1 are shown). RUNX1 promotes lipotoxicity and fibrosis by the regulation of 50.88% and 62.07% of the effector proteins involved in these motives, respectively. In total, 17 proteins specific to lipotoxicity, 24 to fibrosis, and 12 involved in both motives are regulated by RUNX1. The proteins most regulated by RUNX1 involved in lipotoxicity-related processes are: JNK1, CEBPB, and IKBKB, and those involved in fibrosisrelated processes are mothers against decapentaplegic homolog 3 (SMAD3), angiopoietin-2 (ANGPT2), apoptosis regulator BAX, type-1 angiotensin II receptor (AGTR1), and TGF-β. Effector proteins with a role in both pathophysiological processes most activated by RUNX1 are NFκB, NADPH oxidase (NOX)-1, NOX4, CCL2, and TNF. The proteins most inhibited by RUNX1 are SIRT1 (lipotoxicity) and PTEN (fibrosis). Note that the list of proteins in Table 3 is not limited to those shown in the Figure 3.  NAFLD, nonalcoholic fatty liver disease; RUNX1, runt-related transcription factor 1. Green color indicates a positive interaction between the effector protein and RUNX1, while red color indicates a negative one. A more intense color indicates a higher intensity of activation/inhibition.

Overlap between the Mechanistic Pathways Modulated by RUNX1 Activation in IR and Lipotoxicity & Fibrosis Stimulation
The NAFLD motives that have been studied for the generation of the two MoA in this project seem to be pathophysiologically related to each other since there is an overlap of effector proteins from the three motives, as described in Figure 4.
This high relationship prompted us to explore whether an overlap existed in the pathways regulated by the activation of RUNX1 in promoting these motives, and therefore, to be able to relate them. Thus, we have evaluated the similarities that interrelate the motives at the level of common effector proteins and/or pathways modulated by RUNX1 according to our models.
Common NAFLD effector proteins regulated by RUNX1 downstream mechanisms have been recognized by studying the overlap for the three motives together and studying pairs of motives separately. The proteins that we consider to be RUNX1-regulated with an activation value >0.1 are shown in Table 4. In this sense, cannabinoid receptor 1 (CNR1), which was found to be one of the common effector proteins with the three NAFLD analysed motives, presented an activation value lower than 0.1, and it is for this reason that we stop taking this protein into account from now on.
As shown in Table 4, overlapping of RUNX1-regulated proteins is observed in all three motives and in each pair. Despite finding six effector proteins that share the three NAFLD motives, only five presented sufficient signal intensity to be considered downstream effector proteins of RUNX1 inducing NAFLD; these are NFκB1, NFκB2, TNF, ADIPOQ, and IL-6.

Discussion
The novelty of this work lies in the fact that we aimed to perform a high-throughput screening to determine the molecular mechanisms that could establish a link between the activity of RUNX1 and NAFLD pathogenesis.
Until now, the connection between RUNX1 and NAFLD remains uncertain. On one hand, Kaur et al. showed a significant correlation between RUNX1 expression and inflammation, fibrosis, and NASH activity score in patients presenting NASH; they also reported RUNX1 function as a pro-angiogenic factor in SS and NASH [20]. On the other hand, Liu et al. presented low levels of RUNX1 in hepatocellular carcinoma. In this sense, these authors suggested that RUNX1 is a tumour suppressing factor that inhibits angiogenesis [33]. Regarding our previous study, we reported that the mRNA and protein expression of RUNX1 in liver seems to be involved in first steps of NAFLD with a proangiogenicrepairing role; meanwhile, RUNX1 appears to be downregulated in the NASH stage [21]. Since these disagreements, an exhaustive study of the relationship between RUNX1 MoA and NAFLD/NASH pathogenesis need to be performed to clarify this issue. In addition, this study could help to recognize RUNX1 as a potential therapeutic target of NAFLD. Previous reports have suggested that RUNX1 could be a potential therapeutic target of cancers, such as acute myeloid leukaemia, since this protein is an important regulator of haematopoiesis in vertebrates [34,35]. The beneficial effect of the therapeutic amelioration of RUNX1 in patients with nonsmall-cell lung cancer has also been described, since the RUNX1 overexpression is correlated with enhanced metastasis [36]. In addition, RUNX1 have been suggested as a potential therapeutic target to limit the progression of adverse cardiac remodeling and heart failure [37,38]. In this regard, to analyze the potential use of RUNX1 as a therapeutic target of NAFLD should be thoroughly studied. For example, investigating liver targeting through liposomes or bile acids in liver cancer [39] could be possible future strategies to evaluate the role of RUNX1 in the pathogenesis of NAFLD.
In this sense, when we performed an ANN analysis concerning the probability of the relationship between RUNX1 protein and NAFLD motives, our first main finding is that RUNX1 seems to show a medium intensity relationship with both motives-hepatic IR and lipotoxicity-and a high intensity relationship with hepatic injury and liver fibrosis motives, suggesting that this protein probably plays a role in these processes. In this regard, this result matches with Kaur et al., who reported a relevant association between RUNX1 expression and fibrosis and inflammation, two of the main NASH parameters [20]. However, this result contradicts our previous reported hypothesis about the potential protective role of RUNX1 in early stages of NAFLD [21]; in contrast, our current result has shown a low or medium intensity relationship with SS-related parameters. Given that our ANN approach provides the probability of functional relationship-regardless of the activity status (up or downregulation)-and the current conflicting results in the literature, further studies in humans or in vivo are required to clarify these contradictions, although the current available evidence clearly supports an involvement, either by presence or absence, of RUNX1 in NAHLD and NASH.
In the current literature, no direct role of RUNX1 on IR has been described yet. However, as a novelty, we demonstrate in the present study that RUNX1 interacts with proteins involved in this pathophysiological process. IR can be defined as a reduced response of the liver to the effects of insulin, which triggers impaired glucose homeostasis (gluconeogenesis and glucose uptake). IR may exert multiple effects on hepatic metabolism such as increased lipogenesis, increased free fatty acids (FFA) uptake, impaired FFA export, and decreased FFA oxidation. Moreover, outside the liver, IR causes increased serum FFA levels because of failure of insulin to suppress hormone sensitive lipase-mediated lipolysis in adipose tissue [3,[40][41][42]. In this situation, the PKCε, a downstream intermediate of RUNX1 signaling [20], is activated by the accumulation of diacylglycerol and participates in hepatic IR through impairing insulin signaling [43,44]. In addition, it is believed that RUNX1 could also be involved in IR through the transcription of IL-17 [45], a cytokine that leads to neutrophil and monocyte infiltration in the liver, thereby increasing IR [46]. In contrast, RUNX1 has been shown to inhibit the expression of Suppressor of cytokine signaling (SOCS)-3 [47]-an intracellular protein interfering with insulin signaling via ubiquitin-mediated degradation of IRS1 and IRS2 [48]-therefore ameliorating the IR. According to these facts, RUNX1 seems to have a dual role both promoting and/or preventing hepatic IR.
Another crucial event clearly involved in NAFLD progression is the lipotoxicity resulting from an excessive FFA influx to hepatocytes. Hepatic lipotoxicity occurs when the capacity of the hepatocytes to manage and export FFA as triglycerides is overwhelmed [49]. The molecular mechanisms responsible for lipotoxicity in NAFLD include endoplasmic reticulum and oxidative stress and impaired autophagy, processes that in turn activate apoptotic cascades, thus promoting tissue damage and inflammation [49].
Consequently, in conditions of hypoxia induced by steatosis [50] and inflammation, angiogenesis is triggered in chronic liver diseases [51]. It was demonstrated that proangiogenic factors have an early function in NAFLD progression from SS to NASH since proangiogenic treatments reduce not only inflammation but also steatosis [52]. In this regard, RUNX1, a pivotal regulator of hematopoiesis and angiogenesis [12,53], could be activated in order to repair the liver damage in early stages of NAFLD [21,54]. In contrast, RUNX1 activates target genes involved in lipotoxicity [55][56][57][58] such as CEBPB [59] and Cyclic AMP-dependent transcription factor (AT6) in a regulatory feed-back loop with the transcription factor AP-1 and JNK [60]. If exposure of hepatocytes to lipotoxicity and liver injury continues, it can induce apoptosis [61] and trigger inflammation by interacting with toll-like receptors (TLRs). Inflammation is a component of the wound healing process that leads to fibrosis, the deposition of extracellular matrix in liver parenchyma [62]. Additionally, RUNX1 may contribute to fibrosis and inflammation by modulation of proinflammatory cytokines (IL-1β, IL-6, TNF, etc.) [47,63], tissue inhibitor of metalloproteinase 1 (TIMP-1) [64], osteopontin [65] and TLRs [66], among others. Hence, RUNX1 seems to play a dual role, inducing pro-inflammatory cytokines and triggering liver damage, but at the same time having a protective effect by trying to repair the hepatic damage via angiogenesis-related processes.
The second notable finding of this work was obtained because we performed the TPMS technology to identify the key proteins and the most probable paths that link the activation of RUNX1 with the most strongly related motives according to ANNs analysis. In this regard, we wanted to evaluate IR first, since it is one of the main parameters involved in the first stages of NAFLD [67]. Accordingly, IR effector proteins most activated by RUNX1-dependent downstream pathways are NFκB, JNK, PKCε, TNF, IKBKB, and PTGS2, while the most inhibited ones are IRS1, PTEN, IRS2, and SIRT1.
RUNX1 could activate the PTGS2/cyclooxygenase-2 (COX-2) through its interaction with HDAC1 [68,69]. When PTGS2/COX-2 signaling is activated during inflammation in adipose tissue, it can act as a crucial factor for the promotion of obesity-induced IR and fatty liver [70,71]. According to the inflammatory role of PTGS2/COX-2, this enzyme can be induced by growth factors and different cytokines, such as TNF-α, that play a feed-back regulation role [71]. The cytokine TNF-α, produced by adipocytes and macrophages, is highly activated by the downstream mechanisms of RUNX1, particularly via the interaction with the proto-oncogene c-Jun [60,72] or the activation of NFκB [73]. The IκB kinase (IKBK) complex is the master regulator for activation of the NFκB signaling pathway. The kinase complex comprises the two catalytic subunits, IKK1 (IKBKA) and IKK2 (IKBKB), and the regulatory subunit NEMO (IKBKG), which mediates NFκB activation in response to a number of different stimuli such as RUNX1, by phosphorylating IκB proteins [74]. NFκB plays an important role in the regulation of a wide range of proteins/molecular pathways involved in IR. Its activation can be induced by TNF-α and JNK mechanisms [75,76] and can lead to the up-regulation of TNF-α, IL-6 and neutrophil gelatinase-associated lipocalin (LCN-2), contributing to IR-related processes in NAFLD [48,73,77,78]. In addition, the transcription factor AP-1 aggravates IR by inflammation-related processes, inducing the expression of IL-6 [79,80] and TNF-α [72]. TNF-α could be importantly contributing to the development of IR by inhibition and degradation of the IRS mediated by a serine phosphorylation through different mechanisms: (1) SOCS3 is induced by the NFκB/JNK-mediated activation of TNF-α and IL-6 [81], or via CEBPA activation [82], inducing ubiquitin-mediated degradation of IRS1 and IRS2 [83]; (2) MTOR can be activated by TNF-α or PKCβ pathways [84,85] due to hyperglycemia, leading to phosphorylation of multiple serine residues in IRS1 and IRS2 with their further degradation; (3) JNK1 promotes IRS1 and IRS2 serine phosphorylation [86,87]. The inhibitory effects of JNK1 could be also stimulated by PKCβ and PKCε [88,89]. In this regard, some studies have identified associations of PKC activity with disruption of the insulin-induced signal transduction pathway [90][91][92].
In contrast, apart from the degradation/silencing of IRS induced by RUNX1, which was explained above, our analysis has also reported the negative effect of RUNX1 in phosphatase and tensin homologue (PTEN) signaling. Decreased PTEN activity would lead to excessive fat deposition in the liver [40]. PTEN physiological functions negatively regulate the activity of phosphatidylinositol 3-kinase (PI3K)/AKT pathway, which in normal conditions induces lipogenesis in hepatocytes, consequently triggering IR [93,94]. PTEN downregulation has been reported to be carried out by mechanisms involving the sequential activation of MTOR and NFκB [95]. On the other hand, we have also reported a strong repression of SIRT1 by RUNX1 action. SIRT1 is an essential negative regulator of pro-inflammatory pathways, mainly through down-modulating NFκB transcriptional activity, decreasing de novo lipogenesis, and increasing fatty acid β-oxidation [96]. Hence, RUNX1mediated inhibition of SIRT1 interrupts the beneficial effect of this protein, thus promoting IR. In short, the action of all these effector proteins together gives rise to IR mechanisms.
Regarding the second main finding of this work, we wanted to focus the study of the most implicated motives in NASH stage [97]. In this sense, the effector proteins with a role in NASH (lipotoxicity and fibrosis related processes) most activated by RUNX1 are NFκB, TNF, CCL2, NOX1, and NOX4; the proteins most inhibited by RUNX1 are SIRT1 (lipotoxicity) and PTEN (fibrosis).
Additionally, the isoforms NOX1 and NOX4 seem to be upstream regulated by RUNX1. These proteins show a crucial role on both lipotoxicity and fibrosis-related processes, specially by regulating the activation of hepatic stellate cells and apoptosis, which are two important aspects in the fibrogenic process in NASH [108]. Oxidative biomolecular damage and dysregulated redox signaling induce high oxidative stress and thereby liver injury. Moreover, several studies have shown that the inhibition of NOX1 and NOX4 leads to decreased oxidative stress, lipid peroxidation, hepatic injury, inflammation, and fibrosis in NASH [108,109]. RUNX1 could induce NOX4 expression via PKCε [110,111] and NFκB dependent pathways [112], and induce NOX1 only through PKCε activation [111].
Conversely, RUNX1 have shown to inhibit SIRT1 and PTEN. Some studies have reported that liver-specific disruption of SIRT1 not only causes hepatic steatosis but also promotes the progression to an advanced metabolic disorder stage such as lipotoxicity [113]. Additionally, it seems that dysregulations of PTEN expression/activity in hepatocytes represents an important and recurrent molecular mechanism contributing to the development of liver disorders [114], given that further aberrant activation of hepatic stellate and Kupffer cells trigger the development of liver fibrosis and inflammation [95]. In summary, the pathway that constitute these effector proteins gives rise to processes of lipotoxicity and liver damage.
Accordingly, we have reported for the first-time specific MoA that RUNX1 could play a role in NAFLD pathogenesis motives, but this is only an in silico study and needs to be further validated in experimental research.
The last main objective of the present study is to analyze the overlapping proteins between the studied motives involved in NAFLD. In this sense, the shared proteins between IR and lipotoxicity most activated by RUNX1 are JNK1, IKBKB, MTOR, and LCN2, while the most inhibited by RUNX1 is SIRT1. On the other hand, the overlapping proteins observed in lipotoxicity and fibrosis motives that are the most positively modulated by RUNX1 mechanisms are NADPH oxidase NOX1 and NOX4, the chemokine CCL2, the cytochrome CYBB, and the TLRs 2, 4, and 9. The only effector that is shared between IR and fibrosis negatively modulated by RUNX1 is PTEN. Finally, the main contribution of this study is that we found five RUNX1-regulated proteins with a direct involvement in the three main NAFLD motives, which are NFκB 1, NFκB 2, TNF, ADIPOQ, and IL-6. These proteins are indicators of the relevance of their processes in terms of the relationship with RUNX1 mechanisms towards promoting NAFLD. NFκB1 and TNF present a high activation due to RUNX1 activity, as we explained above. In fact, NFκB-dependent pathways seem to definitely be a key element in these MoA due to its high number of up/downstream links, and for its important regulation of a lot of effector proteins of these motives, especially immune response-related proteins that trigger inflammation, fibrosis, or lipotoxicity [96].
On the other hand, in this study, NFκB 2, IL-6, and ADIPOQ present moderated values of activation, which differ from those of NFκB 1 and TNF. NFκB 2 is an important regulator of RUNX1. It was shown that transcription levels of NFκB 2 were increased in RUNX1deficient cells [115]. Furthermore, IL-6, as we already mentioned, is a pro-inflammatory cytokine that takes part in fibrosis and tissue damage induced by RUNX1 [47]. High TNF-α and IL-6 levels have been found in NAFLD patients, indicating an important role of these cytokines in the disease. In fact, IL-6 reduction was significantly correlated with both weight loss and insulin sensitivity [48]. Conversely, ADIPOQ seems to be downregulated by RUNX1 signaling. It has been shown that significantly up-regulated ADIPOQ expression in white adipose tissue leads to increased serum adiponectin concentrations. Low adiponectin levels are closely related to the severity of liver histology in NAFLD [116].
Our approach, as all modelling approaches, is subjected to limitations. First, it is limited by the current knowledge on the key studied elements, in this case RUNX1 functions and interactors and NAFLD molecular pathophysiology; thus, the models and conclusions are susceptible to being updated over time if prospective data and new information are generated. Nevertheless, TPMS models are built by considering the whole human protein network and a wide range of drug-pathology relationships (Supplementary Table S3); not only limited to the key studied elements, or even to hepatic involvement, they present accuracies against the training set above 80% in the case of ANN models, and above 90% in sampling methods-based models [32]. Thus, systems biology and artificial intelligence approaches allowed us to explore and present mechanistic hypotheses that are in agreement with current knowledge, providing a guide for further pre-clinical investigation in the advancement towards defining treatments for NAFLD. Further studies are needed for confirmation and advancement of these data.

Conclusions
NAFLD pathophysiological motives most functionally related to RUNX1, according to an ANNs-based analysis, are hepatic insulin resistance, lipotoxicity, and hepatic injury-liver fibrosis. These three pathophysiological processes are molecularly related, since they share NFκB1, NFκB2, TNF, ADIPOQ, and IL-6 as effector proteins. This connection suggested that RUNX1 could regulate molecular pathways involved in NAFLD pathogenesis, but more studies in this field are needed.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A
The Sampling Methods TPMS sampling-based methods [29] generate models like a Multilayer Perceptron of an Artificial Neural Network over the human protein network (where neurons are the proteins, and the edges of the network are used to transfer the information). This methodology can be used for describing with high capability all plausible relationships between an input (or stimulus, in this case RUNX1) and an output (or response, in this case NAFLD motives: IR, lipotoxicity and fibrosis protein effectors).
Although this type of network would generate many possible mechanistic solutions, it can be limited by constraints and restrictions that must be respected: (a) the topology of the protein network, (b) the functional, medical, and biological information stored in Anaxomics' databases, and (c) the available data about the drug (known effects on the target and target biology). Various different approaches and optimization systems can be used for such a purpose, from those based on randomized systems (such as a Montecarlo-based system [30]) to those based on information derived from the topology of the network, in order to solve each parameter of the equation, i.e., the weights associated to the links between the nodes in the human protein network.
The algorithms construct and analyze the regularities of the sampling of different plausible solutions. This information is used to construct feature vectors descriptive of the most probable protein network interaction structure and network activation signal flow derived from the space of plausible protein interaction solutions. The feature vectors are further used as input to supervise machine learning methods as ensembles of classifiers that allow us to infer new clinical and protein level knowledge. K-Fold and leave-one-out cross-validation methods are employed to assess generalization capability.
The mathematical algorithm can be envisioned as an extremely complex multi-parametric function, where each parameter corresponds to the relative weight of a link connecting nodes (genes/proteins) in a graph (protein map). Mathematical models of biological systems have more variables than restrictions (e.g., the number of entries in the training set is always smaller than the number of parameters-link weights-required by the algorithm), so various sets of parameters are equally valid. Therefore, using TPMS sampling-based methods, we could generate populations of solutions that comply with the biological restrictions of the training set.
From this base set of valid mathematical solutions, this approach can be employed to trace back the biological effects on molecules by analyzing the different populations of solutions. This methodology traces the most probable path (in biological and mathematical terms) that leads from the stimulus to the response through the biological network. In other words, it identifies the most probable MoA that achieve a physiological response when the system is stimulated with a specific stimulus. Not all solutions are used for the analysis, as solutions that comply with the general knowledge collated in the training set are preferred. That is, only MoAs that are plausible from the standpoint of currently accepted scientific knowledge are considered. Accuracy is calculated considering the number of restrictions in the training set that the model complies with, and only models with accuracy above 90% are considered.
In this study, the predicted MoAs were aimed at the elucidation of the mechanisms of RUNX1 activity that leads to the promotion of NAFLD motives: (A) IR and (B) lipotoxicity and fibrosis mechanisms. A set of 250 biologically plausible solutions have been calculated, with a mean accuracy of 94%.