Cytokinopathy with aberrant cytotoxic lymphocytes and profibrotic myeloid response in SARS-CoV-2 mRNA vaccine–associated myocarditis

Our clinical cohort consists of 23 patients with vaccine-associated myocarditis and/or pericarditis. The cohort was predominately male (87%) with an average age of 16.9 ± 2.2 years (ranging from 13 to 21 years), in congruence with prior epidemiological reports (24). Patients had largely noncontributory past medical histories and were generally healthy before vaccination. Most patients had symptom onset 1 to 4 days after the second dose of the BNT162b2 mRNA vaccine (Fig. 1A and tables S1 and S2). Six patients either first experienced symptoms after a delay of >7 days after vaccination (P18, P20, P22, and P23) or were incidentally positive for SARS-CoV-2 by polymerase chain reaction (PCR) testing upon hospital admission (P19 and P21) (fig. S1A); these six patients were thus excluded from further analyses, although they potentially reflect the breadth of clinical presentations of vaccine-associated myopericarditis. Our remaining cohort of patients showed no evidence of recent prior SARS-CoV-2 infection, with antibodies to spike (S) protein but not to nucleocapsid (N) protein and negative nasopharyngeal swab reverse transcription quantitative PCR at hospital admission.

Symptom presentation was consistent with acute myocarditis and/or pericarditis, including chest pain, palpitations, fever, shortness of breath, headaches, myalgia, diaphoresis, fatigue, nausea/emesis, and/or congestion (table S1). Laboratory testing in most patients during hospital admission demonstrated elevations in the maximum observed levels of troponin, C-reactive protein (CRP), and B-type natriuretic peptide (BNP) (Fig. 1, B to D and table S1), indicating acute systemic inflammation with direct myocardial injury. Some patients also had leukocytosis (Fig. 1E; for complete blood count differential, see fig. S1, C to K) and elevated neutrophil-to-lymphocyte ratio (NLR) (Fig. 1F), consistent with systemic inflammation. Levels of eosinophils and basophils were within normal limits (fig. S1, F and G).

Electrocardiogram (ECG) performed at admission revealed abnormal findings ranging from nonspecific abnormalities to ST segment elevations and PR segment depressions consistent with pericarditis (table S1). Echocardiogram findings demonstrated borderline low or abnormally reduced left ventricle ejection fraction in several patients, consistent with the diagnosis of myocarditis and/or pericarditis. Most patients further underwent cardiac magnetic resonance (CMR) imaging during acute hospitalization, revealing imaging abnormalities consistent with acute or subacute myocarditis and/or pericarditis (table S1 and fig. S1B).

Patients were treated with nonsteroidal anti-inflammatory drugs, and some also received steroids and intravenous immunoglobulin (Fig. 1A and fig. S1A). Patients were discharged home with improving symptoms and clinical laboratory findings after 1 to 6 days. Although patients showed rapid resolution of clinical symptoms with improved laboratory findings, most of them maintained some imaging abnormalities. These included late gadolinium enhancement (LGE) at longitudinal clinical follow-up at least 2 months after hospital discharge, as assessed by CMR (table S1).

To first assess whether patients with myopericarditis generate overexuberant humoral responses to vaccination, we performed enzyme-linked immunosorbent assays (ELISAs) for SARS-CoV-2–specific S, S1 subunit, and receptor binding domain (RBD) antibodies. Compared with healthy vaccinated controls (VCs, n = 16), including those closer in age to patients where there is an increased risk for vaccine-associated myopericarditis [younger VCs (YVC), n = 6], patients with myopericarditis (n = 9) had no evidence of enhanced anti–SARS-CoV-2 antibodies (fig. S2, A to C). Given the rather moderately reduced levels of SARS-CoV-2–specific antibodies observed in the patients, we next investigated whether the development of myopericarditis was associated with blunted generation of neutralizing antibody responses in these patients. We performed plaque reduction neutralization tests (PRNT₅₀) and found a reduction in the estimated 50% serum inhibitory concentration (IC₅₀) between patients and healthy VCs with no prior history of SARS-CoV-2 exposure (Fig. 2A), possibly owing to the differences in sample collection time relative to vaccination (3 to 11 days for patients versus 7 days for controls). These results suggest that patients with myopericarditis do not display enhanced SARS-CoV-2–specific and neutralizing antibody responses at the time of acute disease presentation but rather show comparable or potentially blunted responses compared to healthy VCs.

