Metabolic Footprint, towards Understanding Type 2 Diabetes beyond Glycemia

Type 2 diabetes (T2D) heterogeneity is a major determinant of complications risk and treatment response. Using cluster analysis, we aimed to stratify glycemia within metabolic multidimensionality and extract pathophysiological insights out of metabolic profiling. We performed a cluster analysis to stratify 974 subjects (PREVADIAB2 cohort) with normoglycemia, prediabetes, or non-treated diabetes. The algorithm was informed by age, anthropometry, and metabolic milieu (glucose, insulin, C-peptide, and free fatty acid (FFA) levels during the oral glucose tolerance test OGTT). For cluster profiling, we additionally used indexes of metabolism mechanisms (e.g., tissue-specific insulin resistance, insulin clearance, and insulin secretion), non-alcoholic fatty liver disease (NAFLD), and glomerular filtration rate (GFR). We found prominent heterogeneity within two optimal clusters, mainly representing normometabolism (Cluster-I) or insulin resistance and NAFLD (Cluster-II), at higher granularity. This was illustrated by sub-clusters showing similar NAFLD prevalence but differentiated by glycemia, FFA, and GFR (Cluster-II). Sub-clusters with similar glycemia and FFA showed dissimilar insulin clearance and secretion (Cluster-I). This work reveals that T2D heterogeneity can be captured by a thorough metabolic milieu and mechanisms profiling—metabolic footprint. It is expected that deeper phenotyping and increased pathophysiology knowledge will allow to identify subject’s multidimensional profile, predict their progression, and treat them towards precision medicine.


Introduction
It is widely accepted that type 2 diabetes is a multifactorial and heterogeneous condition progressing through multiple pathophysiological mechanisms [1]. Despite recent therapeutic advances and enormous expenditure in treatment [2] type 2 diabetes carries a heavy socioeconomic burden [3] and is a leading cause of renal failure, blindness, cardiovascular events, and major amputations of the lower limbs [4]. Current type 2 diabetes diagnosis criteria and therapeutic targets are mainly focused

Subjects
The study population comprises participants of a diabetes prevalence study performed in Portugal (PREVADIAB2). Between 2008 and 2009, the PREVADIAB1 study recruited 5167 subjects attending the national health care system across the country that were screened for diabetes status [13]. Normoglycemia, prediabetes, and diabetes were defined using the International Diabetes Federation (IDF)/World Health Organization (WHO) glycemia criteria [14,15]. In 2014, we selected within this cohort 1088 subjects that did not have diabetes in PREVADIAB1 and attended 55 health units evenly representing the geographical distribution of the Portuguese population (PREVADIAB2; Figure 1). For all participants, a letter was sent containing information for study participation. . After the application of exclusion criteria and preprocessing the data we performed a cluster analysis of 974 subjects. Finally, we profiled the clusters with several metabolic parameters. T2D-Type 2 Diabetes; SOM-Self-Organizing Map.

