Gene-expression profiles of abdominal perivascular adipose tissue distinguish aortic occlusive from stenotic atherosclerotic lesions and denote different pathogenetic pathways

Perivascular adipose tissue (PVAT) helps regulate arterial homeostasis and plays a role in the pathogenesis of large vessel diseases. In this study, we investigated whether the PVAT of aortic occlusive lesions shows specific gene-expression patterns related to pathophysiology. By a genome-wide approach, we investigated the PVAT transcriptome in patients with aortoiliac occlusive disease. We compared the adipose layer surrounding the distal aorta (atherosclerotic lesion) with the proximal aorta (plaque-free segment), both within and between patients with complete aortoiliac occlusion (Oc) and low-grade aortic stenosis (St). We found that PVAT of the distal versus proximal aorta within both Oc- and St-patients lacks specific, locally restricted gene-expression patterns. Conversely, singular gene-expression profiles distinguished the PVAT between Oc- and St-patients. Functional enrichment analysis revealed that these signatures were associated with pathways related to metabolism of cholesterol, vessel tone regulation, and remodeling, including TGF-β and SMAD signaling. We finally observed that gene-expression profiles in omental-visceral or subcutaneous fat differentiated between Oc- and St-patients, suggesting that the overall adipose component associates with a different atherosclerosis burden. Our work points out the role of PVAT and, likely, other adipose tissues play in the pathophysiological mechanisms underlying atherosclerotic disease, including the abdominal aortic occlusive forms.

Nonetheless, none of these DE transcripts stood correction for multiple testing, in any comparison. Consistently, the histograms of the P-value distributions did not fit that expected for truly DE genes ( Figure IIA,B,C in the online-only Data Supplement 1). These results indicate a lack of substantial locally-restricted differences in the PVAT of the distal aorta compared to the proximal abdominal aorta, both for patients with occlusive and stenotic lesions.
Comparison of PVAT of the abdominal aorta between patients with occlusive versus stenotic abdominal aortic lesions. AT samples clustering by PCA on the global gene expression profiles ( Fig. 1) suggested that the abdominal aortic PVAT could present relevant differences between the two subgroups of PAD patients. Indeed, differential expression analysis comparing PVAT in Oc-vs. St-patients revealed 210 DE transcripts with |log 2 FC | ≥0. 38, ranging from 1.38 to -2.06, that stood correction for multiple testing at an false discovery rate (FDR)-adjusted P-value < 0.05. Among them, 26   www.nature.com/scientificreports www.nature.com/scientificreports/ both the distal (diseased) and proximal (plaque-free) segments of the abdominal aorta (Fig. 3A). Interestingly, unsupervised clustering by applying a PCA to omental-visceral and subcutaneous AT based on these 210 DE transcripts also provided good separation between samples of patients with occlusive lesions and those with stenotic lesions (Fig. 3B).

Functional inferences from genome-wide differential expression analysis of Oc-versus
St-PVAT. We tested which Gene Ontology (GO) Biological Processes (BP)/pathways were associated with the Oc-and St-PVAT phenotypes by Gene Set Enrichment Analysis (GSEA). We found that 7 and 129 GO-BP/ pathways were positively associated with the Oc-PVAT and the St-PVAT samples, respectively (FDR q value < 0.05; see online-only Data Supplement 3A). We further summarized GSEA results into an enrichment network to visualize the relationships among the most relevant GO-BP/pathway gene-sets and help data interpretation (Fig. 4). GO-BP/pathways associated with Oc-PVAT mainly related to cholesterol, sterol, and alcohol biosynthetic processes. Conversely, the most interconnected GO-BP/pathways upregulated in St-PVAT involved those related to muscle and circulatory system processes (including smooth muscle contraction and aorta development), regulation of blood circulation, regulation of cell-substrate adhesion and cell junction assembly, collagen and elastic fiber formation, cell-matrix adhesion, negative regulation of coagulation and platelet degranulation. Among signaling pathways, we observed significant associations of St-PVAT with those linked to platelet-derived growth factor (PDGF), cyclic guanosine monophosphate (cGMP), type I interferon, response to TGF (transforming growth factor) beta, and positive regulation of pathway-restricted SMAD protein phosphorylation, as well as with positive regulation of cytokine-mediated signaling pathway and regulation of calcium and calcineurin-mediated signaling. tively. Significant GO-BP/pathways characterizing Oc samples included eukaryotic translation, granulocyte migration, and response to chemokine, whereas those related to St samples encompassed fatty acid metabolism, low-density lipoprotein receptor particle metabolic process, Golgi associated vesicle biogenesis, and regulation of ERBB/EGFR (epidermal growth factor receptor) signaling pathway.

