Integrative systems biology framework discovers common gene regulatory signatures in multiple mechanistically distinct inflammatory skin diseases

More than 20% of the population across the world is affected by non-communicable inflammatory skin diseases including psoriasis, atopic dermatitis, hidradenitis suppurativa, rosacea, etc. Many of these chronic diseases are painful and debilitating with limited effective therapeutic interventions. However, recent advances in psoriasis treatment have improved the effectiveness and provide better management of the disease. This study aims to identify common regulatory pathways and master regulators that regulate molecular pathogenesis. We designed an integrative systems biology framework to identify the significant regulators across several inflammatory skin diseases. With conventional transcriptome analysis, we identified 55 shared genes, which are enriched in several immune-associated pathways in eight inflammatory skin diseases. Next, we exploited the gene co-expression-, and protein-protein interaction-based networks to identify shared genes and protein components in different diseases with relevant functional implications. Additionally, the network analytics unravels 55 high-value proteins as significant regulators in molecular pathogenesis. We believe that these significant regulators should be explored with critical experimental approaches to identify the putative drug targets for more effective treatments. As an example, we identified IKZF1 as a shared significant master regulator in three inflammatory skin diseases, which can serve as a putative drug target with known disease-derived molecules for developing efficacious combinatorial treatments for hidradenitis suppurativa, atopic dermatitis, and rosacea. The proposed framework is very modular, which can indicate a significant path of molecular mechanism-based drug development from complex transcriptomics data and other multi-omics data.


INTRODUCTION
Nearly 25% of the world's population is affected by non-communicable chronic in ammatory diseases 1 .Some of the in ammatory skin diseases include acne, atopic dermatitis (AD), actinic keratoses (AK), psoriasis (PS), hidradenitis suppurativa (HS), and three types of rosacea (RS) share speci c immune regulator activities with psoriasis [1][2][3] .Most of these diseases manifest autoimmune signatures on a temporal basis.Psoriasis is one of the best-described in ammatory skin diseases that affects approximately 1-3% of global adults 4 .Psoriasis patients manifest comorbidities such as diabetes, metabolic disorders, and severe cardiovascular disease 5 .Previous reports suggested that keratinocytes are major drivers of psoriasis in the skin epidermis, which was later expanded to the involvement of immune cells, speci cally T helper 17 cells (Th17) are the key drivers of disease 6,7 .Eventually, targeting Th17 cell signaling and several interleukin cytokines including IL-23, IL-17A, and IL-17F have improved long-term remissions by 85-100% in psoriasis patients 5,7 .Similarly, AD is an immune-mediated disease that is prevalent in approximately 25% of children and 10% of adults 8 .Likewise, HS is associated with in ammatory bowel disease and has a distinct anatomical etiology including the occurrence of in ammatory nodules, abscesses, and pus-draining sinus tracts at the places where skin rubs against each other like the underarm and groin 6 .The treatment regimens for AD and HS have very little effectiveness and cannot serve the majority of moderate to severe disease patients.Therefore, understanding the shared and unique genetic signatures, activated pathways, and regulators through innovative systems biology and integrative multi-omics can be to our advantage.
The omnipresence of bulk transcriptomics dataset introduces an extremely rewarding line of investigation to systems biologists for integration of transcriptomics to other omics layers including biological networks for a comprehensive and extrapolated exploration of some limited-studied diseases including immune-mediated diseases 9 .Network sciences have been applied to diverse biological domains to identify the associated biomolecules in disease, cancer, infections, and other stress conditions.In simple terms, biological networks are graphical representations containing nodes (genes/proteins) and edges (links) for a speci c condition 10,11 .These networks tend to change their interacting partners in different scenarios including normal/healthy and disease response 12 .Networkbased systems biology analytics has proven a great tool to unravel the structural properties of biological networks and highlight the signi cant contributors to disease pathogenesis and organismal development across different systems 9,[13][14][15][16][17][18] .In this study, we designed an integrative systems biology framework to identify the signi cant regulators across several in ammatory skin diseases.We exploited the gene coexpression network (GCN)-and protein-protein interaction (PPI)-based networks to identify shared genes and protein components in eight skin diseases with relevant functional implications.Further, we identi ed 55 high-priority proteins (HPPs) with increased network indices, which are also associated with immunemediated pathways.Finally, we explored the candidate drug-gene interactions to highlight some therapeutic-relevant target selection for putative treatment strategies.In summary, our integrative framework unravels the data-based regulatory signatures and activated molecular pathways across in ammatory skin diseases.

