Early pathogenesis profiles across SARS-CoV-2 variants in K18-hACE2 mice revealed differential triggers of lung damages

Abstract

The on-going COVID-19 pandemic has given rise to SARS-CoV-2 clades and variants with differing levels of symptoms and severity. To this end, we aim to systematically elucidate the changes in the pathogenesis as SARS-CoV-2 evolved from ancestral to the recent Omicron VOC, on their mechanisms (e.g. cytokine storm) resulting in tissue damage, using the established K18-hACE2 murine model. We reported that among the SARS-CoV-2 viruses tested, infection profiles were initially similar between viruses from early clades but started to differ greatly starting from VOC Delta, where the trend continues in Omicron. VOCs Delta and Omicron both accumulated a significant number of mutations, and when compared to VOCs Alpha, Beta, and earlier predecessors, showed reduced neurotropism and less apparent gene expression in cytokine storm associated pathways. They were shown to leverage on other pathways to cause tissue damage (or lack of in the case of Omicron). Our study highlighted the importance of elucidating the response profiles of individual SARS-CoV-2 iterations, as their propensity of severe infection via pathways like cytokine storm changes as more variant evolves. This will then affect the overall threat assessment of each variant as well as the use of immunomodulatory treatments as management of severe infections of each variant.

Introduction

The emergence of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) has resulted in the coronavirus disease 2019 (COVID-19) pandemic that caused major healthcare and economic burden worldwide (1). The infection SARS-CoV-2 caused is highly heterogeneous ranging from asymptomatic and mild to severe and even death (2). In addition, the heterogeneity is further enhanced by the rapid mutation of the virus giving rise to different clades and variants that caused varying pathologies over the course of the pandemic (3). Hence, it is imperative to characterize the differential pathology of SARS-CoV-2 clades and variants and their mechanisms to improve management of the infection, especially for development of therapies that are based on the modulation of these mechanisms.

Coronaviruses have been shown to recombine extensively due to their proofreading mechanism that can facilitate homologous recombination events, giving rise to novel variants of varying pathology and severity (4). To date, SARS-CoV-2 is clustered into 9 major clades based on the Global Initiative on Sharing All Influenza Data (GISAID) database (5, 6). Clade L is the first novel viral genome sequence serving as the reference ancestral strain (7), and other representative clades with major genetic mutations include clades S, V, G, GK, GH, GR, GV, and GRY. Among the clades, clade G and its derivatives have rapidly become the dominant globally circulating SARS-CoV-2 that gave rise to the major variants of concern (VOCs) (8), including Alpha (B.1.1.7), Beta (B.1.351), Gamma (P.1), Delta (B.1.617.2) and Omicron (B.1.1.529). The VOCs have caused concerns as they were shown to have increased virulence, greater transmissibility, resistance to antibody neutralization, and evasion to host- and vaccine-mediated immunity (9–11). It is understood that the SARS-CoV-2, especially the VOCs, resulted in different degrees of severity and symptoms in humans and animal models (12, 13). However, due to the rapid evolution of each VOCs and them rarely overlapping in emergence, there have yet to be a systematic comparison of the early triggers of airway and lung damage among the major SARS-CoV-2 variants, especially between all VOCs.

Several studies have shown that the expression of pro-inflammatory cytokines is associated with disease progression in COVID-19 patients (14–17). Hypercytokinemia, also known as cytokine storm, is one of the major triggers for tissue damage, multiple organ dysfunction, and lethal outcomes in SARS-CoV-2 infection (18, 19). Further, it is well documented that the build-up of activated neutrophils and macrophages in the airways and lungs following SARS-CoV-2 infection results in increased pro-inflammatory cytokines that contribute to airway pathology (20). On the other hand, the evolution of SARS-CoV-2 VOCs over the course of the pandemic showed great variation in their severity and pathology, which is especially apparent in Omicron VOC, suggesting differential cytokine activation that contributes to their pathogenesis (21). Henceforth, we employed the established K18-hACE2 mouse model to systematically compare the disease phenotypes across clades and variants and correlate them with their respective early response cytokine profiles, to determine involvement of cytokine storm and selected pathways’ contribution to their pathogenesis. We aim to elucidate the evolutionary progression of SARS-CoV-2 pathogenesis to understand the pathology spectrum of the virus, which will help in current and future management and drug discovery effort against the virus.

Materials and methods

Cell line and SARS-CoV-2 variants

African green monkey kidney cells, (Vero E6; ATCC CRL-1586™) were cultured in Dulbecco’s Modified Eagle’s Medium (DMEM; Cytiva) supplemented with 10% heat-inactivated fetal calf serum (FCS).

SARS-CoV-2 viruses were isolated from nasopharyngeal swabs of qRT-PCR confirmed COVID-19 patients and propagated in Vero E6 for infection experiments, where low passages of not more than 4 passages of the virus were used. Genome sequences from the swab samples uploaded to GISAID by the National Public Health Laboratory, National Centre for Infectious Diseases, Singapore, were used to confirm the variants’ identity. A total of eight viruses were used: 1. Clade L (ancestral), Lineage B (EPI_ISL_574502); 2. Clade V, Lineage B.29 (EPI_ISL_493419); 3. Clade G, Lineage B.1.610 (EPI_ISL_443240); 4. Clade GR, Lineage B.1.1 (EPI_ISL_443199); 5. VOC Alpha variant, Clade GRY (thereafter GRY-α), Lineage B.1.1.7, (EPI_ISL_754083); 6. VOC Beta variant, Clade GH/501Y.V2 (thereafter as GH-β), Lineage B.1.351.3 (EPI_ISL_1173248); 7. VOC Delta variant, Clade GK (thereafter GK-δ), Lineage B.1.617.2 (EPI_ISL_2621925); 8. VOC Omicron variant, Clade GRA (thereafter GRA-o), Lineage B.1.1.529 (EPI_ISL_7195620). All virus experiments were performed in NUS Medicine biosafety level 3 (BSL-3) core facility and all protocols were approved by the BSL-3 biosafety committee (BBC) and the institutional biosafety committee (IBC) of the National University of Singapore (NUS).

