Lung IL-13 gene signatures are associated with raised tissue eosinophils in COPD
Chronic Obstructive Pulmonary Disease (COPD) is now the third leading cause of global mortality . Current treatments remain inadequate, with only modest impacts on morbidity and mortality. The prospect of stratifying patients for certain therapies is attractive and previous reports highlight blood eosinophil counts as a marker of steroid efficacy in COPD (reviewed in ). Notably, GOLD recommends using blood eosinophil counts ≥ 300 cells/µl as a biomarker to be used alongside clinical assessment to choose those most likely to benefit from inhaled corticosteroids . Those patients who have blood eosinophil counts ≤ 100 cells/ µl are less likely to benefit from inhaled corticosteroid therapy . The potential of more targeted therapeutics, such as antibodies targeting IL-5 or its receptor, requires further understanding after clinical trials using eosinophil-specific therapies showed disappointing results . However, a positive large scale clinical trial of dupilumab (a monoclonal antibody targeting the IL-4 and IL-13 shared receptor component) has demonstrated reductions in exacerbation frequency . Whilst peripheral blood eosinophil counts appear useful to select patients for dupliumab treatment overall, it remains uncertain how useful this approach was to determine individual clinical response. Indeed the efficacy of this intervention was greater in subjects further stratified by a raised FeNO level, suggesting direct measures of pulmonary inflammation may offer additional benefits. Previous reports suggest that COPD patients show increased numbers of lung eosinophils compared with healthy controls, even when allergy and asthma are excluded [ – ]. Furthermore, recent work has demonstrated an association between raised blood eosinophils and the development of obstructive lung disease . Eosinophilic inflammation may, therefore, play an important role in COPD immune dysregulation. However, the presence and proportion of lung eosinophils varies considerably across COPD patients [ , , ]. This may reflect differences in disease activity, as there is also an association of eosinophils with COPD exacerbations . Sputum has been used to characterise lung eosinophilia in patients with COPD. However, the wider clinical utility of this approach is limited as sputum analysis is largely restricted to research centres . Therefore, blood eosinophils are used as a surrogate marker to enable patient phenotyping more broadly. Whilst there appear to be strong correlations between blood and sputum eosinophils , the association of blood eosinophils with eosinophils in lung tissue and bronchoalveolar lavage (BAL) is unclear. Furthermore, the impact of eosinophils in these different lung compartments on disease measures is not well defined . Understanding the nature of eosinophilic inflammation in the COPD lung itself is a key first step to delivering a step-change in treatment outcomes which is already being achieved in asthma . To address these questions, we investigated the distribution of eosinophils in blood, BAL, and bronchial biopsies from deeply phenotyped COPD patients. Furthermore, we characterised gene expression differences between blood, bronchial biopsies and epithelial brushings and their associations with the presence of eosinophils in these different compartments, aiming to understand the role and possible drivers of increased eosinophils in COPD. 

