Metabolite profiling identifies anandamide as a biomarker of nonalcoholic steatohepatitis

Nonalcoholic fatty liver disease (NAFLD) is characterized by histology or imaging-defined hepatic steatosis, in the absence of secondary causes (1). NAFLD comprises a spectrum of pathology, ranging from steatosis to nonalcoholic steatohepatitis (NASH), which is a progressive form of NAFLD that can advance to cirrhosis and hepatocellular carcinoma (1). The development of NAFLD is linked to elevated BMI, central obesity, and the cardiometabolic syndrome (2, 3). Moreover, NAFLD may be a predictor for the development of diabetes mellitus and is an independent risk factor for cardiovascular disease (4, 5). Although NAFLD is connected to metabolic derangements associated with obesity, the pathways that link NAFLD, obesity, and cardiometabolic disease are poorly understood.

Metabolite profiling is an approach that has previously identified novel associations in diabetes mellitus and other cardiometabolic diseases in large community-based cohorts (6–11). To further elucidate the metabolic links among obesity, cardiometabolic disease, and NAFLD, we developed a metabolite-profiling method that detects both positively and negatively charged analytes from a single sample preparation. Coupling amide chromatography (12) and tandem mass spectrometry (13), we optimized conditions for profiling of human plasma samples.

Using this metabolite-profiling method, we evaluated the association between metabolites and cardiometabolic traits in the Framingham Heart Study Generation 3 (FHS Gen 3) cohort. We identified anandamide (AEA) as a candidate metabolite for cardiometabolic traits and then extended our findings in a hospital-based NAFLD cohort that demonstrated that AEA was associated with NAFLD severity. Cumulatively, these studies identify AEA as a candidate biomarker linking obesity and NAFLD.