Experimental infection

Eight to nine weeks-old female K18-hACE2 transgenic mice (InVivos Ptd Ltd, Singapore) were acclimatized in the ABSL-3 facility for 72 hours prior to infection and were infected with approximately 1 × 10³ PFU of SARS-CoV-2 virus suspension in PBS, via the intranasal route. Baseline body weight were measured prior to virus infection. Body weight, physiological conditions, and survival were monitored daily by two personnel for the duration of the experiment (14-days post infection, dpi) or until the humane endpoint was reached. Each infection group was performed at n=6, except variant GRY-α at n=5, and mock infection at n=3. Scoring of the infected mice physiological conditions were based on five criteria: appearance of mouse coat, level of consciousness, activity level, eye condition, and respiratory quality. Conditions were scored on a scale from 1 to 5, using an observation system adapted from Shrum, et al. (22), with 1 denoting normal physiological state and 5 being the most severe. Additional groups of mice at n=7, except GK-δ at n=5, and GRA-o at n=6 were sacrificed on 4 dpi to assess the viral load in the brain, lung, liver, and spleen; and histological analyses of the lung. Tissues were halved and homogenized in DMEM supplemented with antibiotic-antimycotic and the supernatant was used for plaque assay, while the remaining tissue pellets were kept in -80°C for RNA isolation. The other half of the tissue was fixed with 3.7% formaldehyde (10% formalin) solution for at least 19 h before removal from BSL-3 containment for histological analysis. All animal experiments were conducted in NUS Medicine Biosafety Level 3 (ABSL-3) facility in accordance with NUS IACUC protocol no. R20-0504, using NUS IBC and BBC approved SOPs.

Plaque assay

For virus titre determination, viral supernatants from homogenized tissues were serially diluted in 10-fold increments in DMEM supplemented with antibiotic-antimycotic (Gibco). 250 µl of each serially diluted supernatant were added to confluent Vero E6 cells and incubated for 1 h at 37°C with rocking at 15 min intervals. After 1 h of adsorption, virus inoculum was removed, and cells were washed once with phosphate-buffered saline (PBS). Overlay media containing 1.2% microcrystalline cellulose-DMEM supplemented with antibiotic-antimycotic were added to each well and incubated for 3 days at 37°C at 5% CO₂ for plaque formation. Cells were fixed in 10% formalin overnight and stained with crystal violet for plaque visualization. Number of plaques were determined, and virus titre of individual samples were expressed in logarithm of plaque forming units (PFU) per organ.

Histopathological analyses

Formaldehyde fixed tissues were routinely processed, embedded in paraffin blocks (Leica Surgipath Paraplast), sectioned at 5 µm thickness, and stained with hematoxylin and eosin (H&E; Thermo Scientific) following standard histological procedures. A multiparametric, semiquantitative scoring system was further used to assess the magnitude of histo-morphological and -pathological changes in lung tissues based on six criteria: inflammatory cell infiltrates, hemorrhage, edema, degeneration of alveolar epithelial cells, parenchymal wall expansion and bronchiole epithelial cell damage (23, 24). For the histopathological parameters, a score of 0 – 3 was ordinally assigned, where 0 indicated normal; 1 indicated less than 10%; 2 indicated 10 – 50%; and 3 indicated more than 50% of lung regions affected. The average of histopathological scores of the mice from each group was taken as the final evaluation index.

Immunohistochemistry (IHC)

The paraffin embedded lung sections were deparaffinized and rehydrated, followed by heat-mediated antigen retrieval process. The sections were then treated with hydrogen peroxide blocking and protein blocking for 10 min, respectively. After blocking the non-specific binding sites, sections were incubated with anti-SARS-CoV-2 nucleocapsid protein monoclonal antibody (1:1000; Abcam), followed by the respective horseradish peroxidase (HRP)-conjugated secondary antibody. The sections were subsequently visualized using DAB solution and counterstained with hematoxylin.

RNA isolation and qRT-PCR

Homogenized lung tissue pellets were subjected to total RNA extraction using QIAGEN RNeasy^® mini kit in accordance with manufacturer’s guidelines. RNA purity and concentration were determined using Tecan Infinite^® M200 Pro coupled with NanoQuant Plate™. Real-time reverse transcription polymerase chain reaction (qRT-PCR) was performed using QIAGEN RT² First Strand Kit, QIAGEN RT² SYBR^® Green Mastermix, and QIAGEN RT² Profiler PCR Array (PAMM-150ZE-4) in accordance with manufacturer’s recommendations.

Gene expression analysis

