key: cord-0820825-zbainwlh
authors: Silva, Kushan De; Demmer, Ryan T.; Jönsson, Daniel; Mousa, Aya; Forbes, Andrew; Enticott, Joanne
title: Highly perturbed genes and hub genes associated with type 2 diabetes in different tissues of adult humans: A bioinformatics analytic workflow
date: 2022-02-10
journal: bioRxiv
DOI: 10.1101/2022.02.07.479483
sha: 698ab5d2f4b42f3480012d0c0922391ab84586ac
doc_id: 820825
cord_uid: zbainwlh

Introduction Type 2 diabetes (T2D) has a complex etiology which is not fully elucidated. Identification of gene perturbations and hub genes of T2D may assist in personalizing care. Objectives We aimed to identify highly perturbed genes and hub genes associated with T2D in different tissues of adult humans via an extensive workflow. Methods Workflow comprised five sequential steps: systematic review of NCBI GEO database; identification and classification of differentially expressed genes (DEG); identification of highly perturbed genes via meta-analysis; identification of hub genes via network analysis; downstream analyses. Three meta-analytic strategies: random effects model (REM); vote counting approach (VC); p-value combining approach (CA), were applied. Nodes having above average betweenness, closeness, and degree in the network were defined as hub genes. Downstream analyses included gene ontologies, Kyoto Encyclopedia of Genes and Genomes pathways, metabolomics, COVID-19 related genes, and Genotype-Tissue Expression profiles. Results Analysis of 27 eligible microarrays identified 6284 DEG (4592 down-regulated and 1692 up-regulated) within four tissue types. Tissue-specific gene expression was significantly greater than tissue non-specific (shared) gene expression. Meta-analysis of DEG identified 49, 27, and 8 highly perturbed genes via REM, VC, and CA, respectively, producing a compiled set of 79 highly perturbed (41 down-regulated and 38 up-regulated) genes. The 28 hub genes comprised 13 up-regulated, 9 down-regulated, and 6 predicted genes. Downstream analyses identified enrichments of: shared genes with other diabetes phenotypes; insulin synthesis and action related pathways and metabolomics; mechanistic associations with apoptosis and immunity-related pathways, COVID-19 related gene sets; and cell types demonstrating over- and under-expression of marker genes of T2D. Conclusions We identified highly perturbed genes and hub genes of T2D and revealed their associations with other diabetes phenotypes and COVID-19 as well as pathophysiological manifestations such as those related to insulin, immunity, and apoptosis. Broader utility of the proposed pipeline is envisaged.

According to an analysis of global data through years 1990-2018, diabetes was prevalent in almost half a billion people, a number expected to rise by 25% and 51%, respectively, by 2030 and 2045 [1] . Type 2 diabetes (T2D) is by far the most common, accounting for over 90% of all diabetes cases [2] , a condition found to have affected 6.3% (462 million) of the world population as in year 2017 [3] . Current evidence suggests a complex and multifactorial etiology of T2D characterized by genetic and environmental interactions [4] , although T2D pathogenesis is not yet fully elucidated. Moreover, there is evidence of varying degrees of shared genetic origins in the pathogenesis of T2D and other diabetes phenotypes such as type 1 diabetes (T1D) [5] , latent autoimmune diabetes in adults (LADA) [6] , and maturity-onset diabetes of the young (MODY) [7] . Recent studies also support the deconstruction of T2D heterogeneity to define T2D sub-types [8] and the delineation of a continuum of diabetes subtypes [9] as opposed to the status quo characterized by a few distinct diabetes phenotypes.

Understanding the genetic basis of T2D is fundamental to precision medicine approaches striving for impeccable matching and an individualized level of T2D care [10] .

Hyperglycemia in T2D is triggered by both impaired action and secretion of insulin [11] , unlike T1D the hallmark of which is a lack of insulin synthesis owing to autoimmunemediated pancreatic β -cell destruction [12] . Cardinal tissues of the body impacted by heightened insulin resistance (IR) and diminished insulin secretion in T2D include the pancreas, liver, skeletal muscle, and adipose tissue [13] . Alterations in the pancreas with T2D