Subjects We recruited ex-smoking subjects with mild or moderate COPD (as per GOLD) (n = 31) and healthy ex-smoking volunteers (HV-ES) (n = 20) as the most relevant control group, all with ≥ 10-pack year history and had stopped smoking ≥ 6 months prior to enrolment, as part of the MICAII study (Fig. ). Patients with a history of asthma or atopy were excluded from the cohort. Additional details about this cohort have been reported previously [ – ]. For additional analyses, COPD subjects were split dependent on exacerbation frequency. COPD subjects were classified as either infrequent exacerbators (IE) (≤ 1 exacerbation in the preceding 12 months before enrolment) or frequent exacerbators (FE) (≥ 2 exacerbations in the preceding 12 months before enrolment). All subjects gave written informed consent, and the study was approved by National Research Ethics Service South Central ethical standards – Hampshire A and Oxford C Committees (LREC no: 15/SC/0528). Sampling was undertaken using fibreoptic bronchoscopy and epithelial brushings, bronchial biopsies and BAL were recovered and processed as previously described [ , , ]. RNA isolation and sequencing Total RNA was extracted from epithelial brushing and bronchial biopsy samples using the AllPrep DNA/RNA/miRNA Universal Kit (Qiagen), whilst RNA was extracted from whole blood using the PAXgene Blood RNA Kit (Qiagen). The quantity and quality of RNA samples were determined using the standard RNA analyzer kit on a 96-channel Fragment analyzer (Agilent Technologies). Extracted samples with a yield concentration > 25 ng/µl total RNA, and a DV 200 value (percentage of RNA fragments > 200nucleotides) >  = 30% were deemed to be of sufficient quantity and quality for TotalRNA-seq analysis. Samples were diluted to 25 ng/µl using a Tecan Fluent liquid handling automation system (Tecan). Library preparation was done in four separate runs, one 96 well plate per run. The Kapa RNA HyperPrep Kit with RiboErase (HMR) was used for reverse transcription, generation of double stranded cDNA and subsequent library preparation and indexing to facilitate multiplexing (Roche), all of which was performed through automation on a Tecan fluent. The libraries were quantified with the 96-channel Fragment Analyzer using the standard sensitivity next generation sequencing (NGS) kit (Agilent Technologies). Samples from each preparation plate were pooled and the final pools (4 in total) were quantified using a Qubit instrument for concentration determination with the DNA High Sensitivity kit (ThermoFisher Scientific). Fragment size was determined using the Fragment Analyzer, standard sensitivity NGS kit (Agilent Technologies). Three of four library pools were further diluted to 1 nM and sequenced on a NovaSeq 6000 (Illumina) using NovaSeq 6000 S4 Reagent Kit, 2 × 76 cycles. The remaining library pool was diluted to 1.9 nM and sequenced on NovaSeq 6000 (Illumina) using 2 NovaSeq 6000 SP S1 Reagent Kits, 2 × 51 cyclers. Average reads per sample were 53.3 million. RNASeq analysis Fastq files from 307 paired-end sequencing libraries generated from 120 epithelial brushings, 125 bronchial biopsies, and 62 blood samples were collected and read quality for all libraries was accessed using FastQC (v0.11.9) , Qualimap (v2.2.2d)  and samtools stats (v1.15) . Quality control (QC) metrics for Qualimap were based on a STAR (v2.7.10a)  alignment against the human genome (GRCh38, Gencode v43). Next, QC metrics were summarized using MultiQC (v1.12) . Two libraries were excluded; one due to a low mapping rate (57% vs [79%−97%]) and another due to low sequencing throughput (210 k reads vs [20 M-86 M]), leaving 118 epithelial brushings, 125 bronchial biopsies, and 62 blood samples for further analysis. Sequencing adapters were then trimmed from the remaining libraries using NGmerge (v0.3) . A human transcriptome index consisting of cDNA and ncRNA entries from Gencode (v43) was generated and reads were mapped to the index using Salmon (v1.7.0) . The bioinformatics workflow was organized using Nextflow workflow management system (v20.10)  and Bioconda software management tool . Differential gene expression were assessed with DESeq2 (v 1.34.0) , using ashr (v2.2_54)  for fold change shrinkage, all in R (v4.1.3) . Estimated counts from Salmon were used as input for DESeq2 (v1.34.0) using tximport (v1.22.0)  in R (v4.1.3). In the models used to assess differential expression between subject groups, effects from gender and a technical batch-effect (library prep plate) were taken into account. To ensure that the identified DEGs had robust and biologically meaningful expression, we applied a threshold requiring a median expression value of > 0.5 log2TPM in at least one comparison group, alongside an adjusted p-value < 0.05 determined by the Benjamini–Hochberg multiple testing correction method. DEGs were visualised using the EnhancedVolcano package (v 1.20.0) and DEG overlaps visualised using ggvenn (0.1.10). The clusterProfiler  package (v4.10.0) was used to perform Over Representation Analysis (ORA) and Gene Set Enrichment Analysis (GSEA) on the DEG lists. ORA was performed using enrichGO() for gene ontology (GO) categories Biological Process (BP), Molecular Function (MF) and Cellular Component (CC), enrichKEGG() for KEGG pathways and enrichPathway() for REACTOME pathways. GO categories and KEGG pathways were obtained from Bioconductor org package org.Hs.eg.db (v 3.18.0) , and REACTOME pathways were obtained from ReactomePA (version 1.46.0) . All genes tested for DEG analysis were used as the background gene set, min gene set size was 5 and max gene set size was 500. The Benjamini–Hochberg multiple testing correction method was applied, and significant terms/pathways were filtered according to an adjusted p-value < 0.05. The top 5 significant terms/pathways for the ORA were visualised in R using ggplot2 (v 3.4.4). GSEA was performed using GSEA() on the brushings dataset to further investigate the presence of IL-13 signalling in this sample type. All genes were ranked taking into account both the log2FC (ashr shrunken) and p-value. The INTERLEUKIN_4_AND_INTERLEUKIN_13_SIGNALING curated C2 pathway (R-HSA-6785807) was downloaded using msigdbr (v 7.5.1). Result visualisation was performed using gseaplot2() from the enrichplot package (v 1.22.0). Proteomics Proteins in BAL supernatants were analysed using LC MS/MS as previously described . LC–MS/MS analysis of TMT labelled peptides was carried out on a Q Exactive HF-X (Thermo Fisher Scientific) mass spectrometer interfaced with a Dionex 3000 RSLCnano system. Peptides were captured on a 2 cm × 75 µm C18 trap column (ReproSil-Pur 120 C18-AQ 7um) and samples were separated on a monolithic column (50 cm, cut from a 2 m long column, 100 µm ID, GL Sciences Inc. USA) using a gradient of solvent A (0.2% formic acid) and solvent B (0.2% formic acid in 90% acetonitrile). Peptides were separated using a 90 min gradient of solvent B as follows: 4% to 16.5%B in 2.5—52.5 min; 33.5% B in 73 min followed by a stay at 98% B for 3 min and re-equilibration at 2% B. A flowrate of 0.7 µL/min was used. Peptides were sprayed in an electrospray ionization (ESI) source using a stainless-steel emitter with 2 kV at a capillary temperature of 275 °C. A full-scan MS spectrum was collected at 60,000 resolution at m/z of 200 and scanned at 350–1200 m/z with automatic gain control (AGC) of 3E6. The top 10 precursors were selected, and an MS/MS scan was obtained at 45,000 resolution with 50 ms injection time, isolation window of 0.9 m/z with offset 0.1 m/z, normalized collision energy (NCE) of 29. For MS2, minimum AGC target was set to 1.7E4. Dynamic exclusion duration was set to 15 s. The fixed first mass was set to 100 m/z. Charge state exclusion was set to ignore unassigned, 1, and 7 and greater charges. For internal mass calibration, lock mass of 371.10124 m/z was used. Mass spectrometry data was analysed using Proteome Discoverer 2.3 (Thermo Fisher Scientific) software with search engine Mascot (version 2.6.0). Data was searched using latest UniProt Human protein database. Unfragmented precursor and TMT reporter ions were removed using a non-fragment filter in the PD 2.3 workflow. Search parameters included 3 missed cleavages for trypsin, oxidation (M) and deamidation (N, Q) as variable modifications. Tandem label (229.163 Da) at N-terminus and lysine residues and carbamidomethylation on cysteine residues were set as fixed modifications. The mass tolerances on precursor and fragment masses were set at 20 ppm and 0.05 Da, respectively, for MS2 analysis. Consensus step in PD2.3 included several nodes for spectrum, peptide and protein grouping and FDR calculation. Reporter ions for TMT labelled peptides were quantified using the PD quantitation node and peak integration tolerance was set at 20 ppm by considering most confident centroid peaks. Signal to noise values were calculated in addition to measurement of intensities of the TMT reporter ion for peptide and protein quantitation. The intensities were normalized by total peptide amount in PD 2.3. To account for protein input, the global quantitative proteome data was reviewed before normalization and no samples showed an unexpected pattern of distribution. Albumin and haemoglobin abundances were not significantly different between sub-cohorts. Further normalization of the data across all samples was carried out using Reporter Ion Quantitation in Proteome Discoverer, which calculates the total sum of the abundance values for each TMT channel over all peptides identified within a file. The channel with the highest total abundances served as a reference for correcting abundances across the remaining channels by a constant factor. Immunofluorescence Biopsies were only available from 19 COPD and 8 HV-ES for imaging analysis as samples for RNASeq analysis was prioritised (Fig. ). Biopsies for immunofluorescence were fixed in 4% (w/v) formaldehyde and embedded in paraffin wax, as previously prescribed . Ten serial sections of 5 μm thickness were obtained using a Leica RM2135 microtome (Leica Biosystems, Germany). The sections were then mounted onto APES (3-aminopropyltriethoxysilane) coated microscope slides. (CellPath, Powys). To deparaffinise the tissue, microscope slides containing the tissue were placed into slide racks and into clearene (Leica Biosystems) twice with ten minutes each in a clearene tank. The slides were rehydrated by soaking them in graded alcohols solutions of 100% ethanol, 75% ethanol and 70% ethanol, 5 min each. To perform H&E stains, histology slides were placed into Mayer’s Haematoxylin (CellPath, Powys) solution for 5 min, then placed under running tap water for 5 min. The slides were then placed into eosin stain (CellPath, Powys) for further 5 min. Slides were then briefly placed into 100% ethanol, after which they were dehydrated by taking them into 95% ethanol and 100% ethanol (1 min on each). After which they were then placed into clearene three times (3 min each). The slides were mounted in XFT mounting medium (CellPath, Powys) and cover slipped. Antigen retrieval was performed by pipetting 5% (v/v) pronase solution onto the tissue. These slides were left to incubate at room temperature for 10 min, and then rinsed in 1 × PBS. The tissue slides were incubated at room temperature for a further 60 min with blocking solution (1 × PBS + 1% bovine serum albumin & 2% foetal calf serum). Excess blocking buffer was removed, primary antibodies were added and incubated overnight at 4 °C. Primary antibodies were added to all slides except the negative controls. 0.00125 mg/ml mouse anti-EG2 (Diagnostic Development, Uppsala, Sweden) or 0.005 mg/ml rabbit anti-IL-5Rα (ThermoFisher, UK) primary antibodies incubated overnight at 4 °C. Tissue was also incubated overnight with just the blocking buffer. After overnight incubation, the slides were washed in PBS and then 0.002 mg/ml AlexaFluor647 goat anti-mouse or 0.004 mg/ml AlexaFluor647 goat anti-rabbit secondary antibodies (both ThermoFisher) were incubated at room temperature for 1 h, after which the slides were washed in PBS. The tissue was then stained with DAPI (Roche, Germany) at room temperature for 10 min, after which it was washed off. The slides were mounted onto coverslips with mowiol (Sigma-Aldrich, UK). Statistics Analysis of two groups was performed using a Mann–Whitney U test. Fishers exact test was used for categorical data (GraphPad Prism v9, GraphPad Software Inc., San Diego, USA). Associations were assessed using Spearman's correlation with rho and p values presented. Results were considered significant if p < 0.05. 