Clinical Measurements
Clinical history of the participants in PREVADIAB2 was assessed, including anti-diabetic medication. Participants were weighed on an electronic scale, to the nearest 0.01 kg while wearing minimal clothes and without shoes. Waist circumference was measured with a flexible ruler, at the midline between the top of the iliac crest and the inferior limit of the ribs [16]. Height was measured to the nearest 0.1 cm with a stadiometer according to the standardized WHO procedures. Blood collection was performed, in a seated position from the antecubital vein at rest, after an overnight fast, at 0, 30 and 120 min of the oral glucose tolerance test (OGTT). The blood was collected into dry tubes and into tubes containing ethylenediamine-tetraacetic acid (EDTA) as an anticoagulant. Biological samples were centrifuged at 3500 rpm at 4 °C for 10 min and plasma samples were frozen at −80 °C for posterior analysis. Serum samples were used to analyze the lipid profile of the participants, including total cholesterol, low density lipoprotein (LDL-c), high density lipoprotein (HDL-c) cholesterol, free fatty acids (FFA), and triglycerides (TG) and to quantify creatinine, gamaglutamiltranspeptidase (GGT), alanine aminotransferase (ALT), and aspartate aminotransferase (AST) using colored enzymatic tests in an automated analyzer (Auto Analyzer Olympus AU640, Beckman Coulter, Brea, CA, USA). Glycated hemoglobin was determined by high performance liquid chromatography (HPLC) with boronate affinity (Menarini Premier Hb 9210, Paço de Arcos, Lisboa, Portugal). Plasma glucose levels were measured with a glucose analyzer (glucose oxidase method, Olympus AU640, Beckman Coulter, Brea, California, USA). Plasma insulin and C- Figure 1. Methodology flowchart. In 2014, 1088 subjects from the PREVADIAB1 study without T2D in 2008-2009, were reevaluated 5 years later (PREVADIAB 2 cohort). After the application of exclusion criteria and preprocessing the data we performed a cluster analysis of 974 subjects. Finally, we profiled the clusters with several metabolic parameters. T2D-Type 2 Diabetes; SOM-Self-Organizing Map.

Clinical Measurements
Clinical history of the participants in PREVADIAB2 was assessed, including anti-diabetic medication. Participants were weighed on an electronic scale, to the nearest 0.01 kg while wearing minimal clothes and without shoes. Waist circumference was measured with a flexible ruler, at the midline between the top of the iliac crest and the inferior limit of the ribs [16]. Height was measured to the nearest 0.1 cm with a stadiometer according to the standardized WHO procedures. Blood collection was performed, in a seated position from the antecubital vein at rest, after an overnight fast, at 0, 30 and 120 min of the oral glucose tolerance test (OGTT). The blood was collected into dry tubes and into tubes containing ethylenediamine-tetraacetic acid (EDTA) as an anticoagulant. Biological samples were centrifuged at 3500 rpm at 4 • C for 10 min and plasma samples were frozen at −80 • C for posterior analysis. Serum samples were used to analyze the lipid profile of the participants, including total cholesterol, low density lipoprotein (LDL-c), high density lipoprotein (HDL-c) cholesterol, free fatty acids (FFA), and triglycerides (TG) and to quantify creatinine, gamaglutamiltranspeptidase (GGT), alanine aminotransferase (ALT), and aspartate aminotransferase (AST) using colored enzymatic tests in an automated analyzer (Auto Analyzer Olympus AU640, Beckman Coulter, Brea, CA, USA). Glycated hemoglobin was determined by high performance liquid chromatography (HPLC) with boronate affinity (Menarini Premier Hb 9210, Paço de Arcos, Lisboa, Portugal). Plasma glucose levels were measured with a glucose analyzer (glucose oxidase method, Olympus AU640, Beckman Coulter, Brea, CA, USA). Plasma insulin and C-peptide levels were determined by commercial chemiluminescence assays, Liaison (DiaSorin, Salugia, Pidemont, Italy).

Statistical Analysis
Comparison between genders for parameters informing the cluster algorithm and profiling the clusters was performed with Mann-Whitney test for continuous variables and with chi-square test for categorical variables.