Discussion
PVAT is a well-recognized regulator of vessel homeostasis and its dysfunction may strongly influence the pathogenesis of vascular diseases 21,22 . Herein, we tested on a genome-wide scale whether distinctive gene expression patterns were associated with the PVAT surrounding occlusive and stenotic segments of the abdominal aorta in PAD patients. Our study showed that perilesional aortic PVAT in PAD patients lacks specific gene expression signatures, as we did not observe any significant DE gene compared to the PVAT of the plaque-free segments in subjects with either occlusive or stenotic disease. This indicates that in PAD the entire PVAT of the abdominal www.nature.com/scientificreports www.nature.com/scientificreports/ aorta had uniform gene expression profiles with no evident locally restricted differences in the diseased segment. On the contrary, the PVAT of the abdominal aorta displayed gene expression signatures distinguishing PAD patients with a different burden of atherosclerotic lesions, ie, patients with occlusive versus sub-occlusive lesions. Furthermore, unsupervised clustering analysis of AT samples by PCA showed a clear, not obvious separation of PAD patients with occlusive from those with stenotic lesions for all the AT depots considered. This suggests that subgroups of patients with PAD that differ in severity of abdominal aortic disease may show significant different expression profiles at the systemic level involving other AT depots.
Consistent with this hypothesis, we identified three main biological processes that characterize the AT of these two forms of PAD, with possible pathogenetic implication: (i) metabolism of lipids (ie cholesterol and fatty acids), (ii) maintenance/remodeling of the arterial vessel and (iii) involvement of immune response.
The lipid component may play a role in the pathogenesis of the two forms of PAD through different mechanisms. Accumulation and modification (aggregation, oxidation or enzymatic cleavage) of cholesterol in the sub-endothelial layer of the arteries is a well-known mechanism to produce pro-atherogenic cholesterol particles. Over-representation of cholesterol biosynthetic pathways in the PVAT of Oc-patients suggests that changes in its metabolism 23 may result in a strong pro-atherogenic effect that may produce a severe and faster occlusion of the aorta in these patients. Conversely, the association of St-patients with the "fatty acid metabolism" pathway in both the omental-visceral and subcutaneous AT and the "low-density lipoprotein particle metabolic process" in the subcutaneous AT is intriguing, since both the visceral and subcutaneous AT are energy storage tissues and lipid catabolic regulation can be triggered to reduce the effects of atherogenic stimuli, as demonstrated in in vivo animal models 24 . Moreover, PAD has been suggested to have a unique lipoprotein signature compared to coronary and cerebrovascular diseases, mainly characterized by increased low-density lipoprotein particles rather than LDL cholesterol content 25 . As AT are also known to regulate serum lipids 26 , reducing LDL particles through upregulation of the LDL particle receptor in AT may represent a feedback mechanism to lessen circulating LDL particles in St-patients, similar to what happens with the use of lipid-lowering therapies 27 . An exhaustive assessment of the lipid profile, including standard lipid concentrations, lipoproteins (eg, LDL particles), and proteins affecting lipoprotein homeostasis (eg, PCSK9), could provide a deeper insight into the role of atherogenic lipids and lipoproteins in these two types of PAD.
The second element which may account for a different pathogenesis between Oc-and St-patients concerns the maintenance and remodeling of the vessel structure and the control of the vessel tone. We may hypothesize www.nature.com/scientificreports www.nature.com/scientificreports/ that a long-lasting atherosclerotic process, which characterize our St-patients, may trigger an adaptive response that tries to maintain the balance between arterial function and structure. Indeed, St-patients over-expressed genes related to pathways of elastic fiber and collagen formation, aorta and muscle structure development and contraction, and regulation of cell-substrate adhesion that may play a relevant role to counteract a vessel injury 28 . This is consistent with the knowledge that PVAT has either pro-or anti-atherosclerotic properties 29 . In this context, "TGF (transforming growth factor)-beta response", which was interconnected with the above-mentioned pathways, is known to have athero-protective effects and tissue repair properties 30 . Moreover, TGF-beta also exert its pro-fibrotic effects through SMAD signaling to induce matrix-related genes, such as collagens, fibronectin, plasminogen activator inhibitor, and proteoglycans (see functional network; Fig. 4). Interestingly, recent findings showed that PVAT-derived mesenchymal stem cells contribute to vascular remodeling in vivo through smooth muscle cell differentiation and metabolic reprogramming promoted by TGF-signaling and specific microRNA regulation 31 . Further investigations for histological analysis of the full-thickness aorta to find marker of fibrosis and extracellular matrix remodeling may strengthen the hypothesis that TGF-beta signaling may have indeed a distinct pathogenic role in Oc-versus St-patients.
Finally, the omental-visceral AT of St-patients showed over-represented pathways related to immune functions, including both innate and adaptive response. Both the omentum and subcutaneous adipose tissue are populated with immune cells and are associated with adverse metabolic risk factors, although the omental-visceral AT is recognized to play a major role [32][33][34] . Since the activity of immune cells in AT can affect adjacent tissues and organs 35 , we may hypothesize that the omental AT of St-patients triggers an adaptive response that attempts to counteract an evolving atherosclerotic process. Testing if Oc-and St-patients have a different immune profile, eg by examining their T-and B-cell receptor repertoires in these ATs and/or systemic blood circulation, would provide new insights into the role of various components of the immune system.
To our knowledge, this is the first study that describes the transcriptome of PVAT in PAD patients with AIOD or diffuse stenosis of the aorta and CIA through a genome-wide approach and has several strengths.
First, transcriptome analysis is an effective approach to reveal even subtle changes due to environmental stimuli, genetic and epigenetic background, and a different pathophysiological status. Furthermore, it favors applying a data-driven strategy that often allows identifying unexpected findings and new interpretations 16 . Second, the use of a paired-data approach reduces the effects of heterogeneity among subjects and increases the overall sensitivity and statistical power of the analysis 18 . Third, the state-of-the-art data mining procedure applied herein limited the effect of biases due to data heterogeneity, which commonly affects genome-wide data, and led to a more accurate data interpretation 19,36 . Fourth, functional data interpretation was performed with a 'competitive' , multivariate method, which ensures to capture consistent and not spurious relationships between phenotypes and genes related to GO-BP/pathways with biological meaning, even for small sample size studies [37][38][39][40] . Consistently, we found groups of interconnected pathways that look coherently related to our PAD phenotypes.
This study has also clear limitations. Our PAD patients' cohort is relatively small, and this is mainly due to the paucity of biopsy material available in current clinical settings. Furthermore, PAD is a group of vascular diseases anatomically characterized by stenosis or occlusion of one or more arteries between the aorta and the upper or lower extremity arteries and, thus, presents a heterogeneous population of patients. For this reason, the differences we revealed in our study between Oc-versus St-patients might have been guided by a center-bias selection of the patients. To generalize the findings, our results should be tested on an independent cohort of patients. Finally, we drew our conclusions based on an inferential analysis. Although large-scale and multivariate analyses have been widely shown to be effective in identifying molecular mechanisms related to specific phenotypes, interfering with such mechanisms would likely clarify their putative role in the context of AIOD.
In conclusion, although this is a proof of concept study, our work highlights and confirms the importance of PVAT for the understanding of the underlying pathophysiological mechanisms in abdominal aortic diseases. It would also support the notion that PVAT plays a different pathogenic role in aortic distal atherosclerosis compared to AAA 17 , which are to be considered diseases with different development and evolution mechanisms. Furthermore, it suggests that altered pathways in PVAT, mainly involved in the regulation of lipid metabolism, maintenance and remodeling of the aorta, and immune functions, are functionally associated with different subgroups of PAD patients and may indicate a different natural history of the two conditions. As a clinical exploitation of these findings, peripheral biomarkers could be sought to help distinguish between different PAD patients, e.g. through non-invasive circulating blood testing and/or analysis of easily accessible subcutaneous fat deposits with minimally invasive procedures. This could help to classify (or re-classify) patients showing a different obstructive arteriopathy and characterize their clinical outcomes. Finally, the inferences on specific pathogenic mechanisms presented herein could motivate additional research with potentially clinical significance. For example, genes associated with aorta remodeling, including the TGF-β pathway, can represent possible targets for therapeutic interventions.

