key: cord-0811551-y9tu4g18
authors: Zarinsefat, Arya; Hartoularos, George; Yee, Chun J.; Sarwal, Minnie M.
title: Single-cell RNA sequencing of Tocilizumab-treated peripheral blood mononuclear cells as an in vitro model of inflammation
date: 2020-09-11
journal: bioRxiv
DOI: 10.1101/2020.09.11.281782
sha: c8380bd88b374c675a75e1de4381b012f7d8ff41
doc_id: 811551
cord_uid: y9tu4g18

COVID-19 has posed a significant threat to global health. Early data has revealed that IL-6, a key regulatory cytokine, plays an important role in the cytokine storm of COVID-19. Multiple trials are therefore looking at the effects of Tocilizumab, an IL-6 receptor antibody that inhibits IL-6 activity, on treatment of COVID-19, with promising findings. As part of a clinical trial looking at the effects of Tocilizumab treatment on kidney transplant recipients with subclinical rejection, we performed single-cell RNA sequencing of comparing stimulated PBMCs before and after Tocilizumab treatment. We leveraged this data to create an in vitro cytokine storm model, to better understand the effects of Tocilizumab in the presence of inflammation. Tocilizumab-treated cells had reduced expression of inflammatory-mediated genes and biologic pathways, particularly amongst monocytes. These results support the hypothesis that Tocilizumab may hinder the cytokine storm of COVID-19, through a demonstration of biologic impact at the single-cell level.

COVID-19 has posed a significant threat to global health. Early data has revealed that IL-6, a key 26 regulatory cytokine, plays an important role in the cytokine storm of COVID-19. Multiple trials 27 are therefore looking at the effects of Tocilizumab, an IL-6 receptor antibody that inhibits IL-6 28 activity, on treatment of COVID-19, with promising findings. As part of a clinical trial looking 29 at the effects of Tocilizumab treatment on kidney transplant recipients with subclinical 30 rejection, we performed single-cell RNA sequencing of comparing stimulated PBMCs before 31 and after Tocilizumab treatment. We leveraged this data to create an in vitro cytokine storm 32 model, to better understand the effects of Tocilizumab in the presence of inflammation. infections are not severe 1-3 , the overall global burden of the disease has been significant with 51 up to nearly 20% mortality in certain geographic/demographic groups 4,5 . While notable 52 progress has been made in the understanding the virology and disease process, the abrupt 53 onset and lack of effective vaccination has made treatment of COVID-19 difficult 6,7 . 54 55 Interleukin (IL)-6 is a key regulatory cytokine for the innate and adaptive immune response and 56 is a growth factor for B cell proliferation and differentiation, an inducer of antibody production, 57 and a regulator of CD4+ T cell differentiation 8, 9 . Early data from the COVID-19 outbreak has 58

shown that the complications from the disease are partly due to increases in various cytokines, 59 including IL-6 10-13 , and that elevated IL-6 levels may be associated with worse outcomes 13-15 . 60 Tocilizumab is an IL-6 receptor antibody, which binds to both the membrane-bound and 61 soluble forms of the IL-6 receptor (IL-6R), thereby inhibiting the action of the 62 cytokine/receptor complex and interfering with the cytokine's effects 16 . It is a well-studied 63 and accepted therapy for rheumatoid arthritis 17-19 , and has also been studied in giant cell 64 arteritis 20 and organ transplantation 9,21,22 . As such, multiple global investigators are currently 65 undertaking clinical trials to further assess the efficacy of Tocilizumab in the treatment of 66 COVID-19 and its complications (ClinicalTrials.gov). Thus far, it has been shown that COVID-19 67 patient plasma inhibits the expression of HLA-DR which may be partially restored by 68

Tocilizumab treatment, and that treatment with Tocilizumab may also improve lymphopenia 69 associated with COVID-19 23 . Preliminary data for Tocilizumab treatment on COVID-19 70 outcomes has shown improvement in clinical outcomes 24 . While the clinical effects of 71