We recruited ex-smoking subjects with mild or moderate COPD (as per GOLD) (n = 31) and healthy ex-smoking volunteers (HV-ES) (n = 20) as the most relevant control group, all with ≥ 10-pack year history and had stopped smoking ≥ 6 months prior to enrolment, as part of the MICAII study (Fig. ). Patients with a history of asthma or atopy were excluded from the cohort. Additional details about this cohort have been reported previously [ – ]. For additional analyses, COPD subjects were split dependent on exacerbation frequency. COPD subjects were classified as either infrequent exacerbators (IE) (≤ 1 exacerbation in the preceding 12 months before enrolment) or frequent exacerbators (FE) (≥ 2 exacerbations in the preceding 12 months before enrolment). All subjects gave written informed consent, and the study was approved by National Research Ethics Service South Central ethical standards – Hampshire A and Oxford C Committees (LREC no: 15/SC/0528). Sampling was undertaken using fibreoptic bronchoscopy and epithelial brushings, bronchial biopsies and BAL were recovered and processed as previously described [ , , ]. 

Total RNA was extracted from epithelial brushing and bronchial biopsy samples using the AllPrep DNA/RNA/miRNA Universal Kit (Qiagen), whilst RNA was extracted from whole blood using the PAXgene Blood RNA Kit (Qiagen). The quantity and quality of RNA samples were determined using the standard RNA analyzer kit on a 96-channel Fragment analyzer (Agilent Technologies). Extracted samples with a yield concentration > 25 ng/µl total RNA, and a DV 200 value (percentage of RNA fragments > 200nucleotides) >  = 30% were deemed to be of sufficient quantity and quality for TotalRNA-seq analysis. Samples were diluted to 25 ng/µl using a Tecan Fluent liquid handling automation system (Tecan). Library preparation was done in four separate runs, one 96 well plate per run. The Kapa RNA HyperPrep Kit with RiboErase (HMR) was used for reverse transcription, generation of double stranded cDNA and subsequent library preparation and indexing to facilitate multiplexing (Roche), all of which was performed through automation on a Tecan fluent. The libraries were quantified with the 96-channel Fragment Analyzer using the standard sensitivity next generation sequencing (NGS) kit (Agilent Technologies). Samples from each preparation plate were pooled and the final pools (4 in total) were quantified using a Qubit instrument for concentration determination with the DNA High Sensitivity kit (ThermoFisher Scientific). Fragment size was determined using the Fragment Analyzer, standard sensitivity NGS kit (Agilent Technologies). Three of four library pools were further diluted to 1 nM and sequenced on a NovaSeq 6000 (Illumina) using NovaSeq 6000 S4 Reagent Kit, 2 × 76 cycles. The remaining library pool was diluted to 1.9 nM and sequenced on NovaSeq 6000 (Illumina) using 2 NovaSeq 6000 SP S1 Reagent Kits, 2 × 51 cyclers. Average reads per sample were 53.3 million. 