Methods
Anonymized data and materials have been made publicly available at the NCBI's GEO repository and can be accessed at https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE136822. Exclusion criteria included Marfan syndrome and other genetic disorders of the elastic fiber system, active or recent (5 years) cancer, recent major surgery (6 months), aneurysms or disorders of the immune system such as autoimmune diseases or vasculitis. The Ethical Committee of Centro Cardiologico Monzino approved the study and all the participants signed written informed consent.  17 ) as follows: subcutaneous abdominal fat, omental-visceral fat, PVAT surrounding the distal abdominal aorta, and PVAT surrounding the proximal (upper) abdominal aorta free of plaques or thrombus (by angiography and surgeon's visual inspection at the time of proximal graft anastomosis). Thus, we distinguish the AT surrounding abdominal aorta according to the position of the segment and on the type of lesion, i.e.: DAOc-and PxOc-PVAT, and DASt-and PxSt-PVAT. AT samples were promptly snap-frozen in liquid nitrogen and stored at −80 °C until processing.
Microarray gene expression analysis. The TRIzol Reagent (Thermo Fisher Scientific) was used to extract total RNA from 50-100 mg of frozen samples. To remove genomic contamination, RNAs were treated with TURBO DNase (Thermo Fisher Scientific), following the manufacturer's instructions. RNA yield/purity and integrity were assessed using the Infinite M200 PRO multimode microplate reader (Tecan) and the 2100 Bioanalyzer (Agilent Technologies), respectively. Of the total 44 AT samples, 1 subcutaneous AT was discarded because of its poor RNA yield and quality. Gene expression assays were performed through the HumanHT-12 v4 Expression BeadChips (Illumina). RNA isolation and microarray protocols are described in detail in Piacentini et al. 17 .
Data processing. Array data export and quality control were performed with the Genome Studio Software v2011.1 (Illumina). Raw data were imported into the R software v3.5.0 and normalized with the lumi R/ Bioconductor package 41 . Probes with a detection P value < 0.01 in at least 30% of the total samples were retained. The DaMiRseq R/Bioconductor package was used to identify unwanted sources of variation (aka latent variables) to control the systematic heterogeneity as produced by high-dimensional data 42 . The matrix of expression values was then adjusted for the presence of latent variables by the DaMiR.SVadjust function. Microarray probes were annotated through the lumiHumanIDMapping and biomaRt R/Bioconductor package to retain only those with the most up-to-date annotation 43,44 . Statistical analysis. Statistical analysis was carried out in the R environment v3.5.0. Unsupervised clustering analysis of AT samples was performed by PCA using normalized and adjusted expression data.
The limma R/Bioconductor package was used to perform differential expression analysis 45 . An additive linear model for a multi-level experiment of paired samples, including "tissue-type" (ie, DA-PVAT, Px-PVAT, omental-visceral and subcutaneous AT) and "patient sub-group" (ie, those with occlusive or stenotic lesions of abdominal aorta) as factors, was designed as suggested by Smith et al. 46 The implementation of this statistical model, adjusted for the latent variables, allowed computing paired-samples comparisons within DA-PVAT and Px-PVAT specimens and between the two patient groups.
Transcripts with an absolute log 2 fold-change | (log 2 FC)| ≥0.38 at a FDR-adjusted P-value < 0.05 were reckoned as significantly different. The robustness of the differential expression analysis results was assessed by exploring the histograms of the P-value distribution. For truly DE genes, the shape of the histogram should show a uniformly flat distribution across the unit interval (null P values) with a peak near zero (P values for alternative hypotheses) 47 .
The sizepower R/Bioconductor package 48 was used to assess the power of the differential expression analysis for the two patients' subgroups. A sample size of n = 5 paired-samples was estimated to provide a statistical power of 90% under certain conditions (see the online-only Data Supplement 1 for parameter details).
Functional inferences on genome-wide expression profiles. Biological functions associated with the differences observed by differential expression analysis were inferred by taking advantage of prior biological knowledge on genes grouped by GO-BP and by Reactome pathway database (http://www.reactome.org/) using GSEA, v3.0 38 . To reduce redundancy and visually interpret GSEA data, a network of the most significant GO-BP/Reactome pathway gene-sets was drawn through the Enrichment Map software v.3.2.0 49 , implemented as a plug-in in the Cytoscape v.3.7.0 platform 50 .
Details on data processing and analysis parameters are available in the online-only Data Supplement 1.