To further interrogate whether patients with myopericarditis generate self-targeted autoantibodies after vaccination, we used rapid extracellular antigen profiling (REAP) to screen patients’ plasma for autoantibodies, as has been previously validated (43, 44). In total, autoantibodies recognizing 6183 different human extracellular or secreted proteins/epitopes were assessed by REAP. Of these, we focused on 526 antigens/epitopes relevant to cardiac function, defined either as members of the Gene Ontology (GO) terms “circulatory system process” or “heart contraction” or as “heart tissue enriched/enhanced” in the Human Protein Atlas. Our analysis showed no elevation in the frequency of targeting autoantibodies in patients with myopericarditis compared to YVCs (Fig. 2B and fig. S2D). REAP applied to a separate control cohort of hospitalized patients with COVID-19 (n = 60) showed positive signals for several of these cardiac-related autoantibodies, demonstrating the ability to detect such autoantibodies when present (fig. S2E).

Having found no evidence to support hypersensitivity myocarditis or antibody-mediated processes, we next sought to take a broad and unbiased approach to uncover immune changes in patients with myopericarditis during acute disease. To this end, we first performed serum protein profiling in patients (n = 9) and healthy YVCs (n = 6). Analysis of circulating inflammatory mediators in myopericarditis revealed elevations [false discovery rate (FDR) < 0.05] of the interleukin-1 (IL-1) family cytokines (IL-1β and IL-1RA) (Fig. 3A), consistent with earlier clinical findings suggesting acute, systemic inflammation. In addition, the cytokine IL-15, which can activate NK and T cells, was elevated in myopericarditis. Increases in lymphocyte- and monocyte-directed chemokines were also observed (CCL4, CXCL1, and CXCL10) (Fig. 3A). Given the inflammatory profile and cardiac injury observed in patients with myopericarditis, we also explored levels of matrix metalloproteases and corresponding counter-regulatory proteins and found that they were also elevated [matrix metalloproteinase 1 (MMP1), MMP8, MMP9, and tissue inhibitor of matrix metalloproteinase 1 (TIMP1)]. We next performed exploratory principal components analysis (PCA) involving all measured proteins/cytokines (84 total from Fig. 3A and fig. S3) to identify whether distinct disease immune phenotypes exist among the patients and compared with healthy YVCs. Groups separated well [principal component 1 (PC1): 47.44% of variance], indicating clear differences in immune phenotypes and suggesting ongoing cytokinopathy in patients with myopericarditis (Fig. 3B and fig. S4, A to D). Patients additionally appeared to separate into two groups (PC2: 21.50% of variance), where one was characterized by MMPs, including MMP1 and MMP8 (Fig. 3B and fig. S4, A and B). Upon CMR imaging, patients in this group (P1, P3, P6, and P9) displayed LGE at presentation, whereas the patient in the other group (P2) had no LGE (table S1).

