Integrated lipidomics and transcriptomics characterization upon aging

Subscriber access provided by WEBSTER UNIV

Article

Integrated lipidomics and transcriptomics characterization upon aging related changes of lipid species and pathways in human bone marrow mesenchymal stem cells Xin Lu, Yue Chen, Huiyu Wang, Yunfan Bai, Jianxiang Zhao, Xiaohan Zhang, Li Liang, Yang Chen, Chenfei Ye, Yiqun Li, Yi Zhang, Yu Li, and Ting Ma J. Proteome Res., Just Accepted Manuscript • DOI: 10.1021/acs.jproteome.8b00936 • Publication Date (Web): 04 Mar 2019 Downloaded from http://pubs.acs.org on March 6, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Proteome Research

Integrated lipidomics and transcriptomics characterization upon aging related changes of lipid species and pathways in human bone marrow mesenchymal stem cells Xin Lu1,‡, Yue Chen2,‡, Huiyu Wang3, Yunfan Bai2, Jianxiang Zhao2, Xiaohan Zhang2, Li Liang1, Yang Chen1, Chenfei Ye1, Yiqun Li2, Yi Zhang4, Yu Li2,* and Ting Ma1,* 1. School of Electronic and Information Engineering, Harbin Institute of Technology, Shenzhen, Shenzhen, Guangdong, 518000, China 2. School of Life Science and Technology, Harbin Institute of Technology, Harbin, Heilongjiang, 150080, China 3. School of Pharmacy, Qiqihar Medical University, Qiqihar, Heilongjiang, 161000, China 4. Tian Qing Stem Cell Co. Ltd., Harbin, Heilongjiang, 150080, China

ACS Paragon Plus Environment

1

Journal of Proteome Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 40

CORRESPONDING AUTHOR(S): Ting Ma School of Electronic and Information Engineering, Harbin Institute of Technology, Shenzhen HIT campus of University Town of Shenzhen, Shenzhen, 518055, China Fax/Tel: (+86)-755-26033608 E-mail: [email protected] Yu Li School of Life Science of Technology, Harbin Institute of Technology No. 2, Yikuang Street, Harbin 150081, China Fax/Tel: (+86)-451-86402690 E-mail: [email protected]


2

Page 3 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ABSTRACT

Aberrant differentiations of bone mesenchymal stem cells (BMSCs) have proved to be associated with the occurrence of senile osteoporosis. However, mechanisms of this phenomenon relative to abnormal lipid metabolism remain unclear. This study was conducted to characterize the lipidomics alterations during BMSC passaging, aiming at uncovering the aging-related lipid metabolism that may play an important role in aberrant differentiations of BMSCs. Principal component analysis presented the sequential lipidomics alterations during BMSCs passaging. Majority of glycerophospholipids including phosphatidylcholines, phosphatidylethanolamines, phosphatidylglycerols as well as sphingolipids revealed significant elevations, whereas the others including phosphatidic acids, phosphatidylinositols and phosphatidylserines presented decreases in aged cells. Double bond equivalent versus carbon number plots demonstrated that the changing trends and significances of lipids during passaging were associated with the chain length and degree of unsaturation. In the correlation networks, the scattering patterns of lipid categories suggested the category-related metabolic independences and potential conversions among phosphatidic acids, phosphatidylinositols and phosphatidylserines. The lipid-enzyme integrated pathway analysis indicated the activated metabolic conversions from phosphatidic acids to CDP-diacylglycerol to phosphatidylinositols and from sphingosine to ceramides to sphingomyelins with BMSC passaging. The conversions among lipid species described the lipidomics responses that were potentially the inducement for the aberrant differentiations during BMSC aging.

KEYWORDS lipidomics, transcriptomics, BMSC, aging, glycerophospholipids, sphingolipids


3


Page 4 of 40

INTRODUCTION Mesenchymal stem cells (MSCs), a member of multipotent stem cells with great powers in multipotential differentiation, self-renewal, immunomodulation and hematopoiesis, are formed in various tissues including bone marrow, adipose tissue, umbilical cord blood and so on.1 In particular, bone mesenchymal stem cells (BMSCs) also termed as bone marrow stromal cells are mainly isolated from bone marrow which could differentiate into sorted of functional somatic cells such as osteoblasts, cartilage, adipocytes, myocardial cells and fibroblasts, making it play important roles in the formation and regeneration of bones.2 More importantly, the excessive differentiation of BMSCs into adipose cells rather than osteoblasts, which frequently occurs under exogenous and endogenous factors such as aging, menopause and abuse of hormones, tends to induce the imbalanced bone metabolism resulting in the loss of bone mass and even osteoporosis (OP).3 Therefore, the overall understanding to the cellular functionalities and metabolic alterations of BMSCs along with aging will be of great significance in mechanism exploration and clinical therapy of OP. As a self-organized subset of metabolome, cellular lipids function in variety of biological events

involving

energy

storage,

membrane

construction,

cellular

signaling

and

immunoregulation.4 Moreover, the associations between the disordered lipid metabolism and fat deposition in organs were observed in various diseases such as obesity and hepatic steatosis that demonstrated the active role of altered lipidome played in the abnormal accumulation of adipose during development of chronic metabolic diseases like senile osteoporosis, of which the mechanism remains not fully understood.5,6 In spite of boomingly developed for years, lipidomics, a new metabolomics branch aiming at globally profiling perturbations of various lipid species in organism, was still not fully concerned


4

Page 5 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and served for the investigation on the mechanism of osteoporosis. And even so, majority of the enclosed studies were carried out based on clinical samples such as bone tissues that provided valuable evidence in characterization of lipids significantly varied in patients with osteoporosis but were still insufficient to dynamically illustrate the in-depth mechanism.7,8 Notably, by using a shot-gun lipidomics strategy, Tigistu-Sahle et al. characterized the alterations of membrane glycerophospholipids (GPs), fatty acids (FAs) and mRNA expression of related enzymes in BMSCs with different expansion status and those from donators of different ages to investigate the association between membrane GPs and immunomodulation functionality during BMSCs aging.9 The results showed major lipidomics changes in GP and FA profiles during expansion and senescence, demonstrating the close correlation of the immunological functionality of BMSCs with the changes of membrane lipid compositions. Nevertheless, the specifically acquired GPs and FAs might not uncover the overall fluctuation of lipid species tightly connected and the disturbed metabolic conversions that play key roles during aging. Also, due to severe ion suppression, shot-gun based lipidomics analysis is still inadequate to effectively cover lipid species with low abundances in a complex biological system.10 As thus, to present a more comprehensive characterization upon lipidomics variations and their latent interactions during BMSC aging, a non-targeted lipidomics study was performed based on ultra-performance liquid chromatography coupled to mass spectrometry (UPLC-MS), which is able to reduce the complexity of the matrix and enhance the sensitivity.11 Multivariate analyses involving principal component analysis (PCA) and partial least square discriminant analysis (PLS-DA) were both conducted to reveal the global clustering layout as BMSCs aging until senescence. Univariate analysis was then conducted to evaluate the differential significance of lipids among passaged BMSCs and the correlations of lipid alterations with passages. To observe the diversity of


5


Page 6 of 40

