Predicting Functional Responses of Progenitor Cell Exosome Potential with Computational Modeling.

Congenital heart disease can lead to severe right ventricular heart failure (RVHF). We have shown that aggregated c-kit+ progenitor cells (CPCs) can improve RVHF repair, likely due to exosome-mediated effects. Here, we demonstrate that miRNA content from monolayer (2D) and aggregated (3D) CPC exosomes can be related to in vitro angiogenesis and antifibrosis responses using partial least squares regression (PLSR). PLSR reduced the dimensionality of the data set to the top 40 miRNAs with the highest weighted coefficients for the in vitro biological responses. Target pathway analysis of these top 40 miRNAs demonstrated significant fit to cardiac angiogenesis and fibrosis pathways. Although the model was trained on in vitro data, we demonstrate that the model can predict angiogenesis and fibrosis responses to exosome treatment in vivo with a strong correlation with published in vivo responses. These studies demonstrate that PLSR modeling of exosome miRNA content has the potential to inform preclinical trials and predict new promising CPC therapies. Stem Cells Translational Medicine 2019;8:1212-1221.

With ongoing clinical trials using stem and progenitor cells in children, there is a greater need to better understand the potential mechanisms. Computational modeling has been used to investigate the signals that contribute to repair of the injured right ventricle. This study examines the contents of pediatric progenitor cell exosomes, uses existing computational models to predict functional mechanisms, and then makes a priori predictions on in vitro and in vivo data. This could have far‐reaching implications in personalized precision medicine for children with heart failure.

Introduction

Congenital heart disease (CHD) affects an estimated 40,000 infants each year in the U.S. 1. Improved surgical outcomes and medical management has led to an aging of the CHD population 2. Stem cell therapy, including c‐kit+ progenitor cells (CPCs), is currently being studied as a potential treatment for right ventricular heart failure (RVHF) occurring in complex forms of CHD, such as hypoplastic left heart syndrome (NCT03406884). But older CPCs, starting as early as 1 year old, have been shown to reduce reparative ability 3, 4. In a recent publication, we show that aggregating CPCs into spheres can improve their ability to repair RVHF, likely due to exosome‐mediated effects 5.

Exosomes are 20–150 nm vesicles that carry a variety of lipids, proteins, and noncoding RNAs. Altering the parental cells under different treatment conditions can alter exosome cargo, but the effect is not deterministic. We have previously shown that subjecting CPCs to hypoxic treatment alters exosome content, but not all miRNAs enriched in the exosomes were upregulated in the cell 6. Additionally, while thousands of different miRNAs were identified in these exosomes, miRNAs comprised only ∼1% of exosome noncoding RNA sequences. Therefore, it is important to investigate the combinatorial effect of groups of miRNAs and how they change in response to different treatments as it relates to exosome functional responses. We have previously demonstrated that a systems biology approach using statistical tools can model variations in exosome cargo and predict functional outcomes from exosome therapy 4, 6, 7.

In the present study, we use partial least squares regression (PLSR), a data‐driven computational modeling technique, to establish a mathematical relationship between miRNA levels of exosomes produced from monolayer (2D) and spherically aggregated (3D) CPCs, and biological responses of X and Y in P and Q types of cells after uptake of these exosomes in vitro. Then, using an unbiased approach, we reduce that model to the most important variables (miRNA signals) for the model's prediction of biological responses. We demonstrate the capability of this model to make a priori predictions of in vitro responses from additional biological cues based on the miRNA levels, and then demonstrate the model's ability to predict in vivo responses of angiogenesis and fibrosis in a rat myocardial infarction model, with strong correlation to experimentally observed responses. To the best of our knowledge, this is the first example of using a miRNA computational model trained on in vitro responses to predict in vivo outcomes. Our results represent a step toward bridging gaps between in vitro and animal or human studies.

Materials and Methods

Human Sample Acquisition and Isolation of Human CPCs

This study was approved by the Institutional Review Board at Children's Healthcare of Atlanta and Emory University. Human c‐kit+ CPCs used in this study were isolated from right atrial appendage tissue routinely removed during surgical repair of congenital heart defects as previously described 4. Child CPCs were isolated from three patients aged between 3 and 6 years, identified as 1048, 1063, and 1092.

CPC Culture and Sphere Formation

Cells were used between passages 5 and 9. As previously described, CPC spheres were formed by seeding 1,500 cells per microwell in an AggreWell400 microwell array (Stemcell Technologies, Vancouver, BC, Canada) and culturing overnight 5. CPC spheres have been previously shown to have high cell viability and promote endothelial‐like differentiation 5.

Exosome Generation

Exosomes were generated as previously described 7. Briefly, after 3 or 7 days in culture, 2D and c‐kit+ progenitor cells cultured as spheres (3D CPCs) were quiesced and cultured in fetal bovine serum (FBS)‐free media for 12 hours to generated conditioned media. Exosomes were generated from conditioned media by sequential ultracentrifugation (Optima XPN‐100, SW32Ti rotor, Beckman Coulter, Indianapolis, IN) at 10,000g for 35 minutes to remove cell debris and then at 120,000g for 70 minutes to concentrate exosomes. Successful exosome isolation was confirmed by transmission electron microscopy. Protein content of the exosome suspension was quantified by A280 measurements (NanoDrop 2000, Thermo Fisher Scientific Inc., Waltham, MA) and used to determine treatment loads.

