Comprehensive analysis of m⁶A methylomes in idiopathic pulmonary arterial hypertension

Figures & data

Figure 1. Construction a transcriptome-wide map of m⁶A mRnas identified in the peripheral blood of healthy controls and IPAH patients. (a) the m⁶A content in total RNA form the peripheral blood of CON group and IPAH group. (b) Venn diagram showing the number of methylated m⁶A sites in CON and IPAH. (c) Pie chart showing the proportion of differentially methylated m⁶A sites in the five transcriptional sequences-5’UTR, 3’UTR, stop codon, start codon, CDs. (d) Bar chart showing the distribution of differentially methylated m⁶A sites on chromosomes 1–22 and sex chromosomes. (e) Relative density of DMMSs in each chromosome normalized by length of the respective chromosome. Orange, upmethylated m⁶A sites in IPAH patients; Green, downmethylated m⁶A sites in IPAH patients.

Table 1. Total numbers of differentially methylated N6-methyladenosine sites and associated genes.

Item	upmethylated site	upmethylated gene	downmethylated site	downmethylated gene
mRNA	506	479	2170	1858

Figure 2. Function annotations of genes containing DMMSs. (a) Top 10 biological processes and molecular function of differentially downmethylated genes. (b) Top 10 biological processes and molecular function of differentially upmethylated genes. (c) Top 10 pathway analysis of differentially downmethylated genes. (d) Top 10 pathway analysis of differentially upmethylated genes.

Table 2. Biological process of downregulated m⁶A methylation.

GO.ID	Term	Count	Fold. Enrichment	P-value	FDR
0097201	negative regulation of transcription from RNA polymerase II promoter in response to stress	6	5.431692308	0.000341824	0.011978003
0045945	positive regulation of transcription by RNA polymerase III	5	4.937902098	0.001889292	0.04367133
0007183	SMAD protein complex assembly	6	4.655736264	0.00094497	0.025665989
0097284	hepatocyte apoptotic process	5	4.526410256	0.002994627	0.060818804
0043922	negative regulation by host of viral transcription	5	4.526410256	0.002994627	0.060818804
0010457	centriole-centriole cohesion	5	4.526410256	0.002994627	0.060818804
0031065	positive regulation of histone deacetylation	7	4.473158371	0.000469004	0.015333207
0061157	mRNA destabilization	9	4.444111888	7.6092E–05	0.003477097
0033148	positive regulation of intracellular estrogen receptor signaling pathway	4	4.345353846	0.009544623	0.138817811
0031848	protection from non-homologous end joining at telomere	4	4.345353846	0.009544623	0.138817811
0031053	primary miRNA processing	4	4.345353846	0.009544623	0.138817811
0007016	cytoskeletal anchoring at plasma membrane	4	4.345353846	0.009544623	0.138817811
0002664	regulation of T cell tolerance induction	4	4.345353846	0.009544623	0.138817811

‘Count’ stands for the number of DE genes associated with the listed GOID.

Table 3. Biological process of upregulated m⁶A methylation.

GO.ID	Term	Count	Fold. Enrichment	P-value	FDR
0034453	microtubule anchoring	6	10.23362319	1.94662E–05	0.003761417
0045063	T-helper 1 cell differentiation	4	10.03296391	0.000556815	0.067245369
0034134	toll-like receptor 2 signaling pathway	3	9.840022297	0.003074075	0.189130533
0045625	regulation of T-helper 1 cell differentiation	2	8.528019324	0.021798932	0.580418014
0042908	xenobiotic transport	2	8.528019324	0.021798932	0.580418014
0038003	opioid receptor signaling pathway	2	8.528019324	0.021798932	0.580418014
0034058	endosomal vesicle fusion	2	8.528019324	0.021798932	0.580418014
0032836	glomerular basement membrane development	2	8.528019324	0.021798932	0.580418014
0014894	response to denervation involved in regulation of muscle adaptation	2	8.528019324	0.021798932	0.580418014
0014877	response to muscle inactivity involved in regulation of muscle adaptation	2	8.528019324	0.021798932	0.580418014
0003356	regulation of cilium beat frequency	2	8.528019324	0.021798932	0.580418014

‘Count’ stands for the number of DE genes associated with the listed GOID.

Table 4. Molecular function of downregulated m⁶A methylation.

GO.ID	Term	Count	Fold. Enrichment	P-value	FDR
0098505	G-rich strand telomeric DNA binding	5	5.336982968	0.001208508	0.028476835
0043422	protein kinase B binding	6	5.336982968	0.000376501	0.010843234
0070530	K63-linked polyubiquitin modification-dependent protein binding	9	4.57455683	5.62867E–05	0.002026321
0043047	single-stranded telomeric DNA binding	5	4.447485807	0.003237143	0.060753992
0008641	ubiquitin-like modifier activating enzyme activity	5	4.269586375	0.002994627	0.060818804
0001206	transcriptional repressor activity, RNA polymerase II distal enhancer sequence-specific binding	4	3.881442159	0.014790908	0.19169017
0098847	sequence-specific single stranded DNA binding	5	3.812130692	0.006980259	0.11128415
0001222	transcription corepressor binding	5	3.812130692	0.006980259	0.11128415
0008536	Ran GTPase binding	12	3.767282095	3.39091E–05	0.001292537
0070410	co-SMAD binding	4	3.557988646	0.020567791	0.220954474
0042975	peroxisome proliferator activated receptor binding	4	3.557988646	0.020567791	0.220954474
0008093	cytoskeletal adaptor activity	4	3.557988646	0.020567791	0.220954474

