Study design

All study procedures involving human participants were conducted in accordance with all relevant ethical regulations and the principles of the Declaration of Helsinki. The study protocol was approved by the institutional review boards at Fred Hutchinson Cancer Center, the National Cancer Institute, New York University Langone Medical Center, Johns Hopkins University, Yale University, Mount Sinai School of Medicine, the University of California, San Francisco (Zuckerberg and Parnassus campuses), Louisiana State University Health Sciences Center, the University of Alabama at Birmingham, Roswell Park Cancer Center and the University of Maryland. Written informed consent was obtained from all participants prior to enrollment and sample collection. Participant confidentiality was maintained throughout the study by deidentifying all data and specimens. The research at Case Western Reserve University and Emory University was performed in compliance with institutional, national and international guidelines governing research on human subjects and biospecimen use. Participants in this study were part of the CITN-12 trial (NCT02595866), a multicenter, open-label, non-randomized phase 1 study of participants (n = 30) living with HIV and advanced cancer^27,30. We analyzed samples from a cohort of 30 virally suppressed individuals with HIV and advanced cancer, encompassing 11 distinct cancer types. These participants, who were stratified by their baseline CD4⁺ T cell counts, demonstrated clinical benefit from anti-PD-1 therapy, including objective cancer responses and evidence of HIV reservoir modulation. (Supplementary Data 1, D1 and D2). Participant demographics reflected the underlying composition of the CITN-12 parent trial, in which 30 of 34 participants were male and four were female. The present multiomic substudy was performed on the first 30 participants for whom matched biospecimens were available, comprising 29 males and one female. Sex was self-reported for all participants. Gender identity data were not separately collected. The predominance of male participants reflects the demographics of PLWH and cancer enrolled in US-based trials, which largely recruit from populations with higher HIV prevalence among men who have sex with men. No sex-stratified analyses were performed due to the limited number of female participants. Informed consent was obtained for publication of indirect identifiers and full demographic details^27,30. Additionally, in the context of HIV outcomes, Uldrick et al. previously showed that a reduction in cells harboring HIV DNA, transcription and elevated plasma viremia was observed using the first week of therapy³⁰. The participants received a pembrolizumab regimen (200 mg administered intravenously every 3 weeks for up to 2 years), with a minimum of two treatment cycles. Participants were stratified into three cohorts based on their baseline CD4⁺ T cell counts: cohort 1 (100−199 CD4⁺ T cells per microliter, n = 6), cohort 2 (200−350 CD4⁺ T cells per microliter, n = 12) and cohort 3 (>350 CD4⁺ T cells per microliter, n = 12) (Fig. 1a,b). These participants presented with 11 different cancers, including AIDS-associated cancers such as Kaposi sarcoma and non-Hodgkin lymphoma (Fig. 1a,b and Extended Data Fig. 1a). The distribution of cancer subtypes was not significantly associated with the cohort based on CD4⁺ T cell numbers (Extended Data Fig. 1a). Notably, complete cancer remission was observed in one participant, and partial remission occurred in four participants (Supplementary Data 1, D1). Demographic data, clinical data, cancer profiling data, response to therapy data and all omic data collected are included in Supplementary Data 1, D1.

Sample collection and processing

Peripheral blood was collected longitudinally from study participants at the following timepoints: pretreatment (C01D01), 24 hours after treatment (C01D02), 1 week after treatment (C01D08) and at EOT (Supplementary Data 1, D3). Blood was collected in acid-citrate dextrose (ACD) tubes (BD Biosciences), and PBMCs were isolated using Ficoll-Hypaque-based density centrifugation. The virological and immunological assays run at each timepoint per PLWH (plasma and PBMCs) are summarized in Supplementary Data 1, D3. It should be noted that bulk RNA-seq of PBMCs and plasma cytokine analyses were prioritized during the study, whereas single-cell and flow cytometry-based readouts were assessed on limited samples (depending on sample availability). In brief, fresh blood was diluted in HBSS (without calcium or magnesium), and PBMCs were isolated using Ficoll-Hypaque density gradient centrifugation. PBMCs were resuspended in freezing media (RPMI media with 10% dimethylsulfoxide (DMSO) and 12.5% human serum albumin) and frozen in aliquots of 5–20 million cells per cryovial. Cells were stored in temperature-monitored LN₂ until use. Plasma was isolated from peripheral blood using a high-speed centrifugation procedure⁶³. Prior to use, PBMCs were gradually thawed in a 37 °C water bath and resuspended in prewarmed complete RPMI (RPMI + 10% FBS + 1% penicillin−streptomycin) and then centrifuged at 400g for 7 minutes at room temperature. Cells were then resuspended in complete RPMI and benzonase nuclease (1 µl of benzonase/1 ml of media) for 10 minutes at 37 °C to avoid cell clumping. Cells were then counted using trypan blue dye and a Countess automated cell counter (Invitrogen), washed once more and resuspended in complete RPMI media at distinct concentrations based on desired downstream assay.

