Porcine respiratory coronavirus as a model for acute respiratory disease: mechanisms of different infection outcomes

Abstract

Porcine respiratory coronavirus (PRCV) is a naturally occurring pneumotropic coronavirus in the pig, providing a valuable large animal model to study acute respiratory disease. PRCV pathogenesis and the resulting immune response were investigated in pigs, the natural large animal host. We compared 2 strains, ISU-1 and 135, which induced differing levels of pathology in the respiratory tract to elucidate the mechanisms leading to mild or severe disease. The 135 strain induced greater pathology which was associated with higher viral load and stronger spike-specific antibody and T-cell responses. In contrast, the ISU-1 strain triggered mild pathology with a more balanced immune response and greater abundance of T regulatory cells. A higher frequency of putative T follicular helper cells was observed in animals infected with strain 135 at 11 days postinfection. Single-cell RNA-sequencing of bronchoalveolar lavage revealed differential gene expression in B and T cells between animals infected with 135 and ISU-1 at 1 day postinfection. These genes were associated with cell adhesion, migration, and immune regulation. Along with increased IL-6 and IL-12 production, these data indicate that heightened inflammatory responses to the 135 strain may contribute to pronounced pneumonia. Among bronchoalveolar lavage (BAL) immune cell populations, B cells and plasma cells exhibited the most gene expression divergence between pigs infected with different PRCV strains, highlighting their role in maintaining immune homeostasis in the respiratory tract. These findings indicate the potential of the PRCV model for studying coronavirus-induced respiratory disease and identifying mechanisms that determine infection outcomes.

Graphical Abstract

Open in new tabDownload slide

Introduction

Respiratory coronaviruses (CoVs) are a significant threat to global health and have caused 3 major human epidemics of severe respiratory disease since 2003, including the recent SARS-CoV-2 pandemic. In each case, CoVs emergence in humans has been associated with zoonotic transmissions from animals. Many livestock species are infected with CoVs causing variable morbidity and mortality, associated with lung disease and immunological impairment, leading to substantial economic losses.^1,2 Porcine respiratory CoVs (PRCVs) are globally endemic in pigs. PRCV is a naturally occurring spike deletion variant of the Alphacoronavirus Transmissible Gastroenteritis Virus (TGEV).^3,4 TGEV formerly caused enteric disease with high mortality in neonatal piglets, but since the emergence of PRCV, the impact of TGEV has been greatly reduced owing to cross-protective immunity.^5,6 PRCV infection is usually subclinical with most cases being mild unless secondary infections occur.

Previously it was suggested that PRCV infection with ISU-1 strain induces pulmonary pathology similar to that of SARS-CoV.⁷ However, the lung pathology was mild unless potentiated with corticosteroids.⁷ Given that pigs are large, long-lived natural hosts for CoVs with respiratory anatomy and physiology similar to humans, we investigated the effects of different PRCV strains in vivo. We showed varying degrees of lung pathology: PRCV 135 strain caused a prominent pneumonia similar to human SARS-CoV-2, with viral tropism and damage to the bronchial and alveolar cells; in contrast, the pathological changes driven by ISU-1 infection were very mild.⁸ Although PRCV strains 135 and ISU-1 replicated equally well in continuous porcine cell lines in vitro, only the more pathogenic 135 strain was able to replicate in ex vivo tracheal organ cultures, as well as in vivo in the trachea and lung.⁸ PRCV uses aminopeptidase N (APN) as a receptor, although recent studies suggested another host entry or host restriction factor might be involved.^9–11 Differences in the spike and gene 3 sequences between the strains maybe also responsible for the variation in replication patterns.¹²

Mouse, ferret, hamster, and nonhuman primate models have been developed to study SARS-CoV-2. However, these species are not natural hosts for the virus, and it remains uncertain which model most accurately mimics human disease.^13–15 In contrast, the pig is a natural host for several porcine CoVs including pneumotropic viruses.^1,16,17 As with SARS-CoV-2, PRCV infection can be mild or asymptomatic, but in some instances can lead to pneumonia. Better understanding of the mechanisms that lead to mild or severe disease is essential for developing new control strategies for emerging severe coronavirus diseases of livestock or humans. Here we used infection with the less and more pathogenic PRCV strains ISU-1 and 135, respectively, for in-depth mechanistic evaluation of the pathogenesis, virology, and immune responses of this important family of viruses. We performed high resolution analysis of early and late systemic and mucosal immune responses using multiparameter flow cytometry to reveal how effector functions of the cellular and humoral responses mediate lung pathology. We also performed single-cell RNA sequencing (scRNA-seq) on cells from bronchoalveolar lavage (BAL), as we previously did to characterize changes during influenza infection and respiratory immunization.¹⁸ We used the scRNA-seq data to interrogate early (day 1) and late (day 20) transcriptional changes following infection with ISU-1 and 135 PRCV strains. We propose that the pig is an appropriate and powerful model for understanding immunity to CoVs in a natural large animal host, particularly the early innate and subsequent adaptive immune responses at the site of infection, which are critical for protective immune responses and vaccine efficacy.

Materials and methods

Viruses

The previously characterized PRCV strains 86/135308 (referred to as 135, GenBank accession number OM830318), ISU-1 (GenBank accession number OM830321), 86/137004 (referred to as 137, GenBank accession number OM830320), AR310 (referred to as 310, GenBank accession number OM830319), and TGEV strain FS772/70 were used.¹ All viruses were grown and quantified in swine testis (ST) cells and titers expressed as the number of PFU per milliliter (PFU/ml).

Ethics statement

Animal experiments were approved by the ethical review processes of The Pirbright Institute and the Animal and Plant Health Agency (APHA) according to the UK Animals (Scientific Procedures) Act 1986 under project license PP7764821.

Animal study

Forty outbred female pigs, aged 4 to 6 weeks and genetically composed of one-quarter Large White, one-quarter Landrace, and one-half Hampshire breeds, were included in the study. Prior to the challenge, all animals were assayed for the presence of spike antibodies to ensure their PRCV-free status and randomly assigned into 2 groups of 20 pigs. One week after acclimatization, one group was inoculated with 1 × 10⁷ PFU PRCV 135 strain and the other with 1 × 10⁷ PFU PRCV ISU-1 strain in a total of 4 ml (2 ml per nostril) using a mucosal atomisation device (MAD). Clinical parameters including demeanor, appetite, respiratory signs, sneezing, coughing, nasal and eye discharge, feces consistency, and rectal temperature were assessed. One ISU-1- and five 135-inoculated animals showed mild clinical signs between 2 and 5 days postinfection (DPI), which resolved rapidly and did not exceed a score of 2 out of maximum 19. Animals were humanely culled at 1, 5, 11, and 20 DPI with 5 animals sampled from each group at each time point for assessment of viral load, pathology, and immune responses. Daily nasal swabs were collected for evaluation of viral shedding. At postmortem, BAL, tracheobronchial lymph nodes (TBLN), and spleen were collected, mononuclear cells isolated and cryopreserved as previously described.¹⁹ Because of a technical failure no spleen cells were available to perform intracellular cytokine staining (ICS) at 20 DPI. Viral load was also assessed in BAL, lung, trachea, tracheal swabs, and eyelids.

Lung gross pathology, histopathology, and immunohistochemistry

