Diverse repertoire of human adipocyte subtypes develops from transcriptionally distinct mesenchymal progenitor cells.

Single-cell sequencing technologies have revealed an unexpectedly broad repertoire of cells required to mediate complex functions in multicellular organisms. Despite the multiple roles of adipose tissue in maintaining systemic metabolic homeostasis, adipocytes are thought to be largely homogenous with only 2 major subtypes recognized in humans so far. Here we report the existence and characteristics of 4 distinct human adipocyte subtypes, and of their respective mesenchymal progenitors. The phenotypes of these distinct adipocyte subtypes are differentially associated with key adipose tissue functions, including thermogenesis, lipid storage, and adipokine secretion. The transcriptomic signature of "brite/beige" thermogenic adipocytes reveals mechanisms for iron accumulation and protection from oxidative stress, necessary for mitochondrial biogenesis and respiration upon activation. Importantly, this signature is enriched in human supraclavicular adipose tissue, confirming that these cells comprise thermogenic depots in vivo, and explain previous findings of a rate-limiting role of iron in adipose tissue browning. The mesenchymal progenitors that give rise to beige/brite adipocytes express a unique set of cytokines and transcriptional regulators involved in immune cell modulation of adipose tissue browning. Unexpectedly, we also find adipocyte subtypes specialized for high-level expression of the adipokines adiponectin or leptin, associated with distinct transcription factors previously implicated in adipocyte differentiation. The finding of a broad adipocyte repertoire derived from a distinct set of mesenchymal progenitors, and of the transcriptional regulators that can control their development, provides a framework for understanding human adipose tissue function and role in metabolic disease.

Adipose tissue plays critical roles in systemic metabolism, through multiple functions that include lipid storage and release, as well as secretion of cytokines, such as leptin and adiponectin, which control satiety and energy utilization (1). In addition to its metabolic roles, adipose tissue plays a pivotal role in thermoregulation, both by providing an insulating layer under the skin and by the production of heat (2). Other functions of adipose tissue include mechanical and immune-protective roles in multiple sites, including the mesentery (3). These functions are mediated by adipocytes, of which 3 functionally distinct subtypes have been recognized, best defined in the mouse (4). White adipocytes are lipogenic and responsible for lipid storage and release, brown adipocytes form a defined depot restricted to the interscapular region and contain high basal levels of the mitochondrial uncoupling protein UCP1, and “brite/beige” adipocytes are interspersed within white adipose tissue and express UCP1 in response to stimulation. Lineage-tracing and gene-expression studies point to distinct developmental origins for these adipocyte subtypes (5, 6). In adult humans, no specific depot is solely composed of UCP1-containing adipocytes, but UCP-1+ cells can be found interspersed within supraclavicular, paravertebral, and perivascular adipose tissue (7), resembling the brite/beige phenotype. Thus, in humans, only 2 adipocyte types, UCP-1+ (thermogenic) and UCP-1− (white) have been recognized to exist.

The diversified nature of human adipose tissue is exemplified by differences in lipolytic capacity and cytokine secretion profiles between adipocytes from different depots (8–10). Their numerous functions may be mediated through modulating the 2 recognized adipocyte types, for example via differences in vascularization and innervation between different depots (4). Alternatively, additional yet uncharacterized developmentally and epigenetically distinct adipocyte subtypes that are not readily distinguished by their morphology may exist, similarly to the repertoire of functionally distinct immune cells that share similar morphologies. If so, the contribution of adipose tissue to whole-body homeostasis could result from the relative abundance of distinct adipocyte subtypes, in addition to differences in their functionality. Distinguishing between these possibilities is critical to our understanding of human adipose tissue biology and physiopathology.

Recent advances in the ability to study single-cell transcriptomes can be leveraged to better understand the cellular composition of adipose tissue (11–13). Single-cell transcriptomics of adipose tissue progenitor cells in mouse fat has revealed specific populations primed to expand in response to adrenergic stimulation (13), stromal cells capable of inhibiting adipogenesis (11), and cells with enhanced adipogenic or fibro-inflammatory properties (14). The use of single-cell transcriptomics to define mature adipocyte subtypes is complicated by the large size and high buoyancy of these cells, and the fact that the extensive proteolytic digestion required to isolate them from the adipose tissue compromises their viability and alters their phenotype (15). While elegant solutions to these limitations have been developed (16), they are dependent on transgenic approaches and are thus less amenable to application for human studies.