Complementing our serum analysis, we performed single-cell RNA sequencing (scRNA-seq) on peripheral blood mononuclear cells (PBMCs) of the four patients (P1 to P4) from whom suitable samples could be collected during acute illness to better define and correlate the immune effectors in vaccine-associated myopericarditis. Because inflammatory effects from vaccination may be transient, we integrated these samples with ones from healthy young male donors matched for time after vaccination [early-young vaccinated controls (E-YVCs), n = 4]. These controls received the third (booster) dose of mRNA vaccine, which has been reported to induce greater inflammatory side effects compared with the second dose (17, 45), thus controlling for the benign inflammation early after vaccination in young adult males to help separate the likely pathological signatures in myopericarditis. We additionally included samples from three of the patients at follow-up/recovery (R-P1 191 days, R-P3 195 days, and R-P4 190 days after vaccination), four samples from pediatric male healthy donor (HD) controls, and four samples from multisystem inflammatory syndrome in children (MIS-C) patients after SARS-CoV-2 were analyzed using cellular indexing of transcriptomes and epitopes by sequencing [CITE-seq; (46)] to enable refined annotation of cell subsets using surface protein markers (Fig. 3C and tables S2 and S3). After quality control and processing, we obtained a total of 221,727 cells segregating into 44 distinct cell clusters, including unique subsets of myeloid, NK, B, and T cells (Fig. 3D and fig. S5, A to D). We observed the presence of plasma B cells and dividing plasmablasts in patients with myopericarditis (fig. S6, A and B), which was further validated using flow cytometry (fig. S6C). To further evaluate whether these B cells had features of autoantigen-reactive responses, we performed B cell receptor (BCR) repertoire analysis. We observed no significant reduction in clonal diversity or increase in the proportion of immunoglobulin G (IgG) clones harboring mutated BCR regions (fig. S6, D and E), providing little evidence for clonal expansion or somatic hypermutation and supporting findings from our earlier neutralizing antibody and REAP analyses. Differences were observed in the isotype distribution and proportion of cells expressing IgG1 or IgG3 between myopericarditis and healthy E-YVCs, although these did not reach statistical significance (fig. S6, F to I). Overall, these findings potentially suggest nonspecific expansion of plasma B cells and plasmablasts because of the broad cytokinopathy and/or secondary responses to systemic inflammation in these patients.

Given the cytokinopathy in myopericarditis and specifically the elevation of IL-15 observed in our serum analysis, we first investigated changes in NK cell subsets. The elevation of IL-15 was validated in our scRNA-seq cohort using ELISA (fig. S7A). In addition, we noted increased proportions of a CD16⁺ NK cell subset (P = 0.087 after correction for multiple comparisons) in patients (fig. S7B). This trend was not mirrored in innate-like T cells, including Vδ2 γδ T and mucosal-associated invariant T (MAIT) cells (fig. S7, C and D), possibly indicating unique effects in NK cells. To further investigate whether this NK subset is up-regulating genes downstream of the activating cytokine IL-15, we made use of a published gene set of the IL2RB pathway (GSEA Molecular Signatures Database (MSigDB) M8615) and saw up-regulation of various genes (FDR < 0.05), including both the IL2RG and IL2RB subunits of the IL-15 receptor (fig. S7E). To better define the transcriptional signature in this subset, we used differential gene expression in comparison with the rest of the CD16⁺ NK cells. Top up-regulated genes included NK receptors, including the activating receptor NKG2D (or KLRK1) that recognizes stress-induced cell surface ligands, cytotoxicity-related genes (GZMA and PRF1), degranulation and activation markers (LAMP1, LAMP2, and CD69), as well as various integrins and chemokines, among others (fig. S7F). Top down-regulated genes notably included various alarmins (S100A4, S100A6, S100A10, and S100A11) and interferon (IFN)–stimulated genes (ISG15, ISG20, LY6E, IFITM1, IFITM3, and IFITM2), which may be associated with subset-specific dysregulated NK cell signaling, because this gene expression pattern was not reflected across other cell types. Similarly, consistent with the known internalization of NKG2D after ligand engagement, which acts as potential negative feedback to desensitize NK cell responses (47, 48), we observed significantly reduced levels of NKG2D protein on the surface of NK cells in patients by flow cytometry (fig. S7G). One hypothesis is that these parallel analyses of transcripts and surface protein for NKG2D potentially reflect recent NK cell cytotoxic engagement, consistent with their capacity to mediate tissue damage. Overall, we describe a signature of NK cell dysregulation and activation, which may be linked to the development of vaccine-associated myopericarditis.