The percentages of pulmonary consolidation, including both dorsal and ventral aspects, were assessed blindly by a veterinary pathologist at necropsy.¹⁹ Formalin-fixed tissues were processed by a routine histology method. Hematoxylin and eosin (H&E) staining and immunohistochemistry (IHC) against PRCV nucleoprotein (N) using a monoclonal antibody were performed on serially sectioned formalin-fixed paraffin embedded tissues.¹⁹ Positive control slides containing either TGEV or PRCV-infected or uninfected cells, as well as pig tissues derived from PRCV-challenged experiments, were included in IHC experiments to confirm specificity of immunolabelling.

Virological assays

The evaluation of infectious viral progeny involved the analysis of harvested BAL fluid, tracheal and nasal swabs obtained postmortem, and daily nasal swabs. All samples, stored in Viral Transport Medium (VTM), were subjected to plaque assay in ST cells using a previously established protocol.¹⁹ Tissues collected at post-mortem, including trachea, lung, cheek, and eyelid, were homogenized in phosphate-buffered saline (PBS) using a Tissuelyser II (Qiagen) with a 5-mm stainless steel ball (Qiagen). The resulting tissue-derived supernatant was titrated in a plaque assay in ST cells. Prior to homogenization, all tissues were weighed, with 0.25 g for trachea, cheek, and eyelid tissues, and 0.5 g for lung tissue. The supernatant was clarified through low-speed centrifugation. The trachea was dissected into 3 sections: top, middle, and bottom. The final value was calculated as the average (mean) of the titers generated from each section.

Serological assays

Ninety-six-well Maxisorp ELISA plates (BioLegend, UK) were coated overnight at 4 °C, with pre-optimized concentrations of ISU-1 full-length spike (S) protein (1 µg/ml). Two-fold dilutions of serum, BAL, or nasal swab, starting from a suitable dilution in PBS-Tween plus 4% semi-skimmed milk were applied. Samples were collected either pre-inoculation (day 0) for serum only or at postmortem. Binding of antibodies was detected with goat anti-pig IgG (Fc-specific) (AA141P, Bio-Rad, Watford, UK) or goat anti-pig IgA (AA140P, Bio-Rad, Watford, UK), conjugated to horseradish peroxidase at optimal dilutions. 1-Step Ultra TMB-ELISA Substrate Solution (Thermo Scientific Pierce) was added, and OD values were read for each well at dual wavelengths (450 nm and 630 nm) using an ELx808 Absorbance Microplate Reader (Agilent Technologies). Quantities of antibodies were determined as the reciprocal value of the dilution that yielded the first reading above the cutoff value (endpoint titer). Cutoff values were established as mean blank ODs plus 2-fold standard deviations.

Virus neutralization assay

Virus neutralizing antibodies against all 4 strains of PRCV and TGEV involved were determined in serum and BAL fluid as previously described.¹ Briefly 10³ PFU of the respective virus was incubated with an equal volume of a 2-fold dilution series of serum, starting at 1:5 dilution for serum from 135-infected animals and 1:2 for serum from ISU-1-infected animals. After 30-min incubation at room temperature the virus/serum mix was transferred to ST cells. After 48-h incubation at 37 °C, cells were fixed and stained using 3.3% formaldehyde and 0.1% crystal violet. The presence or absence of CPE was noted, with the absence of CPE indicating neutralization. Neutralization titers were calculated as the quantity required to achieve 50% neutralization per milliliter (NT50/ml) using a Reed and Muench endpoint calculation.²⁰

ELISPOT assay

MultiScreen-HA ELISPOT plates (Millipore, Watford, UK) were coated overnight at 4 °C with 0.5 µg/ml mouse anti-pig IFN-γ (BD Biosciences, San Jose, California, United States) in carbonate buffer (0.1 M Na2CO3/NaHCO3, pH 9.6). Following 3 washes with PBS, the plates were incubated for 2 h at 37 °C with 200 µl per well of blocking buffer (RPMI, 10% FBS, supplemented with 100 U/ml penicillin and 100 µg/ml streptomycin). After removing the blocking buffer, stimuli and cells were plated. Cryopreserved cells (TBLN, BAL, spleen) were thawed, and 100 µl/well at 2 × 10⁶/ml cells were added to the ELISPOT plates and cultured for 48 h at 37 °C, 5% CO2 in triplicate for each condition. Stimuli included medium alone (RPMI-1640 medium with GlutaMax-I and 25 mM HEPES, supplemented with penicillin, streptomycin, and 10% heat-inactivated FCS), live PRCV viruses 135 and ISU-1 at multiplicity of infection 1 (MOI: 1), peptide pools covering the spike (S) of PRCV 135 virus, or 2.5 µg/mL Concanavalin A (Con A, Invitrogen). Peptides (16-mer) overlapping by 12 residues (Mimotopes, Melbourne, Australia) with 3 pools for the S protein (residues 1-100, 101-200, and 201-305) at a final concentration of 2 µg/ml were used.

After washing of the plates with PBS containing 0.05% Tween 20, 100 µl per well biotinylated mouse anti-pig IFN-γ (0.5 µg/ml, BD Biosciences) was added for 2 h at room temperature. The plates were washed 5 times with PBS containing 0.05% Tween 20, and 100 µl per well streptavidin conjugated to alkaline phosphatase (1/1,000, Invitrogen Ltd) was added for 1 h at room temperature. Following another 5 washes with PBS containing 0.05% Tween 20, the plates were developed with 100 µl per well alkaline phosphatase substrate solution (Bio-Rad Laboratories, Hercules, CA) for 20 min at room temperature. After rinsing with tap water, the plates were allowed to dry overnight at room temperature before counting the dark blue-colored immunospots using the AID ELISPOT reader (AID Autoimmun Diagnostika GmbH, Strassberg, Germany). Results were expressed as the number of IFN-γ-producing cells per 10⁶ stimulated cells after subtracting the average number of spots in medium-stimulated control wells.

Intracellular cytokine staining

To assess the production of IL-2, TNF, and IFN-γ by CD8 and CD4 T cells, ICS was performed. Cell suspensions obtained from BAL, TBLN, and spleen were thawed and seeded in duplicate wells of a 96-well plate (1 × 10⁶ cells per well). The cells underwent overnight stimulation with either 135 or ISU-1 (MOI 0.5), or they were treated with media (RPMI containing stable glutamine, 100 IU/ml penicillin, 100 μg/ml streptomycin, and 10% FCS) as a negative control. Incubation occurred at 37 °C with 5% CO₂. Alternatively, some cells were stimulated for 5 h with 2 µg/ml of 100 overlapping peptides (16-mers) spanning the first 100 residues of the PRCV 135 spike protein. Positive control wells were stimulated with a 1/500 dilution of a cocktail containing phorbol 12-myristate 13-acetate (PMA) and ionomycin (BioLegend). Brefeldin A (GolgiPlug, BD Biosciences) was added to all wells and incubated for 4 h. Subsequently, cells were centrifuged for 5 min at 400 × g and washed twice with PBS.

Primary antibodies (Table S9, section A) were applied to the cells for 20 min at 4 °C in the dark. The cells were then permeabilized and fixed using 100 μl of BD Cytofix/Cytoperm (BD Biosciences). Following another centrifugation at 400 × g and 2 PBS washes, the cells were incubated for 30 min at 4 °C in the dark with directly conjugated anti-IFN-γ, anti-TNF, and unconjugated anti-IL2 antibodies (Table S9, section A). After an additional centrifugation and 2 washes with PBS, the cells were treated with fluorochrome-conjugated secondary antibodies. Following 2 more washes and resuspension in PBS, cytokine-producing CD4 and CD8 T cells were analyzed using a MACSquant Analyzer 16 (Miltenyi). Wells containing only cells and medium for each animal were considered as a negative control (unstimulated), and the frequency of cytokine-producing cells was determined by subtracting the values from unstimulated cells. Data analysis was performed using FlowJo software version 10.10.0.

