All data analyses were carried out in R version 3.4.0. Peak areas from positive ion mode, SPE Fraction 2 samples were obtained from MultiQuant (AB SCIEX, Framingham, MA, USA) software. Data for all 5 technical replicates for each condition (control and Meth) of each donor (six total) were combined to give a matrix with columns corresponding to samples and rows corresponding to metabolites. All data were log2 transformed and then normalized by subtracting the area of the internal standard abacavir in each sample. This data was visualized in a heat map presented in Figure 4 using the heatmap.2 function available in the gplots R package 17 . Hierarchical clustering was performed on the matrix rows (metabolites). The peak area values were converted to z-scores to show relative downregulation (blue) and upregulation (red) of metabolites.
For differential expression analyses, we calculated the mean of each metabolite across the technical replicates separated by condition to obtain a matrix of mean peak abundance values. The moderated paired t test within the limma R package was used to determine the differences between control and Meth groups and to identify significantly differentially expressed metabolites 18 (link)–19 . The p-values were adjusted using the Benjamini & Hochberg method, and an adjusted p-value of no more than 0.1 was considered to be significant.