Peripheral leukocyte microRNAs as novel biomarkers for COPD.

COPD is a multifactorial disease caused by environmental determinants as well as genetic risk factors. The prevalence and mortality of COPD continue to increase, and underdiagnosis of COPD remains a critical issue. Previous reports investigated promising microRNAs (miRNAs) to reveal the molecular mechanism for the development of COPD; however, diagnostic and therapeutic markers for COPD have not yet been found. For this study, 20 representative COPD patients were separated into four groups based on increasing severity (A, B, C, and D) and compared to six healthy controls. Small RNA profiles of peripheral leukocytes were differentially expressed miRNAs (analyzed via next-generation sequencing) were validated via quantitative reverse transcriptase-polymerase chain reaction. Compared to healthy controls, 19 differentially expressed miRNAs were found in COPD patients. For all COPD groups, miR-3177-3p was downregulated, while 17 miRNAs were upregulated. Furthermore, the results revealed 21 differentially expressed miRNAs, of which miR-183-5p was continually downregulated from A to B to D. Between respective bronchodilator reversibility positive and negative groups of COPD different groups (A, B, C, and D), 10 miRNAs were differentially expressed, while miR-100-5p was upregulated in the negative groups. In conclusion, miR-106b-5p, miR-125a-5p, miR-183-5p, and miR-100-5p are central for the development of COPD. The severity of COPD was attenuated by miR-106b-5p, thus suggesting this miRNA as potential target for disease treatment.

Introduction

COPD is a major public health problem. Prevalence ranged from 0.2% to 37%, and mortality ranged from 3 to 111 deaths per 100,000 population. Furthermore, mortality steadily increased during the last 30–40 years.1 In 2020, COPD is projected to rank fifth worldwide in terms of disease burden and third in terms of mortality.2 In China, the reported prevalence of COPD varied from 5% to 13% in different provinces and cities, and the prevalence of COPD was significantly higher in men than in women (8.3%–18.9% vs 3.8%–7.1%, respectively). Furthermore, the prevalence and burden of COPD were predicted to increase throughout future decades due to continued exposure to COPD risk factors and the aging of the world’s population.2,3 A further central feature of COPD epidemiology is the elevated risk of underdiagnosis. Mainly for mild to moderate disease processes, 60%–85% of patients are estimated to remain undiagnosed.4,5 Due to the progressive nature of the disease, forced expiratory volume in 1 second (FEV1) typically declines over time.6 Clinical manifestations of COPD are highly variable, and comorbidities are significantly increased in patients, severely affecting survival quality and life.7 To date, the majority of treatments do not prevent the decline of lung function.8 Therefore, we assume to screen the biomarkers and treatment targets that are associated with severity of the disease, aiming to guide individual treatment for COPD.

The markers derived from lung tissue are not routinely available for clinical disease monitoring, whereas blood is readily accessible. Biomarkers in blood, especially inflammatory markers, are associated with exacerbations and mortality in larger COPD cohort studies.9 MicroRNAs (miRNAs) are a subtype of noncoding RNA and have been extensively studied primarily due to their vital role as post-transcriptional gene regulators. Owing to its high stability, strong specificity, high sensitivity, and detection easily in blood, miRNAs have been implicated in a variety of lung diseases,10 and Leidinger et al11 suggested that specific peripheral miRNA profiles could distinguish lung cancer from COPD. The authors12 revealed miR-3620-3p as a promising biomarker, separating COPD and asthma. Kara et al13 demonstrated that miR-29c and miR-126 are essential for the development of COPD.

Alveolar macrophages have been implicated in the pathophysiology of COPD.14 In our present study, we extracted peripheral blood leukocytes from COPD patients and healthy controls, utilizing miRNA screening via next-generation sequencing. We found four miRNAs (miR-106-b-5p, miR-125a-5p, miR-183-5p, and miR-100-5p) to be significant for the development of COPD, including disease transformation and evaluation of bronchodilator reversibility (BDR). Both central miRNAs (miR-106b-5p and miR-125-a-5p) were validated via quantitative real-time polymerase chain reaction (qRT-PCR) in a large sample set and acted as diagnostic and therapeutic markers. Hence, peripheral leukocyte miRNAs could act as novel biomarkers for COPD.