Phenotypic analysis of B and T cells by flow cytometry

Thawed cell suspensions from BAL, TBLN, and spleen were resuspended in 10 ml of medium (RPMI containing stable glutamine, 100 IU/ml penicillin, 100 μg/ml streptomycin, and 10% FCS) and centrifuged at 400 × g for 5 min. Subsequently, 2 × 10⁶ cells per sample from BAL, TBLN, and spleen were seeded into a 96–U-bottom-wells plate and washed twice with PBS. A cocktail of primary antibodies for either B-cell or T follicular helper cell (Tfh)/regulatory T cell (Treg) phenotyping was added to the wells (Table S9, section B and section C, respectively) and incubated for 20 min at 4 °C. The plates were washed twice with PBS, and secondary antibodies added, followed by incubation for 20 min at 4 °C. In cases requiring blocking, mouse serum was added. For fixation and permeabilization, eBioscience Foxp3/Transcription Factor Fixation/Permeabilization Concentrate and Diluent (Invitrogen) were employed for intracellular markers. Following 2 PBS washes, antibodies for intracellular markers were added (Table S9, section B and C) to the wells and incubated for 30 min at 4 °C. Finally, the cells were resuspended in PBS. The samples were analyzed using either the MACSquant Analyzer 16 (Miltenyi) or Cytek Aurora (Cytek Biosciences). FlowJo software was employed for the analysis of the collected data.

Luminex assay

Cells from BAL at 1 and 5 DPI were thawed and resuspended in duplicate wells of a 96-well plate (1 × 10⁶ cells per well). The cells were stimulated overnight with the homologous viruses either PRCV 135 or ISU-1 at MOI 1, or they were treated with medium (RPMI containing stable glutamine, 100 IU/ml penicillin, 100 μg/ml streptomycin, and 10% FCS) as a negative control. The supernatant was collected and used for Luminex assay. Cytokine concentrations for IL-1β, IL12p40, IL-4, IL-6, IL-8, IL-10, TNF-α, IFN-γ, and IFN-α were determined by Swine Cytokine Magnetic 9-plex Panel (ProcartaPlex Multiplex Immunoassay, Thermo Fisher) performed according to the manufacturer’s recommendations.

Separation of BAL cell populations for scRNA-seq

BAL cells were isolated and subsequently cryopreserved as previously described.^18,21 To account for the strong domination of macrophages in BAL cell preparations, we sorted major leukocyte subsets as described in ref. 18. In brief, thawed BAL cells were labelled with antibodies against CD4 (clone 74-12-4, PerCP-Cy5.5 conjugated, BD Biosciences), CD8β (clone PPT23, FITC conjugated, Bio-Rad), CD3 (clone PPT3, APC conjugated, Southern Biotech), CD16 (clone G7, AF647 conjugated, Bio-Rad), CD172a (clone 74-22-15, PE conjugated, Bio-Rad) and Fixable Viability Dye eFluor780 (ThermoFisher) and sorted on a FACSAria UIII (BD Biosciences). This allowed the sorting of CD4⁺, CD8β⁺, CD16⁺CD172a⁺, and CD3^-CD16^- cells (the latter phenotype for enrichment of B cells). Per phenotype, 12,000 cells were collected. Cells from each fraction were processed for partitioning and barcoding using a 10x Genomics Chromium iX Controller.

Library preparation and sequencing

Single cells were processed using the Chromium iX Controller with the Next GEM chip G and the single cell 3’ kit v3.1 library preparation kit from 10x Genomics (Pleasanton, California, United States). For each of the BAL samples, approximately 10,000 cells were targeted for recovery during partitioning. The samples were grouped randomly and run in pairs, using 3 NextSeq High Output 150-cycle reagent kits (Illumina, San Diego, California, United States), with 1% PhiX, aiming for about 200 million reads per sample, with the run format aligning with 10x guidelines. The sequencing was performed on an Illumina NextSeq 550, using 3 high output 200-cycle reagent kits and flow cells. Bcl files were demultiplexed using “cellranger mkfastq” and aligned/counted with “cellranger count” (cellranger-7.0.0), using default parameters and the Sus scrofa genome (genome assembly 11.1, Ensembl release 108).

scRNA-seq quality control and pre-processing

Data analysis was conducted using R (version 4.3.2). To ensure data integrity, we excluded empty droplets and barcode-swapped droplets from the unfiltered CellRanger output using DropletUtils (1.22.0).^22,23 Further filtering involved the removal of cells with total Unique Molecular Identifiers (UMIs) counts exceeding 50,000 or falling below 500. Genes without any counts and cells exhibiting an outlying proportion of mitochondrial UMIs (greater than 5 median absolute deviations from the median) were also excluded.

Normalization was achieved using pooled factors,²⁴ and samples from the same sequencing run were merged to form a single analysis dataset. Specifically, samples 1 to 8 from sequencing run 1 were consolidated to create the 1 DPI analysis dataset, while samples 9 to 16 from sequencing run 2 were combined to form the 20 DPI analysis dataset. Batch effect correction was deemed unnecessary as visualization and clustering revealed no batch effects within each run.

scRNA-seq clustering and annotation

Shared-nearest-neighbors graph-based clustering was performed on each analysis dataset using clusterCells (bluster 1.12.0) with the top 5,000 highly variable genes (HVGs) identified by batchelor (1.18.1).²⁵ The clustering was conducted using the Jaccard index to weight edges and the Louvain method for community detection, with a k value of 20.

For the 1 DPI dataset, the clustering process generated 19 initial clusters. Cluster 18, which exhibited particularly low counts, had a limited number of genes detected, and showed dispersive distribution in uniform manifold approximation and projection (UMAP). Therefore, it was considered as poor quality and excluded from further analysis. Cluster 14 in the 1 DPI dataset displayed a dichotomy of CD4 and CD8 expression. This cluster was further split into CD4, CD8, and NK subclusters. The subclustering process involved identifying the top 1,000 HVGs within the cluster, selecting those genes that correlated best with CD4 and CD8B expression using a zero-inflated Kendall Tau correlation metric implemented in scHOT (v1.14.0), and re-running the clustering command using Jaccard index to weight edges and the fast greedy modularity optimization algorithm for community detection using only those HVGs for cluster 14. Similarly, in the 20 DPI dataset, cluster 9 exhibited a dichotomy of CD4 and CD8 expression. This cluster was split into CD4+ and CD8β+ subclusters using the same strategy as described for cluster 14 in the 1 DPI dataset.

In the 1 DPI dataset, varied expression of genes including KLRB1 and NCR1 within cluster 9 prompted subclustering. This was performed in an unsupervised manner by identifying the top 2,000 HVGs in cluster 9 cells and then re-running the clustering using these HVGs. The Jaccard index was used to weight edges, and the fast greedy modularity optimization algorithm was applied for community detection, resulting in the division of cluster 9 into 2 subclusters. Using the same approach, cluster 12 was split into 3 subclusters, and cluster 3 into 2 subclusters.

The 20 DPI dataset generated 20 initial clusters. In this dataset, clusters 2, 4, 15, and 19 were each split into 2 subclusters using the unsupervised method described above to better distinguish and define cell types within. Subclusters containing fewer than 100 cells were excluded from further analysis owing to their inability to produce significant statistical results.