Fastq files from 307 paired-end sequencing libraries generated from 120 epithelial brushings, 125 bronchial biopsies, and 62 blood samples were collected and read quality for all libraries was accessed using FastQC (v0.11.9) , Qualimap (v2.2.2d)  and samtools stats (v1.15) . Quality control (QC) metrics for Qualimap were based on a STAR (v2.7.10a)  alignment against the human genome (GRCh38, Gencode v43). Next, QC metrics were summarized using MultiQC (v1.12) . Two libraries were excluded; one due to a low mapping rate (57% vs [79%−97%]) and another due to low sequencing throughput (210 k reads vs [20 M-86 M]), leaving 118 epithelial brushings, 125 bronchial biopsies, and 62 blood samples for further analysis. Sequencing adapters were then trimmed from the remaining libraries using NGmerge (v0.3) . A human transcriptome index consisting of cDNA and ncRNA entries from Gencode (v43) was generated and reads were mapped to the index using Salmon (v1.7.0) . The bioinformatics workflow was organized using Nextflow workflow management system (v20.10)  and Bioconda software management tool . Differential gene expression were assessed with DESeq2 (v 1.34.0) , using ashr (v2.2_54)  for fold change shrinkage, all in R (v4.1.3) . Estimated counts from Salmon were used as input for DESeq2 (v1.34.0) using tximport (v1.22.0)  in R (v4.1.3). In the models used to assess differential expression between subject groups, effects from gender and a technical batch-effect (library prep plate) were taken into account. To ensure that the identified DEGs had robust and biologically meaningful expression, we applied a threshold requiring a median expression value of > 0.5 log2TPM in at least one comparison group, alongside an adjusted p-value < 0.05 determined by the Benjamini–Hochberg multiple testing correction method. DEGs were visualised using the EnhancedVolcano package (v 1.20.0) and DEG overlaps visualised using ggvenn (0.1.10). The clusterProfiler  package (v4.10.0) was used to perform Over Representation Analysis (ORA) and Gene Set Enrichment Analysis (GSEA) on the DEG lists. ORA was performed using enrichGO() for gene ontology (GO) categories Biological Process (BP), Molecular Function (MF) and Cellular Component (CC), enrichKEGG() for KEGG pathways and enrichPathway() for REACTOME pathways. GO categories and KEGG pathways were obtained from Bioconductor org package org.Hs.eg.db (v 3.18.0) , and REACTOME pathways were obtained from ReactomePA (version 1.46.0) . All genes tested for DEG analysis were used as the background gene set, min gene set size was 5 and max gene set size was 500. The Benjamini–Hochberg multiple testing correction method was applied, and significant terms/pathways were filtered according to an adjusted p-value < 0.05. The top 5 significant terms/pathways for the ORA were visualised in R using ggplot2 (v 3.4.4). GSEA was performed using GSEA() on the brushings dataset to further investigate the presence of IL-13 signalling in this sample type. All genes were ranked taking into account both the log2FC (ashr shrunken) and p-value. The INTERLEUKIN_4_AND_INTERLEUKIN_13_SIGNALING curated C2 pathway (R-HSA-6785807) was downloaded using msigdbr (v 7.5.1). Result visualisation was performed using gseaplot2() from the enrichplot package (v 1.22.0). 

Proteins in BAL supernatants were analysed using LC MS/MS as previously described . LC–MS/MS analysis of TMT labelled peptides was carried out on a Q Exactive HF-X (Thermo Fisher Scientific) mass spectrometer interfaced with a Dionex 3000 RSLCnano system. Peptides were captured on a 2 cm × 75 µm C18 trap column (ReproSil-Pur 120 C18-AQ 7um) and samples were separated on a monolithic column (50 cm, cut from a 2 m long column, 100 µm ID, GL Sciences Inc. USA) using a gradient of solvent A (0.2% formic acid) and solvent B (0.2% formic acid in 90% acetonitrile). Peptides were separated using a 90 min gradient of solvent B as follows: 4% to 16.5%B in 2.5—52.5 min; 33.5% B in 73 min followed by a stay at 98% B for 3 min and re-equilibration at 2% B. A flowrate of 0.7 µL/min was used. Peptides were sprayed in an electrospray ionization (ESI) source using a stainless-steel emitter with 2 kV at a capillary temperature of 275 °C. A full-scan MS spectrum was collected at 60,000 resolution at m/z of 200 and scanned at 350–1200 m/z with automatic gain control (AGC) of 3E6. The top 10 precursors were selected, and an MS/MS scan was obtained at 45,000 resolution with 50 ms injection time, isolation window of 0.9 m/z with offset 0.1 m/z, normalized collision energy (NCE) of 29. For MS2, minimum AGC target was set to 1.7E4. Dynamic exclusion duration was set to 15 s. The fixed first mass was set to 100 m/z. Charge state exclusion was set to ignore unassigned, 1, and 7 and greater charges. For internal mass calibration, lock mass of 371.10124 m/z was used. Mass spectrometry data was analysed using Proteome Discoverer 2.3 (Thermo Fisher Scientific) software with search engine Mascot (version 2.6.0). Data was searched using latest UniProt Human protein database. Unfragmented precursor and TMT reporter ions were removed using a non-fragment filter in the PD 2.3 workflow. Search parameters included 3 missed cleavages for trypsin, oxidation (M) and deamidation (N, Q) as variable modifications. Tandem label (229.163 Da) at N-terminus and lysine residues and carbamidomethylation on cysteine residues were set as fixed modifications. The mass tolerances on precursor and fragment masses were set at 20 ppm and 0.05 Da, respectively, for MS2 analysis. Consensus step in PD2.3 included several nodes for spectrum, peptide and protein grouping and FDR calculation. Reporter ions for TMT labelled peptides were quantified using the PD quantitation node and peak integration tolerance was set at 20 ppm by considering most confident centroid peaks. Signal to noise values were calculated in addition to measurement of intensities of the TMT reporter ion for peptide and protein quantitation. The intensities were normalized by total peptide amount in PD 2.3. To account for protein input, the global quantitative proteome data was reviewed before normalization and no samples showed an unexpected pattern of distribution. Albumin and haemoglobin abundances were not significantly different between sub-cohorts. Further normalization of the data across all samples was carried out using Reporter Ion Quantitation in Proteome Discoverer, which calculates the total sum of the abundance values for each TMT channel over all peptides identified within a file. The channel with the highest total abundances served as a reference for correcting abundances across the remaining channels by a constant factor. 