Murine cytokines mRNA expression levels from cDNA conversion were analyzed using QIAGEN GeneGlobe Design and Analysis hub, and gene expression were expressed as Log₂ fold regulations. Expression levels were calculated relative to Gapdh and normalized to the mock infected mice using the ΔΔC_T method with the QIAGEN GeneGlobe software. The genes from the QIAGEN panel were further sub-divided into selected pathways that may potentially contribute to tissue damage using DAVID functional annotation tool (25, 26). Pathways were selected from Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases for visualization of enrichment in each pathway across variants.

Statistical analysis

Analyses were calculated using two-tailed Mann-Whitney U test to evaluate the data obtained using Graphpad Prism version 9.0 with * denoting that p < 0.05, ** denoting that p < 0.01, and *** denoting p < 0.001. For gene expression analysis, p-values were obtained from the analysis with QIAGEN GeneGlobe Design and analysis hub. Only genes fulfilling both criteria of >2 fold-regulation compared to uninfected controls and having a p-value of <0.05 were considered significant.

Results

To characterize the pathogenesis of various SARS-CoV-2 viruses from different clades and VOCs, and their respective mechanisms leading to airway and lung damage, we tracked a total of eight major SARS-CoV-2 iterations – L (ancestral), V, G, GR, GRY-α, GH-β, GK-δ, and GRA-o (Figure 1). The repertoire of our coverage includes the major SARS-CoV-2 clades and variants of the pandemic to date. The mutations across the iterations based on their clinical isolates’ sequences obtained from GISAID were listed in Table 1.

Figure 1

Changes in survival, symptoms, and tissue tropism over the course of SARS-CoV-2 evolutionary iterations

We first compared the survival of K18-hACE2 mice following infection with the different SARS-CoV-2 viruses. It was observed that, with the exception of V (34% survival) and GRA-o (83.3% survival), 0% of infected mice survived the infection, where they succumbed between 6 to 8 dpi, with GH-β and GK-δ hitting the end point the earliest at 6 dpi (Figure 2A). When comparing the Kaplan-Meier survival curves, it was shown that only GRA-o, the latest evolutionary iteration, was significantly different in mice survival compared with all other infected groups except for V, another iteration of the virus where there were mice that survived (Figure 2B). The weight loss trends largely followed the survival curves in which rapid weight loss were observed between 4 and 7 dpi for all groups other than V and GRA-o, for which we observed a rebound in weight from surviving mice (Figure 2C). We then compared clinical symptom scores (averaged from appearance, eye condition, consciousness, activity, and respiration) and observed that L and GRY-α showed the most symptom presentation, followed by the others that showed moderate symptoms, while GRA-o was shown to be not presenting any clinical symptom throughout the duration of study (Figure 2D, Figures S1A–E). We then examined viral titre retrieved from the tissues (lung, brain, liver, and spleen) of SARS-CoV-2 infected mice. No infectious viruses were recovered from liver and spleen homogenates using plaque assay (data not shown). For viral titre in the lungs, high titres of SARS-CoV-2 was retrieved from the lungs of L, V, G, GR and GRY-α infected mice, at geometric mean titres of 4.27 × 10⁵, 4.55 × 10⁴, 2.75 × 10⁵, 1.09 × 10⁶, and 3.12 × 10⁵ PFU/organ, respectively. Titres from V are slightly lower, while GR slightly higher when compared with L infected group (p < 0.05). On the other hand, significantly lower titre was observed in GH-β, GK-δ and GRA-o when compared with L, at 6.86 × 10², 7.2 × 10¹ and 1.1 × 10¹ PFU/organ, respectively (p < 0.05) (Figure 2E). As for viral titre in the brain, comparatively high titre was found in the brain tissue of L, GRY-α and GH-β infected groups, at 1.93 × 10⁵, 6.02 × 10³ and 3.4 × 10³ PFU/organ, respectively; while significantly lower titre in the brain was retrieved in brain tissue of G, GR and GK-δ at 1.03 × 10², 4.4 × 10¹, and 8.7 × 10¹ PFU/organ, respectively. No virus was retrieved from the brain of GRA-o infected mice (Figure 2F). Notably, unlike in the lung, viral titre in the brain was more variable where about half of the mice did not have viral titre retrieved from the brain in most groups, regardless of the titre observed.

Figure 2

Early histopathological damage of SARS-CoV-2 viruses were not clearly distinguishable

From the lung tissue collected at 4 dpi, it was observed that histopathological damage of SARS-CoV-2 variants was largely similar, with histopathological index score ranging between 0.63 to 1.69. Tissues from GR and GK-δ infected mice were on the higher end of the scores while those from L and GRA-o were with lower scores (Figure 3A). The factors that likely contributed to the higher histopathological index were inflammatory cell infiltrates, where the trend was consistent with the overall index compared with other factors (Figure 3B;Figures S2A–F). On the other hand, no histological abnormality was observed in the bronchiole epithelial cells following all variants’ inoculation, and neither bronchitis nor bronchiole epithelial denudation was present in the lung tissues (data not shown). Next, we correlated the extent of infected cells and viral presence in the lungs of the ancestral L infected group, with the VOC infected groups. Using IHC of SARS-CoV-2 nucleocapsid protein, we observed that the viral distribution in the lungs at 4 dpi varied greatly (Figures 3C, D). GK-δ, notably, showed major nucleocapsid protein presence in the lung (93.2% coverage), contrary to the low infectious titre detected from the infected murine lungs while GRA-o had almost no nucleocapsid protein detected in the lung tissues (5.83% coverage). Our observation suggested that there was no correlation between nucleocapsid coverage, which are produced in excess during coronavirus infection (27, 28), in the lungs and the lung histopathological damage index between variants.