Clusters were manually annotated with conventional names (cell types) by examining expression of typical defining markers, as listed in Table S1 and illustrated in Fig. S4. Additionally, top marker genes within each cluster were determined using the score Markers function (scran 1.30.2), as detailed in Tables S10 and S11.

scRNA-seq analysis, differential abundance analyses

To compare the cell type proportion (abundance of cells within each cluster) across 2 virus infection conditions, we performed differential abundance analyses using a negative binomial generalized linear model with empirical Bayes quasi-likelihood F-tests,^26,27 as implemented in edgeR (v4.0.16). The model was fitted with 2 conditions (135 vs. ISU-1) to test for differential abundance in condition 135, using ISU-1 as the baseline. Trend estimation was disabled when estimating dispersion and quasi-likelihood dispersion.

To account for the sort and pool strategy used to balance cell types among the sequenced cells, clusters were divided according to their membership in the 4 sorted cell populations (macrophages, CD4 T cells, CD8 T cells, B cells) before running differential abundance analyses within each group (Table S2). Clusters that did not clearly fall into one of these cell types were excluded to ensure the universes for differential abundance analysis were equivalent to or smaller than the sorted cell populations. This approach was intended to avoid sorting biases in intra-cell-type differential abundance analyses.

scRNA-seq analysis, differential expression analyses

Owing to the limited number of subjects in each condition, we conducted pseudobulk differential expression analyses within each cluster using edgeR (Figs. S7, S8, and S10). To compare between 135 and ISU-1 infections, counts from all cells within each cluster were pooled for each sample and clusters with fewer than 10 cells within each sample were removed. A negative binomial generalized linear model was then fitted for each cluster, and differential expression between conditions was tested using empirical Bayes quasi-likelihood F-tests, as implemented in edgeR (v4.0.16). The differential expression results are provided in Tables S12 and S13.

Afterwards, we first combined all genes that were differentially expressed in any cluster within each major cell population (eg CD4 T cells) (Table S5). Then, we created Venn diagrams to display those differentially expressed genes within 4 main cell populations (macrophages, CD4 T cells, CD8 T cells, and B cells). We also checked the direction of effect in clusters within each main population; these directions were all consistent across clusters from the same main population.

Gene ontology and Kyoto Encyclopedia of Genes and Genomes enrichment

Lists of differentially expressed genes were analyzed for Gene Ontology (GO) enrichment using clusterProfiler (4.10.1), with the genome-wide annotation for pig (DOI: 10.18129/B9.bioc.org.Ss.eg.db, release 3.19) focusing on the biological process sub-ontology.

When using clusterProfiler, hypergeometric testing was employed to identify GO terms significantly enriched in the set of DE genes in each cluster compared to a background set of all measured genes. To account for multiple testing, the Benjamini-Hochberg (BH) method was applied to adjust P-values. GO terms with an adjusted P-value (q-value) <0.05 and annotated by at least 3 DE genes were considered statistically significant and are displayed in bubble plots. All results are available in Tables S7 and S8. Since results of 1 DPI only contained very few significant GO terms in per cell type, no bubble plots were created. The same lists of DE genes were analyzed for Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment using clusterProfiler and adjusted for multiple testing by the BH method. The resulting significantly enriched KEGG pathways were manually filtered to concentrate on immune-related pathways and displayed in bubble plots. All results are available in Tables S7, S8, and S14.

Statistical analysis

GraphPad Prism 9.2.0 (GraphPad Software, San Diego, California, United States) was used for statistical analyses. The specific type of analysis used is indicated in each figure legend. Data were first tested for normality. For normally distributed data sets a one-way or two-way ANOVA was applied followed by Bonferroni’s multiple comparisons test. For non-normally distributed data a nonparametric unpaired t test was used. A P-value of <0.05 was considered to be statistically significant. To assess correlation between virological/immune parameters and pathological measures, a nonparametric Spearman correlation coefficient (ρ) was computed for each pair of measures (viral load, antibody response, T-cell responses and pathology) using GraphPad Prism v 10.0.3. All 135- and ISU-1-infected samples at the specified time points were included in the analysis.

Results

Viral load and respiratory pathology following PRCV ISU-1 and 135 infections

Forty pigs were inoculated intranasally/intratracheally with 1 × 10⁷ PFU of either PRCV ISU-1 or 135 strains. Five pigs from each group were humanely culled at 1, 5, 11, and 20 DPI (Fig. 1A). The quantity of infectious virus present in daily nasal swabs was determined by plaque assays. Peak viral shedding occurred between 2 to 6 DPI for both ISU-1- and 135-infected animals, followed by a rapid decline from 6 DPI. Viral shedding of 135-infected animals was significantly higher on 1, 4, 5, and 6 DPI compared to ISU-1-infected animals (Fig. 1B). Total viral shedding over the time course of infection was also significantly greater in the 135-infected animals as determined by the area under the curve (Fig. 1C). The 135-infected animals exhibited higher viral load in tracheal swabs, tracheal tissue, BAL, and lung, with a maximum level detected at 5 DPI and return to baseline at 11 DPI (Fig. 1D–G). Two out of 5 animals showed ISU-1 virus in BAL, 3 out of 5 in tracheal tissues, and 2 out of 5 in lungs at 1 or 5 DPI. Infectious virus was detected in the eyelid of the 135 group (4/5 pigs), which was comparable with studies reported from SARS-CoV-2 (Fig. 1H).²⁸ A section of skin from the cheek (masseter region) was evaluated, but no infectious virus was isolated (data not shown) indicating that the observed replication of 135 is specific for the conjunctiva rather than skin.

At necropsy, pulmonary consolidation was only observed in the 135 group, with an average of 1.2% lung areas affected at 1 DPI, 11.2% at 5 DPI, a reduction by 11 DPI to 5.8%, and absence of overt consolidation at 20 DPI (Fig. 2A). Only the cranial and middle lung lobes were affected. No gross lung lesions were noted in the ISU-1 group throughout the study period. At 5 DPI, the ISU-1-infected upper respiratory tract had mild rhinitis and tracheitis, whereas 135-infected animals exhibited overt pathology in both upper and lower respiratory tracts (Fig. 2B and C). Using IHC, PRCV antigen (nucleocapsid) was detected in the nasal turbinate (respiratory and olfactory aspects of nasal mucosa) of both groups at 1 DPI, with higher levels in 135-compared to ISU-1 infected animals and consistently remained elevated until day 5 DPI. Viral antigen was not detected at 11 and 20 DPI (Fig. 2B and D).

At 5 DPI, animals infected with 135 developed bronchiolar exudation and epithelial damage, and within the alveoli there were frequent neutrophil accumulation (4/5 pigs) and edema (2/5 pigs), all of which these lesions colocalized with virus antigen (Fig. 2B, C). In contrast, the ISU-1 group showed mild epithelial attenuation and rare neutrophilic exocytosis, and virus antigens were absent (Fig. 2B–D). The overt respiratory pathology in 135-infected animals at 5 DPI was positively correlated with higher viral load in tracheal swabs, lung, and BAL compared to the ISU-1 strain (Fig. 2E).

At 11 DPI, virus antigens were absent in the nasal turbinate, trachea, and lung of both groups (Fig. 2B–D). While virus antigens were not detected in the nasal turbinate in either group (Fig. 2B), mucosa epithelial damage with evidence of regeneration persisted in the 135 group. In the trachea, both groups showed mild exudation, and the 135 group exhibited low levels of submucosal lymphoplasmacytic infiltrates. In the lung, both groups displayed low-level bronchiolar epithelial exudation and attenuation. Additionally, there was type II pneumocyte hyperplasia in 135-infected pigs (2/5 pigs).