include the deposition of β -amyloid [14] , reduced β -cell mass and increased β -cell apoptosis [15] , as well as decreased size, irregularities of morphology and increased fat content [16] . A salient manifestation of T2D is the compromised insulin-mediated glucose uptake predominantly by skeletal muscle and to a lesser extent by adipose tissue resulting in hyperglycemia [17] . More granular findings associated with T2D such as modifications of skeletal muscle proteins responsible for work capacity [18] and reduction of skeletal muscle capillary density and microvascular function [19] have also been reported. Adipose tissue inflammation [20] and lipotoxicity [21] have been implicated in T2D pathogenesis. These adipose tissue changes in people with T2D may be caused by numerous transcriptional and epigenetic changes [22] . Similarly, liver fat has been implicated in the pathogenesis of T2D

[23] possibly via oxidized fatty acids [24] .

Exploration of gene perturbation in T2D can deepen our understanding of the molecular etiopathology of T2D. Highly up-and down-regulated genes expressed consistently across different tissue types may uncover potential genome-wide biomarkers of the disease.

Identification of such gene signatures is also an integral component of precision diagnostic, prognostic, monitoring, and treatment approaches. Topologically, hub genes are defined as highly and tightly connected nodes in typically scale-free gene regulatory networks (GRN) [25] . As such, criteria such as high correlation in candidate modules [26] and above-average betweenness, closeness, and degree [27] in GRN have been used to demarcate hub genes in previous studies. Functionally, they perform critical regulatory roles in biological processes interacting with many other genes in associated pathways. Given their crucial structural and functional characteristics, hub genes are highly sought-after in precision medicine approaches as plausible niches for developing drug and treatment targets [28] . In this context, the importance of identifying highly perturbed genes and hub genes associated with T2D for the purpose of individualizing T2D care is unequivocal.

Downstream analyses of gene sets provide invaluable insights into associated core biological functions, pathways, diseases, drugs and many other aspects. Frequently used gene ontology (GO) analysis provide evidence as a snapshot of contemporary biological knowledge related to a given gene including its function at the molecular level, the cellular location(s) it functions at and the biological processes reliant on it [29] . Pathway analyses such as Kyoto Encyclopedia of Genes and Genomes (KEGG) [30] entail derivation of coherent and meaningful biological phenomena attributable to input genes [31] . Besides classic GO and pathway analyses, further downstream analyses such as diseases, drugs, metabolomes, and tissue enrichment are also available and can provide valuable insights into associated genes.

Taken together, downstream analysis of highly perturbed-and hub genes associated with T2D may render valuable information on aspects such as affected biological processes, dysregulated pathways, related diseases, and metabolomic biomarkers.

Advances in high-throughput technologies have generated a wealth of gene expression data whilst the availability of open-source platforms such as the National Center for Biotechnology Information Gene Expression Omnibus (NCBI GEO) [32, 33] and simultaneous advent in big data and bioinformatics analytic tools such as microarray and RNA-seq meta-analysis and gene-gene interaction network analysis strategies have offered unprecedented opportunities for high-level evidence synthesis from a multitude of gene expression datasets. Such approaches are likely to render new knowledge on complex diseases like T2D and acquire adequate statistical power to identify genes associated with a disease that may not have been evoked via prior analysis of a single or a few datasets.

Yet, to date, no comprehensive evidence synthesis study has been performed to identify highly perturbed genes and hub genes associated with T2D in human adults using an extensive bioinformatics analytic pipeline. Prior studies were limited to identifying hub genes in a few (n = 3) microarrays from a single tissue type (pancreatic islets) [34] or performing an ad hoc gene expression meta-analysis of all diabetes phenotypes (n = 13) [35] . In this study, we aimed to identify highly perturbed genes and hub genes associated with T2D in different tissues of adult humans via a pre-defined and extensive bioinformatics analytic workflow consisting of systematic review, meta-analysis, identification and classification of differentially expressed genes (DEG), network analysis and downstream analysis.

