Development and validation of ScreenDMTTo identify which lipid’s association with BMI is replicated, we sought a test that accounted for direction and had the strongest adjusted p-values in terms of the false discovery rate (FDR). Because we cared about directionality it was natural to consider DMT. By comparison, it was recently shown that the MaxP test is the likelihood ratio test (LRT) for the non-directional mediation null hypothesis17, which is mathematically equivalent to testing replication without regard for direction in two studies. The LRT has many beneficial properties such as often being optimal, which is why LRTs so often used in practice, e.g. the t-test for means. For directional mediation and replication, we mathematically prove that DMT is the LRT in Supplementary Note 1. To our knowledge, analyses of replication and mediation had not been previously connected to analysis of qualitative interactions, which occurs when an effect changes sign between groups of subjects, e.g. a treatment benefits one group of patients while it harms another. Surprisingly, we find that directional replication and mediation are connected to qualitative interaction, because applying the DMT in search of effects that have opposite signs between groups is the LRT for qualitative interactions18,19.We considered how to improve the adjusted p-values from our directional MaxP test with screening or filtering. The filtering approach in AdaFilter when applied to two studies uses the minimum p-value (pmin) per analyte as the filtering p-value and the maximum p-value as the selection p-value. We cannot directly use AdaFilter by naively setting the DMT p-value (pDMT) as the selection p-value, since the selection p-value must exceed the filtering p-value, which would not always be true. Instead, we developed a way to apply screening to our directional MaxP test, which we term “ScreenDMT”. ScreenDMT defines its selection p-value per analyte as the maximum of pDMT and pmin. We mathematically prove ScreenDMT to be valid for controlling the FDR and the family-wise error rate (FWER; the rate controlled by the Bonferonni method) in Supplementary Note 1 for large sample sizes. For small sample sizes, the screening-based multiple testing adjustment is very slightly too strong, but this is where the DMT p-value, like similar statistical tests, is weak anyway. Like in AdaFilter14, theoretical control of the FWER requires independence of p-values within each study, while control of the FDR for a large number of analytes allows weak dependence between analytes. Such weak dependence typically holds, for example, between gene loci and expression14. Importantly, not all FDR methods theoretically allow dependence, such as RepFdr16.We compare ScreenDMT’s FDR and power via simulation against competing directional replication methods, including the directional approaches of AdaFilter14, radjust-sym15, RepFdr16, and DMT (Fig. 1). We also compare against a directional alternative to the MaxP test that was adapted from a meta-analysis due to Pearson14,20, which applies MaxP for each analyte to right-sided p-values of both studies, then applies MaxP to left-sided p-values, and doubles the minimum from the two applications. This procedure also underlies AdaFilter’s directional approach.Fig. 1: Comparison of ScreenDMT’s power and false discovery rate to competing approaches in replication of two datasets via simulation.Simulation of two datasets with strong or weak dependence between analytes when the signal strength or effect size per analyte is equal between the two datasets and when it is unequal. The x-axis represents the percent of analytes that do not truly directionally replicate but show some effect in one of the studies or show effect in both studies of opposite direction. a Probability of detecting true associations (power). b False discovery rate. The FDR threshold is 5%.We see in Fig. 1 that RepFdr often does not control its FDR below the nominal level, and apart from RepFdr, ScreenDMT retains the most power across the scenarios. We also see that DMT has more power than the approach adapted from Pearson14,20, which helps explain ScreenDMT’s power advantage over AdaFilter. We perform the analogous simulation for FWER using methods that account for directionality and provide FWER in Supplementary Fig. 1. We see that all methods properly control FWER, DMT is slightly more powerful than the method adapted from Pearson14,20, and ScreenDMT’s power is very similar but slightly stronger than AdaFilter’s.Inverse association of 12,13-diHOME with BMI replicatesWe wanted to test if 12,13-diHOME negatively associates with BMI in the Study 2 cohort of 83 individuals (Table 1). Our previous association of 12,13-diHOME to BMI in Study 1 was assessed using Spearman rank correlation, so we also use a nonparametric test in Study 2 (Figs. 2a, b). Since Study 2 is split into obese and non-obese subjects, we evaluate whether lipids are differentially abundant between groups with a non-parametric t-test (i.e. Wilcoxon rank sum or Mann-Whitney U test) (Figs. 2c, d). We apply this test in a one-sided manner to test if 12,13-diHOME is lower in obese subjects. Since 12,13-diHOME had already been prioritized, its p value does not need to be adjusted for multiple testing. We obtain a one-sided p-value of 0.00093 against a significance p value threshold of 0.05, which demonstrates replication (Fig. 2d).Table 1 Anthropometrics of Study 2 cohortFig. 2: 9,10-diHOME and 12,13-diHOME versus BMI in both cohorts.Plasma levels of 9,10-diHOME and 12,13-diHOME versus BMI in participants from the original cohort (Study 1) and from the new cohort (Study 2). To improve visualization, 9,10-diHOME values above its 99th percentile were winsorized to its 99th percentile. The y-axis represents log2-transformed area ratios, which are in arbitrary units (A.U.). Spearman’s rank correlation coefficient or Wilcoxon rank sum’s z-score are provided. a 12,13-diHOME in Study 1. b 12,13-diHOME in Study 2. c 9,10-diHOME in Study 1. d 9,10-diHOME in Study 2.Application of ScreenDMT to test replication of other lipids’ association with BMIIn addition to 12,13-diHOME, we measured a panel of more than a hundred oxidized lipid species, including metabolites of “parent” lipids: linoleic acid (LA), α-linolenic acid, arachidonic acid, dihomo-γ-linolenic acid, Docosahexaenoic acid and Eicosapentaenoic acid (Fig. 3a). To assess if any lipids other than 12,13-diHOME have an association with BMI that replicates, we applied the non-parametric t-test to all the lipids except 12,13-diHOME and then assessed replication with ScreenDMT (Supplementary Data 1). The only other lipid that replicated was 9,10-dihydroxy-12Z-octadecenoic acid (9,10-diHOME), which is a regioisomer of 12,13-diHOME (Fig. 3b), whose FDR is 5%. Whereas if we had tested replication with DMT, the FDR of 9,10-diHOME would be nearly double at 9%. We show the plots of 9,10-diHOME vs. BMI and of 12,13-diHOME vs. BMI in both cohorts in Fig. 2.Fig. 3: Biosynthesis of linoleic acid diols 9,10-diHOME.a Lipids measured in our lipidomic panel are shown as a product of their precursor fatty acids and the oxidative pathways that are the first step in their biosynthesis. b Biosynthesis of linoleic acid diols.We next assessed lipid sets for enrichment of lipids whose nominal DMT replication p-value was below 5%. We defined one lipid set per parent lipid from the left-most column of Fig. 3a, so we created 6 lipid sets. For example, the LA set contained the eight lipids in the first row of Fig. 3a. We tested if each set was over-represented among the nominally significant lipids using a one-sided Fisher exact test and adjusted the six p-values for multiple hypothesis testing using the Benjamini-Hochberg FDR (Supplementary Data 1). We found that the LA set was the only significant set (FDR < 10−5) and that all of the LA-derived lipids had an inverse (but not necessarily significant) relationship with BMI in both studies.9,10-diHOME and 12,13-diHOME activate adipocyte calcium fluxes12,13-diHOME and its stereoisomer 9,10-diHOME are replicably associated with BMI, and we had previously shown that 12,13-diHOME increases cellular fatty acid uptake by inducing the translocation of fatty acid transporters FATP1 and CD36 to the cell surface; however, the signaling pathway that mediates this effect is unknown. 12,13-diHOME and 9,10-diHOME have both been found to activate ion flux in primary neurons and CHO cells21, so we measured calcium flux in adipocytes treated with each diol. We used in vitro differentiated WT1 mouse brown preadipocytes and in vitro differentiated 3T3-L1 mouse white preadipocytes to model brown and white adipocytes, respectively. Both cell lines were differentiated using a standard adipogenic induction media, and then loaded with Fura-2. After we recorded baseline calcium flux, cells were treated with Hanks Buffered Saline Solution (HBSS) containing 1% fatty acid-free bovine serum albumin (BSA) alone or containing different concentrations of linoleic acid diol, and calcium flux was monitored every 50 seconds for approximately 9 minutes in a plate reader (Fig. 4). Both 9,10-diHOME and 12,13-diHOME triggered calcium influx in brown and white adipocytes, suggesting the signaling pathway activated downstream of these two lipids overlaps.Fig. 4: Linoleic acid diols activate calcium flux in cultured adipocytes.X-axis is timepoint after baseline (baseline is timepoint zero). Each timepoint is approximately 50 seconds apart. Immediately after baseline, vehicle control or different concentrations of linoleic acid diol were injected. Y-axis is FLUOFORTE fluorescence quantification in arbitrary units (A.U.). To normalize, values were divided by their mean baseline value. Each data point is from an independent well; n = 3 independent wells per group. a Differentiated WT1 brown adipocytes treated with 12,13-diHOME or vehicle control. b Differentiated WT1 brown adipocytes treated with 9,10-diHOME or vehicle control. c Differentiated 3T3-L1 white adipocytes treated with 12,13-diHOME or vehicle control. d Differentiated 3T3-L1 white adipocytes treated with 9,10-diHOME or vehicle control. For each panel, normalized post-baseline fluorescence was regressed on continuous timepoint and categorical dose; each dose was compared to vehicle control with a two-sided t-test; and p-values from all panels were Bonferroni corrected. ***Bonferroni-corrected P < 0.001 vs. vehicle control. Statistical results are in Supplementary Data 1.
