RNA-seq 1 and RNA-seq 2 for donor versus PAH: Quality assessment, library preparation, and sequencing

The attached word file has the requested information regarding the versions of packages used to analyze the data are provided in the method description. Please find the attached schematic of the bioinformatics pipeline used for analyzing RNA-seq 1 and RNA-seq 2 datasets. 

RNA-seq 1 and RNA-seq 2 for donor versus PAH: quality assessment, library preparation, and sequencing 

FBs isolated from donor human pulmonary arteries (n = 4) and patients with PAH (n = 4) and cultured ex vivo [RNA-seq 1: early passage (GEO accession: GSE193776) and RNA-seq 2: late passage (GEO accession: GSE144932)] were used for generating transcriptome profiles using RNA-seq and were subjected to bioinformatics analysis. RNA was isolated using the miRNeasy Mini Kit (Qiagen) combined with On-Column DNase I Digestion Set (DNase-Free DNase Set, Qiagen) to avoid contamination by genomic DNA. The RNA and library preparation integrity were verified using a BioAnalyzer 2100 (Agilent Technologies) or LabChip Gx Touch 24 (Perkin Elmer) instrument. The integrity of total RNA was evaluated, and the RNA integrity number ranged between 7.4 and 10. Total RNA (2 µg) was used as input for TruSeq Stranded mRNA Library preparation, according to the low sample protocol (Illumina). Sequencing was performed on an Illumina NextSeq500 instrument using v2 chemistry. 

RNA-seq: general bioinformatics pipeline

Raw reads were assessed for quality, adaptor content, and duplication rates using FastQC 0.10.1; trimmed by Trimmomatic 0.33; and aligned to the Ensemble human genome version hg19/GRCh37 (RNA-Seq 1) or hg38/GRCh38 (RNA-Seq 2-5) by STAR >2.4.0a. The number of reads aligning to genes was counted using the featureCounts >1.4.5-p1 tool from the Subread package. DEGs were identified using DESeq2 version >1.62. For RNA-Seq 1, only genes with a minimum FC of ±2, a maximum Benjamini-Hochberg corrected P value of 0.05, and a minimum combined mean of 5 reads were classified as significant DEGs. For RNA-Seq 2 to 5, only genes with a minimum FC of ±1.5, a maximum Benjamini-Hochberg corrected P value of 0.05, and a minimum combined mean of 5 reads were classified as significant DEGs [mean read count ≥ 5, −0.59 ≤ log2 fold change (FC) ≥ 0.59, false discovery rate (FDR) ≤ 0.05]. Correlations of raw gene counts per sample were assessed using Spearman’s rank correlation algorithm. Heatmaps of DEGs were prepared according to DESeq2 (mean expression ≥ 5, −1 ≤ log₂ FC ≥ +1 and −0.59 ≤ log₂ FC ≥ 0.59, FDR ≤ 0.05). Hierarchical clustering based on Pearson’s distance correlation was performed to achieve clustering of upregulated and downregulated genes.
 

Figure S2: Bioinformatics pipeline used for RNA-seq data analysis.
 

Pathway analysis of RNA-seq datasets

A pathway overrepresentation analysis was performed for RNA-seq DEGs against several KOBAS databases (KEGG pathway, KEGG disease, GO, and Reactome pathways), and we used uncorrected P values to shortlist significantly enriched terms are displayed as bar charts with −log₂ (P value). For network representation, DEGs were imported into the Cytoscape application using ClueGO. An enrichment analyses was performed vs. KEGG pathways at a P value < 0.05. A two-sided hypergeometric test was performed and corrected using Bonferroni step-down method. The kappa value was set at 0.35.