Biopsies were only available from 19 COPD and 8 HV-ES for imaging analysis as samples for RNASeq analysis was prioritised (Fig. ). Biopsies for immunofluorescence were fixed in 4% (w/v) formaldehyde and embedded in paraffin wax, as previously prescribed . Ten serial sections of 5 μm thickness were obtained using a Leica RM2135 microtome (Leica Biosystems, Germany). The sections were then mounted onto APES (3-aminopropyltriethoxysilane) coated microscope slides. (CellPath, Powys). To deparaffinise the tissue, microscope slides containing the tissue were placed into slide racks and into clearene (Leica Biosystems) twice with ten minutes each in a clearene tank. The slides were rehydrated by soaking them in graded alcohols solutions of 100% ethanol, 75% ethanol and 70% ethanol, 5 min each. To perform H&E stains, histology slides were placed into Mayer’s Haematoxylin (CellPath, Powys) solution for 5 min, then placed under running tap water for 5 min. The slides were then placed into eosin stain (CellPath, Powys) for further 5 min. Slides were then briefly placed into 100% ethanol, after which they were dehydrated by taking them into 95% ethanol and 100% ethanol (1 min on each). After which they were then placed into clearene three times (3 min each). The slides were mounted in XFT mounting medium (CellPath, Powys) and cover slipped. Antigen retrieval was performed by pipetting 5% (v/v) pronase solution onto the tissue. These slides were left to incubate at room temperature for 10 min, and then rinsed in 1 × PBS. The tissue slides were incubated at room temperature for a further 60 min with blocking solution (1 × PBS + 1% bovine serum albumin & 2% foetal calf serum). Excess blocking buffer was removed, primary antibodies were added and incubated overnight at 4 °C. Primary antibodies were added to all slides except the negative controls. 0.00125 mg/ml mouse anti-EG2 (Diagnostic Development, Uppsala, Sweden) or 0.005 mg/ml rabbit anti-IL-5Rα (ThermoFisher, UK) primary antibodies incubated overnight at 4 °C. Tissue was also incubated overnight with just the blocking buffer. After overnight incubation, the slides were washed in PBS and then 0.002 mg/ml AlexaFluor647 goat anti-mouse or 0.004 mg/ml AlexaFluor647 goat anti-rabbit secondary antibodies (both ThermoFisher) were incubated at room temperature for 1 h, after which the slides were washed in PBS. The tissue was then stained with DAPI (Roche, Germany) at room temperature for 10 min, after which it was washed off. The slides were mounted onto coverslips with mowiol (Sigma-Aldrich, UK). 

Analysis of two groups was performed using a Mann–Whitney U test. Fishers exact test was used for categorical data (GraphPad Prism v9, GraphPad Software Inc., San Diego, USA). Associations were assessed using Spearman's correlation with rho and p values presented. Results were considered significant if p < 0.05. 