Quantification of plasma cytokines

Custom MesoScale Discovery Human U-Plex assays were used to assess the following cytokines in participant plasma: IL-17A, IL-1β, IL-2, IL-4, IL-6, IL-8, IL-9, IP-10, VEGF, TNF, TGFβ1, TGFβ2, TGFβ3, IFN-α2a, IFNβ, IL-10, IL-15, IL-16, IL-18, IL-22 and IL-7. Plasma cytokines were assessed at the following timepoints: pretreatment (C01D01), 24 hours after treatment (C01D02) and 1 week after treatment (C01D08). Assays were carried out according to the manufacturer’s instructions using 25 µl of undiluted plasma in duplicate wells. TGFβ1/TGFβ2/TGFβ3 cytokines were assayed on a separate plate, and plasma was pretreated with acid, followed by base neutralization to activate the TGF family proteins. Plates were read on the MESO QuickPlex SQ 120MM (MesoScale Discovery) instrument, and cytokine levels were calculated using a standard curve of known cytokine quantities. Prior to univariate analyses, exploratory principal component analysis (PCA) and hierarchical clustering were performed to confirm the inverse correlation between TGFβ levels and inflammatory cytokines, as shown in Supplementary Fig 2a,b.

Quantification of HIV RNA and DNA

Plasma HIV viral load was quantified in study participants⁶³. In brief, a quantitative real-time polymerase chain reaction (qPCR)-based single-copy assay (HMMCgag) was used to target the noncoding 5′ region of the HIV gag gene. Cell-associated HIV DNA (LTR-gag) was measured as described^64,65. Single-copy plasma HIV RNA was measured at multiple timepoints, including baseline (C01D01), 24 hours after treatment (C01D02) and 1 week after treatment (C01D08). In addition, plasma viremia was monitored prior to each of the first four anti-PD-1 cycles and, subsequently, after every three treatment cycles³⁰. Longitudinal measurements of HIV RNA and DNA were used to evaluate virological outcomes and compared between participant groups over the course of therapy.

RNA-seq

Bulk RNA-seq was performed on whole blood collected in PAXgene tubes from participants at the following timepoints: C01D01 (n = 29), C01D02 (n = 28) and EOT (n = 14). These extracted RNA samples were stored at a temperature of −80 °C until they were ready for the preparation of mRNA libraries. The quantity and integrity of the RNA were confirmed through microelectrophoresis using the High Sensitivity RNA ScreenTape Kit (Agilent Technologies) on the 2200 TapeStation system (Agilent Technologies). After this, the construction of mRNA libraries was performed according to the Illumina TruSeq V2 library preparation kit. The validation of these libraries occurred through microelectrophoresis on a 2100 Bioanalyzer system (Agilent Technologies), quantification using Kapa Library Quantification Kits (Roche) and eventual pooling and clustering on Illumina TruSeq V2 flow cells. A paired-end 100-cycle run was performed on the Illumina HiSeq 2500 platform, averaging 43 million paired-end reads per sample.