Metabolites associated with cardiometabolic traits in the Framingham Heart Study. Using a multifunctional approach that detects both positively and negatively charged metabolites, we developed a metabolite-profiling method that assessed 179 analytes across important metabolic pathways (see Supplemental Methods for details; supplemental material available online with this article; https://doi.org/10.1172/jci.insight.92989DS1). We analyzed plasma samples of 997 participants in FHS Gen 3, a longitudinal, community-based study that began enrollment in 2002. The clinical characteristics of the cohort are shown in Table 1. Participants had no overt cardiovascular disease, had a mean age of 40 ± 9 years, a BMI of 26.6 ± 5.3 kg/m², and a homeostatic model assessment of insulin resistance (HOMA-IR) value of 1.04 ± 0.83.

Table 1

Characteristics of Framingham Heart Study Generation 3 participants at exam 1

In regression analyses adjusted for age and sex, we assessed metabolite associations with cardiometabolic traits, including generalized obesity, as measured by BMI; central obesity, as assessed by waist circumference (WC); and glucose metabolism, as assessed by fasting glucose level and HOMA-IR. Associations that exceeded the Bonferroni-corrected P value threshold are presented in Table 2. The full set of associations across cardiometabolic phenotypes, which included systolic and diastolic blood pressure, HDL, and triglycerides, is further presented in Supplemental Tables 2–9.

Table 2

Metabolite-phenotype associations in the Framingham Heart Study

Several known associations with obesity and insulin metabolism were validated, including uric acid (14); 2-aminoadipic acid (15); glutamate (16); the branched-chain amino acids leucine, isoleucine, and valine (17, 18); the short-chain acyl carnitines C3- and C6-carnitine (18); the tryptophan catabolite kynurenic acid (16); and amino acids tyrosine, alanine, and glycine (17). We also identified associations with metabolites that map to pathways previously associated with cardiometabolic disease, including the branched-chain ketoacid 2-ketoisovaleric acid and the tryptophan catabolites quinolinic acid, xanthurenic acid, indole-3-propionic acid, and kynurenine. Intermediate products of purine degradation and uric acid production included inosine, hypoxanthine, and xanthine. Several analytes that spanned glucose metabolism and the tricarboxylic acid cycle were elevated, including glucose, lactic acid, pyruvic acid, α-ketoglutaric acid, aconitic acid, and malic acid. Taken together, the identification of known metabolite-phenotype associations and pathways confirmed the validity of our approach.

AEA and associations with obesity and cardiometabolic traits. In the FHS Gen 3, analytes associated with cardiometabolic disease traits included uridine diphosphate N-acetylglucosamine (UDP-GlcNAc), reduced glutathione, hippuric acid, and AEA (see Table 2). Of these four candidates, AEA was notable for its association across metabolic traits, including generalized and central obesity and glucose metabolism, with the greatest association with markers of obesity. Linear regression with BMI yielded an adjusted β coefficient of 0.22 (P = 3.6 × 10^–13); for WC, the age- and sex-adjusted β coefficient was 0.20 (P = 1.6 × 10^–12). In age- and sex-adjusted logistic regression analyses, for each unit increase in AEA, the odds ratio (OR) for obesity (defined as BMI ≥ 30) was 3.20 (95% CI 2.10–4.88, P = 6.3 × 10^–8) and for excessive WC (defined as ≥88 cm), the OR was 2.79 (95% CI 1.85–4.20, P = 9.2 × 10^–7). AEA was also associated with elevated fasting serum glucose (glucose >110 mg/dl; OR 3.12, 95% CI 1.80–5.40, P = 4.9 × 10^–5) and abnormal HOMA-IR (HOMA-IR ≥2.0, OR 2.16, 95% CI 1.51–3.09, P = 2.5 × 10^–5) after adjustment for age and sex.

We sought to extend our work on AEA by assessing for associations with NAFLD, an obesity-related phenotype. Prior animal studies have suggested that the putative AEA receptor CB1 may modulate the accumulation of hepatic fat (19, 20), and a recent case-control cohort study identified an association between the endocannabinoid system and NAFLD (21). We therefore assessed for an association between AEA and radiographic hepatic steatosis in the FHS Gen 3 cohort and found an age- and sex-adjusted β coefficient of –0.013 (P = 0.03).

AEA levels correlate with the severity of NASH. The observed associations in the FHS Gen 3 cohort were present in a cohort of individuals free of overt disease, and the effect size with hepatic steatosis was accordingly subtle. To verify the specificity and extend our observations, we next analyzed a hospital-based, case-control cohort that had more advanced and precisely phenotyped liver pathology. These patients consisted of obese individuals with biopsy-proven NASH (cases) compared with those with obesity and normal liver histology (controls). The clinical characteristics of the age-, sex-, and BMI-matched cohort are listed in Table 3, and the histologic characteristics of the liver biopsies are summarized in Supplemental Table 10. In univariate analysis, AEA was significantly higher in the NASH group compared with the control group (Table 3, P < 0.002). Multiple logistic regression demonstrated that AEA remained an independent predictor of NASH (OR 5.2, 95% CI 1.6-17.0, P = 0.007), after adjusting for known risk factors, including HDL, triglyceride level, HOMA-IR, and alanine aminotransferase (ALT) (Table 4). AEA was also associated with increased severity of NASH, as assessed by the NAFLD activity score (OR 4.8, 95% CI 1.6–13.8, P = 0.005) and the presence of fibrosis (22) (OR 2.6, 95% CI 1.01–6.5, P = 0.04; see Supplemental Table 10). These findings suggest that AEA is not only associated with biopsy-proven NASH, but also associated with advanced liver histopathology, and that it may serve as a biomarker for NASH.

Table 3

Characteristics of the NASH cohort

Table 4

Multivariable predictors for biopsy-proven NASH

Exogenous AEA alters glucose and lipid metabolism in zebrafish. To take the first preliminary steps toward ascertaining whether AEA may play a functional role in a model organism, we administered exogenous AEA to the zebrafish Danio rerio. We examined glucose and lipid metabolism, both of which were altered in our analyses of the FHS Gen 3 cohort. We conducted dose-response experiments in larval zebrafish that were exposed to increasing concentrations of AEA, which was added to the water to provide a constant reservoir of exogenous compound. Of note, several prior reports and our own preliminary studies showed that acute exogenous AEA administration in mice results in rapid metabolism, most likely by the fatty acid amide hydrolase, most likely by the fatty acid amide hydrolase, preventing the accumulation of AEA in blood and tissues (data not shown) (23, 24). When evaluating glucose homeostasis or lipid deposition, we found that administration of exogenous AEA was associated with increased glucose and lipid deposition (see Supplemental Methods). Furthermore, coadministration with the cannabinoid receptor inhibitor rimonabant ameliorated the AEA-induced elevation in glucose or lipid, underscoring the specificity of the observations (see Supplemental Methods).

In this study, we developed a streamlined metabolomics method to investigate metabolites associated with cardiometabolic disease and NAFLD. Our analysis of a large, well-established, community-based cohort identified AEA as a candidate biomarker linking obesity and hepatic steatosis. Given the overall healthy status of the participants, the findings were subtle and represent associations many years before disease phenotypes develop. For this reason, we also studied AEA in a hospital-based, patient cohort with a histologically defined NAFLD together with carefully selected controls matched for age and BMI. Taken together, our study highlights the endocannabinoid system as a putative link among obesity, the metabolic syndrome, and liver pathology.

Our findings in relation to AEA are also in accord with prior studies of its cognate receptor. A role for the cannabinoid receptor CB1 in obesity and energy homeostasis has been established in relation to food intake and energy balance (25, 26). Activation of the CB1 receptor centrally and peripherally leads to altered energy balance and weight gain in model organisms (19, 27–29). Selective blockade of CB1 with rimonabant in human clinical trials demonstrated an effect on weight loss (30, 31), but further evaluation of cardiovascular outcomes was limited by psychiatric side effects, including an apparent increase in suicide risk (32). On the other hand, there is less data available on the direct role of AEA in metabolic disease. In accordance with our zebrafish data, rodents that are administered a high-fat diet demonstrate elevated AEA associated with weight gain and hyperglycemia (33, 34). However, the association of AEA with human obesity and glucose metabolism has been less certain (35–37) and has never been evaluated in a large human cohort free of overt disease — an important strength of the present study. The confirmation of our finding in two independent cohorts underscores the relevance of this pathway in human disease.

Our analysis also identified AEA as a candidate biomarker for NASH. Obesity and diabetes are risk factors that contribute to the increasing prevalence of hepatic steatosis and NAFLD in the US (2, 38). NASH can progress to cirrhosis or end-stage liver disease and is predicted to be the leading indication for liver transplantation by 2020 (39, 40). As only a subset of those with NAFLD progress to NASH and fibrosis, the ability to risk stratify patients and target only high-risk patients for liver biopsy is essential to limiting biopsy-related morbidity and mortality, while allowing for the accurate identification of high-risk patients (22). The finding that AEA level is associated with NASH severity, independent of traditional NASH risk factors, such as liver enzyme, HOMA-IR, HDL cholesterol, and lipid levels, suggests that this metabolite may serve as a biomarker of the presence and severity of NASH and may aid in elucidating the pathways that drive NASH progression. Although previous studies have implicated the endocannabinoid system in hepatic steatosis, either via altered regulation of the cannabinoid receptors (26) or in association with the related endocannabinoid 2-arachidonoylglycerol (21), our study is the first to our knowledge to implicate AEA specifically. The reasons for the closely related but not identical findings may be due to differences in the type of controls, matching of cases and controls, and differences in outcome measures. Future studies in a prospective cohort may help to clarify which endocannabinoid metabolites are associated with NAFLD and at which stage of severity this association begins.

Strengths of this analysis included a metabolite-profiling method with simplified preanalytical sample preparation, rapid separation, and optimized detection for human plasma metabolites. Maximizing the strengths of targeted methods, we developed a rapid approach that can routinely detect and quantify 179 analytes across key biologic pathways in human plasma. Our platform facilitates analysis in large cohorts, helping to overcome the risk for false-positive associations that can occur in multiple hypothesis testing (41). Whereas other studies have reported methods compatible with positive and negative polarity mass spectrometry (12, 42) or have developed multiplexing separation (43), our study represents both an improvement in the cycle time (~12.5 minutes) and a simplified preparation that facilitates metabolite quantification and replication in large patient cohorts.

Our study has limitations. Our biomarker analyses in the human cohorts identify associations and do not imply a causal relationship between AEA and NAFLD. Additional prospective cohort studies are needed to confirm our findings, adjust for potential confounders, and validate concentration thresholds for differentiating NASH from hepatic steatosis. Although our preliminary studies of exogenous administration of AEA in a vertebrate model demonstrated that excess AEA can lead to accumulation of lipid in the liver, additional work is needed to elucidate mechanisms and definitively demonstrate a causal relationship at more physiologically relevant levels. We also identified other candidate metabolites associated with cardiometabolic risk traits, including UDP-GlcNAc, hippuric acid, and reduced glutathione. It is possible that these biomarkers, which require further validation and study, may also predict liver disease, and future work is needed to determine the relative contributions of these other metabolites.

Taken together, our study supports the notion that metabolites such as AEA are fundamental biomarkers for metabolic disease, including obesity, insulin resistance, and NAFLD.

Subjects. The FHS Gen 3 cohort enrolled 4,095 individuals from 2002 to 2005 in a community-based longitudinal cohort study. Of 1,006 randomly selected participants, 1,000 were free of prevalent myocardial infarction or congestive heart failure at the first examination cycle and were selected for metabolite profiling. Cardiometabolic traits, routine biochemical testing, and metabolite data were available for analysis in 997 of these subjects. 470 participants from FHS Gen 3 participants who had metabolite profiling also underwent an abdominal CT scan with evaluation for the presence of liver fat, as assessed by the liver-to-phantom ratio (44).

The Massachusetts General Hospital NAFLD cohort is a hospital-based cohort of individuals who have undergone clinically indicated liver biopsies and have available liver histology (2011–2016). The initial cohort contains a spectrum of NAFLD phenotypes, including normal liver histology, steatosis, NASH, and NASH with fibrosis. Individuals with other causes of chronic liver disease or with excess alcohol consumption, defined as >21 drinks per week for men and >14 drinks per week for women were excluded from this cohort (45). For the case-control study, patients with diabetes were excluded to remove the confounding effect of diabetes. Subjects with biopsy-proven NASH were matched by age, gender, and BMI to controls that exhibited normal liver histology. Liver biopsies were evaluated by a single blinded hepatopathologist and graded according to the NASH clinical research criteria for NAFLD activity score and fibrosis stage (22). NASH was defined by the presence at least of grade 1 steatosis, lobular inflammation, and hepatocyte ballooning. Normal liver histology had no evidence of steatosis, portal or lobular inflammation, hepatocyte ballooning, or fibrosis. Plasma samples used for analysis were drawn after a 12-hour fast within 3 months prior to liver biopsy.

Plasma samples. EDTA blood samples were collected and immediately centrifuged to separate cellular material from plasma. Aliquots of plasma were frozen on dry ice and stored at –80°C until analysis. 30-μl samples were deproteinized with 70 μl acetonitrile/methanol (75:25; v/v) containing deuterated internal standards (25 μM thymine-d₄ [Sigma-Aldrich], 10 μM inosine-15N₄ [Cambridge Isotope Laboratories], 10 μM citrulline-d₇ [Sigma-Aldrich], 25 μM phenylalanine-d₈ [Cambridge Isotope Laboratories], and 10 μM valine-d₈ [Sigma-Aldrich]). For some experiments in zebrafish lysates, AEA-d₄ [Cayman Chemical] was spiked into the matrix to generate the standard curve used to quantify AEA levels. After vortexing, the samples were centrifuged at 20,000 g at 4°C for 15 minutes, and supernatants were transferred to HPLC quality glass vials with inserts (MicroSolv). Additional sample preparation solvents that were evaluated (see Supplemental Methods) included acetonitrile/methanol/formic acid (75:25:0.1; v/v/v), acetonitrile/water/ammonium hydroxide (85:14:0.7), and methanol/water (80:20).

Chromatography. Deproteinized plasma extracts were subjected to normal phase hydrophilic interaction chromatography using a 2.1 × 100 mm 3.5-μm Xbridge amide column (Waters). Mobile phase A was 95:5 (v/v) water/acetonitrile, with 20 mM ammonium acetate and 20 mM ammonium hydroxide (pH 9.5). Mobile phase B was acetonitrile. Ammonium acetate, ammonium hydroxide, and Optima grade solvents were purchased from Fisher Scientific. For amide-negative mode, the chromatography system consisted of a 1260 Infinity autosampler (Agilent) connected to a 1290 Infinity HLPC binary pump system (Agilent). The initial conditions were 0.25 ml/min of 85% mobile phase B, followed by a linear gradient to 35% mobile phase B over 6 minutes. This was followed by a linear gradient to 2% mobile phase B over 0.5 minutes held for an additional 0.5 minutes and then a 0.5-minute gradient return to 85% mobile phase B. Column equilibration was continued for 4.5 minutes at 0.5 ml/min for a total cycle time of approximately 12.5 minutes. The column compartment was maintained at 30°C.

The chromatography for positive mode utilized an HTS PAL autosampler (Leap Technologies) connected to a 1260 HPLC binary system (Agilent). The same mobile phases and column were used as in the amide-negative mode, except that initial conditions were 0.25 ml/min of 90% mobile phase B. A linear gradient to 10% mobile phase B over 6 minutes was followed by a hold for 1 minute at 10% mobile phase B. The initial conditions were restored over 0.5 minutes, and the column was equilibrated for 4.5 minutes at 0.5 ml/min. Although these methods were developed in parallel on separate LC-MS/MS instruments, the method can be run routinely on a single instrument, alternating a sample injection in the positive and then negative modes.

Mass spectrometry. Using purified reference standards, 157 metabolites were optimized for negative-mode detection on a 6490 QQQ (Agilent) triple-quadrupole mass spectrometer equipped with an electrospray ionization source. Transitions for each compound were assessed for sensitivity, selectivity, and retention time in the pooled plasma matrix with and without the spiked reference standard. 78 of the metabolites were reproducibly detected in human plasma and were selected for dynamic multiple reaction monitoring (MRM), which had a minimum dwell time of 30 milliseconds for each metabolite. Final mass spectrometry settings for the 6490 QQQ were as follows: sheath gas temperature, 400°C; sheath gas flow, 12 l/min; drying gas temperature, 290°C; drying gas flow, 15 l/min; capillary, 4,000 V; nozzle pressure, 30 psi; nozzle voltage, 500 V; and delta EMV, 200 V.

A similar strategy was used to generate a positive ionization mode MRM method, except that a QTrap 4000 (Applied Biosystems/Sciex) mass spectrometer was used. Of 114 metabolites optimized for positive mode, 101 were detectable in human plasma and included in the final MRM method. The final mass spectrometry settings for the QTrap 4000 were as follows: source temperature, 450°C; ion spray voltage, 5,000 V; CAD gas, 10; CUR gas, 20.

Peak integration, normalization, and quality control. Metabolite quantification was determined by integrating peak areas using MassHunter QQQ Quant (Agilent) or Multiquant software (version 2.0; Applied Biosystem/Sciex). All metabolite peaks were manually reviewed for peak quality in a blinded manner. In addition, pooled plasma samples were interspersed within each analytical run at standardized intervals every 10 injections, enabling the monitoring and correction for temporal drift in mass spectrometry performance (46). Each of these samples was prepared, extracted, and processed from a larger pool of normal human plasma. The nearest neighbor flanking pair of pooled plasma was used to normalize experimental samples in a metabolite-by-metabolite manner. In addition, a separate pooled plasma sample was included every 20 injections to assess the performance of the normalization procedure and determine the coefficient of variation (CV) for each metabolite. The CV of each analyte is included in Supplemental Table 1. Any sample with deuterated internal standard values ≥2 SD was excluded from peak integration and further analysis. Validation of peak identity of AEA is shown in the Supplemental Methods.

Statistics. All metabolite values were natural logarithmically transformed because of their nonnormal distribution and then standardized (to mean = 0, SD = 1) within the cohort. Age- and sex-adjusted regression analyses were performed in each study sample to examine the relation of each metabolite (predictor variable) to each clinical metabolic trait (response variables): BMI, fasting glucose, log HOMA-IR, systolic blood pressure, diastolic blood pressure, log triglycerides, and HDL cholesterol. Linear regression analyses were performed with adjustment for age, sex, and batch on the LC-MS run (see Supplemental Tables 2–9). Given that 179 metabolites were analyzed across both methods, we used a Bonferroni-corrected P value threshold of 2.79 × 10^–4 to account for the number of metabolites analyzed. Because the majority of metabolites were correlated within well-defined biological groups (amino acids, nucleotides, Krebs cycle organic acids, etc.), this correction was conservative.

Clinical and laboratory continuous variables were compared between two independent groups using 2 tailed Student’s t test or Mann-Whitney test, as appropriate. Categorical variables were compared using Fisher’s exact test. For the metabolite analysis in the animal samples, we used a P value threshold of 0.05, using Mann-Whitney or Student’s t test, depending on data normality. Statistical analyses were performed using SAS (SAS Institute) (see Supplemental Figures).

Study approval. All human study protocols were approved by the institutional review board (Massachusetts General Hospital), and all study participants provided written informed consent.