Figure 3

SARS-CoV-2 cytokine profiles between evolutionary iterations revealed differences in lung damage triggers

To further elucidate potential mechanisms in each of the SARS-CoV-2 viruses that led to their respective phenotypes, we performed a qRT-PCR array on 84 genes related to inflammatory responses from lung tissues of SARS-CoV-2 infected K18-hACE2 mice. After taking into account fold-regulation and statistical significance (>2 fold-regulation and p< 0.05), the cytokine profiles between the infected groups differed, some significantly, from each other. We noted the groups with the larger significant changes from the ancestral L virus (18 up; 0 down) were G (23 up; 6 down), GR (21 up; 2 down), GRY-α (26 up; 2 down), GH-β (21 up; 5 down) and GRA-o (26 up; 1 down). On the other hand, V (4 up; 0 down) and GK-δ (6 up; 4 down) had fewer significant changes compared with L (Figure 4). The exact changes in fold-regulation of the genes tested were shown in Table S1. In general, the host responses of the SARS-CoV-2 infection of all groups tested deviated significantly from the ancestral L virus, with some having unique significant responses as observed in G (Lta, Tnf, Pf4), GH-β (Il24, Il1a, Mstn, Il22), GK-δ (Il7, Tnfsf11b, Il18), and GRA-o (Ccl22, Bmp4, Tnfsf10, Mif, Hc, Bmp7, Gpi1, Cx3cl1). In view of the differences present, particularly in GRA-o which has the most of them, further dissection on how the differential cytokine changes led to the pathogenesis of current and future variants is thus warranted.

Figure 4

As cytokine storm was perceived to be a major trigger of lung damage in COVID-19 (18, 19), we grouped the panel of genes we tested into pathways associated with cytokine storm using DAVID. We selected pathways involved directly in cytokine storm (19), and associated pathways like MAPK cascade, TNF responses, neutrophil chemotaxis, IL-17 regulation, negative regulation of inflammatory response, and T-cell mediated cytotoxicity to determine how much they played a role in tissue damage for each virus iteration (Figure 5). Firstly, in the cytokine storm pathway, we observed that for iterations that come after L and V, more genes such as Ifng, Il12b and Il10 were found to be significantly differentially expressed in G, GR, GRY-α and GH-β, before being reduced again in GK-δ and GRA-o (Figure 5A). Similar trends were observed in TNF response and neutrophil chemotaxis pathways in genes like Xcl1, Ccl1 and Ccl11 (Figures 5B, C). The MAPK cascade pathway was prominently differentially expressed in G, GR, GRY-α and GH-β, but in this case, L and GK-δ as well (Figure 5D). On the other hand, negative regulation of inflammatory responses early in infection were linked with cytokine storm (29), and we found that G, GR, GRY-α and GH-β further have differentially up-regulated expression of Il10 (Figure 5E). IL-17 pathway, associated with neutrophil infiltration with genes from the pathway like Il15 and Il12b were up-regulated in all groups except V and GK-δ (Figure 5F). Finally, T-cell mediated cytotoxicity pathway was found to be prominently changed for G, GR, GRY-α, GH-β, and GRA-o (Figure 5G). Interestingly, only GRA-o had a further significant up-regulation of Il12a. Overall, genes from cytokine storm and associated pathways (Ccl2, Ccl7 and Ccl12) were found to be differentially regulated in infections of L, G, GR, GRY-α, and GH-β while not as strongly involved in V and GK-δ. GRA-o, while having similar trend of cytokine storm associated pathway gene changes, had more subdued expression, and a slight difference in T-cell mediated cytotoxicity (Il12a).

Figure 5

Discussion

The use of K18-hACE2 murine model to understand SARS-CoV-2 infection has been instrumental to the understanding of SARS-CoV-2 pathogenesis mechanisms since the start of the pandemic, including comparisons of early clades’ differential severity (13, 30–34). Moreover, the K18-hACE2 mice were shown to be simulating severe SARS-CoV-2 infection useful in the identification of immune responses leading to severe outcome (35). Hence, here we reported results using the model to correlate early cytokine profiles (4 dpi) of early SARS-CoV-2 evolutionary iterations (L, V, G, GR), and the VOCs Alpha (GRY-α), Beta (GH-β), Delta (GK-δ) and Omicron (GRA-o), with their infection outcome. Interestingly, the majority of the pathogenesis mechanism and infection outcome remained largely similar up until the Beta VOC, before major changes occur in VOCs Delta and Omicron that emerged later. This may be partially due to the larger number of mutations accumulated (Table 1) for Delta and Omicron VOC, which potentially altered their biology significantly.