Single-cell RNA-seq was conducted on PBMCs from participants across distinct timepoints: C01D01 (n = 12), C01D02 (n = 8) and C01D08 (n = 10). The sequencing protocol adhered to the recommended guidelines of the 3′ V3 chemistry kit by 10x Genomics. Initially, single-cell gel beads-in-emulsion (GEMs) were generated by combining 350,000–400,000 PBMCs, barcoded Single Cell 3′ v3.1 beads and partitioning oil. These GEMs underwent processing via a 10x Chromium Controller machine (10x Genomics). Subsequent steps encompassed gel bead elimination, cell lysis and reverse transcription in a C1000 Touch Thermal Cycler (Bio-Rad) to yield barcoded cDNA. To extract oil residues, the GEM mixture underwent purification using Dynabeads followed by cDNA amplification. After amplification, cDNA was purified using SPRIselect reagent. Evaluation of cDNA quality and quantity occurred using a 2100 Bioanalyzer (Agilent Technologies). Further stages entailed fragmentation, end repair, A-tailing and adaptor ligation in accordance with the manufacturer’s protocol. Library quality was verified using the Bioanalyzer. These libraries were subsequently sent to the Beijing Genomics Institute and sequenced on a DNBseq-T7 (MGI Tech) machine, targeting 20,000 reads per cell. Data analysis involved uploading FASTQ files to the 10x Genomics cloud platform using Cell Ranger, without applying depth normalization. The filtered count matrix from the analysis was processed through the Seurat package⁶⁶ in R for subsequent analysis. Cell annotation was accomplished using SingleR⁶⁷, a reference expression dataset derived from the MonacoImmuneData atlas⁶² in the celldex R package⁶⁷. To enhance analysis quality, the DoubletFinder package⁶⁸ in R was employed to identify and remove doublet cells. Exploration of differential gene expression was facilitated using the MAST R package⁶⁹. Cell clusters were visually represented using uniform manifold approximation and projection (UMAP). PCA was performed to assess global differences among the clusters.

Bioinformatic analysis

Multivariate model of plasma cytokines and transcriptomic modules associated with treatment and HIV RNA

Feature selection was performed using least absolute shrinkage and selection operator (LASSO)⁷⁰, implemented through the glmnet R package. The LASSO technique was used to determine the combination of plasma cytokines that best associated with treatment 24 hours after treatment compared to baseline and also to identify the combination of transcriptomic modules that best predicted the increase in HIV RNA over time. The model was optimized using leave-one-out cross-validation, and the features with the least cross-validated mean square error (MSE) were identified.

Processing of RNA-seq data

An integrated automated pipeline developed in-house was used to preprocess FASTQ files per sample. The raw sequencing reads were trimmed off any adapter sequence contaminants using Trimmomatic⁷¹. The trimmed reads were aligned to the Ensembl version of the human genome (GRCh38) using the STAR aligner⁷². The gene abundance was estimated by counting the number of reads mapping to the exons unique to the transcripts of a gene using HTSeq⁷³. The obtained gene counts were normalized using trimmed mean of M values (TMM) by correction for the library size⁷⁴. Normalized gene expression data from whole-blood RNA-seq samples collected at pretreatment (C01D01), 24 hours after treatment (C01D02) and EOT were analyzed using PCA. Significant transcriptomic shifts were detected at both C01D02 and EOT relative to baseline (C01D01), as assessed by Wilcoxon test (P 1a,b).

Differential gene expression analysis

The differences in transcriptomic profiles between the timepoints were determined by fitting a generalized linear model (GLM) to each gene. Gene expression was treated as the dependent or response variable and timepoint as the independent variable. A likelihood ratio test was then performed to test whether the coefficients or fold changes differed from zero. The obtained P values were then corrected for multiple comparisons using the Benjamini−Hochberg method, and P values less than 5% were reported as significant. All these tests were performed as part of the edgeR package in R⁷⁵.

Pathway enrichment analysis

To test the enrichment of pathways/gene sets among the differentially expressed genes, a GSEA⁷⁶ pre-ranked with 1,000 permutations was performed using the R package fgsea. The Hallmark gene sets and the Immunologic Signature gene sets (C7)⁷⁷ from the MSigDB, along with transcription factor targets from the ChEA database⁷⁸, were used as gene set databases. The differentially expressed genes were pre-ranked by the decreasing order of their −log₁₀(P) times their signed (log fold change) to identify pathways upregulated and downregulated comparing timepoints. The obtained pathway P values were corrected for multiple comparisons using the Benjamini−Hochberg method. Pathways enriched at a P value of less than 5% were considered significant. Enriched pathways were visualized using the R package enrichmentmap. See Supplementary Data 1, D9 for differentially expressed transcription factor targets. Differentially enriched transcription factor targets are listed in Supplementary Data 1, D10. Notably, although pathway modules did not significantly differ across the three pretreatment CD4⁺ T cell-defined cohorts, they were consistently modulated in all groups at 24 hours after infusion (Supplementary Fig. 1c).

SCimilarity analyses to map the modules to publicly available datasets