At 20 DPI, virus antigens were not detected in the respiratory tissues. There were scattered lymphocytic infiltrates in the nasal turbinate submucosa for both 135 and ISU-1, and in the tracheal submucosal more commonly in 135 pigs (3/5) compared to ISU-1 (1/5). In the lung, both groups had BALT formation and rare bronchiolar exudation, and type II pneumocyte hyperplasia in a 135-infected pig (1/5 pigs). In other tissues including submandibular salivary gland, olfactory bulb, duodenum, or colon, no virus antigens or lesions were detected.

Pathology evaluation demonstrated that the PRCV 135 strain induced greater pathology in the respiratory tract than ISU-1, which coincided with peak viral load at 5 DPI, and with minimal to mild histopathological changes persisting until 20 DPI. ISU-1 replication was largely restricted to the upper respiratory tract. Additionally, this was the first report of PRCV detection in eyelids of pigs, suggesting a potential virus dissemination from the respiratory tract to the eye.

Antibody and B-cell responses during ISU-1 and 135 infections

Antibody responses following ISU-1 or 135 inoculation were determined in serum, BAL, and nasal swabs. Spike-specific IgG and IgA were measured by endpoint titer ELISA and were detected in serum, BAL, and nasal swabs at 11 DPI and maintained at 20 DPI with significantly higher titers for IgG in the135-infected group compared to the ISU-1-infected group which were below the level of detection (Fig. 3A and B). IgA titers were lower compared to IgG. Serum from both 135- and ISU-1-infected animals cross-neutralized 135 and ISU-1 virus strains in addition to PRCV strains 137 and 310 and the parent TGEV strain (Fig. 3C). No neutralization of PRCV or TGEV was detected in uninfected serum samples (data not shown). BAL titers were much lower compared to serum, likely owing to the diluting effect of sampling. The histopathology in nasal turbinates and gross lung pathology at 11 DPI correlated with spike-specific antibody responses in BAL, serum, and nasal turbinates on the same day (Fig. 3D).

B-cell responses in BAL were assessed by flow cytometry based on a previously established staining panel²¹ (Fig. S1A). Cells expressing CD79-α, PAX5, IRF4, and Blimp1 were defined as plasmablast/plasma cells. Within the IRF4 and Blimp1 double-negative cells, the following subsets were further defined: proliferating B cells (Ki-67 positive), germinal center B cells (Bcl6 and Ki-67 double positive), and naïve/memory B cells (Bcl-6 and Ki-67 double negative). There was a temporal shift from B cells towards plasma cells in the BAL, suggesting the differentiation of B cells to plasmablast/plasma cells in response to both infections (Fig. 3E). However, this might also be due to the lower starting point in the ISU-1-infected animals, and no differences were detected at 5 and 11 DPI. The shift coincided with a greater abundance of cells with a phenotype characteristic of T follicular helper cells (CD4⁺ICOS⁺Bcl-6⁺) in 135-infected pigs at 11 DPI compared to the ISU-1-infected group (Fig. 3F). Since the effector function of these cells has not been characterized in the porcine model, we refer to them as putative Tfhs. A respiratory resident CD4 T-cell subset with phenotypic and transcriptional characteristics resembling both Tfh and tissue-resident cells has been identified in a murine model and was shown to play a crucial role in inducing and maintaining memory B cells and CD8 T cells specific to influenza antigens.²⁹

In summary a stronger antibody response was detected in the serum, BAL, and nasal swabs of 135-infected animals compared to ISU-1, which was dominated by IgG. Serum antibodies cross-neutralized all PRCV strains and TGEV. For the first time we identified putative Tfh cells in porcine BAL indicating a role in regulating local humoral responses.

Cytokine T-cell responses in BAL, TBLN, and spleen following ISU-1 and 135 infections

As T cells are crucial for control of virus replication, we examined CD4 and CD8β T-cell responses during ISU-1 or 135 infection. We analyzed cytokine-producing cells in BAL, TBLN, and spleen by ELISPOT and ICS (Fig. S1B) following stimulation with 135 or ISU-1 live viruses, or overlapping pools of peptides covering the spike (S) protein as previously described.⁸ An increased number of spike-specific IFN-γ-secreting cells was detected by ELISPOT in BAL of 135-infected animals at 20 DPI, compared to ISU-1-infected animals (Fig. 4A, Fig. S2). Similarly greater proportions of 135- and S-specific CD4 IFN-γ-secreting cells were observed by ICS at 20 DPI, while 135- and S-specific TNF-secreting CD4⁺ T cells were significantly higher at 11 DPI (Fig. 4B). The 135-specific CD8 T cells exhibited the highest IFN- γ production at 11 DPI and were significantly more abundant in 135-infected animals compared to ISU-1-infected animals (Fig. 4C). The percentage of 135-specific IFN-γ and TNF-secreting CD8⁺ T cells, as well as 135- and spike-specific IFN-γ ELISPOT responses, was associated with gross lung pathology and histopathology in nasal turbinates at 11 DPI (Fig. 4D). No differences in proportions of TNF-secreting CD8⁺ or IL-2-secreting CD4 and CD8 cells were detected between groups (Fig. S2A).

Extending ELISPOT and ICS to TBLN and spleen showed greater proportions of cytokine-producing cells in 135-infected animals (Figs S3 and S4). IFN-γ ELISPOT responses revealed significantly higher cytokine secretion in 135-infected animals at 11 DPI in TBLN and spleen (Figs S3A and S4A). ICS in TBLN demonstrated a higher proportion of virus- or spike-specific IFN-γ and TNF-secreting T cells at 11 DPI or 20 DPI, and spike-specific cells at 11 DPI in 135-infected animals (Fig. S3B). CD4 and CD8 IL-2-producing cells were also detected in TBLN with a significantly higher proportion in 135-infected animals at 20 DPI (Fig. S3C).

Taken together these data indicate that 135 infections induced greater IFN-γ, TNF, and IL-2 T-cell responses in BAL, TBLN, and spleen at 11 DPI and 20 DPI compared to ISU-1 infection. The greatest response was in the BAL with high frequencies of CD4 IFN-γ- and TNF-producing cells.

scRNA-seq analyses from BAL cells

Since BAL cell preparations are predominantly composed of alveolar macrophages, we applied cell sorting to obtain matching proportions of macrophages, CD4 T cells, CD8 T cells, and B cells (12,000 cells each) from 4 animals per condition (ISU-1, 135) and 2 time points (1 DPI and 20 DPI), resulting in a total of 16 samples.¹⁸ Cell partitioning and library preparation were performed separately for each sorted and re-pooled BAL sample.

To compare 135 and ISU-1 infection at each time point, samples were sequenced and analysed in DPI-specific batches. Graph-based clustering identified 19 and 20 distinct cell clusters in the 1 DPI and 20 DPI datasets, respectively. In the 1 DPI dataset, cluster 18 was of poor quality and was excluded from further analysis. In the 20 DPI dataset, cluster 19 (containing a mixture of nonconventional T cells and macrophages) was excluded because almost all cells in this cluster originated from sample 12. In all other respects sample 12 showed patterns of cluster membership consistent with the other samples.