alterations associated with the chain length and degree of unsaturation (DOU) of lipid molecules, double bond equivalence versus carbon number plots were established based on lipids of each category. To characterize the conversions between different lipid species biochemically associated, pathway analyses were conducted via the integration of alterations of lipid species and enzymes of which mRNA expressions in passaging BMSCs were measured. EXPERIMENTAL METHODS Ethics Statement This study was approved by the ethical committee at Harbin Institute of Technology and conducted in accordance with the 1964 Declaration of Helsinki and its later amendments or comparable ethical standards. Reagents and Chemicals Methanol and acetonitrile were purchased from Fisher Scientific (Waltham, MA, USA); isopropanol and methyl tert-butyl ether (MTBE) were provided by MREDA (Beijing, China); formic acid and ammonium hydroxide were obtained from Fluka (St. Louis, MO, USA). Deionized water was produced by a Milli-Q ultrapure water system (Millipore, Billerica, MA, USA). Senescence β-Galactosidase Staining Kit was purchased from Cell Signaling Technology (Danvers, MA, USA). All the reagents applied in this study were HPLC grade. Cell Culture Two strains of human BMSCs (hBMSCs) were purchased from Cyagen Biosciences (Suzhou, China). According to the manufacturer’s instructions, hBMSCs were harvested from the bone marrow of two healthy male donors of 37 and 45 years old. Both strains were purified and verified via measurement of surface antigens CD73, CD90, CD105, CD34, CD45, CD14, CD79a and HLA-DR by using Canto-2 flow cytometer (BD Biosciences, San Jose, CA, USA) (Figure


6

Page 7 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


S1). The cells were cultured in Dulbecco’s Modified Eagle Medium (Gibco, Carlsbad, NM, USA) supplemented with 10% fetal bovine serum, 1% L-glutamine and 1% penicillinstreptomycin solution (Invitrogen, Carlsbad, NM, USA). The medium was renewed every three days, and the cells were passaged while the density reached to ~80%. Oil Red O Staining Oil red O staining was conducted with the Oil Red Staining kit (Cyagen Biosciences) based on the manufacturer’s instruction. Briefly, the BMSCs in the plates were washed with phosphate buffer saline (PBS, 1×) for 2 times after removal of the medium. Then, each plate was added by 2 mL of 4% formaldehyde and incubated at room temperature (RT) for 30 min. After washing with 1×PBS, the cells were added by Oil red O solution and incubated for 30 min. The cells were then washed with PBS for 2 times followed by being checked under a microscope. Alizarin Red S Staining The procedure for Alizarin Red S staining was very similar to that for Oil red O staining. In short, the cells were washed with 1×PBS and fixed with 2 mL of 4% formaldehyde for 30 min. Differently, 1 mL of Alizarin Red solution (Cyagen Biosciences) was added into each plate and incubated for 3 min. After removing the Alizarin Red solution, the cells were washed with PBS for another 2 times and checked under a microscope. Senescence β-Galactosidase Staining Senescence β-galactosidase staining was performed to validate the senescence status of BMSCs at old passage. The staining was carried out by using Senescence β-Galactosidase Staining Kit (Cell Signaling Technology, Danvers, MA, USA) based on manufacturer’s instruction. In brief, BMSCs in the plates were washed with PBS for 3 times after removal of media. 1 mL of Fixative Solution (1×) was further added into each plate followed by fixing for 10 min at RT. Then, the


7


Page 8 of 40

plates were washed with PBS for 2 times and added with 1 mL of the β-Galactosidase Staining Solution followed by incubation at 37 °C overnight in a CO2-free incubator. Lastly, the cells were checked under a microscope before the β-Galactosidase dried. Lipid Extraction from Cells The lipid extraction was performed according to the method previously developed with slight modification. In brief, cells reaching density of 90% were washed with PBS for 3 times. After the removal of the residual PBS, the cells were quenched by adding 1 mL of liquid nitrogen followed by being harvested via scraping with 200 L of ice-cooled 75% methanol-water solution. The given suspension was then cracked by sonication at 30 Hz for 1 min in ice-water bath. Next, 500 L of MTBE was added into the homogenate for lipid extraction. The mixture was then oscillated at 1000 rpm, 25 °C for 30 min and added by 125 μL deionized water for stratification. After centrifuged at 12000 g for 15 min, 400 μL of the upper organic layer (lipids mainly dwelling) was collected and dried using rotating-vacuum. The obtained residual was further stored at -80 °C until analysis. UPLC-MS Analysis The chromatographic separation was performed upon an ultrahigh performance liquid chromatography system (Shimadzu, Kyoto, Japan) with a column of HSS T3 2.1 mm × 100 mm, 1.8 μm (Waters, Milford, MA, USA). The mobile phase A consisted of isopropanol/methanol (9/1, v/v) with 10 mM ammonium formate and mobile phase B was 60% acetonitrile-water solution with 10 mM ammonium formate. Prior to LC-MS analysis, the lipid residual was redissolved by 100 μL of 50/50 (v/v) isopropanol-methanol solution and 5 μL of the supernatant after centrifugation was injected. For the lipid separation, a linear gradient was run at flow rate of 0.3 mL/min and compiled as follow: holding 10% A from initial to 1.0 min; 10-80% A from 1.0


8

Page 9 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