To map transcriptional modules from the CITN-12 cohort onto a large, publicly available single-cell atlas, we used the SCimilarity framework. The analysis was conducted using SCimilarity model version 1.1 (cellquery_model_v1.1), which enables efficient cell similarity searches across a reference database of 22.7 million single-cell transcriptomes. The query dataset consisted of scRNA-seq profiles generated from PBMC samples collected in the CITN-12 study at multiple timepoints (C01D01, C01D02 and C01D08). This dataset was formatted for SCimilarity by computing low-dimensional embeddings using the pretrained cellquery model. Transcriptomic modules of interest (as defined in Fig. 1 and Supplementary Data 1, D9 and D10) were scored at the single-cell level, and the resulting module scores were appended to the query scRNA-seq object using the Seurat framework. SCimilarity was used to identify phenotypically similar cells from the precomputed reference database. For each query cell, similarity scores to all database cells were computed. To ensure high-confidence mappings, only the top 1% of most similar cells (that is, 99th percentile matches) were retained for downstream analysis. To focus on biologically relevant contexts, matched cells were filtered to include only those originating from ex vivo PBMC samples derived from healthy individuals, people with cancer or people with infectious diseases. For each module, median similarity scores were calculated per matched sample across the retained database cells. These values were visualized using a row-normalized heatmap to enable comparison of module expression patterns across diverse immune conditions.

Immune cell-type-specific transcriptional signatures module analysis

The significant Immunologic Signature gene sets (of C7 from MSigDB) that were differentially expressed between timepoints were grouped into modules of related gene sets based on a Jaccard coefficient of 0.25—that is, gene sets were grouped into modules if they had a minimum of 25% gene overlap between every pair of gene sets. The minimum module size contained at least three gene sets. The modules were obtained independently for the upregulated and downregulated gene sets (Supplementary Data 1, D9).

SLEA

SLEA³² was used to represent the expression of modules by calculating the z-score of the pathway per sample. The mean expression value of genes enriched in a module was compared to the mean expression of random sets of genes of the same module size for 1,000 permutations for every sample. The difference between the observed and expected mean expression values for each gene set was determined as the SLEA z-score.

Transcriptomic deconvolution

CIBERSORT⁷⁹ was used to infer immune cell subset frequencies from whole-blood transcriptomic data using a reference expression set of distinct hematopoietic cell states that were purified from umbilical cord blood and peripheral blood (Gene Expression Omnibus (GEO) accession ID: GSE24759). Genes of the reference set were assigned to a cell type if they were expressed at a positive fold change at a P value less than 5% in that cell type compared to all the other cell types. This reference expression matrix of genes by cell type was used with CIBERSORT to infer cell subset frequencies at baseline, at 24 hours after treatment and at C01D08 of treatment. The frequencies derived via CIBERSORT were compared to the flow cytometry-based frequencies obtained at C01D01 (that is, pretreatment) and were significantly correlated. The cell subset frequencies obtained via CIBERSORT were further used on the whole-blood RNA-seq data to obtain immune cell subset-specific gene expression. Deconvolution was performed using least squares fit to model the whole-blood gene expression measures across samples as the function of contributions of immune subset-specific gene expression weighted by the corresponding cell frequencies of those subsets obtained via CIBERSORT. Linear regression was applied separately on baseline and 24 hours after treatment. The regression coefficients were used as surrogates for estimated cell subset-specific average gene expression. The difference in the average cell subset-specific gene expression estimates was used as the level of gene expression change per immune cell subset between 24 hours after treatment compared to baseline.

Integration between transcriptomic modules and plasma cytokines

A projection-based approach using mixOmics⁸⁰ was used to integrate the change in plasma cytokine expression with the change in transcriptomic module expression. To minimize the technical effect specific to each data type, a sparse least square regression was used to generate a unique scale for each data type and project all data types on the same scale to assess correlation. Pearsonʼs correlation between the features of the two data types and a P value based on the distribution of the correlation coefficients were calculated. Pearsonʼs correlations below a P value cutoff of 0.05 were considered significant.

In vitro HIV infection assay with ISG-inducing compounds