Tube Formation Assay

Rat cardiac endothelial cells (CECs) were cultured in endothelial cell growth medium‐2 (EGM‐2) (Lonza, Walkersville, MD) supplemented with 2% FBS, 0.4% human fibroblast growth factor‐beta (hFGF‐B), 0.1% vascular endothelial growth factor (VEGF), 0.1% long arginine 3 insulin‐like growth factor (R3‐IGF‐1), 0.1% ascorbic acid, 0.1% human epidermal growth factor (hEGF), 0.1% Gentamicin/Amphotericin‐B (GA‐1000), 0.1% heparin, and 0.04% hydrocortisone per manufacturer's protocol. Following 12‐hour quiescence in EGM‐2 media supplemented with 0.2% FBS only, CECs were treated with 20 μg/ml exosomes for 12 hours. Cells were then lifted and then plated on 75 μl Matrigel (6 μg/ml protein, Corning, Corning, NY) thick gels in 96‐well plates with 15,000 cells per well to allow for tube formation. Cells were stained with 2 μg/ml Calcein AM and imaged 4 hours after plating by using a fluorescent microscope (Olympus IX71, Olympus Corporation, Shinjuku, Tokyo, Japan). Tube length was quantified using ImageJ (Fiji, National Institutes of Health, Bethesda, MD). To account for variability in tube formation between experiments, each experiment was normalized to an untreated negative control that received EGM‐2 media supplemented with 0.2% FBS only. Positive controls were cultured in EGM‐2 media containing all supplements listed above.

Fibroblast TGF‐β Stimulation Assay

Rat cardiac fibroblasts (RCFs) were cultured in Dulbecco's modified Eagle's medium supplemented with 10% FBS, 1% l‐glutamine, 1% penicillin, and 1% streptomycin. RCFs were quiesced for 12 hours in media with reduced FBS (1%), treated with 20 μg/ml exosomes for 12 hours, and then stimulated with transforming growth factor beta (TGF‐β, final concentration of 10 ng/ml) for another 12 hours. Cells were then collected and RNA was extracted using TRIzol reagent (Invitrogen, Carlsbad, CA) per manufacturer's protocol. Real‐time PCR (RT‐PCR) RT‐PCR was performed to evaluate transcript expression of connective tissue growth factor (CTGF), collagen type 1 pro‐α1 chain (COL1A1), collagen type 1 pro‐α2 chain (COL1A2), collagen type 3 pro‐α1 chain (COL3A1), and vimentin (VIM). Data were compared with positive control (TGF‐β stimulation without exosome treatment). Negative control represents RCFs that were not stimulated with TGF‐β and not treated with exosomes.

Real‐Time PCR

Real‐time PCR was then performed on the StepOne System (Applied Biosystems, Foster City, CA) based on SYBR Green fluorescence detection of PCR products. miRNA specific primer sequences were designed by Quantabio (Beverly, MA).

miRNA Sequencing and Analysis

miRNA was purified from exosomes using Qiagen miRNeasy kit (Redwood City, CA). Purified miRNA was analyzed (2100 Bioanalyzer, Agilent Genomics, Santa Clara, CA) for miRNA size, quality, and quantity. miRNA library was prepared using CleanTag Small RNA Library Prep Kit (TriLink Biotechnologies, San Diego, CA). Next generation sequencing was performed using a HiSeq3000 (Illumina, San Diego, CA) by the Emory Yerkes National Primate Research Center. miRNA sequence alignment and differential expression was processed using the Small RNA App in BaseSpace (Illumina, San Diego, CA).

To evaluate miRNA expression in additional exosome samples, miRNA was similarly isolated using Qiagen miRNeasy kit (Redwood City, CA), followed by reverse transcription using Quantabio miRNA cDNA synthesis kit (Beverly, MA) according to manufacturer's protocol. cDNA samples were then subjected to real‐time PCR to determine miRNA levels. Evaluation of miRNA in exosome‐treated cells was similarly performed using total RNA isolated from exosome‐treated cells with TRIzol as starting material.

Principal Component and PLSR Analysis

miRNA sequencing analysis determined mature miRNA hits. Hits values were first transformed to their natural logs. To generate the PLSR model, SIMCA‐P software (Umetrics, Sartorius Stedim Biotech, Umeå, Sweden) was used as previously described to solve the PLSR problem with the nonlinear iterative partial least squares algorithm 8, 9. miRNA expression determined in additional exosome samples by real‐time PCR was normalized to real‐time PCR miRNA expression data from 3D 7‐day exosomes using the fold change (ΔΔCt method). For each matched CPC donor population, fold change values were then used to generate theoretical hits per miRNA by comparing with miRNA sequencing data from 3D 7‐day exosomes. The normoxia/hypoxia data set, previously generated by Affymetrix microarray analysis, was transformed and incorporated into the PLSR model trained on in vitro data to form in vivo predictions of angiogenesis and fibrosis. From the raw expression data, the minimum value for each miRNA was subtracted across the exosome samples. A scaling factor was then created for each miRNA by considering the average miRNA expression from the 2D/3D versus the normoxia/hypoxia exosomes. In vivo angiogenesis and fibrosis were measured in this previously published study as capillaries per mm2 and percent left ventricle fibrotic area 6. Angiogenesis measurements were normalized to the lowest value before inserting into the PLSR model; higher values represent more angiogenesis and lower values represent less angiogenesis. Fibrosis measurements were normalized to the average of all values before inserting into the PLSR model; higher values represent more fibrosis and lower values represent less fibrosis.