Our study demonstrated that early cytokine profiles, viral titre and histopathological damage worked in concert to influence and determine the infection outcome. Overall, the infection outcome remained largely similar from L to GH-β Beta VOC, and then deviated in VOCs that emerged late in the pandemic i.e. GK-δ Delta and GRA-o Omicron VOCs. Among the early clades and VOCs up to Beta, even the milder clade V virus with ~30% survival rate had similar rates of neurotropism based on viral titre from brain tissue, sufficiently high viral titre in the lungs, and cytokine profiles consistent with previous report (36). While clade V virus had fewer significantly differentially expressed genes, the major genes involved in cytokine storm; TNF pathway and neutrophil chemotaxis, Ccl2 and Ccl7 remained, suggesting similar, but milder presentation of pathology for its infection outcome. The subsequent clade and VOCs in the earlier part of the pandemic, G, GR, GRY-α and GH-β had further inflammatory suppression gene Il10 differentially up-regulated that potentially augmented their pathology slightly, as seen in the histopathological index. Furthermore, the viruses from earlier in the pandemic, up to VOC GH-β, showed higher clinical presentation scores, which may be linked to neurotropism, prominent in these groups of infections.

On the other hand, VOCs that emerged later in the pandemic, GK-δ Delta and GRA-o Omicron, followed a different pathogenesis mechanism that resulted in tissue damage. Firstly, neurotropism was observed to be less apparent in Delta, and completely absent in the Omicron variant in the K18-hACE2 mice. It was observed that mechanisms that may contribute to cytokine storm were also less apparent in GK-δ and GRA-o infection. In addition, immune suppression, a mechanism common in SARS-CoV-2 (37), was observed strongly in the GK-δ Delta variant, albeit different from Il10 mediated suppression observed in clades and VOCs that came before. The Delta variant in our study instead showed an almost universal suppression of host responses early post infection at 4dpi (with very little significant differential regulation). This suppression may potentially contribute to the increased nucleocapsid expression and inflammatory infiltrates in the lungs, which may suggest diffusely infected cells in the lungs. On the other hand, the GRA-o variant showed the most varied response among the variants where we observed the least tissue damage, low viral titre and lack of neurotropism. This observation is in line with the finding that Omicron variant has more preference to the upper airway than the lungs (38). In addition, the cytokine profile of our Omicron GRA-o variant infection showed a response largely different than that of its predecessors. It had generally weaker changes in expression of cytokine storm and related pathway genes, coupled with moderately stronger T-cell mediated cytotoxicity gene expression in Il12a and Il12b. The Omicron cytokine profile suggests a different pathogenesis that may contribute to the overall milder tissue damage congruent with other studies to date (21, 39). The effects of the major changes in Omicron cytokine profiles require further studies to validate their correlation with disease outcome, as well as potential interaction with risk factors and chronic conditions (e.g. asthma) due to their major differences and lower lethality.

The finding from our study is crucial in the understanding of SARS-CoV-2 pathogenesis mechanisms. It showed the evolutionary potential of SARS-CoV-2 where its ability to accumulate high number of mutations may result in highly variable pathogenesis and infection severity in the variants tested and may further change in future variants (4, 40). This is particularly evident in our elucidation of Delta and Omicron cytokine profiles suggesting a major change in the biology of the VOCs late in the pandemic (41, 42). Indeed, the emergence of Delta and Omicron variants both caused major outbreaks that dwarfed those that of their predecessors (43, 44). Studies have shown that Delta and Omicron VOCs has evolved into their distinct evolutionary group and therefore may explain the significant differences of their infection profiles and pathogenesis mechanisms (41). In addition, other reasons such as gestation and adaptation in immunocompromised patients (45), spill-over and spill-back to and back from animal reservoirs (45–47), and potential vaccine induced selection pressure (48), may all contribute to the accumulation of these mutations and the differential pathogenesis mechanisms. Overall, our study inferred that the adaptation of the virus may likely lead to emergence of future variants having mutations that drives their pathogenesis differently than existing SARS-CoV-2 viruses. Therefore, this implied that management strategies based on immunomodulation and mediation of host responses have to be revisited in different variants, mirroring the management of influenza, where propensity of cytokine storm differs (49). Finally, our finding may also partially explain the continued efficacy of broad spectrum immunomodulators like corticosteroids; but not more targeted therapies targeting specific pathways due to the potential differences in different SARS-CoV-2 iterations’ pathogenesis mechanisms.

Our study, however, is not without its limitation. We only elucidated the histopathology and cytokine profile at 4 dpi, that only gave us insights of host response early during infection. This also created potential disconnect such as the case between the high Delta nucleocapsid detection in IHC compared to its low lung infectious viral titre. Nevertheless, high viral protein, especially in the case of nucleocapsid, which were produced in excess during infection (28), may not necessarily translate to infectious viral titre. In addition, we only use female mice for our comparative study and thus may not capture any potential sex specific differences between the cytokine profiles of different virus iterations (50). However, a previous study showed that cytokine and chemokine responses between male and female were largely similar, where the major differences lie at the magnitude and timing of the responses (51). Therefore, our study remained representative of the differences in cytokine profiles between variants despite only testing it in female mice. Nevertheless, future studies with male mice can be performed to investigate the extent of differential cytokine expression levels that the different SARS-CoV-2 virus iteration can induce. Finally, we also did not manage to assess the protein levels of the cytokines in our study. Despite the limitations however, our data showed that, by systematically comparing the profiles across the major evolutionary iterations of SARS-CoV-2, a clear heterogeneity of SARS-CoV-2 evolution and adaptation was observed, especially in VOCs Delta and Omicron that emerged later

In conclusion, our study provided insights on the evolution of pathogenesis and its potential mechanisms of SARS-CoV-2 throughout the course of the pandemic. We have observed that adaptation throughout the pandemic can give rise to variants of very different pathogenesis mechanisms that contribute to tissue damage. Our study thus emphasized the importance in differentiating host responses to different SARS-CoV-2 variants. This is especially true as we enter the endemic phase, where the surveillance of differential host responses between variants will become increasingly crucial when assessing future variants, for their threat and impact, as well as when devising immunomodulatory treatments for severe infections.