As a complementary approach to studying human adipocyte subtypes and their development, we leveraged a recent finding by our laboratory that allows us to generate large numbers of mesenchymal progenitor cells from human adipose tissue with little loss of multipotency (17–19). Many of these progenitor cells differentiate into large, long-lived unilocular adipocytes (18), and notably, some correspond to thermogenic adipocytes, as they respond to stimulation by inducing UCP1 and other brown adipose tissue markers (19). To identify the distinguishing features of progenitor cells that give rise to the adipocyte subtypes generated from these progenitors, we obtained transcriptomic profiles before and after adipogenic induction, and before and after thermogenic stimulation of clonal cell populations derived from single mesenchymal progenitor cells.

Our results uncover a minimum of 4 distinct adipocyte subtypes differing in gene-expression profiles associated with distinct, recognizable adipocyte metabolic functions. Progenitor cells that give rise to each adipocyte subtype express distinct complements of transcriptional regulators and cytokines. Our findings can explain numerous existing results of the effects of genetic ablation or overexpression of these factors on adipose tissue function. Thus, our results are predictive of in vivo physiology and suggest further testable models to elucidate human adipose tissue function and systemic metabolism.

Results

Morphological Evidence for Distinct Adipocyte Subtypes.

To explore the mechanisms by which human adipocyte subtypes originate, we developed an approach summarized in Fig. 1 and detailed in . Mesenchymal progenitor cells sprouting from human adipose tissue explants were plated as single cells in individual wells of 384-well multiwell plates, and proliferation was allowed until near confluence, when they were split 1:3. One well was maintained in a nondifferentiated state (“C”); 1 well was subjected to adipose differentiation by exposure to methylisobutylxanthine, dexamethasone and, insulin (MDI, “M”); and 1 was subjected to differentiation and subsequently stimulated with Forskolin (Fsk, “F”) for the last 3 d of culture. Wells were then imaged and RNA extracted. Phase images of clones could be categorized into 3 different classes. One consisted of cells that proliferated to confluence but failed to accumulate lipid droplets in response to MDI (Fig. 1); this class was not considered further in this study. A second class accumulated lipid droplets but was not visibly affected by Fsk (Fig. 1). A third class accumulated lipid droplets and responded to Fsk with clear changes in droplet size and number (Fig. 1). To determine whether these differences were explained by stochastic variation, or whether they reflected the existence of true distinct adipocyte subtypes, we analyzed the frequency distributions of droplet size, measured after image processing as exemplified in Fig. 1 , Middle and Bottom. Aggregated values for lipid droplet size from all clones in the M state (Fig. 1) produced a monophasic histogram representative of the normal variation of adipocyte droplet sizes in these cells collectively. However, when we plotted the mean droplet size per clone, we obtained a nonmonophasic distribution (Fig. 1), revealing the existence of adipocyte subtypes that differ in mean droplet size. We then quantified the responsiveness of each clone by comparing the size of lipid droplets for each clone in the M and F states. While a decrease in droplet size was seen in most clones (Fig. 1), Fsk responsiveness was not represented by a monophasic distribution, as the histogram of the values (droplet size in the F state as percentage of size in the M state) contained at least 2 peaks, a major peak centered around 100% (no response) and another at around 60% (decrease of lipid droplet size to 60% of nonstimulated value) (Fig. 1). These results indicate that at least 2 adipocyte subtypes exist, based on their mean droplet size and differential responsiveness to Fsk.

Fig. 1.

Single mesenchymal progenitors give rise to morphologically different adipocyte subtypes. (A) Mesenchymal progenitor cells from adipose tissue are induced to proliferate under proangiogenic conditions in a 3D hydrogel. Single cells are plated into single wells of 384-well multiwell plates, and proliferation allowed until near confluence, when they are split 1:3. One well is maintained in a nondifferentiated state (C), 2 wells are subjected to adipose differentiation (M), and 1 of the differentiate wells is stimulated with Fsk for the final 3 d of culture (F). Wells are then imaged and RNA extracted. (B) Example of a clone that failed to undergo adipose differentiation. (C) Example of a clone that underwent adipocyte differentiation but was unresponsive to Fsk as assessed by lipid droplet size. (D) Example of a clone that underwent adipocyte differentiation and responded to Fsk with decrease in lipid droplet size. (Scale bars: A, 200 μm; B–D, 50 μm and enlargements, 100 μm.) (E) Frequency distribution of all lipid droplets measured in all clones in the M state. (F) Frequency distribution of the mean droplet size per clone. (G) Measurement of mean droplet sizes per clone comparing M and F states. (H) Frequency distribution of responsiveness to Fsk, where mean lipid droplet size is expressed as percent of droplet size in M for each clone.

Identification of Adipocyte Subtypes through Transcriptional Profiling.