To validate ex vivo findings that activation of TLR3, TLR7/TLR8 and IL-15 pathways induces antiviral responses and restricts HIV infection, we conducted an in vitro infection assay using PBMCs freshly isolated from five healthy human donors. Cells were cultured in AIM-V medium supplemented with 10% serum replacement (Corning, 355500) and 10 mM HEPES. The experimental design consisted of two arms. In arm 1, memory CD4⁺ T cells were directly isolated from PBMCs using a negative selection enrichment kit (STEMCELL Technologies, 19157). In arm 2, CD8⁺ T cells were first depleted from the PBMCs, and the remaining cell population was cultured for stimulation. In both arms, cells were subjected to one of the following stimulation conditions for 24 hours: unstimulated control, IFNγ (5 ng ml⁻¹; PeproTech, 300-02), IFNβ (5 ng ml⁻¹; PeproTech, AF-300-02B), poly I:C (10 ng ml⁻¹; InvivoGen, tlrl-pic-5), R848 (10 ng ml⁻¹; InvivoGen, tlrl-r848-1) or IL-15 (100 ng ml⁻¹; PeproTech, 200-15). To evaluate whether the antiviral effects were mediated by IFN signaling, each stimulation condition was also tested in the presence of neutralizing antibodies against IFNAR (10 µg ml⁻¹; PBL Assay Science, 21385) and IFNγ (10 µg ml⁻¹; Miltenyi Biotec, 130-095-743). At 16 hours after stimulation, activation of antiviral signaling pathways—including phosphorylation of STAT1, IRF3 and IRF7 as well as expression of IFIT1 and APOBEC3G—was assessed. At 24 hours, memory CD4⁺ T cells were isolated from arm 2, and cells from both arms were then infected with the p89.6 dual-tropic HIV-1 strain by spinoculation in the presence of saquinavir to ensure single-round infection. Viral replication was quantified by measuring intracellular p24 levels at day 4 after infection.

Ex vivo HIV-specific response profiling

HIV-specific PBMC responses were assayed using flow cytometry. After the thawing procedure, PBMCs were resuspended at 1.0 × 10⁷ cells per milliliter, and 100 µl (or 1.0 × 10⁶ cells) was added per well of a 96-well U-bottom plate. Cells were unstimulated, stimulated with 2 µg ml⁻¹ HIV-1 PTE Gag peptide or stimulated with 2 µg ml⁻¹ Staphylococcal enterotoxin B (SEB) in the presence of 10 µg ml⁻¹ brefeldin A (BFA). After 6 hours of stimulation, cells were stained with DAPI LIVE/DEAD stain to discriminate between live and dead cells. After a wash with FACS buffer (PBS + 4% human serum), cells were incubated for 20 minutes at room temperature with a surface antibody staining cocktail mix diluted in Brilliant Stain Buffer (BD Biosciences).

Flow cytometry