The methodological approach is summarized in S1 Figure and consisted of five sequential steps: (1) Systematic review of NCBI GEO expression data and related publications (2) Analysis of microarrays to identify DEG (3) Meta-analysis of DEG to identify highly perturbed genes in T2D (4) Network analysis of gene-gene/protein-protein interactions to identify hub genes in T2D (5) Downstream analysis of highly perturbed genes and hub genes associated with T2D.

A preliminary search on the NCBI GEO database was first run on 1 st June 2021 using a pre- 

The 27 microarrays were imported to R using 'getGEO' function of 'GEOquery' package [36] . In each expression set, 'phenoData' component was examined to determine the number of eligible samples and confirm the presence of outcome (T2D vs controls) variable whilst 'featureData' component was explored to verify the presence of gene annotation information.

Where required, non-Normalized gene expression matrices were log 2 transformed in order to alleviate skewness and create symmetric distributions [37] . Samples were assessed for the presence of any batch effects between the two groups by running principal components analysis on transposed expression matrices and were rectified using 'removeBatchEffect' function in 'limma' package [38] . Relevant features were annotated with expression matrices to generate curated data for running DEG analysis. Samples were ascribed to the relevant group (T2D/control) using 'model.matrix' function of 'limma' package [38] producing a binary design matrix. As the detection of DEG can be enhanced by filtering genes with a low expression level, we assumed a median (50%) cut-off for the gene expression level. A uniform analytic pipeline consisting of the following sequential steps were applied to each of the 27 microarrays: (1) Median expression levels were calculated and those above the median were retained. (2) From the resulting genes, those expressed in more than two samples were retained while the others were removed. (3) Model fitting was performed using 'lmFit' function of 'limma' [38] to enumerate expression levels of T2D and control groups. (4) Contrasts were defined as 'T2D -control' and 'empirical Bayes' step was run to derive differential expression results (5) The DEG, defined as those with log 2 FC > 1 & Benjamini-Hochberg (BH)-adjusted p < 0.05 for up-regulated genes and log 2 FC < -0.5 & BH-adjusted p < 0.05 for down-regulated genes, were identified for each microarray.

The number of DEG identified by each dataset is shown in Table 1 . Microarrays with no DEG (n = 11) were excluded and those with non-zero DEG (n = 16) were visualized with a stacked bar chart (Figure 1) . Finally, DEG were classified by tissue types (n = 4; circulatory/ adipose/ digestive/ skeletal muscle) and visualized via a Venn diagram (Figure 2 ). Clinical and other information of the 16 microarrays with non-zero DEG is presented in S3 Table. Details of DEG identified by each dataset are presented in S4 Table, while DEG classified by tissue type are presented in S5 Table. 2

In order to identify highly perturbed genes associated with T2D, we conducted meta-analysis of DEG using 'MetaVolcanoR' package [39] . We implemented all 3 meta-analytic strategies incorporated in this package, namely, random effects model (REM), vote counting approach (VC) and p-value combining approach (CA).

In brief, REM synthesizes a summary fold change of multiple microarrays based on variance, producing a summary p-value which indicates the probability that the summary fold-change is not different to zero. The 'metathr' parameter can be specified to filter the desired percentage of the top-most consistently perturbed genes. Gene perturbation is ranked as per the 'topconfects' approach [40] . The VC algorithm produces DEG according to user-specified p-values and fold change cut-off levels, taking into account both the number of studies in which a DEG appeared and its gene fold change sign consistency. Here also, 'metathr' parameter can be defined to extract the required percentage of highly perturbed genes. Metasynthesis of gene perturbation by CA algorithm is at the mean or median level along with pvalues derived by Fisher method. A required proportion of top-most DEG can be identified by specifying 'metathr' parameter with CA as well [39] .

As required by the package, the 16 datasets each consisting of the columns: gene name (Symbol), fold change (log 2 FC), p-value, and confidence intervals of the fold change (CI.L & CI.R) were merged to build a list item. For all 3 meta-analytic models, 'metathr' was set at 0.01. For VC, p-value and absolute fold change thresholds were set at 0.05 and 0, respectively. Highly perturbed genes identified by each model as well as the compiled list of all highly perturbed genes (n = 79) is presented in Table 2 . Volcano plots illustrating the topmost perturbed genes identified by REM, VC, and CA are given in Figures 3 -5, respectively. Inverse cumulative distribution of consistently differentially expressed genes as per VC was plotted ( Figure 6 ) to demonstrate the number of genes with perturbed expression in ≥ 1 studies. We present detailed meta-analytic output from all 3 approaches in S6 Table   and highly perturbed genes identified by each method in S7 Table. 2