We then asked whether the morphological evidence for distinct human adipocyte subtypes is supported by differences in gene expression. We conducted RNA sequencing from 156 samples comprising the nondifferentiated (“C”), differentiated (M), and differentiated plus thermogenically activated (F) states for each of the 52 clones displaying adipose differentiation. Principal component analysis of aligned reads revealed clear segregation of C, M, and F states (Fig. 2), verifying that differentiation and Fsk stimulation had major effects on gene expression in these clones. Interestingly, a broader spread is seen in the first principal component of M and F compared to C states, indicating generation of divergent transcriptional landscapes following induction of differentiation and thermogenic stimulation.

Fig. 2.

Identification of thermogenic adipocytes within distinct adipocyte subtypes. (A) Principal component analysis of 52 clones in the C, M, or F states in black, green, and red symbols, respectively. (B) Unsupervised hierarchical clustering of 447 genes correlated with lipid-droplet size. Each row is independently scaled from blue (lower) to red (higher) values. (C and D) Transcripts per million (TPM) values for PPARG in each clone (C) and mean and SEM of TPM values in each cluster (D) in the C, M, and F states. Statistical comparison was done using 1-way ANOVA controlling for multiple comparison testing with the Holm–Sidak method. *P = 0.014; ***P < 0.0006; ****P < 0.0001. (E) RT-PCR of thermogenic genes at different times of exposure of heterogeneous adipocyte cultures to Fsk. Ct values for UCP1 were 29 to 31 at t = 0 and decreased to 20 to 22 at t = 6 h. Values represent fold relative to the lowest value. Shown are means and SEM of 2 independent experiments using cells from 2 independent subjects assayed in duplicate. (F) Immunofluorescence staining to mitochondrial Hsp70 (red) and UCP1 (green) without or with exposure to Fsk for 7 d. (Scale bars, 200 μm.) (G–I) Mean and SEM of TPM values from each cluster for the genes indicated in the y axis. Statistical comparison was done as described in D. In G, **P = 0.0097; in H, **P = 0.0067, ****P < 0.0001. (J–L) Volcano plot (J) and enrichment analysis (K and L) of differential expression comparing cluster 2 with all other clusters. (M) Hierarchical clustering of 1,175 genes showing the highest variance within adipocytes from neck (NeckSQ) or abdominal (AbdSQ) subcutaneous depots. (N) Values for genes indicated in the x axis in adipocytes from NeckSQ and AbdSQ. Red and green color text indicates genes decreased or increased, respectively in cluster 2 relative to other clusters. Red downward arrows and green upward arrows indicate whether the gene is decreased or increased respectively in NeckSQ relative to AbdSQ.

To determine whether the divergent transcriptomes of samples in the M state could classify adipocyte subtypes, we employed a clustering approach. We focused on genes with expressions that correlated with lipid-droplet size (as assessed in Fig. 1), a key property of adipocytes, with a Pearson rank-order correlation coefficient higher than 0.5 or lower than −0.5 (Dataset S1). Unsupervised hierarchical clustering of the expression of these 447 genes resulted in 4 major clusters (Fig. 2), indicating 4 distinct adipocyte subtypes. To verify the robustness of these results, we performed an additional, unsupervised clustering analysis using the top 1,500 genes displaying the highest variance. This analysis produced 5 distinct clusters, largely overlapping with the clusters generated by hierarchical clustering of genes correlating with lipid-droplet size ().

To better understand the differences between the adipocyte subtypes defined by the clustering, we looked at the expression levels of the major regulatory gene PPARG (Fig. 2 ). The levels of PPARG were higher in cluster 1 compared to all other clusters in the M state (Fig. 2 , green bars), despite similar levels in the C state (Fig. 2 , black bars). Additionally, in the F state the levels of PPARG in cluster 2 were elevated (Fig. 2 , red bars) but not in cluster 1, and less in other clusters. The different levels and responsiveness of PPARG indicate that these clusters represent functionally distinct adipocyte types.

Identification of Cluster 2 as Corresponding to the Thermogenic Brite/Beige Adipocyte Subtype.

The higher responsiveness of cluster 2 to Fsk suggested these clones could correspond to thermogenic beige/brite adipocytes. In heterogeneous cultures, Fsk strongly induced the expression of thermogenic genes UCP1, DIO2, and CIDEC after acute (6 to 24 h) and chronic (24 h to 7 d) exposure (Fig. 2). Notably, immunofluorescence staining of UCP1 was heterogeneous, consistent with brite/beige adipocytes representing a subtype within the population (Fig. 2). For unclear reasons, despite strong signals by qRT-PCR, very low UCP1 reads were detected by RNA sequencing in individual clones. Nevertheless, DIO2 and CIDEC expression were clearly detected in all clones and were significantly more strongly induced in cluster 2 (Fig. 2 ), consistent with its identity as brite/beige adipocyte subtype. The higher responsiveness to Fsk seen in cluster 2 was specific to thermogenic genes, as no differences were observed on housekeeping genes, such as GAPDH (Fig. 2).