Count’ stands for the number of DE genes associated with the listed GOID.

Table 5. Molecular function of upregulated m⁶A methylation.

GO.ID	Term	Count	Fold. Enrichment	P-value	FDR
0070840	dynein complex binding	5	9.061241351	0.000182742	0.032471479
0032356	oxidized DNA binding	2	8.336342043	0.022748279	0.737044237
0017166	vinculin binding	2	8.336342043	0.022748279	0.737044237
0098748	endocytic adaptor activity	2	7.578492766	0.027367914	0.738933667
0035615	clathrin adaptor activity	2	7.578492766	0.027367914	0.738933667
0000774	adenyl-nucleotide exchange factor activity	2	7.578492766	0.027367914	0.738933667
0043422	protein kinase B binding	2	6.946951702	0.0323281	0.805715735
0008093	cytoskeletal adaptor activity	2	6.946951702	0.0323281	0.805715735
0019789	SUMO transferase activity	3	6.946951702	0.008557751	0.553923826
0003841	1-acylglycerol-3-phosphate O-acyltransferase activity	3	6.581322665	0.009984023	0.553923826

‘Count’ stands for the number of DE genes associated with the listed GOID.

Table 6. KEGG pathway analysis of downregulated m⁶A methylation.

Pathway.ID	Definition	Count	P-value	FDR
hsa04064	NF-κB signaling pathway	20	8.41253E–05	0.007949932
hsa04620	Toll-like receptor signaling pathway	20	0.000536203	0.01729145
hsa04010	MAPK signaling pathway	36	0.001646395	0.03186836
hsa04621	NOD-like receptor signaling pathway	12	0.002627967	0.0394195
hsa04068	Foxo signaling pathway	19	0.01853444	0.1426183
hsa04350	TGFβ signaling pathway	13	0.02468286	0.1609752
hsa04330	Notch signaling pathway	8	0.04808801	0.2326839

‘Count’ stands for the number of differentially methylated genes associated with the listed Pathway.ID.

Table 7. KEGG pathway analysis of upregulated m⁶A methylation.

Pathway.ID	Term	Count	P-value	FDR
hsa04668	TNF signaling pathway	8	0.002409407	0.053333
hsa04064	NF-κB receptor signaling pathway	6	0.01457591	1
hsa04152	AMPK signaling pathway	7	0.01698368	1
hsa04310	Wnt signaling pathway	7	0.03271221	1
hsa04068	NOD-like receptor signaling pathway	4	0.03344811	1

‘Count’ stands for the number of differentially methylated genes associated with the listed Pathway.ID.

Figure 3. Conjoint analysis of differentially methylated and differentially expressed genes. (a) Scatter diagram showing the number of differentially expressed genes (fold change ≥ 2, P<0.05 and FPKM ≥ 0.1 at least one sample). (b) Venn diagram showing the differentially expressed genes containing DMMSs. (c) Ten hub genes identified by cytoHubba in the development of IPAH. (d) m⁶A enrichment of SMAD3, TP53, RPS27A and FoxO3 mRNA in the peripheral blood of CON group and IPAH group. (e) the mRNA expression levels of SMAD3, TP53, RPS27A and FoxO3 were measured in the peripheral blood of CON group and IPAH group. *P<0.05 (different from CON),**P<0.01 (different from CON), ****P<0.0001 (different from CON).

Figure 4. The identification of hypoxia-induced pulmonary arterial hypertension rat model. (a) RVSP and the ratio of RV/(LV + S) were measured in Normoxia group and Hypoxia group. (b) H&E staining was performed in the lung of Normoxia group and Hypoxia group. scale bar = 50 μm. Quantitative measurement of media thickness of pulmonary artery (c) and vessel muscularization (d). Data are present as means±SD (n = 6 rats in per group). *P<0.05 (different from Normoxia group), **P<0.01 (different from Normoxia group), ***P<0.001 (different from Normoxia group), ****P<0.001 (different from Normoxia group).

Figure 5. The identification of MCT-induced pulmonary arterial hypertension rat model. (a) RVSP and the ratio of RV/(LV + S) were measured in Control group and MCT group. (b) H&E staining was performed in the lung of Control group and MCT group. scale bar = 50 μm. Quantitative measurement of media thickness of pulmonary artery (c) and vessel muscularization (d). Data are present as means±SD (n = 3 rats per group). *P<0.05 (different from Control group), **P<0.01 (different from Control group), ***P<0.001 (different from Control group).