Materials and methods

Ethics statement

This study was approved by the Clinical Research Ethics Committee of the Shanxi Medical University. Written informed consent was obtained from participants prior to their enrollment in the study.

Patients

The study is a cross-sectional observation study. Between January 1, 2015, and November 31, 2015, a total of 1,809 cases were tested in the pulmonary function room of the Shanxi Dayi Hospital. According to the Global Initiative for Chronic Obstructive Lung Disease (GOLD) guideline,2 141 patients were diagnosed using spirometry as having COPD (postbronchodilator FEV1/forced vital capacity (FVC) ratio <70%), and using combined assessment, this led to allocation of 16 cases to A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and modified British Medical Research Council (mMRC) grade 0–1 or CAT score <10), 52 cases to B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or COPD Assessment Test (CAT) score ≥10), 25 cases to C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and 48 cases to D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10). We found BDR-positive patients in each group. Our inclusion criteria were as follows: 1) Meet the spirometry-diagnosed criteria. 2) Region of origin, age, and smoking status were consistent with other already diagnosed COPD patients. Our exclusion criteria were as follows: 1) Prior history of other respiratory diseases. 2) Complete blood cell counts were abnormal during the month prior to this study. 3) Having other comorbidities. We recruited 20 male COPD patients aged between 48 and 79 years (average age of 63 years) with a smoking history of at least 20 pack-years. The reason why we selected only male smokers was as follows: it is widely accepted that smoking is the key risk factor in the development of COPD; furthermore, in our diagnosed 141 COPD patients, 92% of patients had smoking history and 76% of patients were males who had a smoking history of 20 pack-years. Furthermore, according to a previous study,3 the ratio of COPD was significantly higher in men than in women (8.3%–18.9% vs 3.8%–7.1%, respectively).

Age-matched, healthy, male smokers (at least 20 packs per year of smoking history) were enrolled in this study, and they possessed normal lung function.

Consequently, the following groups were obtained: A, B, C, D, and control. Each group BDR-N and BDR-P patients: A-N, A-P; B-N, B-P; C-N, C-P; D-N, D-P.

During disease progression, the lung function and symptoms of patients varied frequently, which led to changes regarding the combined assessment and severity of COPD. We followed-up the 20 recruited COPD patients for 1 year with combined assessment per 3 month and clinically observed that the development of COPD severity often followed the progression pattern of A-B-D (progressing from A to B group to D group) or A-C-D (progressing from A group to C group to D group), with only few patients progressing step-by-step from B to C groups. Hence, we proposed two clinical experimental models for further comparison within A-B-D and A-C-D groups.

Sample

Using EDTA anticoagulation tubes, peripheral blood was collected. All samples were numbered, and leukocytes were immediately extracted.

Isolation of peripheral leukocytes

The protocol of leukocyte isolation was as follows:

Fresh whole blood was slowly added to a 15 mL centrifuge tube containing the lymphocyte separation solution or ficoll. Blood was placed over lymphocyte separated liquid, ensuring that both liquids did not mix and that the interface was clear.

The samples were centrifuged at 3,500 rpm for 20 min at 18°C, then delaminated, discarding the supernatant fluid.

The leukocyte layer was carefully separated with a pipette and removed to a new centrifuge tube.

Phosphate-buffered saline (10 mL) was added to the centrifuge tube, samples were centrifuged at 2,000 rpm for 10 min at 18°C, and the upper fluid was discarded.

Phosphate-buffered saline (5 mL) was added again to the centrifuge tube, followed by centrifugation at 2,000 rpm for 5 min at 18°C. The fluid was discarded and stored at −80°C in liquid nitrogen.

RNA extraction

Total RNA was isolated via the Trizol extraction method (Invitrogen, Carlsbad, CA, USA). RNA integrity and concentration were tested using an Agilent 2100 BioAnalyzer (Agilent Technologies, Santa Clara, CA, USA). We obtained an average of 4.4 µg of total RNA. The mean RNA integrity number was 8.9, and the ratio of 28S/18S was >1.5.