Having identified the brite/beige adipocyte subtype on the basis of responsiveness of thermogenic genes to Fsk allowed us to determine other features that distinguish these cells prior to simulation. Differential expression analysis between cluster 2 and the combined clusters 1, 3, and 4 in the M state revealed 12 genes increased and 32 genes decreased in level of expression by 1.5-fold or more with a P-adjusted (P-adj) value <0.05 (Fig. 2 and Dataset S2). Genes underexpressed in cluster 2 were enriched in the ferroptosis and mineral absorption pathways, including the ferroportin gene SLC40A1|FPN1 (Fig. 2 and Dataset S3), which is the only mechanism known to mediate iron egress from cells. Genes increased in cluster 2 were enriched in pathways related to carbohydrate, carbon, and nitrogen metabolism (Fig. 2), and the gene family of carbonic anhydrases was highly enriched (P-adj = 4.505E-4) by the simultaneous presence of CA3 and CA9 (Fig. 2 and Dataset S3). These results suggest that beige/brite adipocytes are primed to accumulate iron and withstand oxidative stress, both necessary to adapt to mitochondrial biogenesis and fatty acid oxidation upon thermogenic stimulation.

A key question emerging from these findings is whether adipocyte subtypes defined by these clusters represent those found in humans in vivo. To address this question, we asked whether genes that are differentially expressed in cluster 2 would be similarly differentially expressed in thermogenic compared to nonthermogenic depots. We leveraged datasets previously generated from human depots using Affymetrix HTA-2 arrays. Adipocytes generated from supraclavicular (NeckSQ) and abdominal (AbdSQ) adipose depots had distinct transcriptomic signatures (Fig. 2), and genes differentially expressed in cluster 2 were consistently differentially expressed in adipocytes from supraclavicular adipose tissue with a similar directionality (Fig. 2). To more comprehensively compare the datasets, we used RankProd, which allows a nonparametric comparison of datasets generated on different platforms (20, 21): 832 genes were significantly enriched, and 699 significantly decreased (percentage of false predictions <0.05) in both cluster 2 and adipocytes from the supraclavicular region (Dataset S4), demonstrating remarkable correspondence between cluster 2 adipocytes and adipocytes from thermogenic depots.

Another key question we addressed is whether thermogenic adipocytes identified through transcriptional profile are functionally distinct. Fatty acids derived from endogenous hydrolysis of lipid-droplet triglycerides represent a major source of fuel during thermogenesis. Reflecting this process, we observed significant glycerol release stimulated by Fsk treatment of adipocytes differentiated from unsorted progenitors (Fig. 3). Notably, images of live adipocytes over time revealed marked heterogeneity in lipolytic responses between cells, as assessed by the size and number of lipid droplets; over time a clear distinction could be seen between adipocytes in which lipid droplets continued to increase and become unilocular (Fig. 3, arrowheads) and those in which lipid droplets became smaller and very numerous (Fig. 3, arrows). The different responses of lipid droplets to Fsk were also evidenced when cells were stained with an antibody to perilipin-1 (PLIN1); in the absence of Fsk stimulation, PLIN1 was evenly distributed around lipid droplets, and there was little heterogeneity between adipocytes (Fig. 3 , Upper). However, in response to Fsk, some cells contained large unilocular lipid droplets with evenly distributed PLIN1 (Fig. 3 , Lower, arrowhead), while others contained numerous small droplets with punctate PLIN1 staining (Fig. 3 , Lower, arrows). To determine whether cells developing the multilocular phenotype correspond to cluster 2 thermogenic adipocytes, we performed in situ hybridization for the long-noncoding RNA LINC00473, which is highly enriched in human thermogenic adipose tissue (22) and highly enriched in cluster 2 (Fig. 3 ). Fluorescence in situ hybridization demonstrated LINC00473 was expressed predominantly in cells containing small lipid droplets (Fig. 3 ). These results demonstrate that cluster 2 adipocytes are functionally distinct, developing multilocular lipid droplets and increased lipolysis in response to cAMP, consistent with their thermogenic phenotype.

Fig. 3.