Next, to inspect changes in the proportions of T cell subsets across the groups while accounting for compositional dependencies between the subsets, we used the Bayesian model scCODA (myopericarditis versus E-YVCs, FDR < 0.05) (49). We observed significantly decreased proportions of a CD4⁺ naïve T cell subset paired with robust expansions of cytotoxic T lymphocytes (CTLs), both CD4⁺ and CD8⁺, and proliferating T cells in myopericarditis (Fig. 4, A to D). Upon examining the transcriptional signature of these expanded subsets, we noted the unique expression of various activation markers, including programmed cell death protein-1 (PD-1) (PDCD1) and CD38/HLA-DR (Fig. 4E), which was further validated by flow cytometry in CD8⁺ T cells (Fig. 4F). These CTLs and proliferating T cells further expressed the chemokine receptors CXCR3 and CCR5 (Fig. 4E), whose ligands—CXCL10 and CCL4, respectively—were shown to be significantly elevated in our earlier serum analysis of the patient cohort (Fig. 3A). These chemokines are known to play central roles in recruitment of T cells to heart tissue, leading to cardiac infiltration of proinflammatory CTLs, and they have previously been described in the context of various inflammatory cardiovascular diseases, including classic myocarditis (50–58). Moreover, these cell subsets express perforin (PRF1) and a diverse set of granzymes (including GZMA, GZMB, GZMH, and GZMK), equipping them to be tissue damaging (Fig. 4E). Further, CD8⁺ T cells in patients showed significantly increased phosphorylated signal transducers and activators of transcription 3 (STAT3) (pSTAT3), indicating ongoing cytokine sensing, compared with E-YVCs by flow cytometry (Fig. 4G). In contrast, published studies of healthy vaccinated donors showed maximum elevation in pSTAT3 1 day after mRNA secondary vaccination (59). To test the hypothesis that these cells are similar to heart CTLs, we trained a logistic regression model using published single-cell data of heart T cells (60) to predict concordant cell subsets in our dataset on the basis of gene expression. Although multiple cell subsets mapped to heart CD4⁺ effector-memory T (CD4⁺ T_tem) cells (suggesting that these cells are most transcriptionally similar to blood circulating T cells) and others did not map well to any single heart population (likely due to tissue adaptability or disease-specific signatures), CD8⁺ CTLs mapped with high frequency to the corresponding heart subset (Fig. 4H), potentially suggesting similarities beyond the cytotoxicity profile.

In addition, to assess the possibility of clonal expansion in these dysregulated T cell subsets, we analyzed the T cell receptor (TCR) repertoire of patients with myopericarditis. An assessment of repertoire diversity, as measured by richness and evenness (Shannon index/richness), revealed that non-naïve CD4⁺ and CD8⁺ T cells were highly diverse across cohorts (fig. S8A). Although CTLs and proliferating T cells were expanded in patients (Fig. 4, B to D), we did not see significant changes in the Vβ usage and found little evidence of highly expanded clones in these cells (fig. S8, B and C), although the possibility of autoreactive cells cannot be fully excluded. Of note, patients with myopericarditis did not display an expansion of TRBV11-2⁺ pathogenic T cells as seen in SARS-CoV-2–associated MIS-C (fig. S8D) (61–64). In summary, T cells in myopericarditis bear an activated signature despite having comparable diversity to healthy VCs, consistent with a cytokine-dependent activation of CTLs after vaccination.