Cluster Analysis
To inform the clustering algorithm and stratify the population we chose age, BMI, waist circumference, glycemia, insulinemia, C-peptide, and FFA measurements at three time points (0 min, 30 min and 120 min of the OGTT). We reasoned that these measurements, available in the clinic, would reflect clinically relevant variations in the metabolic milieu and are sensitive descriptors of dysglycemia heterogeneity.
The methodology workflow used for clustering and cluster profiling is outlined in Figure 1. In the PREVADIAB2 cohort 58 subjects under anti-diabetes medication were excluded from analysis. Likewise, further 59 subjects with missing values on variables informing the clustering algorithm were also excluded. BMI and waist circumference were scaled by min-max standardization by gender. We did not exclude outliers. Odd values of measured parameters were confirmed by a repeated measurement and outliers may be part of the dysglycemia heterogeneity we aim to capture. Furthermore, cluster analysis is an appropriate tool to identify outliers on a multidimensional level.
We used a hierarchical self-organizing map followed by a hierarchical clustering algorithm to stratify the population [28]. This methodology has been used to solve several complex clustering problems, namely in clustering analysis of type 1 diabetes complications [29]. The first algorithm allows data dimensionality reduction. We reduced data dimensionality to 27 units, using an implementation of multilayer self-organizing mapping in R (kohonen package [30]); considering that there are three accepted classes for fasting glycemia and for glucose tolerance, and that we used three time points values of the OGTT (0 min, 30 min and 120 min), 27 are the possible combinations. We used five toroidal-grids, each of which with 27 hexagonal units. The first four included three time point values of the OGTT of glycemia, C-peptide, insulinemia, and FFA, respectively, while the fifth included age, BMI, and waist circumference. After tuning, the remaining algorithm hyperparameters were set as follows: distance metric used was Euclidean distance; the algorithm was randomly initialized, and trained sequentially (500 iterations); default learning rate was used, decreasing linearly from 0.05 to 0.01; Gaussian neighborhood function, starting with default 2/3 of units, and decreasing linearly. The algorithm was run 1000 times and the model with the lowest quantization error was selected.
The resulting 27 units were then clustered and explored with the hierarchical clustering algorithm implementation in R, using the Ward method (hclust function, stats package [31]), which allows to explore the aggregation of SOM units in a hierarchical manner using different granularity (cut levels). Optimal number of clusters for the 27 SOM units was determined with the silhouette width (factoextra package [32]). Cluster stability was assessed by running the described analysis on 100 bootstrapped resamples and calculating the mean Jaccard coefficient [33].
We profiled the clusters for the overall population as well as separately for both genders, by assigning the median of the variable of each cluster to their distribution in the overall population (percentiles 10th, 25th, 75th, and 90th). In profiling the clusters, we used several surrogate indexes of important pathophysiological mechanisms and traits, in addition to the parameters informing the clustering.