Functional differences in lipid-droplet dynamics in the brite/beige adipocyte subtype. (A) Concentration-dependent release of glycerol into the medium by adipocytes incubated with Fsk for 24 h. (B) Lipid droplets visualized in the same field after the indicated times in the presence of Fsk. Arrowhead points to an adipocyte refractory to Fsk containing a large lipid droplet, and arrows to adipocytes containing small multilocular droplets developing in response to Fsk. (C) PLIN1 staining of adipocytes incubated without (Upper) or with (Lower) Fsk for 7 d. Arrowhead points to an adipocyte refractory to Fsk containing a large lipid droplet, and arrows to adipocytes containing small multilocular droplets developing in response to Fsk. (Scale bars, 100 μm.) (D and E) TPM values for LINC00473 in each clone (D) and mean and SEM of TPM values in each cluster (E) in the M state. Statistical comparison was done using 1-way ANOVA controlling for multiple comparison testing with the Holm–Sidak method. (F) PLIN1 staining (red) following in situ hybridization of LINC00473 (green) in adipocytes stimulated with 10 μM Fsk for 6 h. (Scale bar, 50 μm.) (G) Mean number of LINC00473 spots in cells in 10 independent fields; LINC+ cells were defined as those with more than 10 spots. (H) Mean droplet size in LINC− and LINC+ cells. Means and SEM are indicated. Statistical significance of differences between LINC− and LINC+ were calculated using the Wilcoxon matched-pairs signed rank test as implemented in Prism 8.

Having defined clones comprising cluster 2 as corresponding to thermogenic adipocytes, we then asked whether their progenitors could help define mechanisms that control the development of this adipocyte subtype. For this, we compared cluster 2 in the C state with all others. We find 21 genes with lower expression and 18 with 2-fold or higher expression with a P-adj value of 0.05 (Fig. 4 and Dataset S5). Genes decreased in cluster 2 mapped significantly to MAD box and basic helix–loop–helix gene families (MEF2C, ID1, and NPAS1) (Fig. 4 and Dataset S6). This finding is consistent with gene-ablation studies that have identified Id1 as a negative regulator of thermogenic adipocyte development in mouse (23), further supporting the identity of cluster 2 progenitors as corresponding to those for thermogenic adipocytes in vivo. Interestingly, many genes enriched in cluster 2 progenitors corresponded to cytokines, including CSF3, IL1B, IL11, and CXCL8, and pathway enrichment analysis mapped to the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway for IL-17 signaling (Fig. 4 and Dataset S6). Expression of these cytokines was suppressed upon differentiation (Fig. 4 ), suggesting that their function may be important in defining intercellular communications with immune cells known to be important for the establishment of a thermogenic niche (24).

Fig. 4.

Features of progenitor cells that give rise to the brite/beige adipocyte subtype. (A) Volcano plot of differential expression analysis comparing progenitors corresponding to cluster 2 with all others. (B) KEGG pathway analysis of genes underrepresented in cluster 2. (C) KEGG pathway analysis of genes overrepresented in cluster 2. (D–H) Mean and SEM of TPM values from each cluster for the genes indicated in the y axis. Statistical comparison between each column was done using 1-way ANOVA with correction for multiple comparisons using the Holm–Sidak method as implemented in GraphPad Prism 7. *P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001.

To further explore the relevance of cluster 2 progenitors we focused on the finding that the FTO obesity-associated locus controls IRX3 expression during early adipocyte differentiation, ultimately resulting in altered thermogenesis (25). We find that IRX3 is significantly lower in cluster 2 progenitors compared to all others, is increased after differentiation, and suppressed by Fsk (Fig. 4), consistent with its role in suppressing thermogenic induction. These results further highlight the relevance of cluster 2 adipocytes and their progenitors as a class of cells central to the control of whole-body energy metabolism.

Adipocyte Subtypes Specialized for Regulated Expression of Adiponectin and Leptin.

To further explore the differences between adipocyte subtypes we looked at the expression levels of genes central to adipose tissue endocrine function, ADIPOQ and LEP. The values for LEP or ADIPOQ expression in all clones in all states revealed a highly heterogenous distribution (Fig. 5 ). The average of the values across each cluster indicate that ADIPOQ was expressed at significantly higher levels in cluster 1 compared to all others, and its expression was strongly stimulated by Fsk in all clusters (Fig. 5). More surprisingly, LEP expression did not parallel that of ADIPOQ, trending to be highest in cluster 4. Moreover, responsiveness of LEP to Fsk was only seen in cluster 4, where expression was strongly suppressed (Fig. 5). The differential expression of LEP was not associated with adipocyte lipid-droplet content, which was similar between all clusters (Fig. 5 ).

Fig. 5.

