Transcriptomic analysis of circulating extracellular vesicles during the perioperative period of Fontan and Glenn surgery

Introduction

Single-ventricle defects comprise a spectrum of malformations that generally result in one, instead of two functioning ventricles¹. These defects can be palliatively managed through a series of surgical procedures, in what is commonly referred to as Fontan palliation. The initial palliative intervention varies based on the defect type, but the procedures culminate into the Glenn (at ~6 months of age) and Fontan surgeries (2–4 years of age)². The goal is to use the single functioning ventricle to systemically pump blood while rerouting the systemic venous circulation directly to the pulmonary arteries without the use of a sub pulmonary pump. The resulting circulatory physiology is characterized by significantly elevated systemic pressure and low normal cardiac output and will be referred to herein as Fontan physiology³. Although lifesaving, we now appreciate that patients undergoing these procedures suffer from significant end organ dysfunction which leads to increased morbidity and mortality^2,4. While it is believed that this dysfunction results from abnormal hemodynamics, there are currently no known mediators for these adverse outcomes, and surveillance is primarily focused on identifying end-organ damage⁴. Thus, the need to identify novel biomarkers to predict complications early in the Fontan physiology is crucial to help determine the timing of intervention and to better understand the diverse complications which arise over time.

Analyzing the circulating transcriptome is a promising method to identify novel biomarkers in the Fontan physiology. Current efforts to identify non-invasive diagnostic markers use a targeted proteomic approach to correlate serum protein levels to post-Fontan clinical symptoms⁵. Although protein-focused methods have experienced some success, few studies investigate the potential of other molecular species like RNAs as possible biomarkers. These circulating RNAs have been shown to act as strong predictors for clinical outcomes and exist either as free molecules or enclosed within vesicle-like structures such as small extracellular vesicles (sEVs). sEVs are nano-sized lipid membranes released into the interstitial space that contain cargo consisting of DNA, protein, lipid, and RNA⁶. Evidence shows that sEVs can act both as highly sensitive biomarkers and also as potent mediators for adverse outcomes through direct cargo delivery to recipient cells⁷. As a result, sEVs are thought to reflect the physiological status of the body and are a suitable vehicle for identification of novel biomarkers in single-ventricle patients. Thus, elucidating the circulating sEV transcriptome in the Fontan physiology can not only identify novel, highly sensitive biomarkers, but it can also provide further insight into systemic changes that occur through staged palliation⁸.

Given the potential of sEVs in diagnosing post-Fontan complications, we hypothesized that the Fontan operative sequence triggers a change in the circulating sEV transcriptome in patients and such changes can serve as valuable predictors of post-Fontan complications. To investigate our hypothesis, we collected and analyzed RNA from serum and serum-derived sEVs from a cohort of pediatric patients during the perioperative period of the Fontan and Glenn procedure. First, we compared the preoperative and postoperative transcriptome of serum and sEVs from both Glenn and Fontan patients and found distinct transcript signals originating from each. An in vitro wound healing assay revealed that Fontan post-op sEVs have a higher angiogenic potential than pre-op sEVs. Thereafter, we utilized a partial least squares regression model to identify key RNA drivers of angiogenesis, a key molecular driver of venovenous collateral formations. Lastly, we leveraged a tissue deconvolution tool to identify changes in sEV expression in preoperative and postoperative samples in comparison to a previously published sheep Fontan model. This study highlights the potential of serum and sEV RNA-based biomarkers and supports them as mediators for the formation of venovenous collaterals in patients with single-ventricle physiology, paving the way for the discovery of long-term RNA-based biomarkers.

Methods

Ethics

Blood samples were collected under an approved IRB (00005500) prior to surgery and before discharge from patients undergoing Glenn or Fontan surgery. Informed consent was received from parents prior to study enrollment.

Statistics

P-values for differentially expressed genes were corrected for multiple testing using the Benjamini-Hochberg adjustment. Genes with a |log₂FC| greater than or equal to 1 and an adjusted p-value less than or equal to 0.05 were considered differentially expressed.

For the wound healing assay, a two-way repeat measures ANOVA with Tukey posthoc was performed to assess differences between pre-op and post-op (Glenn or Fontan) samples to the negative control. An adjusted p-value less than or equal to 0.05 was considered statistically significant. RT-qPCR mean expression differences were compared using a one-way ANOVA with Tukey posthoc. An adjusted p-value less than or equal to 0.05 was considered statistically significant. GraphPad PRISM 8 software (GraphPad, San Diego, CA) was used for statistical analysis for both wound healing and RT-qPCR experiments. For miRNA inhibition experiments, a paired Student’s t-test was used to assess differences in wound closure between groups with a p-value less than 0.05 being considered significant.