Subject demographics This study included 31 COPD subjects and 20 HV-ES as the most relevant control group; characteristics are summarised in Table . No significant differences were seen in age, sex, or BMI. Expected differences were seen in lung function. Eosinophils and associated proteins are increased in COPD We first characterised whether blood eosinophil levels were different between COPD subjects and HV-ES and observed a significantly higher number of eosinophils in the blood of COPD subjects (p = 0.0264). Furthermore, we found a trend towards increased proportions of eosinophils in the BAL of COPD subjects vs HV-ES (p = 0.0547) (Table ). We subsequently conducted an unbiased proteomic analysis of the BAL supernatant from our cohort and, of the 906 proteins detected, 7 proteins were significantly more highly expressed in the BAL of COPD subjects compared to HV-ES (Fig. A). These 7 proteins all associated with granulocyte activation, and included proteases such as MPO and ELANE. The protein with the greatest fold-change (Log 2 FC 1.81) and significance (adjp = 0.002) was the eosinophil-associated Charcot-Leyden Crystal (CLC/Galectin-10). In addition to numbers of eosinophils in blood and BAL, we characterised the presence of tissue eosinophils in formalin-fixed, paraffin-embedded (FFPE) bronchial biopsies from a sub-cohort of 19 COPD subjects and 8 HV-ES (demographics in Table S1). Using immunofluorescence (IF) targeting eosinophil cationic protein (ECP), we identified the presence of tissue eosinophils in both COPD and HV-ES subjects (Fig. B). Furthermore, we found significantly greater numbers of eosinophils in tissue from COPD subjects vs HV-ES (p = 0.0108, Fig. C). Within this IF subcohort, numbers of blood and BAL eosinophil proportions were also significantly elevated in COPD patients compared with HV-ES (Table S1). Furthermore, the signal of significantly greater CLC/Galectin-10 expression in the COPD vs HV-ES BAL proteome was maintained (Log 2 FC 1.12, adjp = 0.014). BAL and tissue eosinophil levels do not correlate with blood eosinophils in COPD Following demonstration of increased numbers of blood and tissue eosinophils and proportions of BAL eosinophils in COPD vs HV-ES, we next investigated whether levels of eosinophils correlated between each compartment. Within the whole cohort, we found only a weak positive correlation between blood eosinophil numbers and the proportion of eosinophils in BAL (rho = 0.2608, p = 0.0337). There was no correlation between blood eosinophil numbers and the proportion of eosinophils in BAL in COPD subjects alone (rho = 0.01762, p = 0.4632). We found no correlation between blood eosinophil numbers and tissue eosinophil numbers in COPD and HV-ES (rho = 0.1192, p = 0.2768) or in COPD subjects alone in the IF subcohort (rho = −0.09192, p = 0.3541). There was a weak negative correlation between tissue eosinophil numbers and the proportion of eosinophils in BAL in the IF subcohort (rho = −0.3433, p = 0.0430). This was also seen in IF subcohort COPD subjects alone (rho = −0.4961, p = 0.0181). Gene expression changes associated with blood eosinophils within COPD To further investigate differences that may be driving mechanisms underlying eosinophilic inflammation specifically in COPD, we compared COPD subjects who had ≥ 300 cells/µl blood eosinophils or not upon recruitment  (Table ); these subjects are referred to as high and low blood eosinophil COPD subjects respectively. High blood eosinophil COPD subjects had better preserved lung function than low blood eosinophil COPD subjects. However, there were no other significant differences between the groups besides blood eosinophil numbers. To understand the processes that might be contributing to increased eosinophils in COPD, we compared genes that are differentially regulated between high and low blood eosinophil COPD subjects in blood, epithelial brushing and bronchial biopsy samples. In blood, 8 differentially expressed genes (DEGs) were identified (Fig. A, all upregulated), including IL5RA, SIGLEC8 and CLC (Table S2A). In epithelial brushings, 21 DEGs (Fig. B, 10 upregulated, 11 downregulated) were identified (Table S2B). Furthermore, there were 137 DEGs in bronchial biopsies (Fig. c, 116 upregulated, 21 downregulated) (Table S2C). There was no overlap of DEGs between the three compartments sampled. (Fig. D). Over-representation analysis revealed the bronchial biopsy DEGs were enriched in metabolic processes (Fig. E). Gene expression changes associated with tissue eosinophils within COPD IF subcohort To understand which genes may determine an increase in tissue eosinophils, we next investigated both the phenotypic and gene expression associations with tissue eosinophils in COPD subjects within the IF subcohort. As there are no definitive thresholds of high tissue eosinophils in the literature we used the median value of 200 cells/mm 2 as a threshold to define the low and high tissue eosinophil groups. Unlike the previous blood eosinophil analysis, there was no significant difference in lung function between the 2 groups of COPD subjects (Table ). In line with the negative correlation described above, the proportion of BAL eosinophils was significantly lower in the tissue eosinophil high group than in the tissue eosinophil low group (Table ). Analysing gene expression differences between COPD subjects with high and low tissue eosinophils from this smaller IF subcohort revealed 5 DEGs in blood (Fig. A, upregulated, 2 downregulated). In epithelial brushings, 32 DEGs (18 upregulated, 14 downregulated) were identified, including POSTN , encoding a matrix protein induced by IL-13 (Fig. B) . In bronchial biospies, 13 DEGs (4 upregulated, 9 downregulated) were identified (Fig. C). Again, no DEGs were commonly differentially expressed between all 3 compartments (Fig. D). However, IFI6 and LRRC37A were downregulated in both epithelial brushings and biopsies (Table S3). Gene set enrichment anlaysis (GSEA) was performed to further investigate the potential impact of IL-13 signalling in epithelial brushings and identified positive significant enrichment of the REACTOME IL-13/IL-4 signaling pathway (NES = 1.54, p = 0.005, Fig. E). Tissue IL-5Rα expression To further investigate why tissue eosinophils may be increased in COPD, we also probed the biopsy tissue in the IF subcohort for the α-subunit of the receptor for the eosinophil survival factor IL-5 (IL5Rα). IL5Rα appeared to be widely expressed in submucosal glandular tissue in addition to eosinophils (Fig. A). In contrast to tissue eosinophils, there was no significant difference in total expression of IL5Rα between COPD subjects and HV-ES subjects in the IF subcohort (p = 0.1633). There was also no significant differential expression of IL5Rα between COPD patients with high blood eosinophils compared to COPD patients with low blood eosinophils in the IF subcohort (p = 0.2294). In the IF subcohort, we observed a weak but significant correlation between tissue IL5Rα expression and tissue eosinophil numbers in COPD and HV-ES (rho = 0.3148, p = 0.0478) and this correlation strengthened when only COPD subjects from the IF subcohort were included (rho = 0.5628, p = 0.0121, Fig. B). There was significant differential expression of IL5Rα between COPD high tissue eosinophils vs COPD low tissue eosinophils (p = 0.031, Fig. C). In terms of IF subcohort patient phenotype, there was a significant increase in tissue IL5Rα expression in those patients who experienced frequent exacerbations (p = 0.0131, Fig. D) but tissue eosinophils were not different based on exacerbation history (p = 0.1388, Fig. E). 

This study included 31 COPD subjects and 20 HV-ES as the most relevant control group; characteristics are summarised in Table . No significant differences were seen in age, sex, or BMI. Expected differences were seen in lung function. 

We first characterised whether blood eosinophil levels were different between COPD subjects and HV-ES and observed a significantly higher number of eosinophils in the blood of COPD subjects (p = 0.0264). Furthermore, we found a trend towards increased proportions of eosinophils in the BAL of COPD subjects vs HV-ES (p = 0.0547) (Table ). We subsequently conducted an unbiased proteomic analysis of the BAL supernatant from our cohort and, of the 906 proteins detected, 7 proteins were significantly more highly expressed in the BAL of COPD subjects compared to HV-ES (Fig. A). These 7 proteins all associated with granulocyte activation, and included proteases such as MPO and ELANE. The protein with the greatest fold-change (Log 2 FC 1.81) and significance (adjp = 0.002) was the eosinophil-associated Charcot-Leyden Crystal (CLC/Galectin-10). In addition to numbers of eosinophils in blood and BAL, we characterised the presence of tissue eosinophils in formalin-fixed, paraffin-embedded (FFPE) bronchial biopsies from a sub-cohort of 19 COPD subjects and 8 HV-ES (demographics in Table S1). Using immunofluorescence (IF) targeting eosinophil cationic protein (ECP), we identified the presence of tissue eosinophils in both COPD and HV-ES subjects (Fig. B). Furthermore, we found significantly greater numbers of eosinophils in tissue from COPD subjects vs HV-ES (p = 0.0108, Fig. C). Within this IF subcohort, numbers of blood and BAL eosinophil proportions were also significantly elevated in COPD patients compared with HV-ES (Table S1). Furthermore, the signal of significantly greater CLC/Galectin-10 expression in the COPD vs HV-ES BAL proteome was maintained (Log 2 FC 1.12, adjp = 0.014). 