Evidence for adipocyte subtypes specialized for ADIPOQ and LEP expression. (A–F) TPM values for ADIPOQ (A) and LEP (C), and droplet-size values (E) in each clone, separated by clusters as indicate below the x axis. Black, green, and red bars represent values in the C, M, and F states, respectively. Mean and SEM of TPM values for ADIPOQ (B) and LEP (D) and droplet-size values (F) in each cluster in the C, M, and F states. Statistical comparison between columns was done using 1-way ANOVA with correction for multiple comparisons using the Holm–Sidak method as implemented in GraphPad Prism 7. *P < 0.05, **P < 0.01, ***P < 0.001. (G) Immunofluorescence staining of adipocytes from 2 individuals (subjects 37 and 45) after 14 d of differentiation. Arrows indicate adipocytes expressing mostly leptin (green), mostly adiponectin (red), or both (yellow). (Scale bars, 100 μm.) (H and I) Cumulative frequency distribution of staining intensities per cell for ADIPOQ (H) or LEP (I) from 100 cells in fields exemplified in G. (J and K) qRT-PCR for ADIPOQ (J) or LEP (K) in cultures exemplified in G, before or after 14 d of differentiation. Values are expressed as fold-difference relative to the lowest detectable value for each gene. **P < 0.001. (L and M) qRT-PCR for ADIPOQ (L, left y axis) or LEP (M, left y axis), and corresponding concentration of ADIPOQ (L, right axis) and LEP (M, right axis) in medium after 24 h of culture. Values represent the mean and SEM of 2 independent cultures assayed in duplicate for each subject. (N) Immunofluorescence staining of primary human adipocytes isolated by collagenase digestion. The single field shown includes 2 cells, 1 mostly expressing LEP (green) and the other ADIPOQ (red). (Scale bars: single optical section and projection, 100 μm.) Top and Middle are a single optical plane through the cells, and Bottom is the projection of all optical planes comprising the 3D volume of the cells. Nuclei are stained in blue.

To examine whether the observed differential expression of ADIPOQ and LEP between clusters translates into differences in expressed adipokines, we performed immunofluorescence staining on adipocytes differentiated from unsorted progenitors from AbdSQ tissue from 2 randomly selected subjects. In both cases, cells expressing LEP and very little ADIPOQ, and vice versa, could be readily detected (Fig. 5). Cumulative histograms of fluorescence intensity of 100 cells per coverslip revealed differences between subjects in the number of cells expressing very high levels of either ADIPOQ or LEP (Fig. 5 , indicated by ovals). Notably, the differences in ADIPOQ and LEP mRNA levels between the 2 subjects, assessed by qRT-PCR (Fig. 5 ) (differences in mean values 18,077 and 115,516 for ADIPOQ and 146 and 78 for LEP for subject 37 and subject 45, respectively) were not fully explained by the differences in mean expression level per cell measured by immunofluorescence (differences in mean values 242 and 300 for ADIPOQ, and 280 and 190 for LEP for Hu37 and Hu45, respectively), implying that the number of high adipokine-expressing cells contribute significantly to the mean levels of each cytokine within a population. In an additional 2 subjects, we examined whether the detected expression in cells would translate to differences in levels of secreted adipokines. ADIPOQ mRNA levels were induced to a similar level in adipocytes from subjects 49 and 50 (Fig. 5), but LEP mRNA was higher in adipocytes from subject 50 (Fig. 5). These differences were closely reflected by the levels of secreted ADIPOQ and LEP, which could be readily measured by human-specific ELISAs in supernatants 24 h after media replacement.

To verify that the observed differences in adipokine expression were present in mature adipocytes under physiological conditions we tested whether adipocytes preferentially expressing LEP or ADIPOQ can be found in freshly isolated primary cells. We stained cells obtained by collagenase digestion of subcutaneous adipose tissue and found adipocytes expressing high levels of either ADIPOQ or LEP (Fig. 5). Interestingly, staining of adipocytes from mouse epididymal adipose tissue has already been remarked as being similarly heterogeneous, with some cells expressing no detectable leptin (26). Thus, both transcriptional evidence and single-cell immunostaining of primary cultured and freshly isolated cells support the existence of distinct subtypes of white adipocytes specialized for production of ADIPOQ or LEP.

To investigate other differences between adipocyte subtypes expressing high levels of ADIPOQ or LEP we performed differential expression analysis between clusters 1 and 4 in the M state. We find 174 genes were more highly expressed by >5-fold in cluster 1, and 48 genes were more highly expressed >5-fold in cluster 4 (Fig. 6 and Dataset S7). Pathway enrichment analysis of genes more highly expressed in cluster 1 included PPARG signaling and lipogenesis pathways, containing genes such as DGAT2, FABP4, PLIN1, and PLIN4 (Fig. 6 and Dataset S8). In contrast, genes more highly expressed in cluster 4 were associated with extracellular matrix interaction and TGF-β signaling (Fig. 6 and Dataset S8). These results suggest that cluster 1 comprises adipocytes specialized for lipogenesis, while cluster 4 comprises adipocytes that may be more sensitive to extracellular matrix and regulatory signaling.

Fig. 6.