Transcriptomics analysis.
The differentially expressed genes (DEGs) analysis was performed on all eight disease datasets with the same threshold parameters 26 .DESeq2 was used for RNA-Seq experiments and edgeR was used for microarray experiments.

Co-expression network constriction.
We performed the gene co-expression network (GCN) construction and analysis for all eight disease datasets through weighted gene co-expression network analysis WGCNA 27 .We achieved eight diseasespeci c networks ranging from 320 nodes with 264 edges in Acne to 5,465 nodes with 132,849 edges in PS.This analysis demonstrated the challenges in the handling of multivariate transcriptomes.Networks were visualized in Cytoscape 28 .
Protein-protein interaction network construction and analysis.
To establish the disease-speci c protein-protein interaction network we fused the results from coexpression networks, DEGs, and the largest publicly available human protein-protein interaction network from the STRING database 29 .Finally, we had a collection of eight disease-speci c protein-protein interaction networks ranging from 90 nodes with 90 edges in Acne to 3,622 nodes with 33,173 edges in PS.Next, we performed network centrality analyses on all four ltered networks to identify the most important players in each disease 9 .Speci cally, we analyzed the networks by betweenness centrality (bottlenecks), degree (connections), eigenvalue (combined effect of degree and betweenness), information centrality (ability to pass stimuli information), and inner core proteins identi cation by weighted k-sell decomposition method 13 .Networks were visualized in Cytoscape 28 .

Network exposibility analysis.
To check the extent of coverage in a network and inspired by the COVID-19 exposibility scenario, we calculated the network exposibility by dividing the total number of nodes by the total connected components in the network.Based on the network exposibility threshold ≥ 10 to determine the reliability of eight PPI networks, we ltered four disease-speci c interaction networks including HS, PS, AD, and RS for further analysis.

High-priority protein identi cation.
To identify the High-priority proteins (HPPs) or signi cant proteins for each PPI network we calculated the top 5% of centrality value nodes with the degree, betweenness, eigenvector, information centrality, and the inner layer proteins from the weighted k-shell decomposition method.Next, if these signi cant proteins are also signi cantly activated or inhibited regulators characterized by pathway analysis in each of the eight in ammatory skin diseases, then these proteins are termed HPPs or signi cant regulators.Additionally, we identi ed the most important regulators shared among four diseases and mapped them against their regulated pathways.

Drug-gene interaction network.
We extracted the HPPs interacting drugs from the publicly available databases DGIdb 30 .Only the interactions with listed publication information were included in this analysis.These drug-gene (HPPs) interactions can provide signi cant information about the potential therapeutic options for a disease or class of diseases.The drug-gene interaction network was visualized in Cytoscape 28 .
Gene Ontology and Pathway Analysis.
The canonical pathway and regulator analysis was performed by IPA with default parameters.The gene ontology and pathway analysis were performed through metascape 31 with standard parameters.

Signi cance analysis.
The DEG analysis was performed with (log2FC ≥|1|; FDR < 0.05) parameters.The canonical pathway analysis was performed with a signi cance test BH P-value < 0.05.The gene ontology analysis was performed with a signi cance test P-value < 0.05.The network power-law distribution tness threshold is r 2 ≥ 0.5.The network exposibility threshold was 10 for reliable networks.The correlation among values of different centralities is positive if r ≥ 0.5.The signi cance of degree and betweenness among the inner and outer layers of four PPIs were tested by Mann-Whitney test for non-parametric distributed values.The signi cant regulator analysis through IPA was tested by BH P-value < 0.05.