Following demonstration of increased numbers of blood and tissue eosinophils and proportions of BAL eosinophils in COPD vs HV-ES, we next investigated whether levels of eosinophils correlated between each compartment. Within the whole cohort, we found only a weak positive correlation between blood eosinophil numbers and the proportion of eosinophils in BAL (rho = 0.2608, p = 0.0337). There was no correlation between blood eosinophil numbers and the proportion of eosinophils in BAL in COPD subjects alone (rho = 0.01762, p = 0.4632). We found no correlation between blood eosinophil numbers and tissue eosinophil numbers in COPD and HV-ES (rho = 0.1192, p = 0.2768) or in COPD subjects alone in the IF subcohort (rho = −0.09192, p = 0.3541). There was a weak negative correlation between tissue eosinophil numbers and the proportion of eosinophils in BAL in the IF subcohort (rho = −0.3433, p = 0.0430). This was also seen in IF subcohort COPD subjects alone (rho = −0.4961, p = 0.0181). 

To further investigate differences that may be driving mechanisms underlying eosinophilic inflammation specifically in COPD, we compared COPD subjects who had ≥ 300 cells/µl blood eosinophils or not upon recruitment  (Table ); these subjects are referred to as high and low blood eosinophil COPD subjects respectively. High blood eosinophil COPD subjects had better preserved lung function than low blood eosinophil COPD subjects. However, there were no other significant differences between the groups besides blood eosinophil numbers. To understand the processes that might be contributing to increased eosinophils in COPD, we compared genes that are differentially regulated between high and low blood eosinophil COPD subjects in blood, epithelial brushing and bronchial biopsy samples. In blood, 8 differentially expressed genes (DEGs) were identified (Fig. A, all upregulated), including IL5RA, SIGLEC8 and CLC (Table S2A). In epithelial brushings, 21 DEGs (Fig. B, 10 upregulated, 11 downregulated) were identified (Table S2B). Furthermore, there were 137 DEGs in bronchial biopsies (Fig. c, 116 upregulated, 21 downregulated) (Table S2C). There was no overlap of DEGs between the three compartments sampled. (Fig. D). Over-representation analysis revealed the bronchial biopsy DEGs were enriched in metabolic processes (Fig. E). 

To understand which genes may determine an increase in tissue eosinophils, we next investigated both the phenotypic and gene expression associations with tissue eosinophils in COPD subjects within the IF subcohort. As there are no definitive thresholds of high tissue eosinophils in the literature we used the median value of 200 cells/mm 2 as a threshold to define the low and high tissue eosinophil groups. Unlike the previous blood eosinophil analysis, there was no significant difference in lung function between the 2 groups of COPD subjects (Table ). In line with the negative correlation described above, the proportion of BAL eosinophils was significantly lower in the tissue eosinophil high group than in the tissue eosinophil low group (Table ). Analysing gene expression differences between COPD subjects with high and low tissue eosinophils from this smaller IF subcohort revealed 5 DEGs in blood (Fig. A, upregulated, 2 downregulated). In epithelial brushings, 32 DEGs (18 upregulated, 14 downregulated) were identified, including POSTN , encoding a matrix protein induced by IL-13 (Fig. B) . In bronchial biospies, 13 DEGs (4 upregulated, 9 downregulated) were identified (Fig. C). Again, no DEGs were commonly differentially expressed between all 3 compartments (Fig. D). However, IFI6 and LRRC37A were downregulated in both epithelial brushings and biopsies (Table S3). Gene set enrichment anlaysis (GSEA) was performed to further investigate the potential impact of IL-13 signalling in epithelial brushings and identified positive significant enrichment of the REACTOME IL-13/IL-4 signaling pathway (NES = 1.54, p = 0.005, Fig. E). 

To further investigate why tissue eosinophils may be increased in COPD, we also probed the biopsy tissue in the IF subcohort for the α-subunit of the receptor for the eosinophil survival factor IL-5 (IL5Rα). IL5Rα appeared to be widely expressed in submucosal glandular tissue in addition to eosinophils (Fig. A). In contrast to tissue eosinophils, there was no significant difference in total expression of IL5Rα between COPD subjects and HV-ES subjects in the IF subcohort (p = 0.1633). There was also no significant differential expression of IL5Rα between COPD patients with high blood eosinophils compared to COPD patients with low blood eosinophils in the IF subcohort (p = 0.2294). In the IF subcohort, we observed a weak but significant correlation between tissue IL5Rα expression and tissue eosinophil numbers in COPD and HV-ES (rho = 0.3148, p = 0.0478) and this correlation strengthened when only COPD subjects from the IF subcohort were included (rho = 0.5628, p = 0.0121, Fig. B). There was significant differential expression of IL5Rα between COPD high tissue eosinophils vs COPD low tissue eosinophils (p = 0.031, Fig. C). In terms of IF subcohort patient phenotype, there was a significant increase in tissue IL5Rα expression in those patients who experienced frequent exacerbations (p = 0.0131, Fig. D) but tissue eosinophils were not different based on exacerbation history (p = 0.1388, Fig. E). 