The list of highly perturbed genes (n = 79) was fed into GENEMANIA [41] to determine the gene-gene interactions network. The interactions network formulated by the GENEMANIA gene function prediction program, based on the Multiple Association Network Integration Algorithm (MANIA), incorporates a multitude of functional associations including coexpression, pathways, physical interactions, co-localization, genetic interactions, and protein domain similarity. It has been found more accurate and computationally-efficient than other gene function prediction methods [42, 43] . The gene-gene interactions network was visualized in GENEMANIA (Figure 7 ) while the details of the constructed network were given in S8 Table. Next, these interactions were imported to visualize the protein-protein interaction (PPI) network (Figure 8 ) in STRING version 11.0 [44] . We used the 'Centiscape'

application [45] in the 'Cytoscape' software [46] to analyze the PPI network and determine the hub nodes. After removing nodes with a connection number < 2, the network was visualized in 'Cytoscape' (Figure 9 ). Hub genes were defined as those nodes in the network with betweenness, closeness, and degree higher than their mean values. A similar approach has been used in previous studies to demarcate hub genes [27] . Topological features of the 28 hub nodes are summarized in Table 3 while details of the PPI network derived by 'Centiscape' are provided in S9 Table.

Using 'Enrichr' [47] [48] [49] platform, we ran a series of downstream analyses for both highly perturbed genes (n = 79) and hub genes (n = 28) as outlined below.

• Ontologies: GO Biological Process 2018 [29] • Pathways: KEGG 2019 Human [30] • Diseases/Drugs: COVID-19 Related Gene Sets The Genotype-Tissue Expression (GTEx) portal contains tissue-specific gene expression and regulation data [50] , whereas the Human Metabolome Database (HMDB) records human metabolomics data [51] . Findings from downstream analyses of highly perturbed genes are summarized in Tables 4 -6 and illustrated in Figures 10 -12 while details given in S10 Table. Results from downstream analyses of hub genes are presented in Tables 7 -9 and visualized in Figures 13 -15 with details in S11 Table. 3. RESULTS

Of all the DEG identified by different tissues (n = 6284), the proportion of down-regulated genes (n = 4592) were significantly higher than the proportion of upregulated genes (n = Figure 1 ).

Proportions of tissue-specific DEG were significantly higher than the proportions shared with one or more other tissue types in circulatory (242 specific vs 89 non-specific; p < 0.000001), adipose (1989 specific vs 342 non-specific; p < 0.000001), and digestive (1875 specific vs 409 non-specific; p < 0.000001) groups, while skeletal muscle (244 specific vs 217 nonspecific; p = 0.112959) had no such over-expression (Figure 2 ).

There was no significant difference (p = 0.410983) between the proportions of up-(38/79) and down-regulated (41/79) genes that constituted the highly perturbed gene set (n = 79) identified by all 3 meta-analytic algorithms ( Table 2 C19orf33, HLA-DRB4, OR8B12, HLA-DRB5, IFNA7, LOC389286, FUT11,   APOL4, LOC731682, KIAA1984, POPDC3, TAP2, CTSC, DYRK2, MCOLN3, TMEM37,   API5, ARG2, C14orf132, COG2, DHRS2, ENPP2, ENTPD3, HADH, KIAA1279, LARP4,   MARK1, MTRR, NAALAD2, PAAF1, PPM1K, PPP1R1A, RPL14, SLC2A2, SNAP25, STMN2, NMNAT2, ASCL2, and IAPP.

genes.

The 28 hub genes consisted of 13 up-regulated, 9 down-regulated and 6 predicted genes.

There was no significant difference ( (Table 3) .

Ontology analysis of highly perturbed genes identified 89 enriched biological processes (S10 Table) . Top Table 4 ).

Ontological analysis of hub genes identified 179 enriched biological processes (S11 Table) . Table 7) .