Some clusters were additionally split into subclusters to separate CD4 from CD8 T cells, nonconventional T cells from NK cells, or macrophages from monocytes, creating a final dataset with 39,679 and 149,069 cells for annotation in the 1 DPI and 20 DPI datasets, respectively. To label each cluster with conventional cell-type names, we used a combination of marker genes characteristic of each cluster based on a list of commonly defining marker genes (Table S1). From this information, each cluster was manually annotated with an established cell-type name for the 1 DPI and 20 DPI datasets (Fig. 5A, Fig. S5A). In the 1 DPI dataset, subclusters 9b, 14b, and cluster 19 were enriched for NK cell-typical genes. The same applied to subcluster 2a and cluster 16 in the 20 DPI data set. Given that our cell sorting strategy did not include NK cells, these clusters were excluded to avoid a misrepresentation of these otherwise larger cell populations.

Differential abundance testing

We next tested the scRNA-seq data for differential abundance of cells within the sorted cell types (CD4 and CD8 T cells, B cells, macrophages) of ISU-1- versus 135-infected animals at each DPI (Tables S2 and S3, Fig. 5B, C). No differences were detected at 1 DPI, but at 20 DPI pigs infected with the ISU-1 strain had significantly higher proportions of cells in cluster 4b, annotated as Tregs, within the sorted CD4⁺ T-cell population (adjusted P = 1.1 × 10⁻⁴, Fig. 5B top). Tregs were further investigated by flow cytometry, based on a CD4⁺CD25^highFoxp3⁺ phenotype (Fig. S6A, Fig. 5B bottom). There was a significant increase of Tregs among CD4 T cells in ISU-1-infected pigs at 11 DPI compared to 1 DPI, in part driven by a low percentage and low intersample variance at 1 DPI. Of note, differences between 135 and ISU-1 infection were not statistically significant. Lack of precise concordance between scRNA-seq and flow cytometry Treg measurements likely reflects different definition parameters.¹⁸ Nonetheless, our data indicate a greater proportion of CD4 T cells with a regulatory phenotype at later stages of ISU infection.

Differential abundance by scRNA-seq in Treg cluster 4b was mirrored by further differential abundances of the other clusters within the sorted CD4 T-cell population (Fig. S6B). Cluster 4a was also less abundant in 135-infected pigs (Fig. 5C, top), while cluster 8 was more abundant (Fig. 5C, bottom, Fig. S6B). While we cannot deduce whether these differences between PRCV strains result from expansion or reduction of the investigated clusters, it is of interest to look at key transcripts from clusters 4a and 8 to understand the subpopulations they represent. Cells in cluster 8 expressed higher levels of FCG3A, encoding CD16, but also IFITM1 or CD40LG, suggesting that cells in this cluster are capable of a variety of immune functions. Cells in cluster 8 expressed higher levels of RPS genes (eg RPS7, RPS8, RPS15A, RPS20, and RPS27A), indicating high protein synthesis activity. Other highly expressed transcripts of cells in cluster 8, such as TMSB4X, TMSB10, RAC2, MYL6, S100A4, and CORO1A, are involved in cell migration, suggesting recent influx of these cells into BAL. In contrast, cluster 4a was characterized by transcripts involved in T-cell activation, including ICOS, CD28, RIPO2R, FYN, FYB1, CD247, and IL18R1. Together, these results indicate that at 20 DPI, the BAL CD4 T-cell population in 135-infected animals was dominated by a translationally active, migratory population, while that in ISU-1-infected animals had a greater abundance of regulatory and active helper T-cell populations.

Overall, these data indicate that PRCV infection with strains of differing virulence changes the balance of macrophage and CD4 T-cell subsets at the late stages of infection 20 DPI. The higher proportion of Tregs in ISU-1-infected pigs suggest a role in modulating the immune response, potentially contributing to milder lung pathology. In contrast, 135-infected animals showed greater abundances of T cells with high protein synthesis and involvement in cell migration, including inflammatory exudation in the bronchiole and alveoli.

Differential expression analysis

To investigate differential gene expression between the 135 and ISU-1 infections at 1 and 20 DPI, we used edgeR to conduct a pseudobulk differential expression analysis within each cluster. The number of DE genes per cluster are summarized in Fig. S7.

To determine whether DE genes were unique or shared between cell types (macrophages, CD4 T cells, CD8 T cells, B cells), significantly DE genes with a false discovery rate below 0.05 and absolute log fold change higher than 0.6 were visualized using Venn diagrams for 1 and 20 DPI (Figs. 6A and 7A, respectively). At 1 DPI, only a small number of genes were differentially expressed between 135- and ISU-1-infected pigs (Fig. 6A, B, Tables S4, S5). Nearly all of them were lower in 135-infected pigs, and most of them were detected in B cells or plasma cells.

Two genes were expressed at a significantly lower level in 135-infected animals in 2 cell types: NEU1 in CD4 T cells and B cells, and LGALS3 in CD8 T cells and B cells. NEU1 has the capacity to influence cell adhesion, which might alter cell–cell interactions of CD4 T cells and B cells in the BAL of 135-infected pigs.³⁰ Galectin 3, encoded by LGALS3, has been shown to negatively regulate T-cell activation;³¹ lower expression of LGALS3 in CD8 T cells and B cells in 135-infected pigs may reflect increased activation at 1 DPI.³¹ Further DE genes that were lower in 135-infected pigs included DUSP1 and ZFP36, which play a role in negatively regulating inflammatory responses. DUSP1 dephosphorylates MAP kinases, while ZFP36 is a zinc finger protein that inhibits proinflammatory mRNA translation.^32,33 Both DUSP1 and ZFP36 were also found to be downregulated in severe COVID-19, suggesting that reduced expression of these genes in 135-infected B cells and T cells may contribute to higher inflammatory responses (Fig. S10).^34,35

In addition to these shared DE genes, there were several DE genes unique for B cells, the majority of which were lower in the 135-infected pigs. Of note, most of these genes are derived from plasma cell clusters 4 and 15 (Fig. 6B). Similar to NEU1, some of them are involved in cell adhesion and migration: TMSB10 and TMSB4X regulate actin polymerization, while the tetraspanin CD9 is involved in the organization of cell membranes. Additionally, DE genes that were lower in 135-infected pigs included FOS and TYROBP both of which are involved in immune cell signal transduction.^36–38 Of note, porcine BAL preparations contain high proportions of plasma cells within total B cells (Fig. 3D).²¹ Their contribution to antibody production is currently unclear, and some of these cells could also have regulatory functions.³⁹

Given the lower expression of genes related to cell adhesion and cell signalling, this could indicate reduced regulatory capacity of plasma cells in 135-infected pigs. Such postulated lower regulatory function may be linked to another experimental finding. We analyzed early cytokine responses in BAL cell preparations using Luminex assays to measure IFN-α, IL-6, IL-12, IL-4, IL-1β, IL-8, TNF, and IL-10 at 1 and 5 DPI after ex vivo stimulation with virus homologous to the infection strain. BAL cells from 135-infected animals produced significantly higher concentrations of IFN-α, IL-6, and IL-12 in comparison with BAL cells from ISU-1-infected animals at 1 DPI (Fig. 6C and Fig. S9), indicating stronger inflammatory responses in 135-infected pigs. This also aligns with higher expression of the antiviral MX1 gene in Tregs (cluster 3d, Fig. S8) and the plasma cell cluster 1 (Fig. 6B).