Small RNA library construction and Sequencing

The procedure followed a previously published protocol.15–18 According to the manufacturer’s protocol, we prepared the library as follows:

Filtering small RNA: We used 200 ng to 1 µg of RNA samples, then separated RNA segments of different sizes via sodium dodecyl sulfate-polyacrylamide gel electrophoresis (SDS-PAGE), selected 18–30 nt (14–30 ssRNA Ladder Marker, Takara,) stripe, and recycled.

Adaptor ligation: We prepared a connection 3′adaptor system (TruSeq Small RNA Sample Pre Kit, Illumina, San Diego, CA, USA) utilizing the following reaction conditions: 70°C for 2 min; 28°C for 1 h; 28°C for 15 min. Then, we added the 5′adaptor mix system with the following reaction conditions: 70°C for 2 min; 28°C for 1 h.

RT-PCR: We prepared the first strand master mix and super script II (Invitrogen) reverse transcription (reaction condition: 70°C for 2 min; 50°C for 1 h); this was followed by several rounds of PCR amplification with PCR primer cocktail and PCR master mix to enrich the cDNA fragments (reaction conditions: 98°C for 30 s; 11 cycles of 98°C for 10 s, 60°C for 30 s, 72°C for 15 s; 72°C for 10 min; 4°C hold).

Purify PCR products: Subsequently, the PCR products were purified with SDS-PAGE, dissolving the recycled products in ethidium bromide solution.

Validation of the library: The final library was quantified in two ways: determination of the average molecule length using the Agilent 2100 bioanalyzer instrument (Agilent DNA 1000 reagents) and quantifying the library via qRT-PCR (TaqMan Probe, Thermo Fisher Scientific, Waltham, MA, USA).

Sequencing libraries: The qualified libraries were amplified on cBot to generate the cluster on the flowcell. Furthermore, the amplified flowcell was sequenced single end on the HiSeq 2000 System (Illumina); reading length 50 was the most frequently used sequencing strategy.

Small RNA sequencing data analysis

Subsequent to filtering out low quality and adapter reads, the clean reads were aligned to the human reference genome (version 19) utilizing the short oligonucleotide alignment program 2.0 program (Microsoft Corporation, Redmond, WA, USA).19 Occurrences of each of the unique sequence reads were counted as sequence tags (the number of reads for each tag reflects the relative expression level), and only those small RNA sequences of 18 to 30 nt were retained for further analysis.16 The expression amount of every miRNA in each sample was calculated. The DEGseq method20 was used to identify differentially expressed miRNAs within the COPD groups (A/B/C/D and BDR-N and BDR-P groups) as well as between the groups (COPD vs control group) in accordance with the standard of selection: log2 ratio ≥1 and P≤0.001. The miRNAs target prediction was done via DIANA-microT21 (web server v5.0). To obtain dependable prediction, KEGG pathway22 analysis was carried out to identify significantly enriched gene ontology terms and pathways associated with the target genes using miRPath,23 with a cutoff of corrected P-value <0.05.

RT-PCR validation of large samples

The clinical sample size of COPD (A/B/C/D groups) was expanded to 20 cases per group. The differentially expressed miRNAs were validated via qRT-PCR. The procedure included three steps: miRNA isolation (miRNeasy Mini Kit, QIAGEN, Hilden, Germany), cDNA synthesis (miScript II RT Kit, QIAGEN), and qRT-PCR (miScript SYBR Green PCR Kit, QIAGEN). According to the manufacturer’s protocol and previously published protocol,24 a specific experimental process was carried out. The qRT-PCR was run on the Applied Biosystems 7500 Fast Real-Time PCR System (Foster City, CA, USA), using U6 as internal reference gene.

Statistical analysis

Small RNA sequencing data were filtered and screened in accordance with the standard (qvalue ≤0.001 and abs (log2(Y/X)) ≥1). Relative expression level profiles of selected miRNAs from qRT-PCR experiments were used as independent variable. The relative expression of miRNAs was calculated using Microsoft Excel (2007), comparing cycle threshold (Ct) values. The data are depicted with the software Graphpad prism 5 (Graphpad Software Inc., La Jolla, CA, USA).