As per KEGG analysis of highly perturbed genes, 21 pathways were enriched (S10 Table) .

The top 10 pathways were: Epstein-Barr virus infection; Antigen processing and presentation; Table 4 ).

According to KEGG analysis of hub genes, 12 pathways were enriched (S11 Table) . Table 7) .

Downstream analyses revealed 20 COVID-19 related gene sets associated with the highly perturbed genes of T2D (S10 Table) and the top 10 of these are visualized in Figure 11 (a).

There were 23 COVID-19 related gene sets associated with the hub genes of T2D (S11 Table) , the top 10 of which are visualized in Figure 14 (a).

Four HMDB metabolites [Zinc (HMDB01303), Manganese (HMDB01333), Magnesium (HMDB00547), C10H13N2O7P (HMDB01570)] were associated with the highly perturbed genes of T2D (S10 Table) . Of these, Zinc (HMDB01303) was associated with up-regulated genes in T2D while C10H13N2O7P (HMDB01570) was associated with down-regulated genes in T2D ( Table 5) .

There were 45 HMDB Metabolites associated with the hub genes of T2D (S11 Table) . The (Table 8 ).

There were 168 down-regulated GTEx profiles (S10 Table) associated with highly perturbed genes of T2D, the top 10 of which consisted of one brain, two thyroid, and seven blood tissue samples ( Table 6 & Figure 12 ). There were 77 up-regulated GTEx profiles (S10 Table) associated with highly perturbed genes of T2D, the top 10 of which consisted of six adiposeand one each of lung, breast, heart, and blood vessel tissue samples ( Figure 12 ).

There were 223 down-regulated GTEx profiles (S11 Table) associated with hub genes of T2D, the top 10 of which consisted of five nervous system (three brain, one nerve, one pituitary), three blood, one breast and one bladder tissue samples ( Table 9 & Figure 15 ).

There were 178 up-regulated GTEx profiles (S11 Table) associated with hub genes of T2D, the top 10 of which consisted of six adipose-, two blood and two skin tissue samples ( Table 9 & Figure 15 ).

In this study, we identified highly perturbed genes and hub genes associated with T2D in different tissues of adult humans, via an extensive bioinformatics analytic workflow. By doing this, we revealed valuable insights with respect to T2D pathogenesis, including associations with other diabetes phenotypes and COVID-19, patterns of tissue-specific and tissue non-specific differential gene expression as well as pathophysiological manifestations such as those related to insulin action, immunity, and apoptosis. Salient findings of the study which contribute towards the understanding of the genetic basis of T2D are further discussed below. The comprehensive evidence synthesis approach with open-source gene expression data exemplified in this study can be replicated to gain high-level evidence synthesis for other clinical conditions.

Our findings indicate that T2D seems rather a disorder of gene down-regulation than upregulation, when the whole genome is considered. This is consistent with previous studies where a preponderance of gene down-regulation was associated with T2D (52, 53) . Also, hyperglycemia-induced global downregulation of gene expression in adipose and skeletal muscle tissues have been documented previously (54) . A similar pattern has been observed in T1D (55) as well as other endocrine disorders such as polycystic ovary syndrome (56) . In contrast, highly perturbed genes and hub genes associated with T2D, which might together constitute the candidate gene set critical for pathogenicity of T2D, were found to contain both up-and down-regulated genes. This presentation suggests a more complex dysregulation at the crux of the GRN of T2D, involving actions and interactions between both repressed and augmented genes.

Results of the present study indicate the predominance of tissue-specific DEG in T2D which may have important implications for guiding biomarker discovery process. It supports the use of target tissue gene expression analysis as a viable avenue for identifying tissue-specific T2D biomarkers. A previous analysis integrating multiple tissue transcriptomics and PPI data to explore molecular biomarkers of T2D confirmed the presence of common signatures (57) .

We also observed common DEG across different tissue types which can act as confluent molecular signatures of T2D. Identification of tissue-specific and non-specific molecular gene signatures of T2D facilitates downstream exploration of key pathways amenable to therapeutic targeting and drug repurposing efforts.