In contrast to 1 DPI, at 20 DPI the number of genes differentially expressed between ISU-1- and 135-infected pigs was much higher (634 in total compared to 23 at 1 DPI), and comprised all cell types, with several DE genes shared between cell types (Fig. 7A, Table S6) and cell clusters (Fig. S10). Overall, these DE genes showed a very heterogeneous picture. For example, among the genes with lower expression in CD4 and CD8 T cells from 135-infected pigs were CD69 and CD28, involved in T-cell activation and costimulation, respectively. Transcripts of the anti-apoptotic BIRC3 gene were also lower, suggesting a reduced capacity for stimulation or survival, but transcripts of the pro-apoptotic FAS gene were also lower. The dichotomy of simultaneously lower BIRC3 and FAS expression was also found in B cells and macrophages. Similarly, genes like ERBB4 and LRRK2, which are involved in cell activation,⁴⁰ were both more highly expressed in CD4 and CD8 T cells of 135-infected pigs and contrasted the lower expression of genes encoding CD69 and CD28. The anti-apoptotic BCL2 gene was also expressed to a greater level in CD4 T cells of 135-infected pigs. This suggests that the 2 PRCV strains drive separate immune activation pathways, which are, however, often similar in their overall function.

We also performed GO enrichment to investigate broader biological functions of the DE genes at 20 DPI (Tables S7 and S8). Several GO terms were shared between B-cell clusters (including plasma cells), CD4 and CD8 T cells (Fig. 7B, Fig. S11). Nearly all of those terms represented high-level biological functions like “multicellular organism development,” “transcription by RNA polymerase II,” “cellular response to chemical stimulus,” “regulation of molecular function,” or “cell differentiation.” In a complementary analysis, we performed KEGG pathway enrichment analysis. To focus on immunological pathways, we post hoc excluded pathways related to cancer, bacterial, and parasitic diseases (Fig. 7C). One of the most prominent findings was enrichment of apoptosis-related pathways across all plasma cell and CD4 T-cell clusters, as well as in macrophage cluster 20. Another pathway broadly enriched was TNF signalling, which may be linked to apoptosis. Some plasma cell and CD4 T-cell clusters also showed differences in MAPK and NFκB pathways, as well as in the parathyroid hormone signalling pathway, the latter potentially associated with calcium regulation and active immunoglobulin production. Therefore, both GO and KEGG analysis indicate that even 20 days after infection, when PRCV has already been cleared for 10 days (Fig. 1) the different viral strains continue to differentially impact gene regulation and cell metabolism in immune cells in the BAL. Of note, as on 1 DPI, B cells and the various plasma cell clusters showed the highest number of DE genes compared to other cell types (Fig. 7A, Fig. S7B), indicating a hitherto unrecognized impact of viral strain pathogenicity on gene expression in local B cells and plasma cells.

In summary, these data suggest that the strong inflammatory processes triggered as early as 1 DPI have lasting consequences on the transcriptional landscape of immune cells in BAL, even after the virus has been cleared. Our data also indicate that differences in gene expression between pigs infected with the ISU-1 or 135 PRCV strains are most prominent in B cells during the early stage of infection. The absence of DE genes in macrophages suggests that this cell type responded similarly to both infections. Additionally, higher IL-6 and IL-12 cytokines in 135-infected animals may be associated with greater pulmonary lesions in these animals.

Discussion

The spectrum of outcomes following SARS-CoV-2 exposure in humans is wide and considered to be largely determined by the host response. Most infections are asymptomatic or mild, indicating that effective host responses can successfully dampen viral replication and symptoms. However, the factors that precipitate severe disease are not fully understood. There is therefore a need for models that closely resemble human SARS-CoV-2. Currently mouse, ferret, hamster, and nonhuman primate models have been developed to study SARS-CoV-2, but these are not natural hosts. In contrast pigs are natural hosts for PRCV, but little is known about the mechanism of pathogenicity and protective immunity in this species. Here we used infection with the weakly and strongly pathogenic PRCV strains ISU-1 and 135, respectively, for in-depth mechanistic evaluation of the pathogenesis and immune control of respiratory CoVs.

Pathological examination revealed prominent respiratory lesions in 135-infected animals at 5 and 11 DPI which were still present at 20 DPI. This was characterized by a marked upper and lower respiratory pathology and with an associated high viral shedding in nasal swabs and viral load in trachea, lung, and BAL at 5 DPI. Further, infectious virus was detected in the eyelids of a large proportion of 135-infected animals at 5 DPI, which suggest a similar pathogenesis as that reported in humans with SARS-CoV-2 infection.^28,41,42 While the virus tropism was not examined by microscopy, the ability of 135 to replicate in eyelid tissue, likely conjunctiva epithelium, may represent an understudied site of viral replication and potential route of both infection and transmission for respiratory CoVs.

As well as being more virulent and pathogenic, 135 infection induced greater spike-specific IgG and IgA titers in serum, BAL, and nasal swabs compared to ISU-1 infection. The increased antibody response was associated with a corresponding increase in plasmablast/plasma cells in the BAL. Additionally, we identified, for the first time, cells with a Tfh phenotype that likely play a crucial role in maintaining local B cell and CD8 T cell memory.²⁹

Infection with 135 induced a higher proportion of virus-specific IFN-γ and TNF secreting CD4 and CD8 T cells. The frequencies of cytokine-secreting CD4 and CD8 T cells were greater in the BAL compared to TBLN and spleen in both groups. The highest frequencies of IFN-γ and TNF-secreting cells and the highest antibody titers were detected at 11 and 20 DPI, while the peak of lung pathology and viral load was at 5 DPI, implying that, along with the greater ability of 135 virus to infect more cells, the infection also resulted in a dysregulated and/or more vigorous innate immune response. Consequently, the presence of more virus antigen and this robust innate response may further drive a more vigorous adaptive immune response that may additionally contribute to the lung damage. This could be consistent with findings in SARS-CoV-2 patients where higher viral load was associated with increased disease severity along with robust activated CD4 T cells, activated CD8 TEMRA, hyperactivated or exhausted CD8 T cells and plasmablasts.^43–45

We therefore examined the early (1 DPI) and late (20 DPI) events in the BAL of animals infected with 135 and ISU-1 to identify gene expression changes that may lead to more severe or less severe inflammatory and pulmonary pathology. We analyzed BAL as this is an easily accessible site for sampling in both humans and pigs, and we have previously established the feasibility of performing single-cell transcriptomics with these samples in pigs.¹⁸ Single-cell RNA-seq analysis of BAL cells revealed that the more pathogenic 135 strain shifts the balance of CD4 T cells away from Tregs (cluster 4b) toward translationally active CD4 T cells (cluster 8) at 20 DPI. The higher proportion of biosynthetic CD4 T cells in the 135-infected animals may have contributed to the increased antibody production, while the more abundant Tregs in ISU-1-infected pigs were associated with lower pulmonary lesions.

Differential expression analysis also indicated variation in gene expression between 135- and ISU-1-infected pigs at 1 DPI among genes related to cell migration and immune system processes. In 135-infected animals, higher expression of the antiviral gene MX1 in plasma cells, along with increased cytokine production of IL-6 and IL-12 in PRCV-stimulated BAL cultures, indicates a heightened inflammatory response, which might contribute to more severe pulmonary pathology. The higher concentration of IFN-α in 135-infected pigs may be attributed to higher viral load but is unlikely to be associated with the increased pulmonary pathology.⁴⁶ However, elevated IL-6 production has been associated with severe disease in COVID-19 patients, while IL-12 promotes Th1 differentiation and NK cell activation.^47–49