Results

Between January 1, 2015, and November 31, 2015, 141 COPD patients were diagnosed and invited to participate in this study. The ratio of male to female was 8–2. The age distribution is 38–82 years. A number of male COPD patients were 112, and 88 were agreed to take part in the study. Sixty-eight male patients, with a smoking history of at least 20 pack-years, agreed to fill in the questionnaire survey and agreed for collection of peripheral blood samples. However, 48 patients had other comorbidities. In the end, we selected 20 stable COPD patients (Figure 1), 5 patients per group (A/B/C/D groups, Table 1). The control group consisted of six patients with a mean age of 61.5 years.

Figure 1

Flowchart of selected 20 COPD patients.

Abbreviation: pts, patients.

Table 1

General clinical technique of grouping COPD patients and controls (male, with a smoking history of at least 20 pack-years)

Group	Patient number	Age, years	Height (cm)	Weight (kg)	FEV1/FVC (%)	FEV1	CAT	mMRC
A	6	66	158	60	68.92	2.46 (108%)	4	0
	2	79	171	71.7	53.65	1.87 (64.7%)	6	1
	3	76	168	75	59.79	1.86 (73.5%)	6	1
	4	65	178	88	51.18	2.39 (72.3%)	0	1
	5	62	167	63	57.66	2.32 (79.2%)	9	1
B	20	72	165	63	45.45	1.42 (56.5%)	16	1
	7	63	177	101	68.7	2.7 (81.2%)	12	2
	8	63	172	64	49.39	2.29 (74.2%)	15	1
	9	61	160	67	58.03	2.13 (80.3%)	15	2
	10	56	158	69	65	2 (74.1%)	19	2
C	11	65	165	60	40.77	0.95 (38.7%)	6	1
	12	48	159	45	40.13	0.76 (25.3%)	8	1
	13	52	175	73	39.45	1.51 (43%)	7	1
	14	79	174	69	50.42	1.35 (50%)	9	1
	15	66	170	61	41.59	1 (34.1%)	3	1
D	16	57	169	78	41.38	1.08 (34.4%)	28	3
	17	65	167	52	40.19	1.01 (33.9%)	20	2
	18	67	165	49	50.77	0.87 (32.2%)	12	2
	19	53	163	67	54.85	1.52 (49.4%)	11	2
	1	66	169	76	48.34	1.36 (47.1%)	6	1
Control	13	60	165	54	74.33a	2.56 (92.3%)a	–	–
	19	57	163	60	89.15a	2.45 (84.7%)a	–	–
	21	68	164	67	80.88a	3.32 (103.8%)a	–	–
	22	60	178	59	84.11a	3.74 (104.5%)a	–	–
	24	60	175	89	72.28a	2.9 (84.9%)a	–	–
	27	58	163	65	78.6a	3.25 (113.4%)a	–	–

Notes:

FEV1/FVC, FEV1: ABCD group postbronchodilator. A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symp toms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symp toms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviations: CAT, COPD assessment test; FEV1, forced expiratory volume in 1 second; FVC, forced vital capacity; mMRC, modified British Medical Research Council; GOLD, Global Initiative for Chronic Obstructive Lung Disease.

Small RNA sequencing data description and annotation

The average of the raw sequencing data (NCBI SRX2502315) for the next-generation sequencing reads of 26 samples was 14.3 M. Following filtering out low quality reads and adapter reads, the mean value of final clean reads was 12.9 M. The total percentage of sequence tag length between 20 and 24 (nt) was 74, and the mean value of detected miRNAs in each sample was 890 (Table S1). The expression level of every miRNA was calculated for each sample.

Differentially expressed miRNAs between COPD patients and control

Compared with the healthy controls, there were 1,715 differentially expressed miRNAs, including 51 upregulated and 18 downregulated, in COPD patients.

The comparison of differentially expressed miRNAs among ABCD groups

Compared to control, all COPD groups had 19 differentially expressed miRNAs in common. Of these, miR-3177-3p was downregulated and 17 miRNAs were upregulated in all groups (Table 2). A Venn diagram presents the 19 miRNAs in detail in Figure 2.