With the extensive inflammatory profile observed in myopericarditis and signs of cardiac injury and LGE on CMR imaging, we lastly investigated changes in the innate myeloid compartment. Compositional analysis (scCODA, FDR < 0.05) revealed a significant decrease in CD14^dim CD16⁺ nonclassical monocytes, which are commonly anti-inflammatory (65, 66), paired with an increase in inflammatory CD14⁺ CD16⁻ classical monocytes (Fig. 5, A and B). These classical monocytes further showed increased expression of genes from the S100A family of alarmins in myopericarditis (Fig. 5C), supporting their role in the inflammatory response and consistent with the elevations we saw in IL-1β (39, 67, 68). Such classical monocytes can further differentiate into tissue macrophages and contribute to chronic disease (65). To better define this possibility in myopericarditis, we used a published dataset of 238 genes (GSEA MSigDB M3468) encoding enzymes/proteins and regulators involved in extracellular matrix remodeling. We observed an up-regulation of this remodeling and profibrotic signature in classical monocytes of patients with myopericarditis (Fig. 5D), consistent with our earlier serum protein/enzyme analysis in these patients (Fig. 3A). In addition, we performed differential gene expression analysis (FDR < 0.05, logFC >0.1), revealing the up-regulation of various specific genes shared across patients (Fig. 5E). These genes included CD163, which marks tissue-resident macrophages (69) and has been associated with a profibrotic phenotype (70). Positive CD163 immunostaining was previously reported on myocardial inflammatory infiltrates in myocarditis after COVID-19 (71, 72) and adenovirus-vectored (Ad26.COV2.S) SARS-CoV-2 vaccination (73). These monocytes also expressed CCR2 (Fig. 5E), which enables migration to sites of tissue injury and differentiation into inflammatory cardiac macrophages, including in myocarditis (66, 74–76). The abundance of CCR2⁺ macrophages has been associated with cardiac remodeling and fibrosis (74, 77). These results are in agreement with the elevations in various MMPs shown by our earlier serum analysis of the patient cohort (Fig. 3A). Other up-regulated genes included ones involved in procollagen processing (ADAMTS2), adhesion and migration (ITGAM, also known as MAC-1 or CD11b), as well as likely inflammation and oxidative stress (CLEC4E, RNASE2, FKBP5, GSR, MARC1, and SULT1A1) (Fig. 5E). Upon monocyte/macrophage activation, CD163 is shed from the cell surface, and its serum levels can be further associated with fibrosis (78, 79). In our patients, soluble CD163 (sCD163) was confirmed to be significantly elevated in serum (Fig. 5F). Together, our findings potentially provide a mechanistic link to the LGE findings in the majority of patients in our cohort (Fig. 5G), which can persist for months after vaccination (Fig. 5H). These results may further suggest the need for careful long-term monitoring of patients with myocarditis and/or pericarditis exhibiting signs of cardiac fibrosis. In summary, our results point to the effector immune populations in vaccine-associated myopericarditis, consistent with published biopsy reports demonstrating a predominate macrophage and lymphocytic immune infiltrate of affected heart tissue (30, 34–37, 80).

The primary aim of this study was to investigate the systemic immunopathological landscape underlying rare cases of myopericarditis occurring after SARS-CoV-2 LNP-mRNA vaccination. A total cohort of 23 patients developing myocarditis and/or pericarditis with elevated troponin, BNP, and CRP levels as well as cardiac imaging abnormalities shortly after vaccination was compared with healthy VCs, among other groups. Different subsets of these patients were profiled using multimodal approaches, including SARS-CoV-2 antibodies and neutralization, high-throughput exoproteome cardiac autoantibody profiling, serum proteomics, flow cytometry, ELISA, and peripheral blood single-cell transcriptome and repertoire sequencing to define immunopathological alterations. Last, clinical follow-up including CMR imaging was conducted between 2 to 9 months after vaccination to assess potential long-term disease effects.

Human participants in this study provided informed consent to use their samples for research and to publish de-identified data in accordance with Helsinki principles for enrollment in research protocols that were approved by the Institutional Review Board (IRB) of Yale University (protocols 2000028924 and 1605017838). For a subset of patients, whose inclusion was solely for medical record review for population surveillance and vaccine investigations, the IRB approved the study and waived the requirement for informed consent (protocol 2000028093). All relevant ethical regulations for work with human participants were followed.

For the antibody, neutralization, REAP, and cytokine assays, whole blood was collected in heparinized CPT blood vacutainers (BDAM362780, BD) and processed on the same day of collection. Plasma samples were collected after centrifugation of whole blood at 600g for 20 min at room temperature. Undiluted plasma was transferred to 1.8-ml Eppendorf polypropylene tubes and stored at −80°C for subsequent analysis. For blood clinical parameters, significance was determined through comparison against normal reference ranges used at Yale New Haven Hospital for each test.

