Training-Associated Alterations In Equine Respiratory Immunity – A Multi-Omics Comparative Approach.

Background: Neutrophilic airway inammation is highly prevalent in racehorses in training, with the term mild to moderate equine asthma (MMEA) being applied to the majority of such cases. Our proposed study is largely derived from the strong association between MMEA in racehorses and their entry into a race training program; this has led to our primary aim of measuring the effect of race training on pulmonary immune cell function. The objectives of this study are to characterise the effect of training on the local pulmonary immune system and quantify the magnitude of effect by dening (a) the gene expression of tracheal wash (TW) derived cells and (b) the protein expression of TW samples derived from Thoroughbred racehorses prior to (T0) and following commencement of race training (T1). Results: Tracheal wash samples were obtained at both time points (T0, T1) from the same Thoroughbred horses (n=16). Gene expression of TW derived cells, determined by RNAseq and analysed by DEseq2 detected 2,138 differentially expressed (DE) genes during the training period; 1,122 of these were upregulated. The proteome of TW samples was also evaluated and contained 260 DE proteins during the training period; 103 of these were upregulated and the rest downregulated. Gene and protein sets were enriched for biological processes related to acute phase response and oxidative stress, as well as to immune response and inammation. Many genes, including ISG15, ISG20, IFI35, SOCS1 and TRIM21, were highly enriched for IFN signaling. Interestingly, pathway analysis also highlighted genes and proteins related to haemopoietic processes. Conclusions: This study demonstrated TW samples to represent a rich source of airway cells, protein and RNA to study airway immunity in the horse and highlighted the benets of a multi-omics methodological approach to studying the dynamics of equine airway immunity. Intense training induced quantiable alterations in both gene and protein expression of airway derived cells is consistent with deregulation of airway immunity and haemopoietic abnormalities. Respectively, these ndings likely reect the known associations between race training and both airway inammation and bleeding, in particular offering further insight into the potential mechanisms which underpin training associated airway inammation.


Introduction
During the early phase of training, racehorses commonly develop airway in ammation. Prevalence rates as high as 70-80% have been reported (1)(2)(3)(4), and generally reduce as training progresses (5,6). This airway in ammatory response, when associated with poor athletic performance and/or chronic coughing has been termed Mild to Moderate Equine Asthma (MMEA) (4,7). Although the precise aetiopathogenesis of MMEA remains unde ned, proposed aetiological candidates have included both infectious (bacteria and viruses) and non-infectious (organic dust) agents; indeed, some cases may be attributable to a combination of infectious and non-infectious factors acting in concert (4). In relation to infectious factors, the reported association with the early training period may re ect an initial period of susceptibility with the subsequent temporal establishment of immunity. Although entry into race training may result in increased exposures to airborne infectious and non-infectious agents, it is also associated with a signi cant increase in exercise intensity and frequency. The well recognized association between high intensity exercise and symptoms of respiratory infection amongst human athletes has fueled interest in the impact of training on immune function (8,9). Despite localization of in ammation to the airways in MMEA (7), few studies have assessed the in uence of training on immune cell function speci cally at this anatomical site, a key consideration in light of the recently reported disassociation between the response of blood derived equine monocytes and alveolar macrophages (AMs) to training (10,11). Both an exercise-associated reduction in equine AM phagocytic capacity (12) and a training-associated derangement in the responsiveness of equine AMs to various TLR ligands have been reported, theoretically re ecting increased susceptibility to opportunistic infection (10). In agreement, we recently demonstrated a training-associated alteration in equine AM basal gene expression, and con rmed a degree of immunomodulation at the level of the airways (11). Consequently, there appears to be su cient justi cation for further de ning the link between training and airway immunity.
Exercise-induced pulmonary haemorrhage (EIPH) is also highly prevalent in racehorses, with almost 100% of racing Standardbreds having direct (erythrocytes) or indirect (haemosiderophages) cytological evidence of recent bleeding within the airways (13). Despite the high prevalence and co-existence of both MMEA and EIPH within this population, there is little evidence for any aetiological or pathogenetic link.
This study was designed to further expand upon and adapt our previously published observations, through targeting TW-derived cells from a larger number of horses and adopting a more rigorous study design and extensive analytical approach (RNAseq and proteomics). The adopted analytical approach was selected primarily to further de ne the mechanisms underpinning training-associated immunomodulation in racehorses, with the added potential of allowing further characterization of training-associated airway diseases and identify disease-and/or disease susceptibility-associated biomarkers and/or novel therapeutic targets.
Thus, the principle objective of the study was to investigate the impact of entry into race training on airway immune function. The "open window" theory, re ecting a temporal association between intense training and increased susceptibility to opportunistic infection, is well recognized within the eld of human exercise immunology (8, 14) and may explain the association between increased tracheal mucus (a marker of MMEA) and both time in training and the isolation of Streptococcus zooepidemicus and nonhaemolytic streptococci from respiratory tract secretions of racehorses in training (6). Such supportive information has the capacity to inform management practices which, in light of the inherent necessity for episodes of high speed exercise within a training regimen, are likely to focus on either opportunistic pathogen avoidance schemes and/or novel immunomodulation strategies. Evidently, the latter approach would be highly dependent on a more thorough understanding of the mechanisms and pathways underpinning MMEA with a view to identifying potential therapeutic and/or prophylactic targets; the adopted transcriptomic and proteomic analyses were considered as appropriate methodological approaches to further our understanding of this highly prevalent disease.
We hypothesised that intense training has a signi cant impact on the respiratory immune response of racehorses and that a combined transcriptomics and proteomics approach would further elucidate the speci c immunological downstream molecular cascades impacted by this activity. Moreover, we hypothesised that the data generated may further inform future efforts aimed at establishing the potential role of infection in the aetiopathogensis of MMEA. In order to address the above hypotheses, we studied the gene expression pro le of tracheal wash (TW)-derived cells and the proteomic pro le of whole TW samples derived from UK national hunt Thoroughbreds, both before and after entry into a training schedule. Furthermore, as the racing Thoroughbred offers an ideal model for the study of the impact of exercise it was considered that the results of this study may have signi cant translational application to the increasing body of evidence relating to human exercise immunology.