As revealed by KEGG pathway analyses, both MODY and T1D were enriched pathways associated with highly perturbed genes of T2D. In addition, MODY was an enriched pathway associated with hub genes of T2D as well. Specifically, SLC2A2 (GLUT2) and IAPP genes were commonly enriched in MODY, while HLA-DRB4 and HLA-DRB5 genes were underlying the enrichment with T1D. A gene expression meta-analysis also revealed the existence of possible pleiotropic mechanisms manifest via common gene signatures (PGRMC1 and HADH) across different diabetes phenotypes [35] .

Down-regulation of SLC2A2 is associated with not only T2D [58] but also neonatal diabetes [59] and early childhood diabetes [60] suggesting a likely role in insulin secretion. Amylin (IAPP), a gluco-modulatory hormone co-expressed with insulin by pancreatic β cells, is down-regulated in both T1D and advanced T2D [61] while amylin agonists are considered as novel therapeutic agents for treating diabetes [62] . It has also been found that human amylin plays a protective role against autoimmune diabetes inducing CD4+Foxp3+ regulatory T cells [63] . Downregulation of HLA-DRB4 in peripheral blood mononuclear cells (PBMC) is associated with T2D as well as dyslipidemia and periodontitis [64] , while a meta-analysis revealed that the lack of HLA-DRB5 increased T2D risk [65] . Intriguingly, both HLA-DRB4

and HLA-DRB5 are associated with β cell autoantibodies and T1D [66] . In fact, previous studies reported that both T1D and T2D share HLA class II locus components [65] .

Interestingly, two of the hub genes of T2D found in our study (MMP9, ARG1) have been found as hub genes of T1D in a previous analysis [67] . Collectively, these findings support some degree of shared genetic architecture between T2D and other diabetes pathologies.

Downstream analyses provided insights into the characterization of T2D as a disorder of insulin secretion and action. We found that insulin secretion was the most significant KEGG pathway associated with hub genes, whereby two down-regulated hub genes (SLC2A2, SNAP25) in T2D were underlying this enrichment. Zinc, the most significant HMDB metabolite associated with both highly perturbed genes and hub genes of T2D, is an essential element with key regulatory roles in insulin synthesis, storage, and secretion [68] . Other metabolomic biomarkers associated with highly perturbed genes included magnesium which is necessary for insulin signaling [69] as well as manganese which is involved in insulin synthesis and secretion [70] . Together, these findings underscore the effects on insulin production and action as pivotal to T2D pathogenesis.

Downstream analyses revealed that multiple GO and KEGG pathways associated with apoptosis, including intrinsic apoptotic signaling pathway, were enriched in T2D. It is known that hyperglycemia-induced β cell apoptosis, a hallmark in T2D progression, occurs via intrinsic pathways causing reduced islet mass and metabolic abnormalities [71] . Hyperglycemia-induced apoptosis has been reported to occur in other sites such as renal cells [72] and coronary arteries [73] , indicating a possible role in disease progression and the onset of complications.

Downstream analyses also revealed that multiple immunity-related GO and KEGG pathways, encompassing both innate and humoral immune responses, were enriched in T2D. These included multiple ontologies involving neutrophils, antigen processing and presentation as well as severe immune reactions such as graft-versus-host disease and allograft rejection.

Impaired immunity in T2D and consequent susceptibility to infections and complications is frequently observed [74] . A deeper understanding of the genomics underlying impaired immunity in T2D might provide opportunities to personalize the management of comorbidities and pharmacotherapy.

Epidemiological studies strongly suggest poorer prognosis of COVID-19 among people with T2D [75] , although underlying mechanisms are not well-understood [76] . Downstream analysis of highly perturbed genes and hub genes of T2D in the present study revealed a large number of enriched COVID-19 related gene sets, providing support for this putative link at a more granular level.

Downstream analysis of GTEx profiles identified tissues that are likely to demonstrate underand over-expression of DEG associated with T2D. Findings indicate that adipose tissue tends to over-express marker genes of T2D, while these might be under-expressed in other tissues such as those of the nervous system. These findings have implications for biomarker discovery and can guide further research on tissues which should be explored for DEG identification.