miRTarBase was used to identify miRNAs with known targets (validated by at least three assays; http://mirtarbase.mbc.nctu.edu.tw). miRNA target pathways were analyzed by Ingenuity Pathway Analysis (Qiagen, Redwood City, CA). Top ranking networks were selected with scores defined as the −log(Fischer's Exact test result).

Statistical Analysis

The OriginPro 2016 software (OriginLab Corporation, Northampton, MA) or GraphPad (Graphpad Software, San Diego, CA) was used to perform statistical tests as described in the figure legends. All data sets were significantly drawn from a normally distributed population at the 0.05 level using Kolmogorov–Smirnov's test for normality.

Chemicals and Reagents

Unless otherwise stated, all reagents were purchased from Sigma–Aldrich (St. Louis, MO).

Results

Effect of Exosomes on Tube Formation and Fibrotic Gene Expression

Different child CPC donor populations (1048, 1063, and 1092) were cultured as monolayers (2D) or spherical aggregates (3D) produced using a microwell array as previously described 5. Exosomes were generated from 2D and 3D CPCs at days 3 and 7 in culture (Supporting Information Fig. S1). The angiogenic effects of exosomes on CECs were investigated using a tube formation assay. 3D 7‐day CPC exosomes significantly improved total tube length compared with 2D CPC exosomes and 3D 3‐day CPC exosomes (Fig. 1A, 1B). Average tube length was not significantly improved. The antifibrotic effects of exosomes on rat cardiac fibroblasts (RCFs) was investigated using a fibroblast TGF‐β stimulation assay. Treatment of TGF‐β stimulated fibroblasts with exosomes from 3D 3‐day CPCs significantly reduced COL3A1 and VIM transcript expression compared with TGF‐β only control (+; Fig. 1C). Exosomes from 3D 7‐day CPCs significantly reduced CTGF, COL1A1, COL1A2, and VIM transcript expression compared with TGF‐β only control (+), as previously published 5. RCFs had greater uptake of exosomes than CECs, with no significant differences in uptake between exosome treatment groups (Supporting Information Fig. S2).

Figure 1

Effect of c‐kit progenitor cell (CPC) exosomes on tube formation and fibrotic gene expression. (A): Representative images of blood vessel‐like structures formed by rat cardiac endothelial cells (CECs) stained with Calcein on Matrigel. CECs were treated with positive growth factor‐rich control media (+), negative growth factor‐devoid control media (−), or exosomes from 2D CPCs at day 3, 2D CPCs at day 7, 3D CPCs at day 3, and 3D CPCs at day 7. Scale bar = 100 μm. (B): Total tube length (top panel) and average tube length (bottom panel) formed by CECs were quantified in ImageJ and normalized to negative control. Data represented as individual samples and average ± SD. Dotted lines represent positive control value (+) and negative control value (−). *, p ≤ .05 with one‐way analysis of variance (ANOVA) with Tukey's post hoc analysis. (C): Fibrotic gene expression in TGF‐β stimulated rat cardiac fibroblasts (RCFs) after treatment with exosomes from 2D or 3D CPCs at day 3. Data were normalized to TGF‐β stimulation only (+). Negative control (−) represents no TGF‐β stimulation. *, p ≤ .05 versus TGF‐β only with one‐way ANOVA with Tukey's post hoc analysis. Fibrotic gene expression in TGF‐β stimulated RCFs after treatment with exosomes from 2D or 3D CPCs at day 7 were recently published in reference 5. 1048, 1063, and 1092 represent different child CPC patients.

PLSR Modeling of CPC Exosome miRNA Signals and In Vitro Responses to Exosome Treatment

PLSR was used to establish a relationship between the cues (2D or 3D and 3‐day or 7‐day) and the responses of in vitro angiogenesis and antifibrosis responses, using the signals (miRNA levels) to derive mathematical relationship in an unbiased approach (cue‐signal‐response paradigm). Of the 2,588 mature miRNAs sequenced, 218 miRNAs had nonzero hits and 107 miRNAs had nonzero values for at least three samples. PLSR identified the most important signals for the biological responses by calculating variable importance for projection (VIP) scores using a weighted sum of squares for a 2‐component model. The signals projecting strongly (either positively or negatively) with a response were highly ranked. The model was reduced to the top 40 highest ranking miRNAs (Supporting Information Fig. S3). The top 40 VIP miRNAs were then fit to a 3‐principal component model. In a 3‐principal component model, these 40 miRNAs captured the variance in the signals (X‐variable matrix) and responses (Y‐variable matrix) with high coefficients of determination (R 2 X = 0.856, R 2 Y = 0.514). The cumulative variance captured by each component is given in Supporting Information Figure S4A. The variance captured by each component for individual in vitro responses is given in Supporting Information Figure S4B. As expected, principal components 1 and 2 captured the majority of in vitro response variance. However, the inclusion of principal component 3 greatly improved the predictability of COL3A1, VIM, and average tube length (Supporting Information Fig. S4B). We generated separate models with varying numbers of top scoring VIP miRNA signals included and found that a model with the top 40 VIPs maximized R 2 X and R 2 Y values while minimizing data input required. For example, a three‐principal component model with more VIPs (top 60) had only marginally higher R 2 X (0.858) and a model with less VIPs (top 20) had reduced overall R 2 Y (0.412).

Scores plot of the 3‐principal component fit of the data shows separation of 2D and 3D CPC exosomes across principal component 1 (Fig. 2A). 3D 7‐day exosomes were further separated by principal component 2 from 3D 3‐day exosomes. These data suggest that component 1 likely represents dimension of culture condition and component 2 likely represents dimension of culture time. The poor separation of 2D 3‐day and 7‐day exosomes suggests no significant difference in miRNA expression between these conditions. Plotting the miRNA signals and the in vitro responses in the same component space in a loadings plot illustrates the influence that each miRNA signal has on a response in relation to all the other signals. The loadings plot shows separation of angiogenesis measures (total tube formation length and average tube formation length) and fibrosis measures (CTGF, COL1A1, COL1A2, COL3A1, VIM) across principal components 1 and 2 (Fig. 2B). Angiogenesis responses had negative loadings values in components 1 and 2, corresponding with 3D 7‐day cues in the scores plot. This affirms our observed data showing improved tube formation with 3D 7‐day exosome treatment. CTGF, COL1A1, and COL1A2 had positive loadings values in principal components 1 and 2. This suggests that these responses are negatively correlated with 3D 7‐day exosome cues, as affirmed by our previously published data showing reduced expression of CTGF, COL1A1, and COL1A2 with 3D 7‐day exosome treatment in vitro 5. Similarly, VIM and COL3A1 were negatively correlated with 3D 3‐day exosome cues. The top 40 VIP miRNAs predicted in vitro responses with high coefficients of determination, with the exceptions of COL1A1 and COL1A2 (Fig. 2C). The predictive ability of the model was calculated based on cross‐validation boot‐strapping method, performed by omitting an observation from training the model and then using that model to predict responses without the withheld observation. This procedure was repeated until every observation had been excluded exactly once. The cumulative value for the predictive ability of the model was evaluated after addition of each component and is given as Q 2(cum) in Supporting Information Figure S4A. Scores and loadings plots showing principal components 1 versus 3 and components 2 versus 3 are given in Supporting Information Figure S5.

Figure 2

Partial least squares regression modeling of c‐kit progenitor cell (CPC) exosome miRNA signals and in vitro responses to exosome treatment.A three‐component model was trained using the top 40 miRNA variables of importance for the model projection (variable importance for projections [VIPs]). (A): Scores plot of component 1 versus 2 from partial least squares regression (PLSR) analysis trained with 2D or 3D, 3‐day or 7‐day treatment with top 40 VIPs. 1048, 1063, and 1092 represent different child CPC donors. (B): Loadings plots of components 1 versus 2 show VIP miRNA signals covarying with in vitro responses of angiogenesis (total tube formation length, average tube formation length) and fibrosis (mRNA expression of connective tissue growth factor, COL1A1, COL1A2, COL3A1, and VIMENTIN). miRNAs with the highest weighted coefficients for angiogenesis (light blue) and fibrosis (light red) in component 3 are highlighted. (C): PLSR model predictions of in vitro responses correlated with observed measurements.

Validation of miRNA Expression and Covariance

To validate the covariance of miRNAs with responses, the expression of miRNAs miR‐28‐3p, miR‐320a, miR‐423‐5p, and miR‐323a‐3p was quantified in CECs and RCFs treated with 2D and 3D 7‐day exosomes (Fig. 3). These miRNAs were chosen because they had the highest weighted coefficients for angiogenesis and fibrosis responses in principal component 3. Exosome treatment increased expression of miRNAs compared with untreated cells, suggesting that exosomes successfully targeted and delivered miRNA cargo to CECs and RCFs. These data agree with previously published data demonstrating CPC exosome uptake by endothelial cells and fibroblasts 4, 6. Because angiogenesis and fibrosis responses were segregated across principal components in our model, we expect that miRNAs that covary with angiogenesis measures would lead to opposite changes in exosome‐treated fibroblasts versus exosome‐treated endothelial cells. Indeed, we found that 3D CPC exosome treatment increased the expression of miRNAs that covaried with tube formation (miR‐28‐3p, miR‐320a, and miR‐423‐5p) in CECs but not in RCFs. miR‐323a‐3p, which covaried with fibrotic gene expression responses, was increased in RCFs with 3D exosome treatment, but not in CECs. These data demonstrate that cumulative, small nonsignificant changes in groups of miRNAs can produce significant biological responses. miRNAs yield themselves to amplified responses because they can target several pathways and can be agonists or antagonists to multiple genes. Transfection of miRNAs with high‐positive weighted component 3 coefficients for tube formation (miR‐320a and miR‐423‐5p) into CECs increased tube formation compared with nontransfected CECs (Supporting Information Fig. S6).

Figure 3

miRNA expression in exosome‐treated cardiac endothelial cells and rat cardiac fibroblasts. Cardiac endothelial cells and rat cardiac fibroblasts were treated with 2D or 3D 7‐day exosomes. miR‐28‐3p, mir‐320a, and mir‐423‐5p covaried with tube formation measures (bottom left quadrant of loadings plot in Fig. 2B) and mir‐323a‐3p covaried with fibrotic gene expression (top right quadrant). The expression of these miRNAs was quantified using real‐time PCR in 2D and 3D exosome‐treated cardiac endothelial cells (A) and cardiac fibroblasts (B). Data represent cycle threshold (Ct) ± SD. Lower Ct values represent higher expression. Dotted line represents untreated control. Significance was tested using Student's paired t test. 1048, 1063, and 1092 represent different child CPC donors.

miRNA Target Analysis

The top 40 VIP miRNAs were queried for experimentally validated targets using miRTarBase. Thirty‐seven of the 40 miRNAs had confirmed targets by at least three assays. Pathway analysis of these targets showed significant fit to tissue development networks and important cardiovascular signaling pathways (Fig. 4). Angiogenesis is a complex process involving multiple signaling pathways. Several pathways that have been shown to play key roles in angiogenesis were highlighted by ingenuity pathway analysis including: PTEN signaling 10, JAK/stat signaling 11, P2Y purigenic receptor signaling 12, angiopoietin signaling, HIF1α signaling 13, 14, and nitric oxide signaling 15. Pathway analysis also identified apelin signaling, which has been shown to regulate endothelial cell proliferation and angiogenesis in the heart and prevent cardiac fibroblast activation and collagen production 16, 17. Additionally, 15 of the 218 miRNAs with nonzero hits were verified to target one or more of our fibrosis measures: CTGF, COL1A1, COL1A2, COL3A1, and VIM (Supporting Information Table S1). Of these 15, three miRNAs (miR‐378a‐3p, miR‐134‐5p, and miR‐320a) were included in our top 40 VIPs and all target VIM. Although the effect of exosomes on cardiac hypertrophy was not measured in vitro, pathway analysis highlighted the calcineurin‐NFAT signaling pathway, which has been shown to regulate the cardiac hypertrophic response 18. This suggests that 3D CPC exosomes could have played a key role in attenuating cardiac hypertrophy in our previously published study 5.

Figure 4

Target pathways of important miRNAs identified by partial least squares regression model. miRNA targets for the top 40 variables of importance for the model projection (variable importance for projection [VIPs]) were determined using miRTarBase. 37/40 VIPs had validated targets. Ingenuity pathway analysis was used to examine target network pathways and significant cardiovascular signaling pathways. Higher scores represent greater fit of the network or signaling pathway to miRNA targets.

A Priori Prediction of the Effect of Notch Inhibition on Exosome Function in 3D CPCs

Notch signaling plays an important role in regulating CPC differentiation 19, 20, 21, 22. Moreover, we have recently shown that 3D CPCs improve RV vessel density in a RVHF rat model in a Notch‐dependent manner 5. However, the effect of Notch1 inhibition on 3D CPC exosome‐mediated angiogenic function was not explored. To determine whether our model could be used to make a priori predictions of angiogenic function in response to exosomes from 3D 7‐day CPCs with Notch1 inhibition (3D + Ni), we used real‐time PCR to measure the expression of the top 40 VIP miRNAs in 3D + Ni exosomes. These newly generated data were normalized to 3D exosome miRNA expression as detailed in the methods and inserted into the PLSR model described in Figure 2. Regression algorithms calculated predicted total tube length with 3D + Ni exosome treatment as lower than predicted 3D exosome treatment (Fig. 5A). To verify model predictions, tube formation assay was then performed using 3D + Ni exosomes (Fig. 5B). Total tube formation length with CECs was significantly lower with 3D + Ni 7‐day exosome treatment compared with 3D 7‐day exosome treatment. Predicted total tube length values correlated strongly with observed responses (Fig. 5C). Correspondingly, one of the miRNAs with high positive weighted component 3 coefficients for tube formation (miR‐423‐5p) was decreased in CECs treated with 3D + Ni 7‐day exosomes compared with 3D exosomes (Supporting Information Fig. S7).

Figure 5

PLSR model predicts in vitro tube formation in response to 3D + Ni exosome treatment. Top 40 variable importance for projection (VIP) miRNAs were quantified in 3D + Ni 7‐day exosome samples using real‐time PCR and plugged into partial least squares regression model. (A): Model predictions of in vitro responses to exosomes from 3D 7‐day c‐kit+ progenitor cells (CPCs) transduced with Notch1‐shRNA (3D + Ni) compared with 3D 7‐day CPCs. (B): Representative image of blood vessel‐like structures formed by rat cardiac endothelial cells (CECs) stained with Calcein on Matrigel. CECs were treated with 3D + Ni exosomes. Scale bar = 100 μm. Total tube length formed by CECs were quantified in ImageJ and normalized to negative control. Data represented as individual samples and average ± SD. Dotted lines represent positive (+) control value and negative (−) control values. *, p ≤ .05 with Student's t test. (C): Predicted total tube length based on VIP miRNA expression in 3D + Ni samples correlated with new observations. 1048, 1063, and 1092 represent different child CPC donors.

Using In Vitro Trained Computational Model to Predict In Vivo Functional Outcomes

We then tested the capability of the PLSR model, trained with in vitro functional data, to predict in vivo angiogenesis and fibrosis results from a previously published study 5. We have previously shown that exosomes from CPCs isolated from neonate patients (0–1 month old) have greater ability to repair the heart after myocardial ischemia and reperfusion injury than exosomes from CPCs isolated from child patients (1–5 years old), and that hypoxic treatment enhances the ability of both neonate and child CPC exosomes to repair the heart 6. In these prior studies, angiogenesis was measured as the number of vessels in the peri‐infarct zone 28 days postexosome delivery. Reductions in cardiac fibrosis 28 days postexosome delivery were observed using picrosirius red staining, which stains for collagen I and collagen III fibers, and quantified as percent fibrotic tissue area. In vivo measurements were normalized as detailed in the methods. miRNA levels were previously quantified in neonate and child CPC exosomes after hypoxic and normoxic treatment using Affymetrix GeneChip analysis. These miRNA expression values were transformed to match our current model's VIP miRNA values analyzed using RNAseq, as detailed in the methods, and inserted into the model. Thirty out of 40 VIP signals were analyzed in the prior publication. Average tube formation length and COL3A1 expression responses were modeled with high percent predictability (i.e., coefficient of correlation) and were therefore chosen as the analogs for in vivo angiogenesis and fibrosis data, respectively (Fig. 2C). Without retraining the model against these new data, regression algorithms calculated predicted in vivo outcomes based on neonate and child, hypoxia, and normoxia exosome miRNA values. Predicted in vivo angiogenesis and fibrosis values correlated well with observed responses (Fig. 6). Using COL1A1, COL1A2, or VIM as analogs of in vivo fibrosis produced predicted values with reduced, but appreciable correlation with observed responses (COL1A1: y = 0.71x − 0.02, R 2 = 0.6742; COL1A2: y = 1.26x − 0.09, R 2 = 0.6594, VIM: y = 0.6129x + 0.35, R 2 = 0.838). Using total tube length as an analog of in vivo angiogenesis produced values with reduced predictability (y = 0.6459x + 0.527, R 2 = 0.2552). These results suggest that in vivo angiogenesis and fibrosis responses to exosomes from hypoxia and normoxia‐treated neonate and child CPCs could have been accurately predicted prior to performing the animal experiments previously published in reference 6.

Figure 6

In vitro model successfully predicts in vivo angiogenesis and fibrosis in response to c‐kit progenitor cell (CPC) exosomal therapy. Exosomes were isolated from CPCs donated by neonate (0–1 month old) and child (1–5 years old) patients cultured under normoxic and hypoxic conditions. Exosome miRNA content, determined by Affymetrix GeneChip analysis, and in vivo angiogenesis and fibrosis responses to CPC exosome therapy in a myocardial ischemia/reperfusion injury rat model were previously published in reference 6. In vitro average tube formation length and COL3A1 expression responses served as the analogs for in vivo angiogenesis and fibrosis data, respectively, in the PLS model described in Figure 2. Without additional model training, the in vitro PLS model predictions correlated with observed in vivo responses.

Discussion

Exosomes are generated by nearly all cell types and have emerged as an attractive therapy to improve cardiac repair. Herein, we use a systems biology approach to create a predictive model using exosome miRNA signals. Our results demonstrate that, by sorting for top VIP miRNA signals, a reduced model can be generated that strongly predicts angiogenesis and fibrosis changes in response to exosome treatment both in vitro and in vivo. Many studies have shown that stem cell exosomes mimic the reparative properties of their parental stem cells and that blockage of exosome secretion reduces the effectiveness of stem cell therapy 23, 24, 25, 26. Exosome therapies have low immunogenicity and have quickly progressed to clinical trials for intractable cutaneous ulcers (NCT02565264) and type‐I diabetes mellitus (NCT02138331) 27. Although these studies demonstrate promise for exosome therapy in patients with cardiovascular diseases, more comprehensive characterization of exosome cargo and better understanding of the mechanisms governing the loading of exosome cargo from the parent cell are required.

In recent years, data gathering methods have shifted focus to high‐throughput techniques capable of measuring thousands of variables. In our study, we used miRNA sequencing and evaluated the expression of 2,488 mature miRNAs under 12 different conditions. The large quantity of data variables generated by high‐throughput techniques, such as RNA sequencing, places emphasis on data analysis tools that can reduce the dimensionality of the data to form meaningful conclusions. PLSR is helpful in this regard by identifying linear relationships within the variables and matching them to outputs, which the user can associate with biological information. An important advantage of PLSR is that the fully trained model can be used to calculate quantitative predictions of responses of signal data sets from new cues that were not included in the training set. Because PLSR considers contributing vectors and not each individual data point for every cue, PLSR has the added benefit of accommodating unknown coefficients and incomplete data sets, which may have been key in using in vitro data collections to predict in vivo outcomes during this study.

In the present study, we trained a PLSR model using miRNA expression data from CPC exosomes from children undergoing CHD repair. Although patient‐to‐patient variability in responses were captured by exosomal miRNA content in our PLSR model, exosomes produced by 2D 3‐day CPCs from donor 1048 generated higher COL1A1 and lower COL1A2 fold change expression in treated RCFs than expected (COL1A1: predicted: 1.97, observed: 5.09; COL1A2: predicted: 0.92, observed: 0.26; Fig. 2C). Considering exosomes contain many other factors outside of mature miRNAs, these data suggest that exosomes from CPCs generated from patient 1048 may contain non‐miRNA factors that are responsible for the observed functional response. Studying other noncoding RNAs or proteomics analyses in this CPC population may reveal novel targets to modulate fibroblast activity. Our data additionally showed differences in miRNA expression between endothelial cells and fibroblasts treated with 3D 7‐day exosomes, suggesting the presence of two subpopulations of exosomes (Fig. 3). Quantifying exosome size also demonstrated two subpopulations of 3D 7‐day exosomes separated by size (Supporting Information Fig. S1). In light of studies demonstrating the cellular release of multiple subpopulations of exosomes with distinct biological properties upon activation or stimulation, 3D aggregation may promote the release of exosomes designed specifically to target and alter the function of endothelial cells or fibroblasts 28, 29. Additionally, there was clearly a time‐dependent expression of miRNAs within the exosomes as exosomes derived from earlier time points only induced changes in fibrotic markers, whereas exosomes derived from later time points induced changes in both fibrotic markers and angiogenesis. As several groups are now creating designer exosomes, these data could facilitate the rational design of exosome cargo and vesicle membrane to maximize exosome function and control exosome uptake 30.

Out of the 40 miRNAs selected for the model described in Figure 2, 30 miRNAs were matched with data from the study by Agarwal et al. 6. Inputting matching data for just 30 out of the top 40 miRNAs enabled the model to predict in vivo responses with strong correlation (Fig. 6). Comparing the top 40 VIP signals identified in this model with the top 30 VIPs published by Agarwal et al. reveals four miRNAs in common: miR‐125b‐1‐3p, miR‐148a‐3p, miR‐335‐3p, and miR‐486‐5p. The presence of these miRNAs in both models suggests that they may play well‐conserved roles in angiogenesis and antifibrosis across CPC donor ages (neonate and child) and CPC treatment conditions (aggregation and hypoxia). Using average tube length and COL3A1 as analogs generated the highest predictability of in vivo responses to child and neonate, normoxia, and hypoxia exosomes. Interestingly, the variance of average tube length values and COL3A1 expression captured within the model improved drastically with the inclusion of component 3, which captures an additional vector of covariance among the signals, distinct from principle components 1 and 2 (Supporting Information Fig. S4B). Together, these findings suggest that component 3 may represent the dimension of CPC donor patient age. Although we have previously published predictive PLSR models trained with exosome miRNA signals, this is the first example, to the best of our knowledge, of using a miRNA model trained on in vitro responses (cells) to predict in vivo outcomes (rodent models).

Although the model in our current study utilizes only a handful of cues (2D vs. 3D, age, oxygen tension, Notch activity), our results demonstrate that our data driven, statistical model may help to predict novel combinations of CPC therapies that maximize exosome therapeutic functionality and provide a computational platform to test them prior to costly experiments. Studies using CPCs and CPC‐related cardiac stem cells have used multiple different strategies to improve CPC function, including environmental and chemical cues, hydrogels, biomimetic scaffolds, electrical stimulation, gene‐editing tools, and lentiviral transduction 31, 32, 33, 34, 35. As more exosome miRNA data sets from these, and other, CPC and CPC‐related stem cell modifications become available, our PLSR model can be expanded to generate more robust predictions for a larger breadth of responses and cellular modifications. It should also be noted that our sample size was low for generating a model. This was due to patient availability and was consistent with our published studies using three patients as a representative cohort. Despite this low number, the patients were similar in age and diagnosis (septal defects) and we were still able to capture the variability. As more samples are obtained, we will be able to add to the model and refine our conclusions.

PLSR analysis allows the user to evaluate different types of signals and their contributions to given responses, not just limited to miRNA levels. For example, large amounts of clinical data are collected on patients currently enrolled in clinical trials for cardiac stem cell therapies, such as blood pressure, oxygen saturation, body weight, and so forth. Each of these parameters alone may not be predictive of patient clinical improvement. However, we have shown that small changes in groups of miRNAs can be correlated to significant biological responses in vitro (Figs. 1 and 3). Similarly, small changes in groups of clinical parameters may be correlated to significant clinical responses in patients. Incorporation of patient characteristics and clinical data into the PLSR model could create a clinically useful tool to predict which patients may respond well to cell therapy and which patients may respond poorly to cell therapy and may benefit from greater intervention. As exosome therapies transition into clinical trials, identifying patient characteristics associated with highly therapeutic “responder” exosomes can facilitate the development of more potent therapeutics.

Conclusion

The strong prediction of in vivo exosome function from in vitro analyses represents a step toward producing patient‐specific exosome therapeutics. These studies can identify optimal donors and stem cell manipulations predicted to improve exosome function before performing in vivo experiments. Combinations of high‐performing donor cell exosomes or stem cell manipulations can produce maximally functional exosome therapies. Further studies may use exosome miRNA data in combination with patient characteristics to predict clinical outcomes in patients. Although our approach does not identify causative mechanisms, the unbiased and quantitative selection of cues and signals allows for improved understanding of CPC exosome effects and has the potential to be extended to predict exosome function for similar outcomes in other diseases.

Author Contributions

D.T.: conception and design, financial support, administrative support, provision of study material or patients, collection and/or assembly of data, data analysis and interpretation, manuscript writing, final approval of manuscript; J.R.H.: conception and design, collection and/or assembly of data, data analysis and interpretation, manuscript writing, final approval of manuscript; S.B.: collection and/or assembly of data, data analysis and interpretation, final approval of manuscript; J.T.M.: conception and design, collection and/or assembly of data, data analysis and interpretation, final approval of manuscript; M.O.P.: data analysis and interpretation, manuscript writing, final approval of manuscript; M.E.D.: conception and design, financial support, administrative support, manuscript writing, final approval of manuscript.

Disclosure of Potential Conflicts of Interest

The authors indicated no potential conflicts of interest.

Supplemental Figure 1 Exosome size was measured using electron microscopy. (A) Exosomes were visualized using transmission electron microscopy. (B) Exosome vesicle diameter was quantified under each condition.

Click here for additional data file.

Supplemental Figure 2 CECs and RCFs uptake 2D, 3D, and 3D + Ni exosomes. 7‐day exosomes were stained with Calcein and dialyzed overnight in ddH2O to remove free dye. RCFs and CECs were treated with 20 μg/ml Calcein‐stained exosomes for 12 hours at 37°C. Before imaging, cells were washed with PBS and counter‐stained with Hoechst and CMTPX. Dialyzed free Calcein was used as a control. Representative 2D‐confocal images from a central focal plane with orthogonal images on the bottom and right for RCF (A) and CEC (B) with internalized Calcein‐stained exosomes. Blue = Hoechst, Green = Calcein, Red = CMTPX. Scale bars = 5 μm. Quantification of exosome uptake (C) for RCF and CEC cells was determined from individual cells (n = 12–20). Dotted line represents dialyzed, exosome‐free Calcein treatment.

Click here for additional data file.

Supplemental Figure 3 Top 40 miRNAs of importance for model projection. PLSR analysis identified the most important signals for the biological responses by calculating variable importance for projection (VIP) scores using a weighted sum of squares for a 2‐component model. The signals projecting strongly (either positively or negatively) with a response were highly ranked. The top 40 VIP miRNAs were chosen to model the relationship between the cues and responses in a 3‐component model.

Click here for additional data file.

Supplemental Figure 4 Contribution of model components to variable fit. (A) Sum of squares for the variation in the X matrix (R2X), sum of squares for the variation in the Y matrix (R2Y), and fraction of the total variation in the Y matrix that can be predicted (Q2[cum]) with the addition of each component in a 3‐component model of the cues and responses. (B) Proportion of in vitro response variance captured by each component in a 3‐component model of the cues and responses.

Click here for additional data file.

Supplemental Figure 5 Scores and loadings plots of principal components 1 vs 3 and 2 vs 3. Scores and loadings plots of principal components 1 vs 3 (A and B) and principal components 2 vs 3 (C and D).

Click here for additional data file.

Supplemental Figure 6 Transfection of top VIP miRNAs into CECs affects tube formation. Quiescent rat cardiac endothelial cells (CECs) underwent transfection, followed by a tube formation assay. Briefly, CECs were grown in complete EGM‐2 media and prior to transfection, were quiesced for 12 hours. For the transfection, 10 nM miRNA with 3% HiPerFect Transfection Reagent (Qiagen) in bare EGM‐2 media was added to CEC cultures and incubated for 8 hours. Low serum, quiescence media was then added to the cultures for 10 hours and followed by a tube formation assay as described above. Successful transfection was confirmed through RT‐qPCR. (A) Representative images of blood vessel‐like structures formed by transfected CECs stained with Calcein on Matrigel. Top 40 VIP miRNAs with high‐positive weighted component 3 coefficients for tube formation (miR‐320a, miR 423‐5p, miR‐28‐3p) were used for transfection. miRNA free, growth factor rich media with 3% HiPerFect was the positive control (+) and miRNA/HiPerFect free, growth factor‐devoid media was the negative control (−). Scale bar = 100 μm. (B) Total tube length formed by the CECs was quantified using ImageJ and normalized to the positive control. Dotted lines represent positive (+) and negative (−) control values. *p ≤ .05 with Mann–Whitney U test. n = 5‐6. miR‐320a and miR‐423‐5p had significantly higher levels of tube formation compared to the non‐transfected CECs, as suggested by the model. (C) Table shows the high‐positive miRNA weighted coefficients for total tube formation from our model.

Click here for additional data file.

Supplemental Figure 7 miRNA expression in exosome‐treated cardiac endothelial cells and rat cardiac fibroblasts. Cardiac endothelial cells (CECs) were treated with 3D or 3D + Ni 7‐day exosomes. miR‐28‐3p, mir‐320a, and mir‐423‐5p covaried with tube formation measures. The expression of these miRNAs were quantified using real‐time PCR in 3D and 3D + Ni exosome‐treated cardiac endothelial cells. Data represents cycle threshold (Ct) ± standard deviation. Lower Ct values represent higher expression. n = 3. *p ≤ .05 vs 3D using Student's t‐test.

Click here for additional data file.

Supplemental Table 1 miRNAs that target fibrosis genes. 15 of the 218 miRNAs with non‐zero hits were verified to target one or more fibrosis measures. 10 out of 15 of these miRNAs were included in the first reduction of the model with 107 miRNAs that had at least 3 measurements different from the median. 3 out of 15 of the miRNAs were included in our top 40 VIPs. miRNA targets were determined using miRTarBase.

Click here for additional data file.