Table 2

The expression profiles of selected miRNAs in ABCD group

miRNA_id	Log2 ratio (cases/controls)
	A	B	C	D
hsa-miR-3177-3p	−0.7266	−0.5033	−0.6503	−0.5537
hsa-miR-342-3p	1.5549	0.6783	0.9226	0.8622
hsa-miR-19b-3p	1.1053	0.5859	0.5666	0.6701
hsa-miR-17-5p	1.1454	0.6625	1.1415	0.734
hsa-miR-3613-5p	1.0345	0.7642	0.7434	1.0095
hsa-miR-93-5p	1.2347	0.7427	1.3273	0.9713
hsa-miR-29b-3p	1.3183	0.667	0.8991	0.7947
hsa-miR-125a-5p	1.4198	0.7671	0.5207	0.9874
hsa-miR-374a-5p	2.161	1.3767	1.2994	1.7128
hsa-miR-30b-5p	1.5249	0.5128	0.8939	0.5213
hsa-miR-29a-3p	1.0511	0.5396	0.6523	0.915
hsa-miR-106b-5p	0.9186	0.7551	0.7154	0.5881
hsa-miR-26b-5p	1.3534	0.6498	0.8861	1.3839
hsa-miR-223-3p	1.7663	0.584	1.273	0.8297
hsa-miR-142-5p	1.2301	0.7453	0.8384	1.1657
hsa-miR-142-3p	1.3586	0.937	0.6522	1.1217
hsa-miR-16-5p	1.27	0.8286	1.1467	1.0301
hsa-miR-19a-3p	0.7737	0.668	0.5932	0.671
hsa-miR-451a	0.7738	1.4441	1.5878	−2.3274

Notes: A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviation: miRNAs, microRNAs.

Figure 2

Intersections of differentially expressed genes among the A, B, C, D groups of COPD.

Notes: A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Compared with the healthy controls, there were 33 downregulated miRNAs and 87 upregulated miRNAs in A group; there were 23 downregulated miRNAs and 27 upregulated miRNAs in B group; there were 44 downregulated miRNAs and 80 upregulated miRNAs in C group; and there were 24 downregulated miRNAs and 65 upregulated miRNAs in D group. The profile of differentially expressed miRNAs is given in Figure 3.

Figure 3

The profile of differentially expressed miRNAs between the A, B, C, D groups and control group.

Notes: A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviation: miRNAs, microRNAs.

Surprisingly, we detected that the expression level of miR-106b-5p in A group was the highest and in D group the lowest. In other words, miR-106b-5p may negatively correlate with disease severity of COPD. miR-125a-5p was upregulated in all groups and was inversely correlated with ABC, but the expression level was elevated again in D-group. miR-451a was upregulated in ABC and positively correlated with ABC; however, expression was downregulated in D-group. Thus, these three miRNAs (miR-106b-5p, miR-125a-5p, and miR-451a) are noteworthy, especially the common upregulated miRNAs (miR-106b-5p and miR-125a-5p).

In our study, enrichment analysis was carried out to reveal relevant pathways. Via KEGG analysis, the representative signaling pathways of these targeted genes were found to be involved in the regulation of autophagy, toll-like receptor signaling pathway, non-small-cell lung cancer, and arrhythmogenic right ventricular cardiomyopathy (Table 3), all of which are involved as complications of COPD.

Table 3

Top 20 of KEGG pathways enrichment with miRNA target-predicted genes

