By varying experimental parameters, we found that the activity of Tn5 is a limiting factor for scATAC-seq. We then devised a qPCR assay with lambda DNA as substrate to measure relative Tn5 activities against a standard curve (Supplementary Tables S1, S2). It was used to compare the activity of different concentrations of an in-house Tn5 preparation (Tn5-H) loaded with the adapters listed in Table 1 to Tn5 enzymes from 10× Genomics (Tn5-TXG) as well as the Illumina TDE1 enzyme (Tn5-ILMN) (Figure 1C).

The highest concentration of in-house Tn5 (Tn5-H100, 83 µg/ml or 1.6 µM) used showed a four-fold higher activity per volume in terms of product yield as compared to the Tn5 enzyme provided in the 10x Genomics scATAC-seq version 2 (v2) kits (Tn5-TXGv2) and an about 1.3-fold higher activity per volume than Tn5-TXGv1.1 (Figure 1C, Supplementary Table S2). We also tested different buffer compositions and found that the buffer provided in the 10x Genomics kit performed best. The Tn5-ILMN enzyme activity was similar to that of Tn5-TXGv1.1 (Figure 1C).

We then compared our Tn5-H100 and Tn5-H30 (25 µg/ml or 0.48 µM) preparations to the Tn5-TXGv1.1 enzyme using mouse embryonic fibroblasts (MEFs) treated with interferon β (IFNβ) for 6 hours. We found that Tn5-H100 yielded the best results with respect to different QC parameters (Figure 2A, Supplementary Figure S2A; Supplementary Figures S1A, B, Supplementary Table S3): Recovered cell barcodes showed 1) higher TSS enrichment scores, 2) higher number of unique fragments, and 3) smaller fragments. In the following, we refer to the use of Tn5-H100 with the buffers provided by 10x Genomics as scTurboATAC vs scATAC for the use with Tn5-TXG.

For further analysis, we selected high-quality cell populations of 5,000 to 6,000 cells for the three samples (Figure 2A, red rectangle, Supplementary Table S4). The scTurboATAC protocol yielded an about 2-fold higher number of unique fragments normalized to sequencing depth than cells from the regular scATAC workflow (Figure 2B). Similar to the mouse data, a higher number of unique fragments (absolute and normalized to sequencing depth) was detected in cryopreserved peripheral blood mononuclear cells (PBMC) as a test case when using Tn5-H100 as compared to Tn5-TXGv2 (Figure 2C, Supplementary Figure S1C). Embedding and clustering of both data sets revealed an increased resolution of the scTurboATAC data for both MEF (Figure 2D) and PBMC data (Figure 2E).

Two distinct MEF subtypes (epithelial- and mesenchymal-like MEFs) were annotated with marker genes from our previous work (Muckenhuber et al., 2023) after removing minor clusters 1 and 2 (scATAC) and cluster 1 (scTurboATAC) that showed increased apoptosis module scores (Supplementary Figure S1D). Cluster 4 of scATAC and cluster 2 of scTurboATAC showed increased mesenchymal-like MEF marker module scores (Figure 2D, Supplementary Figure S2A). The epithelial-like MEF markers were enriched in cluster 3 of scATAC and clusters 3 and 4 of scTurboATAC (Figure 2D, Supplementary Figure S2B). A differential analysis of gene activity scores yielded 552 differentially active genes between epithelial clusters C3 and C4 of scTurboATAC (FDR ≤ 0.05 and log2FC ≥ 0.5). The most enriched GO term of these differentially active genes was GO:0061844 “antimicrobial humoral immune response” (Fisher p-value of 3.6e-06), likely arising from the IFNβ treatment. When computing module scores of previously identified interferon stimulated genes (ISGs) in MEFs for clusters C2, C3 and C4 of scTurboATAC, cluster C4 showed the lowest activity which points to a weaker response (Supplementary Figure S2C).