Flow cytometry staining was performed using surface, intracellular and phospho-specific protocols depending on the experimental design. Surface staining was carried out at 37 °C for 20 minutes in 25 µl of PBS supplemented with 2% FBS (PBS + 2% FBS). After incubation, cells were washed with 150 µl of PBS + 2% FBS buffer. For intracellular staining, cells were washed and then incubated in fixation/permeabilization solution (BD Biosciences, FoxP3 Transcription Factor Kit) for 30 minutes at 4 °C. After fixation, cells were washed with the accompanying fixation/permeabilization buffer and incubated with intracellular antibody cocktails—prepared in the same buffer—for 45 minutes at 4 °C. For phospho-flow staining, cells were fixed in 50 µl of BD Phosflow Fix Buffer I (cat. no. 557870) for 30 minutes at 4 °C, washed with 150 µl of PBS + 2% FBS and permeabilized using 50 µl of cold BD Phosflow Perm Buffer III (cat. no. 558050) on ice for 15 minutes. Cells were washed twice with 150 µl of 1× permeabilization buffer (Invitrogen, 00-8333), and all residual buffer was carefully removed. Intracellular staining for phospho-proteins was performed by adding antibodies diluted in 25 µl of 1× permeabilization buffer (Invitrogen) and gently mixing. Cells were incubated for 60 minutes at 4 °C and subsequently washed with 150 µl of 1× permeabilization buffer. After the final wash, cells were resuspended in 100 µl of PBS for acquisition. All centrifugation steps were carried out at room temperature at 500g for 5 minutes. Data were acquired using a BD FACSymphony A5 flow cytometer operated via BD FACSDiva software, collecting a minimum of 100,000 live cells per sample. Flow cytometry data were analyzed using FlowJo software (version 10.8.1). High-dimensional flow cytometry data were analyzed using unsupervised methods to visualize and identify distinct immune cell populations. After compensation, quality control and gating to exclude debris, doublets and dead cells, events corresponding to viable single cells were exported for downstream analysis. Live CD3⁺CD4⁺ T cells or live CD3⁺CD8⁺ T cells were selected and exported for clustering analysis. Marker expression values were biexponentially transformed, and data from all samples within each panel were concatenated for joint analysis. Dimensionality reduction was performed using UMAP for two-dimensional visualization of cellular heterogeneity. UMAP projections were used to display marker density and distribution across cell subsets. Unbiased clustering was then performed using the PhenoGraph algorithm (https://github.com/jacoblevine/PhenoGraph) with the number of nearest neighbors (k) set to 60, enabling the identification of phenotypically distinct cellular clusters. Cluster-level marker expression was summarized using median fluorescence intensity (MFI) and visualized via heatmaps in R to compare marker profiles across identified populations. UMAP plots and cluster overlays were generated using FlowJo and R, and representative density plots were used to highlight markers of interest. Antibodies used for HIV-specific profiling are as follows: BUV496 anti-human CD3, BUV563 anti-human CD4, BUV615 anti-human CD45RA, BUV661 anti-human CD27, BUV737 anti-human CD8, BUV805 anti-human CD45RO, BV480 anti-human CD14, BV480 anti-human CD19, BV650 anti-human CD69, BV750 anti-human CCR7, BB515 anti-human CCR5, BB630 anti-human SLAMF6, BB660 anti-human CD101, PE-Cy7 anti-human PD-1, PE-Cy5 anti-human CD95, APC-R700 anti-human CD28, APC-Cy7 anti-human CD39, BUV395 anti-human BCL2, BV421 anti-human TCF7, BV570 anti-human Ki-67, BV605 anti-human IRF4, BV711 anti-human EOMES, BV786 anti-human TBET, BB700 anti-human IFNγ, BB750 anti-human NR4A1, PE anti-human TOX, PE-CF594 anti-human TNF and APC anti-human ID2. Antibodies for in vitro experiments are as follows. The staining panel included: LIVE/DEAD BV510, Life Technologies, L34957; CD3 BUV615, BD Biosciences, 612992, clone UCHT1; CD4 BV605, BioLegend, 317438, clone OKT4; CD45RA BV650, BioLegend, 304136, clone HI100; CD27 AF647, BD Biosciences, 567026, clone O323; CCR7 BUV563, BD Biosciences, 741317, clone 3D12; pIRF3 PerCP, Novus Biologicals, BS-3195R-PERCP, clone polyclonal; pIRF7 BUV496, BD Bioscience, custom, clone K47-671; pSTAT1 (p701) AF488, BD Biosciences, AB_2737715, clone 4a; IFIT1 APC, Novus Biologicals, NBP2-71005APC, clone OTI3G8; APOBEC3G AF700, Novus Biologicals, NBP1-77206AF700, clone polyclonal; and p24 RD1, Beckman Coulter, 6604667, clone KC57. All antibodies were commercially purchased and validated by the vendor, and the details of the clone, fluorophore and dilution used per stain are provided in Supplementary Data 1, D11.

Statistical analyses and reproducibility

Whole-blood, PBMC and plasma samples were collected longitudinally from participants and varied based on availability (see Supplementary Data 1, D3 for sample availability per omic). No replicates were analyzed or added at any timepoint. The samples were collected based on the clinical trial design (see ‘Study design’ section above), and no further statistical tests were done to predetermine the sample size. We were not blinded to the outcomes tested. Samples were run in a single batch for plasma cytokine and single-cell analyses. Whole-blood RNA-seq was done in two batches, and transcript count matrices were corrected using the ComBat package in R. The uncorrected and corrected count matrices will be made available on the GEO upon finalization and before publication of this paper. FlowJo (version 10.10.0), GraphPad Prism (version 10.2.3) and R (version 4.4.0 in RStudio 2024.04.1 Build 748) were used to generate figures and for statistical analyses. The Benjamini−Hochberg test was used to correct for multiple testing during whole-blood RNA-seq data analyses. All figure legends indicate the P values and statistical tests used. The figure panels were assembled using Adobe Illustrator 2024 (version 28.5).

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.