Funding

This research was supported by the following grants: NUHSRO/2020/066/NUSMedCovid/01/BSL3 Covid Research Work, NUHSRO/2020/050/RO5+5/NUHS-COVID/4, Ministry of Education, Singapore MOE2017-T2-2-014, Singapore NMRC Centre Grant Program – Diabetes, Tuberculosis and Neuroscience CGAug16M009, Ministry of Health MOH-COVID19RF2-0001.

Acknowledgments

We are grateful to the National University of Singapore, Yong Loo Lin School of Medicine BSL-3 Core Facility for their support with this work.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

• 1

  HuBGuoHZhouPShiZL. Characteristics of SARS-CoV-2 and COVID-19. Nat Rev Microbiol (2021) 19:141–54. doi: 10.1038/s41579-020-00459-7

  • CrossRef
  • Google Scholar

• 2

  MizrahiBShiloSRossmanHKalksteinNMarcusKBarerYet al. Longitudinal symptom dynamics of COVID-19 infection. Nat Commun (2020) 11:6208. doi: 10.1038/s41467-020-20053-y

  • CrossRef
  • Google Scholar

• 3

  TaoKTzouPLNouhinJGuptaRKde OliveiraTKosakovsky PondSLet al. The biological and clinical significance of emerging SARS-CoV-2 variants. Nat Rev Genet (2021) 22:757–73. doi: 10.1038/s41576-021-00408-x

  • CrossRef
  • Google Scholar

• 4

  GribbleJStevensLJAgostiniMLAnderson-DanielsJChappellJDLuXet al. The coronavirus proofreading exoribonuclease mediates extensive viral recombination. PloS Pathog (2021) 17:e1009226. doi: 10.1371/journal.ppat.1009226

  • CrossRef
  • Google Scholar

• 5

  MaxmenA. One million coronavirus sequences: Popular genome site hits mega milestone. Nature (2021) 589:21. doi: 10.1038/d41586-021-01069-w

  • CrossRef
  • Google Scholar

• 6

  ElbeSBuckland-MerrettG. Data, disease and diplomacy: GISAID’s innovative contribution to global health. Glob Chall (2017) 1:33–46. doi: 10.1002/gch2.1018

  • CrossRef
  • Google Scholar

• 7

  WuFZhaoSYuBChenYMWangWSongZGet al. A new coronavirus associated with human respiratory disease in China. Nature (2020) 579:265–9. doi: 10.1038/s41586-020-2008-3

  • CrossRef
  • Google Scholar

• 8

  World Health, O. COVID-19 weekly epidemiological update. edition 68. Geneva: World Health Organization (2021).

  • Google Scholar

• 9

  YangTJYuPYChangYCLiangKHTsoHCHoMRet al. Effect of SARS-CoV-2 B.1.1.7 mutations on spike protein structure and function. Nat Struct Mol Biol (2021) 28:731–9. doi: 10.1038/s41594-021-00652-z

  • CrossRef
  • Google Scholar

• 10

  TegallyHWilkinsonEGiovanettiMIranzadehAFonsecaVGiandhariJet al. Detection of a SARS-CoV-2 variant of concern in south Africa. Nature (2021) 592:438–43. doi: 10.1038/s41586-021-03402-9

  • CrossRef
  • Google Scholar

• 11

  WangPNairMSLiuLIketaniSLuoYGuoYet al. Antibody resistance of SARS-CoV-2 variants B.1.351 and B.1.1.7. Nature (2021) 593:130–5. doi: 10.1038/s41586-021-03398-2

  • CrossRef
  • Google Scholar

• 12

  BoehmEKronigINeherRAEckerleIVetterPKaiserLet al. Novel SARS-CoV-2 variants: the pandemics within the pandemic. Clin Microbiol Infect (2021) 27:1109–17. doi: 10.1016/j.cmi.2021.05.022

  • CrossRef
  • Google Scholar

• 13

  LeeTYLeeHKimNJeonPKimJWLimHYet al. Comparison of SARS-CoV-2 variant lethality in human angiotensin-converting enzyme 2 transgenic mice. Virus Res (2021) 305:198563. doi: 10.1016/j.virusres.2021.198563

  • CrossRef
  • Google Scholar

• 14

  SuCMWangLYooD. Activation of NF-kappaB and induction of proinflammatory cytokine expressions mediated by ORF7a protein of SARS-CoV-2. Sci Rep (2021) 11:13464. doi: 10.1038/s41598-021-92941-2

  • CrossRef
  • Google Scholar

• 15

  ChenGWuDGuoWCaoYHuangDWangHet al. Clinical and immunological features of severe and moderate coronavirus disease 2019. J Clin Invest (2020) 130:2620–9. doi: 10.1172/JCI137244

  • CrossRef
  • Google Scholar

• 16

  HuangCWangYLiXRenLZhaoJHuYet al. Clinical features of patients infected with 2019 novel coronavirus in wuhan, China. Lancet (2020) 395:497–506. doi: 10.1016/s0140-6736(20)30183-5

  • CrossRef
  • Google Scholar