For the scRNA-seq analysis, PBMCs were first isolated by Ficoll-Paque Plus (GE Healthcare) or Lymphoprep (STEMCELL Technologies) density gradient centrifugation. Cells were subsequently washed twice in phosphate-buffered saline (PBS) and resuspended in complete RPMI 1640 (cRPMI) medium (Lonza) containing 10% fetal bovine serum (FBS), 2 mM glutamine, and penicillin and streptomycin (100 U/ml each; Invitrogen). PBMCs were then stored at 10⁶ cells per ml in 10% dimethyl sulfoxide in FBS and stored in −80°C overnight before further storage in liquid nitrogen. Serum was isolated by centrifugation of serum tubes and saving the supernatant in aliquots, which were flash-frozen in liquid nitrogen before cryopreservation in −80°C. All patient PBMCs were processed and cryopreserved within 1 to 6 hours of blood draw. HD PBMCs were all processed within 24 hours to account for overnight shipping.

ELISA was performed as previously described (141). In short, Triton X-100 and ribonuclease A were added to serum samples at final concentrations of 0.5% and 0.5 mg/ml, respectively, and incubated at room temperature for 30 min before use to reduce risk from any potential virus in serum. MaxiSorp plates (96 wells; 442404, Thermo Fisher Scientific) were coated with recombinant SARS-CoV-2 S total (50 μl per well; 100 μg; SPN-C52H9, ACROBiosystems), S1 (100 μg; S1N-C52H3, ACROBiosystems), and RBD (100 μg; SPD-C52H3, ACROBiosystems) at a concentration of 2 μg/ml in PBS and were incubated overnight at 4°C. The coating buffer was removed, and plates were incubated for 1 hour at room temperature with 200 μl of blocking solution [PBS with 0.1% Tween 20 (PBS-T) and 3% milk powder]. Plasma was diluted serially at 1:100, 1:200, 1:400, and 1:800 in dilution solution (PBS-T and 1% milk powder), and 100 μl of diluted serum were added for 2 hours at room temperature. Human anti-spike (SARS-CoV-2 human anti-spike clone AM006415; 91351, Active Motif) and anti-nucleocapsid (SARS-CoV-2 human anti-nucleocapsid clone 1A6; MA5-35941, Invitrogen) were serially diluted to generate a standard curve. Plates were washed three times with PBS-T, and 50 μl of horseradish peroxidase anti-human IgG antibody (1:5000; A00166, GenScript) diluted in dilution solution were added to each well. After 1 hour of incubation at room temperature, plates were washed six times with PBS-T. Plates were developed with 100 μl of 3,3’,5,5'-tetramethylbenzidine (TMB) substrate reagent set (555214, BD Biosciences), and the reaction was stopped after 5 min by the addition of 2 N sulfuric acid. Plates were then read at wavelengths of 450 nm and 570 nm.

Vero E6 kidney epithelial cells (C1008, CRL-1586) were cultured in Dulbecco’s modified Eagle’s medium (DMEM) supplemented with 1% sodium pyruvate, nonessential amino acids, and 5% FBS at 37°C and 5% CO₂. The cell line was obtained from the American Type Culture Collection and was negative for contamination with Mycoplasma using commercially available kits. SARS-CoV-2 (ancestral strain, D614G) USA-WA1/2020 was obtained from BEI Resources (NR-52281) and minimally amplified in Vero E6 cells to generate working stocks of virus. All experiments were performed in a Biosafety Level 3 facility with approval from the Yale Environmental Health and Safety office.

Patient and VC sera were isolated as above and heat-treated for 30 min at 56°C. Plasma was serially diluted from 1:10 to 1:2430 and incubated with SARS-CoV-2 (ancestral strain, D614G) for 1 hour at 37°C. Inoculum was subsequently incubated with Vero E6 cells in a six-well plate for 1 hour for adsorption. Next, infected cells were overlayed with DMEM supplemented with NaHCO₃, 2% FBS, and 0.6% Avicel mixture. Plaques were resolved at 40 hours after infection by fixing in 4% formaldehyde for 1 hour followed by staining in 0.5% crystal violet. All experiments were performed in parallel with negative control sera with an established viral concentration sufficient to generate 60 to 120 plaques in control wells.

Antibody purification, yeast adsorption and library selections, NGS library preparation, and REAP score calculation were performed as previously described (43, 44, 86). A more detailed protocol is also included in Supplementary Methods.