Our study used a deeply-phenotyped cohort to demonstrate a clear increase in eosinophils in the blood, BAL and tissue of mild-moderate COPD patients compared to healthy, ex-smokers. These eosinophils appear to be activated in the lung due to increased levels of CLC/Galectin-10 in the BAL of COPD patients. We further found that blood eosinophil levels did not correlate with eosinophil levels in BAL or lung tissue in COPD subjects, raising new questions about the utility of this measure alone for understanding lung tissue inflammation. We observed a significant correlation between tissue eosinophils and tissue IL5Rα expression. However, lung IL5Rα expression was not limited to eosinophils, with substantial IL5Rα expression in submucosal glands, further indicating a possible role for IL-5 receptor signalling in the epithelium . Of note, in the context of recent trials, tissue eosinophils were associated with DEGs known to be regulated by IL-13. Our study highlights the complexity of eosinophilic inflammation in COPD with the impact of eosinophils on disease characteristics being subtle and related to the compartment in which eosinophils are measured. Current GOLD guidelines recommend using blood eosinophil counts ≥ 300 cells/µl as a biomarker to identify those with the greatest likelihood of treatment benefit with inhaled corticosteroids . Using RNASeq, bronchial biopsies had the greatest number of DEGs that associated with this measure of increased blood eosinophils with enrichment of genes in numerous metabolic pathways. However, there were very few bronchial biopsy genes associated with tissue eosinophilia. These data highlight the impact that measuring eosinophils in different compartments has on defining patient lung relevant endotypes. Furthermore, they suggest a complex relationship between raised eosinophils and T2-gene signatures, independent of IL-5, which has implications for using eosinophils as a treatable trait in COPD. Our data also highlight the complexity in defining raised eosinophils in COPD as whilst previous data supports a correlation between blood and sputum [ , , ], blood does not appear to be a good biomarker of tissue or BAL eosinophils in patients and we observed marked compartmental differences. Adding further complexity, there was a negative correlation between BAL and tissue eosinophils. Eltoboli et al. (2015) previously demonstrated an association of increased tissue eosinophils with reticular basement membrane thickening of COPD patients, suggesting that the presence of these cells in tissue is important for disease . Whilst there is some evidence for the stability of eosinophil numbers in bronchial biopsies from COPD patients , our observation of a lack of correlation between compartments agrees with a study including 294 COPD patients that also demonstrated no correlation between compartments . The weak correlations between blood and BAL suggests vascular leakage may be playing a role in the detection of eosinophils in the BAL but not residency of these cells in the tissue. Indeed, as there is a negative correlation between the tissue and the airway lumen, this observation might suggest that opposing mechanisms are involved in retaining eosinophils in the tissue or trafficking to the lumen. Despite the lack of correlation between blood and tissue eosinophils, we did observe a gene expression signature associated with blood eosinophils in the bronchial biopsy tissue, many of which were associated with metabolic processes. This increased expression of metabolic genes could be a result of the ongoing energy demands of eosinophil-driven inflammatory processes in the tissue . Alternatively, given that high blood eosinophils were associated with more preserved lung function, this increase in metabolic processes may represent ongoing repair processes . Further work will be required to either confirm or refute these speculations. To understand the drivers of increased tissue eosinophils in COPD, we investigated the gene expression differences associated with this trait. In brushings, there is further evidence of increased T2 inflammatory processes, with increased expression of POSTN , encoding periostin. Periostin expression in the epithelium is increased by the T2-cytokine, IL-13, and in asthma has been associated with both increased airway eosinophils and mucus secretion . In biopsies, there was also an increase in CCL26 , encoding the eosinophil chemokine eotaxin-3, expression of which is also known to be increased by IL-4 and IL-13 . Our GSEA further supports the presence of T2 inflammation through identification of positive enrichment of IL-13/IL-4 signalling in the brushings dataset. This evidence of increased IL-13 signalling in patients with increased tissue eosinophils may provide some explanation for the observed efficacy of dupliumab in COPD . Taken together, these data suggest that elevated eosinophil counts are common in COPD but that tissue eosinophils may be more relevant to continuing disease processes than blood or luminal eosinophils. This observation may explain the limited efficacy of anti-IL-5 treatments in preventing COPD exacerbations and disease progression as these studies used blood eosinophils to stratify these patients . These anti-IL-5 treatments are effective at reducing blood eosinophils, but our data suggest that IL-13-driven pathways may be responsible for the maintenance and survival of these cells in COPD lung tissue. Further support for this observation comes from the BOREAS trail which demonstrated no significant effect of 52-week dupliumab treatment on blood eosinophil levels despite a significant effect on exacerbations and prebronchodilator FEV1, although FeNO was reduced . We recognise that this study is not without its limitations. Due to the deep characterisation of subjects and intensive sampling, the cohort is small and the study captured only cross-sectional measures, providing no insight into eosinophil stability in the different compartments over time and disease states (e.g. stable vs exacerbation). We thus cannot fully rule out that lack of correlations are not due to cohort size. In particular, the size of the IF cohort was limited due to capture of tissue for RNASeq being prioritised. Additionally, sputum data was not available from all subjects and thus we have no data as to the correlation between tissue and sputum eosinophils. Our study included mild-moderate COPD patients and gives insights about earlier disease. Comparison of our findings with those in more severe disease in future studies is merited. Our study demonstrates a clear increase in eosinophils in COPD compared to health in both blood, BAL and lung tissue. Furthermore, we demonstrate that blood eosinophil levels did not correlate with eosinophils nor IL5Rα expression in tissue. We have identified differentially expressed genes that associated with eosinophils in different compartments. Blood eosinophils do define an inflammatory endotype that can be detected in lung tissue but do not reflect the expected IL-5-mediated pathways. Further delineating the complex signalling pathways driving tissue eosinophilic inflammation in COPD with mechanistic studies could provide information on optimal targeting of existing and novel therapies. 

There is active eosinophilic inflammation in the lungs of COPD patients compared to controls but blood eosinophils alone do not reflect tissue eosinophils or gene expression. Understanding of lung eosinophil biology is needed to tailor new therapy. 

Additional file 1 Additional file 2 Additional file 3 