Serum sample collection

Human blood samples were collected under an approved IRB before discharge from patients undergoing Glenn or Fontan surgery. Sex and race were determined based on clinical registry. Samples were clotted using a serum blood collection tube (BD Vacutainer, Cat # 02-683-94) for 15 min at room temperature. Samples were then spun at 1500 × g for 15 min at 4 °C followed by immediate serum collection from the supernatant. Serum samples were stored at –80 °C for long-term use.

Sheep animal model

Details on the development of the sheep Fontan model are available in our group’s previous work^8,9. In brief, the sheep SVC was directly anastomosed to the main pulmonary artery and the IVC was connected to the pulmonary artery with an ePTFE graft. Venous blood samples were drawn from the internal jugular vein of the sheep preoperatively and up to one week postoperatively. Sheep subjects ranged between one to two years of age, an equivalent of early adolescence and early adulthood in humans. Sheep serum samples were stored at –80 °C for downstream applications.

Extracellular vesicle isolation and characterization

Serum samples were thawed on ice in preparation for sEV isolation. 1000 µL of serum was centrifuged at 10,000 × g for 10 min to remove large particles and debris. 500 µL of supernatant was loaded into a qEV 35 nm original column (Izon Science, Christchurch, New Zealand, Cat # ICO-35) to isolate sEVs. sEV fractions were eluted in 1 mL of PBS and stored at –80 °C for downstream applications. sEV size and concentration were measured via nanoparticle tracking analysis using a NanoSight NS300 (Salisbury, UK). In brief, sEV samples were diluted 1:50 in PBS, and 1 mL of sample was injected into the NanoSight at a flow rate of 500 µl/min. Three 60 s videos at different sample locations were gathered for size and concentration analysis. Additionally, sEV morphology was also characterized by negative stain and transmission electron microscopy imaging through the Emory Integrated Electron Microscopy Core.

RNA isolation and sequencing

In total, 300 µL of serum or 700 µL of EVs were used to isolate RNAs using a Plasma/Serum RNA Purification Mini Kit (Norgen Biotek, Canada) and an Exosomal RNA Isolation Kit (Norgen Biotek, Canada, Cat # 55000) per the manufacturer’s instructions. Quality control of the isolated RNAs was performed using an Agilent Bioanalyzer (TapeStation 4200, Agilent Technologies). miRNA and total RNA library preparation and next-generation sequencing were done by the Yerkes Nonhuman Primate Genomics Core at Emory University.

RNA alignment

Raw data from Illumina was converted to fastq format. Total RNA sequences were cleaned using the TrimGalore wrapper package. Illumina adapters from 3’ and 5’ ends were trimmed using the TrimGalore subpackage, cutadapt, and reads with less than 15 base pairs were removed. Additionally, reads with average quality score below 20 were removed. Cleaned reads were then further processed using the bbduk.sh package in the BBMap suite. Reads were analyzed for polybase sequences and ribosomal RNA contaminants and removed. The subsequent reads underwent alignment using the STAR aligner to GRCh38 for human total RNA and OviAri4 for sheep total RNA¹⁰. miRNA sequences were processed using a similar pipeline. QIAseq miRNA adapters from 3’ and 5’ ends were trimmed using the TrimGalore wrapper package, followed by quality trimming. Alignment of the miRNA was performed using miRDeep2¹¹.

Differential gene expression analysis

Aligned reads were filtered using the filterByExpr function from the edgeR package in R to remove lowly expressed genes and normalized with the trimmed mean of M-values (TMM) method¹². Principal component analysis was performed using the prcomp function from the stats package in R. Heatmaps were made using the pheatmap function with a euclidean distance calculation and a ward.D2 clustering method. Differential gene expression of serum and serum-EVs was done using the dream (differential expression for repeated measures) linear mixed modeling approach to account for paired sequencing data¹³. Weights were estimated using the voomwithDreamWeights function using default parameters. Sample grouping between pre-op vs post-op was treated as a fixed effect while sex and race were treated as random effects.

Data mined GEO data were filtered using the filterByExpr function and normalized with the TMM method from the edgeR package. Differential gene expressions were done using the limma-voom workflow with sample group (pre-op vs post-op) treated as a fixed effect. |log₂FC| greater than or equal to 1 and an adjusted p-value less than or equal to 0.05 were considered differentially expressed.

Pathway and protein-protein interaction analysis

For total RNA pathway analysis, differentially expressed genes were input on Metascape as a multiple gene list to compare pre-op to post-op¹⁴. Only GO Biological Processes were considered for downstream analysis. Target genes from VIP miRNAs were identified using miRTarBase¹⁵. Targets were defined as those having been validated by at least three different experimental methods. Protein-protein interactions were determined using the STRING (v.12.0) database¹⁶. Searches were limited to Homo Sapiens protein hits and the required score was set to 0.900 (highest confidence). Textmining was excluded as one of the active interaction sources.

sEV tissue origin

sEV tissue origin was determined by the EV-Origin algorithm published by Li et al.¹⁷. Briefly, sEV reads were TPM (transcripts per million) normalized. Differences in estimated tissue fractions were used to compare pre-op and post-op human and sheep samples. To compare our results with previously published data, we mined circulating sEV RNA data with a similar experimental design to ours from the Gene Expression Omnibus (GEO). When needed, gene names were changed to human gene names using the biomaRt package in R¹⁸.

HUVEC culture

HUVECs were cultured in endothelial growth medium (EGM-2, Lonza, Bend, OR, Cat # CC-4176) supplemented with 2% FBS, 1% penicillin-streptomycin, 0.4% hFGF-β, 0.1% vascular endothelial growth factor (VEGF), 0.1% ascorbic acid, 0.1% long arginine 3 insulin-like growth factor (R3-IGF-1), 0.1% heparin, 0.1% human epidermal growth factor (hEGF), 0.04% hydrocortisone, and 0.1% Gentamicin/Amphotericin-B (GA-1000), per the manufacturer’s instructions.

Wound healing assay

HUVECs were cultured to 80% confluency with endothelial cell growth medium and then seeded on 12-well plates with 1 × 10⁵ cells per well. Cells were incubated for 48 h until they formed a monolayer. A P1000 pipette tip was used to gently form a scratch on the monolayer followed by two PBS washes. Brightfield images were taken as an initial assessment of the wound size using a phase contrast microscope (Olympus IX71, Olympus Center Valley, PA). Cells were then treated with endothelial basal medium (EBM-2, Lonza, Bend, OR, Cat # 00190860) with 20 µg of pre-operative sEVs (n = 8 biological replicates with 3 technical replicates) or post-operative sEVs (n = 8 biological replicates with 3 technical replicates) and incubated for 16 hours. PBS was used as a negative control. After 16 hours, plates were imaged again to quantify wound closure.

miRNA inhibition experiments were performed using similar methods. Lyophilized miRCURY LNA miRNA inhibitors (Qiagen, Hilden, Germany) for hsa-miR-340-5p and hsa-miR-199b-5p were reconstituted with 1X TE Buffer following manufacturer’s instructions. HUVECs were cultured to 80% confluency and seeded on 12-well plates with 1 × 10⁵ cells per well. Cells were incubated for 48 h until they formed a monolayer. A P1000 pipette tip was used to form a scratch on the monolayer followed by two PBS washes. Brightfield images were taken as an initial assessment of wound size. For transfection, a 1:1 mixture of re-suspended inhibitor and EBM-diluted Lipofectamine 3000 (Thermo Fisher Scientific) was used. Cells were treated with endothelial basal medium containing either 20 µg of sEVs, or 20 µg of sEVs and 10 nM miRNA inhibitor. After 16 hours, plates were imaged again to quantify wound closure. Images for all wound healing experiments were analyzed with ImageJ (Fiji, National Institutes of Health, Bethesda, MD) using the wound healing size tool with manual input¹⁹.

sEV RNA uptake and RT-qPCR

HUVECs were seeded in a 6-well plate at 5 × 10⁵ cells per well and incubated overnight to 80% confluency. A P1000 pipette tip was used to form a scratch followed by two PBS washes. HUVECs were incubated with EBM and one of the following: 100 µgs of pre-op, post-op sEVs or a PBS equivalent. After eight hours, cells were washed to remove unbound vesicles and collected via trypsinization. RNA extraction was performed using a PureLink RNA Mini Kit (Invitrogen, Waltham, MA). RNA concentration and purity were assessed using a NanoDrop One (Thermo Fisher Scientific). Reverse transcription of 200 ng total RNA was performed using a High-Capacity cDNA Reverse Transcription Kit (Applied Biosystems, Waltham, MA) following the manufacturer’s protocol. Quantitative PCR was performed using PowerUp SYBR Green Master Mix (Thermo Fisher Scientific) in a StepOnePlus system (Applied Biosystems, Waltham, CA). Primers were designed using PrimerBank software and specific primer sequences can be found in Supplementary Table 7²⁰.

Partial least squares regression

Partial least squares regression (PLSR) analysis was performed using SIMCA-P software (version 16, UMetrics, Sartorius Stedim Biotech, Umeå, Sweden) as previously described to solve the PLSR problem with a nonlinear iterative partial least squares algorithm^21,22. mRNA and miRNA data were normalized for library size, log-transformed, scaled and centered before being incorporated into the PLSR model.

Results

Study design and small extracellular vesicle characterization

Blood serum was collected from patients with single-ventricle defects undergoing either Fontan or Glenn surgery during the perioperative period. A complete list of patient characteristics and clinical diagnosis can be found in Supplementary Table 1 and 2. A total of 25 patients were enrolled in the study with 13 out of the 25 patients undergoing the Fontan procedure and 12 out of the 25 patients undergoing the Glenn. A combined 40 serum samples and 35 sEV samples were used for analysis due to low sEV and serum RNA yields on four and nine samples, respectively (Fig. 1). 3 samples (Patients 24, 25, and 28) were used exclusively for sEV assays. Samples taken before surgery were labeled as preoperative (pre-op) and samples taken after surgery were designated as postoperative (post-op).While serum RNA alignment metrics were lower than anticipated, this is in-line with previous serum clinical reports (Supplementary Table 3, 4)²³.

Transcriptomic analysis of circulating extracellular vesicles during the perioperative period of Fontan and Glenn surgery — **Fig. 1: Graphical representation of the study design.**

We characterized the size and shape of isolated sEVs via negative stain transmission electron microscopy (TEM). TEM images showed a canonical sEV double-membrane and globular shape for both sample types (Fig. 2a). Additionally, we measured the sEV size distribution and concentration via nanoparticle tracking analysis. Both sample types had similar size ranges with pre-op sEVs having a mean size of 110 ± 10.7 nm and post-op sEVs having a mean size of 107 ± 9.8 nm (Fig. 2b).

**Fig. 2: sEV morphological characterization.**

Unsupervised transcriptomic analysis of serum and circulating extracellular vesicles

We first examined transcriptome variations between pre-op and post-op samples from Fontan and Glenn patients through clustering and differential gene expression analysis. Hierarchical clustering of Glenn sEV and serum samples revealed a clustering pattern separating pre-op and post-op samples irrespective of sex or race with larger variability present in sEV samples (Fig. 3a). Principal component analysis of our datasets showed similar separation patterns with the top two principal components accounting for 25.6% and 35.8% of the variance in sEV and serum datasets, respectively (Fig. 3b). sEV variance was primarily captured by the top six principal components while serum variance was primarily captured by the top three principal components. Clustering patterns for both serum and sEV samples were consistent across the top principal components when analyzed together (Supplementary Fig. 1a). Differential gene expression analysis between pre-op and post-op sEV Glenn samples identified 76 and three differentially expressed genes, respectively. Serum samples had fewer differentially expressed genes, with two in pre-op and three in post-op samples (Fig. 3c).

**Fig. 3: Unsupervised and differential gene expression analysis of Fontan and Glenn serum and sEV samples.**

We analyzed our Fontan samples using the same pipeline. An unbiased heatmap of top variable genes demonstrated distinct clustering between pre-op and post-op samples, independent of sex, in both sEVs and serum. (Fig. 3d). A principal component analysis revealed a visible separation between pre-op and post-op samples in serum but with more overlap in sEVs (Fig. 3e). For both, the first two principal components accounted for approximately 40% of the total variance in the dataset. Like Glenn samples, clustering patterns for sEVs and serum were consistent across the top principal components (Supplementary Fig. 1b). However, unlike the Glenn samples, we detected a larger number of differentially expressed genes in Fontan sEVs with 330 and 132 genes found in pre-op and post-op samples, respectively. Fontan serum samples, however, had six differentially expressed pre-op genes and no differentially expressed genes in post-op (Fig. 3f). Plotting the top differentially expressed pre-op and post-op sEV genes for both Fontan and Glenn showed significant variability in the pre-op samples (Fig. 3g, h). Less variance was observed in the top differentially expressed genes for serum samples in both Glenn and Fontan samples (Supplementary Fig. 2a, b). To account for the possibility that changes in the sEV transcriptome could result from surgery rather than the Fontan or Glenn procedures specifically, we mined four published sEV total RNA-sequencing datasets from the Gene Expression Omnibus with similar experimental designs as our study^24,25,26,27. Our comparison found that our top differentially expressed genes are uniquely up-and-down-regulated for both Fontan and Glenn samples (Supplementary Fig. 2c, d).

Post-Op sEV cargo enhances endothelial cell migration

We next sought to identify differences in the transcriptome between serum and their sEVs postoperatively. After combining both our Fontan and Glenn datasets, our analysis revealed that 8967 RNA species were uniquely found in sEVs, 54 were found in the serum, and 476 were found in both (Fig. 4a, b). A comparative pathway analysis of differentially expressed genes in serum and sEV samples showed that the serum transcriptome is enriched for pathways involved in regulating RNA-protein complex assemblies and cytoplasmic translation. In contrast, sEV pathways are primarily enriched for extracellular matrix organization and circulatory system processes (Fig. 4c).

**Fig. 4: Validation of post-op sEV pro-angiogenic cargo via wound healing assay.**

We then conducted a pathway interaction network analysis in our pre-op and post-op sEV samples. Pre-op sEVs showed a wider diversity of pathways ranging from cytoskeleton organization, ALK signaling, aerobic respiration, and apoptosis (Fig. 4d). Post-op sEVs, however, had more restricted pathways. Notably, post-op sEVs were enriched for motility-focused pathways including microtubule assembly, cytoskeleton regulation, collagen proteins, and angiogenesis (Fig. 4e). Our pathway analysis indicated an enrichment of motility and angiogenesis pathways in post-op sEVs, prompting us to validate these findings. A wound healing assay in HUVEC cells treated with pre-op and post-op sEVs revealed a significant improvement in wound healing (p < 0.01) after a 16-hour treatment of Fontan post-op sEVs compared to Fontan pre-op sEVs, with a mean wound closure of 21.1% ± 5.9% in pre-op Fontan sEVs and 35% ± 4.1% in post-op Fontan sEVs (Fig. 4f). Fontan post-op sEVs, however, narrowly did not meet the statistically significant threshold when compared to the negative control (p = 0.0555) (Fig. 4g). Interestingly, there was no significant wound healing improvement (p = 0.506) between Glenn pre-op and post-op sEVs, who had a mean wound closure of 29.4% ± 8.5% and 37.9% ± 8.6% for pre-op and post-op samples, respectively (Fig. 4g). To investigate this phenotypic difference, we performed a pathway analysis between post-op sEV and serum in Fontan and Glenn samples. Our analysis revealed that, although both sample groups had overlapping pathways, Fontan sEVs are more enriched for pathways involving angiogenesis – vasculature development and tube morphogenesis – than Glenn samples (Fig. 4h).

Identification of key post-op sEV RNAs in enhancing wound healing

After having identified the pro-angiogenic phenotype in our sEVs, we next sought to pinpoint specific RNAs that may drive the observed response. First, we trained a PLSR model with both pre-op and post-op sEV total RNA to create a predictive model for our wound healing results. Our model accurately matched our experimental results (R² = 0.965) but with a mild predictability score (Q² = 0.438) (Fig. 5a). We next sought to improve the predictability of our model by re-training it and selecting the top 100 variables of importance for projection (VIPs) (Fig. 5b, Supplementary Table 5). Restricting our model to the top 100 VIPs significantly improved its predictability score (Q² = 0.936) and better matched our experimental results (R² = 0.997), suggesting that our selected VIP genes may be sufficient to build a robust PLSR model (Supplementary Fig. 3a-b). To further validate the contribution of selected VIP genes, we performed RT-qPCR to assess whether direct RNA transfer via pre-op and post-op Fontan sEVs would enhance or decrease the presence of VIP transcripts in endothelial cells. Among the genes tested, GLE1 (p < 0.05) and MAPK6 (p < 0.01) were significantly upregulated following post-op sEV treatment when compared to both pre-op and negative control (mean difference of –0.69 ± 0.18 and –1.09 ± 0.17 for GLE1 and MAPK6, respectively) (Fig. 5e). Additionally, genes such as ENAH and NFATC2 exhibited significant differences (p < 0.05) between post-op and the negative control, though no significant changes were observed between pre-op and post-op sEV treatments (Fig. 5e). To better elucidate the functional role of our VIP genes in the cell, we performed a protein-protein interaction network analysis using the STRING database¹⁶. Using our top 100 VIP genes as input, we identified shared protein interactions among them and selected for the top number of interactors under the premise that overlapping interactions would have magnified influence on cellular processes. Interestingly, our network revealed that MAPK6 and GLE1 had one of the fewest interactors while NFKB1 and GFRA1 had the highest number of protein interactors among all VIP genes. Conversely, CDH1, CTNNB1, and SRC had the highest number of VIP gene interactions (Supplementary Fig. 3e).

**Fig. 5: Identification of key sEV RNA cargo in promoting angiogenesis.**

We next turned our attention to studying whether miRNAs contribute to the observed post-op sEV phenotypes. In parallel to our previous model, we built a PLSR model to identify key miRNA drivers in our sEVs. This initial model exhibited a low predictability score (Q² = 0.063) but accurately matched our experimental results (R² = 0.730) (Fig. 5c). Restricting the model to the top 50 miRNAs significantly improved its predictability (Q² = 0.447) (Fig. 5d, Supplementary Fig. 3c, d). GO pathway analysis of the top 50 miRNAs revealed enrichment for angiogenesis-associated processes: tube morphogenesis, positive regulation of cell migration, positive regulation of cell motility, and positive regulation of phosphorylation (Fig. 5f). A full table of VIP miRNAs can be found on Supplementary Table 6. To assess the functional significance of our top VIP miRNAs in angiogenesis, we transfected miRNA inhibitors following sEV treatment. Out of our tested miRNAs, we observed that the inhibition of hsa-miR-340-5p significantly diminished the wound healing response in post-op sEVs (mean difference: –24.7 ± 9.4, p < 0.05) while inhibition of hsa-miR-199b-5p improved the wound healing response in pre-op sEVs (mean difference: 14.2 ± 5.1, p < 0.05) (Supplementary Fig. 4a, b).

sEV tissue origin shifts following a Fontan or Glenn procedure

Finally, after having pinpointed the role of specific sEV cargo in driving functional responses, we investigated whether the circulating sEV transcriptome had the potential to act as a robust tool to monitor tissue response and adaptation post-Fontan or Glenn surgery. We used the EV-Origin deconvolution algorithm to quantify the fraction of sEVs coming from 16 different tissue types¹⁷. EV-Origin deconvolution revealed a dramatic shift in sEV tissue origin following surgery. For Glenn samples, we observed an increase in estimated tissue origin fraction from nerve (0.05), brain (0.05), and pancreatic (0.04) tissues. Interestingly, we also observed a decrease in cardiac (0.06) fraction postoperatively (Fig. 6a). Fontan samples had relatively similar fraction patterns. Notably, we observed an increase in the origin fraction from the brain (0.10), esophagus (0.09), and a slight increase in hepatic fractions (0.01). Cardiac (0.02), nerve (0.08), and adipose (0.1) fractions significantly decreased following the Fontan surgery (Fig. 6b).

**Fig. 6: EV-origin deconvolution of Fontan and Glenn pre-op and post-op sEV samples.**

To strengthen our findings, we employed two validation strategies: first, we used EV-Origin to analyze sEV transcriptomes from our group’s previously published sheep Fontan animal model. Second, we utilized the datasets previously mined to control our differential gene expression analysis to evaluate if the tissue responses we observed were a result of a general surgical response rather than Fontan or Glenn-specific outcomes. Like our human samples, we observed a significant shift in sEV origin fraction on the sheep model following surgery. The fraction of sheep sEVs originating from nervous, hepatic, lung, and dermis tissues increased by 0.14, 0.04, 0.03, and 0.21, respectively. We also detected a decrease in tissue origin fraction of 0.06 from cardiac tissue, 0.12 from the stomach, 0.11 from the pancreas, and 0.14 from the colon (Fig. 6c).

A direct comparison of Fontan, Glenn, sheep, and other publicly available circulating sEV datasets revealed that some of our EV-Origin results are unique to the Fontan and Glenn procedures (Fig. 6c). Our results showed that the Fontan and Glenn procedure significantly decreases cardiac sEV fractions and increases sEVs originating from the brain. Interestingly, we observed that all datasets except for the Glenn showed a significant increase in hepatic sEVs.

Discussion

To our knowledge, we report the first transcriptomic analysis of circulating RNA species in serum and serum-derived small extracellular vesicles (sEVs) throughout the perioperative period of the Fontan and Glenn procedures. This study aimed to identify transcriptomic changes in response to the Fontan and Glenn procedures, pinpoint viable circulating RNA-based biomarkers that could serve as predictive indicators of adverse outcomes in patients with the Fontan physiology and elucidate the roles of circulating sEVs as mediators in the development of systemic complications thereby enhancing our understanding of the underlying mechanisms and improving prognostic capabilities. Through this study, we gained valuable insight into the potential role of sEVs in promoting collateral development within the Fontan physiology, while identifying specific RNA cargo that could serve as biomarkers for this process.

We first chose to conduct an untargeted characterization of the serum and circulating sEV transcriptome. We gathered serum and circulating sEV samples from patients with a variety of single-ventricle defects before and after a Glenn or Fontan procedure. After isolation of these sEVs, we assessed their morphology using negative staining and nanoparticle tracking analysis. Neither method detected any meaningful change in sEV size or shape—an observation consistent with previous serum-derived sEV reports – with both methods reporting a globular sEV shape and a mean size range of around 100 nm^25,28. Following sEV physical characterization, we compared the pre- and postoperative serum and sEV transcriptomes in Fontan and Glenn patients. Interestingly, we identified a substantial amount of RNAs differentially expressed in Fontan sEVs while all other samples showed almost no increase in post-op RNAs when compared to their pre-op counterparts. This implies a more significant systemic response to the Fontan hemodynamics than to the Glenn, which could underline the stronger sEV phenotypic outcomes seen in our in vitro assays. This systemic response was reinforced in our hierarchical clustering of sEVs and serum where we showed that the pre- and postoperative samples cluster separately from one another. In terms of random effects, no clustering patterns were detected based on a patient’s race and sex in Glenn samples; however, we could not clearly distinguish an effect of race in our Fontan samples as many samples came from a predominantly white patient population which limited our race-based clustering. Strikingly, we observed a better clustering pattern separating pre-op and post-op serum and sEV samples in heatmaps rather than the principal component analysis, suggesting that much of the biological differences between the two states lies in the top variable genes in the datasets. The overlap in pre-op and post-op principal component analysis groups can be due to technical or patient-to-patient differences, as shown in previous clinical reports⁸. These technical differences can be seen from our low alignment rates in some samples – especially in serum samples – which may have restricted our ability to completely capture the biological variability for each group. To be mindful of this, we employed a paired-matched design and normalized the library-size to account for differences in total reads and focused our study primarily on the more robust sEV transcriptome. Nonetheless, for future studies, greater sequencing depth through a larger cohort collection and volume of sEVs or serum RNAs will be required.

We have previously reported on the transcriptome differences between serum-derived sEVs and their respective serum in a sheep Fontan model⁸. The same was seen in our cohort where serum RNAs were predominately involved in regulating translation processes while pre-op and post-op sEV RNAs from both surgical groups appeared to be involved in pathways involved in regulating angiogenesis. Circulating sEVs have been shown to modulate angiogenesis, so we reasoned that a shift in the transcriptome profile – as observed in our pre-op and post-op sEVs—would invariably result in significant alterations in the capacity of these sEVs to modulate angiogenesis²⁹. Additionally, many reports highlight the formation of venovenous collaterals following a Fontan procedure^9,30,31,32. The formation of these structures is mechanically attributed to the balancing of the abnormally high central venous pressure caused by a patient’s remodeled single-ventricle physiology³². Our transcriptomics findings led us to hypothesize that post-op sEVs may act as mediators in collateral formation. However, it is well-established that chronic hypoxia promotes vessel growth via the secretion of pro-angiogenic factors, and thus, due to the perioperative nature of our experimental design, we could have been observing remnants of a chronic hypoxic sEV response in pre-Fontan physiology⁹. A scratch wound healing assay directly comparing pre-op and post-op sEVs revealed a statistically significant difference in wound healing between Fontan pre-op and post-op sEVs but not in Glenn samples, suggesting that there is further upregulation of pro-angiogenic cargo independent of chronic hypoxia in Fontan post-op sEVs. Pathway analysis supports our findings with more pro-angiogenic cargo in post-op Fontan sEVs than pre-op Fontan sEVs. While the p-value for post-op Fontan samples (p = 0.0555) was just above the threshold for significance when compared to negative control, we believe that the magnitude of the biological response supports our hypothesis. Nevertheless, additional sEV RNA data from a distant period after Fontan operation will be needed to validate our results. Aside from hypoxia, the activation of regenerative and healing mechanisms from a large post-operative wound may also be the impetus and origin for our observed phenotypes. A clear separation between a short-term regenerative response and long-term collateral formation via sEVs should also be the focus in future longitudinal studies.

While sequencing technologies allow for thorough profiling of circulating sEVs, pinpointing specific sEV cargo that best predicts our observed response would allow for quicker and more cost-efficient methods to monitor collateral formations. Our group has previously shown that computational modeling via partial least squares regression (PLSR) can accurately pinpoint specific sEV RNA species that are vital for promoting an observed phenotype³³. Initially, we built a PLSR model using the complete total RNA and miRNA datasets for pre-op and post-op sEVs, however, while both models had high correlation scores to our experimental data, both models also presented with low predictability scores. We reasoned that such scores were a result of unnecessary noise feeding into our model, and by restricting the model to the top VIP RNAs or miRNAs, we observed much higher predictability scores. These more refined models indicate that the gene set critical for monitoring collateral formation may be significantly smaller than we previously anticipated. We tested our model by quantifying RNA transcripts following pre-op or post-op Fontan sEV treatment in endothelial cells. Out of all the tested RNAs, MAPK6 and GLE1 showed a significant upregulation following post-op sEV treatment. GLE1 is an RNA exporter with unclear roles in angiogenesis^34,35. However, MAPK6 is a critical integration point for multiple signaling pathways, and its role in promoting angiogenesis is well-documented^36,37. In addition to MAPK6 and GLE1, our inhibition experiments indicate that miR-199b-5p and miR-340-5p play opposing roles in promoting wound healing, with miR-199b-5p acting as an inhibitor and miR-340-5p as a stimulator. These findings align with existing literature and suggest that miR-199b-5p is upregulated in pre-op sEVs but decreases post-Fontan, allowing a shift toward a more pro-angiogenic profile^38,39. miR-340-5p’s pro-angiogenic role has also been documented though it is important to note that our findings are based on an in vitro HUVEC model⁴⁰. While this model provides valuable insights, it may not fully replicate the complexity of collateral formation in vivo. Further research should explore whether the up- or down-regulation of these transcripts in circulating sEVs persists in patients with developed collaterals and whether they can be used for diagnostic purposes, potentially incorporating more physiologically relevant models for validation.

Lastly, we wanted to determine whether circulating sEVs could also be used to monitor tissue responses induced by surgery and the unique physiology introduced by staged surgical palliation. A key aspect of sEVs is their unique property to reflect their parent cell via their cargo⁶. This allows sEVs to serve as non-invasive indicators of tissue activity, with their secretion levels potentially reflecting responses to the altered hemodynamics of the Fontan or Glenn procedures. We decided to use EV-Origin, a powerful deconvolution tool, to assess systemic tissue responses postoperatively¹⁷. Our EV-Origin analysis revealed that a decrease in cardiac sEVs and an increase in brain sEVs appear to be a distinct response to the Fontan or Glenn procedures. Interestingly, despite the sheep animal model being closer to the Fontan physiology, we were not able to detect a complete similarity in sEV origin between human and sheep Fontan samples, however, the decrease in cardiac sEV remains consistent for the sheep Fontan. The observed gradual decrease in cardiac sEVs is interesting and could be attributed to the more favorable conditions in the single ventricle, such as reduced volume load and improved oxygenation, that follow the Fontan procedure⁴. Similarly, the increase in brain sEVs coincides well with the known complications in the brain following the Fontan palliation⁴¹. Nevertheless, a more detailed study should be conducted by isolating and characterizing the transcriptome of these cardiac and brain sEVs. The isolation of these sEVs – via cardiac and brain specific sEV markers – would be an effective approach to further understand the systemic implications of this response. Further, more longitudinal studies will be needed in both human and animal models to better correlate the increased or decreased secretion of circulating sEVs and the physiological state of their tissue of origin.

In summary, our study revealed a clear transcriptomic difference between both Fontan and Glenn preoperative and postoperative sEV samples with milder differences found in the circulating serum transcriptome. Fontan post-op sEVs showed a stronger pro-angiogenic function when compared to sEVs originating from other groups, and we identified MAPK6, GLE1, hsa-miR-340-5p, and hsa-miR-199b-5p as key RNA drivers for the observed wound healing response. Further, sEV deconvolution revealed that cardiac and brain sEVs change in response to either Fontan or Glenn surgery in humans. The outcome of this study provides evidence for sEVs as potent biomarkers in the formation of venovenous collaterals and systemic patient monitoring. Lastly, our analysis proposes key RNAs and miRNAs that can be used as tools to be used in prognostic applications.