Most of the DE genes at 1 DPI were found in B cells, primarily in plasma cells. While the total number of DE genes was much higher at 20 DPI, with a high number in B and plasma cells, many were also shared between B cells and CD4 and CD8 T cells. A notable finding was that most of those DE genes were lower in cells from 135-infected pigs compared to ISU-1-infected pigs with some identified as negative regulators of inflammatory responses.^32–35 This may also indicate a hitherto unnoticed role of B cells and plasma cells in local immune homeostasis in the lung. Regulatory functions of B cells and plasma cells have been described in a wide range of immune responses in mice and humans.^39,50 Since B cells with regulatory functions do not form a defined cell lineage, no universal phenotype or master transcription factor is available for their identification. However, CD9, which was among the genes with lower expression in B cells and plasma cells of 135-infected pigs, has been described as a molecule associated with IL-10-competent regulatory B cells in mice.⁵¹

GO analyses indicated that at 20 DPI genes defferentially expressed between 135- and ISU-1-infected pigs are associated with high level cellular functions, such as gene transcription and cell metabolism. This applied to B cells (including plasma cells) but also CD4 and CD8 T cells. Similarly, KEGG pathway investigations showed differential apoptosis, TNF, and innate signalling in most of those cell types. This suggests that, despite clearance of virus and a reduction in pathology at this time point, longer lasting effects caused by the 2 different PRCV strains are maintained. This may have significant implications for immune responses to subsequent viral or bacterial infections in the respiratory tract, commonly seen with “porcine respiratory disease complex.”^52,53 This phenomenon may arise not only from the pathogens involved, but also by the underlying priming of the immune system by previous infections, despite successful clearance. This is reminiscent of post-viral syndrome reported in humans, such as the long-term damage in COVID-19 patients, commonly known as “long COVID.”⁵⁴ The pig model may offer an opportunity to experimentally explore these effects by following PRCV infection with another infection, providing a platform to study post-viral syndrome changes.

A limitation of this study is the lack of uninfected control pigs, which would have indicated whether the observed changes in gene expression represented up- or downregulation. Therefore, we were careful to refer to higher and lower expression in one viral strain versus the other. Nevertheless, the data suggest that the more pathogenic 135 strain leads to increased virus replication and triggers stronger, potentially damaging innate and adaptive immune responses, characterized by heightened inflammation, leading to more prominent respiratory pathology. In contrast, the less pathogenic ISU-1 strain induces a more regulated immune response, with more Tregs and controlled cytokine production, resulting in milder pathology. These results agree with similar studies in humans highlighting the importance of timely and appropriately scaled immune responses to SARS-CoV-2 and that early intervention with therapies that dampen immune responses, such as steroids, has proved crucial in preventing progression to severe disease.^55–58

Our data provide a reference map for changes in gene expression early in PRCV infection. Novel vaccine and therapeutics could be tested in the pig PRCV model, which more accurately reflects pathological inflammation in the respiratory tract of humans. Our model will also contribute to determining how to achieve robust virus-specific immunity to clear virus infection without causing exacerbated inflammation and immunopathology, which is crucial for limiting life-threatening disease.

Acknowledgments

We thank the flow cytometry, sequencing, and bioinformatics facilities at The Pirbright Institute for enabling this research. We are grateful to the staff at Animal and Plant Health Agency for providing excellent animal care. C.R. is supported by the Equal Opportunities Foundation (Hong Kong) and the Braithwaite Family Foundation. The visual abstract was created using BioRender.com.

Author contributions

Ehsan Sedaghat-Rostami (Data curation [Equal], Formal analysis [Equal], Investigation [Equal], Methodology [Equal], Validation [Equal], Visualization [Equal], Writing—original draft [Equal], Writing—review & editing [Equal]), Brigid Veronica Carr (Formal analysis [Equal], Investigation [Equal], Visualization [Equal], Writing—review & editing [Supporting]), Liy Yang (Data curation [Equal], Formal analysis [Equal], Investigation [Equal], Methodology [Equal], Software [Equal], Validation [Equal], Visualization [Equal], Writing—original draft [Equal], Writing—review & editing [Equal]), Sarah Keep (Conceptualization [Supporting], Formal analysis [Equal], Funding acquisition [Supporting], Methodology [Equal], Visualization [Supporting], Writing—review & editing [Supporting]), Fabian Z.X. Lean (Formal analysis [Equal], Investigation [Equal], Methodology [Equal], Visualization [Supporting], Writing—review & editing [Supporting]), Isabella Atkinson (Formal analysis [Supporting]), Albert Fones (Formal analysis [Supporting]), Basudev Paudyal (Formal analysis [Supporting], Investigation [Supporting], Visualization [Supporting], Writing—review & editing [Supporting]), James Krik (Formal analysis [Supporting]), Eleni Vatzia (Formal analysis [Supporting]), Simon Gubbins (Formal analysis [Supporting], Methodology [Equal], Validation [Supporting], Writing—review & editing [Supporting]), Erica Bickerton (Conceptualization [Supporting], Funding acquisition [Supporting], Writing—review & editing [Supporting]), Emily Briggs (Formal analysis [Supporting], Investigation [Supporting], Writing—review & editing [Supporting]), Alejandro Nunez (Investigation [Supporting], Writing—review & editing [Supporting]), Adam MacNee (Formal analysis [Supporting], Investigation [Supporting]), Katy Moffat (Formal analysis [Supporting], Investigation [Supporting], Methodology [Supporting], Writing—review & editing [Supporting]), Graham Freimanis (Data curation [Supporting], Formal analysis [Supporting], Investigation [Supporting], Methodology [Supporting], Writing—review & editing [Supporting]), Christine Rollier (Writing—review & editing [Supporting]), Andrew Muir (Formal analysis [Supporting], Investigation [Supporting], Methodology [Supporting], Software [Equal], Validation [Equal], Visualization [Supporting], Writing—review & editing [Supporting]), Arianne C. Richard (Formal analysis [Supporting], Investigation [Supporting], Resources [Equal], Software [Lead], Supervision [Equal], Validation [Supporting], Visualization [Supporting], Writing—review & editing [Supporting]), Nicos Angelopoulos (Data curation [Equal], Formal analysis [Supporting], Investigation [Supporting], Resources [Equal], Software [Equal], Supervision [Equal], Validation [Equal], Visualization [Supporting], Writing—review & editing [Supporting]), Wilhelm Gerner (Conceptualization [Supporting], Investigation [Supporting], Methodology [Equal], Supervision [Equal], Validation [Equal], Visualization [Equal], Writing—original draft [Equal], Writing—review & editing [Equal]), and Elma Tchilian (Conceptualization [Lead], Data curation [Lead], Formal analysis [Equal], Funding acquisition [Lead], Investigation [Lead], Methodology [Lead], Project administration [Lead], Resources [Lead], Supervision [Lead], Validation [Equal], Visualization [Equal], Writing—original draft [Lead], Writing—review & editing [Lead])

Supplementary material

Supplementary material is available at The Journal of Immunology online.

Funding

This work was supported by UKRI Biotechnology and Biological Sciences Research Council (BBSRC) grant BB/X014266/1 and Institute Strategic Programme and Core Capability Grants to The Pirbright Institute (BBS/E/PI/230001A, BBS/E/PI/230001B, BBS/E/PI/230001C and BBS/E/PI/23NB0003), BBS/E/B/000C0545, MR/W016303/1.

Conflicts of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Data availability

All R scripts used for processing of the scRNA-seq data are available at: https://github.com/LynnGoodnight/bsp185. Fastq and h5 files from sequencing are available at: GEO (GSE288145) https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE288145.

References

Author notes