Abstract
Coronavirus Disease 19 (COVID-19) is a respiratory disease caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), which has grown to a worldwide pandemic with substantial mortality. Immune mediated damage has been proposed as a pathogenic factor, but immune responses in lungs of COVID-19 patients remain poorly characterized. Therefore we conducted transcriptomic, histologic and cellular profiling of post mortem COVID-19 (n=34 tissues from 16 patients) and normal lung tissues (n=9 tissues from 6 patients). Two distinct immunopathological reaction patterns of lethal COVID-19 were identified. One pattern showed high local expression of interferon stimulated genes (ISGhigh) and cytokines, high viral loads and limited pulmonary damage, the other pattern showed severely damaged lungs, low ISGs (ISGlow), low viral loads and abundant infiltrating activated CD8+ T cells and macrophages. ISGhigh patients died significantly earlier after hospitalization than ISGlow patients. Our study may point to distinct stages of progression of COVID-19 lung disease and highlights the need for peripheral blood biomarkers that inform about patient lung status and guide treatment.
Introduction
COVID-19 is a pandemic respiratory disease with 2-3% lethality and a particularly severe course in males, patients with cardiovascular comorbidities and in the elderly1,2. Lymphopenia, high levels of pro-inflammatory cytokines in the circulation3, and phenotypic changes of pro-inflammatory macrophages in bronchoalveolar lavages4 in severe patients have led to the notion that the immune response against the causative virus SARS-CoV-2 may contribute to devastating end organ damage5. Since patients with severe COVID-19 may develop acute respiratory distress syndrome (ARDS) and many patients die from respiratory failure with diffuse alveolar damage6, it is critical to understand the immunological profiles in the lungs of these patients.
To better understand the molecular and cellular underpinnings of COVID-19 lung disease, we used histologic and transcriptional analyses of post mortem lung tissues in a cohort of patients where the cause of death was respiratory failure. We detected two distinct immunological and cellular profiles in the lungs of these patients, defined by their differential expression of interferon stimulated genes (ISGs) and immune infiltration patterns, which we termed ISGhigh and ISGlow. ISG subgroups strongly differed in regards to the characteristics and the extent of pulmonary damage, pulmonary viral loads, disease course and time to death from hospitalization. These data highlight two distinct patterns of immune pathology of pulmonary COVID-19 and may give insight into the natural progression of COVID-19 in lungs.
Results
Two patterns of gene expression in COVID-19 autopsy lungs
Here we analyzed 34 post mortem lung samples from 16 deceased COVID-19 patients and 9 post mortem lung samples from 6 patients, who died from non-infectious causes (Table 1). The primary cause of death in all patients of this cohort was respiratory failure, sometimes multi-organ failure including failure of the respiratory system. Since lung samples from the same patients did not always appear morphologically uniform, all lung specimens were subjected to differential gene expression analysis based on a commercially available targeted next generation sequencing (NGS) assay (OIRRA) designed for quantification of immune cell and inflammatory transcripts (Supplementary Table 2). Among the 398 genes investigated, we identified 68 up-regulated and 30 down-regulated genes in COVID-19 infected lungs compared to normal tissue (Figure 1a,b; Supplementary Table 3), and a PCA analysis showed segregation of COVID-19 patients in two well-defined clusters that showed distinct association with viral load (Figure 1a,c).
Using a consensus of 30 different indices7 we identified three groups of samples defined by distinct expression levels of the deregulated genes by K-means clustering (Figure 1a). Clusters 1 (50% of samples) and 2 (41%) contained COVID-19 samples while cluster 3 contained all normal lung samples as well as three COVID-19 samples (9%). To understand why the majority of COVID-19 lung tissues segregated into defined clusters 1 or 2, we undertook a gene ontology analysis. We identified ISGs as a key upregulated pathway in COVID-19 autopsy lungs (Table 2), which was differentially represented in cluster 1 and 2, respectively (Figure 1d). Identification of an ISGhigh cluster (Cluster 1, ISGhigh) was surprising, since SARS-CoV-2 was recently proposed to lead to limited ISG induction, yet only in comparison to other respiratory viruses8. Our data suggest that autopsy lungs of COVID-19 patients, who died from respiratory failure, showed two different gene expression patterns with different levels of ISG activation (ISGhigh and ISGlow).
Clinical differences of COVID-19 patients with an ISGhigh versus ISGlow lung profile
Patients, whose lung samples all segregated in the ISGhigh or in the ISGlow subgroups, were called ISGhigh and ISGlow patients. To investigate whether there were clinical differences between these two patient groups, we compared their clinical and epidemiological information. The majority of COVID-19 patients in our cohort (81%) were male and the average body mass index (BMI) was 31.4. kg/m2. Patient level analysis revealed no correlation of sex or BMI with the ISG patterns (Figure 1a). All patients in our cohort died from respiratory failure or multi-organ failure including failure of the respiratory system, independent of ISG subgrouping. When we analysed comorbidities and autopsy findings, we found that 5 out of 7 (71%) ISGlow patients, but none of the ISGhigh patients had an autoptic finding of a thromboembolic event in the lungs and/or disseminated intravascular dissemination (DIC) indicating an abnormally activated blood coagulation (hypercoagulability) exclusively in ISGlow patients. COVID-19 associated coagulopathy was unlikely to contribute to rapid exacerbation of pulmonary COVID-19 since this diagnostic criterion was only observed in patients displaying the ISGlow phenotype. Most notably, when we compared the disease course and hospitalization time between the patients with different ISG lung profiles, we found a significantly longer hospitalization time for ISGlow patients from admission to death compared to ISGhigh patients (Figure 1e). Of note, we did not find differences in the anamnestic onset of disease symptoms and time from positive testing for COVID-19 infection by nasopharyngeal swab to hospitalization or death between ISGhigh and ISGlow patients. To exclude bacterial and viral superinfections as a confounder of clinical course, we performed whole genome sequencing on all samples to detect bacterial and/or viral DNA. Bacterial superinfections were found in three lung tissue samples, in 3/16 COVID-19 patients, that were equally distributed among the different groups (Figure 2a-e). Based on the limited number of samples with evidence of bacterial superinfection, there was no correlation with clinical subgroups or the duration of the disease.
Taken together, expression of the ISGhigh profile in COVID-19 lungs is associated with accelerated disease course and early lethal outcome. The ISGlow profile is associated with coagulopathies and later lethal outcome. This points to two distinct clinical courses of fatal COVID-19 pulmonary disease. Since autopsy studies always focus on disease endpoints, our data do not allow us to draw direct conclusions about natural stages of disease progression in COVID-19 lungs. However, it is very suggestive that the ISGhigh profile precedes the ISGlow profile, consistent with a longitudinal study in peripheral blood showing that ISG expression was high in early COVID-19 and declined later9.
Immune microenvironment characteristics of the ISGhigh and ISGlow COVID-19 lung profiles
In line with a recent study showing a correlation of ISG expression and viral load in nasopharyngeal swabs10, expression of ISGs was positively correlated with pulmonary viral load (Figure 3a), and immunohistochemical staining confirmed the presence of SARS-CoV-2 nucleocapsid protein in ISGhigh lungs, mainly localized to pneumocytes (Figure 3b).
SARS-CoV-2 induces a strong antiviral immune response. Therefore we analysed frequencies of specific immune cells in lungs by computational image analysis. T cells (CD3+) of the CD4+ and CD8+ lineages, B cells (CD20+), and macrophages (CD68+) were selectively enriched in lung tissues from ISGlow patients (Figure 3c,d; Figure 4a,b). A strong enrichment for CD68+ and CD163+ monocytes in lung tissue was observed with spatial correlation of stains for both markers indicating co-expression. Since circulating monocytes in COVID-19 patients co-express CD68 and CD163, it was not surprising that CD68 and CD163 expression in lungs followed a similar pattern11 (Figure 3c,d). Surprisingly, CD123+ plasmacytoid dendritic cells (pDCs) did not show elevated frequencies in ISGhigh lungs (Figure 4a,b), and the type 1 IFN (IFN-I) genes that we measured (IFNa17, IFNb1) were not higher expressed in ISGhigh lungs with high viral load compared to ISGlow lungs (data not shown), potentially because SARS-CoV-2 inhibits IFN-I production8. Therefore our analysis did not allow us to identify the upstream trigger of ISGs in lungs.
A cytokine storm has been proposed to cause adverse outcome of COVID-1912. It has been suggested that peripheral monocytes do not contribute to it13, yet cytokines are highly expressed in bronchoalveolar lavages (BALs) of COVID-19 patients4. Therefore we investigated expression of a pro-inflammatory cytokine signature (TNF, IL1B, IL6, CCL2, IFNA17, IFNB1, CXCL9, CXCL10, CXCL11) in lung samples from lethal COVID-19. The proinflammatory gene signature was significantly enriched in the ISGhigh subset (p=0.0061) (Figure 5a). Activated CD8+ T cells are essential for elimination of coronaviruses14,15. Therefore we defined and investigated an activated cytotoxic T cell signature (CD38, GZMA, GZMB, CCR5) and found that it was inversely correlated to viral counts, particularly in ISGlow cases (Figure 5b). This suggests that activated CD8+ T cells may indeed contribute to the elimination of SARS-CoV-2 in lungs.
Taken together these data show that COVID-19 autopsy lungs with an ISGhigh profile show high virus titers, high local expression of innate cytokines and weak immune cell infiltration, while COVID-19 lungs with an ISGlow profile show low virus titers, lower expression of innate cytokines and strong immune cell infiltration. This pattern could indicate expression of the lung ISGhigh profile at an earlier, innate disease stage, i.e. at a time when the virus is not yet controlled, and expression of the ISGlow profile at a later disease stage, i.e. after T cell priming.
The ISGhigh and ISGlow lung immunoprofiles correlate with morphological changes
To investigate the potential immunological causes for lung damage in COVID-19, we studied whether ISG profiles in COVID-19 post mortem lungs were associated with specific histomorphological features of fatal COVID-19. Diffuse alveolar damage (DAD) was mostly found in ISGlow patients (Figure 5c), but intra-alveolar hemorrhage (IAH) was not associated with lung ISG status.
As the cytokine storm was implicated in decline of COVID-19 patients, we analysed whether expression of the above defined pro-inflammatory cytokine signature was associated with IAH but this was not the case (Figure 5d). However, within this cytokine signature we identified co-regulated subgroups (IL1B/IL6/TNF, IFNB1/IFNA17, CCL2/CXCL9/CXCL19/CXCL11) (Figure 5e). Of these, the CXCL9/10/11 sub-signature was positively associated with IAH (Figure 5f-i). This is in line with observations that these chemokines compromise endothelial integrity via CXCR316, and that CXCL10 is a key determinant of severe COVID-1917. Interestingly, basal levels of CXCL9 or CXCL10 are elevated in elderly, hypertensive and obese individuals, who were strongly represented in our autopsy cohort (Table 1) and are predisposed to severe COVID-1918,19.
It has been proposed that infiltrating monocytes and macrophages play a role in lung damage4,20. In support of this data, we found CD68+ macrophage infiltrates to be positively associated with DAD (Figure 5j). In addition, DAD was associated with the activated cytotoxic T cell signature (p=0.0022) (Figure 5k), yet not with the overall numbers of pulmonary CD8+ T-cells (Figure 5l). This raises the possibility that activated CD8+ T cells contribute to DAD as they eliminate virus from infected lungs. None of the above pulmonary cytokine sub-signatures, however, was positively associated with DAD (Figure 5m-p), suggesting that none of these cytokines drives lung pathology directly.
In summary, we did not find distinct features of lung damage in ISGhigh patients, suggesting that extra-pulmonary factors may contribute to mortality in these patients. However, ISGlow patients show prominent DAD, associated with peri-alveolar foci of CD68+ macrophages and an activated T cell signature. Local expression of most cytokines did not correlate to lung damage, except for CXCL9/10/11, which correlated to IAH (p=0.018) (Figure 5i). Based on studies that associated serum CXCL10 levels with general disease severity17, it will be interesting to investigate whether serum levels of these cytokines predict a specific pattern of lung damage.
The ISGlow lung profile shows signs of tissue regeneration and T cell exhaustion
Since ISGlow lung samples were derived from patients with a longer disease course, we investigated specific pathways of local immune regulation and tissue regeneration. ISGlow lung samples expressed elevated p53 and Ki67 (Figure 1a, Figure 6a), i.e. reactive markers indicating lung remodeling in DAD21.
Since we found local upregulation of C1QA (p=0.017) and C1QB (p=0.0012) specifically in ISGlow lungs (Figure 1a, Figure 6b), we hypothesized that complement activation may further contribute to lung damage in these patients. Consistently, we found strong staining for C3d and C5b-9 complex deposition in lung tissue indicating complement activation in lungs of ISGlow patients (Figure 6c). Since C1Q also restrains antiviral CD8+ effector T cell responses22, it may contribute to the local regulation of effector T cells. In line with previous observations23, we found a higher frequency of CD8+PD1+ T-cells in the ISGlow subgroup (p=0.001), potentially indicative of T cell exhaustion (Figure 4a,b).
Overall our results identify two patterns of pulmonary COVID-19 disease that lead to death from respiratory failure. Patients of the ISGhigh subgroup die early with high viral loads and high cytokine and ISG expression levels in lungs. Their lungs are morphologically relatively intact, and our data do not identify a uniform pathomechanism underlying lethal outcome, although some show CXCL9/10/11-associated IAH (Figure 5i). The distinct ISGlow group of patients dies later, with low viral loads in the lungs, low local expression of cytokines and ISGs, yet strong infiltration of pulmonary tissue by CD8+ T cells and macrophages, which both correlate to severity of DAD and local complement activation. Some of these patients show IAH in addition to DAD, and many of them suffer from coagulopathies. Altogether this patient group appears to suffer from severe pulmonary immunopathology. The design of our autopsy study does not allow to directly conclude that ISGlow lungs might have undergone a previous ISGhigh phase, though circumstantial evidence about the general course of coronavirus infections may suggest this possibility.
Discussion
Here we describe two immunopathological patterns in lungs of fatal COVID-19 patients based on ISG expression. The ISGhigh pattern is observed in patients, who die early after hospitalization and is characterized by high viral load and high levels of pro-inflammatory cytokines, yet relatively intact lung morphology, while the ISGlow pattern is characterized by low viral load, massive lung damage, marked lung immune cell infiltrates and late death. Our findings are consistent with epidemiological data showing two peaks of mortality24, and another study of four COVID-19 autopsies, where one patient died early after hospital admission, with striking upregulation of pulmonary IL-1b/IL-6 in lungs and little lung damage, while three patients expressed low levels of pulmonary cytokines, massive DAD and delayed death25. Therefore our study allows us to propose two immunopathological stages of pulmonary COVID-19.
The segregation of autopsy lung samples from COVID-19 patients in two groups based on ISG expression contributes to our understanding of the interferon response against SARS-CoV-2. Like other coronaviruses, SARS-CoV-2 is particularly sensitive to type I interferons8,26. Therefore, and similar to other coronaviruses27, it has evolved strategies to evade the interferon response, and SARS-CoV-2 leads to relatively weak IFN-I/III release in host cells at low multiplicities of infection8. This initial delay of IFN-I/III production may facilitate initial virus replication in lungs, as studies with SARS-CoV in mice have shown, and enable an eventual increase of the IFN-I response and death15. A similar observation was made for fatal SARS-CoV infections in humans, which were accompanied by elevated expression of ISGs28. Since the SARS-CoV-2 receptor ACE2 is itself an ISG on lung epithelial cells29, virus infection and the interferon response may promote each other in this phase of the infection. This may explain the observed correlation of high ISG expression and high viral load in lungs and widespread presence of SARS-CoV-2 in lung epithelial cells. Together this may contribute to fatal outcome of SARS-CoV-2 infections in the ISGhigh group.
The observation of the ISGhigh pattern in COVID-19 autopsy lungs seems to be at odds with initial observations that critical COVID-19 patients express on average lower ISGs in blood than patients with a milder course of disease30. One possible explanation is that the blood ISG status is not directly reflective of the lung. In support of this idea, bronchoalveolar lavages (BALs) from severe/critical patients showed high proportions of ISGhigh macrophages and high expression of CXCL9/10/11, IL6, IL-1b, TNF, CCL24, which is reflective of our ISGhigh phenotype in lung autopsies. An alternative explanation for the apparent disconnect of lung and blood ISG status may come from the overall frequency of the ISGhigh subtype of critical/fatal COVID-19: ISGs are highly expressed in blood of COVID-19 patients during the early, innate phase of COVID-199, and we show that patients with an ISGhigh status in lungs die early upon hospitalization. While the early-mortality subset accounted for 44% of all deaths in our study, epidemiological data from France24 suggest that this early critical subset is actually smaller: only 15% of patients died early after hospitalization. This percentage is consistent with the study by Hadjadj et al., which detects high ISG expression in blood of 3/17 (18%) critical COVID-19 patients, yet this signal gets diluted in the majority of ISGlow critical cases30. These two alternative explanations show how critical it will be to compare gene expression in blood and lungs of individual patients at different times of the infection and to identify peripheral biomarkers for COVID-19 lung status.
ISGlow COVID-19 patients in our study die with classical features of DAD31, on average 9.1 days after hospitalization. Later death compared to patients with an ISGhigh pattern and progressive decline of systemic ISG expression during COVID-199,13 led us to infer that the ISGlow pattern in lungs reflects a later phase of pulmonary COVID-19. ISGlow lungs show higher frequencies of T and B lymphocytes, compared to ISGhigh lungs. None of our fatal cases showed lung lymphocyte counts below control levels. Therefore COVID-19 associated lymphopenia in blood23,32 or spleens3 does not translate into lymphocyte depletion in infected lungs. Potential reasons are that the infected lung acts as a potent sink for circulating lymphocytes and that local proliferation limited recruitment from the blood, as was shown for CD8+ T cells in BAL of severe patients4. Consistent with previous observations33 we describe an activated CD8+ T cell signature in lungs of ISGlow patients that contain low viral counts. This suggests that CD8+ T cells are critical for antiviral protection, and may transition into a protective memory pool, as observed for SARS-CoV34,35. In addition, we found elevated frequencies of CD8+PD1+ cells in ISGlow lungs compared to ISGhigh lungs, but not above control levels. The observation that PD-1 levels are elevated in peripheral CD8+ T cells of severe COVID-19 infection, and whether this indicates exhaustion, remains controversial13,23. Overall, although we did not have paired serum antibody levels available, the infiltration pattern of ISGlow lungs suggested adaptive immune activation.
While our study sheds further light on COVID-19 lung disease, conclusions on therapy must be drawn with caution. We found that early after hospitalization, ISGhigh autopsy lungs had uniformly high titers of SARS-CoV-2, and others found that viral loads in swabs and sputum are highest in early COVID-1936. This could indicate that treatment with compounds that directly interfere with the SARS-CoV-2 replication cycle, e.g. protease or polymerase inhibitors, should start early. However, high expression of ISGs in some lung autopsies raises caution about the use of IFN-I/III as therapeutics, at least as long as causes and consequences of interferon signaling in COVID-19 lungs remain unclear. We found reduced viral counts in ISGlow patients, but did not identify the moment at which the body is cleared of the virus. Therefore our findings of reduced viral loads in ISGlow patients should not be taken as justification to withhold compounds that directly interfere with the SARS-CoV-2 replication cycle from patients. Extending previous work37, we found signs of elevated complement activation specifically in ISGlow lungs. However, it is not known whether complement is synchronously activated in patient lungs and plasma38. Therefore, our results do not provide a further step towards personalized patient care but strengthen the hypothesis that complement inhibitors may show therapeutic benefit, at least in some COVID-19 patients.
Our study has several limitations. The fact that patients with the innate / ISGhigh stage die early while patients with the lymphocytic / ISGlow profile die late after hospitalization, together with knowledge about the immune reaction against other coronaviruses, strongly suggests that COVID-19 lung disease progresses from an ISGhigh to an ISGlow stage. However, an autopsy study is, by design, not longitudinal. Therefore we do not have formal proof that all COVID-19 infected ISGlow lungs have undergone a previous ISGhigh stage. Also, it is unknown why some patients die early and others late. In addition, our focus on the lung only allowed us to investigate pulmonary factors of patient mortality, i.e. an overshooting innate immune activation with IAH in both ISGhigh and ISGlow cases and DAD in ISGlow cases. No reported cause of death was enriched in ISGhigh or ISGlow patients, and multi-organ failure was reported as a cause of death only for two of our patients. Another limitation is that we lack gene expression data from blood of autopsy patients or serological data at the time of death. Therefore we were not able to identify peripheral biomarkers predicting specific immunological profiles in the lung. Finally, we analyzed our lungs with a focused gene expression set since quality and quantity of autopsy-derived RNA is often insufficient for unbiased methods. In spite of this technical limitation, which restricted our analysis, we were able to uncover two novel and distinct immunopathological profiles in lungs of fatal COVID-19.
Taken together, our autopsy study sheds light on two distinct courses of lethal COVID-19 in lungs. It remains to be seen whether interferon signaling is only associated with or causally involved in these disease courses. However, our study strengthens the notion that interferon signaling is a central determinant of the pulmonary immune response against SARS-CoV-2.
Methods
Ethics statement
This study was conducted according to the principles expressed in the Declaration of Helsinki. Ethics approval was obtained from the Ethics Committee of Northwestern and Central Switzerland (Project-ID 2020-00629). For all patients, either personal and/or family consent was obtained for autopsy and sample collection.
Patients and sample collection
The study is based on the analysis of 16 out of 21 consecutive COVID-19 autopsies performed between March 9th and April 14th 2020 at the Institute of Pathology Liestal and Institute of Medical Genetics and Pathology Basel, Switzerland. Clinical features including symptoms, course of disease, comorbidities, laboratory results and therapy are listed in Table 1a. Detailed autopsy findings for each patient were recently published, and the identifiers (with the prefix “C”) for each COVID-19 patient are consistent with the description of this Swiss COVID-19 autopsy cohort6. In this study, we analysed formalin fixed and paraffin embedded (FFPE) lung tissue of distinct areas of the lungs of 16 of these COVID-19 patients. All 16 COVID-19 patients had positive nasopharyngeal swabs collected while alive. In all COVID-19 patients, diagnosis was confirmed by detection of SARS-CoV-2 in postmortal lung tissues. 5/16 patients were additionally tested by postmortal nasopharyngeal swabs which were positive for SARS-CoV-2 in all 5 cases.
As a control cohort, we selected 6 autopsies performed between January 2019 and March 2020 at the Institute of Pathology Liestal (“normal” patients N1 – N6). These control patients died of other, non-infectious causes and had a similar age, gender and cardiovascular risk profile. Patients with infections were excluded from this control cohort. Another control cohort consisted of 4 autopsies of patients suffering from various infections mainly with bacteria affecting the lung (patients with lung pathology, P1 – P4). Details for both control cohorts are listed in Table 1b,c. SARS-CoV-2 was ruled out for each control patient by PCR-examination of lung tissue samples.
Nucleic acid extraction
RNA was extracted from up to six sections of FFPE tissue blocks using RecoverAll Total Nucleic Acid Isolation Kit (Cat No. AM1975, Thermo Fisher Scientific, Waltham, MA, USA). Extraction of DNA from up to 10 sections of FFPE tissue samples was automated by EZ1 Advanced XL (Qiagen, Hilden, Germany) using the EZ1 DNA Tissue Kit (Cat No. 953034, Qiagen, Hilden, Germany). Concentration of DNA and RNA were measured with Qubit 2.0 Fluorometer and Qubit dsDNA HS Assay or Qubit RNA HS Assay (Cat No. Q33230 & Q32852, Thermo Fisher Scientific, Waltham, MA, USA), respectively.
Quantification of SARS-CoV-2 in FFPE tissue samples
Post mortem viral load was individually measured in all lung tissue blocks from all patients included in this study. SARS-CoV-2 was detected in 15ng of human total RNA using the TaqMan 2019-nCoV Assay Kit v1 (Cat No. A47532, Thermo Fisher Scientific, Waltham, MA, USA), which targets three genomic regions (ORFab1, S Protein, N Protein) specific for SARS-CoV-2 and the human RNase P gene (RPPH1). The copy numbers of the SARS-CoV-2 viral genome was determined by utilizing the TaqMan 2019-nCoV Control Kit v1 (Cat No. A47533, Thermo Fisher Scientific, Waltham, MA, USA) and a comparative “ΔΔC?” method. The control kit contains a synthetic sample with a defined amount of target molecules for the human RPPH1 and the three SARS-CoV-2 assays, and was re-analyzed in parallel with patient samples. For each patient sample, this method resulted in individual copy numbers of the human RPPH1 and the three SARS-CoV-2 targets. Finally, the mean copy number of the SARS-CoV-2 targets was normalized to 1 × 106 RPPH1 transcripts.
Profiling of immune response by targeted RNAseq
The expression levels of 398 genes, including genes relevant in innate and adaptive immune response and housekeeping genes for normalization, were analyzed with the Oncomine Immune Response Research Assay (OIRRA, Cat No. A32881, Thermo Fisher Scientific, Waltham, MA, USA). The OIRRA is a targeted gene expression assay designed for the Ion™ next-generation sequencing (NGS) platform. Our study focused on the analysis of rare autopsy tissue samples from COVID-19 patients collected in clinical routine during the COVID-19 pandemic. An inherent problem for transcriptomic studies of autopsy tissues is that it is often not possible to extract high quality RNA in sufficient amounts. To avoid sample dropout due to these reasons, we decided to use a robust and straightforward targeted gene expression assay (OIRRA) rather than whole transcriptome analysis. Since the focus of our study was to investigate the immune profile of lungs, an immunoprofiling assay was deemed most appropriate. The OIRRA gene expression assay was originally designed to interrogate the tumor microenvironment to enable mechanistic studies and identification of predictive biomarkers for immunotherapy in cancer. The assay is optimized to measure the expression of genes involved in immune cell interactions and signaling, including genes expressed at low levels and involved in inflammatory signaling. The 398 genes covered by this assay are listed in Supplementary Table 2. The accessibility of such commercially available assays could be an encouragement to hospitals around the world to conduct similar molecular profiling studies of diagnostic tissue samples from COVID-19 patients, allowing relatively fast and easy stratification of patients into distinct biological groups as a starting point for targeted intervention strategies.
The NGS libraries were prepared as recommended by the supplier. In brief, 30ng of total RNA were used for reverse transcription (SuperScript VILO, Cat No. 11754250, Thermo Fisher Scientific, Waltham, MA, USA) and subsequent library preparation. The libraries were quantified (Ion Library TaqMan Quantitation Kit, Cat No. 4468802, Thermo Fisher Scientific, Waltham, MA, USA), equimolarly pooled and sequenced utilizing the Ion GeneStudio S5xl (Thermo Fisher Scientific, Waltham, MA, USA). De-multiplexing and gene expression level quantification were performed with the standard setting of the ImmuneResponseRNA plugin (version 5.12.0.1) within the Torrent Suite (version 5.12.1), provided as part of the OIRRA by Thermo Fisher Scientific, Waltham, MA, USA.
Detection of co-infections by whole genome sequencing
To identify potential pathogens accompanying an infection with SARS-CoV-2, we analyzed the DNA of the same tissue samples used for detection and profiling of the SARS-CoV-2-specific immune response. First, 250ng of genomic DNA was enzymatically sheared (15 minutes at 37°C) and barcoded using the Ion Xpress Plus Fragment Library Kit (Cat No. 4471269, Thermo Fisher Scientific, Waltham, MA, USA). Subsequently, the libraries were quantified (Ion Library TaqMan Quantitation Kit, Cat No. 4468802, Thermo Fisher Scientific, Waltham, MA, USA) and up to three libraries were pooled at equimolar levels for analysis with Ion GeneStudio S5xl (Thermo Fisher Scientific, Waltham, MA, USA). Sequencing data for each sample was analysed using the CLC genomics workbench (version 20.0.3, Qiagen, Hilden, Germany) in combination with the microbial genomics module (version 20.0.1, Qiagen, Hilden, Germany): The raw reads were trimmed by quality (Mott algorithm with limit 0.05 and a maximum of 2 ambiguous bases per read) and mapped to the human genome (GRCh37 hg19, match score: 1, mismatch cost: 2, indel opening cost: 6, indel extension cost: 1). Unmapped reads were analysed by taxonomic profiling to identify reads of viral or bacterial origin. The profiling utilized an index of 11’540 viral genomes with a minimum length of 1’000 bp and 2’715 bacterial reference genomes with a minimum length of 500’000 bp, retrieved from the NCBI Reference Sequence Database (date of download: 2020-04-02).
Immunohistochemistry
Immunohistochemical analyses for CD3, CD4, CD8, CD15, CD20, CD68, CD123, CD163, PD-1, MPO, p53, Ki67, C3d and C5b-9 were performed on all lung tissue blocks used in this study. Antibodies, staining protocols and conditions are detailed in Supplementary Table 4.
Qualitative and semiquantitative assessment of histopathological lung damage and neutrophilic infiltration
Hematoxylin and eosin (H&E) and Elastica van Gieson (EvG) stained sections of all lung tissues used in this study were independently evaluated by two experienced and board certified pathologists (VZ and KDM) (Supplementary Table 5). Both pathologists evaluated the presence of diffuse alveolar damage (DAD), and if present, its stage, intra-alveolar edema and hemorrhage. The characteristic three phases or stages of DAD - exudative (1), proliferative / organizing (2), fibrotic (3) - were assessed as described39. In our cohort of COVID-19 lungs, we observed only DAD stages 1 and/or 2, and the fibrotic phase (3) was not observed. In addition, both pathologists evaluated the severity of histopathological changes in COVID-19 lungs (1 = mild / discrete alterations, 2 = moderate, 3 = severe changes) based on resemblance between normal and pathologically altered lung tissues. Parameters that were taken into account included reduction of alveolar air-filled spaces, typical histologic features of DAD with hyaline membrane formation, infiltration of lymphocytes, monocytes and neutrophils into interstitial and alveolar spaces, type 2 pneumocyte hyperplasia, desquamation of pneumocytes, histologic features of organizing pneumonia including intra-alveolar fibrin deposition and fibrosis (acute fibrinous and organizing pneumonia, AFOP)40,41. The number of neutrophils per lung tissue section was estimated on H&E stained sections and by immunohistochemical stains for CD15 and MPO using a three tiered system (1 = few or no neutrophils, 2 = moderate number of neutrophils, 3 = high number of neutrophils). Assessment of the two pathologists was concordant in the vast majority of cases. Discrepant cases were reviewed by a third pathologist (NW) to reach consent.
Digital image analysis
Slides were digitalized on a 3DHistech™ P1000 slide scanner at 400x magnification (3DHISTECH Ltd. Budapest, Hungary). Digital slide review and quality control was performed by a board-certified pathologist (VHK). Tissue regions with staining artefacts, folds or other technical artefacts were excluded from analysis. A deep neural network (DNN) algorithm (Simoyan and Zisserman VGG, HALO AI™ on HALO™ 3.0.311.167, Indica Labs, Corrales, NM) was trained using pathologist annotations to automatically localize and measure the area of each lung tissue sample on the digital slides. Background regions and glass were excluded from analysis. Mark-up images for tissue classification were generated and classification accuracy was confirmed through pathology review. For cell-level analysis, color deconvolution for DAB, AP and hematoxylin channels was performed and nuclear segmentation was optimized using cell-morphometric parameters. Marker-positive cells in stromal and epithelial regions were quantified. For CD3, CD4, CD8, CD20, CD68, CD123, CD163 and PD1, staining detection was optimized for the cytoplasmic / membranous compartment and marker expression was measured on a continuous scale at single cell resolution. For assessment of CD8/PD1 double stains, color deconvolution was optimized for separation of DAB (PD1) and AP (CD8) staining products. Internal controls (non-immune cells) and external controls (tonsil) were used to calibrate the detection limits and cross-validated by visual review. For each tissue sample, the total area of lung tissue in mm2, the absolute number of marker-positive cells, cell morphometric parameters and staining intensity were recorded.
Identification of SARS-CoV-2 immune response pattern
Gene expression analysis
Samples were included in the study based on quality of libraries and alignment performance. Applied inclusion criteria are: >1 million of mapped reads, good concentration of libraries, average read length >100bp, > 300 target genes with more than 10 reads. One sample with > 1 Mio reads was excluded from the study because of shorter read length and a low library concentration. Notably this sample had the longest time between death and autopsy (72h) before analysis. Differential expression analysis was performed using the edgeR package comparing normal lung samples, COVID-19 samples and samples from patients with other infections. Genes were selected for downstream analyses by fdr <0.05 and |logFC| >1 for clustering analysis. Clustering analysis was performed using k-means algorithm and complete linkage. Ideal number of clusters (n=3) was chosen based on 30 different algorithms42 and the final clustering derives from the consensus of 2000 iterations. Expression of gene signatures was calculated as median of log2(cpm + 1) of selected genes.
Functional enrichment analysis
Biological processes enrichment was performed using the enrichGO function of the package clusterProfiler43 setting all the genes included in the assay as universe.
Statistical analysis
All the analyses and graphical representations were performed using the R statistical environment software44 and the following packages: ggplot245, circlize46, ComplexHeatmap47, ggfortify48, reshape249 and factoextra50. Correlation between transcripts and viral counts was performed using Pearson’s correlation. Association between continuous and categorical data were tested using Wilcoxon rank sum test.
Based on the often non-uniform histopathological appearance of lung samples from the same patient, transcriptomic, morphologic or histopathological analyses were performed at the tissue sample level. Analyses involving patients’ clinical or demographical data were performed at the patient level, and patients in which all analyzed lung samples expressed an ISGhigh or an ISGlow profile were called ISGhigh or ISGlow patients.
Box-plots elements indicate the median (center line), upper and lower quartiles (box limits) and show all the data points. Whiskers extend to the most extreme value included in 1.5x interquartile range.
Data Availability
The datasets generated and analysed during this study can be accessed in GEO (GEO Submission (GSE151764) [NCBI tracking system #20999475]) and are available from the corresponding author upon request.
Data availability
The datasets generated and analysed within this study can be accessed in GEO (GSE151764) and are available from the corresponding author upon request.
Author contributions
RN, YC, VHK, FD, TJ and KDM jointly conceived the study, performed data interpretation and prepared the manuscript. AT, MB, HM, MT and CA provided intellectual input, provided critical resources and critically reviewed the manuscript. RN, YC, VHK, MH, TH and TJ performed bioinformatic and statistical analysis. AT, JDH, TM, NS, AF collected autopsy specimens, patient data and performed experiments. VZ, NW, WK and CA performed histomorphological evaluation. All authors approved the final manuscript.
Competing interests
VHK has served as an invited speaker on behalf of Indica Labs. TH and TJ are employees of Novartis. The other authors declare no competing interests.
Acknowledgements
VHK gratefully acknowledges funding by the Promedica Foundation (F-87701-41-01). AT, JDH, TM and KDM are supported by and gratefully acknowledge the Botnar Research Centre for Child Health. We would like to thank Christian Tosch, Beat Béni, Daniel Turek, Melanie Sachs, Anne Graber, Christina Herz, Arbeneshe Berisha, Norbert Wey and the USZ pathology IT team, André Fitsche, Marcel Glönkler, Christiane Mittman and the USZ pathology laboratory team for expert technical support and scanning of slides.