It is noted that differences in sequencing depth can affect the resolution and might contribute to the differences observed. However, when sub-sampling the ATAC data to obtain the same read numbers for each Tn5 preparation, we still obtained higher number of fragments per cell with an additional cluster for scTurboATAC (Supplementary Figures S2D, E, Supplementary Table S5). Furthermore, the scTurboATAC PBMC data also displayed an increased number of clusters (14 instead of 10) as compared to scATAC at comparable read depth. Thus, the scTurboATAC protocol resolves two epithelial-like clusters, while only one cluster was detected with the standard scATAC protocol. The resolution improvements obtained with the scTurboATAC protocol are likely to reflect a higher complexity of the data that arises from the higher number of unique fragments.

scTurboATAC was further evaluated for the analysis of transcription factor (TF) binding and gene regulatory mechanisms in MEFs stimulated with IFNβ for 6 hours. Calling peaks from the pseudo-bulk data yielded a higher number of accessible peaks per cell with the scTurboATAC compared to the scATAC protocol (Supplementary Figures 2F, G). Here, scATAC recovered roughly 15,000 and scTurboATAC 26,000 accessible peaks per cell on average. Furthermore, accessibility footprints at transcription start sites (TSSs) and at STAT1 and CTCF binding motifs showed similar or better enrichment in pseudo-bulk scTurboATAC compared to scATAC (Figure 3A).

Next, we applied our previously described approach for co-accessibility analysis with the RWire software to the MEF data set (Muckenhuber et al., 2023). RWire reveals enhancer-promoter links from correlation between open sites of scATAC-seq data. For epithelial-like MEFs of scTurboATAC, we only analyzed cluster 3 as it represents a homogeneous cell population. We observed twice as many co-accessible links above background co-accessibility and minimal percent accessible cells in epithelial- as well as mesenchymal-like MEFs of scTurboATAC than scATAC (Figure 3B). Co-accessible links from scTurboATAC showed higher correlation coefficients and percent accessible cells than co-accessible links from scATAC (Supplementary Figure S2; Supplementary Figure S2H).

Transcription of interferon stimulated genes (ISGs) is partly regulated by simultaneous STAT1 and STAT2 (STAT1/2) binding to ISG promoters as well as to distal regulatory sites (Muckenhuber et al., 2023). We were able to link roughly 10 % of ISGs without STAT1/2 promoter binding to a distal STAT1/2 binding event by co-accessibility analysis of 6h IFNβ-stimulated epithelial- and mesenchymal-like MEFs from scTurboATAC (Figure 3C). In contrast, only 5–8% of ISGs were linked to a distal STAT1/2 binding event by co-accessibility analysis of 6h IFNβ-stimulated epithelial- as well as mesenchymal-like MEFs from scATAC. For both epithelial- and mesenchymal-like MEFs newly linked promoters by TurboATAC belonged to genes associated with “cellular response to interferon-beta” (p-values 7x10⁻¹¹ and 1.9x10⁻¹² respectively) highlighting the relevance of the identified connections. An example for the improved detection of ISG regulation by distal STAT1/2 binding events in the scTurboATAC data is shown for the Ifit gene cluster in Figure 3D. The 150 kb locus depicted contains six ISGs (Ifit1, Ifit1bl1, Ifit1bl2, Ifit2, Ifit3, and Ifit3b) (Muckenhuber et al., 2023). For Ifit1, Ifit1bl2 and Ifit3b promoters, direct STAT1/2 binding was observed in IFNβ-stimulated MEFs. Three co-accessible links were detected between ISG promoters and distal STAT1/2 binding events in epithelial- and mesenchymal-like MEFs of scATAC. These co-accessible links were also detected by scTurboATAC where the percent accessible cells was higher. In addition, six co-accessible links between ISG promoters and STAT1/2 binding events were only detected with scTurboATAC. In mesenchymal-like MEFs of scATAC, we detected one unspecific link between the Ifit2 promoter and a distal peak without a STAT1/2 binding event or other ISG promoter. This link was not observed in mesenchymal-like MEFs of scTurboATAC.