Term	Count	P-value
hsa04140:Regulation of autophagy	13	6.24E-06
hsa04510:Focal adhesion	32	1.30E-04
hsa05410:Hypertrophic cardiomyopathy	17	6.56E-04
hsa04620:Toll-like receptor signaling pathway	18	0.001668041
hsa04622:RIG-I-like receptor signaling pathway	14	0.002677985
hsa04960:Aldosterone-regulated sodium reabsorption	9	0.011756963
hsa04650:Natural killer cell-mediated cytotoxicity	19	0.013152554
hsa04012:ErbB signaling pathway	14	0.015295956
hsa05223:Non-small-cell lung cancer	10	0.021165637
hsa05414:Dilated cardiomyopathy	14	0.023529701
hsa04623:Cytosolic DNA-sensing pathway	10	0.023644926
hsa04810:Regulation of actin cytoskeleton	26	0.024474467
hsa04512:Extracellular matrix–receptor interaction	13	0.026785185
hsa04722:Neurotrophin signaling pathway	17	0.028234004
hsa05200:Pathways in cancer	36	0.029159074
hsa04910:Insulin signaling pathway	18	0.02986907
hsa05412:Arrhythmogenic right ventricular cardiomyopathy	12	0.030460335
hsa04115:p53 signaling pathway	11	0.034589008
hsa05320:Autoimmune thyroid disease	9	0.040026109
hsa05218:Melanoma	11	0.044781153

Abbreviation: miRNA, microRNA.

Differentially expressed miRNAs in the progression of experimental groups

In the progression of ABD, we detected 21 common differentially expressed miRNAs, of which only miR-183-5p was found to be continually downregulated. In the progression of ACD, we found 33 common differentially expressed miRNAs (including miR-183-5p); however, there were no continuously upregulated or downregulated miRNAs (Table 4).

Table 4

The status of differentially expressed miRNAs among all groups

Control group	Comparison group	UP	Down
Control	A	71	39
Control	B	35	14
Control	C	68	42
Control	D	48	16
A	B	41	38
B	D	26	20
A	C	34	25
C	D	38	36
A_P	A_N	40	30
B_P	B_N	5	81
C_P	C_N	89	42
D_P	D_N	73	110

Notes: A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviation: miRNAs, microRNAs.

Detected miRNAs between BDR-P and BDR-N groups

Within each group, an average of 51.75 miRNAs were upregulated and 65.75 miRNAs were downregulated in BDR-N groups, compared to BDR-P groups. Respective BDR-P and -N groups of ABCD had 10 differentially expressed known miRNAs in common. Only miR-100-5p was upregulated in all negative groups. These 10 miRNA sets are presented in detail in the Venn diagram (Figure 4). Consequently, miR-100-5p may be used as a biological marker for the evaluation of bronchial obstruction reversibility.

Figure 4

Intersections of differentially expressed genes between BDR positive vs negative groups.

Notes: A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviation: BDR, bronchodilator reversibility.

Discussion

COPD is a multifactorial disorder caused by environmental determinants, most commonly cigarette smoking, and genetic risk factors.6 Unfortunately, few studies have successfully unraveled the genetic predisposition and pathogenetic mechanisms for COPD.25 Dysregulation of miRNAs has been reported to be involved in the pathogenesis of pulmonary diseases, including COPD.26 In previous studies, miR-146a27,28 has been considered closely related to inflammatory responses of COPD, via regulation of its targets COX2, TRAF6, and IRAK1 and involvement in NF-kB and toll-like receptor signaling. Moreover, Wang et al29 observed that miR-34 could inhibit the expression of CD80 and CD86 as well as the secretion of IFN-α to restrain dendritic cell dysfunction in COPD patients. Furthermore, miR-34c has been reported to be correlated with the severity of emphysema in COPD by modulating its predicted targets.30 Consequently, via regulation of mRNA degradation and translational repression, miRNAs affect different biological process of COPD, such as inflammatory responses, emphysematous lesions, and tissue repair.31

It was reported that 234 miRNAs were detected by miRNA microarray in current smokers with COPD, and miR-145-5p could be a promising biomarker of COPD.12 Because of the different characteristics of selected COPD patients and the different method used, our results might be different. By next-generation sequencing, we found 1,715 differentially expressed miRNAs, and 7 miRNAs (hsa-miR-3177-3p, hsa-miR-19b-3p, hsa-miR-3613-5p, hsa-miR-93-5p, hsa-miR-29a-3p, hsa-miR-142-5p, and hsa-miR-19a-3p) might be novel in COPD. 12 miRNAs (miR-342-3p,32 miR-17-5p,13 miR-16-5p,13 miR-29b-3p,32,33 miR-26b-5p,33 miR-125a-5p,13,34 miR-106b-5p,33 miR-223-3p,35 miR-142-3p,36 miR-30b-5p,36 miR-374a-5p,37 and miR-451a37) had been reported in previous COPD studies. It has been reported that miR-29b38 and miR-106b39 are involved in TGF-β1 signaling. Previous studies33–39 employed the spirometric-based staging of the degree of airflow obstruction and examined airway epithelial cells, serum, plasma, and lung tissue using microarray or RT-PCR.