RESULTS
A Network Science-based Framework to Prioritize Genes/Pathways from Transcriptome Datasets.
To study the shared and unique genetic signatures, regulators, and prioritizing genes in different chronic in ammatory skin diseases, we implemented an integrative multi-omics framework utilizing the transcriptome datasets and extrapolating the network-based analysis pipeline.In this framework, 11 different steps are followed to identify the most appropriate drug-gene interactions (Fig. 1). 1. First, we extracted publicly available transcriptome datasets in eight different in ammatory skin diseases (namely: acne, atopic dermatitis (AD), actinic keratoses (AK), psoriasis (PS), hidradenitis suppurativa (HS), and three types of rosacea (RS)) as described in the methods.2. Subsequently, we performed the gene coexpression network construction and analysis for all eight disease datasets through WGCNA 27 .We achieved eight disease-speci c networks ranging from 320 nodes with 264 edges in Acne to 5,465 nodes with 132,849 edges in PS.Which demonstrated the challenges in the handling of multivariate transcriptomes.3. The third step was to perform DEGs analysis on all eight disease datasets with the same threshold parameters 26 .We observed some similarities and differences in the expression pro les among eight diseases.4. To establish the disease-speci c protein-protein interaction network, we next combined the results from co-expression networks, DEGs, and the largest publicly available human protein-protein interaction network from the STRING database 29 .Finally, we had a collection of eight disease-speci c protein-protein interaction networks ranging from 90 nodes with 90 edges in Acne to 3,622 nodes with 33,173 edges in Psoriasis.5. Based on the network explosibility, a comparable number of nodes and interactions, as well as power-law distribution of eight networks, we selected four diseasespeci c interaction networks including atopic dermatitis (AD), psoriasis (PS), hidradenitis suppurativa (HS), and rosacea (RS) for further analysis.6. Afterwards, we performed network centrality analyses on all four ltered networks to identify the most important players in each disease 9 .Speci cally, we analyzed the networks by betweenness centrality (bottlenecks), degree (hubs), eigenvalue (combined effect of degree and betweenness), information centrality (ability to pass stimuli information), and inner core proteins identi cation by the weighted k-sell decomposition method.7.With that, we obtained a collection of high-priority proteins (HPPs) in these four in ammatory skin diseases.When comparing the appearance of these HPPs in all four diseases, we identi ed shared and unique HPPs in each disease.8.Then, we mapped the HPPs in signi cantly activated/inhibited signaling pathways in these four skin diseases.Interestingly, we found some highly activated signaling pathways are shared among all four diseases.Additionally, we identi ed the most important regulators shared among four diseases and mapped them against their regulated pathways.9. Moving forward, we identi ed the HHPs interacting with drugs from the publicly available databases DGIdb 30 .These drug-gene (HPPs) interactions can provide signi cant information about the potential therapeutic options for a disease or class of diseases.10.From these data, we propose that these drugs can be repurposed for treating several in ammatory skin diseases either alone or in combination.11.All the predictive outcomes including the role of genes, pathways, HPPs, drug-HPP interactions, and other predictions are required to be validated experimentally for each disease condition.Taken together, we proposed a modular framework that can be applied to any polygenic chronic in ammatory disease to re ne the regulator identi cation, pathway association, and drug-gene prediction for the alternative therapeutic targets in any immunogenetic-related disease.
Different chronic in ammatory skin diseases display unique and conserved transcriptome activity.
The conventional transcriptome-based analysis demonstrates the expression behavior in disease for the corresponding cohort of samples.To begin, we identi ed the total number of DEGs from the transcriptomics (RNA-Seq and microarray) datasets in in ammatory skin diseases including acne, AD, AK, CD, HS, ICD, PS, and three types of RS; erythematotelangiectatic rosacea (ER), phymatous rosacea (PhyR), and papulopustular rosacea (PapR) with a threshold (log2FC ≥|1|; FDR < 0.05).Based on the transcriptomics platform and diseases, we captured a diverse number of DEGs per disease (Fig. 2A, Table S1).For example, ICD has 65 up-and 47 down-regulated genes, whereas most DEGs were discovered in PS disease with 1763 up-and 3243 down-regulated genes, respectively.Whereas, acne, CD, and HS have a similar breakdown of DEGs (Table S1).Afterward, we performed the functional analysis for DEGs in each disease to explore the activated/inhibited pathways in those datasets.We report that the Th1 pathway, dendritic cell maturation, Th17 activation pathway, interferon signaling, Th2 pathway, IL-8 signaling, GP6 signaling pathway, CD28 signaling in T helper cells, HOTAIR regulatory pathway, NF-κB signaling are some of the signi cant activated canonical pathways in these chronic in ammatory skin diseases (BH P-value < 0.05; Fig. 2B, Table S1).Interestingly, we found that the senescence pathway, 1L-7 signaling, and IL-13 signaling are activated only in acne and AK, whereas B Cell receptor signaling was activated in HS, PapR, and CD.
These initial results suggest that there is a huge complexity among chronic in ammatory skin diseases with only a few shared and distinct genes and associated pathways.Thus, we explored the system biology-based methods to unravel intricate relationships, novel gene regulators, and master regulators of pathways altered in these diseases.
The correlation method predicted novel genes involved in four (AD, HS, PS, and RS) chronic skin diseases.
The correlation-based gene clustering is used to reduce the high-dimensional transcriptomics datasets for easier human interpretations by clustering genes into several groups based on their expression pro le.GCN-based clustering is also a subtype of conventional clustering with reduced complexity through network representation.The resultant GCN can be a source point to predict the relationships among genes in a biological pathway through the guilt-by-association concept 32 .However, the main challenge is to restrict the spurious correlations, which can be achieved to a great extent by using datasets with a higher number of samples/conditions, robust GCN network construction algorithms, and standard parameters for each experiment.Thus, we performed our correlation-based analysis through WGCNA with the same parameters of power-cutoff, weight-cutoff, and clustering algorithm for eight in ammatory skin disease transcriptomes 27 .There are eight resultant GCNs of which acne (320 nodes and 264 edges) is the smallest, HS is the largest by an edge (4,162 nodes and 420,422 edges) and, PS (5,465 nodes and 132,849 edges) is the largest by nodes.AK, CD, and ICD have a comparable number of nodes and an edge.(Fig. 3A, Table S2).
To begin our GCNs analysis, we performed several network centrality analyses including the degree to check the scale-free properties of each network 33 .As expected, we found that all eight GCNs; Acne, AD, AK, CD, HS, ICD, PS, and RS follow scale-free properties calculated by the power-law distribution with r 2 values 0.94 for AD and 0.64 for ICD (Fig. 3B, Table S2).Though most of these GCNs are biologically relevant through power-law distribution, they don't seem comparable to each other.To check the extent of coverage in a network, we calculated the network exposibility as described in the methods section.We found that four GCNs corresponding to acne, AK, CD, and ICD have a small number of nodes and a high number of connected components, whereas the remaining four other GCNs related to AD, HS, PS, and RS have a high number of nodes but a smaller number of connected components comparably (Fig. 3C, Table S2).To determine the exposibility-based reliability of eight GCNs, we set a network exposibility threshold ≥ 10 and ltered our four disease-speci c GCNs including HS, PS, AD, and RS for further analysis (Fig. 3D).
It is well reported that some of the immune-mediated diseases including psoriasis, HS, and AD share most of their disease-associated pathways and comorbidities with each other 34,35 .Therefore, to understand the shared gene patterns among these diseases, we explored the core co-expression network of four diseases.Interestingly, we report a core co-expression network among PS, HS, and AD with 145 genes and 161 associations (Fig. 3E, Table S2).Interestingly, some of these genes are RPS-, RPL-, MRPL-, EIF-, COX-family genes, ILF2, PPIA, NDUFAB1, NDUFA4, ATP5C1, and FH.Furthermore, to verify the biological associations of 145 core genes, we performed gene ontology analysis.The ontology analysis identi ed the enriched pathways including translation initiation complex formation, diabetic cardiomyopathy, ribosome assembly, regulation of proteolysis, TNF-alpha/NF-kappa B signaling complex 6, telomerase RNA localization to Cajal body, translation factors, amoebiasis, Interferon type I signaling, Emerin complex-52, protein targeting, and DNA repair (P-value < 0.05; Fig. 3F, Table S2).Similarly, we explored the core co-expression network among PS, HS, and RS.In this analysis, we found a small core network encompassing 30 genes and 30 associations (Fig. 3G, Table S2).Interestingly, some of these genes are IKAROS Family Zinc Finger 1 (IKZF1), CD2 Molecule (CD2), CD53 Molecule (CD53), SAM And SH3 Domain Containing 3 (SASH3), Lysosomal Protein Transmembrane 5 (LAPTM5), Interleukin 10 Receptor Subunit Alpha (IL10RA), Protein Tyrosine Phosphatase Receptor Type C (PTPRC), and Pleckstrin (PLEK).Furthermore, to verify the biological associations of 30 core genes in HS, AD, and RS, we performed gene ontology analysis.We found that most of the enriched pathways are cell activation, signaling by Interleukins, CSI at the vascular wall, downstream TCR signaling, heterotypic cell-cell adhesion, neutrophil degranulation, B cell receptor signaling, cytokine in immune response, and cellular homeostasis (P-value < 0.05; Fig. 3H, Table S2).These results suggest the signi cance of correlationbased GCNs construction and analyses to identify the emerging and shared players in associated pathways among in ammatory skin diseases.
Interactome analysis identi ed the central proteins in psoriasis, hidradenitis suppurativa, atopic dermatitis, and rosacea.
The molecular mechanisms and intricacies of biological processes associated with immune-mediated skin diseases are mostly determined by protein interactions in certain pathways 36 .Given that the coverage of DEGs from transcriptome to differentially expressed proteins is a maximum of only ~ 60%, we used the disease-speci c co-expressed gene list to extract the PPI network.These interactions can be more relevant in any stress condition if they are also part of corresponding GCNs by satisfying the guiltby-association assumption 37 .To explore these possibilities, we constricted the disease-speci c PPI network by integrating the co-expressed genes into the known interactions from the STRING database 29 .
The resultant PPI networks/interactomes for these four diseases include PS with the greatest number of nodes (3,622 nodes and 33,173 edges) and RS with the least number of nodes (851 nodes and 3,234 edges) (Fig. 4A, Table S3).It is reported that some of the molecular pathways and disease comorbidities are shared among psoriasis, hidradenitis suppurativa, atopic dermatitis, and rosacea 1 .To identify the core proteins in four PPIs, we performed the commonality analysis and found that 53 proteins are shared among these four in ammatory skin diseases (Fig. 4B, Table S3).Some of these proteins are CD2, CD3D, CD53, CD96, CDK1, COL1A1, COL1A2, COL3A1, COX5A, COX7B, CXCR4, DCN, FH, IKZF3, IL10RA, ITGAL, KRT6A, KRT6B, LAPTM5, and LCK, which participate in lymphocyte activation, diabetic cardiomyopathy, neutrophil degranulation, membrane raft distribution, lymphocyte-mediated immunity, TYROBP causal network in microglia, response to molecule of bacterial origin, microglia pathogen phagocytosis pathway, homeostasis of several cells, formation of the corni ed envelope, negative regulation of defense response, cell cycle phase transition, positive regulation of endopeptidase activity, regulation of glial cell differentiation, acute viral myocarditis, human cytomegalovirus infection, and immune effector process (Figure S1, Table S3).
The disease-speci c interactomes position their important proteins in a particular arrangement to communicate throughout the network most e ciently.To explore these structural arrangements, we performed network analysis on four disease-speci c PPIs including PS, HS, AD, and RS.To begin our analysis rst, we tested the biological relevance of interactomes by performing the power-law distribution analysis of four PPIs to verify their scale-freeness 33,38 .Interestingly, we report that all four PPIs for PS, HS, AD, and RS follow scale-free properties calculated by the power-law distribution and the r 2 values are similar (Fig. 4C, Table S3).Previously it has been reported that the network-based method has successfully highlighted the central and core proteins associated with disease pathogenesis 9,36,38 .Therefore, we leveraged a part of this analysis in our framework and analyzed four PPIs with a degree, betweenness centrality, eigenvector, information centrality, and weighted k-sell decomposition method.First, we explored the degree distribution of expressed genes (Exp; FDR < 0.05) and non-expressed (Not-Exp; FDR > 0.05) for each disease.Interestingly, we found that expressed genes of all four diseases namely HS, PS, AD, and RS encompass a high degree distribution than their non-expressed counterparts (Fig. 4D, Table S3), Mann-Whitney test P-value (HS < 0.05, PS < 0.0001, AD < 0.0001, and RS > 0.66).
Previously, we have demonstrated that some of the highly connected genes also possess other network properties, which make them extremely vulnerable during disease pathogenesis 9,18 .Therefore, we calculated the correlation between degree centrality and other centralities including information, eigenvector, and betweenness centralities of all four PPI networks.Interestingly, we found that the degree of PS and HS are strongly correlated with information centrality with r > 0.6 for both, whereas this correlation is not strong in AD and RS with r < 0.5 (Fig. 4E).Similarly, we computed the correlation among degree and eigenvector centrality and reported that PS, HS, and AD have the strongest correlation (r < 0.5), whereas RS has the weakest correlation (r < 0.5, Fig. 4F).Differently, the correlation between degree and betweenness is very poor in PS, AD, and RS (r < 0.5) and strong in HS (r > 0.5) (Figure S1B).Further, we decomposed the PPI networks into k-shells through the Weighted k-shell decomposition method to identify the hidden genes, which the conventional centralities fail to identify 13 .Next, we identi ed the inner layer and peripheral (outer) layer proteins in each PPI as mentioned by Ahmed et al., 2018 13 .We report that there are 461, 573, 393, 95 proteins in the inner layer and 3135, 1982, 1466, and 754 proteins in the outer layer of PS, HS, AD, and RS, respectively (Table S3).Moving forward, we explored the degree and betweenness centrality distributions of the inner layer and outer layer proteins in four PPIs.We report that inner-layer proteins possess a signi cantly higher degree compared to outer-layer proteins (Fig. 4G; Mann-Whitney test P-value < 0.0001, Table S3).However, we do not see much difference in the betweenness centrality distribution of inner and outer layer proteins (Fig. 4H; Mann-Whitney test P-value < 0.0001).These analyses con rm that some of the proteins of PPIs share different types of high centrality values, which can be central in the disease pathogenesis.Additionally, we investigated the high centrality proteins (signi cant proteins) by selecting the top5% of high centralities i.e Hub (degree), bottleneck (Betweenness, BW), information centrality (IC), and eigenvector centrality (EV) for individual diseasespeci c PPI network.We demonstrate that 46, 51, and 23 proteins are shared by all four centralities in PS, HS, and AD, respectively (Fig. 4I-K, Table S3).Whereas there was not a single protein shared by all four centralities in RS (Fig. 4L, Table S3).Interestingly, we also found that most of the high betweenness centrality (bottleneck) proteins are not shared by any other centrality.These analyses indicate that most of these bottlenecks are part of smaller subnetworks in four PPIs.
Network centrality-based prioritization of proteins and pathways in chronic in ammatory skin diseases.
It has been previously described in several instances that in the disease-speci c interactomes, the most signi cant contributing proteins have a high degree (connections) and high betweenness (bottlenecks) against other proteins throughout the PPI network [14][15][16]36,39 . Furthr, the analysis was expanded to other centralities including information and eigenvector centralities, which improved the identi cation of signi cantly contributing proteins in the interactome 9 .Inspired by these studies, we identi ed the signi cant proteins for disease-speci c PPI networks of PS.HS, AD, and RS individually.The top 5% of centrality value nodes with the degree, betweenness, eigenvector, information centrality, and the inner layer proteins from the weighted k-shell decomposition method were identi ed as signi cant proteins.Next, we put our effort into classifying these signi cant proteins as regulators in four skin diseases in this study and termed them as HPPs or signi cant regulators (Fig. 5A, Table S4).As a result, we identi ed 55 HPPs that contribute signi cantly to the pathogenesis of four selected diseases (Table S4).Further, we explored the regulator activity of these 55 HPPs and found that most of these proteins are signi cantly activated in at least seven in ammatory skin diseases (Fig. 5B, BH P-value < 0.05; Activation z-score > |1|; Table S4).Interestingly, some of these 55 HPPs are involved in RTKs Signaling, proteoglycans in cancer, cell morphogenesis, hemopoiesis, PID CXCR4 pathway, epithelial cell differentiation, Cytokine Signaling in the Immune system, adherens junction, regulation of kinase activity, cellular response to lipid, head, and neck SCC, response to wounding, Hippo signaling regulation, Interferon type I signaling, PI3K-Akt signaling, leukocyte activation, response to estradiol, lymphocyte activation, Th1, and Th2 cell differentiation, and NK cell-mediated cytotoxicity (Fig. 5C, Table S4; P-value < 0.05).Most of these HPPs and associated pathways are known for their major contribution to several immune-mediated diseases 5- 7,20,25,36,38,40,41 . Mving forward, we investigated the reported potential therapeutics for candidate 55 HPPs, which will help in designing highly effective drug repurposing strategies.We extracted the Drug-Gene interaction network pairs from DGIdb 30 for 55 HPPs with available literature publications.We found a total of 32 HPPs have a known compendium of 199 drug compounds with 237 interactions (Fig. 5D, Table S4).Interestingly, we report that some 55 HPPs i.e., CD2, LCK, STAT1, TNFRSF1B, IKZF1, APP, and BMP7 can be a high-priority drug target for several of these chronic skin diseases.
IKZF1 is a potential therapeutic target for several in ammatory skin diseases.
Given the overlap in the molecular pro les of multiple in ammatory skin diseases, we explored the shared signi cant regulator proteins in HS, PS, AD, and RS.For this, we complied all the signi cant proteins based on ve centrality methods (degree, betweenness, information, eigenvector centrality, and weighted k-shell decomposition) for each of these diseases.As a result, we have identi ed a list of 236, 195, 168, and 72 signi cant proteins respectively in PS, HS, AD, and RS (Fig. 6A, Table S4).Afterward, we studied the distribution of shared and unique signi cant proteins in each PPI network and found that six proteins (NDUFAB1, MRPL3, DDX1, EIF2S1, SLIRP, and ATP5C1) are contributing signi cantly to these four diseases.Whereas there is only one protein IKZF1 contributing signi cantly to HS, AD, and RS with high network centrality values, we explored the drug-gene interaction of IKZF1.In this regard, six drugs including lenalidomide, daunorubicin, methotrexate, cytarabine, imatinib, and udarabine are known to target IKZF1 (Fig. 6B).Interestingly, we also found that daunorubicin also targets APP and BMP7 which are also part of predicted 55 HPPs.Similarly, methotrexate also targets CCND1 and BMP7, and cytarabine targets BMP7 too (Fig. 5D).Thereafter, we extracted all the possible interactions of IKZF1 from HS-PPI, HS-GCN, and known TF-target relationships to explore the extrapolating effect of a drug on the IKZF1 and its partners.We found that IKZF1 can be a really good candidate as it acts as a master regulator by interacting with 19 targets, 62 proteins, and 88 co-expressed genes (Fig. 6C, Table S4).We hypothesized that the activities of these proteins and co-expressed genes can be altered by targeting IKZF1 with a speci c drug, such as lenalidomide.As IKZF1 is shared among three of these skin diseases, we surveyed the expression pattern of IKZF1 interacting and co-expressing genes in HS, PS, AD, and RS.Interestingly, we found that most of the IKZF1 interacting partners; 19 targets, 62 proteins, and 88 co-expressed genes follow a similar expression pattern across four in ammatory skin diseases (Fig. 6D-F, Table S4).
Taken together, our study identi ed high-value proteins involved in the pathogenesis of chronic skin diseases including HS, PS, AD, and RS.Additionally, we designed a robust framework to identify signi cant contributors utilizing transcriptome datasets and integrative multi-omics approaches supported by network systems biology for more accurate predictions of regulators in most chronic skin diseases.