Results
We analyzed 974 subjects of the PREVADIAB2 Portuguese cohort with median age of 62 years and 60% women. According to IDF/WHO guidelines for dysglycemia, 72% had normoglycemia, 22% prediabetes, and 6% diabetes, with no gender difference (p = 0.08). Table 1 summarizes descriptive statistics of the variables informing cluster analysis and cluster profiling. BMI and waist circumference scaled by gender, age, and the three time points (0 min, 30 min and 120 min of OGTT) of glycemia, insulinemia, C-peptide, and FFA levels were chosen to inform the clustering procedure. We reasoned that these measurements can be obtained in clinical practice, reflect clinically relevant variations in the metabolic milieu and together are sensitive descriptors of dysglycemia heterogeneity. In addition, cluster profiling incorporated surrogate indices of relevant dysmetabolism mechanisms. We assigned the clusters median of the profiling variables to the corresponding percentile of their distribution in the overall cohort (10th, 25th, 75th, and 90th percentile). Overall, the 25th-75th interval encompassed the accepted normal values of the profiling variables. Firstly, we reduced data dimensionality to 27 units, using a SOM algorithm (Supplemental Figure  S1). To resolve the phenotypic variation of the 27 SOM units we used a hierarchical clustering algorithm, capturing subject's similarity at different hierarchical aggregation levels. After the training procedure we found that the optimal number of clusters was two (silhouette index) with cluster stability (mean Jaccard coefficient 0.8).
The first cutting level of hierarchical clustering of the 27 SOM units identified two clusters ( Figure 2; K = 2), where Cluster-I (C-I) included the majority of subjects under analysis (77%) with high representation of normoglycemia. In C-I, the median of all profiling variables fell within the range of the 25th to 75th percentiles of the global population. In contrast, Cluster-II (C-II) showed higher BMI, waist circumference, glycemia, insulinemia, and C-peptide levels, and therefore is associated with insulin resistance and fatty liver (median in 75th-90th percentiles). An analysis by gender showed similar profiling patterns, though expected differences in some variables such as BMI and waist circumference were detected (Supplemental Figures S2 and S3). We noticed that C-I had a lower proportion of subjects with glycemia but, interestingly their absolute count was higher than C-II (145 subjects in C-I vs. 122 in C-II). This indicates that both clusters contain a significant number of individuals with glycemia suggesting an underlying phenotypic heterogeneity which is not captured by the optimal number of clusters. Our analysis proceeded towards higher granularity to an intermediate hierarchical cut level (K = 11) and full cluster discrimination (K = 27). The initial breakdown of C-I in three clusters (C-A, C-B, and C-C) reveals very different profiles ( Figure 2; K = 11). While C-C represents ~60% of C-I with large representation of normoglycemia, subjects with glycemia were overrepresented in C-B and showed a distinct profile with impaired insulin secretion and increased insulin clearance. Further C-B discrimination (K = 27, C-1 to C-27) distinguished Clusters 6 (C-6) and 7 (C-7) that retained the combined insulin secretion and insulin clearance unbalances but are differentiated by the increased free fatty acids levels in C-7 but not in C-6 ( Figure 3; Supplementary Figures S4 and S5 for males and females, respectively). These two clusters represent 9% of the all cohort and uncovered a distinct pathophysiological trajectory entailing impaired insulin secretion and increased insulin clearance that can be associated or not with lipid dysmetabolism. Cluster 8 (C-8) did not show unbalances of evaluated dysmetabolism mechanisms, nonetheless subjects were exposed to increased median levels of both glycemia and FFA. This profile may represent individuals that develop glycemia and dyslypidemia by mechanisms that are uncommonly assessed. In C-A ( Figure 2; K = 11), the median levels of the evaluated mechanisms fell within the 25th and 75th percentiles of the overall population, except for insulin secretion. Nevertheless, C-A harbored a cluster (Cluster 4) representing 3% of the cohort with normal insulin secretion but with hyperinsulinemia, which is associated with insulin clearance suppression and fatty liver. Noteworthy, in this cluster, median glycemic levels were increased at 30 minutes  The initial breakdown of C-I in three clusters (C-A, C-B, and C-C) reveals very different profiles ( Figure 2; K = 11). While C-C represents~60% of C-I with large representation of normoglycemia, subjects with dysglycemia were overrepresented in C-B and showed a distinct profile with impaired insulin secretion and increased insulin clearance. Further C-B discrimination (K = 27, C-1 to C-27) distinguished Clusters 6 (C-6) and 7 (C-7) that retained the combined insulin secretion and insulin clearance unbalances but are differentiated by the increased free fatty acids levels in C-7 but not in C-6 ( Figure 3; Supplementary Figures S4 and S5 for males and females, respectively). These two clusters represent 9% of the all cohort and uncovered a distinct pathophysiological trajectory entailing impaired insulin secretion and increased insulin clearance that can be associated or not with lipid dysmetabolism. Cluster 8 (C-8) did not show unbalances of evaluated dysmetabolism mechanisms, nonetheless subjects were exposed to increased median levels of both glycemia and FFA. This profile may represent individuals that develop dysglycemia and dyslypidemia by mechanisms that are uncommonly assessed. In C-A ( Figure 2; K = 11), the median levels of the evaluated mechanisms fell within the 25th and 75th percentiles of the overall population, except for insulin secretion. Nevertheless, C-A harbored a cluster (Cluster 4) representing 3% of the cohort with normal insulin secretion but with hyperinsulinemia, which is associated with insulin clearance suppression and fatty liver. Noteworthy, in this cluster, median glycemic levels were increased at 30 minutes ( Figure 3; Supplemental Figures S4 and S5 for males and females, respectively). In C-II we focused in three sub-clusters of hierarchical discrimination for K = 27 (Clusters 18, 22, and 24) that covered distinct glycemia/FFA profiles and considering their representativeness in the population (Figure 3; Supplemental Figures S4 and S5 for males and females, respectively). In Cluster 18 all subjects had normoglycemia and the median glycemia and FFA values fell within the middle range of the overall population. However, they were exposed to hyperinsulinemia, plausibly explained by decreased insulin clearance as well as increased beta-cell function, in agreement with increased DI. This profile may represent a dysmetabolic compensation of the observed insulin resistance, namely in liver, that associates with NAFLD.
Clusters 22 and 24 are similarly associated with NAFLD and decreased eGFR (Figure 3; Supplemental Figures S4 and S5 for males and females respectively). However, both clusters were exposed to hyperglycemia and hyperinsulinemia, nevertheless only C-24 was exposed to higher FFA levels while showing lower adipose tissue insulin resistance than C-22. Strikingly, higher adipose insulin resistance was uncoupled from increased FFA levels in C-22.