• 17

  ZhouZRenLZhangLZhongJXiaoYJiaZet al. Heightened innate immune responses in the respiratory tract of COVID-19 patients. Cell Host Microbe (2020) 27:883–890 e882. doi: 10.1016/j.chom.2020.04.017

  • CrossRef
  • Google Scholar

• 18

  CoperchiniFChiovatoLRicciGCroceLMagriFRotondiMet al. The cytokine storm in COVID-19: Further advances in our understanding the role of specific chemokines involved. Cytokine Growth Factor Rev (2021) 58:82–91. doi: 10.1016/j.cytogfr.2020.12.005

  • CrossRef
  • Google Scholar

• 19

  YangLXieXTuZFuJXuDZhouYet al. The signal pathways and treatment of cytokine storm in COVID-19. Signal Transduct Target Ther (2021) 6:255. doi: 10.1038/s41392-021-00679-0

  • CrossRef
  • Google Scholar

• 20

  ChannappanavarRFehrARVijayRMackMZhaoJMeyerholzDKet al. Dysregulated type I interferon and inflammatory monocyte-macrophage responses cause lethal pneumonia in SARS-CoV-Infected mice. Cell Host Microbe (2016) 19:181–93. doi: 10.1016/j.chom.2016.01.007

  • CrossRef
  • Google Scholar

• 21

  HalfmannPJIidaSIwatsuki-HorimotoKMaemuraTKisoMScheafferSMet al. SARS-CoV-2 omicron virus causes attenuated disease in mice and hamsters. Nature (2022) 603:687–92. doi: 10.1038/s41586-022-04441-6

  • CrossRef
  • Google Scholar

• 22

  ShrumBAnanthaRVXuSXDonnellyMHaeryfarSMMcCormickJKet al. A robust scoring system to evaluate sepsis severity in an animal model. BMC Res Notes (2014) 7:1–11. doi: 10.1186/1756-0500-7-233

  • CrossRef
  • Google Scholar

• 23

  LiFHanMDaiPXuWHeJTaoXet al. Distinct mechanisms for TMPRSS2 expression explain organ-specific inhibition of SARS-CoV-2 infection by enzalutamide. Nat Commun (2021) 12:866. doi: 10.1038/s41467-021-21171-x

  • CrossRef
  • Google Scholar

• 24

  LeistSRDinnonKH3rdSchaferATseLVOkudaKHouYJet al. A mouse-adapted SARS-CoV-2 induces acute lung injury and mortality in standard laboratory mice. Cell (2020) 183:1070–85 e1012. doi: 10.1016/j.cell.2020.09.050

  • CrossRef
  • Google Scholar

• 25

  Huang daWShermanBTLempickiRA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc (2009) 4:44–57. doi: 10.1038/nprot.2008.211

  • CrossRef
  • Google Scholar

• 26

  ShermanBTHaoMQiuJJiaoXBaselerMWLaneHCet al. DAVID: a web server for functional enrichment analysis and functional annotation of gene lists (2021 update). Nucl Acids Res (2022) 50:W216–21. doi: 10.1093/nar/gkac194

  • CrossRef
  • Google Scholar

• 27

  SurjitMLalSK. The SARS-CoV nucleocapsid protein: a protein with multifarious activities. Infect Genet Evol (2008) 8:397–405. doi: 10.1016/j.meegid.2007.07.004

  • CrossRef
  • Google Scholar

• 28

  SavastanoAIbanez de OpakuaARankovicMZweckstetterM. Nucleocapsid protein of SARS-CoV-2 phase separates into RNA-rich polymerase-containing condensates. Nat Commun (2020) 11:6041. doi: 10.1038/s41467-020-19843-1

  • CrossRef
  • Google Scholar

• 29

  LuLZhangHDaupharsDJ. & he, y. w. a potential role of interleukin 10 in COVID-19 pathogenesis. Trends Immunol (2021) 42:3–5. doi: 10.1016/j.it.2020.10.012

  • CrossRef
  • Google Scholar

• 30

  OladunniFSParkJGPinoPAGonzalezOAkhterAAllue-GuardiaAet al. Lethality of SARS-CoV-2 infection in K18 human angiotensin-converting enzyme 2 transgenic mice. Nat Commun (2020) 11:6122. doi: 10.1038/s41467-020-19891-7

  • CrossRef
  • Google Scholar

• 31

  YindaCKPortJRBushmakerTOffei OwusuIPurushothamJNAvanzatoVAet al. K18-hACE2 mice develop respiratory disease resembling severe COVID-19. PloS Pathog (2021) 17:e1009195. doi: 10.1371/journal.ppat.1009195

  • CrossRef
  • Google Scholar

• 32

  ZhengJWongLRLiKVermaAKOrtizMEWohlford-LenaneCet al. COVID-19 treatments and pathogenesis including anosmia in K18-hACE2 mice. Nature (2021) 589:603–7. doi: 10.1038/s41586-020-2943-z

  • CrossRef
  • Google Scholar

• 33

  SunSHChenQGuHJYangGWangYXHuangXYet al. A mouse model of SARS-CoV-2 infection and pathogenesis. Cell Host Microbe (2020) 28:124–133 e124. doi: 10.1016/j.chom.2020.05.020

  • CrossRef
  • Google Scholar

• 34

  DongWMeadHTianLParkJGGarciaJIJaramilloSet al. The K18-human ACE2 transgenic mouse model recapitulates non-severe and severe COVID-19 in response to an infectious dose of the SARS-CoV-2 virus. J Virol (2022) 96:e0096421. doi: 10.1128/JVI.00964-21

  • CrossRef
  • Google Scholar