Furthermore, we found that miR-106b-5p and miR-125a-5p might be potential markers of COPD. And it was an exciting discovery that miR-106b-5p was negatively correlated to severity of COPD; this was validated via qRT-PCR in large samples (Table 5 and Figure 5A). Recent studies about miR-106b-5p40–43 and miR-125a-5p44–48 are mainly related to other diseases, and few studies relate to COPD. Soeda et al33 reported that miR-106b plasma levels are evidently downregulated in COPD and are negatively correlated with disease duration. In this study, miR-106b-5p was upregulated. The cause for this inconsistency might due to different samples used, different technology followed, and different regions studied. However, the correlation was consistent; the previous study merely did not classify the patients with COPD and clarify regulatory mechanisms. We speculate the upregulation of miR-106b-5p in our study might be involved in signaling pathways, including TGF-β signal pathway, cell cycle, p53 signal pathway, apoptosis, and lung cancer. These pathways have been shown to be related to COPD development. Kosacka et al49 assessed serum levels of apoptosis-related markers (sFasL, TNFα, p53, and bcl-2) in COPD patients, the study described the vital role of apoptosis in the pathogenesis of COPD. Morty et al50 covered the impact of TGF-β signaling on COPD. Cigarette smoke may affect the reparative potential of lung fibroblasts, altering the progression of the cell cycle to S phase, which contributes to the development of pulmonary emphysema.51 And, it is widely known that lung cancer is one of the comorbidities of COPD. Previous studies13,34,52,53 have suggested miR-125 to be an important regulatory factor related to inflammation and further proposed that miR-125 might be connected with COPD and involved in disease pathogenesis; however, the authors did not expound the specific relationship. We found miR-125a-5p to be upregulated in COPD in all ABCD groups, which was confirmed via qRT-PCR in large samples (Table 5 and Figure 5B), and it might also be involved in the Wnt signaling pathway via its targeted gene SMAD4, thus affecting the process of COPD. Consequently, we proposed miR-106b-5p and miR-125a-5p to be potential targets, especially the graded marker miR-106b-5p due to its relevancy to COPD severity.

Table 5

Relative expression (2−ΔΔCT) of two validated miRNAs in different groups (mean ± SD)

miRNA ID	Control	A	B	C	D
miR-106b	1.03±0.09	557.54±35.13	495.96±77.12	268.67±36.26	157.17±27.94
miR-125a	1.01±0.17	28.54±0.85	23.2±4.55	19.8±3.62	27.33±2.41

Notes: A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviations: miRNAs, microRNAs; SD, standard deviation.

Figure 5

Box-plot diagrams of relative expression of 2 validated miRNAs in control and A, B, C, D groups.

Notes: (A) miR-106b-5p; (B) miR-125a-5p. Relative miRNAs expression level is shown in y axis. A-group (low risk, less symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade 0–1 or CAT score <10), B-group (low risk, more symptoms GOLD 1–2 mild or moderate airflow limitation, 0–1 exacerbations per year, and mMRC grade ≥2 or CAT score ≥10), C-group (high risk, less symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year, and/or ≥1 hospitalized exacerbation per year, and mMRC grade 0–1 or CAT score <10), and D-group (high risk, more symptoms GOLD 3–4 severe or very severe airflow limitation, and/or ≥2 exacerbations per year/≥1 hospitalized exacerbation per year, and mMRC grade ≥2 or CAT score ≥10).

Abbreviation: miRNAs, microRNAs.