Plasma was isolated from patients as described above. Levels of cytokines were assessed using commercially available vendors, described at length previously (142). Briefly, sera were shipped to Eve Technologies (Calgary, Alberta, Canada) on dry ice, and levels of cytokines and chemokines were measured using the Human Cytokine Array/Chemokine Array 71-403 Plex Panel (HD71) and the Human Matrix Metalloprotease Panel (HDMYO-3-12). All samples were measured upon first thaw. Samples below the range of detection were replaced with the lowest observed value for quantitation.

PCA and hierarchical clustering were performed using MATLAB 2020b (MathWorks Inc.) using standard built-in functions. Gene enrichment was performed on identified clusters using GO analysis to identify significance, and the top 10 hierarchical, biological processes were reported (143, 144).

Serum was isolated as above alongside PBMC isolation, flash-frozen, and stored at −80°C until further use. Serum was then thawed at 37°C and plated in duplicate or triplicate. Standards were made fresh for each assay. The assays were run using the ELISA MAX Deluxe Set Human IL-15 (BioLegend, catalog no. 435104) and the CD163 Human ELISA (Thermo Fisher Scientific, catalog no. EHCD163) kits following the manufacturer’s instructions.

Cryopreserved PBMCs were thawed in a water bath at 37°C for ~2 min and removed from the water bath when a tiny ice crystal remained. Cells were transferred to a 15-ml conical tube of prewarmed growth medium. The cryovial was rinsed with additional growth medium (10% FBS in DMEM) to recover leftover cells, and the rinse medium was added to the 15-ml conical tube while gently shaking the tube.

Thawed PBMCs were centrifuged at 400g for 8 min at room temperature, and the supernatant was removed without disrupting the cell pellet. The pellet was resuspended in 1× PBS with 2% FBS, and cells were filtered with a 30-μM cell strainer. The cellular concentration was adjusted to 1000 cells per μl on the basis of the cell count, and cells were immediately loaded onto the 10x Chromium Next GEM Chip G according to the manufacturer’s user guide [Chromium Next GEM SingleCell V(D)J Reagent Kits v1.1]. We aimed to obtain a yield of ~10,000 cells per lane.

cDNA libraries for gene expression and TCR/BCR sequencing were generated according to the manufacturer’s instructions (Chromium Next GEM SingleCell V(D)J Reagent Kits v1.1). Each library was then sequenced on an Illumina NovaSeq 6000 platform. The single-cell transcriptome sequencing data were processed using CellRanger v5.0.1 and aligned to the GRCh38 reference genome (145). CITE-seq surface protein reads were quantified using a provided tagged antibodies dictionary.

After data alignment and quantification, we performed quality control assessment before proceeding with data processing and downstream analysis. First, genes expressed in fewer than five cells and cells with fewer than 200 genes or more than 10% mitochondrial gene fraction were removed.

The filtered single-cell data were thereafter integrated, and batch effects were corrected using single-cell variational inference (scVI) (146) with a generative model of 64 latent variables and 500 iterations. More specifically, scVI’s negative binomial model was used on raw counts, selecting 5000 highly variable genes identified by the tool’s native method in the combined datasets to produce latent variables. To minimize the dependence of clustering on cell cycle effects, previously defined cell cycle phase–specific genes in the Seurat package (147) were excluded from this list of highly variable genes. Data from patients with acute myopericarditis, E-YVCs, recovered patients, and HDs were integrated for downstream analysis. Four MIS-C CITE-seq samples were also included to aid with cell annotation using surface protein markers.

The latent representation generated by scVI was used to compute the neighborhood graph (scanpy.pp.neighbors), which was used for Louvain clustering (scanpy.tl.louvain) and Uniform Manifold Approximation and Projection (UMAP) visualization (scanpy.tl.umap).

Before clustering, to aid with cell annotation, Single-Cell Remover of Doublets (Scrublet) (148) was used to preliminary label doublets by computing a doublet score for each cell. In particular, a Student’s t test (P < 0.01) and Bonferroni correction were used within fine-grained subclustering of each cluster initially identified by the Louvain algorithm.