Materials And Methods
Animals used in the study Tracheal secretions were collected from 16 male Thoroughbred horses [6.9 ± 0.5 years (mean ± SEM)]. All horses were housed at the same racing yard and were trained by the same trainer. Samples were obtained as part of a routine assessment of respiratory health (differential cytology); however, residual sample was retained for an ongoing proteomic and transcriptomic study approved by the Veterinary Ethical Review Committee of the School of Veterinary Medicine, University of Edinburgh and informed consent was obtained for the use of residual sample material. Enrolled horses underwent physical examination, endoscopy of the respiratory tract, and TW. All horses were clinically healthy and treated according to standard welfare procedures. All methods reported are in accordance with ARRIVE guidelines (https://arriveguidelines.org) (15).

Sample collection and cell isolation
Field work -Tracheal wash sample collection Samples were collected at two different time points: T0, when the horses were at rest between racing seasons and T1, when the horses were in active race training. T0 sampled were collected following 59 ± 5.7 days at rest and T1 samples were collected when the horses were considered "racing t", having entered into the training program for approximately 2.5 months (76 ± 11.4 days) and having already ran in at least 2 competitive races (average number of starts 3.9 ± 0.4) that season. All TW samples were collected using a trans-endoscopic technique, as previously described (16). Horses were manually restrained and a nose twitch applied when necessary. A 1500mm working length, 9.2 mm outer diameter video endoscope (2.8mm biopsy channel) (Aohua, China), was passed via the ventral nasal meatus, via the pharynx and into the tracheal lumen via the rima glottidis. Following assessment of the amount and nature of mucus deposits within the tracheal lumen, approximately 20 mL of sterile 0.9% saline at room temperature was instilled via a catheter passed via the biopsy channel of the endoscope at the proximal aspect of the trachea. The endoscope was then further advanced to the level of the thoracic inlet where the pool of instilled uid had gravitated. As much uid as possible was subsequently aspirated via the transendoscopic catheter. Samples were stored on ice and submitted for laboratory analysis at the Roslin Institute and the Royal (Dick) School of Veterinary Studies and processed within 4 hours of collection.

Differential Cell Count of Tracheal Wash
An aliquot of 0.5 ml was submitted to the pathology lab at the Royal (Dick) School of Veterinary Studies for differential cell count analysis (DCC). Differential leucocyte count (minimum of 200 cells) was performed on dithiothreitol (DTT) treated samples and expressed as a percentage of total non-squamous and non-epithelial nucleated cells, as previously described (17). Horses were considered free of in ammatory airway disorders based on the following DCC cut off values: neutrophils < 20%, eosinophils (< 1%) (16, 18). Samples containing haemosiderophages were considered consistent with a diagnosis of EIPH (13).

Total cell Isolation from Tracheal Wash
A 1ml aliquot of sample was stored at -80°C for future proteomic analyses and the remaining sample incubated for 15 minutes at room temperature in 0.1% DTT to depolymerize secreted mucin, as previously described (17). Dithiothreitol has been demonstrated to cause no deleterious effects on human sputum derived cells or interfere with surface marker measurements using ow cytometry (19,20). Following gravity ltration through a 100µm pore mesh lter, the sample was centrifuged at 400g for 10min at 4 o C.
Supernatant was carefully removed and the cell pellet re-suspended in Dubelcco's PBS, from which a total cell count (excluding squamous epithelial cells) and cell viability (Trypan Blue exclusion staining) was performed using a haemocytometer. Afterwards, the cell pellet was resuspended in 1ml of Trizol and stored at -80oC for future RNA analysis.

RNA Analysis
Total RNA extraction Total RNA was extracted using a combination of Trizol reagent (Thermo Scienti c™, 15596026) and an RNAeasy plus micro kit (Qiagen, cat no 74034), according to manufacturer's instructions. Following removal of the supernatant, the remaining cell pellet was lysed by the addition of 1mL Trizol Reagent. Subsequently, 0.2 mL 1-Bromo-3-chloropropane (BCP) (Sigma Aldrich, cat no B9673-200ML) was added per 1 mL of Trizol. The sample was shaken vigorously for approximately 30 sec and left at RT for 5 min to allow complete dissociation of nuclear-protein complexes. The homogenate was then centrifuged at 18,000g for 15 min at 4 0 C resulting in the formation of a lower red phenol-chloroform phase, an interphase, and an upper colourless aqueous phase. The aqueous phase contained the RNA and had almost 50% of the volume of the Trizol used, plus the volume of the sample. Following transfer to a clean tube for the precipitation step, 0.5 mL of 70% ethanol was added, the sample stored for 2 hours at -20 0 C and then transferred to an RNeasy spin column and centrifuged at 18,000g for 5 min at 4 0 C. Following centrifugation, the ow through was removed, the RNA washed once with RW1 buffer and DNA treatment was performed using the RNase-Free DNase Set (Qiagen, cat no 79254) according to the manufacturer's instruction. Due to DNA contamination, this step was performed twice or samples were run through gDNA Eliminator Spin Columns twice after the elution step. Afterwards, the RNA membrane was washed with RW1, RPE and 80% ethanol. Finally, RNA was eluted in 15 µl RNase-free water and RNA samples were stored at -80 o C until further use.
RNA quality assessment RNA concentration and purity was measured using ND-1000 Nanodrop spectrophotometer (Thermo Scienti c, Wilmington, USA) by measuring absorbance at 260 and 280 nm (A260, and A280 respectively). Purity of RNA was determined using the A260/A280 ratio; a ratio approximating 2 was considered to be indicative of pure RNA. RNA integrity was con rmed with the High Sensitivity RNA ScreenTape system (Agilent Technologies). A RNA integrity number (RIN) greater than 7 was considered appropriate for qPCR and RNAseq analysis.

RNA sequencing analysis
Total RNA was processed using NEBNext Ultra Directional RNA Library Prep Kit by Novogene (Novogene Europe) to generate cDNA libraries according to the manufacturer's instructions and subsequently sequenced with an Illumina NovaSeq 6000 platform (Illumina, San Diego, USA) at a depth of 35M reads strand-speci c 150-bp paired-end per sample. Ribosomal RNA (rRNA) was depleted from samples for total RNAseq, using biotinylated, target-speci c oligonucleotides with Ribo-Zero rRNA removal beads.
Processing of next-generation sequencing data and differential expression analysis Raw data is deposited in the Gene Expression Omnibus under the study accession number GSE181221 (https://www.ncbi.nlm.nih.gov/geo/). For each sample, a set of expression estimates, as transcripts per million (TPM), were obtained using the high speed transcript quanti cation tool Kallisto (21) as previously described (22,23). Raw sequencing reads were processed to remove adaptors and poor-quality sequences (Q25 and below) using the Cutadapt program (24). Quality control of trimmed data was assessed using FASTQC (https://www.bioinformatics.babraham.ac.uk/projects/fastqc/) and Multiqc (25). Kallisto (0.46.0.1) was ran on all samples, using as its index the complete set of 59,087 predicted transcripts for the EquCab3.0 genome (ftp://ftp.ensembl.org/pub/release-99/fasta/equus_caballus/cdna/Equus_caballus.EquCab3.0.cdna.all.fa.gz) and the quanti cation of gene expression was calculated. The horse annotation (GTF) was obtained from Ensembl (EquCab3.0) software. For the differential analysis of count data, we used the DESeq2 method which integrates methodological advances with novel features to enable a more quantitative analysis of comparative RNAseq data using shrinkage estimators for dispersion and fold change (26). Differentially expressed genes were then identi ed by applying an FDR cut off of 0.05. Analysis was performed in R v 3.5.0 (https://www.r-project.org/ ).

Quantitative PCR (qPCR)
One microgram total RNA was converted to complementary DNA (cDNA) using the prescription NanoScript reverse transcription kit (SuperScript III First-Strand Synthesis System, Invitrogen, Cat No 18080051), according to the manufacturer's instruction. cDNA was stored at − 20°C until further use.
Transcript levels were quanti ed in triplicate using an MX3005P quantitative polymerase chain reaction (qPCR) system (Stratagene) with primers listed in Table 1 and qPCRBIO SyGreen Mix Lo-ROX kit (PCRBIO). Primer e ciency was validated with a standard curve of ve serial dilution points and SDHA was used as a reference gene. SDHA was used as a reference gene because it remained stable in the RNAseq data and it has previously been determined to be the most stable housekeeping gene to study equine exercise-induced stress (27,28). Reverse transcriptase negative and "no template" control samples were included in each run as negative controls. Data were analysed using Stratagene MxPro software and relative gene expression was calculated using the 2^−ΔΔCT method (29).

Sample preparation for LC-MS
Samples of animals collected at T1 and T0 were pooled and analysed in the Mass Spectrometry Facility at the Roslin Institute (https://www.ed.ac.uk/roslin/facilities-resources/proteomics-and-metabolomicsfacility). Note samples from Horses 3, 7 and 12 were not included in the proteomic analysis for time point T0 and T1, because insu cient protein material was isolated. The protein samples (T0 and T1) were then reduced with dithiothreitol and alkylated with iodoacetamide prior to tryptic digestion on S-TRAP (Proti , USA) cartridges following the manufacturer's protocol. The resulting peptides were cleaned up and labelled by reductive di methylation on Solid Phase Extraction cartridges.
Stable isotopic reductive dimethylation of the resulting peptides was performed on OASIS HLB columns (Waters) following standard methods (31). The peptides were acidi ed by addition of tri uoroacetic acid to 0.1 % (v/v) before loading onto conditioned OASIS HLB columns. Peptides from T0 and T1 samples were labelled with 'light' and 'heavy' labelling reagents, respectively. The 'light' buffer consisted of 0.8 %

LC-MS analysis
The peptides eluted by reversed-phase chromatography were analysed by a micrOTOF-Q II mass spectrometer (Bruker) following earlier method with minor modi cations (33). Spectral data were processed by PoteinScape server 3.1 (Bruker) software and the resulting peak lists were searched using Mascot 2.4 server (Matrix Science) against Uniprot Horse (Equus Caballus) sequence database containing 44,472 entries (https://www.uniprot.org/proteomes/UP000002281). Carbamidomethylation of cysteine was used as a xed modi cation and oxidation of methionine and light and heavy dimethylation of N-terminus and lysine were chosen as variable modi cations. False discovery rate was limited to < 1% for peptide IDs after searching decoy databases. Dimethyl quanti cation was performed by using WARPLC plugin (Bruker) on Proteinscape 3.1 server to integrate extracted ion chromatogram of every precursor. Peptide ratios were normalized based on setting overall peptide median ratio at 1, which corrects for unequal protein sampling and a coe cient of variability of peptide ratios were also determined for each quanti ed proteins.

Gene ontology and pathway analysis
Identi cation of enriched KEGG pathways in the upregulated and downregulated gene lists was performed with DAVID (v6.8) (34). Pathway analysis was performed using Ingenuity Pathway Analysis (35) to infer the functional roles and relationships of the differentially expressed genes based on the log2 fold change value of each gene.

Statistical analysis
Descriptive statistics were performed using R v. 3.5.0 (https://www.r-project.org/ ). After testing for normality, a one-sample t-test was employed to identify differences in cells between time point T0 (rest) and time point T1 (training). Statistical signi cance was assumed at p < 0.05. Numeric results are presented as mean ± SEM.

Cell recovery and populations
The present study used a series of 32 horse-derived TW samples from 16 Thoroughbred racehorses at two time points (T0-resting period, T1-training period) and an average of 10.8x10 6 (± 1.9) cells per sample were isolated, consistent with previously reported data (17,36). None of the animals had clinical signs of respiratory disease or a history of poor performance. Results of the differential cell counts (DCC) are shown in Fig. 1. Although a slightly greater neutrophil and lymphocyte percentage was detected in samples obtained from horses during the training period, there was no signi cant difference between groups for all immune cells in the DCC. As expected, a remarkable increase of haemosiderophages was reported during training. Since all samples were collected from racehorses, the presence of haemosiderophages in these samples was predictable, as almost all racing Thoroughbred horses will bleed to varying degrees into the airways (13).
Intense training-induced gene expression in equine tracheal wash derived cells Whole-transcriptome pro ling (RNAseq) was performed on TW derived cells from 16 Thoroughbred racehorses prior to (T0) and approximately 2.5 months following entry into an intense training program (T1). A total of 21,357 equine genes were identi ed, and the reads were quanti ed to identify those that were differentially expressed between the two timepoints ( Fig. 2). Intense training produced a radical change in equine gene expression, with 2,138 differentially expressed genes (1,122 upregulated and 1,016 downregulated) being identi ed between the two time points. We de ned differentially expressed genes to be those showing up-or downregulation following intense training with a false discovery rate (FDR) below 0.05. The complete lists of mapped genes and differentially expressed genes are presented in the supplemental material (Additional le 1). Figure 3 summarizes the differentially expressed genes using a false discovery rate (FDR) of < 0.05 and a fold change of 1.5. Selected genes which were evaluated by qPCR con rmed their differential expression (Fig. 4).
Intense training-induced protein expression in equine tracheal wash samples In order to expand on the RNA studies, we also de ned the airway protein pro les (total proteome) of TW samples, thus revealing the mechanisms which may underpin any alterations in immune function.
Protein extraction was successfully performed on TW samples from 13 out of 16 horses. Insu cient protein was isolated from TW samples derived from Horses 7 and 12 at T0 and Horse 3 at T1; thus, these animals were excluded from the proteomic analysis at both time points. An average of 1.6 ± 0.2 (+ SEM) mg protein was isolated from 500µl of TW per animal. Twenty micrograms per sample were pooled for each time point and samples were processed as described above. In order to visualise the total protein load of all the samples, samples were run on gradient gels and stained with instant blue protein stain, as previously described (17). Subsequently, protein samples were submitted for proteomic analysis. We identi ed distinctive protein signatures in the Thoroughbred TW samples in response to training. In total,

Functional analysis of RNAseq and proteomic datasets
In line with previous observations (37), the correlation between individual molecular immune transcripts and proteins was low. Interestingly, although analysis of both the RNAseq and proteomic datasets revealed a limited number of gene/proteins with a similar differential expression in response to training (Table 3), a clear correlation between the two complete datasets was apparent at both the pathway and cascade level. KEGG pathway and gene ontology enrichment analysis of the upregulated and downregulated genes (Additional le 3) and proteins (Additional le 4) was performed using the Database for Annotation, Visualization, and Integrated Discovery (DAVID) annotation software [https://david.ncifcrf.gov/ (34)]. These analyses revealed a training associated activation (at the level of both gene and protein expression) of pathways and biological processes involved in acute in ammatory response, immune defence and immune cell activation and differentiation. Additional pathways which were activated in response to training included metabolic pathways, cellular and oxidative stressassociated pathways and pathways related to heme degradation or haematopoietic abnormalities; these latter pathways likely re ected the higher prevalence of EIPH in horses in training. Importantly, in the RNAseq dataset, we detected pathways related to viral infection and interferon signalling. Ingenuity Pathway Analysis (IPA) was also used to functionally assess the list of differentially expressed genes and proteins (35). Canonical pathway network analysis based on the 2,138 differentially expressed genes indicated activation of IFN signalling (Fig. 5). In line with DAVID analysis, cellular stress-related functions, such as oxidative phosphorylation and mitochondrial dysfunction, as well as those related to the immune response and in ammation were identi ed as being highly enriched. Interestingly, despite the small overlap between the transcriptomic and proteomic datasets the majority of the Diseases and Biofunctions determined to be enriched by IPA were common between the two datasets; 76.5% of the upregulated and 56.6% of the downregulated molecules (Fig. 6). Full database content for each IPA category and the genes related to this are presented in Additional les 5, 6.
Training-associated immune response on equine respiratory immunity In line with previous literature in different species, including humans (8) and horses (10,11), the current study also supports the principle that the immune system is very responsive to exercise. Interferons have diverse functional properties and are involved in numerous biological processes related to host response to infection, in ammation, cancer, autoimmunity and metabolic disorders. Based on the publicly available interferome database (http://www.interferome.org/interferome/home.jspx), we detected the differential expression of several interferon regulated genes and proteins during intense training (38). Figure 7 shows the Venn diagrams created based on the number of genes or proteins regulated by one or more IFNs (Type I, II or III).
Many of the recently described panel of 90 highly conserved mammalian type 1 interferon-stimulated genes (ISG), were differentially expressed in our data (39); these genes and a few encoding proteins, most of which were upregulated are presented in Table 4. They are broadly involved in antigen presentation, antiviral response, IFN activation and repression, apoptosis and cell signalling.
Interestingly, the classical sensor for RNA PAMPs, DHX58 was upregulated, whereas, the DNA sensor AIM2 was downregulated. The fragile X mental retardation (FMR1) gene that encodes an RNA-binding protein that plays a role in intracellular RNA transport and in the regulation of translation of target mRNAs, was downregulated in response to training. FMR1 has been recently linked to the IFN response and has been shown to be a proviral factor for in uenza virus and was found to induce mild restriction of HIV-1 (40,41). Other host proviral genes (TOR1B, IFI44, DNAJA1, AREG) and LCN2 protein involved in in uenza A virus replication in mice were also upregulated during training in TW derived cells (41).
Vital IFN-induced antiviral factors such as MX1 and RSAD2, as well as IRF9 (a key transcription factor involved in IFN induction and response) were also upregulated in relation to intense training. Moreover, a link between local synthesis of early components of the complement system and the type I IFN response is now well documented (39). Complement components such as C1r and C2 proteins were also upregulated in TW samples during the training period. In line with this, protein SERPING1 a negative regulator of C1r was downregulated in the same samples. On the other hand, other members of ISG that play a role in the suppression of the IFN system, such as IFI35, TRAFD1, TRIM21 and SOCS1 were also upregulated. The negative regulation of the IFN response is multifaceted and necessary in order to avoid excessive/perpetual activation of IFN-induced pathways (39).
Other pathways involved in immune processes, such as immunological response, respiratory disease and in ammation were downregulated (Additional le 3, 4, 6). Based on the interferome dataset the majority of the interferon regulated genes seem to be downregulated in both datasets (1,155 in RNAseq and 180 in proteomics) during training (Fig. 7). These results revealed key signalling networks involved in Natural killer cell mediated cytotoxicity, T and B cell receptor signalling pathway and TNF signalling. For example, several components such as CD69, NFAT5, NFKBIZ, NKAP, TAGAP, involved in T cell activation and differentiations, as well as the antimicrobial factor LYZ and the pro-in ammatory cytokine TNF were signi cantly downregulated. Besides, the anti-in ammatory elements TIGIT, a T cell suppressor, as well as SOCS1, a negative regulator of cytokine signalling were upregulated. Thus, there seems to be some evidence of training associated immunosuppression on equine TW derived cells. clinical signs, the forms of both conditions were considered to be subclinical. In an attempt to further dissect the effect of the neutrophil in ltration on the gene expression of these animals, whole genome analysis was performed on TW samples derived from Horses with high neutrophil numbers (High N, n = 3) versus healthy individuals / controls (n = 13). Based on RNA sequencing data we identi ed 57 equine genes ( Fig. 8-A) differentially expressed in the High N group compared to the controls (Additional le 7) (following multiple-testing correction, false discovery rate (FDR) p-value < 0.05).
To gain a global view of the results, functional analysis was performed as previously described using IPA and DAVID in both datasets. As expected, in samples derived from animals affected with high neutrophil counts various relevant biological processes related to neutrophil in ltration were observed. Examples are neutrophil and granulocyte migration, neutrophil and granulocyte chemotaxis and neutrophil extravasation. Other immune related terms were also detected; response to biotic stimulus, defence response, response to external stimulus, immune system process, in ammatory disease/response, response to stress and regulation of response to stimulus (Fig. 8-B, Additional le 7).

Discussion
To our knowledge, this is the rst study to report on the "effect of training" on UK Thoroughbred national hunt racehorses, determined by both transcriptomic and proteomic analysis of TW samples harvested before and during a period of training, thus permitting an assessment on the basal gene expression of principal innate immune cells of the equine airway. TW samples are routinely collected in contrast to the less readily available and more technically invasive bronchoalveolar lavage (BAL)-derived samples, which have previously been favoured for transcriptomic and proteomic analyses (42)(43)(44). The transcriptome was derived from high throughput sequencing, and results were analysed with a paired design to account for individual variation. Animals were of similar age, same sex and breed and lived in a similar environment. All animals recruited in this study were from the same yard and subject to the same general management structure. RNA-seq data were analysed using DEseq2 software, a widely used tool in highthroughput RNA sequencing studies using raw counts (26). Compared with a single-end strategy, the paired-end strategy used in this study has been recommended as it minimises the likelihood of false positives (45); moreover, the validity of the data, in terms of differential expression, was supported by the RT-qPCR results. However, irrespective of the methodology used, considering the genetic variability among individuals and the limited sample size, the potential for false-positive results remains, necessitating caution when interpreting individual gene expression data in RNAseq-based studies. Similar caution should be exercised when interpreting the differential expression of individual proteins in the proteomic dataset. Thus, the consideration of linked, rather than individual, genes and proteins conveys more con dence when interpreting these kind of data (45).
Our experimental protocol allowed us to investigate the effect of training on airway innate immune responses at both the gene and protein level. Although this revealed signi cantly limited overlap between the RNAseq and proteomic datasets, this was not unexpected as such a lack of concordance between gene and protein expression has been well documented (46, 47) and can be attributed to a variety of factors such as post-transcriptional machinery, variable half-lives, molecular degradation, or even sampling bias (37). Additionally, certain factors speci c to the current study design may also have contributed to this discordance in gene and protein expression. Firstly, unlike the transcriptomic analysis, proteomic analysis was not restricted to the cellular component of the TW sample, thus permitting the detection of secreted proteins. This approach was considered more likely to provide a more holistic assessment of the airway immune status as well as offering a greater potential to identify biomarkers with potential clinical and/or training applications. Secondly, gene expression data was obtained from each individual animal; whereas, protein analysis was conducted on pooled samples at each time point.
This approach was justi ed on nancial grounds due to the high cost associated with the application of these technologies. It should also be noted that, due to insu cient protein yield, samples from Horses 3, 7 and 12 were not included in the proteomic analysis for both time points. Thirdly, the lack of overlap in gene/protein expression in the current comparison can also be explained by the high ratio (approx. 27) between the expressed genes (21,357) and proteins (802) detected using the applied methodologies.
Lastly, with such transcriptomic studies, the conclusions that can be derived from the data generated are highly dependent on the quality of the annotations of the genome used, which for the horse remains challenging (48). Overall, consideration of these study-related factors and the recognised discordance between gene and protein expression does signi cantly question the validity of inferring changes in protein expression from transcriptomic data alone. Rather, this study highlights the value of applying a combined transcriptomic and proteomic approach when studying cellular mechanisms in both health and disease. Despite the lack of agreement at the individual gene and protein level, signi cant similarities were observed at a pathway and cascade level.
We demonstrated a clear modi cation of both mRNA and protein expression in TW derived samples during the training period. Based on the multi-omics methodologies applied, we identi ed 2,138 differentially expressed equine genes and 260 proteins during this period. Our ndings are consistent with those derived from comparable studies in humans and horses (44,49,50), and support a clear association between intense training and immune system deregulations, haemopoietic and metabolic abnormalities and cellular stress at the level of the airway. It is conceivable that such training-associated changes may play a role in increasing susceptibility to opportunistic infection and airway in ammation in racehorses.

Local immune responses in tracheal airway immune cells during training
Using a multi -omics approach (RNAseq and proteomics) we observed an alteration in a number of gene and protein pathways in association with training (time point T1); many of these pathways were associated with important aspects of airway immunity, such as antigen presentation, complement activation, immune cell chemotaxis, cellular stress, response to bacterial/viral infection and regulation of IFN signalling (Figs. 5 and 7, Additional le 3, 4). Analysis also revealed the increased expression of several myeloid cell chemoattractants (e.g., CCL2, CSF3R, S100A9 and other related S100 genes) at time point T1. Various relevant gene ontology terms were observed, such as neutrophil and monocyte chemotaxis, highlighting a level of immune cell activation in the airways. Additionally, acute phase response, cellular stress and oxidative phosphorylation also characterise the response to intense training.
Overall, both the transcriptomic and proteomic data were consistent with a level of in ammation in the airways during training, a phenomenon very well documented in racehorses. Based on the TW differential cytology results, 3 out of 16 animals developed airway neutrophilic in ammation (as de ned by neutrophil ratio > 20) at time point T1, consistent with a low level of airway in ammation during training, despite the absence of clinical signs or evidence of poor performance. In light of the well-recognised association between intense exercise and an acute phase response in both humans and horses (51), it is quite conceivable that repetitive periods of arduous exercise over a prolonged period of time may also result in low level airway in ammation. Indeed, recent equine studies have shown that an intense training programme will induce an "in ammation-like state" based on the measurement of acute phase protein SAA, which is greater than that induced by a lighter physical activity program (52). Serum amyloid A was also signi cantly induced during the training period (T1) in our study. Despite the relatively low number of horses with a signi cant elevation in TW neutrophil ratio at T1, we did detect a training-associated upregulation of myeloperoxidase (MPO) protein expression, potentially re ective of a low level of airway neutrophil presence/activity below the threshold of detection by TW cytological examination. Indeed, MPO has previously been proposed as a more sensitive marker of airway neutrophil presence compared to cytology (53). Despite its bene cial properties in host defence, overproduction of MPO may exert detrimental effects during in ammation (54), including tissue damage and vascular dysfunction through the generation of potent ROS and activation of MMPs (54,55). Consistent with these downstream effects, we also detected the upregulation of several MMPs (including MMP9) and ROS in association with training (56-58).
In addition to in ammation, training was associated with evidence of tissue remodelling. Chitinase 3 like 1, ARG2 and MMP9, all upregulated during training, have previously been linked to airway remodelling and/or declining lung function (59). MMP9 has been shown to contribute to airway remodelling in severe equine asthma (49). Training was also associated with an upregulation of ORMDL3 and HLA, both of which are highly associated with human and murine asthma (59), with ORMDL3 playing an important role in airway remodelling (increased airway smooth muscle, sub-epithelial brosis and mucus) (60). Interleukin 4 induced 1, another gene which was highly expressed during training, has been shown to inhibit T cell proliferation, regulate the programming of macrophages towards a polarised M2 phenotype and downregulate LPS-induced TNF expression in murine macrophages (61). Indeed, the data suggests that training may induce an overall immunomodulatory milieu at the level of the airway, characterised by increased expression of IL4I1 (61) and ROS (62). This may explain the signi cant downregulation of TNF during the training period and the identi cation of gene ontology terms related with apoptosis and cell death, as ROS have been shown to regulate both necrotic and apoptotic T activated cell death via caspase activation (62).

Training-associated IFN response on equine airways
A large number of IFN regulated and stimulated genes were differentially expressed during the training period. The IFN response is considered a key driver of in ammation in the lung by stimulating the recruitment and activation of immune cells, often in response to viral infection. Despite these bene cial effects, the negative regulation of the IFN response is necessary to avoid excessive and continuous activation of IFN-induced pathways with associated collateral tissue damage (39). Moreover, there is recently increasing evidence that type I IFNs are not only induced by viral infections, but also by bacterial and fungal infections. Indeed, in the case of bacterial lung infections, the type I IFN mediated signalling may be detrimental (63).
Training resulted in a greater differential expression of genes related to the type 1 IFN response, compared to type II and III IFN responses. Many of the recently described and highly conserved mammalian type I IFN ISGs were differentially expressed during training; these were related with antiviral activity, antigen presentation, PAMP detection and apoptosis. Interestingly, we identi ed several of the ISGs upregulated during training to be associated with ubiquitination; these included the ring nger proteins RNF19A, RNF25, RNF181, RNF183 and highlighted protein reformation as part of the IFN response. This was also re ected in the pathway analysis, as protein synthesis was one of the top Biofunctions detected in our datasets. Recent ndings derived from different mammalian species (including the horse), showed the type I IFN response to positively bias the sensitivity of RNA virus surveillance initiatives, recognising the problems associated with differentiating between "self" and "exogenous" cytoplasmic RNA (39). Indeed, it is feasible that the type 1 IFN signature identi ed during training may be attributable to sub-clinical viral infection in this cohort of horses. Further work on a more geographically diverse population would be required to address this possibility. Despite the transcriptomic and proteomic evidence of altered immunity during the training period, the causal factors which underpin these changes remain speculative. Such causal factors may act in isolation or in concert and may include the following: repetitive periods of intense exercise, changes in airborne environment and subclinical infection. Indeed, both infectious and environmental causes have been proposed for the increased risk of airway in ammation within the early training period, potentially re ecting the comingling of horses from diverse locations and increased housing associated with the transition into this period (6, 64-66).
Previous studies have identi ed an association between the bacterial load within the trachea and both cytological and endoscopic evidence of airway in ammation (6, 64, 66). Furthermore, in addition to the early training period, the detection of bacteria (streptococcal species) within the trachea was associated with evidence of airway in ammation; however, this latter association was not con rmed as causal, potentially re ecting a more indirect association (e.g. increased colonisation of mucus following compromised clearance). In our study, the timing of the second sample collection was such that this "early training" association should have subsided, which, together with the limited cytological or endoscopic evidence of signi cant airway in ammation offers more assurance that the transcriptomic and proteomic data was not solely re ective of the early transition to the training environment. Viral infections should also be considered as a potential explanation for the association between training and airway in ammation, particularly in light of the speci c IFN pathways altered by training. Indeed, previous studies have shown that viral infections are common among young horses after entry into training (67, 68). However, attempts to identify viral RNA based on proteomic analysis, failed to identify any appropriate candidates. If a viral cause was indeed present, it was not associated with overt clinical signs, nor did it appear to impact athletic performance.
Previous studies have demonstrated a clear association between housing and airway in ammation, re ective of the increased exposure to organic dust (69, 70). As horses were housed during the training period and largely at pasture during the T0 sampling, it is quite feasible that some of the immunological changes observed, at both the individual gene and protein level and the pathway level, were consistent with an increased exposure to organic dust; however, this was not universally or consistently re ected in the TW differential cytology data. The ventilation within the stable was considered to be good and all horses were bedded on low dust shavings and fed haylage from the ground. Additional data relating to dust exposures would have assisted the interpretation of the data in light of the change in housing; however, this was not possible within the constraints of the study design. Lastly, in line with previous studies (10), the immunological alterations associated with training may have been attributable to the repeatable episodes of high intensity exercise. Frellstedt et al (2014) previously demonstrated trainingassociated changes in immunity at the level of the lower airways in treadmill exercised horses (10).

Training induced cellular stress
As well as the potentially direct effect of training on airway immunity, various transcription modulation processes also seem to play an important role. Cappelli and others (2020) recently demonstrated the dominant role of transcription modulation in orchestrating the genomic response of equine PBMCs to exercise induced stress (71). Muscle ATP demand increases with increasing exercise intensity with the rate of production dependent on the availability of oxygen, carbon substrates, Ca 2+ .and other molecules including, nitric oxide and reactive oxygen and nitrogen species (RONS) (72). During intense exercise, increased RONS production may cause oxidative stress and damage to cellular structures and reduce mitochondrial e ciency, resulting in in ammation and transient immune dysfunction (73).
In the current study, pathway analysis detected several differentially expressed genes related to oxidative stress, oxidative phosphorylation, mitochondrial dysfunction and protein synthesis (Fig. 4, 5; Additional le 3, 4). Reactive oxygen species modulator 1 gene, responsible for increasing the level of reactive oxygen species (ROS) in cells, was signi cantly upregulated and superoxide dismutases (SOD), the rst line of defence against superoxide radicals (73), were downregulated during the training period. At the protein level, both SOD1, with the capacity to limit the detrimental effects of ROS and apoptotic signalling, and SOD2 were signi cantly downregulated during the training period. Proteomic analysis also revealed downregulation of other antioxidants; these included thioredoxin (TXN) and various heat shock proteins (HSP90B1, HSPA5, ST13, HSPD1, HSPE1, HSPB1) (74). TXN has the capacity to directly interact with several transcription factors including NFKB (74), also regulated by ROS and which itself regulates the transcription of a host of acute phase, cytokine and cell surface receptor genes; in some systems, antioxidants have been shown to reduce or block NFKB activation (74). Our data showed the trainingassociated differential expression of several genes (NFKBIZ, CARD19, NKAP) involved in the NFKB pathway. ROS can also activate MAP kinases via a Ras-dependent mechanism (74)

Haemopoietic abnormalities during training
In light of the very high prevalence of EIPH in racing Thoroughbreds (approaching 100%) (13), it is unsurprisingly that most of the horses (~ 69%, 11 out of 16) developed cytological evidence of EIPH during the training period. The pulmonary pathology associated with EIPH re ects changes in the pulmonary interstitium and vasculature, resulting in interstitial oedema and brosis (75,76); yet, to date, none of the EIPH studies have examined the transcriptomics of airway derived cells harvested from cases. Analysis of gene expression data derived from the current study revealed minimal differences between EIPH affected (n = 11) and unaffected (n = 5) horses during the training period (data not shown); however, this has re ected the small group sizes following subdivision of the T1 cohort; however, several genes related to haemopoietic abnormalities were detected in both sub-groups, suggesting that bleeding may have occurred even in those cases with no clear cytological evidence of such. Our data does however support the potential value in expanding the transcriptomic analyses to a larger sample population in order to more critically assess the gene expression pro le of EIPH cases.
Transcriptomic pro le of racehorses with high neutrophil count on TW samples Although only a small number of horses developed a signi cant airway neutrophilia during the training period, supplementary comparative analysis was applied to the data derived from the three horses with the greatest magnitude of TW neutrophilia, compared with the other horses at time point T1. This revealed the differential expression of 57 genes, a number similar to that reported in circulating leukocytes derived from both Standardbred and endurance horses and BALF samples from horses with airway in ammation (28, 77). Also consistent with previous reports, almost all differentially expressed genes were upregulated, with only three being downregulated. The differential expression of IL1RN in the group of horses with airway in ammation (Additional le 7), was consistent with the ndings of Hansen et al (2020) who reported an association between IL1RN expression in BALF derived cells and the neutrophilic form of MMEA (78). Similarly, genes cysteine rich secretory protein LCCL domain containing 2 (CRISPLD2) and transglutaminase 3 (TGM3) were also present in BALF cells derived from horses with neutrophilic form of MMEA (28). Interestingly, the neutrophilic group also showed increased expression of the mucin gene MUC5AC, shown previously to contribute to mucus hypersecretion in severe equine asthma (79) and also highly expressed in bronchial biopsies from asthmatic humans (80). Overall, the transcriptomic changes in this small group of horses were consistent with the airway cytological pro les observed, re ecting pathways associated with neutrophil chemotaxis, response to stress and immune defence (49). Interestingly, prior comparative transcriptomic analyses of endotracheal biopsies derived from asthmatic horses also revealed the differential expression (compared with non-asthmatic horses) of a similar gene list to that identi ed in the current study. This gene set included CSF3R, CXCR2, PLEK, IL1RN, RETN, VSTM1, TREM1 and PLAUR, and is related to neutrophil chemotaxis and immune responses. Colony stimulating factor 3 is responsible for neutrophil production and proliferation and promotes neutrophil tra cking by modulating chemokine and adhesion receptors, such as CXCR2 that mediates neutrophil migration to sites of in ammation. Although di cult to detect in blood, CSF3 is secreted to a greater extent during infectious or in ammatory conditions (81-83) including human asthma where it was positively associated with BALF neutrophil numbers and sputum from asthmatic patients (80). Moreover, CSF3 is considered as a potential target against pathological in ammation and tissue remodelling occurring in human asthma (80). Other genes detected in this group, such as TREM1 and PLAUR have also been reported in human asthma (59,84).
A number of Th1 immune response agents were also enriched in the neutrophilic group of horses (CXCR2, IDO1, CSF3R, ISG20, IL18RAP and TREM1). Thus, it appears that the equine airway neutrophilia detected in these cases was associated with activation of the innate immune system, with a likely Th1 polarised response. Similar observations have been previously reported, although, within the context of MMEA, suggesting that immunological pathways vary according to the cytological pro le of the airway in ammatory response and potentially the phase of the disease (85).
Despite limitations of the current study, this is the rst whole genome sequencing study performed on equine TW derived cells. Results shared important similarities with previous equine and human studies and although the aim of the study was focussed on the impact of training, transcriptomic analysis detected a number of novel genes potentially related to disease pathogenesis which may provide new insights for future studies. It remains feasible that the assessment of gene expression in TW-derived cells may represent a relatively non-invasive method to identify molecular endotypes of equine asthma in larger studies. However, detailed phenotypic characterisation of the horses included in these studies would be crucial. Furthermore, data derived from such analyses, combining proteomic and transcriptional changes has a real capacity to improve our general understanding of the pulmonary immune system of the horse.

Conclusions And Future Perspectives
Although this study focussed on the temporal changes associated with entry into a training program, the combined -omics approach could also provide a more de nitive and comprehensive characterization of lower airway diseases. Although, in this respect, further work comparing TW and BALF samples, similar to that already reported in relation to cytological correlations, is warranted. Finally, as previously suggested, many similarities were observed with human-derived data in relation to exercise immunology, further supporting the use of the horse as an attractive animal model from which translational application of ndings to humans may be justi ed.