Figure 6. The expression of major m⁶A-related proteins in Hypoxia group and Normoxia group. (a) Western blot image of METTL3, FTO, ALKBH5, IGF2BP1, IGF2BP2, YTHDF1, YTHDF2. β-actin was used as an internal reference. (b) Statistical bar graph of METTL3, FTO, ALKBH5, IGF2BP1, IGF2BP2, YTHDF1, YTHDF2. Data are shown as means±SD (n = 6 rats per group). ^n.s.P>0.05 (different from Normoxia group), *P<0.05 (different from Normoxia group), **P<0.01 (different from Normoxia group), ***P<0.001 (different from Normoxia group).

Figure 7. Immunofluorescence staining of YTHDF1 (a), IGF2BP1 (b), ALKBH5 (c) in pulmonary vessels. α-SMA (red), objective protein (green) and nucleus were dyed with DAPI (blue), scale bar = 40 μm.

Figure 8. The expression of major m⁶A-related proteins in PASMCs exposed to hypoxia and PDGF-BB. (a) Western blot image of METTL3, FTO, ALKBH5, IGF2BP1, IGF2BP2, YTHDF1, YTHDF2, PCNA under hypoxia. β-actin was used as an internal reference. (b) Statistical bar graph of METTL3, FTO, ALKBH5, IGF2BP1, IGF2BP2, YTHDF1, YTHDF2, PCNA. (c) Western blot image of METTL3, YTHDF1, PCNA under PDGF-BB stimulation. β-actin was used as an internal reference. (d) Statistical bar graph of METTL3, YTHDF1, PCNA. Data are shown as means±SD (n = 3 each). ^n.s.P>0.05, *P<0.05, **P<0.01.

Figure 9. The protein expression of m⁶A-regulated genes in lungs of experimental PH rat models and PDGF-BB induced PASMCs. Representative immunoblots of TP53, SMAD3, FoxO3 and RPS27A protein expression in (a) lung tissues of hypoxia-induced pulmonary arterial hypertension rats, (b) lung tissues of MCT-induced pulmonary arterial hypertension rats and (c) PDGF-BB induced PASMCs (n = 3 per group) normalized to β-actin. Relative densitometric analysis of TP53, SMAD3, FoxO3 and RPS27A protein expression in (d) lung tissues of hypoxia-induced pulmonary arterial hypertension rats, (e) lung tissues of MCT-induced pulmonary arterial hypertension rats and (f) PDGF-BB induced PASMCs (n = 3 per group) normalized to β-actin. Data are shown as means±SD. ^n.s.P>0.05, *P<0.05.

Table

	CON	IPAH
Age, years	27.00 ± 3.54	31.25 ± 12.70
Gender, male/female	2/3	4/4
BMI, kg/m2	22.92 ± 4.94	23.15 ± 3.29
WHO-Fc, I-II/III-IV	/	3/5
NT-proBNP, pg/ml	/	1090.47 ± 602.77
6MWT, m	/	460.88 ± 141.37
SPAP, mmHg	/	93.25 ± 23.04
mPAP, mmHg	/	59.13 ± 11.86
PAWP, mmHg	/	7.50 ± 4.12
PVR, Wood units	/	20.31 ± 7.42

BMI, body mass index; WHO-Fc, Word Health Organization functional class; NT-proBNP, n-terminal natriuretic peptide type-B; 6MWT, six-minute walking distance; SPAP, systolic pulmonary pressure; mPAP, mean pulmonary artery pressure; PAWP, pulmonary arterial wedge pressure; PVR, pulmonary vascular resistance.

Table A2. Primers used for MeRIP-qPCR.

Gene Name		Sequence (5’−3’)
SMAD3	F R	CAGTCCAACCAGAAACACCC ACTGTCACAAGCGTCAGAGA
TP53	F R	CAGTCAGATCCTAGCGTCGA CAGAACGTTGTTTTCAGGAAGT
RPS27A	F R	TTTCAACAAACCAGAAGACAAGT GCACTACTAAAAGGAAATGGTGT
FoxO3	F R	TTGCGATGGTTTATGGGACG GCTAGGGGAAGGACGGTTAA

Table A3. Primers used for Qrt-PCR.

Name	Primer sequence (5′-3′)	GenBank accession number	Product size (bp)
TP53	CCCTCCTCAGCATCTTATCCG GCACAAACACGCACCTCAAA	NM_001126115.2	261
RPS27A	TGAACCCTCGGATACGATAGAAA GCACCACCACGAAGTCTCAAC	NM_001135592.2	180
SMAD3	CTACCAGTTGACCCGAATGTGC TCTGTCTCCTGTACTCCGCTCC	NM_001145102.1	74
FoxO3	CGGACAAACGGCTCACTCT GCTCTTGCCAGTTCCCTCATT	NM_001455.4	173
GAPDH	GGAAGCTTGTCATCAATGGAAATC TGATGACCCTTTTGGCTCCC	NM_002046	168