Adipocyte subtypes specialized for ADIPOQ and LEP expression differ in metabolic and transcription factor gene expression. (A) Volcano plot comparing cluster 1 and cluster 4. Orange and blue symbols are genes enriched in cluster 1 or cluster 4, respectively. (B and C) KEGG pathway enrichment of enriched genes in cluster 1 or cluster 4, respectively. (D–I) Correlation of transcription factors showing highest correlation with ADIPOQ, indicated in the y axis, with ADIPOQ (D–F) or LEP (G–I). (J–M) Correlation of transcription factors showing highest correlation with LEP, indicated in the y axis, with LEP (J and K) or ADIPOQ (L and M).

We next searched for mechanisms that might differentially control ADIPOQ or LEP production by defining genes encoding for transcriptional regulators that correlate with expression of these adipokines across all clones in the M and F states. Genes correlating with Spearman correlation coefficients higher or lower than 0.7 or −0.7, respectively, were considered, based on ranges for spurious correlation defined in (Statistics). Within these genes we found transcriptional regulators strongly linearly correlated with ADIPOQ expression, including known adipogenic genes PPARG, CEBPA, and RXRA (Fig. 6 ). These correlated weakly or not at all with LEP expression (Fig. 6 ). However, 2 transcription factors, ZNF395 and DDX41, correlated strongly with LEP (Fig. 6 ) but weakly or not at all with ADIPOQ expression (Fig. 6 ). Importantly, ZNF395 has been previously identified in unbiased screens as a zinc-finger protein that promotes adipocyte differentiation in vitro (27, 28). Further experiments will determine whether this transcription factor specifically supports development of LEP producing adipocytes, as suggested by these results.

Discussion

The major finding of this study is the identification of a minimum of 4 distinct human adipocyte subtypes that can differentiate from mesenchymal progenitor cells. Moreover, we identify the specific progenitor populations for each adipocyte subtype, providing insights into their developmental mechanisms. These findings expand our understanding of adipocyte heterogeneity from white vs. brown/brite/beige to additional subtypes that can collectively mediate the multiple diverse functions of adipose tissue. The approach we used in this work was to analyze the transcriptomes of cell populations derived from single human mesenchymal multipotent progenitor cells before and after adipogenic differentiation and before and after thermogenic activation. A significant advantage of this approach is that stochastic variations displayed by single cells, such as in degree and timing of differentiation, are averaged out in each clone. Thus, each adipocyte subtype is defined by multiple clones in which the stochastic responsiveness of individual cells to adipogenic and thermogenic stimulation is eliminated. The input of RNA is also higher than that possible from single cells, thus allowing deeper sequencing and phenotypic characterization.

Our finding of distinct adipocyte subtypes is consistent with prior observations of heterogeneity in populations of 3T3-L1 adipocytes (29–33), and within adipocytes from human and mouse adipose tissue (34). While these prior observations could be potentially attributed to cell-to-cell stochastic variations, our results show that underlying gene expression architecture indeed reflects the existence of human adipocyte subtypes. Moreover, our findings trace these differences to progenitor cells that give rise to each subtype, providing further confidence in the findings.

As expected from prior analyses of adipocyte populations (19), one of the subtypes we identified corresponds to thermogenic brite/beige human adipocytes, which are evidenced by the induction of the thermogenic phenotype in response to elevation in cAMP. The features of these adipocytes prior to stimulation are predicted to enable rapid increases in mitochondrial biogenesis and oxidative stress following stimulation. For example, SLC40A1/FPN1, the only known mechanisms for iron efflux (35), is significantly lower in thermogenic compared to nonthermogenic adipocytes prior to stimulation. These results can explain previous observations of increased levels of genes leading to iron efflux in response to overfeeding, and their reversal with weight loss (36, 37), both conditions associated with decreased and increased content of thermogenic adipose tissue in humans (38). Iron is a rate-limiting and regulatory factor in adipose tissue browning and mitochondrial respiration (39, 40), and is tightly associated with redox stress. It is interesting that predicted increased iron accumulation is accompanied by increased expression of carbonic anhydrases CA3 and CA9, which can potently protect cells from oxidative stress (41–43). Moreover, CA3 is typically expressed at much higher levels in human skeletal muscle (43, 44), consistent with the notion that thermogenic adipose tissue and muscle share developmental features (6). The parallel down-regulation of SLC40A1/FPN1 to increase cellular iron, and up-regulation of CA3 and CA9 to protect from redox stress, provide a permissive environment for rapid mitochondrial biogenesis and enhanced respiratory flux.