miR-183-5p might play a central role in the development of COPD, especially for the evaluation of disease progression (from A to B to D). This seems to explain the clinical phenomenon of different developmental directions for disease severity in patients with COPD. It was reported54 that miR-183-5p was related to immune responses, and so while chronic inflammatory responses were the significant pathophysiology in COPD, miR-183-5p might affect development of COPD by regulating inflammatory mediators.

Moreover, miR-100-5p was commonly upregulated in the BDR negative groups. Akbas et al55 reported miR-100 to be related to COPD. However, there has currently been no correlative reference about miRNAs of BDR. Zanini et al56 reported bronchial reversibility as a more frequently present characteristic in COPD patients, although it could be one single phenotype of COPD. The ECLIPSE study showed that BDR was related to a decline of FEV1.57 Rhee58 speculated on whether pure patients with COPD with BDR can be considered having asthma–COPD overlap syndrome (ACOS). Consequently, we specifically analyzed the sequencing data and found promising viability for miR-100-5p as biological marker to evaluate BDR in COPD, further distinguishing between COPD, asthma, and ACOS.

The limitations of this study include that it focused on male smokers with COPD and did not account for variability in the BD response. Further studies are needed to confirm the role of miR-183-5p and miR-100-5p in l patient with long-term follow-up. Furthermore, the regulatory mechanism of miR-106b-5p and miR-125a-5p might need to be clarified in subsequent research.

Conclusion

Here, we have listed the evidence for miRNA as viable diagnostic biomarkers for COPD owing to its function in regulating apoptosis, NF-kB, Toll-like receptor, TGF-β, or Wnt signaling pathways, thereby affecting COPD progression. In particular, we found that miR-106b-5p might become a graded marker for the evaluation of COPD severity. Therefore, the expression of miRNA and its targeted mRNA may provide a promising possibility for therapeutic intervention at an early stage. However, it is still a challenge to fully expose the function and molecular interactions of miRNA in cell and animal models, which would contribute to formulate specific diagnosis as well as treatment strategies.

Table S1

The description of small RNA sequencing data

Sample name	BDR	Raw tags	Clean tags	Ratio of tags 20–24 nt	Ratio of miRNA	Species of miRNA
A1_P	+	14,595,822	13,556,065	78.26	75.29	975
A2_P	+	13,969,245	12,882,282	52.61	40.58	808
A3_N	−	14,587,837	12,614,597	74.84	69.65	933
A4_N	−	14,998,340	14,216,530	78.04	66.92	957
A5_N	−	14,266,859	12,451,088	79.01	64.59	837
B1_P	+	13,784,913	12,547,380	78.92	77.96	983
B2_P	+	14,307,584	13,717,056	67.22	53.46	843
B3_P	+	13,917,470	11,579,946	73.71	64.24	867
B4_N	−	13,842,456	12,824,493	79.92	69.39	921
B5_N	−	14,186,522	13,286,444	84.48	66.91	864
C1_P	+	14,845,131	13,290,576	70.66	64.25	897
C2_N	−	13,506,361	12,403,551	78.39	76.31	981
C3_N	−	13,787,867	12,590,034	81.54	75.85	962
C4_N	−	13,527,597	11,730,930	76.3	62.74	833
C5_N	−	14,542,514	13,150,839	81.68	66.54	835
D1_P	+	14,646,677	13,358,418	75.5	66.56	991
D2_N	−	14,419,732	12,719,972	75.4	72.72	931
D3_N	−	14,996,182	13,514,779	80.76	66.99	826
D4_N	−	14,656,089	12,247,648	68.13	60.07	819
D5_N	−	14,800,956	13,193,073	78.64	67.28	862
Control_A2	NA	13,682,036	12,429,610	77.12	67.91	968
Control_A3	NA	14,897,096	13,076,502	77.11	63.35	869
Control_A4	NA	14,703,228	13,138,111	79.09	67.85	909
Control_B2	NA	13,720,996	11,777,107	44.65	35.35	814
Control_B3	NA	13,555,103	12,133,258	80.14	63.37	840
Control_B4	NA	14,163,024	13,768,274	54.68	42.91	818

Abbreviations: BDR, bronchodilator reversibility; miRNA, microRNA; NA, not applicable.