Discussion
In this work we analyzed a population-based cohort, comprising subjects with normoglycemia non-treated pre-diabetes, and previously undiagnosed diabetes, and we identified a substructure of two optimal clusters (C-I and C-II). C-II was characterized by high insulin resistance, hyperinsulinemia, and fatty liver, recapitulating an expected profile of type 2 diabetes (insulin In C-II we focused in three sub-clusters of hierarchical discrimination for K = 27 (Clusters 18, 22, and 24) that covered distinct glycemia/FFA profiles and considering their representativeness in the population (Figure 3; Supplemental Figures S4 and S5 for males and females, respectively). In Cluster 18 all subjects had normoglycemia and the median glycemia and FFA values fell within the middle range of the overall population. However, they were exposed to hyperinsulinemia, plausibly explained by decreased insulin clearance as well as increased beta-cell function, in agreement with increased DI. This profile may represent a dysmetabolic compensation of the observed insulin resistance, namely in liver, that associates with NAFLD.
Clusters 22 and 24 are similarly associated with NAFLD and decreased eGFR (Figure 3; Supplemental Figures S4 and S5 for males and females respectively). However, both clusters were exposed to hyperglycemia and hyperinsulinemia, nevertheless only C-24 was exposed to higher FFA levels while showing lower adipose tissue insulin resistance than C-22. Strikingly, higher adipose insulin resistance was uncoupled from increased FFA levels in C-22.

Discussion
In this work we analyzed a population-based cohort, comprising subjects with normoglycemia non-treated pre-diabetes, and previously undiagnosed diabetes, and we identified a substructure of two optimal clusters (C-I and C-II). C-II was characterized by high insulin resistance, hyperinsulinemia, and fatty liver, recapitulating an expected profile of type 2 diabetes (insulin resistance) that, unexpectedly, included already a significant proportion of normoglycemic subjects. Within C-I, we identified the other well-known subgroups of insulin deficiency, that likewise, included individuals with normoand dysglycemia. Overall, this confirms that insulin resistance and insulin deficiency profiles occur across all currently accepted classes of dysglycemia and also in normoglycemic subjects.
Our analysis further revealed that the two distinct clusters are made up of heterogenous groups regarding both, mechanisms and the milieu, that importantly have been independently associated with the development of diabetes' complications [34,35]. In fact, the optimal number of clusters is defined based on subjects' distances in a hyperspace, that are not necessarily clinically significant. However, cluster analysis is also a useful exploratory tool and, therefore, we proceed the analysis to higher levels of granularity.
We found that subjects with similar glycemia disturbances fall in clusters that differ not only in respect to insulin secretion and overall insulin resistance, but also regarding tissue-specific insulin resistance and insulin clearance, as illustrated in Clusters 6 vs. 22 and 8 vs. 24 ( Figure 4). Conversely, within many clusters we found subjects that, despite being currently classified as having different glycemic status including normoglycemia, show a similar pattern of pathophysiological mechanisms. Whether such subjects with normoglycemia are protected by some non-evaluated mechanisms or traits, or on contrary, they are already progressing towards diabetes complications, cannot be discerned at this point. Nevertheless, these results indicate that the standard clinical evaluation of fasting and post-prandial glycemia does not reflect the overall exposure to glucose or express the metabolic milieu heterogeneity and that it is not informative about mechanisms governing the metabolic status. This may partly explain why some people classified with prediabetes might develop complications at early stages and others don't [12,36]. Likewise, we found clusters that, despite having the same fatty liver prevalence, can have or not dysglycemia ( Clusters 18,22,and 24). This might indicate diverse fatty liver etiology [37]. Curiously, the clusters with dysglycemia also had decreased glomerular filtration rate, supporting a role of the kidney in disease progression.
J. Clin. Med. 2020, 9, x FOR PEER REVIEW 10 of 14 resistance) that, unexpectedly, included already a significant proportion of normoglycemic subjects. Within C-I, we identified the other well-known subgroups of insulin deficiency, that likewise, included individuals with normo-and glycemia. Overall, this confirms that insulin resistance and insulin deficiency profiles occur across all currently accepted classes of glycemia and also in normoglycemic subjects.
Our analysis further revealed that the two distinct clusters are made up of heterogenous groups regarding both, mechanisms and the milieu, that importantly have been independently associated with the development of diabetes' complications [34,35]. In fact, the optimal number of clusters is defined based on subjects' distances in a hyperspace, that are not necessarily clinically significant. However, cluster analysis is also a useful exploratory tool and, therefore, we proceed the analysis to higher levels of granularity.
We found that subjects with similar glycemia disturbances fall in clusters that differ not only in respect to insulin secretion and overall insulin resistance, but also regarding tissue-specific insulin resistance and insulin clearance, as illustrated in Clusters 6 vs. 22 and 8 vs. 24 ( Figure 4). Conversely, within many clusters we found subjects that, despite being currently classified as having different glycemic status including normoglycemia, show a similar pattern of pathophysiological mechanisms. Whether such subjects with normoglycemia are protected by some non-evaluated mechanisms or traits, or on contrary, they are already progressing towards diabetes complications, cannot be discerned at this point. Nevertheless, these results indicate that the standard clinical evaluation of fasting and post-prandial glycemia does not reflect the overall exposure to glucose or express the metabolic milieu heterogeneity and that it is not informative about mechanisms governing the metabolic status. This may partly explain why some people classified with prediabetes might develop complications at early stages and others don't [12,36]. Likewise, we found clusters that, despite having the same fatty liver prevalence, can have or not glycemia ( Clusters 18,22,and 24). This might indicate diverse fatty liver etiology [37]. Curiously, the clusters with glycemia also had decreased glomerular filtration rate, supporting a role of the kidney in disease progression.  Clusters 4,6,8,18,22, and 24 is shown in radar plots. The considered percentile ranges to profile the clusters (10th, 25th, 75th, and 90th) are represented in doted circumferences (light grey). Represented parameters are grouped in milieu (yellow), mechanisms (blue), and pathologic conditions (red). In the footprint higher deviations from the center in the milieu correspond to higher exposure, in mechanisms to more severe affection and to the presence and severity of pathological conditions. IQR-Interquartile range.  Clusters 4,6,8,18,22, and 24 is shown in radar plots. The considered percentile ranges to profile the clusters (10th, 25th, 75th, and 90th) are represented in doted circumferences (light grey). Represented parameters are grouped in milieu (yellow), mechanisms (blue), and pathologic conditions (red). In the footprint higher deviations from the center in the milieu correspond to higher exposure, in mechanisms to more severe affection and to the presence and severity of pathological conditions. IQR-Interquartile range.
A previous cluster analysis of subjects with diabetes revealed five clusters, four of which were related to type 2 diabetes [6]. Beyond severe insulin resistance and insulin deficiency, Alqhvist et al. [6] identified one cluster related with age at onset of diabetes and another with BMI. We found a cluster (C-8), in which dysglycemia was not associated with insulin resistance or insulin deficiency, but also could not be explained by age or BMI. Therefore, this cluster may represent a dysglycemic trajectory that operates by mechanism yet to be uncovered.
In our study we focused on clusters with higher prevalence. However, we found that less prevalent clusters comprised mainly subjects with dysglycemia, despite the lower proportion of prediabetes and diabetes in our sample. The occurrence of these metabolic profiles, if validated in larger studies, may deserve attention in the context of the high prevalence of dysglycemia in the world. A limitation of our study is that we could only use cross-sectional data, which does not allow addressing the existence of dynamic paths across the clusters throughout time that could suggest effective preventive interventions. This should be assessed in future longitudinal studies. Furthermore, the high number of included parameters can be also considered a limitation. Analysis with fewer parameters might be more perceptible and thus more appealing. Nonetheless, our work suggests that such approach is not sufficient to reveal the metabolic complexity.
McCarthy proposed that different phenotypes of type 2 diabetes result from the combination of several mechanisms and traits [11]. Pathophysiological mechanisms that drive dysglycemia are affected in a continuous manner and can be differently combined. Of note, this model places the subject in a path from normalcy to disease. Our work, to the best of our knowledge, is the first to stratify a population including subjects with type 2 diabetes, prediabetes, and normoglycemia, not relying in the current diagnostic criteria. We found that metabolic heterogeneity can be captured by taking in account not only glycemia levels but also other components of the milieu, disclosing underlying pathophysiological mechanisms.
We could not address whether the founded clusters are associated with specific type 2 diabetes complications. We speculate that it might be the case, given that they are highly stratified in both milieu and mechanisms and we found clusters with higher fatty liver and others with lower estimated glomerular filtration.
Radar plots in Figure 4 illustrate the metabolic footprint of the clusters discussed above, where we can easily assess each milieu component (yellow), mechanisms affection (blue), as well as presence of pathological conditions (red). This representation supports the notion that glycemic imbalances is only one axis of a broader multidimensional profile. The most recent type 2 diabetes treatment guidelines involve a multifactorial treatment targeted to glycemic and other risk factors control, aiming to prevent chronic complications. The pharmacological armamentarium in diabetes includes new classes of drugs with different mechanisms of action in various organs. Still, many of the choices are based on trial and error methods or on the assessment of cardiovascular complications or kidney disease. A deeper understanding of the underlying mechanisms, taking in account the milieu beyond glycemia, and therefore better forecasting complications, will led to a more efficient therapeutic plan in the concept of precision care.

Conclusions
The metabolic footprint revealed alterations in the metabolic milieu and mechanisms impairments independently associated with diabetes and its complications, detected throughout all classes of dysglycemia and in normoglycemia. This work is a proof of concept and future developments further integrating omics data will aid discovery of additional pathophysiological mechanisms in specific clusters. It is expected that deeper phenotyping and increased pathophysiology knowledge will allow to more precisely identify subjects affected in this multidimensional profile, predict their progression, and treat them accordingly towards precision medicine. We envisage that the metabolic footprint will assist a precision approach to type 2 diabetes Supplementary Materials: The following are available online at http://www.mdpi.com/2077-0383/9/8/2588/s1. Figure S1: Self Organizing Map (SOM) neighbor distance plot (a) and quality plot (b) of the 27 SOM units; Figure  S2: Clusters profiling for the male population; Figure S3: Clusters profiling for the female population; Figure S4: Clusters profiling for K = 27, males; Figure S5: Clusters profiling for K = 27, females; Table S1: Indexes formulas and calculations of the parameters used to profile the clusters.