Data processing was performed using the Scanpy toolkit (149) following published recommended standard practices. Before downstream analysis, data were normalized (scanpy.pp.normalize_per_cell, scaling factor = 10⁴) and log-transformed (scanpy.pp.log1p). For heatmaps of gene expression, the data were further scaled (scanpy.pp.scale, max value = 10).

First, a logistic regression model was used to transfer preliminary cell annotations from our previous published MIS-C dataset (63) using highly variable genes (scanpy.pp.highly_variable_genes) shared between the two datasets to guide cluster identification. Thereafter, clustering was performed using the Louvain algorithm with an initial resolution of three, and the resulting clusters were defined on the basis of the expression of published cell-specific gene markers. For ambiguous clusters, differential gene expression analysis compared with all other clusters was performed (scanpy.tl.rank_genes_groups), and the resulting top genes were used for unbiased determination of cell identity.

The initial doublet predications by Scrublet were revised here, and doublets were confirmed on the basis of cluster coexpression of heterogeneous lineage gene markers (e.g., CD3, CD19, and CD14) along with the n_counts and n_genes distributions. In addition to removing doublets, clusters displaying patterns of dying cells with high mitochondrial content were further removed. Refined subclusters were used for downstream data analysis and visualization, including cell type frequency, differential gene expression, and repertoire analyses.

To determine significant changes in the proportions of identified subsets across the groups while accounting for compositional dependencies between the subsets in the scRNA-seq data, we used the Bayesian model scCODA (49). Specifically, for each subset, the myopericarditis group was compared with the E-YVC group, and default parameters were used. Changes classified by scCODA as credible after correction for multiple comparisons (FDR < 0.05) were considered statistically significant. In some cases where changes were not classified as significant by scCODA’s method, we used the nonparametric unpaired two-sided Wilcoxon rank-sum test with Benjamini-Hochberg FDR correction for multiple comparisons as previously described (150), and we report the exact adjusted P value, given that such changes may still be biologically meaningful.

The limma package (151) was used to obtain differentially expressed genes (FDR < 0.05). The analysis was performed either comparing different cell subsets or in the same cell subset comparing cells from different groups. Top differentially expressed genes by logFC (up-regulated and down-regulated) and published pathway gene sets were used for plotting and heatmap visualization after data scaling.

Published gene sets from the MSigDB were used to calculate the average enrichment scores across cells of defined subsets as specified in the figure legends using the scanpy.tl.score_genes function with default parameters.

For mapping single-cell immune populations between datasets, a logistic regression model was used using scikit-learn with default parameters. For heart cells, T lymphocyte subsets from published single-cell data (60) were used to train the model using expression data of highly variable genes (scanpy.pp.highly_variable_genes) shared between the two datasets. The model was then applied to predict matching cell subsets in our dataset (mean prediction probability = 0.9), and the percentage of cells mapping between subsets was visualized.

TCR diversity, gene usage, and clonotype analyses were done as previously described (63). A more detailed protocol is also included in Supplementary Methods.

Analysis of the BCR repertoire was performed as previously described (63) using the Immcantation framework as listed above in the TCR analysis, and annotations were incorporated from overall PBMC clustering for B cell subsets.

Patient and donor PBMCs were thawed as above and rested in cRPMI for 2 hours at 37°C before staining. The cells were counted, resuspended in PBS, and stained with e780 Fixable Viability Dye (eBioscience, catalog no. 65-0865-14) for 15 min on ice. The cells were then washed and resuspended in PBS with 10% FBS, with Human TruStain FcX (BioLegend, catalog no. 422302) and True-Stain Monocyte Blocker (BioLegend, catalog no. 426102) for 15 min on ice. The cells were then stained with surface antibody cocktails for 25 min on ice. For phospho-staining, cells were first treated with Fc blocker and then fixed and permeabilized with methanol before staining for intracellular proteins. All flow cytometry data were acquired on the LSR Fortessa cytometer, and gating strategies are shown in figs. S9 and S10. Further details about flow cytometry antibodies are included in Supplementary Methods.

Data were analyzed by various statistical tests as described in detail in the corresponding individual figure captions. For all comparisons, error bars, exact P values, and testing levels were shown whenever possible, and corrections for multiple comparisons were done when appropriate.