• 35

  YangSCaoLXuWXuTZhengBJiYet al. Comparison of model-specific histopathology in mouse models of COVID-19. J Med Virol (2022) 94:3605–12. doi: 10.1002/jmv.27747

  • CrossRef
  • Google Scholar

• 36

  McCrayPBJr.PeweLWohlford-LenaneCHickeyMManzelLShiLet al. Lethal infection of K18-hACE2 mice infected with severe acute respiratory syndrome coronavirus. J Virol (2007) 81:813–21. doi: 10.1128/JVI.02012-06

  • CrossRef
  • Google Scholar

• 37

  GamageAMTanKSChanWOYLiuJTanCWOngYKet al. Infection of human nasal epithelial cells with SARS-CoV-2 and a 382-nt deletion isolate lacking ORF8 reveals similar viral kinetics and host transcriptional profiles. PloS Pathog (2020) 16:e1009130. doi: 10.1371/journal.ppat.1009130

  • CrossRef
  • Google Scholar

• 38

  HuiKPYHoJCWCheungMCNgKCChingRHHLaiKLet al. SARS-CoV-2 omicron variant replication in human bronchus and lung ex vivo. Nature (2022) 603:715–20. doi: 10.1038/s41586-022-04479-6

  • CrossRef
  • Google Scholar

• 39

  ShuaiHChanJFHuBChaiYYuenTTYinFet al. Attenuated replication and pathogenicity of SARS-CoV-2 B.1.1.529 omicron. Nature (2022) 603:693–9. doi: 10.1038/s41586-022-04442-5

  • CrossRef
  • Google Scholar

• 40

  PlanteJAMitchellBMPlanteKSDebbinkKWeaverSCMenacheryVDet al. The variant gambit: COVID-19’s next move. Cell Host Microbe (2021) 29:508–15. doi: 10.1016/j.chom.2021.02.020

  • CrossRef
  • Google Scholar

• 41

  BansalKKumarS. Mutational cascade of SARS-CoV-2 leading to evolution and emergence of omicron variant. Virus Res (2022) 315:198765. doi: 10.1016/j.virusres.2022.198765

  • CrossRef
  • Google Scholar

• 42

  ZhanYYinH. & yin, j. y. B.1.617.2 (Delta) variant of SARS-CoV-2: features, transmission and potential strategies. Int J Biol Sci (2022) 18:1844–51. doi: 10.7150/ijbs.66881

  • CrossRef
  • Google Scholar

• 43

  TareqAMEmranTBDhamaKDhawanMTalleiTE. Impact of SARS-CoV-2 delta variant (B.1.617.2) in surging second wave of COVID-19 and efficacy of vaccines in tackling the ongoing pandemic. Hum Vaccin Immunother (2021) 17:4126–7. doi: 10.1080/21645515.2021.1963601

  • CrossRef
  • Google Scholar

• 44

  PazMAldunateFArceRFerreiroICristinaJ. An evolutionary insight into severe acute respiratory syndrome coronavirus 2 omicron variant of concern. Virus Res (2022) 314:198753. doi: 10.1016/j.virusres.2022.198753

  • CrossRef
  • Google Scholar

• 45

  KupferschmidtK. Where did ‘weird’ omicron come from? Science (2021) 374:1179. doi: 10.1126/science.acx9738

  • CrossRef
  • Google Scholar

• 46

  IslamAFerdousJIslamSSayeedMADutta ChoudhurySSahaOet al. Evolutionary dynamics and epidemiology of endemic and emerging coronaviruses in humans, domestic animals, and wildlife. Viruses (2021) 13:1908. doi: 10.3390/v13101908

  • CrossRef
  • Google Scholar

• 47

  ChenQHuangXYLiuYSunMXJiBZhouCet al. Comparative characterization of SARS-CoV-2 variants of concern and mouse-adapted strains in mice. J Med Virol (2022) 94:3223–32. doi: 10.1002/jmv.27735

  • CrossRef
  • Google Scholar

• 48

  KoyamaTMiyakawaKTokumasuRSSJKudoMRyoA. Evasion of vaccine-induced humoral immunity by emerging sub-variants of SARS-CoV-2. Future Microbiol (2022) 17:417–24. doi: 10.2217/fmb-2022-0025

  • CrossRef
  • Google Scholar

• 49

  RyabkovaVAChurilovLPShoenfeldY. Influenza infection, SARS, MERS and COVID-19: Cytokine storm - the common denominator and the lessons to be learned. Clin Immunol (2020) 223:108652. doi: 10.1016/j.clim.2020.108652

  • CrossRef
  • Google Scholar

• 50

  PeckhamHde GruijterNMRaineCRadziszewskaACiurtinCWedderburnLRet al. Male Sex identified by global COVID-19 meta-analysis as a risk factor for death and ITU admission. Nat Commun (2020) 11:6317. doi: 10.1038/s41467-020-19741-6

  • CrossRef
  • Google Scholar

• 51

  QiSNgwaCMorales ScheihingDAAl MamunAAhnstedtHWFingerCEet al. Sex differences in the immune response to acute COVID-19 respiratory tract infection. Biol Sex Differ (2021) 12:66. doi: 10.1186/s13293-021-00410-2

  • CrossRef
  • Google Scholar