Tocilizumab in inflammatory and autoimmune disease has been well-studied, there is a 72 paucity of data on the mechanistic/biologic impact of the drug on our immune system. give the final annotated clusters ( Figure 1B) . 104 105 Feature plots showing the expression of "cytokine storm" 43 related pro-inflammatory genes are 106 cell-type specific, with predominance for expression in T cell and monocyte clusters ( Figure 1C) . 107

Although many genes are known to be involved in the cytokine storm of COVID-19 37,38 , we 108 demonstrate that some of the key pro-inflammatory genes (cytokines, interferons, and tumor 109 necrosis factor) are also noted as part of the inflammatory profile in control (no Tocilizumab) 110 patients ( Figure 1C , control cells). Overall, stimulated PBMCs not exposed to Tociluzimab show 111 a dominant signal for T cell activation. After 6 months of treatment with Tociluzimab there is a 112 shift in peripheral blood subset frequencies observed across no treatment (control) vs. In addition to the effect of Tociluzimab on T cells, we also observed an unexpected polarization 178 of monocytes after Tocilizumab treatment ( Figure 1E) . Notably, the Tocilizumab monocyte 179 cluster was enriched for CD14, suggestive of an increased presence of classical monocytes 47 , 180

All Cells

while CD16/FCGR3A expression was more evenly expressed between the two clusters ( Figure  181 3A). We then performed cell trajectory analysis of these monocytes for Tociluzimab treatment 182 effect, utilizing Monocle. This revealed six distinct cell trajectory branches, with two of the 183 branches containing nearly all control cells not exposed to Tocilizumab, and the other four 184 branches containing nearly all Tocilizumab-exposed PBMCs ( Figure 3B) , supporting the 185 presence of unique PBMC trajectories after patient exposure to IL6-R blockade. We utilized 186

Monocle's BEAM function to perform branched expression analysis modeling of the distinct cell 187 trajectory branches for Tociluzimab-exposed PBMCs (circled branch, Figure 3B ), which showed 188 distinct clusters of cells based on treatment status ( Figure 3C) . Tocilizumab for transplant rejection recipients and the first scRNA-seq analysis for such a study. 217

We show a separation of cell clustering based on treatment status, reduced enrichment of 218 inflammatory pathways in Tocilizumab patients, and relatively reduced expression of IL-6R 219 pathway genes in Tocilizumab-treated cells. As would be expected, we did not observe any 220 differences in IL-6 gene expression between control and Tocilizumab cells (as Tocilizumab is an 221 IL-6R blocker), but rather only effects on the subsequent function of that cytokine's pathways. 222

We also show an enrichment of CD14 expression (associated with classical monocytes) in 223