Taken together, these findings contribute towards the understanding of the genetic basis of T2D, which is fundamental for precision medicine approaches for individualized T2D care.

Further research is warranted to explore and substantiate the molecular mechanisms underlying these findings which would be critical for establishing precision T2D medicine 

The funder was not involved in the design of the study; the collection, analysis, and interpretation of data; writing the report; and did not impose any restrictions regarding the publication of the report.

Authors declare that there are no conflicts of interest. C19orf33, HLA-DRB4, OR8B12, HLA-DRB5, IFNA7, LOC389286, FUT11,  APOL4,  LOC731682,  KIAA1984,  POPDC3, TAP2, CTSC, DYRK2,  MCOLN3, TMEM37, API5, ARG2,   ALDOB, PRIMA1, PVRL2, CNTNAP2,  RASL11B, POMZP3, ELFN1, ESPNL,  PI3,  SCN1B,  EGR2,  MGRN1,  SLC9A3R2,  LOC644422,  MCL1,  ZBTB16,  LOC650885,  UGT2B7,  C19orf33, HLA-DRB4, OR8B12, HLA-DRB5, IFNA7, LOC389286, FUT11,  APOL4,  LOC731682,  KIAA1984,  POPDC3, TAP2, CTSC, DYRK2,  MCOLN3, TMEM37, API5, ARG2 TNFAIP6, TAPBP, PHACTR3,  PHLDA1, ALDOB, PRIMA1,  PVRL2, CNTNAP2, RASL11B,  POMZP3, ELFN1, ESPNL, PI3,  SCN1B,  EGR2,  MGRN1,  SLC9A3R2, LOC644422, MCL1,  ZBTB16   C14orf132, COG2, DHRS2, ENPP2,  ENTPD3, HADH, KIAA1279, LARP4,  MARK1, MTRR, NAALAD2, PAAF1,  PPM1K, PPP1R1A, RPL14, SLC2A2,  SNAP25 , STMN2, NMNAT2, ASCL2, IAPP 

Global and regional diabetes prevalence estimates for 2019 and projections for 2030 and 2045: Results from the International Diabetes Federation Diabetes Atlas

Global aetiology and epidemiology of type 2 diabetes mellitus and its complications

Epidemiology of type 2 diabetes-global burden of disease and forecasted trends

tRNA Biology in the Pathogenesis of Diabetes: Role of Genetic and Environmental Factors

Shared genetic risk contributes to type 1 and type 2 diabetes etiology

Overlap of genetic susceptibility to type 1 diabetes, type 2 diabetes, and latent autoimmune diabetes in adults

Pathogenic variants in actionable MODY genes are associated with type 2 diabetes

Type 2 diabetes genetic loci informed by multi-trait associations point to disease mechanisms and subtypes: A soft clustering analysis

Common and rare forms of diabetes mellitus: towards a continuum of diabetes subtypes

Precision medicine in type 2 diabetes

Mechanisms controlling pancreatic islet cell function in insulin secretion

A machine learning approach to predict response to immunotherapy in type 1 diabetes

Defining the underlying defect in insulin action in type 2 diabetes

Pancreastatin induces islet amyloid peptide aggregation in the pancreas, liver, and skeletal muscle: An implication for type 2 diabetes

Conversion of the death inhibitor ARC to a killer activates pancreatic β cell death in diabetes

Altered volume, morphology and composition of the pancreas in type 2 diabetes

Glucose transporters in adipose tissue, liver, and skeletal muscle in metabolic health and disease

Variation in pleiotropic hub gene expression is associated with interspecific differences in head shape and eye size in Drosophila

Identification of Hub Genes and Key Pathways Associated With Bipolar Disorder Based on Weighted Gene Coexpression Network Analysis

Bioinformatics prediction and analysis of hub genes and pathways of three types of gynecological cancer

Bioinformatic analysis of potential hub genes in gastric adenocarcinoma

Gene Ontology annotations: what they mean and where they come from

KEGG: kyoto encyclopedia of genes and genomes