Next, we compared scTurboATAC and scATAC with the Tn5-TXGv2 enzyme with respect to cell type resolution by clustering and marker gene activity in PBMCs (Figure 4A, Supplementary Table S6). The scTurboATAC protocol yielded a higher number of accessible pseudo-bulk peaks per cell (24,500 on average) compared to the scATAC protocol (19,000 on average, Supplementary Figure S3A). This difference resulted in an improved detection of progenitor cells, classical monocytes/dendritic cells and provided a higher resolution of the B cell cluster (Figure 4A, B). We identified a distinct B cell cluster with scTurboATAC (cluster 1) that was not resolved in scATAC data by integrating B cells from scATAC and scTurboATAC (Supplementary Figure S3B). When investigating the B cell clusters further, we found differential TF activities for PU.1-IRF, OCT2 and ATF3 among others (Figures 4B, C, Supplementary Figure S3C). While OCT2 has been observed in late state B cell differentiation (Corcoran et al., 1993), PU.1-IRF plays a role in earlier lymphoid maturation (Scott et al., 1994). This could account for their contrary activity scores. Elevated ATF3 activity was mainly found in cluster 1 of scTurboATAC data. Since ATF3 modulates metabolism and immunity, especially immune-responsiveness, this population might represent a specific cellular state (Ku and Cheng, 2020). Previous studies suggested that ATF3 is a repressor of IRF targets (Garber et al., 2012; Labzin et al., 2015). In line with this view, the population with high ATF3 activity has low PU.1-IRF scores.

Multi-omics approaches allow to directly link epigenetic heterogeneity with transcriptional profiles. Accordingly, we adapted our scTurboATAC protocol to this type of application and devised a Tn5 loading scheme that is compatible with the 10x Genomics multiome chemistry. The latter utilizes a spacer oligonucleotide to add cell barcodes to the Tn5 fragments possibly by ligation. To account for this difference, we changed the oligonucleotides for loading (Table 1) and additionally used a Protein A-tagged Tn5 for a single-cell C&T-seq assay with RNA information from the same cell (Figure 1B). To avoid self-dimerization for the multiome-ATAC workflow a complementary oligonucleotide with a phosphorothioate (PTO) backbone was used. This modification can be omitted for multiome-C&T due to the Tn5 capture and washing of unbound Tn5. Phosphorylation of Tn5_read1 was found to be essential to facilitate the barcode ligation step in the droplets.

Our protocol yielded a higher number of unique fragments normalized to sequencing depth when using self-loaded Tn5-H100 in comparison to Tn5-TXGv2 with the multiome chemistry in PBMCs (Figure 5A left, Supplementary Table S7) while the number of reads in TSSs was similar (Supplementary Figures S4A, B). UMI counts of transcripts normalized to sequencing depth showed the same yield in comparison to the standard workflow (Figure 5A right, Supplementary Table S7). The number of clusters was higher for scMultiome-ATAC with Tn5-H100 (Figure 5B) compared to the standard scMultiome-ATAC with Tn5-TXGv2 (Figure 5C) for both for the ATAC readout alone as well as co-embedded ATAC and RNA data. When using only the gene expression data, the number of clusters remained the same (Figures 5B, C). Thus, the improved ATAC reaction did not compromise the RNA readout in the scMultiome-ATAC protocol.

The scC&T analysis was performed in isogenic TF1 wildtype and mutant cell lines that carry the IDH2 R140Q mutation (Stein, 2016) with histone H3 trimethylation at lysine 27 (H3K27me3) as C&T target. For both samples the number of H3K27me3 fragments was similar (Figure 5D left, Supplementary Table S8) with a bias for fragments with higher nucleosome numbers (Supplementaty Figure S4C), a high fraction of reads in peaks (Supplementaty Figure S4D) and comparable UMI transcript counts (Figure 5D right, Supplementary Table S8). UMAP of C&T and RNA data from both cell lines showed separate embedding with some overlap according to IDH2 mutation status (Figure 5D far right).

Very broad regions of H3K27me3 enrichment were observed as depicted for an exemplary region of the HOXC gene cluster in Figure 5E. In both cell lines all genes within this H3K27me3-rich region were silenced while adjacent genes lacking H3K27me3 signal were highly expressed.