to 15.0 min; 80-95% A from 15.0 to 20.0 min followed by isocratic elution at 95% A for 4 min; changing the gradient back to 10% A from 24.0 to 25.0 min followed by an equilibration for 2 min. During the analyses, the column temperature was maintained at 40 °C. Data acquisitions including full scan and MS/MS scan were performed using a TripleTOF 5600 plus series quadrupole time-of-flight mass spectrometer equipped with a dual electrospray ion source (AB Sciex, Redwood City, CA, USA) in both positive and negative mode. MS conditions were set as follows: spray voltages were set at 5.5 kV for positive mode and 5.0 kV for negative; both GS1 and GS2 were set as 55 psi while TEM was 550 °C. Mass scan was performed ranging from m/z 300 to 1500 at a rate of 1.5 spectra/s. For MS/MS analysis, the collision energy (CE) used was ranged from 20 to 70 eV as a function of molecular weight (MW). Data Processing Data pre-processing including peak detection, noise filtering, feature alignment and data normalization were all performed by using the off-line xcms package integrated in R-platform. The parameters were adopted with slight adjustment depending on the raw data. In detail, “matchedFilter” was applied as algorithm method for feature detection; the peak width was set as 5 to 30 s; band width was set at 6 s. The other parameters were listed in Table S1. Identification of Lipids Potential lipid features were firstly selected via matching the experimental m/z values with the theoretical m/z values obtained from the LIPIDMAPS structure database (LMSD, http://www.lipidmaps.org/data/structure/index.html) based on the cutoff of smaller than 5 ppm. Furthermore, MS/MS spectra were applied for the identification of lipid species according to their specific collision patterns.12


9


Page 10 of 40

RNA-Seq Total RNA extraction was performed by using RNA extraction kit (Trizol reagent, Invitrogen, Carlsbad, USA), followed by concentration measurement with Agilent 2100 Bioanalyzer and Agilent RNA 6000 Nano Kit (Agilent, Santa Clare, CA, USA). During the library preparation, the purification of the poly-A containing mRNAs was firstly performed using poly-T oligoattached magnetic beads. Following purification, the mRNA was fragmented into small pieces, which were further used for synthesizing the first strand cDNA with reverse transcriptase and random primers. Synthesis of the second strand cDNA was then conducted by using DNA polymerase I and RNase H. After addition of a single “A” base and the ligation of the adapter, the fragmented cDNA products were purified and amplified via PCR. The yielded PCR products were then quantified and collected into a single strand DNA (ssDNA) circle, giving the final library. Subsequently, DNA nanoballs (DNBs) were produced via rolling circle replication (RCR) based on the ssDNA circle, followed by loaded into the patterned nanoarrays. The single-end read was conducted by 50 bp on the BGISEQ-500 platform based on the combinatorial probeanchor synthesis (cPAS) sequencing method.13 Quantitative Real Time PCR Reverse transcription was accomplished with kit based on the official instructions (PrimeScript RT Reagent kit, Takara, Japan). Furthermore, the qRT-PCR quantification was performed in triplicates using SYBR Green Master mix kit (Roche, Indianapolis, USA) on an ABI Prism 7900 HT sequence detection system (Applied Biosystems, Foster City, USA). Primer sequences of the genes anchored in pathway analysis were listed in Table S2. Pathway Analysis


10

Page 11 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


For the visualization of gradually changed reactions, alterations of lipid species and mRNA expressions of enzymes in glycerophospholipid and sphingolipid metabolism were imported into Cytoscape software attached with KEGG XML (KGML) pathway mapping files.14 The significantly changed genes in the RNA-seq dataset were determined based on the criteria including: fold-change values of p13- to p7-cells (FC > 1.5 and FC < 0.7)15 and adjusted pvalues in Student’s t-tests (p < 0.05, Benjamini-Hochberg method). The expressions of anchored enzymes were further validated by using qRT-PCR. Statistical Analysis PCA was conducted upon the lipidomics dataset to reveal the global alterations of BMSCs along with passaging. Hierarchy clustering analysis (HCA) was applied to evaluate the similarity of lipidomics patterns between two-sourced BMSCs. The differential significances and altering trends of lipid species were featured with Student’s t-test (p < 0.05) and Pearson correlation analysis, respectively. To investigate the altered significance specific to the chain length and degree of unsaturation, DVCs were established based on lipids of each category. To investigate the latent affections between lipid species, correlation networks were constructed in GPs and SPs based on Pearson correlation coefficients. To explore the metabolic fluctuations during BMSC passaging, pathway analysis was carried out based on the integration of lipidomics and transcriptomics dataset. The expression of GAPDH was used to do the data normalization of anchored genes in qRT-PCR analysis. The relative expressions of anchored genes were further calculated based on the 2-ΔΔCt method. PCA were performed with SIMCA-p v. 11.5 (Umetrics AB, Umea, Sweden). HCA-heatmap was conducted using MetaboAnalyst (http://www.metaboanalyst.ca/MetaboAnalyst/). DVCs were visualized by using OriginPro 8 software. Correlation network construction was realized


11


Page 12 of 40

with Cytoscape v.3.1.0 (www.cytoscape.org). Calculations of correlation coefficients and Student’s t test were conducted on the R platform v.2.14.2. Other calculations and statistical methods were completed with R platform v.2.14.2. RESULTS Validation and Morphological Alterations of BMSCs Prior to the passaging inducement, the osteogenic and adipogenic differentiation potencies of BMSCs from two donors with matched characteristics (D1 and D2) at early generation were validated via Alizarin red staining and oil red staining, respectively. BMSCs from both sources could be effectively stained by both validation methods as typically presented in Figure 1A and 1B (D1), indicating their initial pluripotency at a younger generation. Gradually, as the development of passaging inducement, both sources of BMSCs revealed significant variations in cell morphology. As shown in Figure 1C and 1D, the p13-BMSCs from D1 presented significant lower cell density than that of p7 by the same initial quantity during the same culturing time, suggesting the slower proliferation of the former than the latter cells. Meanwhile, the cells became morphologically longer as they went elderly (Figure 1E and 1F) that was in accordance with the increased proportion of fibroblast-like BMSCs observed in late passage (p13). More importantly, the p13-BMSCs were significantly green-stained by the βgalactosidase (Figure 1H) compared with p7 cells (Figure 1G). These results presented a more typical morphology of senescence in aged BMSCs.


12

Page 13 of 40


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


13


Page 14 of 40

Figure 1. Validation of the differentiation potentials of BMSCs at early passage and the morphological alterations of BMSCs between early and old passages. (A) Alizarin red staining of p7-BMSCs. (B) Oil red staining of p7-BMSCs. (C) Morphology of BMSCs at p7 (bar=25 m). (D) Morphology of BMSCs at p13 (bar=25 μm). (E) Morphology of BMSCs at p7 with regional magnification (10×). (F) Morphology of BMSCs at p13 with regional magnification (10×). (G) β-Galactosidase staining of p7-BMSCs (bar=12.5 μm). (H) β-Galactosidase staining of p13BMSCs (bar=12.5 μm). Global Profiles of Lipid Species during BMSC Passaging As shown in Table S3, a total of 295 individual lipid molecules were identified, including 204 glycerophospholipids (GPs), 26 glycerolipids (GLs) and 65 sphingolipids (SPs). The typical MS/MS spectra of a portion of lipids were listed in Figure S2. To globally observe the developing patterns of lipids as BMSC passaging, PCA model was established based on all the lipids identified. As shown in the PCA score plot (Figure 2A), the BMSCs from both donors (D1 and D2) at different passages were located in different quadrants and sequentially scattered along with principal component 1 (PC1) which devoted the largest grouping potential in the model demonstrating the significant successional derivations of lipids during BMSC aging. Meanwhile, the plot also displayed a clear separation between the cells of “both ends of age” (p7 and p13) and those of “medium ages” (p9 and p11) by principal component 2 (PC2), demonstrating that some specific lipid species might present call-back-like tendencies along with BMSC aging.


14

Page 15 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 2. The PCA loading-bi plots constructed with various lipid species in BMSCs of sequential passages. (A) Scattering pattern of BMSCs of D1 and D2 at p7, p9, p11 and p13. (B) Scatting patterns of GPs and GLs in the loading-bi plot. (C) Scattering pattern of SPs in the loading-bi plot At the meantime, the high correlations between diverse lipid species and categories with principal components were also observed. As shown in Figure 2B, majority of GPs, including phosphatidylcholines (PCs), phosphatidylethanolamines (PEs), phosphatidylglycerols (PGs) and a certain amount of GLs including diacylglycerols (DGs) and a part of triacylglycerols (TGs) presented general positive correlations with PC1 whereas phosphatidic acids (PAs), most of the phosphatidylserines (PSs), phosphatidylinositols (PIs) and remained TGs negatively correlated


15


Page 16 of 40

with PC1. Monoacylglycerols (MGs) were the only ones that generally gathered in the centre of the plot, suggesting their non-significance in grouping contribution. Similar to GPs, nearly all categories of SPs showed dramatically positive correlations with PC1 (Figure 2C). The varying trends of lipid species and categories were presented in Figure S3, in which (Ly)PCs, (Ly)PEs and (Ly)PGs revealed global increases during BMSC aging in both cell strains while PAs reversely changed (Figure S3A). For PSs and PIs, “hollow”-like patterns could be observed in both plots between observations of p7 (youngest) and p13 (oldest), which were formed by the decreases from p7 to p9 and increases from p11 to p13, reflecting the call-backs and the decreases in these lipids as BMSCs got elderly. Notably, except globosides which revealed an obvious “upheaval”-like call-back from p7- to p13-BMSCs, aging-associated elevations could be globally seen in majority of SP categories including sphingomyelins (SMs), ceramides (Cers), cerebrosides and gangliosides (Figure S3B). Especially, some certain SPs, i.e. ceramides (Cers) and cerebrosides showed the most significant increases than other lipids. For GLs, no obvious alterations could be seen in MGs while DGs were generally up-regulated in the process of passaging (Figure S3C). Interestingly, due to the “bipolar” distribution scattered in the loading plot, some TGs did show significant increases during cell aging whereas other decreased. To sum up, the alterations of lipids revealed in trend plots were preferably matched with the correlations of lipid species with principle components presented in score plot, which demonstrated that the changes of lipids during BMSC aging were not monotonous but diverse. Compared with the alterations between BMSCs of different passages, the cells from the same passage but different donors were closely located in the same quadrant that suggested the similar lipidomics patterns of BMSCs from the two donors during passaging inducement. As shown in the score plot (Figure 2A), the difference between two sourced BMSCs at the same passage


16

Page 17 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


became progressively smaller with the increases in passages. This behaviour was consistent with the increasing similarities and decreasing Euclidean distances between the paired observations displayed in heatmap plot (Figure S4), which to some extent suggested the potential of cell passaging in reducing individual differences of BMSCs under the same culturing conditions. Altering Tendencies and Significances of Lipids Associated with Cell Passaging To more specifically characterize the differential tendencies and significances of lipid species between BMSCs of early and late passage, Pearson correlation analysis and Student’s t-test were conducted upon each lipid specie and category. First of all, due to the high similarity of lipidomics profiles between the two donors, correlation and univariate analyses were carried out based on the combined dataset of both sourced BMSCs to select the characteristics with less individual differences. The Student’s ttests were performed between the youngest (p7) and oldest (p13) cells to evaluate the differential significances of lipids after long-term passaging. As shown in Figure 3A, majority of the lipid categories except MGs revealed overall significant changes (p < 0.05) between p7 and p13 observations. In detail, GPs including 21 of 22 LyPCs, 29 of 38 PCs, 14 of 16 LyPEs, 41 of 55 PEs, 9 of 15 PSs, 7 of 7 PAs, 18 of 26 PIs and 14 of 21 PGs, GLs including 8 of 10 DGs and 9 of 14 TGs, SPs including 21 of 23 SMs, 10 of 10 Cers, 4 of 4 cerebrosides, 10 of 11 globosides and 5 of 8 gangliosides were significantly varied between early (p7) and late passages (p13) (Table S3). According to the fold-change plot and correlation coefficients (Figure 3B and 3C), PAs, PIs and PSs were the only ones presenting obvious decreases in p13-BMSCs while nearly all the others revealed elevations. Specifically, among the GPs with significant alterations (p < 0.05), all of the 7 PAs, 15 of 18 PIs (p < 0.05) and 6 of 9 PSs showed decreased levels in aged BMSCs


17


Page 18 of 40

while all of the 21 LyPCs, 24 of 29 PCs, 12 of 14 LyPEs, all of the 41 PEs and 13 of 14 PGs were up-regulated (Table S3). Cerebrosides, Cers and PAs revealed the most significant alterations between p7- and p13BMSCs and the highest correlations with cell passages. Notably, Cers and PAs were both the simplest lipid species as well as the precursors for the biosynthesis of other SPs and GPs that their remarkable alterations might be the starts of the perturbed metabolic flows. Being the most abundant lipids, (Ly)PCs, (Ly)PEs and SMs revealed the second most significant alterations in both fold changes and correlation coefficients. Diacyl-GPs (PCs and PEs) displayed less significant increases than their monoacyl-counterparts (LyPCs and LyPEs). In line with the patterns shown in PCA, TGs revealed a binary-directional metabolic manner during cell passaging that six TGs increased while seven decreased in p13 cells (Table S3), suggesting the functionalities of TGs varied probably depending on the structural diversity of the acyl-chains.


18

Page 19 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3. Univariate characterizations of alterations of lipid species between BMSCs at p7 and p13. (A) –log(p) values of Student’s t-test of lipid species between BMSCs at p7 and p13. (B)


19


Page 20 of 40

Fold-change values of lipid species in p13-BMSCs relative to p7-BMSCs. (C) Pearson correlation coefficients of lipid species with passaging In addition, to explore the call-back tendencies observed in PCA, the Mann-Whitney rank sum tests were launched between observations of the “both ends of age” (p7+p13) and those of “medium ages” (p9+p11). As shown in Figure S5A, it could be seen that the lipid species that revealed call-backs in PCA, i.e. PSs, PIs and Globosides, presented the most significant variations between the two clusters. Besides, both PSs and PIs showed relatively high foldchanges while Globosides possessed lowest values that were also in accordance with the “hollow” and the “upheaval” pattern in PCA trend plots (Figure S5B). Therefore, in addition to the general decreased patterns of PSs and PIs as well as the increases of Globosides observed previously, the obvious call-backs were also confirmed in some of them which needed to be further classified and characterized. Lipidomics Alterations Specific to Chain Length and Degree of Unsaturation during BMSC Passaging To characterize the lipidomics alterations specific to lipid structures in each category, carbon numbers (CNs) and double bond equivalence values (DBEs) of hydrophobic chains in each lipid molecule were calculated and imported into the DBE versus CN plots (DVCs), in which foldchange (p13/p7) value was used to visualize the altered degree of each individual lipid during BMSC passaging. As shown in Figure 4, several monoacyl-GPs in PCs, PEs, PGs and PIs (commonly CN < 24) showed more significant fold-change (FC) increases than corresponding diacyl-GPs (or 1-alkyl2-acyl-GPs). Meanwhile, the most abundant GPs, i.e. PCs and PEs both revealed diverse foldchange performances in different located areas of the DVC plots. Specifically, the diacyl-PCs


20

Page 21 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and diacyl-PEs locating in the bottom-left regions revealed dramatically lower FC levels than those in the up-right regions (Figure S6), in another word, diacyl-GPs with longer chain lengths and more double bonds presented more significant increases than those with lower CN and DBE values. However, this manner was not clearly observed in the plasmalogens (PC-Os and PE-Os). For monoacyl-GPs, most of the monounsaturated LyPCs revealed less significant increases than saturated and polyunsaturated ones. In contrast, LyPEs revealed a clear CN-associated changing pattern that the ones with CN from 16 to 18 showed significantly higher FC levels than those from 20 to 24.


21


Page 22 of 40

Figure 4. DVC plots showing fold-change alterations of categories of lipid species in p13BMSCs relative to p7-BMSCs For other categories, PIs with DBEs of 1 and 4 as well as SMs with CN of 36 and DBE of 0 showed significant higher FCs than the other ones in corresponding DVC plots. Similar to diacyl-PCs, PSs also showed general increasing FC values as DBE increased. For Cers with DBE


22

Page 23 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of 1, an FC increase-decrease manner was observed from CN of 33 to 42. Especially, TGs with CN of 54 as well as those with DBE of 2 presented much lower FC levels than other ones. As shown in Figure S7, the TGs positioning in these particular regions (CN=54/DBE=2) showed FC levels globally lower than 1 whereas those of other TGs were larger than 1, which were just corresponding to the decreases and increases TGs showed in the trend plots. Correlation Network Analysis To investigate the metabolic independence and latent associations among lipid species, correlation networks were separately constructed based on GPs and SPs by the criteria of |Pearson correlation coefficient| larger than 0.9 which threshold ensured the stability and robustness of the network construction. Consequently, a total of 1213 edges and 150 nodes including 38 PCs, 41 PEs, 11 PGs, 17 PIs, 12 PSs, 8 PAs, 14 LyPCs and 9 LyPEs were enclosed in the GP network diagram (Figure 5A). The positive and negative correlations between lipids were shown by the edges in red and green, respectively. The network, on the whole, could be generally identified into two connected subnets. One was formed by majorities of PAs, PSs, PIs and partial PCs, PEs (subnet-a) while the other subnet comprised of most of PCs, PEs, LyPEs, LyPCs, PGs and minorities of PSs, PIs and PAs (subnet-b). Nearly all sorts of lipids, except PCs, PEs and LyPCs, were generally gathered by their own categories, demonstrating their regular metabolic independences during passaging. In contrast, PCs and PEs were intertwined each other through the whole network bridging the two subnets, which suggested their central roles in extensive conversions with other lipid categories during BMSC aging.


23


Page 24 of 40

Figure 5. Correlation networks constructed with GPs based on |Pearson CC| larger than 0.9. (A) Network visualized by GP categories. (B) Network visualized by fold-change scale


24

Page 25 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


From another point of view (Figure 5B), the lipids in subnet-a and subnet-b showed general decreased and increased alterations in aged BMSCs, respectively. Majorities of PAs, PSs and PIs which proved to present significant decreases in p13 cells were found gathering in the outermost of subnet-a. The positive correlations between PSs and PIs were clearly observed. In addition to these lipids observed in subnet-a, there remained three PAs, six PIs and four PSs staying in subnet-b by their own clusters. As shown in Figure S8, the PAs, PSs as well as PIs in subnet-a and subnet-b were mainly differed by the altering tendencies along with passaging. The PAs in subnet-a were structurally with more carbons in their hydrophobic chains (CN=33-37) and smaller FCs (FC=0.064-0.156, n=5) than those in subnet-b (CN=30-32, FC=0.323-0.397, n=3). Especially, PAs in subnet-a revealed more rapid decreases from p7- to p9-cells than those in subnet-b, reflecting more significant affections they received from early passaging. PSs in subnet-a and subnet-b showed decreases and increases in p13-cells, respectively. Besides, the PSs in subnet-b revealed greater DBE values than the corresponding PSs with the same CNs in subnet-a. PIs in subnet-a and subnet-b revealed three diverse tendencies. Specifically, two PIs in subnet-b with CN of 36 presented increases in p13 whereas other PIs in both subnets generally decreased. Similar to PAs, the PIs in subnet-a showed more rapid declines from p7- to p9-cells than the decreased PIs in subnet-b. Also, the formers were structures with smaller DBEs than corresponding latter ones with the same CNs (0-6 versus 5-9). Notably, nearly all PIs in both subnets displayed obvious call-back manners along with passaging that were agreed with the layout observed in PCA score plot. Unlike GPs, nearly all sorts of SPs were gathering well in the network as shown in Figure S9. A total of 11 Cers, 30 SMs, 8 Gangliosides, 4 Cerebrosides, 10 Globosides, 4 sphingolid bases (SPBs) and 4 phosphosphingolipids (MIPCs) were enrolled in the network construction. As the


25


Page 26 of 40

biochemical precursors for SP biosynthesis, Cers mainly located in the centre of the clusters. The glycosphingolipids such as Cerebrosides, Gangliosides, Globosides were closely gathered, suggesting the potential biochemical relationship between glycosphingolipid categories. The SPBs such as sphinganine and sphingosine were the starting points for Cer biosynthesis. As observed in the network diagram, sphinganine (SPB_2) revealed negative correlations with Cers whereas sphingosine (SPB_1) and its forebody sphingosine-1-phosphate (S-1-P, SPB_3) presented reversed correlations that implied the fluctuated metabolic routes for Cer biosynthesis during BMSC aging. Lipidomics and Transcriptomics Based Pathway Analysis To anchor the metabolic fluctuations relevant to BMSC aging, pathway analysis was conducted via the integration of lipidomics and transcriptomics dataset (Dataset S1). The anchored enzymes were validated by using quantitative real time polymerase chain reaction (qRT-PCR) as shown in Figure S10. As shown in Figure 6A, a total of 43 enzymes (rectangles) were matched in glycerophospholipid metabolism including 5 up-regulated and 10 down-regulated in late passage (p13). For the compound mapping, a total of 13 positions (dots) involving 10 increased and 3 decreased were mapped in the pathway diagram. One of the most eye-catching features in this diagram was the activated metabolic flow of PA to CDP-DG to PI (PA/CDP-DG/PI) conversion. As shown in Table S4, the enzyme of 1-acyl-sn-glycerol-3-phosphate acyltransferase 1 (AGPAT1) for PA production showed significant increased expression in p13-cells. Meanwhile, phospholipid phosphatase 4 (PLPP4) mediating PA/DG conversion was significantly inhibited. Therein the former conversion (PA/CDP-DG) was the key segment for PI de novo biosynthesis while the latter path (PA/DG) was a sideway for PC and PE supplement. Meanwhile, the


26

Page 27 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


reductions in PIs were also observed accompanied with the up-regulation of CDP-diacylglycerolinositol 3-phosphatidyltransferase (CDIPT) which is responsible for CDP-DG/PI conversion. Comprehensively, both decreases in PAs and PIs as well as up-regulations of the enzymes mediating the transformation between them suggested the enhanced de novo biosynthesis of PIs. In contrast, the most abundant membrane GPs, i.e. PCs and PEs were all accumulated in p13BMSCs, which might be due partially to the increased secretions of DGs, the main substrate for the de novo biosynthesis of PCs and PEs. Meanwhile, lysophosphatidylcholine acyltransferases (LPCAT1, LPCAT3 and LPCAT4) and cytosolic phospholipase A2 beta (PLA2G4B) for transformations between PCs/PEs and their monoacyl-counterparts (LyPCs and LyPEs) were also inhibited that might also result in the reductions of PSs. The conversions between diacylGPs and monoacyl-GPs were generally associated with the activity of fatty acid metabolism which plays key roles in energy supplement for cellular events. Therefore, the increases of both diacyl- and monoacyl-glycerophospholipids and the decreases of PLA2G4B and LPCATs suggested the low efficiencies of material conversion and energy metabolism between GPs and LyGPs.


27


Page 28 of 40


28

Page 29 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 6. Transcriptomics and lipidomics integrated pathway diagrams presenting alterations of lipid categories and relevant enzymes in p13-BMSCs relative to p7-BMSCs. (A) Pathway diagram of glycerophospholipid metabolism. (B) Pathway diagram of sphingolipid metabolism. For sphingolipid metabolism, a total of 26 enzymes including 6 up-regulated and 8 downregulated were matched in the pathway diagram (Figure 6B). Correspondingly, 10 metabolitepositions including 7 lipid species and 3 SPBs (S-1-P, sphingosine and sphinganine) were mapped, in which 6 SPs presented increases while 4 positions involving 3 SPBs decreased in BMSCs of late passage (Table S4). Notably, as a starting point for biosynthesis of sphingomyelins and glycosphingolipids, Cers showed clear increases at the expanses of upstream substrates and downstream products during cell passaging. Ceramide synthetase 1 (CERS1) mediating sphingosine/Cer conversion was significantly up-regulated, which might further result in the reductions in S-1-P and sphingosine. For the downstream flow, CerPs also revealed significant reductions with ceramide kinase (CERK) decreased, which mediate Cer/CerP conversion. As the primary destination of Cers, SMs showed significant increases in p13BMSCs, which were just consistent with the enhanced Cer/SM conversion by the up-regulated sphingomyelin synthase 2 (SGMS2). On the other side, no significant alterations could be observed at the enzymes mediating the linkages between Cers and the globally increased glycosphingolipids. Therefore, it could be seen that the increases in Cers were realized at the expanses of S-1-P, sphingosine and CerPs, probably for the excessive generation of SMs during BMSC aging. DISCUSSION As the forebody of osteoblast, BMSC directly affects the quantity and quality of osteoblasts making it play a key role in the formation of bone tissue. Especially, one of the main causes of


29


Page 30 of 40

senile osteoporosis is the adipopexis in bone marrow resulted from the inhibited osteogenic and the enhanced adipogenic differentiation of aged BMSCs.16 These processes have been known to be strongly associated with the aberrant regulations of cell signalling pathways such as BMP/Smad, Wnt, Notch and so on.17,18 Post-transcriptional regulations mediated by microRNA as well as transcription factors such as PPARγ and Osterix were also demonstrated to be with great influences on the adipogenic differentiation capable of leading to further adipose deposition in bone marrow.19,20 Although the downstream metabolites especially various lipid species have also been reported to remarkably influence the adipose metabolism during diverse of biological events, very few studies have been carried out to investigate variations of lipids as BMSCs aging. Therefore, this present study characterized the global disturbance of lipid profiles and relevant metabolic flows during BMSC aging to address the lipidomics responses that may affect differentiation potentials of BMSCs. In a study focusing on the metabolic switch during adipogenesis, majority of GPs, mainly diacyl-PCs ， PC-Os as well as LyPCs were found to be up-regulated during the adipogenic differentiation of murine preadipocytes.21 Similarly, the increased levels of GPs such as LyPC(18:1(9Z)) were also observed in the rosiglitazone treated OP9-DL1 cells during adipogenic differentiation.22 More directly, Barretto et al. demonstrated that the adipose-derived stem cells suffering adipogenic differentiation revealed significant increases in PCs, PEs, PSs, PGs as well as DGs compared to those cultured in osteogenic differentiation and control media.23 All of these findings illustrated the high correlations between the up-regulations of GPs, especially PCs, and the adipogenic differentiation in various cells. As the most abundant GPs in eukaryotic cellular membranes, PCs have been recognized as the membrane lipids indispensable in the assembly and secretion of lipoproteins. The impaired PC biosynthesis proved to trigger the


30

Page 31 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


reduction of very low density lipoproteins (VLDLs) and high density lipoproteins (HDLs), of which increased formation also potentially requests the increased production of PCs.24 PE is another membrane amino-GP with the second highest abundance in eukaryotic cells. It has been reported as a key participant into several biological processes such as regulation of cell cycle, autophagy, modulation of membrane structures and properties.25,26 Both PCs and PEs are mainly de novo synthesized via Kennedy pathway mediated by CEPT1 from CDP-choline and CDPethanolamine, respectively.27 Also, PE is meanwhile an important methyl donor available for PC biosynthesis, which is frequently rebooted while PC being in urgent need.28 PS is well known due in part to its activities in apoptosis and immune response.29,30 Besides, it has also been proven to participate in the construction of the lipid-calcium-phosphate complexes which is capable of promoting bone formation via enhancing mineral deposition.31 In the investigation upon the effect of PSs to MSC differentiation, the addition of PSs into the culture medium presented great capacity in enhancing osteogenic differentiation of hMSCs via ERK signal pathways.32 Therefore, the reduction in PSs as BMSC aging potentially down-regulates the activity of osteogenic differentiation that might to some extent echo the increased potential of adipogenic differentiation based on their theoretical inverse relationship.33 In contrast to the GPs widely distributed in the body, SPs were mainly enriched in neutral cells and tissues presenting great activities in multiple cellular events including cell proliferation, differentiation, autophagy, apoptosis, senescence and so on.34 Particularly, as the conversion hub in sphingolipid metabolism, Cers in WI-38 human diploid fibroblasts (HDF) were found to be elevated dramatically and specifically at senescent phase and that also exhibited inhibition towards the growth of young HDF.35 These alterations were in line with the results in our study. Differently, the increases in Cers reported were mainly induced by the blocked Cer/SM


31


Page 32 of 40

conversion and the increased SM/Cer conversion whereas both reactions were opposite with corresponding results observed in the current study (Cer/SM increased and SM/Cer decreased). Instead, the accumulations of Cers in the present study proved to be associated with the activated conversions of S-1-P/sphingosine/Cer and CerP/Cer. As the main substrate for Cer biosynthesis, sphingosine was also found being with lower level in adipogenically differentiated 3T3-L1 cells than those in preadipocytes.36 Also, S-1-P has also been demonstrated to inhibit the adipogenic differentiation and improve the osteogenic differentiation in C3H10T1/2 multipotent stem cells probably via activation of PI3K/Akt signalling pathway.37 The intriguing CN- and DBE-associated alterations of lipid species during aging were of great interests. The bioconversions among homologous lipids with different chain lengths and DOUs play key roles in the regulation of a wide range of biological events and are capable of reflecting the status of cellular functionalities.38 In this study, with the employment of DVCs, it was observed that the PCs/PEs with high DBE (from 4 to 6) and CN values (from 36 to 40) presented much higher FC values than those with low DBEs (from 0 to 3) and CNs (from 28 to 35). This layout suggested that the PCs/PEs with C16:1, C18:1 and C18:2 fatty acyls (top-right area in DVC) showed relatively higher FC values than those with saturated acyl chains (C16:0 and C18:0) (bottom-left area in DVC). Similarly, the FC values of diacyl-PEs with CN of 20 were highly correlated with their DBE values. These manners were both consistent with the corresponding results revealed in the reported study, in which FA(16:1), FA(18:1), FA(18:2) showed decreases in the BMSCs of late passage while their saturated counterparts i.e. FA(16:0), FA(18:0) increased. Also, the FAs with CN of 20 revealed an increasing trend as their DBE increased.9 The long-chain saturated and mono-unsaturated fatty acids with 16 to 18 carbons are the most abundant fatty acids in animals and plants. As the starting point of the n-3 and n-6


32

Page 33 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


family fatty acids, they also play key roles in the biosynthesis of other fatty acids with shorter chain length and more double bonds. Therefore, the more significantly increased lipids with mono- and diunsaturated fatty acyls were probably synthesized at the expanse of lipids in these particular CN and DBE regions that was to some extent in accordance with the excessive consume of C16 and C18 based mono- and diunsaturated FAs in the culture medium, and the increased formation of n-6 and n-3 containing GPs in the former study.9 CONCLUSION In this study an integrated lipidomics profiling analysis was carried out upon BMSCs of sequential passaging to investigate the metabolic alterations of various lipid species during the senescence process of BMSCs. Multivariate analysis showed the BMSCs presented both a successive changing pattern and a call-back manner along with cell passaging. Majority of the lipids including GPs, SPs and GLs were up-regulated whereas the minority of GPs (PAs, PIs and PSs) and SPBs decreased in the late passage BMSCs. The metabolic alterations of some lipid categories were specific to carbon numbers and DBE values. The correlation network and pathway analysis presented the disturbed metabolic flows as well as the latent routines of metabolic conversions among lipid categories. Comprehensively, various lipid species revealed systematic variations during BMSC aging that were interpreted as potential inducing factors for aberrant differentiations of BMSCs.

Figure Legend Figure 1. Validation of the differentiation potentials of BMSCs at early passage and the morphological alterations of BMSCs between early and old passages. (A) Alizarin red staining of


33


Page 34 of 40

p7-BMSCs. (B) Oil red staining of p7-BMSCs. (C) Morphology of BMSCs at p7 (bar=25 m). (D) Morphology of BMSCs at p13 (bar=25 m). (E) Morphology of BMSCs at p7 with regional magnification (10×). (F) Morphology of BMSCs at p13 with regional magnification (10×). (G) β-Galactosidase staining of p7-BMSCs (bar=12.5 μm). (H) β-Galactosidase staining of p13BMSCs (bar=12.5 μm). Figure 2. The PCA loading-bi plots constructed with various lipid species in BMSCs of sequential passages. (A) Scattering pattern of BMSCs of D1 and D2 at p7, p9, p11 and p13. (B) Scatting patterns of GPs and GLs in the loading-bi plot. (C) Scattering pattern of SPs in the loading-bi plot Figure 3. Univariate characterizations of alterations of lipid species between BMSCs at p7 and p13. (A) –log(p) values of Student’s t-test of lipid species between BMSCs at p7 and p13. (B) Fold-change values of lipid species in p13-BMSCs relative to p7-BMSCs. (C) Pearson correlation coefficients of lipid species with passaging Figure 4. DVC plots showing fold-change alterations of categories of lipid species in p13BMSCs relative to p7-BMSCs Figure 5. Correlation networks constructed with GPs based on |Pearson CC| larger than 0.9. (A) Network visualized by GP categories. (B) Network visualized by fold-change scale Figure 6. Transcriptomics and lipidomics integrated pathway diagrams presenting alterations of lipid categories and relevant enzymes in p13-BMSCs relative to p7-BMSCs. (A) Pathway diagram of glycerophospholipid metabolism. (B) Pathway diagram of sphingolipid metabolism Supporting Information


34

Page 35 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The supporting information is available free of charge. Figure S1. The expressions of surface antigens in p7-BMSCs (PDF). Figure S2. The typical MS/MS spectra of a portion of lipids identified (PDF). Figure S3. Metabolic trend plots of various categories of lipids in BMSCs of D1 and D2 during BMSC aging (PDF). Figure S4. HCA-heatmap diagram showing the metabolic similarities of lipidomics patterns between BMSCs from D1 and D2 during passaging (PDF). Figure S5. Univariate characterizations of alterations of lipid species between BMSCs of p7/p13 and p9/p11 (PDF). Figure S6. DVC plots of LyPCs, LyPEs, diacyl-PCs, diacyl-PEs, PC-Os and PE-Os showing fold-change alterations in BMSCs of p13 relative to those of p7 (PDF). Figure S7. Fold-change values (p13 to p7) of TGs with CN of 54 and DBE of 2 compared with those at other coordinates (PDF). Figure S8. Changing trends of PAs, PSs and PIs in subnet-a and subnet-b at sequential passages (PDF). Figure S9. Correlation networks constructed with SPs of different categories based on |Pearson CC| larger than 0.9 (PDF). Figure S10. Relative mRNA expressions of enzymes anchored in pathway analysis in p7- and p13-BMSCs (PDF). Table S1. Detailed methods and parameters applied for ion feature selection based on R-platform (PDF). Table S2. The primer sequences of genes anchored in pathway analysis (PDF). Table S3. The characterization information of lipid species identified in BMSC lipidomics profiles (PDF). Table S4. The enzymes in glycerophospholipid and sphingolipid metabolism pathways with mRNA significantly altered between p7- and p13-BMSCs (PDF). Dataset S1. Transcriptomics dataset of BMSCs of different passages (txt)

Author Contributions


35


Page 36 of 40

Xin Lu took charge of the experimental design, sample preparation, data analysis and paper writing; Yue Chen took charge of the experimental design, cell culture, image capturing and article modification; Huiyu Wang took charge of sample analysis; Yunfan Bai participated into data analysis; Jianxiang Zhao and Xiaohan Zhang took parts in cell culture and RT-PCR analysis; Li Liang, Yang Chen and Chenfei Ye participated into the construction of figures; Yi Zhang took charge of the quality control for culturing BMSCs; Ting Ma and Yu Li provided overall guidance towards experimental sections. This study was financially supported by Ting Ma. ‡These authors contributed equally. Notes The authors declare no competing financial interest. ACKNOWLEDGMENT We would like to thank Biomedical Engineering Center, Harbin Institute of Technology for the help on cell culturing and sample preparation. This study was supported by the National Key Research and Development Program of China (2018YFC1312000), The Basic Research Foundation

Key

Project

Track

of

Shenzhen

Science

and

Technology

Program

(JCYJ20160509162237418, JCYJ20170413110656460).

ABBREVIATIONS BMSC, bone mesenchymal stem cell; OP, osteoporosis; GP, glycerophospholipid; FA, fatty acid; UPLC-MS, ultra-performance liquid chromatography coupled to mass spectrometry; PCA, principal component analysis; PLS-DA, partial least square discriminant analysis; DOU, degree of unsaturation; GL, glycerolipid; SP, sphingolipid; PC1, principal component 1; PC2, principal


36

Page 37 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


component

2;

PC,

phosphatidylcholine;

PE,

phosphatidylethanolamine;

PG,

phosphatidylglycerol; DG, diacylglycerol; TG, triacylglycerol; PA, phosphatidic acid; PS, phosphatidylserine; PI, phosphatidylinositol; MG, monoacylglycerol; LyPC, monoacyl phosphatidylcholine;

LyPE,

monoacylphosphatidylethanolamine;

LyPG,

monoacylphosphatidylglycerol; Cer, ceramide; SM, sphingomyelin; CN, carbon number; DBE, double bond equivalence; DVC, double bond equivalence versus carbon number plot; FC, foldchange; SPB, sphingolid base; MIPC, phosphosphingolipid; S-1-P, sphingosine-1-phosphate; qRT-PCR, quantitative real time polymerase chain reaction; AGPAT1, 1-acyl-sn-glycerol-3phosphate acyltransferase 1; PLPP4, phospholipid phosphatase 4; CDIPT, CDP-diacylglycerolinositol 3-phosphatidyltransferase; LPCAT, lysophosphatidylcholine acyltransferase; PLA2G4B, cytosolic phospholipase A2 beta; LyGP, monoacylglycerophospholipid; CERS1, ceramide synthetase 1; CERK, ceramide kinase; SGMS2, sphingomyelin synthase 2; VLDL, very low density

lipoprotein;

HDL,

high

density

lipoprotein;

AA,

arachidonic

acid;

PI3K,

phosphoinositide 3-kinase; PPARγ, peroxisome proliferator-activated receptor γ; HDF, human diploid fibroblast; EGFR, epidermal growth factor receptor; MTBE, methyl tert-butyl ether; PBS, phosphate buffer saline; MW, molecular weight; KGML, KEGG XML; ssDNA, single strand DNA; DNB, DNA nanoball; RCR, rolling circle replication; cPAS, combinatorial probe-anchor synthesis; HCA, hierarchy clustering analysis. REFERENCES (1) Ankrum, J. A.; Ong, J. F.; Karp, J. M. Mesenchymal stem cells: Immune evasive, not immune privileged. Nat. Biotechnol. 2014, 32 (3), 252-60. (2) Gregory, C. A.; Prockop, D. J.; Spees, J. L. Non-hematopoietic bone marrow stem cells: Molecular control of expansion and differentiation. Exp. Cell Res. 2005, 306 (2), 330-5. (3) Liu, H.; Xia, X.; Li, B. Mesenchymal stem cell aging: Mechanisms and influences on skeletal and non-skeletal tissues. Exp. Biol. Med. 2015, 240 (8), 1099-106. (4) Wenk, M. R. Lipidomics: New tools and applications. Cell 2010, 143 (6), 888-95. (5) Perreault, L.; Starling, A. P.; Glueck, D. H.; Brozinick, J. T.; Sanders, P.; Siddall, P.; Kuo, M. S.; Dabelea, D.; Bergman, B. C. Biomarkers of ectopic fat deposition: The next frontier in serum lipidomics. J. Clin. Endocrinol. Metab. 2016, 101 (1), 176-82.


37


Page 38 of 40

(6) Scherer, M.; Montoliu, I.; Qanadli, S. D.; Collino, S.; Rezzi, S.; Kussmann, M.; Giusti, V.; Martin, F. Blood plasma lipidomic signature of epicardial fat in healthy obese women. Obesity (Silver Spring) 2015, 23 (1), 130-37. (7) You, Y.; Lin, C.; Liang, H.; Lee, S.; Tsai, K.; Chiou, J.; Chen, Y.; Tsao, C.; Chen, J. Association between the metabolome and low bone mineral density in Taiwanese women determined by 1 H NMR spectroscopy. J. Bone Miner. Res. 2014, 29 (1), 212-22. (8) Lv, H.; Jiang, F.; Guan, D.; Lu, C.; Guo, B.; Chan, C.; Peng, S.; Liu, B.; Guo, W.; Zhu, H. Metabolomics and Its application in the development of discovering biomarkers for osteoporosis research. Int. J. Mol. Sci. 2016, 17 (12), 2018. (9) Kilpinen, L.; Tigistu-Sahle, F.; Oja, S.; Greco, D.; Parmar, A.; Saavalainen, P.; Nikkila, J.; Korhonen, M.; Lehenkari, P.; Kakela, R.; Laitinen, S. Aging bone marrow mesenchymal stromal cells have altered membrane glycerophospholipid composition and functionality. J. Lipid Res. 2013, 54 (3), 622-35. (10) Ejsing, C. S.; Sampaio, J. L.; Surendranath, V.; Duchoslav, E.; Ekroos, K.; Klemm, R. W.; Simons, K.; Shevchenko, A. Global analysis of the yeast lipidome by quantitative shotgun mass spectrometry. Proc Natl Acad Sci U S A 2009, 106 (7), 2136-41. (11) Astarita, G.; Kendall, A. C.; Dennis, E. A.; Nicolaou, A. Targeted lipidomic strategies for oxygenated metabolites of polyunsaturated fatty acids. Biochim Biophys Acta 2015, 1851 (4), 456-68. (12) Han, X.; Gross, R. W. Shotgun lipidomics: Electrospray ionization mass spectrometric analysis and quantitation of cellular lipidomes directly from crude extracts of biological samples. Mass Spectrom. Rev. 2005, 24 (3), 367-412. (13) Fehlmann, T.; Reinheimer, S.; Geng, C.; Su, X.; Drmanac, S.; Alexeev, A.; Zhang, C.; Backes, C.; Ludwig, N.; Hart, M.; An, D.; Zhu, Z.; Xu, C.; Chen, A.; Ni, M.; Liu, J.; Li, Y.; Poulter, M.; Li, Y.; Stahler, C.; Drmanac, R.; Xu, X.; Meese, E.; Keller, A. cPAS-based sequencing on the BGISEQ-500 to explore small non-coding RNAs. Clin. Epigenetics 2016, 8, 123. (14) Jang, M.; Whang, J.; Lewis, C. S.; Park, H. S. An interpretation of biological metabolites and their reactions based on relation degree of compound pairs in KEGG XML files. J. Softw. 2010, 5 (2), 187-94. (15) Black, M. B.; Parks, B. B.; Pluta, L.; Chu, T. M.; Allen, B. C.; Wolfinger, R. D.; Thomas, R. S. Comparison of microarrays and RNA-seq for gene expression analyses of dose-response experiments. Toxicol. Sci. 2014, 137 (2), 385-403. (16) You, L.; Pan, L.; Chen, L.; Gu, W.; Chen, J. MiR-27a is essential for the shift from osteogenic differentiation to adipogenic differentiation of mesenchymal stem cells in postmenopausal osteoporosis. Cell Physiol. Biochem. 2016, 39 (1), 25365. (17) Jacobsen, C. M.; Schwartz, M. A.; Roberts, H. J.; Lim, K. E.; Spevak, L.; Boskey, A. L.; Zurakowski, D.; Robling, A. G.; Warman, M. L. Enhanced Wnt signaling improves bone mass and strength, but not brittleness, in the Col1a1(+/mov13) mouse model of type I Osteogenesis Imperfecta. Bone 2016, 90, 127-32. (18) Song, B. Q.; Chi, Y.; Li, X.; Du, W. J.; Han, Z. B.; Tian, J. J.; Li, J. J.; Chen, F.; Wu, H. H.; Han, L. X.; Lu, S. H.; Zheng, Y. Z.; Han, Z. C. Inhibition of Notch signaling promotes the adipogenic differentiation of mesenchymal stem cells through autophagy activation and PTEN-PI3K/AKT/mTOR pathway. Cell Physiol. Biochem. 2015, 36 (5), 1991-2002. (19) Lee, J.; Ge, K. Transcriptional and epigenetic regulation of PPARγ expression during adipogenesis. Cell Biosci. 2014, 4 (1), 29-29. (20) Guan, X.; Gao, Y.; Zhou, J.; Wang, J.; Zheng, F.; Guo, F.; Chang, A.; Li, X.; Wang, B. miR‐223 regulates adipogenic and osteogenic differentiation of mesenchymal stem cells through a C/EBPs/miR‐223/FGFR2 regulatory feedback loop. Stem Cells 2015, 33 (5), 1589-600. (21) Halama, A.; Horsch, M.; Kastenmuller, G.; Moller, G.; Kumar, P.; Prehn, C.; Laumen, H.; Hauner, H.; Hrabe de Angelis, M.; Beckers, J.; Suhre, K.; Adamski, J. Metabolic switch during adipogenesis: From branched chain amino acid catabolism to lipid synthesis. Arch. Biochem. Biophys. 2016, 589, 93-107. (22) Tan, J.; Wang, Y.; Wang, S.; Zhang, N.; Wu, S.; Yuan, Z.; Zhu, X. Untargeted metabolomics analysis of adipogenic transformation in OP9-DL1 cells using liquid chromatography-mass spectrometry: Implications for thymic adipogenesis. Cell Biol. Int. 2017, 41 (4), 447-56. (23) Barretto, L. S.; Lessio, C.; Sawaki e Nakamura, A. N.; Lo Turco, E. G.; da Silva, C. G.; Zambon, J. P.; Gozzo, F. C.; Pilau, E. J.; de Almeida, F. G. Cell kinetics, DNA integrity, differentiation, and lipid fingerprinting analysis of rabbit adiposederived stem cells. In Vitro Cell Dev. Biol. Anim. 2014, 50 (9), 831-9. (24) Cole, L. K.; Vance, J. E.; Vance, D. E. Phosphatidylcholine biosynthesis and lipoprotein metabolism. Biochim. Biophys. Acta. 2012, 1821 (5), 754-61. (25) Emoto, K.; Kobayashi, T.; Yamaji, A.; Aizawa, H.; Yahara, I.; Inoue, K.; Umeda, M. Redistribution of phosphatidylethanolamine at the cleavage furrow of dividing cells during cytokinesis. Proc Natl Acad Sci U S A 1996, 93 (23), 12867-72. (26) Nebauer, R.; Rosenberger, S.; Daum, G. Phosphatidylethanolamine, a limiting factor of autophagy in yeast strains bearing a defect in the carboxypeptidase Y pathway of vacuolar targeting. J. Biol. Chem. 2007, 282 (23), 16736-43. (27) Payne, F.; Lim, K.; Girousse, A.; Brown, R. J.; Kory, N.; Robbins, A.; Xue, Y.; Sleigh, A.; Cochran, E.; Adams, C.; Dev Borman, A.; Russel-Jones, D.; Gorden, P.; Semple, R. K.; Saudek, V.; O'Rahilly, S.; Walther, T. C.; Barroso, I.; Savage, D. B. Mutations disrupting the Kennedy phosphatidylcholine pathway in humans with congenital lipodystrophy and fatty liver disease. Proc. Natl. Acad. Sci. U S A 2014, 111 (24), 8901-6. (28) Watkins, S. M.; Zhu, X.; Zeisel, S. H. Phosphatidylethanolamine-N-methyltransferase activity and dietary choline regulate liver-plasma lipid flux and essential fatty acid metabolism in mice. J. Nutr. 2003, 133 (11), 3386-91.


38

Page 39 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(29) Verhoven, B.; Schlegel, R. A.; Williamson, P. Mechanisms of phosphatidylserine exposure, a phagocyte recognition signal, on apoptotic T lymphocytes. J. Exp. Med. 1995, 182 (5), 1597-601. (30) Vance, J. E.; Tasseva, G. Formation and function of phosphatidylserine and phosphatidylethanolamine in mammalian cells. Biochim. Biophys. Acta. 2013, 1831 (3), 543-54. (31) Merolli, A.; Santin, M. Role of phosphatidyl-serine in bone repair and its technological exploitation. Molecules 2009, 14 (12), 5367-81. (32) Xu, C.; Zheng, Z.; Fang, L.; Zhao, N.; Lin, Z.; Liang, T.; Zhai, Z.; Zhou, J. Phosphatidylserine enhances osteogenic differentiation in human mesenchymal stem cells via ERK signal pathways. Mater. Sci. Eng. C Mater. Biol. Appl. 2013, 33 (3), 1783-8. (33) Yuan, Z.; Li, Q.; Luo, S.; Liu, Z.; Luo, D.; Zhang, B.; Zhang, D.; Rao, P.; Xiao, J. PPARγ and Wnt signaling in adipogenic and osteogenic differentiation of mesenchymal stem cells. Curr. Stem Cell Res. T. 2016, 11 (3), 216-25. (34) Young, M. M.; Kester, M.; Wang, H. Sphingolipids: Regulators of crosstalk between apoptosis and autophagy. J. Lipid Res. 2013, 54 (1), 5-19. (35) Venable, M. E.; Lee, J. Y.; Smyth, M. J.; Bielawska, A.; Obeid, L. M. Role of ceramide in cellular senescence. J. Biol. Chem. 1995, 270 (51), 30701-8. (36) Choi, K.; Lee, Y.; Choi, M.; Sin, D.; Lee, S.; Ji, S.; Lee, M. K.; Lee, Y.; Yun, Y.; Hong, J. T. Inverse relationship between adipocyte differentiation and ceramide level in 3T3-L1 cells. Biol. Pharm. Bull. 2011, 34 (6), 912-16. (37) Hashimoto, Y.; Matsuzaki, E.; Higashi, K.; Takahashi-Yanaga, F.; Takano, A.; Hirata, M.; Nishimura, F. Sphingosine-1phosphate inhibits differentiation of C3H10T1/2 cells into adipocyte. Mol. Cell Biochem. 2015, 401 (1-2), 39-47. (38) Li, M.; Yang, L.; Bai, Y.; Liu, H. Analytical methods in lipidomics and their applications. Anal Chem 2014, 86 (1), 161-75.


39


Page 40 of 40

For TOC only


40

Integrated lipidomics and transcriptomics characterization upon aging

Recommend Documents