We evaluated the presence/absence of proteins encoded by 14 077 genes in adipocytes obtained from different tissue samples using immunohistochemistry. By combining this with previously published adipocyte‐specific proteome data, we identified proteins associated with 7340 genes in human adipocytes. This information was used to reconstruct a comprehensive and functional genome‐scale metabolic model of adipocyte metabolism. The resulting metabolic model, iAdipocytes1809, enables mechanistic insights into adipocyte metabolism on a genome‐wide level, and can serve as a scaffold for integration of omics data to understand the genotype–phenotype relationship in obese subjects. By integrating human transcriptome and fluxome data, we found an increase in the metabolic activity around androsterone, ganglioside GM2 and degradation products of heparan sulfate and keratan sulfate, and a decrease in mitochondrial metabolic activities in obese subjects compared with lean subjects. Our study hereby shows a path to identify new therapeutic targets for treating obesity through combination of high throughput patient data and metabolic modeling.
Combining large‐scale immunohistochemical analysis and proteomics data, 7340 gene products are identified in human adipocytes. Based on this data, a genome‐scale metabolic model is reconstructed and used to integrate clinical and transcriptome data from lean and obese subjects.
We simulated the metabolic differences between the individuals with different body mass indexes (BMIs) using transcriptome and fluxome data.
An increase in the metabolic activity around androsterone, ganglioside GM2 and degradation products of heparan sulfate and keratan sulfate, and a decrease in mitochondrial metabolic activities are found in obese subjects compared with lean subjects.
We simulated the change in lipid droplet (LD) size and found that lean subjects have large dynamic changes in LD formation compared with obese subjects.
Besides enabling patient stratification, our study allows the identification of novel therapeutic targets for obesity.
White adipose tissue (WAT) serves as an important buffer for handling daily lipid flux by uptake of fatty acids (FAs) in the post‐prandial state and release in the post‐absorptive state (Frayn, 2002). Furthermore, studies in recent years have revealed that WAT is a major endocrine organ sending out hormones and signaling molecules that regulate and coordinate energy homeostasis, insulin sensitivity, lipid metabolism, substrate selection, satiety and appetite (Cristancho and Lazar, 2011). WAT also cooperates with other tissues including liver, muscle, pancreas, heart and brain through the release of FAs (lipokines), secretory factors (adipokines) and pro‐inflammatory cytokines. WAT dysfunction or overload of its lipid storage capacity can lead to wide range of diseases (e.g., immunological and inflammatory diseases), including metabolic diseases such as obesity and its adverse outcomes (Lago et al, 2007; Auffray et al, 2009). Obesity can lead to preventable cause of death and increases the likelihood of coronary heart disease, diabetes and several forms of cancer (Cao, 2010) and ∼33% of the adults older than 20 years living in the United States and 20% of European adult population are obese. Obesity is therefore considered to be one of the greatest threats to global human health (Caveney et al, 2011). However, obesity is not directly related with particular etiological factors and it is therefore essential to enable stratification and prediction of disease risks among obese subjects to ensure early treatment either through diet intervention, exercise or surgery.
An increased understanding of the mechanisms behind obesity and related diseases will provide valuable insights into their etiology, pathogenesis and may lead to new treatment strategies. It is inherently difficult to find the exact disease onset, since such systemic diseases are caused by a combination of different genetic and environmental factors and often result in similar disease phenotypes. Therefore, it remains challenging to identify the cellular and molecular mechanisms associated with obesity, in particular as the metabolism of the single cell involves thousands of metabolites and interconnected chemical reactions that occur simultaneously. Proper understanding of such a complex system requires a holistic approach. For this purpose, the so‐called genome‐scale metabolic models (GEMs) are suitable as they allow for analysis of metabolism at the genome scale but at the same time ensures identification of specific pathways, enzymes and metabolites associated with specific phenotypes (Nielsen, 2009; Thiele and Palsson, 2010).
GEMs integrate biochemical and physiological data on genes and enzymatic reactions and allow the study of relationships between networks, functions, diseases and patients at a systems level. Reconstructed large‐scale GEMs can serve as scaffolds for data analysis and identify biologically meaningful correlations and mechanistic relationships between components of different sub‐systems (Patil and Nielsen, 2005). Hereby, GEMs can be employed for understanding the underlying mechanisms of complex diseases, and hence be used for identification of novel therapeutic and drug targets and discovery of new biomarkers. The diseases that we are discussing are not difficult to diagnose (e.g., obesity and diabetes) but efficient biomarkers would be of interest to predict prognosis and outcome of the disease and thereby enabled patient stratification which will be the basis for developing personalized medicine (Mardinoglu and Nielsen, 2012; Nielsen, 2012; Figure 1).
Two generic literature‐based GEMs of human metabolism, Recon 1 (Duarte et al, 2007) and the Edinburgh human metabolic network (EHMN; Ma et al, 2007) have been reconstructed previously. Although these first reconstructions represent a major advancement, human metabolism is specialized in different cell types, and hence there is a need for reconstruction of cell type or tissue‐specific GEMs. In this context, tissue‐specific GEMs have been developed for liver (hepatocytes) (Gille et al, 2010; Jerby et al, 2010), cardiomyocytes (Karlstaedt et al, 2012), kidney (Chang et al, 2010), brain (Lewis et al, 2010) and alveolar macrophage (Bordbar et al, 2010) and recently, three small cell type‐specific GEMs of hepatocytes, myocytes and adipocytes (iAB586) were developed (Bordbar et al, 2011). Furthermore, using the Integrative Network Inference for Tissues (INIT) algorithm, we previously reconstructed cell type‐specific draft GEMs for 69 different cell types and 16 cancer types (Agren et al, 2012).
With the objective of gaining new insight into adipocyte metabolism at the genome level, we first used human antibodies to evaluate the presence/absence of 14 077 proteins in adipocytes found in breast and two different soft tissue samples and checked their presence call with previously published adipocyte proteome data. Second, we manually reconstructed a high‐quality, simulation ready GEM for adipocytes by using all adipocyte‐specific proteome data. The model is based on previously published GEMs but also on publicly available databases on metabolism. Third, we employed the functional GEM for the analysis of microarrays that profile the gene expression from subcutaneous adipose tissue (SAT) of subjects from the Swedish Obese Subjects (SOS) Sib Pair Study, which includes nuclear families with body mass index (BMI)–discordant sibling pairs (BMI difference⩾10 kg/m2). The male and female participants of the study were divided into three different groups based on their BMIs: lean, overweight and obese. Analysis of transcription data from this study recently demonstrated that there are differences in mitochondrial function between men and women (Nookaew et al, 2012), but there was not performed any analysis on the effect of obesity. We incorporated differentially expressed genes between obese and lean subjects into the GEM. Besides the gene expression data from the SOS Sib Pair Study, additional clinical data (e.g., plasma and WAT lipid concentrations) were also incorporated into the model. By integrating gene expression data and adipose tissue uptake/secretion rates with the reconstructed GEM, we identified metabolic differences between individuals with different BMIs by using the concept of Reporter Metabolites (Patil and Nielsen, 2005) and transcriptionally controlled reaction fluxes (Bordel et al, 2010; Figure 1).
Immunohistochemistry‐based proteomics of human adipocytes
Several studies have reported proteomics data of human WAT that include not only adipocytes but also the connective tissue matrix, nerve tissue, stromal vascular cells and immune cells (Peinado et al, 2012). Recently, Xie et al (2010) characterized the proteome of human adipocytes and reported existence of proteins encoded by 1574 genes in subcutaneous abdominal adipocytes taken from three healthy lean subjects. Although this first adipocyte‐specific proteome study is promising to understand the unique characteristics of the adipocytes, it is still necessary to expand the protein coverage to be able to study adipocyte biology at the genome‐wide level.
The Human Protein Atlas (HPA) is a knowledge‐based portal that cover the annotated protein expression feature for protein targets analyzed with one or more antibodies and the main subcellular localization of protein targets in all major human cell types and cell lines (Uhlen et al, 2010). Here, we expanded the coverage of the HPA so as to define the protein profiles of adipocytes found in breast and two different soft tissues and examined the spatial distribution and the relative abundance of proteins encoded by 14 077 genes in adipocytes (Supplementary Dataset 1a). Soft tissue cores with 1 mm diameter were sampled from different locations of the human body and included in tissue microarrays (TMAs) together with breast tissue core as previously described (Kampf et al, 2012). A total of 17 296 high‐throughput generated affinity‐purified antibodies (Supplementary Dataset 1a) against 14 077 proteins were arranged according to abundance of the corresponding protein target. Immunohistochemically stained sections from TMAs blocks were scanned in high‐resolution scanners and separated to individual spot images to represent each core. Here, proteins with strong, moderate and weak relative expression were included in the reconstruction process of the GEM for adipocytes. It is important to point out that the absolute levels of each protein have not been determined and could vary by orders of magnitude. The high‐resolution images together with annotation of the presence or absence of a particular protein target in adipocytes are publically available through the HPA (http://www.proteinatlas.org).
The proteome data generated for adipocytes found in breast and two different soft tissues were merged with previously published proteome data as presented in Figure 2A (Supplementary Dataset 1b). In total, we have proteome evidence for the presence/absence of proteins associated with 14 337 genes in adipocytes and together with the genes in our Human Metabolic Reaction (HMR) database (see later), we cover 98% of the enzymes (Supplementary Dataset 1c) reported in the HPA database (Figure 2B). Our proteomics analysis resulted in provision of evidence for presence of proteins associated with 7340 genes in adipocytes that could be used for our GEM reconstruction. Each protein that has evidence for presence/absence in adipocytes was annotated with UniProt id (Apweiler et al, 2011) and the corresponding encoding genes were annotated with Ensembl (Flicek et al, 2011) id for standardization.
Since adipocytes obtained from breast and two different soft tissues were used for the protein profiling, there was some variation between the samples. To estimate the effect of these variations on the functionality of the adipocytes, we used the functional annotation tool DAVID to calculate the enrichment in KEGG pathways (Huang et al, 2009). The results are presented in Figure 2C for breast, the two types of soft tissues, as well as for all proteomics data used for the GEM reconstruction. The analysis demonstrated that for all the adipocyte‐specific proteome data, there is enrichment in terms of metabolic pathways including FA metabolism, elongation and biosynthesis as well as major signaling pathways including adipocytokine, insulin, neurotrophin and PPAR signaling pathways. The figure shows that the samples from different tissues exhibit some differences, but that the overall pattern is similar. We therefore decided to reconstruct a general GEM for adipocytes by incorporating all proteins that were expressed in any of the three tissues.
Subcellular localization of the proteins and associated reactions
Metabolism is distributed across subcellular compartments of the cell and precise coordination is important for overall cellular function and maintenance of metabolic homeostasis. For instance, β‐oxidation of very long chain FAs (VLCFAs) is a process that starts in the peroxisome and ends in the mitochondria. The subcellular localization of reactions has large implications on the functionality of GEMs, as only a portion of metabolites can be transported between compartments. Furthermore, compartments can be individually redox and/or energy balanced. It is therefore important to understand the locations of individual enzymes as well as their associated reactions when metabolism is going to be studied at the genome scale. The HPA includes subcellular profiling data using immunofluorescence‐based confocal microscopy in three human cancer cell lines of different origin to be used for further in‐depth functional studies.
Here, proteins were classified into eight different compartments following our HMR database standard (Agren et al, 2012): cytosol, nucleus, endoplasmic reticulum (ER), Golgi apparatus (GA), peroxisome, lysosome, mitochondria and extracellular space. Reactions were assigned to compartments through their association with proteins in these different compartments. We assigned a confidence score ranging from one to three for each protein due the availability of knowledge in HPA and Uniprot (Supplementary Dataset 2a). Proteins localized in the aggresome, centrosome and cytoskeletons were assigned to the cytoplasm, whereas proteins in cell junctions and focal adhesions were assigned to the extracellular compartment (Supplementary Dataset 2b). Lysosome and peroxisome were merged into vesicles in the HPA data, but proteins associated with vesicles in our model were separated based on Uniprot data. In the case a protein has different locations in the two resources, we chose to assign the protein in the location reported based on the HPA data. In cases where the subcellular localization information was absent in HPA, Uniprot and literature (indirect physiological evidence), proteins and their associated reactions were assigned to be present in the cytosol.
Transcriptome data for SAT and its correlation with proteome data
Gene expression in SAT of 304 subjects involved in the SOS Sib Pair Study was analyzed by using Affymetrix U133 Plus 2.0 microarrays (Affymetrix, Santa Clara, CA, USA). The participants of the study, 209 female and 95 male subjects, were divided into three different groups according to their BMI: lean (18.5<BMI <25), overweight (25⩽ BMI<30) and obese (30⩽BMI) (Supplementary Table S1). Our study considers the differences in male and female obese subjects compared with lean male and female subjects independently. Differentially expressed genes between obese, overweight and lean subjects for male and female subjects were identified and differentially expressed genes on male and female obese subjects were incorporated during the reconstruction process of GEM for adipocytes to analyze gene expression data from the SOS Sib Pair Study.
Although the generated proteome data for adipocytes cover the entire set of cellular processes (e.g., signaling, metabolism and cell cycle), GEMs are applicable only for the study of metabolism. To get a general overview of the global changes between obese, overweight and lean subjects, the enrichment of differentially expressed genes was calculated for KEGG pathways (Supplementary Figures S1 and S2) and for biological process Gene Ontology (BP:GO) terms (Supplementary Figures S3 and S4). This was done using DAVID and for male and female subjects (Huang et al, 2009). To check the correlation of the genome‐wide transcription data of SAT with the proteome data, the enrichment of differentially expressed genes in male and female obese subjects was calculated for the most significant KEGG pathways from the analysis of the proteome data (Figure 2C). This was done as the transcriptome data for SAT represent not only adipocytes but also other cell types; including immune cells and preadipocytes linked with different BMIs.
The analysis demonstrates that some of the metabolic and signaling pathways found to be enriched in adipocytes based on the proteome data also show significant changes in gene expression between lean and obese subjects, both in males and in females. Similarly, we find that BP:GO terms that are enriched based on the proteome data (assessed with DAVID; Huang et al, 2009) also show enrichment based on the transcriptome data (Supplementary Dataset 3). Thus, enriched BP:GO terms such as post‐translational protein modification, cellular protein metabolic process, lipid metabolic process, cellular lipid metabolic process and FA metabolic process are found both from the adipocyte‐specific proteome data and from comparison of expression data for lean and obese subjects.
Reconstruction of the adipocytes GEM iAdipocytes1809
GEMs are reconstructed based on high‐throughput data such as genome, transcriptome, proteome, metabolome and fluxome and they provide an excellent scaffold for integrative analysis of this kind of data (Patil and Nielsen, 2005; Cakir et al, 2006). To reconstruct a large‐scale comprehensive GEM for adipocytes, biochemical and genetic evidence were combined with data on protein expression and localization. Besides the proteins associated with 7340 genes identified from our proteomics analysis, we further used differentially expressed genes in SAT of obese male and female subjects for our reconstruction process. However, SAT contains not only adipocytes but also other cell types and it is therefore necessary to check the functionality of the genes before including into the GEM.
HepatoNet1, a GEM for hepatocytes which is reconstructed based on the manual evaluation of the original scientific literature (Gille et al, 2010), was used as a starting point for our reconstruction process and used to generate an initial candidate list of network components (Figure 3A). First, metabolism of lipids and lipoproteins in Reactome, a manually curated and peer‐reviewed pathway database (Croft et al, 2011), was merged into HepatoNet1. Second, the resulting network was combined with the evidence‐based generic human models Recon1 (Duarte et al, 2007) and the compartmentalized EHMN (Hao et al, 2010). This combined reaction list resulted in an updated version of our HMR database that is available at http://www.metabolicatlas.org. The HMR database contains 6000 metabolites in 8 different compartments (3160 unique metabolites), 8100 reactions and 3668 genes associated with those reactions. HMR includes all of the genes and gene‐associated reactions in HepatoNet1, Recon1 and EHMN and other isolated gene‐associated reactions, and represents the most comprehensive human reaction database for genome‐scale modeling (Figure 3B). Third, the existence of each protein coding genes associated with a reaction in HMR was assessed for the presence or absence in adipocytes using previously published and the here generated adipocyte‐specific proteome data. Differentially expressed genes and associated reactions in obese subjects compared with lean subjects were also included in the model by checking their functionality. This process provided us with a list of reactions that occur in adipocytes. Gaps in the resulting network were filled using our recently developed iHuman1512 metabolic network (Agren et al, 2012), public databases such as KEGG (Kanehisa et al, 2010) and HumanCyc (Romero et al, 2005) and manual evaluation of the literature about adipocyte metabolism. This gap filling resulted in generation of iAdipocytes1809 that is a functional and fully connected GEM for adipocytes (Figure 3A). Reactions were included in the GEM depending on evidence from previously published models and databases (Supplementary Table S2) or on the availability of specific experimental evidence, for example, enzyme assay or protein identification, for the occurrence of the corresponding reaction in adipocytes. Reaction directionality in the model was treated as in the original resource and the directionality was only changed if it was necessary to perform successful simulations of known biological functions of adipocytes. iAdipocytes1809 contains 6160 reactions and 4550 metabolites in 8 different compartments (2497 unique metabolites) (Supplementary Table S3) and it is available in the Systems Biology Mark‐up Language (SBML) format at http://www.metabolicatlas.org. There are 1809 genes in the model and 80% of the reactions are associated with one or more genes. In iAdipocytes1809, individual metabolites rather than generic pool metabolites for 59 FAs (Supplementary Table S4) have been used and this allowed us to incorporate measured concentrations of different FAs in human plasma and adipocytes into the model.
Adipocyte‐specific metabolome data from the Human Metabolome Database (HMDB) (Wishart et al, 2009) and a comprehensive database for lipid biology, Lipidomics Gateway (Harkewicz and Dennis, 2010), were used in the reconstruction process. To ensure standardization, each metabolite in the GEM was assigned at least one of the following: HMDB ids, Lipidomics Gateway ids, KEGG ids (Kanehisa et al, 2010), Chemical Entities of Biological Interest (ChEBI) (Degtyarenko et al, 2008), or International Chemical Identifiers (InChI). In the model, new genes were assigned to the reactions in iAdipocytes1809 using EC numbers from UniProt (Apweiler et al, 2011) and the Lipid Map proteome database (Cotter et al, 2006).
The reconstruction of iAdipocytes1809 involved a comprehensive review of lipid metabolism and gene products related to the lipid metabolism in publicly available databases were closely examined for their existence in adipocytes. In all, 1235 protein encoding genes with proteome evidence in adipocytes and 244 genes in iAB586 were included in iAdipocytes1809 (Figure 3C). Another 137 genes were included in the model due to positive evidence for the presence of these proteins from the transcriptome data. These 137 genes were significantly higher expressed both in male and in female obese subjects compared with lean subjects, and are hence relevant for the analysis of gene expression data from the SOS Sib Pair Study. Another 193 genes and their associated reactions were included in the model due to connectivity constraints and known function of adipocytes and these genes are mainly associated with transport across membranes. In all, 73 of these genes do not have any negative evidence whereas 120 of these genes have negative HPA scores. These 120 genes are likely to be false negatives in the HPA, and they can be explained by either poor antibody hybridization to adipocytes or condition‐dependent protein expression. We also compared iAdipocytes1809 with iAB586 and the previously published generic human network, iHuman1512, and we found that iAdipocytes1809 contains all of the genes and associated reactions in iAB586 and adipocyte‐specific reactions and associated genes of iHuman1512 (Figure 3D).
In iAdipocytes1809, 59 different common long and very long chain FAs (Supplementary Table S4) in human plasma can be taken up as NEFAs and lipoproteins (Supplementary Figures S5–S7). Cholesterol and its 59 different CEs can also be taken up from low‐density lipoproteins (LDLs) and high‐density lipoproteins (HDLs). In the post‐prandial state, FAs and CEs can be incorporated into LD structures and in the post‐absorptive state LDs can be broken down to FAs and CEs. The comprehensive iAdipocytes1809 covers all metabolic pathways known to exist in adipocytes as well as extensive knowledge of lipid metabolism in human cells and adipocytes in particular. Lipid metabolism in adipocytes involve uptake of FAs from two potential sources: non‐esterified FAs (NEFAs) and lipoproteins (Supplementary Table S5) including chylomicrons, very low‐density lipoprotein (VLDL), LDL and HDL through LPL. iAdipocytes1809 also cover the transport of FAs, cholesterol and cholesterol esters (CEs) across the plasma membrane, small amount of de novo FA and cholesterol synthesis, FA transport into mitochondria and peroxisomes for β‐oxidation, FA esterification into triacylglycerols (TAGs), lipolysis and lipid droplet (LD) formation and maintenance (see Supplementary information).
Model validation and formation of lipid droplets
Due to the large size of GEMs (typically thousands of metabolic reactions), proper validation is essential to ensure that a reconstructed model has predictive ability (Mardinoglu and Nielsen, 2012). A common problem is lack of mass balancing of the individual reactions, and this typically involves metabolites with ill‐defined formulas, such as polymers or metabolite pools, but different protonation states can also be problematic. This type of problem results in a situation where the model can take up or excrete metabolites in an unbalanced manner. The model was tested so that all reactions except pool reactions were mass balanced. A second issue is that reactions can violate thermodynamic constraints by being written in the wrong direction or where irreversible reactions are written as reversible. This was tested by making sure that the model could not generate high‐energy compounds from low‐energy compounds (such as ATP and water from ADP and phosphate or any organic compound from CO2 and water). A third issue is that reactions may be dead‐ends, meaning that they cannot carry flux, either because the model is not able to produce/take up some of the substrates or because it is not able to consume/excrete some of the products. Dead‐end reactions tend to propagate as one unconnected reaction can give rise to many more. Therefore, many GEMs contain a large proportion of dead‐end reactions. Considerable efforts were spent in making sure that the number of dead‐end reactions was kept to a minimum. Furthermore, with the simulation ready iAdipocytes1809, the production of all metabolites in the model was checked with minimum input to the model (Supplementary Dataset 4) to ensure the connectivity. Lastly, even a well‐connected, thermodynamically correct and balanced model may not be able to perform all relevant metabolic functions, or it may be able to perform functions that it should not do (such as synthesis of essential amino acids or FAs). The model was therefore validated for 250 known metabolic functions of adipocytes, adapted from the definitions provided in connection with setting up HepatoNet1 (Gille et al, 2010; Supplementary Dataset 5). All validation steps and simulations were performed using the RAVEN Toolbox (Agren et al, 2013).
Formation of LDs (Figure 4) is included in iAdipocytes1809. LDs protect the cell from the lipotoxic effects of unesterified lipids, but they are also associated with the development of metabolic diseases, such as obesity, type 2 diabetes, atherosclerosis and liver steatosis that are characterized by the extensive accumulation of LDs. LDs represent a significant energy storage that can be mobilized by catabolism and this process is highly regulated by hormones and signaling pathways. Although LDs are recognized as an organelle called adiposomes (Farese and Walther, 2009), we represented it as a large composite metabolite in our model. LDs have polar surfaces that contain a variety of proteins, PLs and sterols and have a non‐polar core that contain CEs, DAGs and TAGs (Gross et al, 2011; Supplementary Table S6). The diameter of LDs varies from 50 nm to 200 μm and LDs may grow to enormous sizes with different mechanisms including localized synthesis of lipids, transport of lipids and coalescence of LDs.
The function of iAdipocytes1809 was examined by estimating the formation of LDs and by means of available physiological, biochemical and genetic evidence in lean and obese subjects and clinically observed data. Recently, McQuaid et al (2011) measured the delivery and transport of FAs in adipose tissue using multiple and simultaneous stable‐isotope FA tracers in lean and obese subject groups over 24‐h period. Even though abdominally obese subjects have greater adipose tissue mass than control lean subjects, the rates of delivery of NEFAs were downregulated in obese subjects. Clinical data on Supplementary Table S1, NEFA release, TAG extraction, glucose uptake (Supplementary Dataset 6; McQuaid et al, 2011) and amino‐acid uptake rates (Supplementary Dataset 7; Patterson et al, 2002) and FA concentrations in plasma and adipocytes (Supplementary Dataset 8) were incorporated into the model to predict the formation of LDs. Based on measurements of the uptake of glucose and TAG and the release of NEFAs over a 24‐h period (Figure 5A and 5B), we simulated the change in LD size (Figure 5C; Supplementary Dataset 9a). We found from our simulations that lean subjects have large dynamic changes in LD formation compared with obese subjects, which is in agreement with experimental data (Arner et al, 2011). Furthermore, we predicted a lower acetyl‐CoA production in obese subjects, as shown in Figure 5D (Supplementary Dataset 9b).
Identification of obesity‐specific metabolic features
Modeling using iAdipocytes1809 can be applied to predict metabolic states under various perturbations, study regulation of the adipocytes, identify potential therapeutic targets and discover novel biomarkers for the development of more effective therapies. Here, we used the model to identify Reporter Metabolites (Patil and Nielsen, 2005) of male and female obese subjects compared with lean subjects using gene expression data obtained from the SOS Sib Pair Study. Reporter Metabolites are metabolite nodes in the metabolic network around which there are significant transcriptional changes. Here, 20 statistically significant Reporter Metabolites are presented for upregulated and downregulated genes in male and female obese subjects through the employment of iAdipocytes1809 (Figure 6). The most significant results from our Reporter Metabolites analysis for upregulated and downregulated genes are correlated with the KEGG pathways enrichment results of significantly expressed genes in the obese subject groups (Supplementary Figures S1 and S2).
To illustrate the improvement of iAdipocytes1809 over the published iAB586, the Reporter Metabolites were also calculated for male and female obese subjects by using iAB586 (Supplementary Figures S8). Reporter Metabolites involved in the mitochondrial dysfunction as well as different amino acids were identified to be similar to the Reporter Metabolite analysis using iAdipocytes1809. However, iAB586 could not detect several of the most significant and in our view most interesting Reporter Metabolites identified when using iAdipocytes1809 due to the increase in number of reactions, metabolites and genes (see Discussion). Thus, the Reporter Metabolite analysis with iAdipocytes1809 and iAB586 provides an unbiased confirmation that iAdipocytes1809 represent significant advancement of the adipocyte metabolic network compared with iAB586.
Identification of obesity‐specific transcriptionally regulated reactions by random sampling
We identified changes in metabolic fluxes in response to obesity and which of these changes are likely to be associated with transcriptional changes using iAdipocytes1809. We defined a region of feasible flux distributions using uptake rates for TAGs and glucose and NEFA release rates for lean and obese subjects, and used these to calculate a set of possible flux distributions using a random sampling algorithm (Bordel et al, 2010). The average values and standard deviations for each of the fluxes in iAdipocytes1809 were calculated and the changing fluxes were compared with the significance of change in gene transcription for the corresponding enzymes. This allowed us to identify specific reactions for which flux changes are likely to be transcriptionally regulated.
The results from this analysis showed that the following pathway fluxes were transcriptionally downregulated in obese subjects: uptake of glucose, uptake of FAs, oxidative phosphorylation, mitochondrial and perixomal β‐oxidation, FA metabolism and tricarboxylic acid (TCA) cycle. This analysis is consistent with the findings from the Reporter Metabolites analysis and KEGG enrichment pathway analysis (Supplementary Table S7). Furthermore, fluxes associated with beta‐alanine metabolism were found to be transcriptionally downregulated in obese subjects (Figure 7). Previously, it has been reported that blood flow, glucose uptake, release of NEFA and the extraction of TAG from plasma were significantly lower in abdominally obese subjects compared with lean subjects (McQuaid et al, 2011). Our random sampling results clearly indicate that all metabolic pathways in mitochondria are downregulated (Kusminski and Scherer, 2012) similar to the mitochondrial decline that is a hallmark of different diseases associated with aging.
We evaluated the presence/absence of 14 077 proteins in adipocytes using human antibodies and presented a high‐quality, simulation ready and functional GEM for adipocytes, iAdipocytes1809, based on proteome, metabolome and transcriptome data. iAdipocytes1809 was used to analyze clinically observed transcriptome and fluxome data to understand the mechanistic changes in adipocyte metabolism in response to obesity. Gene expression and associated metabolites in obese and lean subjects were studied using a systems biology approach, which allowed us to understand how the transcriptome data in SAT represent different metabolic functions in obesity. An increased understanding of metabolic pathways altered in complex diseases may contribute to the identification of novel therapeutic targets (Cheong et al, 2012) and here potential therapeutic targets for obesity were discovered through Reporter Metabolite analysis and random sampling of flux states. Moreover, protein coding genes related to the therapeutic targets may be used as potential drug targets.
Through our analysis, it was observed that expression of genes involved in metabolic pathways including mitochondrial and perixomal β‐oxidation, FA synthesis, amino‐acid metabolism, pyruvate metabolism, oxidative phosphorylation and TCA cycle are downregulated in obese subjects. Most of these pathways are linked with the mitochondrial dysfunction and different therapeutic interventions including antioxidants and chemical uncoupler treatments have proven to improve the mitochondrial dysfunction (Kusminski and Scherer, 2012). Recently, Canto et al (2012) have reported that increasing the NAD+ levels, identified as therapeutic targets in our study, through supplements (NAD+ precursors) enhances the oxidative metabolism and protects the cells against high‐fat diet‐induced obesity.
Mitochondrial acetyl‐CoA has a central role in different pathways in the mitochondria and it reacts with oxaloacetate to form citrate, which can be transported from the mitochondria to the cytosol where it is participating in FA synthesis (Dean et al, 2009). Acetyl‐CoA derived through other principal sources, including degradation of amino‐acid and ketone bodies, and FA oxidation processes are insufficient for FA synthesis. Increasing the acetyl‐CoA concentration and eventually FA synthesis in adipose tissue of obese subjects results in whole body regulation of metabolism, including stimulation of muscle insulin action and suppression of hepatosteatosis, as reported by Cao et al (2008). Here, we propose to boost the metabolic activity of mitochondria in the adipocytes of obese subjects by increasing the mitochondrial acetyl‐CoA which is another therapeutic target identified in our study. Acetyl‐CoA formation in lean and obese subjects was predicted through use of iAdipocytes1809 in lean and obese subjects and its concentration can be increased through the uptake of beta‐alanine. The effect of beta‐alanine as a dietary supplement was previously examined in football players and it is reported that it has effect on lean tissue accruement and body fat composition (Hoffman et al, 2006). Furthermore, rat studies reported that beta‐alanine decreases the LPL enzyme activity in adipose tissue (Prabha et al, 1988) that may help to decrease the uptake of FAs to be stored in adipocytes. Our results suggest that increasing the level of beta‐alanine in obese subjects may help to decrease the fat composition in obese subjects.
Furthermore, we have identified new candidates for potential therapeutic targets for obesity and their association with the obesity has been reported in different studies. One of the most significant results from the Reporter Metabolite analysis for upregulated genes was androsterone and its precursor 5alpha‐androstane‐3,17‐dione. In obese subjects, increased secretion of androsterone was reported in serum concentrations using various indices calculated from urinary steroid excretion rates (Vierhapper et al, 2004). The measurement of androsterone‐to‐etiocholanolone ratio in urine is commonly used as an indicator of activities of 5α‐ and 5β‐reductases in both male and female subjects and the ratio increased with indexes of insulin resistance, which is an adverse outcome of obesity (Tomlinson et al, 2008). Our results indicated that metabolism of obese subjects around androsterone increases and the detection of androsterone level in plasma may represent a marker for altered metabolism in obese subjects. 5alpha‐pregnane‐3,20‐dione and 3alpha‐hydroxy‐5alpha‐pregnan‐20‐one obtained from Reporter Metabolites analysis are other metabolites involved in steroid hormone synthesis specifically in pathway of progesterone metabolism. Several studies reported that sex steroids affect the adipose tissue deposition through the regulation of preadipocyte proliferation and/or differentiation as well as lipogenesis and/or lipolysis of differentiated adipocytes (Anderson et al, 2001). Zhang et al (2009) previously studied the comparison of progesterone metabolite formation in preadipocytes and lipid‐storing adipocytes and they reported that the production of progesterone metabolites was significantly increased. Our study signified that metabolism around these two metabolites increase in obese subjects relative to metabolism in lean subjects.
Another high‐ranking target for upregulated genes is ganglioside GM2 (GM2). Gangliosides, one of the major glycosphingolipids in mammals, have major roles as mediators for cell to cell or cell to matrix recognition and regulate the transmembrane signal transducers and cell proliferation. Gangliosides in adipose tissues are also associated with insulin signaling mechanisms and it is reported that series of gangliosides GM2, GM1 and GD1a are dramatically increased in adipose tissues of obese mice (Tanabe et al, 2009). The action of GM2 synthase also directs the metabolic flow of a‐series gangliosides toward GD1a in mouse studies.
A third prominent group among the Reporter Metabolites for upregulated genes is the degradation products of heparan sulfate proteoglycans (HSPG) and keratan sulfate. These compounds are classified as glycosaminoglycans and attach to cell surface or extracellular matrix proteins. It has been reported that catalytically active adipose tissue LPL attaches to HSPG at the luminal surface of vascular endothelium (Olivecrona and Beisiegel, 1997; Lafontan, 2008) and hydrolyze the TAGs for uptake of FAs into the cell. The LPL moves between individual HSPG chains within the layer and this creates a high concentration of LPL along the surface layer of HSPG chains (Lookene et al, 1996). In the presence of heparin, more LPL is secreted and increased secretion was balanced by decreased degradation of LPL. There are special mechanisms that inhibit LPL and one mechanism is that LPL forms complexes with FAs (Bengtsson and Olivecrona, 1980). During the LPL hydrolysis and accumulation of FAs in the cells, the LPL is sequestered into enzyme FA complexes, lipolysis is reduced and eventually the binding of LPL to heparan sulfate is broken. If a high‐affinity ligand (e.g., FAs, heparin and apoCII) is available, then the LPL detaches from the cell surface to heparan sulfate chains and without ligand in the medium, the LPL recycles into the cells where it is degraded. Furthermore, several studies have reported that more sulfated polysaccharide chains increase the affinity for binding of LPL (Olivecrona and Olivecrona, 2009). Keratan sulfate, a biomarker of proteoglycan degradation, can be expressed from stem cells in human SAT and its relevance with obesity has been reported earlier. Messier et al (2000) studied the effect of exercise and diet in weight loss in older obese adults with knee osteoarthritis and analyze the levels of total proteoglycan, keratan sulfate and interleukin‐1 beta in their synovial fluid. It is reported that all participants of the study lost weight by changing their diet and exercise regimens and the level of keratan sulfate in synovial fluid decreased.
We envisage that iAdipocytes1809 is a key step to better enable links between molecular processes and patient phenotypes and hereby enable patient stratification through identification of specific molecular mechanisms in adipocyte metabolism. Here, we demonstrated this by predicting differences in the formation of LDs and acetyl‐CoA in lean and obese subjects. Besides enabling patient stratification, this may also lead to identification of novel therapeutic targets for obesity through combination of our model with high throughput patient data. iAdipocytes1809 can also be used as a scaffold for a comprehensive whole adipocyte model that accounts for all of the annotated gene functions identified in adipocytes (Karr et al, 2012). Compared with the previously described adipose model, the here presented model is significantly larger in terms of metabolites/genes and reactions and this allows for identification of metabolic biomarkers that cannot be identified with iAB586. Furthermore, our model has undergone a thorough validation process where 250 metabolic functions were simulated and as our model also included a description of individual FAs and sterolesters it can much better simulate lipid metabolism, including simulation of lipid droplet formation. In conclusion, we demonstrated the high quality adipocyte GEM iAdipocytes1809 is very well suited for integration of omics data and hereby result in a comprehensive understanding of adipocytes biology in response to obesity.
Materials and methods
Proteome data for adipocytes
The proteomic profiling of adipocytes in breast and two different soft tissues using immunohistochemistry was performed as previously described (Uhlen et al, 2005). Representative formalin paraffin‐embedded material from donor blocks was punched (1 mm in diameter) and placed in a recipient block TMA that includes 46 normal tissues, 20 types of cancer and 47 cell lines (Kampf et al, 2012). Thereafter, 4‐μm TMA sections were cut using a microtome and placed on super frost glass slides. Breast and two different soft tissue cores together with 708 previously prepared tissue cores were analyzed for protein expression using IHC (Kampf et al, 2004). The variability introduced by the individual experimental staining protocol, including the choice of antibody dilution and antigen retrieval methods, was addressed by the use of TMAs (Uhlen et al, 2005). Information from databases, for example, Uniprot, ENSEMBL, and published data from PubMed were used as guides to establish if immunohistochemical results represent expected protein expression of intended target protein. In addition, other validation strategies were used to ensure antibody validity, for example, protein arrays, western blot and immunofluorescence experiments. When applicable, paired antibodies raised against separate, non‐overlapping epitopes on the same target protein were used for extended validation (Uhlen et al, 2010). IHC‐stained tissues were scanned and digitalized at × 20 magnification. Annotation of high‐resolution images was manually performed by certified pathologists. Relative expression was indicated with four different color codes ranging from strong (red), moderate (orange), weak (yellow) and no expression (white) (Kampf et al, 2004). Localization information of proteins was inferred from manually curated Uniprot data (Apweiler et al, 2011) and recently generated HPA data on intracellular localization of proteins (Lundberg and Uhlen, 2010).
Transcriptome data for SAT
Gene expression in SAT of subjects involved in SOS Sib Pair Study which includes nuclear families with BMI–discordant sibling pairs (BMI difference ⩾10 kg/m2) was measured by microarray. The study group consisted of 304 subjects (209 female and 95 male) divided into three different groups according to their BMI: lean (18.5<BMI <25), overweight (25⩽ BMI<30) and obese (30⩽BMI) (Supplementary Table S1). Total RNA was prepared as previously described (Carlsson et al, 2009) and gene expression in human SAT was measured using Affymetrix U133 Plus 2.0 microarrays (Affymetrix). Probe hybridization intensity values were summarized using to Probe Logarithmic Intensity Error (PLIER) and the robust quantile method was applied to get normalized expression values. The normalization of the microarrays was carried out using the Expression Console software from Affymetrix and the quality assessment was carried out using R Statistical and Computing language and the Bioconductor software (Gentleman et al, 2004).
The annotation of the presence or absence of protein targets in adipocytes together with the high‐resolution images is publically available through the HPA (http://www.proteinatlas.org).
GEM for adipocytes, iAdipocytes1809 and HMR database is publically available in the Systems Biology Mark‐up Language (SBML) format at Human Metabolic Atlas (http://www.metabolicatlas.org).
Gene expression from SAT of subjects from the SOS Sib Pair Study is publically available in Gene Expression Omnibus (GEO) database with the accession number GSE27916 http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?token=phadhqigowykote&acc=GSE27916.
We acknowledge Prof Fredrik Karpe and Dr Barbara A. Fielding, Oxford Centre for Diabetes, Endocrinology and Metabolism for providing adipose tissue uptake rates for lean and obese subjects. This work was supported by grants from the Knut and Alice Wallenberg Foundation, the Chalmers Foundation, the Wellcome Trust (GR079534) and the Swedish Research Council (K2010‐55X‐11285‐13).
Author contributions: AM reconstructed iAdipocytes1809 and validated the model together with RA. AM and RA performed the analysis of clinical data. CK and AA generated the proteome data for adipocytes. PJ, AJW, PF and LMC generated gene expression data for adipose tissue and IN normalized the gene expression data. JN and MU conceived the project. AM, RA and JN wrote the paper and all authors were involved in editing the paper.
Conflict of Interest
The authors declare that they have no conflict of interest.
Supplementary Tables S1–7, Supplementary Figures S1–8
Supplementary Dataset 1a
The presence/absence of proteins encoded by 14 077 genes in adipocytes is examined based on immunohistochemistry and using antibodies generated within the Human Protein Atlas (HPA) project.
Supplementary Dataset 1b
Xie et al. (Xie et al., 2010) characterized the proteome of human adipocytes and reported existence of proteins encoded by 1 574 genes in subcutaneous abdominal adipocytes taken from three healthy lean subjects.
Supplementary Dataset 1c
The enzymes reported in the Human Protein Atlas database (http://www.proteinatlas.org) are listed.
Supplementary Dataset 2a
Subcellular localization of the proteins based on the Human Protein Atlas (HPA) and Uniprot data and their assignments in the genome‐scale metabolic model for adipocytes, iAdipocytes1809.
Supplementary Dataset 2b
Subcellular localization of the proteins and their confidence scores in the genome‐scale metabolic model for adipocytes, iAdipocytes1809 based on the Human Protein Atlas (HPA) and Uniprot data. Confidence score 2 for HPA, 1 for Uniprot and 3 for both HPA and Uniprot were assigned for each protein.
Supplementary Dataset 3
Enrichment of biological process Gene Ontology terms for encoded genes in here generated and previously published proteome data (assessed with DAVID (Huang et al., 2009)) are compared with the biological process Gene Ontology terms enrichment results of encoded genes in transcriptome data in male and female subjects.
Supplementary Dataset 4
List of the input metabolites and inferred transport reactions in the genome‐scale metabolic model for adipocytes, iAdipocytes1809, to ensure the connectivity.
Supplementary Dataset 5
Known biological function of the adipocytes and the metabolic capacity was demonstrated by the simulation of 250 metabolic functions based on the definition of functions in HepatoNet 1 (Gille et al., 2010).
Supplementary Dataset 6
McQuaid et al. (McQuaid et al., 2011) measured the delivery and transport of fatty acids in adipose tissue using multiple and simultaneous stable‐isotope fatty acid tracers in lean and obese subjects groups over a 24 hour period. Even though abdominally obese subjects have greater adipose tissue mass than lean control subjects, the rates of delivery of NEFAs are down regulated in obese subjects. In their study, uptake/secretion rates for NEFAs, TAGs and glucose in obese and lean subjects are reported. In adipose tissue, blood flow, glucose uptake, release of NEFAs and the extraction of TAGs from plasma was significantly lower in the abdominally obese subjects compared to lean subjects. Values are given as means ± SD.
Supplementary Dataset 7
Paterson et. al (Patterson et al., 2002) have studied the forearm and adipose tissue amino acid metabolism in lean and obese human subjects after 22 hours of fasting. They hypothesized that greater conservation of body protein observed during fasting in obese than in lean subjects and a combination of stable isotope tracer infusion and arteriovenous balance techniques was used to quantify amino acid kinetics. In their study, local net amino acid arteriovenous differences were calculated as the amino acid concentration in arterial plasma minus the concentration in venous plasma and local net fluxes were calculated as the arteriovenous difference multiplied by local plasma flow. Plasma arterial amino acid concentrations, regional subcutaneous abdominal adipose tissue arteriovenous concentration differences are adapted from the study of Paterson et. al (Patterson et al., 2002) and local net fluxes for lean and obese subjects were calculated.
Supplementary Dataset 8
Fatty acid (FA) composition in human plasma (Quehenberger et al., 2010), adipose tissue (Hodson et al., 2008; Raclot et al., 1997) and liver tissue (Shorten and Upreti, 2005) where there is not any information for adipocytes are incorporated in to the model for generation of pool reactions. Lipidomics analysis of a pooled human plasma obtained from healthy individuals after overnight fasting has been used in the model (Quehenberger et al., 2010). In their study, over 500 distinct molecular species distributed among the main lipid categories including fatty acyls, glycerolipids, glycerophospholipids, sphingolipids, sterols, and prenols were quantified. The reported 31 FAs in human plasma used in our model and normalized values are shown. In their study, some plasma FA known to exist in the structure of adipose tissue triacylglycerols and phospholipids have not been reported. These unreported FAs were included in to the structure of pool reactions by assuming their concentrations as half of the minimum concentration reported in the measurement study. Assumptions for the concentration of non reported FAs in lipidomics analysis allow us to include them into the genome‐scale metabolic model for adipocytes, iAdipocytes1809. Raclot et. al. (Raclot et al., 1997) measured the FA composition of TAG in adipose tissue and originated non‐esterified FA (NEFA) concentration from adipose TAGs through lipolysis. The mobilization of 34 individual fatty acids is quantified in adipose tissue of eight healthy non‐obese women and used in our model to generate the pool reactions. The reported FA content was normalized to 1 mol/mol of TAG and shown. Hodson et. al. (Hodson et al., 2008) collected the data for measurements of FA composition in adipose tissue and blood lipids and their changes with different diet. The reported FA composition of NEFA is compared with the lipidomics analysis of a pooled human plasma study (Quehenberger et al., 2010) and an agreement has been seen. The FA molar composition of plasma TAGs and CEs reported in the study is normalized to 1 and used in iAdipocytes1809.
Supplementary Dataset 9a
Lipid droplet formation in adipose tissue of lean and obese subjects have been predicted with the genome‐scale metabolic model for adipocytes, iAdipocytes1809 and exchange reactions have been presented. Uptake/secretion rates for NEFA, TAG and glucose in obese and lean subjects are adapted from McQuaid et al. (McQuaid et al., 2011).
Supplementary Dataset 9b
Acetyl‐CoA formation in the mitochondria of adipose tissue in lean and obese subjects have been predicted through the usage of our genome‐scale metabolic model for adipocytes, iAdipocytes1809 and exchange reactions have been presented. Uptake/secretion rates for NEFA, TAG and glucose in obese and lean subjects are adapted from McQuaid et al. (McQuaid et al., 2011)
This is an open‐access article distributed under the terms of the Creative Commons Attribution License, which permits distribution, and reproduction in any medium, provided the original author and source are credited. This license does not permit commercial exploitation without specific permission.
- Copyright © 2013 EMBO and Macmillan Publishers Limited