To evaluate the data quality further we investigated the heterogeneity of gene expression in the TF1 wildtype cell line and its correlation to promoter-associated H3K27me3. The number of RNA UMIs was >2500 on average per cell and provided a good coverage of the single cell transcriptome profile. It revealed heterogeneous expression of the transcription factor KLF6 in TF1 wildtype cells. While the number of unique fragments per cell obtained from the H3K27me3 C&T analysis yields a relatively sparse coverage of the rather broad genome-wide H3K27me3 peaks, the H3K9me3 signal was sufficient to reveal a medium to high negative correlation of KLF6 expression with the H3K27me3 modification at the KLF6 promoter (Figure 5F). Interestingly, the epigenetic deregulation of KLF6 has been linked to AML (Adelman et al., 2019), where recurrent mutations of IDH1/2 are drivers of tumorigenesis and therapeutic target (Tyner et al., 2018). Taken together, these data illustrate the versatility and specificity of the scTurboATAC and scMultiome-C&T techniques, which open up new ways for data analysis and experimental readouts.

Purified unloaded Tn5 was obtained from the EMBL Protein Expression and Purification Core Facility, which was produced and stored as described previously (Hennig et al., 2018) with some adjustments resulting in the following workflow: The pETM11-Sumo3-Tn5 plasmid encoding His₆-Sumo3-tagged Tn5 E54K L372P was freshly transformed into E. coli BL21(DE3) codon + RIL cells (Stratagene). Precultures were grown overnight at 37°C in LB medium supplemented with 30 µg/ml kanamycin and 34 µg/ml chloramphenicol and used to inoculate the large-scale expression cultures. A 10 ml preculture was added to 1 liter of TB-FB medium supplemented with 2 mM MgSO₄, 0.05% glucose, 1.5% lactose, 30 µg/ml kanamycin and 34 µg/ml chloramphenicol. The cultures were grown at 37°C until the OD₆₀₀ reached −0.6, after which the temperature was reduced to at 18°C. After overnight expression at 18°C, the cultures were harvested by centrifugation (30 min, 5000 x g, 4°C) and the pellets were flash-frozen in liquid nitrogen and stored at −80°C until the start of the protein purification.

For the protein purification, all chromatography steps were performed at 4°C and samples were kept on ice at all times. The cell pellet was resuspended in cold running buffer (20 mM Hepes-NaOH pH 7.2, 800 mM NaCl, 5 mM imidazole, 1 mM EDTA, 2 mM DTT and 10% glycerol) supplemented with cOmplete EDTA-free protease inhibitors (Roche). The cells were lysed by 5 passages through a microfluidizer, followed by centrifugation (30 min, 140000 x g, 4°C). After centrifugation, a 10% polyethyleneimine (PEI) solution at pH 7.2 was added dropwise to the cleared lysate to a final concentration of 0.5% while constantly stirring the solution. To remove the PEI-precipitated fraction, the lysate was subjected to another centrifugation step at 20000 x g for 30 minutes at 4°C. The cleared lysate containing the recombinant His₆-Sumo3-tagged Tn5 protein was then loaded on a 5 ml prepacked cOmplete His-tag purification column (Roche) equilibrated in running buffer. After loading the sample, the cOmplete column was washed with the same buffer until the absorbance signal at 280 nm returned to baseline. The column was eluted with a 1-step elution in running buffer supplemented with 300 mM imidazole. The elution fractions were analysed via SDS-PAGE and the fractions containing the His₆-Sumo3-tagged Tn5 were pooled and dialysed overnight to running buffer in the presence of His₆-tagged SenP2 protease to remove the His₆-Sumo3 fusion tag. The next morning, the dialysed sample was loaded again on a 5 ml cOmplete His-tag purification column equilibrated with running buffer. The untagged Tn5 was collected in the flow through of the column and further purified using size exclusion chromatography. The Tn5 sample was loaded on a HiLoad Superdex 200 16/600 pg column (GE Healthcare) equilibrated with 50 mM Tris pH 7.5, 800 mM NaCl, 0.2 mM EDTA, 2 mM DTT and 10% glycerol. The elution fractions were analysed via SDS-PAGE and the dimer peak fractions were pooled. The final Tn5 sample was aliquoted and stored at −20°C in 25 mM Tris pH 7.5, 800 mM NaCl, 0.1 mM EDTA, 1 mM DTT and 50% glycerol at a concentration of 1 mg/ml.