Our experimental paradigm allows us to explore the potential mechanisms by which each adipocyte can be generated from progenitor cells. Brite/beige adipocyte progenitors express high levels of cytokine genes, which could modulate the local concentrations of immune and stromal cells, which in turn could affect fate determination (Fig. 7). Indeed, multiple reports exist of paracrine interactions with immune cells that enhance thermogenic activation (45–47). In addition to cytokines, progenitors that give rise to thermogenic adipocytes display enhanced expression of PTGS2, the rate-limiting step in prostaglandin synthesis. PTGS2 has been directly associated with formation of brown adipocytes in mice (48), and partial PTGS2 deficiency leads to obesity (49). The 21 genes decreased in thermogenic adipocyte progenitors mapped significantly to MAD box and basic helix–loop–helix gene families (MEF2C, ID1, and NPAS1). Id1 depletion results in increased generation of thermogenic adipose tissue in mice, consistent with its lower levels in human progenitors fated for the brite/beige phenotype (23, 50).

Fig. 7.

Models derived from clustering and functional data. (A) Model for mechanism by which cytokine secretion by progenitors of beige/brite adipocytes create an immunological niche conductive to thermogenesis, different from that created by nonthermogenic adipocyte progenitor subtypes. (B) Model for mechanisms of regulation of adipokine levels, by modulating expression levels in all cells (Left), or by modulating the number of cells expressing each adipokine (Right).

While we expected to find 2 adipocyte subtypes, corresponding to white and brite/beige adipocytes, we were surprised to find additional subtypes differing in key adipocyte functions, such as leptin and adiponectin expression. Our data indicate that levels of LEP expression in mixed adipocyte populations are explained by both the number of adipocytes expressing high LEP levels, and the mean levels of expression per cell (Fig. 5 ). Adipocytes specialized for LEP expression are highly responsive to Fsk (Fig. 5), consistent with these cells being capable of acute regulation of LEP levels in response to fasting and catecholamines. These results imply that expression levels of adipokines in human adipose tissue and in the circulation may be determined by 2 potential mechanisms: 1 consisting in variation in the number of adipocytes expressing high levels of each adipokine, and another in variation in the mean levels of expression per adipocyte (Fig. 7). This mode of regulation may explain why systemic leptin, which varies with degree of adiposity, correlates but is not fully explained by adipocyte size or lipid content (51, 52). Further experiments will test the model described in Fig. 7, where adiponectin and leptin levels in vivo are defined by both number of adipocytes expressing high adipokine levels and the amount of adipokine within each cell.

The specific functions and relevance to human physiology of the diverse adipocyte subtypes identified by this work requires prospective isolation of these cells from adipose tissue, and investigation of differences in their number and properties between human populations. Up to now cell surface markers used in studies of mouse and human adipose tissues are not sufficiently differentially expressed among the 4 subtypes we have identified. Thus, new strategies are required to generate antibodies and antibody combinations to cell surface proteins that are sufficiently differentially expressed among the adipocyte and progenitor subtypes identified here. The depth of sequencing of adipocyte progenitor subtypes obtained by the methods used in this report provides necessary information about less abundant transcripts that can guide the identification of cell-specific marker combinations.

In summary, our results provide evidence for multiple human adipocyte subtypes specialized for key functions associated with adipose tissue, arising from progenitors that can be distinguished by their complement of expressed cytokines and transcriptional regulatory factors. Many of our results are consistent with published observations of the role of these regulatory factors, and many additional results provide the basis for exploration of additional currently unknown mechanisms of adipocyte differentiation and adipose tissue function. The concept that adipose tissue function results from the ensemble of multiple subtypes of cells, rather than from modulation of function of a single cell type, will help deconvolve the complex relationships between adipose tissue biology and metabolic disease.

Methods

We collected neck subcutaneous adipose tissues from carotid endarterectomies with no a priori selection of individual donors. Abdominal subcutaneous adipose tissue was obtained from discarded tissue following panniculectomy. All subjects consented to the use of tissue and all procedures were approved by the University of Massachusetts Institutional Review Board.

Cells sprouting from adipose tissue explants were obtained as previously described (19), and cells sorted as singlets into 384-well plates. Clones that grew to confluency were further split 1:3, as described in . Cells in the control group (C condition) were maintained in DMEM + 10% FBS. Adipogenic differentiation of 2 of the 3 wells was induced by providing confluent cells with DMEM + 10% FBS supplemented with 0.5 mM 3-isobutyl-1-methylxanthine, 1 μM dexamethasone, and 1 μg/mL insulin (MDI). After 72 h the media was replaced with DMEM + 10% FBS. Fifty percent of the media was replaced with fresh media every 48 h until day 14. At that point, 1 of the 2 wells were stimulated by adding 10 μM Fsk every 24 h for 72 h. RNA extraction, library preparation, sequencing, expression data analysis, and other nonspecialized methods are described in .