DISCUSSION
Almost 1/5 th of the human population across the world is affected by some type of non-communicable chronic in ammatory disease.The disease manifestation and immune regulatory signatures are shared among some of these in ammatory skin diseases.The high-throughput 'omics technology and robust multi-omics network integration techniques showed immense potential to accelerate more comprehensive and system-wide discoveries in the pathogenesis of complex in ammatory skin diseases.Here, our designed robust framework identi ed signi cant contributors from different transcriptome datasets and integrative multi-omics supported the network systems biology for more accurate predictions of signi cant regulators of these in ammatory skin diseases (namely: acne, atopic dermatitis (AD), actinic keratoses (AK), contact dermatitis (CD), irritant contact dermatitis (ICD), psoriasis (PS), hidradenitis suppurativa (HS), and three types of rosacea (RS)).Our conventional and unconventional biological data analysis approach extracted the shared and unique biomolecular (gene/protein/regulator) signatures and biological processes in four (PS, PS, AD, and RS) of these eight in ammatory skin diseases.Further, systems genetics and network biology-driven analyses were instrumental in prioritizing the most relevant proteins (55 high-priority proteins; HPPs) across four of these in ammatory skin diseases.These HPPs can serve as a template to study the impact and response of disease pathogenesis of these in ammatory skin diseases, and putative drug targets strategies either as standalone or combinatorial approaches.Overall, this study employs advanced and powerful network systems-based analyses that combine integrative multi-omics to establish a mutual understanding of in ammatory mechanisms in four skin diseases including PS, PS, AD, and RS.
In the last two-decade, PS has been one of the most studied in ammatory skin diseases, which shares a signi cant number of immune-mediated pathways with other in ammatory diseases including acne, AD, AK, HS, and RS.For example, similar to psoriasis, targeting Th2 cells, IL-4, and IL-13 have been the most e cacious (50-70%), and advanced therapy for AD 40 and Janus kinase (JAK) inhibitors has been demonstrated to suppress the cytokine responses more effectively and are also a powerful treatment strategy for AD 40 .HS is also an in ammatory skin disease implicated by the pathogenesis of neutrophilic in ammation, dysbiosis, TNF, interferon responses, hair-and skin-gland abnormalities, autoantibodies, and plasma cells 41 .Our multi-omics analysis identi ed the shared molecules, regulators, and pathways at transcriptome and interactome layers in four diseases including TNF-alpha/NF-kappa B signaling, translation factors, amoebiasis, Interferon type I signaling, Cytokine signaling, and interleukin signaling.
These pathways have been under investigation to manage/treat several immune system diseases.
Our network analytics utilizing different conventional and improved centralities appreciably identi ed the central proteins of these skin disease-speci c interactomes as well as substantially improved the classi cation of signi cant proteins (55 HPPs) in the pathogenesis of PS, HS, AD, and RS.The functional analysis of these HPPs further highlighted the Cytokine Signaling, adherens junction, response to wounding, Hippo signaling regulation, Interferon type I signaling, PI3K-Akt signaling, leukocyte activation, Th1, and Th2 cell differentiation, NK cell-mediated cytotoxicity to name a few.Many of these pathways have been validated in various earlier investigations 36,40,42 .These ndings demonstrate the effectiveness of a network centrality-based framework to untangle the underlying players of these diseases.
Our integrative multi-omics framework has identi ed IKZF1 as a signi cant shared regulator among HS, AD, and RS.Remarkably, IKZF1 is associated with chromatin remodeling and lymphocyte differentiation which can in uence the disease pathogenesis.A considerable amount of IKZF1 interacting partners in HS interactome and co-expression networks suggest the involvement of this transcription factor in the immune environment modulation in the disease pathogenesis.Interestingly, we also found most immunesystem pathways including cellular response to cytokine stimulus, leukocyte activation, neutrophil degranulation, in ammatory response, TCR pathway, CXCR4 pathway, hemostasis, and microglia pathogen phagocytosis pathway is enriched among IKZF1 interactors.The speci c role of neutrophils, macrophages, B cells, and plasma cells has been identi ed in HS pathogenesis [42][43][44] .Similarly, CXCR4 pathways is considered important in regulating B cell functions, homing of plasma cells, and in the migration of T cells 45,46 .Interestingly, many small molecule drugs are known modulators of IKZF1 response.It is likely that these or similar molecules may nd a place in developing therapeutic innervations of these complex in ammatory skin diseases.Overall, our integrative network-based multiomics framework identi ed the underlying regulators, and pathways common among these diverse proin ammatory diseases.
This study also has limitations as available datasets for bulk transcriptome-driven/disease integrative multi-omics in ammatory skin diseases, will not generate the most reliable predictions of signi cant regulators, yet this framework as suggested here can be applied for multiple transcriptome datasets of a single disease for more accurate predictions.Furthermore, expanding the current version of transcriptome-driven integrative systeomics to the next level by implementing other multi-layered datasets including systems genomics, epigenomics and the emergence of new technologies like singlecell RNA sequencing, and spatial transcriptomics will uncover the high-resolution molecular pathogenesis, pathways, and cell-cell communications in chronic in ammatory skin disease.Figure 6 .M., B.M., and M.A. conceived the project.B.M. performed all the computational, network-based, and statistical analyses.B.M. wrote the rst draft and nal version of the manuscript.All the authors discussed the results, critically reviewed the manuscript, and provided valuable comments/edits.

Figures Figure 1 Framework
Figures

Figure 5 The
Figure 5