Lyophilized adapter oligonucleotides were resuspended in annealing buffer (50 mM NaCl, 40 mM Tris·HCl, pH 8) at 100 µM concentration. Adapters (Table 1) were preannealed (read 1 with the mosaic end (ME) sequence and read 2 with ME) separately on a thermocycler for 2 min at 85°C followed by cooling down to 20°C with 0.1°C per minute. Glycerol was then added to a final concentration of 50% (v/v) and the adapters were stored at −20°C. For Tn5 assembly, transposase, annealed primer pairs and dilution buffer (DB, 50 mM Tris·HCl, pH 7.5, 100 mM NaCl, 0.1 mM EDTA, 1 mM DTT, 0.1 % NP-40, 50 % glycerol) were mixed with equal volumes (e.g. 2.5 µl each, 10 µl in total) and incubated for 30 min at room temperature. The loaded Tn5 was then diluted with 20 µl DB yielding the Tn5-H100 (83 µg/ml) preparation and stored at −20°C.

Activity was assesed by qPCR measuring the efficiency in generating fragmented lambda DNA. Different Tn5 preparations were used at the same volume with 10 ng of lambda DNA, 2.5 µl tagmentation buffer (20 mM Tris-HCl pH 7.5, 20 mM MgCl2, 50% v/v DMF) or 10× Genomics tagmentation buffer as indicated. For the reaction, 1.25 µl phage DNA (10 ng), and 1.25 µl Tn5 enzyme were mixed and incubated for 10 min at 55°C. The reaction was stopped by adding 1.25 µl 0.2% SDS. The mix was diluted 1:100 in nuclease free water. A 10 µl qPCR mix was used that contained 5 µl 2x KAPA SYBR Fast, 0.2 µl 10 µM i7 adapter (CAA​GCA​GAA​GAC​GGC​ATA​CGA​GAT), 0.2 µl 10 µM i5 adapter (AAT​GAT​ACG​GCG​ACC​ACC​GAG​ATC​TAC​AC), 1 µl diluted template DNA and 3.6 µl nuclease free water. The following program was run on a StepOnePlus (Applied Biosystems) machine (72°C for 5 min, 95°C for 30 s, followed by 25 cycles of 95°C for 10 s, 60°C for 30 s and 72°C for 60 s. Relative activities were calculated from the measured Ct values (see Supplementary Tables S1, S2). Linear performance was verified using a diluted reference standard.

For scTurboATAC, the transposase provided with the ATAC kit from 10 × Genomics was replaced with Tn5-H100. Purified nuclei were processed using the standard ATAC protocol according to the user guides for Chromium Next GEM Single Cell ATAC v1.1 or v2 kits (CG000496, 10 × Genomics) as indicated. MEFs were cultured and treated as previously described (Muckenhuber et al., 2023). PBMCs were isolated as previously described and stored viably frozen in serum containing DMSO (Mallm et al., 2019). Prior to the experiments they were thawed in a water bath and immediately processed. For transposition, 10,000 nuclei were used and loaded on to the Chromium Next GEM Chip H (PN-1000161). Libraries were generated using standard reagents. Quality control and molarity calculations of the final libraries were performed with the Qubit 3.0 Fluorometer (Q33216, Invitrogen) and the 4200 tapestation system (Agilent Technologies). For multiome analysis the appropriately loaded Tn5 was used (Table 1), while all downstream reactions were performed with standard multiome reagents from 10x Genomics.

TF1-IDH2 wildtype and mutant cell lines were acquired form ATCC (Manassas, Virginia, United States, wildtype, CRL-2003; IDH2 mutant, CRL-2003). They were cultured following the providers instructions. Cells were harvested, washed once in PBS and 1 × 10⁶ cells were fixed with 1 ml 0.2% PFA for five min at room temperature. Fixation was stopped with 125 mM glycine for 5 min followed by one wash with PBS. Cells were then resuspended in 500 µl ice cold nuclei buffer and incubated for 10 min on ice. Cells were pelleted and resuspended in 500 µl antibody buffer (2 mM EDTA, 0.1% BSA, 0.04% digitonin, 20 mM HEPES pH 7.5, 150 mM NaCl and 0.5 mM spermidine).

A 10 µl volume of αH3K27me3 (cell signalling, 9733) and 12.5 µl RNase inhibitor (40 U/µl) were added and incubated over night at 4°C. Cells were washed with wash buffer (0.04% digitonin, 20 mM HEPES pH 7.5, 150 mM NaCl and 0.5 mM spermidine), resuspended in 500 µl wash buffer containing 10 µl 1 µg/µl anti-rabbit secondary antibody raised in guinea pig (Sigma-Aldrich, SAB3700889) and 1 U/µl RNase inhibitor. After 30 min at room temperature cells were washed three times with wash buffer with five min incubation time per wash and resuspended in Tn5 buffer (0.01% digitonin, 20 mM HEPES pH 7.5, 150 mM NaCl, 0.5 mM spermidine and 1 × protease inhibitor). Cells were counted and 500,000 cells were used from this point forward. 10 µl of loaded proteinA-Tn5 was added and incubated at room temperature for 30 min. After three washes with wash buffer, cells were resuspended in tagmentation buffer (10 mM Tris pH 7.5, 10 mM MgCl, 5 mM HEPES, 50 mM NaCl, 25% DMF, 0.01% digitonin, 1 U/µl RNase inhibitor), transferred to a PCR tube and incubated for one hour at 37°C. A volume of 15 µl was then used for GEM formation with the multiome kit from 10x Genomics.

Multiplexed library pools were generated at 10 nM concentrations of each library and sequenced on NovaSeq 6000 system (Illumina, v1.5). Sequencing was paired end with Read1N and Read2N sequenced with 50 cycles each, i7 (index) sequenced with 8 cycles and i5 using 10 × Genomics barcodes with 16 cycles for ATAC and 24 cycles for multiome-ATAC and multiome-C&T. Sequencing depth was at least 25,000 read pairs per nucleus. RNA libraries were sequenced paired end with 28 bp and 90 bp for read 1 and read 2, respectively.

Data were demultiplexed and aligned with Cell Ranger ATAC count v2.0.0 (10× Genomics) using the provided mouse mm10 reference. Further processing of data was conducted with ArchR (v1.0.3) in R (v4.2.2) (Granja et al., 2021). Non-cell barcodes were removed by applying a minimal threshold to log10 number of fragments (3.8 for conventional scATAC, 3.9 for scATAC with Tn5-H30, and 4.3 for TurboATAC) and TSS enrichment score (5 for conventional ATAC, 12 for scATAC with Tn5-H30 and TurboATAC). This yielded 5153 cells for conventional scATAC, 5809 cells for scATAC with Tn5-H30, and 5612 cells for TurboATAC. Cell doublets were removed using Amulet (scDblFinder v1.12.0) (Thibodeau et al., 2021) and single cells were embedded in two-dimensional space using an accessibility matrix of 500 bp tiles, Iterative LSI (LSI components 2–10 for conventional ATAC, 2–12 for TurboATAC) and UMAP. Cells were clustered by SNN modularity optimization with a resolution of 0.1. Apoptotic cell clusters were identified by marker gene modules of single-cell gene activity scores (Casp3, Parp1, Bax, Bak1, Bid, Bbc3, Pmaip1, Bcl2, Bcl2l1, Mcl1, Trp53, Cycs, Apaf1, Tnfrsf10b, and Fas; (Galluzzi et al., 2018)) and removed for further analyses.

Sample-specific and merged peak sets were obtained from pseudo-bulk accessibility data using MACS2 (Zhang et al., 2008) in ArchR (extendSummits = 1000; reproducibility = 1). Transcription factor binding motifs were annotated in sample-specific peaks with Homer mm10 motifs from chromVARmotifs (v0.2.099) (Schep et al., 2017). Transcription factor motif footprints were computed using ArchR (no Tn5 bias normalization, 5 bp smooth window). MEF subtype clusters were identified as epithelial- and mesenchymal-like MEFs by marker gene modules of single-cell gene activity scores of epithelial- and mesenchymal-like MEF markers from our previous study (Muckenhuber et al., 2023).

Simultaneously accessible regions in the same cell were determined by co-accessibility analysis with RWire (https://github.com/RippeLab/RWire-IFN) based on our previously described RWire approach (Mallm et al., 2019; Muckenhuber et al., 2023). Cell numbers between the samples and subtypes were matched to 1,617 cells. The consensus peak set was used for all samples and subtypes. Pearson correlation coefficients were calculated for peaks within 1 Mb for each sample and subtype without aggregation of single cells. Only co-accessible links with positive correlation above background co-accessibility threshold r (0.088 for epithelial-like MEFs from conventional ATAC, 0.079 for mesenchymal-like MEFs from conventional ATAC, 0.066 for epithelial-like MEFs from TurboATAC, 0.066 for mesenchymal-like MEFs from TurboATAC) were considered. The background co-accessibility threshold was determined by the 99th co-accessibility percentile from randomly shuffled accessibility values over cells and peaks. Co-accessible links were further filtered by a lower threshold of 5 % accessible cells in the linked peaks. STAT1/2 bound sites and TSSs of interferon-stimulated genes from our previous study (Muckenhuber et al., 2023) at co-accessibly linked peaks were annotated.

Data were demultiplexed and aligned with Cell Ranger ATAC count v2.1.0 (10 × Genomics) using the provided human GRCh38 reference. Further processing of data was conducted with ArchR (v1.0.3) in R (v4.2.2) (Granja et al., 2021). Non-cell barcodes were removed by applying a minimal threshold to log10 number of fragments (4.1 for conventional ATAC, and 4.4 for TurboATAC) and TSS enrichment score (10 for conventional ATAC, 8 for TurboATAC). This yielded 7658 cells for conventional ATAC, and 8309 cells for TurboATAC. Cell doublets were removed using Amulet (scDblFinder v1.12.0) (Thibodeau et al., 2021) and single cells were embedded in two-dimensional space using an accessibility matrix of 500 bp tiles, Iterative LSI (LSI components 2-21) and UMAP. Cells were clustered by SNN modularity optimization with a resolution of 0.2.

Clusters were assigned to cell types by marker gene modules of single-cell gene activity scores of established markers (Supplementary Table S4). Sample-specific and merged peak sets were obtained from pseudo-bulk accessibility data using MACS2 (Zhang et al., 2008) in ArchR (extendSummits = 1000; reproducibility = 1). Transcription factor binding motifs were annotated in sample-specific peaks with Homer hg38 motifs from chromVARmotifs (v0.2.0) (Schep et al., 2017). A subset of only B cells was embedded in two-dimensional space using an accessibility matrix of 500 bp tiles, Iterative LSI (LSI components 2–8) and UMAP. Cells were clustered by SNN modularity optimization with a resolution of 0.25. B cells from conventional ATAC and TurboATAC were integrated using Harmony (v0.1.1) (Korsunsky et al., 2019). ChromVAR deviations were calculated for 200 bp windows around top 1000 transcription factor binding motifs in sample-specific peak sets by decreasing scores (v1.20.0) (Schep et al., 2017). ChromVAR deviations were imputed using MAGIC (van Dijk et al., 2018) in ArchR.

Data were demultiplexed and aligned with Cell Ranger ARC count v2.0.2 (10× Genomics) using the provided human GRCh38 reference. Further processing of data was conducted with ArchR (v1.0.3) (Granja et al., 2021) and Seurat (v4.3.0) (Stuart et al., 2019) in R (v4.2.2). Non-cell barcodes in ATAC were removed by applying a minimal threshold of 3.8 and 1000 to log10 number of fragments and reads in TSS, respectively. This yielded 4563 cells for conventional ATAC, 7405 cells for ATAC Tn5-H50, and 6690 cells for TurboATAC. ATAC cell doublets were removed using Amulet (scDblFinder v1.12.0) (Thibodeau et al., 2021), cells with a blacklist ratio above 0.01 and/or nucleosome ratio above 3 were excluded, and single cells were embedded in two-dimensional space using an accessibility matrix of 500 bp tiles, Iterative LSI (LSI components 2–20) and UMAP. Cells were clustered by SNN modularity optimization with a resolution of 0.3.

Non-cell barcodes in RNA were removed by applying a minimal threshold of 100 and 500 to number of detected genes and number of UMI counts, respectively, as well as a maximal threshold of 40% to percentage of mitochondrial UMI counts. This yielded 4587 cells for RNA of conventional ATAC, 7916 cells for RNA of ATAC Tn5-H50, and 7101 cells for RNA of TurboATAC. RNA cell doublets were removed using Scrublet with a cutoff of 0.15 (v0.2.3) (Wolock et al., 2019) in Python v3.10.4 (pandas v1.4.3, scipy v1.8.1) and single cells were embedded in two-dimensional space using SCTransform (v0.3.5) (Hafemeister and Satija, 2019), PCA (PCs 1–20) and UMAP. Cells were clustered by SNN modularity optimization with a resolution of 0.5. Co-embedding of single cells using ATAC and RNA information was conducted using Signac (v1.9.0) (Stuart et al., 2021). Combining filtered RNA and ATAC cell information yielded 3620 cells for RNA/ATAC of conventional ATAC, 5999 cells for RNA/ATAC of ATAC Tn5-H50, and 5077 cells for RNA/ATAC of TurboATAC. Single cells were embedded using LSI (LSI components 2–25) for ATAC, SCTransform and PCA (PCs 1–30) for RNA, weighted nearest neighbor graph and UMAP. Cells were clustered by SNN modularity optimization with a resolution of 0.8.

Data were demultiplexed and aligned with Cell Ranger ARC count v2.1.0 (10× Genomics) using the provided human GRCh38 reference. Further processing of data was conducted with ArchR (v1.0.3) (Granja et al., 2021) and Seurat (v4.3.0) (Stuart et al., 2019) in R (v4.2.2). Browser tracks of pseudo-bulk RNA and H3K27me3 signal were generated with plotgardener (v1.4.1) (Kramer et al., 2022). For H3K27me3 peak calling, MACS2 was used with the broad setting and keeping one read duplicate resulting in 57,533 peaks that were used in downstream analysis with ArchR. Non-cell barcodes in H3K27me3 were removed by applying a minimal threshold of 500 and 0.1 to number of fragments and fraction of reads in peaks, respectively. This yielded 2786 cells for H3K27me3 of TF1-IDH2wt, and 2653 cells for H3K27me3 of TF1-IDH2mut. Non-cell barcodes in RNA were removed by applying a minimal threshold of 100 and 500 to number of detected genes and number of UMI counts, respectively, as well as a maximal threshold of 25% to percentage of mitochondrial UMI counts. This yielded 2928 cells for RNA of TF1-IDH2wt, and 2956 cells for RNA of TF1-IDH2mut. Co-embedding of single cells using H3K27me3 and RNA information was conducted using Signac (v1.9.0) (Stuart et al., 2021). Combining filtered RNA and H3K27me3 cell information yielded 1803 cells for RNA/H3K27me3 of TF1-IDH2wt, and 1800 cells for RNA/H3K27me3 of TF1-IDH2mut. Single cells were embedded using LSI (LSI components 3–20) for H3K27me3, SCTransform (v0.3.5) (Hafemeister and Satija, 2019) and PCA (PCs 1–20) for RNA, weighted nearest neighbor graph and UMAP. Differential gene expression analysis was performed on log10-transformed and scaled data using Wilcoxon test in Seurat (only.pos = FALSE, min.pct = 0.1, logfc.threshold = 0).