Identifying significantly impacted pathways: a comprehensive review and assessment

Gene Expression Omnibus: NCBI gene expression and hybridization array data repository

NCBI GEO: archive for functional genomics data sets--update

Identification of Hub Genes in Type 2

Tissue Non-Specific Genes and Pathways Associated with Diabetes: An Expression Meta-Analysis

GEOquery: a bridge between the Gene Expression Omnibus (GEO) and BioConductor

Identification of gene expression signatures for psoriasis classification using machine learning techniques

limma powers differential expression analyses for RNA-sequencing and microarray studies

Gene Expression Meta-analysis Visualization Tool

Topconfects: a package for confident effect sizes in differential expression analysis provides a more biologically useful ranked gene list

The GeneMANIA prediction server: biological network integration for gene prioritization and predicting gene function

GeneMANIA: a realtime multiple association network integration algorithm for predicting gene function

S4

A critical assessment of Mus musculus gene function prediction using integrated genomic evidence

STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets

Analyzing biological network parameters with CentiScaPe

Cytoscape: a software environment for integrated models of biomolecular interaction networks

Enrichr: interactive and collaborative HTML5 gene list enrichment analysis tool

Enrichr: a comprehensive gene set enrichment analysis web server 2016 update

The Genotype-Tissue Expression (GTEx) project

HMDB 4.0: the human metabolome database for 2018

Genome wide analysis of gene expression changes in skin from patients with type 2 diabetes

Gene expression in skeletal muscle biopsies from people with type 2 diabetes and relatives: differential regulation of insulin signaling pathways

Acute hyperglycemia induces a global downregulation of gene expression in adipose tissue and skeletal muscle of healthy subjects

Gene Expression Analysis of the Pre-Diabetic Pancreas to Identify Pathogenic Mechanisms and Biomarkers of Type 1 Diabetes

Meta-analysis of gene expression profiles of lean and obese PCOS to identify differentially regulated pathways and risk of comorbidities

Tissue-Specific Molecular Biomarker Signatures of Type 2 Diabetes: An Integrative Analysis of Transcriptomics and Protein-Protein Interaction Data

Systems biology of the IMIDIA biobank from organ donors and pancreatectomised patients defines a novel transcriptomic signature of islets from individuals with type 2 diabetes

SLC2A2 mutations can cause neonatal diabetes, suggesting GLUT2 may have a role in human insulin secretion

A novel gene in early childhood diabetes: EDEM2 silencing decreases SLC2A2 and PXD1 expression, leading to impaired insulin secretion

Mechanisms of islet amyloidosis toxicity in type 2 diabetes

Mono and dual agonists of the amylin, calcitonin, and CGRP receptors and their potential in metabolic diseases

Human amylin induces CD4+Foxp3+ regulatory T cells in the protection from autoimmune diabetes

Circulating lymphocytes and monocytes transcriptomic analysis of patients with type 2 diabetes mellitus, dyslipidemia and periodontitis

HLA Class II Allele Analyses Implicate Common Genetic Components in Type 1 and Non-Insulin-Treated Type 2 Diabetes

Next-Generation Sequencing Reveals That HLA-DRB3, -DRB4, and -DRB5 May Be Associated With Islet Autoantibodies and Risk for Childhood Type 1 Diabetes

Analysis of potential hub genes involved in the pathogenesis of Chinese type 1 diabetic patients

Effect of zinc supplementation on insulin resistance and metabolic risk factors in obese Korean women

Magnesium in Obesity, Metabolic Syndrome, and Type 2 Diabetes

Mechanisms and Disease Pathogenesis Underlying Metal-Induced Oxidative Stress

Linking metabolic abnormalities to apoptotic pathways in Beta cells in type 2 diabetes

Apoptosis occurs differentially according to glomerular size in diabetic kidney disease

High glucose-induced apoptosis in human coronary artery endothelial cells involves up-regulation of death receptors

Type 2 Diabetes and its Impact on the Immune System

Diabetes Epidemiology in the COVID-19 Pandemic

COVID-19 in people with diabetes: understanding the reasons for worse outcomes