Tocilizumab-treated monocytes, which are believed to be phagocytic, but with reduced 224 inflammatory attributes 47 . This is consistent with our PA described above that shows 225 enrichment of inflammatory pathways in control cells, but not Tocilizumab-treated cells 226 To assign cells to donors of origin in our multiplexed design, we leveraged the genetic 302 demultiplexing tools demuxlet 30 and freemuxlet, both a part of the popscle suite of population 303 genetics tools (https://github.com/statgen/popscle). These tools leverage the genetic 304 polymorphisms present in transcripts to assign the cells found in each droplet to their donor of 305 origin. Demuxlet uses the genotype calls from a genotyping SNP array to classify cells in 306 droplets according to their donor of origin, while freemuxlet "learns" the genotypes of a pre-307 defined number of donors from the transcripts themselves, and assigns the droplets to a 308 respective anonymous donor according to those learned genotypes. Upon first receiving 309 sequencing data, demuxlet was run with input genotypes from all the patients in the cohort. 310

While demuxlet was able to assign most droplets to donors of origin, it revealed that two 311 patients in the genotyping SNP array appeared to have identical genotypes (likely due to human 312 error) and that cells from some patients seemed to drop out (likely due to low viability cells or 313 inaccurate cell counting or mixing). Therefore, to validate demuxlet results, freemuxlet was run 314 using an independent list of SNP sites: exonic SNPs with a minor allele frequency > 0.05 as 315 observed in the 1000 Genomes Project. In order to leverage the droplets across multiple 316 microfluidic reactions, which may enable higher confidence in the learned genotypes, we merged the BAMs from multiple experiments containing the same patients into a single BAM 318 and input this merged BAM into freemuxlet. The droplet assignments from the anonymous 319 donors output by freemuxlet were then compared to those from demuxlet, showing very high 320 concordance. Moreover, comparing the VCF generated from freemuxlet (using the SNPs present 321 in the droplets) to the VCF generated from the SNP genotyping array yielded a 1:1 322 correspondence of anonymous individuals to patients, barring those few problematic patients. 323

Through comparing VCFs and the presence/absence of individuals in each multiplexed 324 experiment, we were able to definitively assign a detected genotype to all detected individuals. 325

Droplet barcodes were then filtered to remove heterotypic droplets containing cells from 326 multiple individuals, and the remaining homotypic droplets were analyzed downstream. 327 328

Raw FASTQ files were processed using CellRanger (v 3.0.1) to map reads against human genome 330 38 as a reference, filter out unexpressed genes, and count barcodes and unique molecular 331 identifiers (UMIs). Subsequent analyses were conducted with Seurat (v 3.1.2) 31 in R (v 3.6.2). 332

We compared PBMCs from all anti-CD3/CD28 stimulated cells from the study baseline, to 333 unstimulated Tocilizumab-treated cells from 3 to 6 months post-treatment with Tocilizumab. 334

Utilizing Seurat, we first filtered cells to only keep those that had less than 10% mitochondrial 335 genes and cells with numbers of features greater than 200 and less than 2,500. Cells were 336 assigned patient identification based on the demuxlet/freemuxlet output described above, and 337 once patients were identified, additional treatment/stimulation/time metadata could be 338 applied. Given that our experiment was divided over 2 days given the high number of samples/cells, we applied Seurat's SCTransform function for data integration to account for any 340 possible batch effects from experiment days 32,33 . Once the data was integrated, we continued 341 downstream data processing. We first determined the principal components (PCA), then 342 constructed a shared nearest neighbor graph (SNN), identified clusters with a resolution of 343 0.75, and finally visualized the cells using uniform manifold approximate and projection 344 (UMAP), per the typical Seurat workflow 31 . Clustering was achieved by using 15 components 345 from the PCA dimensionality reduction. 346

To identify cluster-specific markers following the creation of UMAP plots, we utilized 348 normalized RNA counts of all clusters, scaled the data, and performed differential gene 

A familial cluster of pneumonia associated with the 407 2019 novel coronavirus indicating person-to-person transmission: a study of a family 408 cluster

Clinical features of patients infected with 2019 novel 410 coronavirus in Wuhan

Clinical Characteristics of 138 Hospitalized Patients With

Novel Coronavirus-Infected Pneumonia in Wuhan, China

Case-Fatality Rate and Characteristics of Patients Dying 414 in Relation to COVID-19 in Italy

Characteristics of and Important Lessons From the Coronavirus 416 Disease 2019 (COVID-19) Outbreak in China: Summary of a Report of 72314 Cases From 417 the Chinese Center for Disease Control and Prevention

Preliminary Identification of Potential Vaccine 419 Targets for the COVID-19 Coronavirus (SARS-CoV-2) Based on SARS-CoV Immunological

Use of antiviral drugs to reduce COVID-19 transmission

IL-6 as a keystone cytokine in health and disease

Interleukin-6, A Cytokine Critical to Mediation of 426

Autoimmunity and Allograft Rejection: Therapeutic Implications of IL-6

Exuberant elevation of IP-10

CoV-2 infection is associated with disease severity and fatal outcome. medRxiv

In the eye of the COVID-19 cytokine storm

COVID-19: consider 433 cytokine storm syndromes and immunosuppression

Detectable serum SARS-CoV-2 viral load (RNAaemia) is closely 436 correlated with drastically elevated interleukin 6 (IL-6) level in critically ill COVID-19 437 patients

Risk factors for severity and mortality in adult COVID-19 inpatients 439 in Wuhan

Relationships among lymphocyte subsets, cytokines, and the 441 pulmonary inflammation index in coronavirus (COVID-19) infected patients

Current concepts in the diagnosis and management 444 of cytokine release syndrome

Interleukin-6 receptor inhibition with tocilizumab and attainment 446 of disease remission in rheumatoid arthritis: the role of acute-phase reactants

43-52 systematic review

Risk of adverse events including 451 serious infections in rheumatoid arthritis patients treated with tocilizumab: a systematic 452 literature review and meta-analysis of randomized controlled trials

Tocilizumab for Giant Cell Arteritis

Impact of Tocilizumab (Anti-IL-6R) Treatment on

Antibody-mediated Rejection

Assessment of Tocilizumab (Anti-Interleukin-6 Receptor

Monoclonal) as a Potential Treatment for Chronic Antibody-Mediated Rejection

Transplant Glomerulopathy in HLA-Sensitized Renal Allograft Recipients

tissue and activation signatures in health and disease

Single-Cell Transcriptomics of Regulatory T Cells

A single-cell map for the transcriptomic signatures of 475 peripheral blood mononuclear cells in end-stage renal disease

Single-cell transcriptomics of blood reveals a natural killer cell 478 subset depletion in tuberculosis

Multiplexed droplet single-cell RNA-sequencing 480 using natural genetic variation

Integrating single-cell transcriptomic 482 data across different conditions, technologies, and species

Comprehensive Integration of Single-Cell Data

Normalization and variance stabilization of single-cell RNA-seq 487 data using regularized negative binomial regression. bioRxiv

Reactome pathway analysis: a high-489 performance in-memory approach

clusterProfiler: an R package for comparing biological 491 themes among gene clusters

Single-cell mRNA quantification and differential analysis with Census

The pathogenesis and treatment of the `Cytokine Storm

COVID-19 severity correlates with airway epithelium-497 immune cell interactions identified by single-cell analysis

Human CD56bright NK Cells: An Update

Regulation of CD4 T cell memory by OX40 (CD134)

Defining Memory CD8 T Cell

Utility of flow cytometric detection of CD69 expression as a rapid 504 method for determining poly-and oligoclonal lymphocyte activation

Into the eye of the 507 cytokine storm

SOCS1 regulates 509 the IFN but not NFkappaB pathway in TLR-stimulated human monocytes and 510 macrophages

The IL-6/gp130/STAT3 signaling axis: recent 512 advances towards specific inhibition

JAK3/STAT3 oncogenic pathway and PRDM1 expression stratify 514 clinicopathologic features of extranodal NK/T-cell lymphoma, nasal type

Non-Classical monocytes display inflammatory features: Validation in Sepsis and 518

Capturing Differential Allele-Level Expression 520 and Genotypes of All Classical HLA Loci and Haplotypes by a New Capture RNA-Seq

Signaling in Mouse Hematopoietic Stem and Progenitor Cells Impacts the Ability of the

Molecular heterogeneity in acute renal allograft 526 rejection identified by DNA microarray profiling

Pathological inflammation in patients with COVID-19: a key role for 528 monocytes and macrophages

Pathogenic T-cells and inflammatory monocytes incite 530 inflammatory storms in severe COVID-19 patients

Immune cell profiling of COVID-19 patients in the recovery 532 stage by single-cell sequencing

The authors thank the many individuals without whose enthusiastic participation and help this 385 study would never have been accomplished. We would like to acknowledge the following: TA 386Sigdel who contributed to the study design, JA Liberto who contributed to study design and cell 387 culture/isolation, P Rashmi who contributed to the study design, and AA Da Silva who 388 contributed to cell culture/isolation. A Zarinsefat is funded by the NIH: 5 T32 AI 125222.