TeaCNV: decoding tumor somatic absolute copy number and clonal architecture from single cell chromatin accessibility data

TeaCNV: decoding tumor somatic absolute copy number and clonal architecture from single cell chromatin accessibility data | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Article TeaCNV: decoding tumor somatic absolute copy number and clonal architecture from single cell chromatin accessibility data Shaojun Zhang, Ying Wang, Yuhao Deng, Xinbao Yin, Yanru Zhang, and 6 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-6609843/v1 This work is licensed under a CC BY 4.0 License Status: Posted Version 1 posted You are reading this latest preprint version Abstract Accurate inference of absolute copy numbers beyond simple gains and losses from single-cell chromatin accessibility (scATAC-seq) data remains challenging, thereby obscuring the distinction between genetic and epigenetically driven oncogenic dependencies. Here, we present TeaCNV, a computational framework that reconstructs clonal absolute copy number profiles and tumor clonal architectures from scATAC-seq data without matched DNA baselines. Validated against bulk whole-genome sequencing in renal cell carcinomas, TeaCNV resolved subclonal absolute copy number profiles with less than 10% error and detected copy number variations (CNVs) with 98.6% accuracy, outperforming existing methods. Applied to six cancer types including renal, breast, pancreatic, head and neck, colorectal, and ovarian cancers, TeaCNV delineated polyclonal architectures and revealed distinct chromatin accessibility patterns driven by CNVs in key driver genes, including AKT2, ZNF217 and SOX2. By enabling absolute copy number profiling and clonal deconvolution from epigenomic assays, TeaCNV bridges critical gaps in studying oncogenic dependencies and genotype-phenotype relationships at single-cell resolution. Biological sciences/Computational biology and bioinformatics/Computational models Biological sciences/Computational biology and bioinformatics/Genome informatics Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Introduction Chromosomal copy number variations (CNVs), which drive ploidy changes in tumor genomes, constitute a hallmark of human cancer and critically shape tumor transcriptional landscapes 1 , 2 . Aneuploidy, defined by the presence of abnormal chromosome numbers, is strongly associated with tumor aggressiveness and poor clinical outcomes 3 , 4 , 5 , 6 . Precise quantification of chromosomal absolute copy numbers (CNs) therefore represents a fundamental requirement for understanding cancer evolution and guiding clinical decision-making 7 , 8 . While bulk and single cell whole-genome sequencing (scWGS) remain gold standards for CN profiling 9 , 10 , 11 , their clinical implementation is hindered by technical complexity, high cost, and limited scalability 12 . Furthermore, scWGS fails to capture tumor microenvironment dynamics. Although single-cell RNA sequencing (scRNA-seq) enables microenvironment characterization and tumor subpopulation analysis, existing CNV inference tools, such as inferCNV 13 , 14 , 15 , 16 , 17 , CopyKAT 18 , Numbat 19 and HoneyBADGER 20 , suffer from inherent limitations including transcriptional noise and post-transcriptional regulation frequently confound DNA level CNV estimates 21 . Single cell ATAC sequencing (scATAC-seq), which profiles chromatin accessibility through transposase-accessible DNA sequencing, provides a more direct readout of genomic content than transcriptome-based approach 22 , 23 . Current scATAC-based CNV detection methods face two limitations: dependence on matched normal DNA data and a focus on relative gain/loss detection rather than absolute CN solution 24 , 25 , 26 . These limit their utility in polyploid tumor analysis. To address these limitations, we developed TeaCNV, a computational framework that reconstructs absolute CN profiles and resolves clonal substructures directly from scATAC-seq data. Validated against bulk WGS in renal cell carcinomas, TeaCNV achieved 98.6% accuracy of CNVs calling and less than 10% error in absolute CN profiling while detecting subclonal CNVs undetectable by conventional sequencing. Application to multi-cancer types (renal, breast, pancreatic, head and neck, ovarian and colorectal cancers) revealed lineage-defining amplifications of oncogenic drivers and polyclonal architectures. Results Methodology to infer clonal absolute copy number profile from scATAC-seq We developed TeaCNV, a quantitative framework for reconstructing clonal absolute CN profiles and clonal substructure from scATAC-seq data without requiring matched DNA sequencing references. The algorithm processes an input peak-by-cell matrix with cells pre-annotated as reference (normal) or inferred (tumor or epithelial) cells. Peaks are sorted by chromosomal coordinates. The outputs include clonal partitions of single cells and corresponding clonal absolute CN states for chromosomal segments (Fig. 1 ). The inference pipeline operates through four steps ( Methods ). First is initial cell clustering. For each inferred cell, relative copy ratios are computed by normalizing peak accessibility signals against the median accessibility in matched reference cells. Dimension reduction via principal component analysis (PCA) on these ratios enables initial clustering of inferred cells into subpopulations (Fig. 1 a). Second is chromosomal breakpoint detection. Cells within each subpopulation are aggregated. Multi-scale chromosome breakpoints are identified using the pruned exact linear time (PELT) algorithm 27 . This generates segmentation boundaries for each subpopulation while accommodating scATAC-seq data sparsity (Fig. 1 b). Third is joint ploidy-state optimization. For each subpopulation, absolute CNs of genomic segments are inferred through a combined optimization, including Nelder-Mead optimization to model expected relative copy ratios and CN state relationships, maximum likelihood estimation weighted by segment length and ratio variance, optimal integer CN determination via minimum Akaike information criterion (AIC, Fig. 1 c). Fourth is clonal architecture refinement. Subpopulations undergo pairwise comparison of their absolute and relative CN profiles. Subpopulations with concordant profiles are merged, followed by breakpoint recalculation and integer CN estimation (Fig. 1 d). Last is the consensus profile generation including clone-specific absolute CN matrices across chromosomal segments and single-cell clonal assignments (Fig. 1 e). Identifying Tumor clonal substructures with TeaCNV in clear cell renal cell carcinoma To evaluate TeaCNV’s ability to reconstruct clonal architecture, we performed high-throughput 10X Genomics sequencing on four clear cell renal cell carcinoma (ccRCC) patients, generating two scATAC&RNA-seq co-assayed datasets and two scATAC-seq-only datasets (Fig. 2 a and b , Supplementary Fig. 1) . Cell type annotation, based on chromatin accessibility signatures of canonical markers ( PTPRC , PECAM1 , CD34 , COL1A2 , and EPCAM ) and subsequent UMAP visualization, identified 8,367 epithelial cells for CNV analysis ( Supplementary Figs. 2 and 3 ). Immune, endothelial and stromal cells were used as references. Using these non-epithelial cells as a reference, we applied TeaCNV to calculate single cell relative copy number ratios and reconstruct clonal architecture, while simultaneously inferring genome-wide absolute CNs (Fig. 2 c-i). TeaCNV demonstrated robust reference compatibility, yielding consistent clonal architectures and integer CN profiles regardless of the reference cell type used (immune vs. endothelial, Supplementary Figs. 4 and 5 ). In both scATAC&RNA-seq co-assayed and scATAC-seq-only datasets, TeaCNV accurately identified established ccRCC driver CNVs, including chr3p loss and chr5q gain 28 (Fig. 2 c-i). TeaCNV revealed distinct complexities in clonal substructures across these four ccRCC samples. While ccRCC1 displayed a monoclonal architecture (Fig. 2 c ) , the other three cases displayed polyclonal expansion with a dominant subclone (Fig. 2 d, f and h ). In ccRCC2, two clones shared similar epithelial scores but diverged in chromatin accessibility patterns driven by different absolute CNs at chr14 and chr16q (Fig. 2 d and e ). TeaCNV resolved three clones in ccRCC3, with losses of chr8p and chr16q defining two minor clones that collectively accounted for 32.59% of the tumor cells (Fig. 2 f and g ). A rare subclone (clone3: 5.19%) was further distinguished from clone 2 by co-amplification of chr8q and chr16p, accompanied by corresponding changes of chromatin accessibility (Fig. 2 g). ccRCC4 contained four subclones stratified by absolute CNs of chr20q and chr22 (Fig. 2 h). The dominant clones (clone 1 and 4: 56.9%) exhibited a 4-copy on chr20q, while a minor clone (clone 3: 10.5%) exhibiting co-occurrence of a 4-copy of chr20q and a 1-copy of chr22 (Fig. 2 i). Performance evaluation of absolute copy number estimation using matched bulk WGS data To validate the performance of absolute CNs inference using TeaCNV, we further analyzed three ccRCC patients (ccRCC3, ccRCC4 and ccRCC1) with sample-matched bulk WGS. The absolute CN profiles from bulk WGS data were considered as the ground truth. In ccRCC3, the integer CN profiles estimated by TeaCNV were highly concordant with the WGS profiles, particularly for the trunk events chr3p loss and chr7 gain (Fig. 3 a and b ). As previously described, the subclonal chr8p and chr16q losses in clone 2 and 3, identified by TeaCNV, were further validated by WGS data, showing absolute CNs between 1 and 2 (Fig. 2 f, g, Fig. 3 a). Notably, the rare clonal events of 1-copy gain on both chr8q and chr16p (clone 3) were undetectable in bulk WGS data due to insufficient cell representation (Fig. 3 b and Supplementary Fig. 6a ). Using the scATAC&RNA-seq co-assayed data of ccRCC3, we applied epiAneufinder 24 and Copy-scAT 25 to chromatin accessibility data and inferCNV 13 to gene expression (Fig. 3 c-e). These approaches either diluted signals from aberrant regions or misidentified neutral regions, leading to ambiguous distribution of inferred CNV scores for genomic segments with distinct CN states. Validation on another scATAC&RNA-seq co-assayed sample ccRCC4 (Fig. 3 f-j, Supplementary Figs. 6b) , and the scATAC-seq-only sample ccRCC1 (Fig. 4 a-d ) demonstrated TeaCNV’s better performance in distinguishing chromosomal segments with distinct CN states. Overall, segments with differing CN states showed non-overlapping distribution of estimated integer CNs (average dispersion score: 0.97) in TeaCNV, compared to ambiguous distributions in alternative methods (average dispersion score: 0.06–0.27, Fig. 4 e). Across all samples, TeaCNV achieved 98.6% overall accuracy (precision: 93.0%, recall: 96.1%, F1 score: 0.92) for CNV event identification compared to the bulk WGS ground truth, outperforming existing methods (Fig. 4 f and g ). Quantitative analysis revealed the smallest error (average RMSE: 0.07) in absolute CN profiles estimated by TeaCNV (Fig. 4 h, Methods ). Identifying Tumor clonal substructures in solid tumors We applied TeaCNV to published scATAC-seq datasets spanning pancreatic (PDAC), breast (BRCA), head and neck (HNSCC), colorectal (CRC), and ovarian (OV) cancers ( Supplementary Fig. 7 ) 29 . In the 1,694 epithelial cells from the PDAC sample, TeaCNV identified two clones, one of which was a diploid subpopulation exhibiting a reduced epithelial score (Fig. 5 a and b ). TeaCNV successfully identified driver events of PDAC, including AKT2 amplifications and SMAD4 deletions, which showed different chromatin accessibility patterns between the diploid and aneuploid subpopulations 30 . In the BRCA, HNSCC, CRC and OV samples, all epithelial cells were identified as aneuploid and exhibited polyclonal architectures (Fig. 5 c-j). Driver genes of BRCA, including MCL1 , MYC , CSMD1 , RB1 and WWOX were identified as truncal CNV events across clones in the BRCA sample (Fig. 5 c and d ) 31 . Although these clones showed similar epithelial scores, they were distinguished by subclonal CNV events, such as chr11q and chr20q13.2. In particular, the oncogene ZNF217 showed different CN states, driving changes in ATAC signals across clones (Fig. 5 d). In HNSCC, driver genes including LRP1B , FGFR1 and MYC were identified as truncal CNVs. However, SOX2, PDE4D and KLF1 were identified as subclonal events, showing distinct chromatin accessibility patterns across clones ( Fig. 5 e and f ) 32 . Similarly, driver genes of CRC ( APC , MYC , CCND2 , CDK8 , KLF5 , HNF4A , SMAD4 and MACROD2 , Fig. 5 g and h ) 33 and OV ( MCL1 , MECOM , KRAS , METTL17 and CCNE1 , Fig. 5 i and j ) 34 were successfully identified as truncal CNVs in the corresponding samples. While clones across most samples exhibited conserved epithelial scores, divergent chromatin accessibility patterns at subclonal loci underscored the genomic heterogeneity. Discussion We present TeaCNV, a computational framework that enables robust estimation of absolute CN profiles and clonal architecture directly from scATAC-seq data, without requiring matched bulk DNA sequencing. By aggregating sparse chromatin accessibility signals across epigenetically homogeneous cells, TeaCNV resolves genomic heterogeneity in single cell epigenomic datasets. Validation against bulk WGS in renal cell carcinoma demonstrated 98.6% accuracy in CNV detection and less than 10% error of integer CN profiling. Applied to six cancer types, TeaCNV revealed polyclonal architectures in 77.78% of cases, demonstrating its broad utility across solid tumors. Compared to scRNA-based CNV inference tools, TeaCNV leverages scATAC-seq’s genome-wide coverage to avoid transcriptional confounding, while overcoming data sparsity through two innovations: iterative clustering of cells to aggregate sparse signal and likelihood-based modeling integrating segment length and variance for ploidy estimation. Unlike existing scATAC-based methods 24 , 25 , TeaCNV quantifies absolute CN states rather than relative gains/losses. This precision enables accurate functional dissection of large-scale CNV regions 10 , 35 , 36 , 37 , such as distinguishing the distinct oncogenic roles of co-altered genes in single cell multiomics sequencing data. Three considerations guide the application of TeaCNV. First, it performs optimally on subclones composed of sufficiently homogeneous cells, as sparse data limits the ability to call single-cell CNVs —a limitation shared by all scATAC-based approaches. Tumor cell lines or samples with high subclonal diversity (e.g., those containing numerous rare subpopulations) may yield unstable estimates. The minimal subclonal size we identified is 30 cells. Second, the current resolution of TeaCNV captures arm-level alterations. Focal events may be missed by sparse chromatin accessibility signals. Third, TeaCNV resolves total CN but not allele-specific CNV estimation. Efforts are underway to extend its functionality to enable allele-specific copy number estimation, which would further enhance its utility for analyzing complex tumor genomes. Tumor plasticity arises from aberrant activation of transcriptional programs caused by genetic and non-genetic mechanisms 38 , 39 , 40 . With advancements in single cell sequencing technologies, simultaneous profiling of the epigenome and transcriptome at single-cell resolution is now possible 41 , 42 . TeaCNV bridges this gap in single cell multiomics by connecting clonal genotypes to epigenomic and transcriptomic phenotypes. By enabling absolute CNV analysis in any scATAC-seq dataset, TeaCNV empowers systematic exploration of how genomic instability shapes epigenetic diversity—a critical factor underlying therapeutic resistance and metastatic progression. Methods Human specimens A total of four ccRCC patients were collected with informed consent, following approval from the Institutional Ethics Committee for Clinical Research approval at the Qilu Hospital of Shandong University, including two patients for scATAC-seq, two for scATAC&RNA-seq co-assays, and three patients for bulk WGS. Experimental methods Nuclei isolation from tissues Tissues were dissected, snap-frozen in liquid nitrogen and stored at -80°C. Nuclei were isolated from the frozen tissues using a protocol designed for scATAC-seq (10X Genomics, CG000212 Rev B). For scATAC&RNA-seq co-assays, nuclei isolation was performed according to the manufacturer's recommended protocol (10X Genomics, CG000375 Rev C), which included RNase inhibitors in the buffer to prevent mRNA degradation during cell lysis. Droplet library preparation and sequencing Nuclei and barcoded beads were pooled together and loaded into the 10X Genomics system for single cell ATAC sequencing following the manufacturer’s instructions (10X Genomics, CG000496 Rev A for Chromium Next GEM Single Cell ATAC Library Kit v2). scATAC&RNA-seq co-assays libraries were prepared using the Chromium Next GEM Single Cell Multiome ATAC + Gene Expression kit, following a separate protocol (10X Genomics, CG000338 Rev F). The barcoded libraries were then pooled and sequenced on the Illumina NovaSeq 6000 system with the associated cells. Library Preparation and Sequencing for bulk WGS About 0.6 µg high-quality genomic DNA was sheared with Covaris LE220 Sonicator (Covaris) to about 350 bp. The library was constructed according to the protocol of KAPA Hyper Prep kit (Roche). First the fragmented DNA was purified using sample purification beads, and the product was repaired by the end and A base was added to the 3' end. Then, the adapters are ligated with the specific barcode sequence. The CleanNGS magnetic beads (CleanNA) were used to screen out incomplete connections and self-connecting products. Sequencing libraries were formed by PCR amplification using universal primers complementary to the adaptor sequences. Paired-end sequencing was performed using the NovaSeq 6000 S4 Reagent Kit v1.5 (300 cycles) on Illumina NovaSeq 6000 platform (Illumina, San Diego, USA) by Sequanta technologies (Shanghai, China). scATAC-seq and scATAC&RNA-seq data processing Reads from scATAC-seq and scATAC&RNA-seq datasets were aligned to the GRCh38 (hg38) reference genome and quantified using the cellranger-atac count (v.1.2.0) and the cellranger-arc count (v.2.0) pipelines (10x Genomics), respectively. Peaks were identified using the MACS3 tool (v.3.0.0) 43 through the CallPeaks function in the Signac package (v.1.14.0, https://github.com/timoast/signac ) 44 . Peaks on chromosomes X and Y, and those within the ENCODE Unified GRCh38 Blacklist regions, were removed using ‘blacklist_hg38_unified’ in the subsetByOverlaps function in Signac. The resulting sample-specific peak set was used to generate the peak-count matrix using FeatureMatrix function in Signac package for downstream analyses. Quality control of single cell sequencing data Quality-control filtering of the scATAC-seq and scATAC&RNA-seq data was performed using functions from the Signac package. Filters applied for the cell inclusion were as follows: number of fragments in peaks > 1,000; number of peaks in cell > 2000; and enrichment-score for Tn5-integration events at transcriptional start sites > 3. Normalization, dimensionality reduction, clustering and cell tying The filtered peak-count matrix from scATAC-seq or scATAC&RNA-seq data was normalized using term frequency-inverse document frequency (TF-IDF) normalization implemented in the Signac package. This normalization accounts for variations in coverage across cells and peaks. The top 95% of peaks were selected as features for dimensionality reduction. We used the RunSVD function to perform singular value decomposition (SVD) on the normalized TF-IDF matrix, a method that is also known as latent semantic indexing (LSI) dimension reduction. The resulting 2:30 LSI components were used for nonlinear dimensionality reduction using the RunUMAP function from the Seurat package. Nuclei were clustered using a graph-based clustering approach implemented in Seurat using the 2:30 LSI components. For scATAC-seq data, we annotated cell types based on the activity of canonical cell type-specific markers, including epithelial cell ( EPCAM and KRT family genes), immune cell ( PTPRC ), endothelial cell ( PECAM1 and CD34 ), and stromal cell ( COL1A2 and COL1A3 ). For scATAC&RNA-seq data, we applied Seurat package (v.4.4.0) 45 to the gene-count matrix for scaling, normalization and identification of highly variable genes for unsupervised cell clustering with default parameters. The elbow plot was generated with the ElbowPlot function of Seurat, and based on this, the number of significant principal components (PCs) was determined. In this study, the top 2,000 highly variable genes and the first 30 PCs identified by Seurat were used for unsupervised clustering analysis. We annotated cell types based on the expression of canonical cell type-specific markers. Bulk WGS data processing and CNV analysis Sequencing reads from bulk tumor tissue and matched normal tissues were aligned to the reference human genome (GRCh38) using Burrows–Wheeler Aligner (BWA v0.7.17) software 46 to obtain the original mapping results stored in BAM format. SAMtools 47 were used to sort and index BAM files. The ‘runVarbin’ module of copykit (v0.1.2) 48 was used to count the number of reads in 220-kb genomic bins defined by the GRCh38 genome assembly, with GC correction applied to the counts. We calculated the log ratio of tumor to normal tissue for each genomic bin using the bin counts, enabling CNV calling, segmentation, and the estimation of absolute CNs for each segment. TeaCNV Algorithm Preprocessing and transformation of scATAC-seq data To balance data sparsity and genomic coverage in the input peak-cell matrix, we exclude peaks detected in fewer than 5% of cells, ensuring that the retained number of peaks remains above 10,000. If this threshold is not met, the detection criterion is relaxed by progressively lowering the minimum proportion of cells required for peak detection. In solid tumors, immune cells or confident non-malignant cells serve as reference cells, while epithelial cells or candidate malignant cells are considered inferred cells. High-confidence reference cells are defined as those whose total peak read counts and the number of detected peaks fall within the 5th to 95th percentile range of all reference cells. The same percentile-based filtering is applied to inferred cells using metrics specific to this group. To mitigate bias in ploidy estimation due to chromosomal dropout, cells are excluded if more than 60% of its measurements for any chromosome have a value of zero. We construct two peak-cell matrices: \(\:{X}_{n\times\:m}\) , where \(\:n\) represents the number of peaks ( \(\:n\:\ge\:\:\text{10,000}\) ) and \(\:m\) represents the number of reference cells, and \(\:{Y}_{n\times\:t}\) where \(\:t\) is the number of inferred cells. To reduce the impact of extreme peak values, we cap the values in matrices X and Y to the range [0,4] by replacing all values exceeding 4 with 4. To correct for variations in sequencing depth across cells, we respectively normalize matrices \(\:X\) and \(\:Y\) as follow: $$\:{{X}^{{\prime\:}}}_{i,j}=\frac{{X}_{i,j}}{\sum\:_{k}{X}_{k,j}}\bullet\:\frac{1}{m}\bullet\:\sum\:_{j=1}^{m}\sum\:_{k=1}^{n}{X}_{k,j}$$ $$\:{{Y}^{{\prime\:}}}_{i,j}=\frac{{Y}_{i,j}}{\sum\:_{k}{Y}_{k,j}}\bullet\:\frac{1}{t}\bullet\:\sum\:_{j=1}^{t}\sum\:_{k=1}^{n}{Y}_{k,j}$$ For the i -th peak, we calculate the average value across reference cells using the normalized matrix \(\:{X}^{{\prime\:}}\) : $$\:{\stackrel{-}{X{\prime\:}}}_{i·}=average\left[{X{\prime\:}}_{i,1\dots\:m}\right]\:$$ For each peak in each inferred cell, we calculate the ratio relative to reference average signal to generate the ratio matrix \(\:{R}_{n\times\:t}\) , where: $$\:{R}_{i,j}=\frac{{Y{\prime\:}}_{i,j}}{{\stackrel{-}{X{\prime\:}}}_{i·}}\:(1\le\:i\le\:n,\:1\le\:j\le\:t)$$ Aggregating cell sub-populations The ratio matrix \(\:{R}_{n\times\:t}\) is sorted by the chromosomal location of the analyzed peaks. To capture variations derived from chromosomal segments rather than particular peaks, we calculate the average ratio values within a genomic window of 5 peaks for each chromosome. Subsequently, we identify highly variable features and perform PCA based on the merged matrix using the Seurat package. The top 2,000 highly variable features and the first 50 PCs identified by Seurat are used for unsupervised clustering analysis, with the resolution parameter set to 1. This approach enables the preliminary classification of observed cells into distinct subgroups. For cells within the same subgroup, we aggregate ratios by averaging across cells to obtain subgroup-specific ratios as: $$\:{R}^{s}=[{\stackrel{-}{R}}_{1·}^{s},\dots\:,{\stackrel{-}{R}}_{n·}^{s}]\:$$ $$\:{\stackrel{-}{R}}_{i·}^{s}=average\left[{R}_{i,kϵS}\right]\:$$ Here, \(\:{\stackrel{-}{R}}_{i·}^{s}\) represents the average ratio for the i -th peak corresponding to subgroup S. Genome segmentation of cell subpopulation For the s -th subgroup, genome segmentation is performed on \(\:{R}^{s}\) organized according to chromosomal locations. Breakpoints are identified using the pruned exact linear time (PELT) algorithm from changepoint R package 49 , along with an alternative method employing FPOP algorithm from robseg R package 50 . To ensure the reliability of the segments, the chromosomal segments shorter than 2MB are excluded. The ratio for each chromosomal segment is estimated by calculating the median value of the peaks within the same segment. This results in the segmental ratio, which is considered as relative copy ratio for the s -th subgroup: $$\:{R}^{s}=\left[{R}_{1}^{s},\dots\:,{R}_{l}^{s}\right]$$ $$\:\:{R}_{j}^{s}=median\left({\stackrel{-}{R}}_{{j}_{k\in\:{segment}_{j}},\bullet\:}^{s}\right)$$ Where l is the number of segments. Inferring absolute copy number The algorithm for absolute CN estimation is based on the ABSOLUTE algorithm 8 , modified to accommodate single-cell sequencing data. For the estimation of each subgroup, the input segmental ratio R with standard error \(\:\sigma\:\) consists of \(\:{R}_{j},\:j\in\:\{1,\dots\:,l\}\) , corresponding to a genomic fraction \(\:{w}_{j}\) . Each \(\:{R}_{j}\) is assumed to arise from one of the integer CN states in the set \(\:Q=\{\text{1,2},\dots\:,Q\}\) with probabilities \(\:p\left({q}_{j}\right),\) where \(\:q\in\:Q\) . The observed \(\:{R}_{j}\) is modeled as a mixture of \(\:Q\) Gaussian components located at \(\:\mu\:=\left\{{\mu\:}_{q\in\:Q}\right\}\) representing expected ratio of integer CN state. $$\:P\left({R}_{j}|\mu\:,{\sigma\:}^{2},\theta\:,{w}_{j}\right)={\sum\:}_{q\in\:Q}{w}_{j}\bullet\:P\left({q}_{j}|{\theta\:}_{q}\right)\bullet\:N\left({\mu\:}_{q},{\sigma\:}^{2}\right)$$ Where, \(\:\theta\:={\{\theta\:}_{q\in\:Q}\}\) reprents the mixture weights, indicating the expected genomic fraction allocated to each CN state. Given limited knowledge about copy-state \(\:q\) , the distribution is chosen to have maximum entropy: $$\:P\left({q}_{j}|{\theta\:}_{q}\right)=\frac{{e}^{-{\theta\:}_{q}\bullet\:{q}^{\#}\bullet\:{w}_{j}}}{\sum\:_{k\in\:Q}{e}^{-{\theta\:}_{k}\bullet\:{k}^{\#}\bullet\:{w}_{j}}}$$ Where \(\:{q}^{\#}\) indicates the order of \(\:q\) in the CN state set \(\:Q,\) beginning with 1. The unknown parameters \(\:\theta\:={\{\theta\:}_{q\in\:Q}\}\) are estimated by minimizing the loss function: $$\:\underset{\theta\:}{\widehat{\theta\:}={arg}\text{min}}{\left\{\sum\:_{q\in\:Q}{\left(\sum\:_{j=1}^{l}{w}_{j}\bullet\:P\left({q}_{j}|{\theta\:}_{q}\right)-{\theta\:}_{q}\right)}^{2}\right\}}^{\frac{1}{2}}$$ The full log-likelihood of the input data \(\:R\) is then computed as: $$\:logL\left(R|\mu\:,\sigma\:,\widehat{\theta\:},w\right)=\sum\:_{j=1}^{l}logL\left({R}_{j}|\mu\:,{\sigma\:}^{2},\widehat{\theta\:},{w}_{j}\right)$$ where $$\:L\left({R}_{j}|\mu\:,{\sigma\:}^{2},\widehat{\theta\:},{w}_{j}\right)={\sum\:}_{q\in\:Q}{w}_{j}\bullet\:P\left({q}_{j}|{\widehat{\theta\:}}_{q}\right)\bullet\:N\left({\mu\:}_{q},{\sigma\:}^{2}\right)$$ The unknown parameters \(\:\mu\:={\{\mu\:}_{q\in\:Q}\}\) and \(\:{\sigma\:}^{2}\) are estimated using a combination of the Nelder-Mead optimization algorithm and maximum likelihood estimation. Posterior probabilities are used to infer the copy-state indicators for each segment, with the absolute CN state of genomic segment j is defined as the state corresponding to the maximum posterior probability: $$\:P\left({\widehat{q}}_{j}\right)=p\left({q}_{j}|{\widehat{\theta\:}}_{q}\right)\bullet\:\frac{N\left({R}_{j}|{\widehat{\mu\:}}_{q},{\widehat{\sigma\:}}^{2}\right)}{L\left({R}_{j}|\widehat{\mu\:},{\widehat{\sigma\:}}^{2},\widehat{\theta\:},{w}_{j}\right)}$$ Optimizing copy number estimation Each subgroup corresponds to an estimated expected copy-ratio \(\:\mu\:={\{\mu\:}_{q\in\:Q}\}\) associated with the absolute CN state set \(\:Q\) . Ideally, the interval between two consecutive CN states, defined as \(\:{\Delta\:}={\mu\:}_{q}-{\mu\:}_{q-1}\:(\forall\:q\in\:Q)\) , should remain consistent across the values derived from the same subgroup. However, bias in the estimation of \(\:\mu\:\) may arise due to inhomogeneous observed values. To address this potential bias, we iterate through all possible values of \(\:{\Delta\:}\) for each subgroup, and include \(\:{k}_{1}\bullet\:{\Delta\:}\:({where\:k}_{1}=2)\) for the \(\:{\Delta\:}0.6\) to update expected copy-ratio \(\:\mu\:={\{\mu\:}_{q\in\:Q}\}\) . We calculate the Akaike Information Criterion (AIC) value for all candidate \(\:\mu\:\) and determine the optimal CN estimation for each corresponding subgroup based on the minimum AIC value. Subsequently, the overall ploidy of the subgroup is then calculated based on the estimated integer CN states of segments weighted by their genomic fraction: $$\:ploidy=\sum\:_{j=1}^{l}{w}_{j}\bullet\:{\widehat{q}}_{j}$$ Scoring the confidence of cell subgroup with estimated copy number profile To assess the confidence of identified homogeneous subgroups, we score each subgroup using a combination of mean squared error (MSE) between observed and expected copy ratio, and the proportion of genomic segments explained by the expected copy ratios. First, the MSE between the observed and expected copy-ratios is calculated as: $$\:MSE=\frac{1}{l}\sum\:_{i=1}^{l}\left({R}_{i}-{\mu\:}_{i}\right)$$ where l is the total number of segments and \(\:{\mu\:}_{i}\) corresponds to the expected ratio of CN state for segment \(\:i\) . Next, we calculate the proportion of genomic segments explained by the expected copy ratios: $$\:F=\sum\:_{q\in\:Q}\left(\sum\:_{i\in\:\left\{l\right|\left|{R}_{l}-{\mu\:}_{q}\right|<d\}}\frac{{length}_{i}}{L}\right)\:$$ where \(\:L\) is the total length of the genome. \(\:d\) is the allowed maximal difference between observed and expected ratio (default is 0.1), and \(\:{length}_{i}\) is the length of segment i . The overall score, weighted by a parameter \(\:\beta\:\) , is defined as: $$\:score=F\bullet\:\left(-logMSE\right)\bullet\:\beta\:$$ where: $$\:\beta\:=\frac{{\delta\:}_{{\Delta\:}}}{log\left(1+l\right)\bullet\:{e}^{ploidy-2}}$$ and: $$\:{\delta\:}_{{\Delta\:}}=\left\{\begin{array}{c}\varDelta\:,\:if\:\varDelta\:<0.3\\\:1,\:otherwise\end{array}\right.$$ The weighted parameter \(\:\beta\:\) decreases with smaller interval between consecutive CN states, a larger number of segments, and greater difference in ploidy between estimated and diploidy. Therefore, a higher score indicates a more reliable and highly homogeneous subgroup, reflecting that more genomic regions can be explained by the integer CN states with smaller MSE, appropriate interval between consecutive CN states and consistent genomic segmentation. Subclone partitioning We define an initial adjacency matrix \(\:A\) to represent whether two subgroups can be merged into one clone. $$\:{A}_{i,j}=\left\{\begin{array}{c}1,\:i=j\\\:0,i\ne\:j\end{array}\right.$$ here, a value of 1 indicates that subgroup \(\:i\) and \(\:j\) can be merged, while a value of 0 indicates they cannot. To update the adjacency matrix \(\:A\) , we compare the peak signals of subgroup in a pair-wise manner using the peak value at the group level, defined as follow: $$\:{{Y}^{{\prime\:}}}_{{s}_{i}}={\left[\stackrel{-}{Y}{{\prime\:}}_{1},\dots\:,\stackrel{-}{Y}{{\prime\:}}_{n}\right]}^{{s}_{i}}$$ $$\:{{Y}^{{\prime\:}}}_{{s}_{j}}={\left[\stackrel{-}{Y}{{\prime\:}}_{1},\dots\:,\stackrel{-}{Y}{{\prime\:}}_{n}\right]}^{{s}_{j}}$$ Here, \(\:\:{{Y}^{{\prime\:}}}_{{s}_{i}}\) and \(\:{{Y}^{{\prime\:}}}_{{s}_{j}}\) represents the average values of peaks across cells from subgroup \(\:{s}_{i}\) and \(\:{s}_{j}\) , respectively. We then calculate the odds ratios between subgroup \(\:{s}_{i}\) and \(\:{s}_{j}\) : $$\:{R}_{{s}_{i},{s}_{j}}=\frac{{{Y}^{{\prime\:}}}_{{s}_{i}}}{{{Y}^{{\prime\:}}}_{{s}_{j}}}\:$$ Next, genome segmentation is performed based on \(\:{R}_{{s}_{i},{s}_{j}}\) using the same approach described previously. For each segment, we compare the distribution of peak values (equivalent to copy-ratio) between subgroup \(\:{s}_{i}\) and \(\:{s}_{j}\) using a two-sided t-test, and adjust the resulting P values for multiple hypotheses through the Benjamini-Hochberg method. If the difference in any one segment is significant (adjusted p-value < 0.05) and the estimated integer CN is different between \(\:{s}_{i}\) and \(\:{s}_{j}\) , we set: $$\:{A}_{i,j}=\:{A}_{j,i}=0$$ Otherwise, we update the matrix as: $$\:{A}_{i,j}=\:{A}_{j,i}=1$$ Based on the updated adjacent matrix, we obtain the final subclonal partitioning. For the final subclone substructure, we update the corresponding integer CN profile using the same method as previously described. Application of other CNV inference approaches We detected CNVs using epiAneufinder 24 and Copy-scAT from scATAC-seq of three ccRCC samples which have sample-matched bulk WGS data. The analysis was performed with each method’s default parameters to identify CNVs at single cell resolution. In addition, for two ccRCC samples with scATAC&RNA-seq, we employed inferCNV 13 on the gene expression data, using the same reference as TeaCNV. We followed the recommended parameters for 10X (denoise = TRUE, cutoff = 0.1) and performed CNV calling using the ‘consensus’ i6 HMM mode. Evaluating performance of copy number estimation To evaluate the performance of CNV detection, we used CNVs identified by bulk WGS as the ground truth and calculated precision, recall, ACC and F1 score based on the overlap between the predicted and true CNVs. The genome was divided into 100-kb bins, excluding those overlapping breakpoints in either ground truth or inferred results. We defined true positive (TP) as bins with CNVs present in both the group truth and the inferred results, false positive (FP) as bins with CNVs present in the inferred results but absent in the ground truth, and false negative (FN) as bins with CNVs in the ground truth but not detected in the inferred results. Precision, recall, ACC and F1 score were then calculated using the following formulas: $$\:precision=\:\frac{TP}{TP+FP}$$ $$\:recall=\:\frac{TP}{TP+FN}$$ $$\:ACC=\:\frac{TP+TN}{TP+TN+FP+FN}$$ $$\:F1=\:\frac{2\bullet\:precision\bullet\:recall}{precision+recall}$$ TeaCNV reports CNVs at clonal level, we averaged the precision, recall and F1 score across the subclones derived from each sample. For epiAneufinder and Copy-scAT, CNV events were considered detected at the bulk level if they were identified in a specific proportion of cells, ranging from 10–100% in 10% increments. The average precision, recall, and F1 score were then calculated across these varying cell proportion thresholds. To evaluate the accuracy of inferred CN profiles, we calculated the deviation of estimations from WGS using the root mean square error (RMSE) metric: $$\:RMSD=\:\sqrt{\frac{1}{n}\bullet\:\sum\:_{i=1}^{n}{\left({CN}_{i}-{CN}_{i}^{true}\right)}^{2}}$$ Here, \(\:{CN}_{i}\) represents the inferred integer CN for bin i and \(\:{CN}_{i}^{true}\) is the true integer CN derived from WGS. Because epiAneuFinder and Copy-scAT do not report integer CNs directly, we centered the inferred CNV scores to align with integer CNs for bins with distinct absolute CN states in WGS data. The definition of \(\:{CN}_{i}\) was as follow: $$\:{CN}_{i}=\left\{\begin{array}{c}{CN}_{i}\:\:\:for\:TeaCNV\\\:{CNV\:score}_{i}-\underset{j\in\:\left\{j|{{CN}}_{j}^{{true}}={{CN}}_{i}^{{true}}\right\}}{\text{average}}\left({CNV\:score}_{j}\right)+{CN}_{i}^{true}\:\:\:for\:other\:methods\end{array}\right.$$ For TeaCNV, \(\:{CN}_{i}\) is taken directly from the inferred integer CNs. For other methods, \(\:{CN}_{i}\) is adjusted by subtracting the average of the CNV scores from bins with the same true CN state, followed by adding the true absolute CN. This adjustment allows for a comparison of inferred profiles against the ground truth. In theory, the inferred values should be confined to discrete ranges defined by distinct absolute CN states. Ideally, the density peaks in the distribution of inferred CNs for genomic regions with the same CN state should not overlap with those for regions with different CN states. The dispersion in the distribution of inferred value from the genomic regions with distinct CN states reflects the confidence of estimation. To quantify the dispersion of CN estimations, we use the following formula: $$\:dispersion=\:\sum\:_{s}\left[{\int\:}_{A}^{}max\left({p}_{s}\left(x\right),{p}_{s+1}\left(x\right)\right)dx\right]$$ $$\:A=\left\{x\right|{p}_{s}\left(x\right)\ne\:0\:\&\:{p}_{s+1}\left(x\right)\ne\:0\}$$ Where \(\:s\) represents unique absolute CN states from WGS data. \(\:{p}_{s}\left(x\right)\) is the probability density of the inferred values for the genomic regions with integer CN s in WGS data. For TeaCNV, the inferred value is output integer CN profiles. For other methods, the inferred value is the estimated CNV score. The interval A denotes the range of integration that includes all x values where both \(\:{p}_{s}\left(x\right)\) and \(\:{p}_{s+1}\left(x\right)\) are non-zero. A larger dispersion indicates a greater difference in inferred values between distinct absolute CN states, suggesting stronger confidence in the separation of these states. Declarations Data availability The scATAC-seq, scATAC&RNA-seq aco-assays and bulk WGS validation data of ccRCC can be accessed through the link with token (https://zenodo.org/records/14190637?preview=1&token=eyJhbGciOiJIUzUxMiJ9.eyJpZCI6Ijg1YTczZDU4LTlkNjQtNGZkYS05ZjY3LTA3M2JhMjcwOTY2MiIsImRhdGEiOnt9LCJyYW5kb20iOiJlYWRjZmU3ZDU1M2I1ZjMxYTNlNGQ2MThiYzFiNDhhNCJ9.i4vMudpQCvXLnj3kkcBhJpLoS9GD_6rb-tQW0MjAq7AG1H2RHG84ZWWm5jhc2cmd5fMl-QOcmCXUk_7MVjib_A). scATAC-seq data of BRCA, PDAC, HNSCC, CRC and OV were downloaded from Gene Expression Omnibus (GEO) with accession number GSE240822. Code availability The TeaCNV algorithm is available at https://github.com/ShaojunLab/TeaCNV. The analysis scripts used to reproduce results included in the paper are available at https://github.com/ShaojunLab/TeaCNVanalysis. Competing interests The authors declare no competing interests. Author contributions Y. W., S. Z. and F. W. formulated the study and the overall approach. Y. W. developed and implemented the computational algorithm with contribution form Y. Z.. X. Y. provide clinical samples. Y. D., Y. C. and Z. C. performed single-cell sequencing experiments and bulk WGS. Y.W. performed the analysis with the help from M. Z., X. W. and H. L.. Y. W., F. W. and S. Z. drafted the manuscript. All authors provided suggestions and corrections on the manuscript text. Acknowledgements This work was supported by the National Natural Science Foundation of China (grants 32170666 to S. Zhang, grants 32400528 to Y. Wang; grant 32270699 to F. Wang), Guangdong Pearl River Program (grant 2021QN02Y180 to F. Wang). References Esteller M, Dawson MA, Kadoch C, Rassool FV, Jones PA, Baylin SB (2024) The Epigenetic Hallmarks of Cancer. Cancer Discov 14:1783–1809 Chakravarthi BV, Nepal S, Varambally S (2016) Genomic and Epigenomic Alterations in Cancer. Am J Pathol 186:1724–1735 Friedlander ML, Hedley DW, Taylor IW (1984) Clinical and biological significance of aneuploidy in human tumours. J Clin Pathol 37:961–974 Clark GM, Dressler LG, Owens MA, Pounds G, Oldaker T, McGuire WL (1989) Prediction of relapse or survival in patients with node-negative breast cancer by DNA flow cytometry. N Engl J Med 320:627–633 Kallioniemi OP, Punnonen R, Mattila J, Lehtinen M, Koivula T (1988) Prognostic significance of DNA index, multiploidy, and S-phase fraction in ovarian cancer. Cancer 61:334–339 Choma D, Daures JP, Quantin X, Pujol JL (2001) Aneuploidy and prognosis of non-small-cell lung cancer: a meta-analysis of published data. Br J Cancer 85:14–22 Bielski CM et al (2018) Genome doubling shapes the evolution and prognosis of advanced cancers. Nat Genet 50:1189–1195 Carter SL et al (2012) Absolute quantification of somatic DNA alterations in human cancer. Nat Biotechnol 30:413–421 Salehi S et al (2021) Clonal fitness inferred from time-series modelling of single-cell cancer genomes. Nature 595:585–590 Minussi DC et al (2021) Breast tumours maintain a reservoir of subclonal diversity during expansion. Nature 592:302–308 Casasent AK et al (2018) Multiclonal Invasion in Breast Tumors Identified by Topographic Single Cell Sequencing. Cell 172:205–217e212 Lei Y et al (2021) Applications of single-cell sequencing in cancer research: progress and perspectives. J Hematol Oncol 14:91 Patel AP et al (2014) Single-cell RNA-seq highlights intratumoral heterogeneity in primary glioblastoma. Science 344:1396–1401 Tirosh I et al (2016) Dissecting the multicellular ecosystem of metastatic melanoma by single-cell RNA-seq. Science 352:189–196 Tirosh I et al (2016) Single-cell RNA-seq supports a developmental hierarchy in human oligodendroglioma. Nature 539:309–313 Venteicher AS et al (2017) Decoupling genetics, lineages, and microenvironment in IDH-mutant gliomas by single-cell RNA-seq. Science 355 Puram SV et al (2017) Single-Cell Transcriptomic Analysis of Primary and Metastatic Tumor Ecosystems in Head and Neck Cancer. Cell 171:1611–1624e1624 Gao R et al (2021) Delineating copy number and clonal substructure in human tumors from single-cell transcriptomes. Nat Biotechnol 39:599–608 Gao T et al (2023) Haplotype-aware analysis of somatic copy number variations from single-cell transcriptomes. Nat Biotechnol 41:417–426 Fan J et al (2018) Linking transcriptional and genetic tumor heterogeneity through allele analysis of single-cell RNA-seq data. Genome Res 28:1217–1227 Hou Y et al (2016) Single-cell triple omics sequencing reveals genetic, epigenetic, and transcriptomic heterogeneity in hepatocellular carcinomas. Cell Res 26:304–319 Buenrostro JD et al (2015) Single-cell chromatin accessibility reveals principles of regulatory variation. Nature 523:486–490 Grandi FC, Modi H, Kampman L, Corces MR (2022) Chromatin accessibility profiling by ATAC-seq. Nat Protoc 17:1518–1552 Ramakrishnan A et al (2023) epiAneufinder identifies copy number alterations from single-cell ATAC-seq data. Nat Commun 14:5846 Nikolic A et al (2021) Copy-scAT: Deconvoluting single-cell chromatin accessibility of genetic subclones in cancer. Sci Adv 7:eabg6045 Wu CY et al (2021) Integrative single-cell analysis of allele-specific copy number alterations and chromatin accessibility in cancer. Nat Biotechnol 39:1259–1269 Killick R, Fearnhead P, Eckley IA (2012) Optimal Detection of Changepoints With a Linear Computational Cost. J Am Stat Assoc 107:1590–1598 Ricketts CJ et al (2018) The Cancer Genome Atlas Comprehensive Molecular Characterization of Renal Cell Carcinoma. Cell Rep 23:313–326e315 Terekhanova NV et al (2023) Epigenetic regulation during cancer transitions across 11 tumour types. Nature 623:432–441 Cancer Genome Atlas Research Network. Electronic address aadhe, Cancer Genome Atlas Research N. Integrated Genomic Characterization of Pancreatic Ductal Adenocarcinoma. Cancer Cell 32, 185–203 e113 (2017) Cancer Genome Atlas N (2012) Comprehensive molecular portraits of human breast tumours. Nature 490:61–70 Cancer Genome Atlas N (2015) Comprehensive genomic characterization of head and neck squamous cell carcinomas. Nature 517:576–582 Cancer Genome Atlas N (2012) Comprehensive molecular characterization of human colon and rectal cancer. Nature 487:330–337 Cancer Genome Atlas Research N (2011) Integrated genomic analyses of ovarian carcinoma. Nature 474:609–615 Rice AM, McLysaght A (2017) Dosage sensitivity is a major determinant of human copy number variant pathogenicity. Nat Commun 8 Wang F et al (2019) Integrated transcriptomic-genomic tool Texomer profiles cancer tissues. Nat Methods 16:401–404 Shi H, Williams MJ, Satas G, Weiner AC, McPherson A, Shah SP (2024) Allele-specific transcriptional effects of subclonal copy number alterations enable genotype-phenotype mapping in cancer cells. Nat Commun 15:2482 Shlyakhtina Y, Moran KL, Portal MM (2021) Genetic and Non-Genetic Mechanisms Underlying Cancer Evolution. Cancers (Basel) 13 Marine JC, Dawson SJ, Dawson MA (2020) Non-genetic mechanisms of therapeutic resistance in cancer. Nat Rev Cancer 20:743–756 Tellez-Gabriel M, Ory B, Lamoureux F, Heymann MF, Heymann D (2016) Tumour Heterogeneity: The Key Advantages of Single-Cell Analysis. Int J Mol Sci 17 Lee J, Hyeon DY, Hwang D (2020) Single-cell multiomics: technologies and data analysis methods. Exp Mol Med 52:1428–1442 Baysoy A, Bai Z, Satija R, Fan R (2023) The technological landscape and applications of single-cell multi-omics. Nat Rev Mol Cell Biol 24:695–713 Zhang Y et al (2008) Model-based analysis of ChIP-Seq (MACS). Genome Biol 9:R137 Stuart T, Srivastava A, Madad S, Lareau CA, Satija R (2021) Single-cell chromatin state analysis with Signac. Nat Methods 18:1333–1341 Hao Y et al (2021) Integrated analysis of multimodal single-cell data. Cell 184:3573–3587e3529 Li H, Durbin R (2009) Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics 25:1754–1760 Danecek P et al (2021) Twelve years of SAMtools and BCFtools. GigaScience 10 Minussi DC et al (2022) Resolving clonal substructure from single cell genomic data using CopyKit. bioRxiv , 2022.2003.2009.483497 Killick R, Eckley IA (2014) changepoint: AnRPackage for Changepoint Analysis. J Stat Softw 58 Fearnhead P, Rigaill G (2018) Changepoint Detection in the Presence of Outliers. J Am Stat Assoc 114:169–183 Additional Declarations There is NO Competing Interest. Supplementary Files SupplementaryFigures.pdf Supplementary Figures Cite Share Download PDF Status: Posted Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-6609843","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":458967751,"identity":"f6108751-cc3e-4d25-999a-67a6fa645696","order_by":0,"name":"Shaojun Zhang","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAAz0lEQVRIiWNgGAWjYFACHhDxTw5IMB5gYEggWssBYwYGZgbStCQ2EK1Fvv/swc8Fv+6k90vkHzjwoSKNgb+9G78+gwPnkqVn9j3LnTkjmeHgjDM5DBJnzm7Ar4Wxx0Cat4c5d8ONZIbDvG0VDAYSufi1yDfzGP8Gakk3IFoLwzEeM2meH4cToFpyCGsxOMNjZs3bkGY4s+exAdAvaTwE/SLff8b4Ns8fG3l+9sSHDz5UJMvxt/cScBgIMLYh2DyElYPBHyLVjYJRMApGwcgEAJnuSfyNBTziAAAAAElFTkSuQmCC","orcid":"https://orcid.org/0000-0001-5760-5797","institution":"Guangdong Academy of Medical Sciences and Medical Research Institute, Guangdong Provincial People's Hospital (Guangdong Academy of Medical Sciences), Southern Medical University","correspondingAuthor":true,"prefix":"","firstName":"Shaojun","middleName":"","lastName":"Zhang","suffix":""},{"id":458967752,"identity":"a7fe9fba-659e-417d-9756-06914908b783","order_by":1,"name":"Ying Wang","email":"","orcid":"https://orcid.org/0000-0003-4868-6147","institution":"Guangdong Academy of Medical Sciences and Medical Research Institute, Guangdong Provincial People's Hospital","correspondingAuthor":false,"prefix":"","firstName":"Ying","middleName":"","lastName":"Wang","suffix":""},{"id":458967753,"identity":"65ef18d6-1278-4fc6-a812-809fa02f506e","order_by":2,"name":"Yuhao Deng","email":"","orcid":"","institution":"Institute of Precision Medicine, The First Affiliated Hospital, Sun Yat-Sen University","correspondingAuthor":false,"prefix":"","firstName":"Yuhao","middleName":"","lastName":"Deng","suffix":""},{"id":458967754,"identity":"4cdca4a1-92b0-4c27-b0cd-fb7ff22d8aca","order_by":3,"name":"Xinbao Yin","email":"","orcid":"","institution":"Department of Urology, The Affiliated Hospital of Qingdao University","correspondingAuthor":false,"prefix":"","firstName":"Xinbao","middleName":"","lastName":"Yin","suffix":""},{"id":458967755,"identity":"1db9a9c7-28f4-433c-9415-a39c1d1311ee","order_by":4,"name":"Yanru Zhang","email":"","orcid":"","institution":"Institute of Precision Medicine, The First Affiliated Hospital, Sun Yat-Sen University","correspondingAuthor":false,"prefix":"","firstName":"Yanru","middleName":"","lastName":"Zhang","suffix":""},{"id":458967756,"identity":"688f8e00-2097-436d-befa-5b3bd1fd03a3","order_by":5,"name":"Yurong Chen","email":"","orcid":"","institution":"Department of Gynecological Radiotherapy, Harbin Medical University Cancer Hospital","correspondingAuthor":false,"prefix":"","firstName":"Yurong","middleName":"","lastName":"Chen","suffix":""},{"id":458967757,"identity":"09a7779c-2f72-4f9a-91f7-42bc24919cff","order_by":6,"name":"Min Zhang","email":"","orcid":"","institution":"Institute of Precision Medicine, The First Affiliated Hospital, Sun Yat-Sen University","correspondingAuthor":false,"prefix":"","firstName":"Min","middleName":"","lastName":"Zhang","suffix":""},{"id":458967758,"identity":"26904107-9057-4f3b-86ae-8dc706af4111","order_by":7,"name":"Xin Wang","email":"","orcid":"https://orcid.org/0000-0002-6692-8165","institution":"The First Affiliated Hospital of Sun Yat-Sen University","correspondingAuthor":false,"prefix":"","firstName":"Xin","middleName":"","lastName":"Wang","suffix":""},{"id":458967759,"identity":"f4390a93-0d57-4865-accc-a65af686f852","order_by":8,"name":"Hang Li","email":"","orcid":"","institution":"Guangdong Academy of Medical Sciences and Medical Research Institute, Guangdong Provincial People's Hospital","correspondingAuthor":false,"prefix":"","firstName":"Hang","middleName":"","lastName":"Li","suffix":""},{"id":458967760,"identity":"75939b6c-6111-499d-ac82-cc65866b8447","order_by":9,"name":"Zhizhuo Cao","email":"","orcid":"","institution":"Institute of Precision Medicine, The First Affiliated Hospital, Sun Yat-Sen University","correspondingAuthor":false,"prefix":"","firstName":"Zhizhuo","middleName":"","lastName":"Cao","suffix":""},{"id":458967761,"identity":"fff91269-1344-4c6b-9b9e-15886f10939b","order_by":10,"name":"Fang Wang","email":"","orcid":"https://orcid.org/0000-0002-3510-4550","institution":"The First Affiliated Hospital, Sun Yat-sen University","correspondingAuthor":false,"prefix":"","firstName":"Fang","middleName":"","lastName":"Wang","suffix":""}],"badges":[],"createdAt":"2025-05-07 08:36:48","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-6609843/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-6609843/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":83220586,"identity":"a3db33ad-f030-45c0-b835-233de7a31634","added_by":"auto","created_at":"2025-05-21 10:18:26","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":61916,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eOverview of absolute copy number estimation from scATAC-seq data with TeaCNV\u003c/strong\u003e. \u003cstrong\u003e(a)\u003c/strong\u003e Workflow for preprocessing scATAC-seq data. TeaCNV processes the cell-peak matrix to calculate copy ratios relative to reference cells (Ref) and merges the ratios of adjacent peaks. \u003cstrong\u003e(b)\u003c/strong\u003e Cell clustering and genome-wide segmentation. For each cluster, the average copy ratio of cells is computed and used for PELT segmentation. The average copy ratio of bins within each segment is calculated to generate the copy ratio histogram of genomic segments, with bar lengths representing the genomic fraction (x-axis) for each segment. \u003cstrong\u003e(c)\u003c/strong\u003eAbsolute CN estimation. Three potential interpretations of the copy ratio histogram (corresponding to cluster C1 in (b) are presented in terms of absolute CNs. Dotted lines indicate copy ratios corresponding to specific absolute CNs. The optimal absolute CN estimation based on the minimum AIC value. \u003cstrong\u003e(d)\u003c/strong\u003e Optimizing subclone partitioning. Clusters are merged if no significant differences in absolute and relative CNs are observed across all genomic segments. (\u003cstrong\u003ee)\u003c/strong\u003e Output from TeaCNV showing the subclonal absolute CN profiles.\u003c/p\u003e","description":"","filename":"Figure1.png","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/ffa8c07ef86d270056aaa77e.png"},{"id":83220589,"identity":"e0c93432-3639-4e66-8686-25ef1d03fe88","added_by":"auto","created_at":"2025-05-21 10:18:26","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":897608,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eResolving clonal substructures in ccRCC patients using scATAC-seq and scATAC\u0026amp;RNA-seq co-assayed\u003c/strong\u003e. \u003cstrong\u003e(a, b)\u003c/strong\u003e UMAP plots of cells from ccRCC patients, colored by sample and cell type. (a) scATAC\u0026amp;RNA-seq co-assayed data and (\u003cstrong\u003eb\u003c/strong\u003e) scATAC-seq-only data. \u003cstrong\u003e(c)\u003c/strong\u003e Heatmap of copy number ratios in epithelial cells relative to non-epithelial cells (top) and inferred clonal absolute CN profiles (bottom) for ccRCC1. Rows are annotated by clone. \u003cstrong\u003e(d)\u003c/strong\u003e Heatmap of copy number ratios in epithelial cells from ccRCC2. \u003cstrong\u003e(e) \u003c/strong\u003eInferred clonal absolute CN profiles (top) and UMAP plots of epithelial cells colored by activity scores of tumor epithelial markers (\u003cem\u003eEPCAM\u003c/em\u003e, \u003cem\u003eKRT18\u003c/em\u003e, \u003cem\u003eKRT8\u003c/em\u003e, \u003cem\u003eKRT19)\u003c/em\u003e and ATAC signals from representative subclonal regions (bottom) for ccRCC2. \u003cstrong\u003e(f)\u003c/strong\u003e Heatmap of copy number ratios in epithelial cells from ccRCC3.\u003cstrong\u003e (g) \u003c/strong\u003eInferred clonal absolute CN profiles (top) and UMAP plots of epithelial cells colored by activity scores of tumor epithelial markers and ATAC signals from representative subclonal regions (bottom) for ccRCC3. \u003cstrong\u003e(h)\u003c/strong\u003eHeatmap of copy number ratios in epithelial cells from ccRCC4. \u003cstrong\u003e(i) \u003c/strong\u003eInferred clonal absolute CN profiles (top) and UMAP plots of epithelial cells colored by activity scores of tumor epithelial markers and ATAC signals from representative subclonal regions (bottom) for ccRCC4.\u003c/p\u003e","description":"","filename":"Figure2.png","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/71c6482760762df3c6432b27.png"},{"id":83221603,"identity":"7b8167be-41c0-4b0f-8d2c-33fe5c9f9508","added_by":"auto","created_at":"2025-05-21 10:26:26","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":645305,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eComparison of copy number profiles from bulk WGS and TeaCNV inference\u003c/strong\u003e. \u003cstrong\u003e(a) \u003c/strong\u003eInteger\u003cstrong\u003e \u003c/strong\u003eCN profiles from bulk WGS (left) and genome density distribution of absolute CN states (right) for ccRCC3. Significantly altered regions are highlighted in blue (loss) and orange (gain).\u003cstrong\u003e (b) \u003c/strong\u003eClonal integer CN profiles estimated by TeaCNV based on chromatin accessibility (left) and genome density distribution of estimated integer CN states (right) for ccRCC3. In the left panel, bars of each genomic segment correspond to different clones. \u003cstrong\u003e(c - e)\u003c/strong\u003e Heatmap of relative CN profiles (top) and genome density distributions of relative CNV scores (bottom) estimated by (c) epiAneuFinder, (d) Copy-scAT and (e) inferCNV for ccRCC3. Blue and orange highlight regions identified by bulk WGS. \u003cstrong\u003e(f) \u003c/strong\u003eInteger CN profiles from bulk WGS (left) and genome density distribution of absolute CN states (right) for ccRCC4. \u003cstrong\u003e(g) \u003c/strong\u003eClonal integer CN profiles estimated by TeaCNV based on chromatin accessibility (left) and genome density distribution of estimated integer CN states (right) for ccRCC4. \u003cstrong\u003e(h - j)\u003c/strong\u003e Heatmap of relative CN profiles (top) and genome density distributions of relative CNV scores (bottom) estimated by (h) epiAneuFinder, (i) Copy-scAT and (g) \u0026nbsp;inferCNV for ccRCC4.\u003c/p\u003e","description":"","filename":"Figure3.png","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/d0f1d1445d4048c8f68ba7df.png"},{"id":83220591,"identity":"d8aa6395-a3a9-45ad-9173-feb8473da580","added_by":"auto","created_at":"2025-05-21 10:18:26","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":202100,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eEvaluation of performance\u003c/strong\u003e. \u003cstrong\u003e(a) \u003c/strong\u003eInteger\u003cstrong\u003e \u003c/strong\u003eCN profiles from bulk WGS (left) the genome density distribution of absolute CN states (right) for ccRCC1. Significantly altered regions are highlighted in blue (loss) and orange (gain).\u003cstrong\u003e (b) \u003c/strong\u003eClonal integer CN profiles estimated by TeaCNV based on chromatin accessibility (left) and genome density distribution of estimated integer CN states (right) for ccRCC1. \u003cstrong\u003e(c, d)\u003c/strong\u003eHeatmaps of relative CN profiles (top) and genome density distributions of relative CNV scores (bottom) estimated by (c) epiAneuFinder and (d) Copy-scAT for ccRCC1. Blue and orange highlight regions identified by bulk WGS. \u003cstrong\u003e(e)\u003c/strong\u003eDispersion of estimated CNs for genomic segments with distinct absolute CNs in the bulk WGS profile. RMSE: root mean square error. \u003cstrong\u003e(f)\u003c/strong\u003e Precision, recall and \u003cstrong\u003e(g)\u003c/strong\u003e F1 score of CNVs detection by different methods. Each dot represents a distinct sample. P-values were calculated using Wilcoxon test. \u003cstrong\u003e(h)\u003c/strong\u003eDeviation of inferred CN profiles from bulk WGS.\u003c/p\u003e","description":"","filename":"Figure4.png","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/93e545bd8ca8b3aa4e4428c1.png"},{"id":83220587,"identity":"636769b9-3b33-4b37-a120-d17511f05335","added_by":"auto","created_at":"2025-05-21 10:18:26","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":1260369,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eClonal integer copy number profile estimation in solid tumors\u003c/strong\u003e. TeaCNV estimation for PDAC \u003cstrong\u003e(a and b),\u003c/strong\u003e BRCA \u003cstrong\u003e(c and d),\u003c/strong\u003e HNSCC \u003cstrong\u003e(e and g)\u003c/strong\u003e, CRC \u003cstrong\u003e(g and h)\u003c/strong\u003e and OV \u003cstrong\u003e(i and j)\u003c/strong\u003e samples. \u0026nbsp;(a, c, e, g, i) Heatmaps of relative CNs at single cell resolution annotated by clone. (b, d, f, h, j) Top: Clonal integer CN profiles. Bottom: UMAP plots of epithelial cells colored by activity scores of tumor epithelial markers and ATAC-seq signals from representative subclonal regions. Driver genes are highlighted.\u003c/p\u003e","description":"","filename":"Figure5.png","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/aba049b0f10ee293a793194e.png"},{"id":93477847,"identity":"8b935b2b-5236-442f-918f-81ea61679e75","added_by":"auto","created_at":"2025-10-14 09:34:19","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":3817836,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/f61d659c-7371-40c5-806d-f03fc64ca9c9.pdf"},{"id":83220592,"identity":"a3d15b52-5327-45fd-a2f8-b5a906e8dc82","added_by":"auto","created_at":"2025-05-21 10:18:26","extension":"pdf","order_by":6,"title":"","display":"","copyAsset":false,"role":"supplement","size":6734527,"visible":true,"origin":"","legend":"Supplementary Figures","description":"","filename":"SupplementaryFigures.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6609843/v1/bb3abd6f8ea6a0c6e88b302a.pdf"}],"financialInterests":"There is \u003cb\u003eNO\u003c/b\u003e Competing Interest.","formattedTitle":"TeaCNV: decoding tumor somatic absolute copy number and clonal architecture from single cell chromatin accessibility data","fulltext":[{"header":"Introduction","content":"\u003cp\u003eChromosomal copy number variations (CNVs), which drive ploidy changes in tumor genomes, constitute a hallmark of human cancer and critically shape tumor transcriptional landscapes\u003csup\u003e\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e\u003c/sup\u003e. Aneuploidy, defined by the presence of abnormal chromosome numbers, is strongly associated with tumor aggressiveness and poor clinical outcomes\u003csup\u003e\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e, \u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e, \u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e, \u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e\u003c/sup\u003e. Precise quantification of chromosomal absolute copy numbers (CNs) therefore represents a fundamental requirement for understanding cancer evolution and guiding clinical decision-making\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e, \u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003eWhile bulk and single cell whole-genome sequencing (scWGS) remain gold standards for CN profiling\u003csup\u003e\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e, \u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e, \u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e\u003c/sup\u003e, their clinical implementation is hindered by technical complexity, high cost, and limited scalability\u003csup\u003e\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e\u003c/sup\u003e. Furthermore, scWGS fails to capture tumor microenvironment dynamics. Although single-cell RNA sequencing (scRNA-seq) enables microenvironment characterization and tumor subpopulation analysis, existing CNV inference tools, such as inferCNV\u003csup\u003e\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e, \u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e, \u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e, \u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e, \u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e\u003c/sup\u003e, CopyKAT\u003csup\u003e\u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e\u003c/sup\u003e, Numbat\u003csup\u003e\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e\u003c/sup\u003e and HoneyBADGER\u003csup\u003e\u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e\u003c/sup\u003e, suffer from inherent limitations including transcriptional noise and post-transcriptional regulation frequently confound DNA level CNV estimates\u003csup\u003e\u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003eSingle cell ATAC sequencing (scATAC-seq), which profiles chromatin accessibility through transposase-accessible DNA sequencing, provides a more direct readout of genomic content than transcriptome-based approach\u003csup\u003e\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e, \u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e\u003c/sup\u003e. Current scATAC-based CNV detection methods face two limitations: dependence on matched normal DNA data and a focus on relative gain/loss detection rather than absolute CN solution\u003csup\u003e\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e, \u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e, \u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e\u003c/sup\u003e. These limit their utility in polyploid tumor analysis.\u003c/p\u003e \u003cp\u003eTo address these limitations, we developed TeaCNV, a computational framework that reconstructs absolute CN profiles and resolves clonal substructures directly from scATAC-seq data. Validated against bulk WGS in renal cell carcinomas, TeaCNV achieved 98.6% accuracy of CNVs calling and less than 10% error in absolute CN profiling while detecting subclonal CNVs undetectable by conventional sequencing. Application to multi-cancer types (renal, breast, pancreatic, head and neck, ovarian and colorectal cancers) revealed lineage-defining amplifications of oncogenic drivers and polyclonal architectures.\u003c/p\u003e"},{"header":"Results","content":"\u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003eMethodology to infer clonal absolute copy number profile from scATAC-seq\u003c/h2\u003e \u003cp\u003eWe developed TeaCNV, a quantitative framework for reconstructing clonal absolute CN profiles and clonal substructure from scATAC-seq data without requiring matched DNA sequencing references. The algorithm processes an input peak-by-cell matrix with cells pre-annotated as reference (normal) or inferred (tumor or epithelial) cells. Peaks are sorted by chromosomal coordinates. The outputs include clonal partitions of single cells and corresponding clonal absolute CN states for chromosomal segments (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eThe inference pipeline operates through four steps (\u003cb\u003eMethods\u003c/b\u003e). First is initial cell clustering. For each inferred cell, relative copy ratios are computed by normalizing peak accessibility signals against the median accessibility in matched reference cells. Dimension reduction via principal component analysis (PCA) on these ratios enables initial clustering of inferred cells into subpopulations (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea). Second is chromosomal breakpoint detection. Cells within each subpopulation are aggregated. Multi-scale chromosome breakpoints are identified using the pruned exact linear time (PELT) algorithm\u003csup\u003e\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e\u003c/sup\u003e. This generates segmentation boundaries for each subpopulation while accommodating scATAC-seq data sparsity (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eb). Third is joint ploidy-state optimization. For each subpopulation, absolute CNs of genomic segments are inferred through a combined optimization, including Nelder-Mead optimization to model expected relative copy ratios and CN state relationships, maximum likelihood estimation weighted by segment length and ratio variance, optimal integer CN determination via minimum Akaike information criterion (AIC, Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ec). Fourth is clonal architecture refinement. Subpopulations undergo pairwise comparison of their absolute and relative CN profiles. Subpopulations with concordant profiles are merged, followed by breakpoint recalculation and integer CN estimation (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ed). Last is the consensus profile generation including clone-specific absolute CN matrices across chromosomal segments and single-cell clonal assignments (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ee).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003eIdentifying Tumor clonal substructures with TeaCNV in clear cell renal cell carcinoma\u003c/h3\u003e\n\u003cp\u003eTo evaluate TeaCNV\u0026rsquo;s ability to reconstruct clonal architecture, we performed high-throughput 10X Genomics sequencing on four clear cell renal cell carcinoma (ccRCC) patients, generating two scATAC\u0026amp;RNA-seq co-assayed datasets and two scATAC-seq-only datasets (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ea \u003cb\u003eand b\u003c/b\u003e, \u003cb\u003eSupplementary Fig.\u0026nbsp;1)\u003c/b\u003e. Cell type annotation, based on chromatin accessibility signatures of canonical markers (\u003cem\u003ePTPRC\u003c/em\u003e, \u003cem\u003ePECAM1\u003c/em\u003e, \u003cem\u003eCD34\u003c/em\u003e, \u003cem\u003eCOL1A2\u003c/em\u003e, and \u003cem\u003eEPCAM\u003c/em\u003e) and subsequent UMAP visualization, identified 8,367 epithelial cells for CNV analysis (\u003cb\u003eSupplementary Figs.\u0026nbsp;2 and 3\u003c/b\u003e). Immune, endothelial and stromal cells were used as references.\u003c/p\u003e \u003cp\u003eUsing these non-epithelial cells as a reference, we applied TeaCNV to calculate single cell relative copy number ratios and reconstruct clonal architecture, while simultaneously inferring genome-wide absolute CNs (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ec-i). TeaCNV demonstrated robust reference compatibility, yielding consistent clonal architectures and integer CN profiles regardless of the reference cell type used (immune vs. endothelial, \u003cb\u003eSupplementary Figs.\u0026nbsp;4 and 5\u003c/b\u003e). In both scATAC\u0026amp;RNA-seq co-assayed and scATAC-seq-only datasets, TeaCNV accurately identified established ccRCC driver CNVs, including chr3p loss and chr5q gain\u003csup\u003e\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e\u003c/sup\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ec-i).\u003c/p\u003e \u003cp\u003eTeaCNV revealed distinct complexities in clonal substructures across these four ccRCC samples. While ccRCC1 displayed a monoclonal architecture (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ec\u003cb\u003e)\u003c/b\u003e, the other three cases displayed polyclonal expansion with a dominant subclone (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ed, f \u003cb\u003eand h\u003c/b\u003e). In ccRCC2, two clones shared similar epithelial scores but diverged in chromatin accessibility patterns driven by different absolute CNs at chr14 and chr16q (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ed \u003cb\u003eand e\u003c/b\u003e). TeaCNV resolved three clones in ccRCC3, with losses of chr8p and chr16q defining two minor clones that collectively accounted for 32.59% of the tumor cells (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ef \u003cb\u003eand g\u003c/b\u003e). A rare subclone (clone3: 5.19%) was further distinguished from clone 2 by co-amplification of chr8q and chr16p, accompanied by corresponding changes of chromatin accessibility (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eg). ccRCC4 contained four subclones stratified by absolute CNs of chr20q and chr22 (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003eh). The dominant clones (clone 1 and 4: 56.9%) exhibited a 4-copy on chr20q, while a minor clone (clone 3: 10.5%) exhibiting co-occurrence of a 4-copy of chr20q and a 1-copy of chr22 (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ei).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e\n\u003ch3\u003ePerformance evaluation of absolute copy number estimation using matched bulk WGS data\u003c/h3\u003e\n\u003cp\u003eTo validate the performance of absolute CNs inference using TeaCNV, we further analyzed three ccRCC patients (ccRCC3, ccRCC4 and ccRCC1) with sample-matched bulk WGS. The absolute CN profiles from bulk WGS data were considered as the ground truth. In ccRCC3, the integer CN profiles estimated by TeaCNV were highly concordant with the WGS profiles, particularly for the trunk events chr3p loss and chr7 gain (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ea \u003cb\u003eand b\u003c/b\u003e). As previously described, the subclonal chr8p and chr16q losses in clone 2 and 3, identified by TeaCNV, were further validated by WGS data, showing absolute CNs between 1 and 2 (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003ef, g, Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ea). Notably, the rare clonal events of 1-copy gain on both chr8q and chr16p (clone 3) were undetectable in bulk WGS data due to insufficient cell representation (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003eb \u003cb\u003eand Supplementary Fig.\u0026nbsp;6a\u003c/b\u003e).\u003c/p\u003e \u003cp\u003eUsing the scATAC\u0026amp;RNA-seq co-assayed data of ccRCC3, we applied epiAneufinder\u003csup\u003e\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e\u003c/sup\u003e and Copy-scAT\u003csup\u003e\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e\u003c/sup\u003e to chromatin accessibility data and inferCNV\u003csup\u003e\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e\u003c/sup\u003e to gene expression (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ec-e). These approaches either diluted signals from aberrant regions or misidentified neutral regions, leading to ambiguous distribution of inferred CNV scores for genomic segments with distinct CN states. Validation on another scATAC\u0026amp;RNA-seq co-assayed sample ccRCC4 (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003ef-j, \u003cb\u003eSupplementary Figs.\u0026nbsp;6b)\u003c/b\u003e, and the scATAC-seq-only sample ccRCC1 (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003ea-d\u003cb\u003e)\u003c/b\u003e demonstrated TeaCNV\u0026rsquo;s better performance in distinguishing chromosomal segments with distinct CN states.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eOverall, segments with differing CN states showed non-overlapping distribution of estimated integer CNs (average dispersion score: 0.97) in TeaCNV, compared to ambiguous distributions in alternative methods (average dispersion score: 0.06\u0026ndash;0.27, Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003ee). Across all samples, TeaCNV achieved 98.6% overall accuracy (precision: 93.0%, recall: 96.1%, F1 score: 0.92) for CNV event identification compared to the bulk WGS ground truth, outperforming existing methods (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003ef \u003cb\u003eand g\u003c/b\u003e). Quantitative analysis revealed the smallest error (average RMSE: 0.07) in absolute CN profiles estimated by TeaCNV (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eh, \u003cb\u003eMethods\u003c/b\u003e).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e\n\u003ch3\u003eIdentifying Tumor clonal substructures in solid tumors\u003c/h3\u003e\n\u003cp\u003eWe applied TeaCNV to published scATAC-seq datasets spanning pancreatic (PDAC), breast (BRCA), head and neck (HNSCC), colorectal (CRC), and ovarian (OV) cancers (\u003cb\u003eSupplementary Fig.\u0026nbsp;7\u003c/b\u003e)\u003csup\u003e\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e\u003c/sup\u003e. In the 1,694 epithelial cells from the PDAC sample, TeaCNV identified two clones, one of which was a diploid subpopulation exhibiting a reduced epithelial score (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea \u003cb\u003eand b\u003c/b\u003e). TeaCNV successfully identified driver events of PDAC, including \u003cem\u003eAKT2\u003c/em\u003e amplifications and \u003cem\u003eSMAD4\u003c/em\u003e deletions, which showed different chromatin accessibility patterns between the diploid and aneuploid subpopulations\u003csup\u003e\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003eIn the BRCA, HNSCC, CRC and OV samples, all epithelial cells were identified as aneuploid and exhibited polyclonal architectures (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ec-j). Driver genes of BRCA, including \u003cem\u003eMCL1\u003c/em\u003e, \u003cem\u003eMYC\u003c/em\u003e, \u003cem\u003eCSMD1\u003c/em\u003e, \u003cem\u003eRB1\u003c/em\u003e and \u003cem\u003eWWOX\u003c/em\u003e were identified as truncal CNV events across clones in the BRCA sample (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ec \u003cb\u003eand d\u003c/b\u003e)\u003csup\u003e\u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e\u003c/sup\u003e. Although these clones showed similar epithelial scores, they were distinguished by subclonal CNV events, such as chr11q and chr20q13.2. In particular, the oncogene \u003cem\u003eZNF217\u003c/em\u003e showed different CN states, driving changes in ATAC signals across clones (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ed). In HNSCC, driver genes including \u003cem\u003eLRP1B\u003c/em\u003e, \u003cem\u003eFGFR1\u003c/em\u003e and \u003cem\u003eMYC\u003c/em\u003e were identified as truncal CNVs. However, \u003cem\u003eSOX2, PDE4D\u003c/em\u003e and \u003cem\u003eKLF1\u003c/em\u003e were identified as subclonal events, showing distinct chromatin accessibility patterns across clones \u003cb\u003e(\u003c/b\u003eFig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ee \u003cb\u003eand f\u003c/b\u003e)\u003csup\u003e\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e\u003c/sup\u003e. Similarly, driver genes of CRC (\u003cem\u003eAPC\u003c/em\u003e, \u003cem\u003eMYC\u003c/em\u003e, \u003cem\u003eCCND2\u003c/em\u003e, \u003cem\u003eCDK8\u003c/em\u003e, \u003cem\u003eKLF5\u003c/em\u003e, \u003cem\u003eHNF4A\u003c/em\u003e, \u003cem\u003eSMAD4\u003c/em\u003e and \u003cem\u003eMACROD2\u003c/em\u003e, Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eg \u003cb\u003eand h\u003c/b\u003e)\u003csup\u003e\u003cspan citationid=\"CR33\" class=\"CitationRef\"\u003e33\u003c/span\u003e\u003c/sup\u003e and OV (\u003cem\u003eMCL1\u003c/em\u003e, \u003cem\u003eMECOM\u003c/em\u003e, \u003cem\u003eKRAS\u003c/em\u003e, \u003cem\u003eMETTL17\u003c/em\u003e and \u003cem\u003eCCNE1\u003c/em\u003e, Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ei \u003cb\u003eand j\u003c/b\u003e)\u003csup\u003e\u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e\u003c/sup\u003e were successfully identified as truncal CNVs in the corresponding samples. While clones across most samples exhibited conserved epithelial scores, divergent chromatin accessibility patterns at subclonal loci underscored the genomic heterogeneity.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e"},{"header":"Discussion","content":"\u003cp\u003eWe present TeaCNV, a computational framework that enables robust estimation of absolute CN profiles and clonal architecture directly from scATAC-seq data, without requiring matched bulk DNA sequencing. By aggregating sparse chromatin accessibility signals across epigenetically homogeneous cells, TeaCNV resolves genomic heterogeneity in single cell epigenomic datasets. Validation against bulk WGS in renal cell carcinoma demonstrated 98.6% accuracy in CNV detection and less than 10% error of integer CN profiling. Applied to six cancer types, TeaCNV revealed polyclonal architectures in 77.78% of cases, demonstrating its broad utility across solid tumors.\u003c/p\u003e \u003cp\u003eCompared to scRNA-based CNV inference tools, TeaCNV leverages scATAC-seq’s genome-wide coverage to avoid transcriptional confounding, while overcoming data sparsity through two innovations: iterative clustering of cells to aggregate sparse signal and likelihood-based modeling integrating segment length and variance for ploidy estimation. Unlike existing scATAC-based methods\u003csup\u003e\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e, \u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e\u003c/sup\u003e, TeaCNV quantifies absolute CN states rather than relative gains/losses. This precision enables accurate functional dissection of large-scale CNV regions\u003csup\u003e\u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e, \u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e, \u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e, \u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e\u003c/sup\u003e, such as distinguishing the distinct oncogenic roles of co-altered genes in single cell multiomics sequencing data.\u003c/p\u003e \u003cp\u003eThree considerations guide the application of TeaCNV. First, it performs optimally on subclones composed of sufficiently homogeneous cells, as sparse data limits the ability to call single-cell CNVs —a limitation shared by all scATAC-based approaches. Tumor cell lines or samples with high subclonal diversity (e.g., those containing numerous rare subpopulations) may yield unstable estimates. The minimal subclonal size we identified is 30 cells. Second, the current resolution of TeaCNV captures arm-level alterations. Focal events may be missed by sparse chromatin accessibility signals. Third, TeaCNV resolves total CN but not allele-specific CNV estimation. Efforts are underway to extend its functionality to enable allele-specific copy number estimation, which would further enhance its utility for analyzing complex tumor genomes.\u003c/p\u003e \u003cp\u003eTumor plasticity arises from aberrant activation of transcriptional programs caused by genetic and non-genetic mechanisms\u003csup\u003e\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e, \u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e, \u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e\u003c/sup\u003e. With advancements in single cell sequencing technologies, simultaneous profiling of the epigenome and transcriptome at single-cell resolution is now possible\u003csup\u003e\u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e, \u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e\u003c/sup\u003e. TeaCNV bridges this gap in single cell multiomics by connecting clonal genotypes to epigenomic and transcriptomic phenotypes. By enabling absolute CNV analysis in any scATAC-seq dataset, TeaCNV empowers systematic exploration of how genomic instability shapes epigenetic diversity—a critical factor underlying therapeutic resistance and metastatic progression.\u003c/p\u003e \u003cdiv id=\"Sec8\" class=\"Section2\"\u003e \u003cdiv id=\"Sec9\" class=\"Section3\"\u003e \u003c/div\u003e \u003c/div\u003e"},{"header":"Methods","content":"\u003ch2\u003eHuman specimens\u003c/h2\u003e\u003cp\u003eA total of four ccRCC patients were collected with informed consent, following approval from the Institutional Ethics Committee for Clinical Research approval at the Qilu Hospital of Shandong University, including two patients for scATAC-seq, two for scATAC\u0026amp;RNA-seq co-assays, and three patients for bulk WGS.\u003c/p\u003e\n\u003ch3\u003eExperimental methods\u003c/h3\u003e\n\u003cdiv id=\"Sec11\" class=\"Section2\"\u003e \u003ch2\u003eNuclei isolation from tissues\u003c/h2\u003e \u003cp\u003eTissues were dissected, snap-frozen in liquid nitrogen and stored at -80\u0026deg;C. Nuclei were isolated from the frozen tissues using a protocol designed for scATAC-seq (10X Genomics, CG000212 Rev B). For scATAC\u0026amp;RNA-seq co-assays, nuclei isolation was performed according to the manufacturer's recommended protocol (10X Genomics, CG000375 Rev C), which included RNase inhibitors in the buffer to prevent mRNA degradation during cell lysis.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003eDroplet library preparation and sequencing\u003c/h2\u003e \u003cp\u003eNuclei and barcoded beads were pooled together and loaded into the 10X Genomics system for single cell ATAC sequencing following the manufacturer\u0026rsquo;s instructions (10X Genomics, CG000496 Rev A for Chromium Next GEM Single Cell ATAC Library Kit v2). scATAC\u0026amp;RNA-seq co-assays libraries were prepared using the Chromium Next GEM Single Cell Multiome ATAC\u0026thinsp;+\u0026thinsp;Gene Expression kit, following a separate protocol (10X Genomics, CG000338 Rev F). The barcoded libraries were then pooled and sequenced on the Illumina NovaSeq 6000 system with the associated cells.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec13\" class=\"Section2\"\u003e \u003ch2\u003eLibrary Preparation and Sequencing for bulk WGS\u003c/h2\u003e \u003cp\u003eAbout 0.6 \u0026micro;g high-quality genomic DNA was sheared with Covaris LE220 Sonicator (Covaris) to about 350 bp. The library was constructed according to the protocol of KAPA Hyper Prep kit (Roche). First the fragmented DNA was purified using sample purification beads, and the product was repaired by the end and A base was added to the 3' end. Then, the adapters are ligated with the specific barcode sequence. The CleanNGS magnetic beads (CleanNA) were used to screen out incomplete connections and self-connecting products. Sequencing libraries were formed by PCR amplification using universal primers complementary to the adaptor sequences. Paired-end sequencing was performed using the NovaSeq 6000 S4 Reagent Kit v1.5 (300 cycles) on Illumina NovaSeq 6000 platform (Illumina, San Diego, USA) by Sequanta technologies (Shanghai, China).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003escATAC-seq and scATAC\u0026amp;RNA-seq data processing\u003c/h2\u003e \u003cp\u003eReads from scATAC-seq and scATAC\u0026amp;RNA-seq datasets were aligned to the GRCh38 (hg38) reference genome and quantified using the \u003cem\u003ecellranger-atac count\u003c/em\u003e (v.1.2.0) and the \u003cem\u003ecellranger-arc count\u003c/em\u003e (v.2.0) pipelines (10x Genomics), respectively. Peaks were identified using the MACS3 tool (v.3.0.0)\u003csup\u003e\u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e\u003c/sup\u003e through the \u003cem\u003eCallPeaks\u003c/em\u003e function in the Signac package (v.1.14.0, \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://github.com/timoast/signac\u003c/span\u003e\u003cspan address=\"https://github.com/timoast/signac\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e)\u003csup\u003e44\u003c/sup\u003e. Peaks on chromosomes X and Y, and those within the ENCODE Unified GRCh38 Blacklist regions, were removed using \u0026lsquo;blacklist_hg38_unified\u0026rsquo; in the \u003cem\u003esubsetByOverlaps\u003c/em\u003e function in Signac. The resulting sample-specific peak set was used to generate the peak-count matrix using \u003cem\u003eFeatureMatrix\u003c/em\u003e function in Signac package for downstream analyses.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec15\" class=\"Section2\"\u003e \u003ch2\u003eQuality control of single cell sequencing data\u003c/h2\u003e \u003cp\u003eQuality-control filtering of the scATAC-seq and scATAC\u0026amp;RNA-seq data was performed using functions from the Signac package. Filters applied for the cell inclusion were as follows: number of fragments in peaks\u0026thinsp;\u0026gt;\u0026thinsp;1,000; number of peaks in cell\u0026thinsp;\u0026gt;\u0026thinsp;2000; and enrichment-score for Tn5-integration events at transcriptional start sites\u0026thinsp;\u0026gt;\u0026thinsp;3.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec16\" class=\"Section2\"\u003e \u003ch2\u003eNormalization, dimensionality reduction, clustering and cell tying\u003c/h2\u003e \u003cp\u003eThe filtered peak-count matrix from scATAC-seq or scATAC\u0026amp;RNA-seq data was normalized using term frequency-inverse document frequency (TF-IDF) normalization implemented in the Signac package. This normalization accounts for variations in coverage across cells and peaks. The top 95% of peaks were selected as features for dimensionality reduction. We used the RunSVD function to perform singular value decomposition (SVD) on the normalized TF-IDF matrix, a method that is also known as latent semantic indexing (LSI) dimension reduction. The resulting 2:30 LSI components were used for nonlinear dimensionality reduction using the RunUMAP function from the Seurat package. Nuclei were clustered using a graph-based clustering approach implemented in Seurat using the 2:30 LSI components.\u003c/p\u003e \u003cp\u003eFor scATAC-seq data, we annotated cell types based on the activity of canonical cell type-specific markers, including epithelial cell (\u003cem\u003eEPCAM\u003c/em\u003e and \u003cem\u003eKRT\u003c/em\u003e family genes), immune cell (\u003cem\u003ePTPRC\u003c/em\u003e), endothelial cell (\u003cem\u003ePECAM1\u003c/em\u003e and \u003cem\u003eCD34\u003c/em\u003e), and stromal cell (\u003cem\u003eCOL1A2\u003c/em\u003e and \u003cem\u003eCOL1A3\u003c/em\u003e). For scATAC\u0026amp;RNA-seq data, we applied Seurat package (v.4.4.0)\u003csup\u003e\u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e\u003c/sup\u003e to the gene-count matrix for scaling, normalization and identification of highly variable genes for unsupervised cell clustering with default parameters. The elbow plot was generated with the ElbowPlot function of Seurat, and based on this, the number of significant principal components (PCs) was determined. In this study, the top 2,000 highly variable genes and the first 30 PCs identified by Seurat were used for unsupervised clustering analysis. We annotated cell types based on the expression of canonical cell type-specific markers.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec17\" class=\"Section2\"\u003e \u003ch2\u003eBulk WGS data processing and CNV analysis\u003c/h2\u003e \u003cp\u003eSequencing reads from bulk tumor tissue and matched normal tissues were aligned to the reference human genome (GRCh38) using Burrows\u0026ndash;Wheeler Aligner (BWA v0.7.17) software\u003csup\u003e\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e\u003c/sup\u003e to obtain the original mapping results stored in BAM format. SAMtools\u003csup\u003e\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e\u003c/sup\u003e were used to sort and index BAM files. The \u0026lsquo;runVarbin\u0026rsquo; module of copykit (v0.1.2)\u003csup\u003e\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e\u003c/sup\u003e was used to count the number of reads in 220-kb genomic bins defined by the GRCh38 genome assembly, with GC correction applied to the counts. We calculated the log ratio of tumor to normal tissue for each genomic bin using the bin counts, enabling CNV calling, segmentation, and the estimation of absolute CNs for each segment.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec18\" class=\"Section2\"\u003e \u003ch2\u003eTeaCNV Algorithm\u003c/h2\u003e \u003cdiv id=\"Sec19\" class=\"Section3\"\u003e \u003ch2\u003ePreprocessing and transformation of scATAC-seq data\u003c/h2\u003e \u003cp\u003eTo balance data sparsity and genomic coverage in the input peak-cell matrix, we exclude peaks detected in fewer than 5% of cells, ensuring that the retained number of peaks remains above 10,000. If this threshold is not met, the detection criterion is relaxed by progressively lowering the minimum proportion of cells required for peak detection. In solid tumors, immune cells or confident non-malignant cells serve as reference cells, while epithelial cells or candidate malignant cells are considered inferred cells. High-confidence reference cells are defined as those whose total peak read counts and the number of detected peaks fall within the 5th to 95th percentile range of all reference cells. The same percentile-based filtering is applied to inferred cells using metrics specific to this group. To mitigate bias in ploidy estimation due to chromosomal dropout, cells are excluded if more than 60% of its measurements for any chromosome have a value of zero.\u003c/p\u003e \u003cp\u003eWe construct two peak-cell matrices: \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{X}_{n\\times\\:m}\\)\u003c/span\u003e\u003c/span\u003e, where \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:n\\)\u003c/span\u003e\u003c/span\u003e represents the number of peaks (\u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:n\\:\\ge\\:\\:\\text{10,000}\\)\u003c/span\u003e\u003c/span\u003e) and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:m\\)\u003c/span\u003e\u003c/span\u003e represents the number of reference cells, and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{Y}_{n\\times\\:t}\\)\u003c/span\u003e\u003c/span\u003e where \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:t\\)\u003c/span\u003e\u003c/span\u003e is the number of inferred cells. To reduce the impact of extreme peak values, we cap the values in matrices X and Y to the range [0,4] by replacing all values exceeding 4 with 4. To correct for variations in sequencing depth across cells, we respectively normalize matrices \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:X\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:Y\\)\u003c/span\u003e\u003c/span\u003e as follow:\u003cdiv id=\"Equa\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equa\" name=\"EquationSource\"\u003e\n$$\\:{{X}^{{\\prime\\:}}}_{i,j}=\\frac{{X}_{i,j}}{\\sum\\:_{k}{X}_{k,j}}\\bullet\\:\\frac{1}{m}\\bullet\\:\\sum\\:_{j=1}^{m}\\sum\\:_{k=1}^{n}{X}_{k,j}$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equb\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equb\" name=\"EquationSource\"\u003e\n$$\\:{{Y}^{{\\prime\\:}}}_{i,j}=\\frac{{Y}_{i,j}}{\\sum\\:_{k}{Y}_{k,j}}\\bullet\\:\\frac{1}{t}\\bullet\\:\\sum\\:_{j=1}^{t}\\sum\\:_{k=1}^{n}{Y}_{k,j}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eFor the \u003cem\u003ei\u003c/em\u003e-th peak, we calculate the average value across reference cells using the normalized matrix \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{X}^{{\\prime\\:}}\\)\u003c/span\u003e\u003c/span\u003e:\u003cdiv id=\"Equc\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equc\" name=\"EquationSource\"\u003e\n$$\\:{\\stackrel{-}{X{\\prime\\:}}}_{i\u0026middot;}=average\\left[{X{\\prime\\:}}_{i,1\\dots\\:m}\\right]\\:$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eFor each peak in each inferred cell, we calculate the ratio relative to reference average signal to generate the ratio matrix \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}_{n\\times\\:t}\\)\u003c/span\u003e\u003c/span\u003e, where:\u003cdiv id=\"Equd\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equd\" name=\"EquationSource\"\u003e\n$$\\:{R}_{i,j}=\\frac{{Y{\\prime\\:}}_{i,j}}{{\\stackrel{-}{X{\\prime\\:}}}_{i\u0026middot;}}\\:(1\\le\\:i\\le\\:n,\\:1\\le\\:j\\le\\:t)$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003c/div\u003e \u003c/div\u003e \u003cdiv id=\"Sec20\" class=\"Section2\"\u003e \u003ch2\u003eAggregating cell sub-populations\u003c/h2\u003e \u003cp\u003eThe ratio matrix \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}_{n\\times\\:t}\\)\u003c/span\u003e\u003c/span\u003e is sorted by the chromosomal location of the analyzed peaks. To capture variations derived from chromosomal segments rather than particular peaks, we calculate the average ratio values within a genomic window of 5 peaks for each chromosome. Subsequently, we identify highly variable features and perform PCA based on the merged matrix using the Seurat package. The top 2,000 highly variable features and the first 50 PCs identified by Seurat are used for unsupervised clustering analysis, with the resolution parameter set to 1. This approach enables the preliminary classification of observed cells into distinct subgroups. For cells within the same subgroup, we aggregate ratios by averaging across cells to obtain subgroup-specific ratios as:\u003cdiv id=\"Eque\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Eque\" name=\"EquationSource\"\u003e\n$$\\:{R}^{s}=[{\\stackrel{-}{R}}_{1\u0026middot;}^{s},\\dots\\:,{\\stackrel{-}{R}}_{n\u0026middot;}^{s}]\\:$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equf\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equf\" name=\"EquationSource\"\u003e\n$$\\:{\\stackrel{-}{R}}_{i\u0026middot;}^{s}=average\\left[{R}_{i,kϵS}\\right]\\:$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eHere, \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\stackrel{-}{R}}_{i\u0026middot;}^{s}\\)\u003c/span\u003e\u003c/span\u003e represents the average ratio for the \u003cem\u003ei\u003c/em\u003e-th peak corresponding to subgroup S.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec21\" class=\"Section2\"\u003e \u003ch2\u003eGenome segmentation of cell subpopulation\u003c/h2\u003e \u003cp\u003eFor the \u003cem\u003es\u003c/em\u003e-th subgroup, genome segmentation is performed on \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}^{s}\\)\u003c/span\u003e\u003c/span\u003e organized according to chromosomal locations. Breakpoints are identified using the pruned exact linear time (PELT) algorithm from \u003cem\u003echangepoint\u003c/em\u003e R package\u003csup\u003e\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e\u003c/sup\u003e, along with an alternative method employing FPOP algorithm from \u003cem\u003erobseg\u003c/em\u003e R package\u003csup\u003e\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e\u003c/sup\u003e. To ensure the reliability of the segments, the chromosomal segments shorter than 2MB are excluded. The ratio for each chromosomal segment is estimated by calculating the median value of the peaks within the same segment. This results in the segmental ratio, which is considered as relative copy ratio for the \u003cem\u003es\u003c/em\u003e-th subgroup:\u003cdiv id=\"Equg\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equg\" name=\"EquationSource\"\u003e\n$$\\:{R}^{s}=\\left[{R}_{1}^{s},\\dots\\:,{R}_{l}^{s}\\right]$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equh\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equh\" name=\"EquationSource\"\u003e\n$$\\:\\:{R}_{j}^{s}=median\\left({\\stackrel{-}{R}}_{{j}_{k\\in\\:{segment}_{j}},\\bullet\\:}^{s}\\right)$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eWhere \u003cem\u003el\u003c/em\u003e is the number of segments.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec22\" class=\"Section2\"\u003e \u003ch2\u003eInferring absolute copy number\u003c/h2\u003e \u003cp\u003eThe algorithm for absolute CN estimation is based on the ABSOLUTE algorithm\u003csup\u003e\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e\u003c/sup\u003e, modified to accommodate single-cell sequencing data.\u003c/p\u003e \u003cp\u003eFor the estimation of each subgroup, the input segmental ratio \u003cem\u003eR\u003c/em\u003e with standard error \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\sigma\\:\\)\u003c/span\u003e\u003c/span\u003e consists of \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}_{j},\\:j\\in\\:\\{1,\\dots\\:,l\\}\\)\u003c/span\u003e\u003c/span\u003e, corresponding to a genomic fraction \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{w}_{j}\\)\u003c/span\u003e\u003c/span\u003e. Each \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}_{j}\\)\u003c/span\u003e\u003c/span\u003e is assumed to arise from one of the integer CN states in the set \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:Q=\\{\\text{1,2},\\dots\\:,Q\\}\\)\u003c/span\u003e\u003c/span\u003e with probabilities \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:p\\left({q}_{j}\\right),\\)\u003c/span\u003e\u003c/span\u003e where \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:q\\in\\:Q\\)\u003c/span\u003e\u003c/span\u003e. The observed \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}_{j}\\)\u003c/span\u003e\u003c/span\u003e is modeled as a mixture of \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:Q\\)\u003c/span\u003e\u003c/span\u003e Gaussian components located at \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\mu\\:=\\left\\{{\\mu\\:}_{q\\in\\:Q}\\right\\}\\)\u003c/span\u003e\u003c/span\u003e representing expected ratio of integer CN state.\u003cdiv id=\"Equi\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equi\" name=\"EquationSource\"\u003e\n$$\\:P\\left({R}_{j}|\\mu\\:,{\\sigma\\:}^{2},\\theta\\:,{w}_{j}\\right)={\\sum\\:}_{q\\in\\:Q}{w}_{j}\\bullet\\:P\\left({q}_{j}|{\\theta\\:}_{q}\\right)\\bullet\\:N\\left({\\mu\\:}_{q},{\\sigma\\:}^{2}\\right)$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eWhere, \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\theta\\:={\\{\\theta\\:}_{q\\in\\:Q}\\}\\)\u003c/span\u003e\u003c/span\u003e reprents the mixture weights, indicating the expected genomic fraction allocated to each CN state. Given limited knowledge about copy-state \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:q\\)\u003c/span\u003e\u003c/span\u003e, the distribution is chosen to have maximum entropy:\u003cdiv id=\"Equj\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equj\" name=\"EquationSource\"\u003e\n$$\\:P\\left({q}_{j}|{\\theta\\:}_{q}\\right)=\\frac{{e}^{-{\\theta\\:}_{q}\\bullet\\:{q}^{\\#}\\bullet\\:{w}_{j}}}{\\sum\\:_{k\\in\\:Q}{e}^{-{\\theta\\:}_{k}\\bullet\\:{k}^{\\#}\\bullet\\:{w}_{j}}}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eWhere \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{q}^{\\#}\\)\u003c/span\u003e\u003c/span\u003e indicates the order of \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:q\\)\u003c/span\u003e\u003c/span\u003e in the CN state set \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:Q,\\)\u003c/span\u003e\u003c/span\u003e beginning with 1. The unknown parameters \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\theta\\:={\\{\\theta\\:}_{q\\in\\:Q}\\}\\)\u003c/span\u003e\u003c/span\u003e are estimated by minimizing the loss function:\u003cdiv id=\"Equk\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equk\" name=\"EquationSource\"\u003e\n$$\\:\\underset{\\theta\\:}{\\widehat{\\theta\\:}={arg}\\text{min}}{\\left\\{\\sum\\:_{q\\in\\:Q}{\\left(\\sum\\:_{j=1}^{l}{w}_{j}\\bullet\\:P\\left({q}_{j}|{\\theta\\:}_{q}\\right)-{\\theta\\:}_{q}\\right)}^{2}\\right\\}}^{\\frac{1}{2}}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eThe full log-likelihood of the input data \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:R\\)\u003c/span\u003e\u003c/span\u003e is then computed as:\u003cdiv id=\"Equl\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equl\" name=\"EquationSource\"\u003e\n$$\\:logL\\left(R|\\mu\\:,\\sigma\\:,\\widehat{\\theta\\:},w\\right)=\\sum\\:_{j=1}^{l}logL\\left({R}_{j}|\\mu\\:,{\\sigma\\:}^{2},\\widehat{\\theta\\:},{w}_{j}\\right)$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere\u003cdiv id=\"Equm\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equm\" name=\"EquationSource\"\u003e\n$$\\:L\\left({R}_{j}|\\mu\\:,{\\sigma\\:}^{2},\\widehat{\\theta\\:},{w}_{j}\\right)={\\sum\\:}_{q\\in\\:Q}{w}_{j}\\bullet\\:P\\left({q}_{j}|{\\widehat{\\theta\\:}}_{q}\\right)\\bullet\\:N\\left({\\mu\\:}_{q},{\\sigma\\:}^{2}\\right)$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eThe unknown parameters \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\mu\\:={\\{\\mu\\:}_{q\\in\\:Q}\\}\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\sigma\\:}^{2}\\)\u003c/span\u003e\u003c/span\u003e are estimated using a combination of the Nelder-Mead optimization algorithm and maximum likelihood estimation.\u003c/p\u003e \u003cp\u003ePosterior probabilities are used to infer the copy-state indicators for each segment, with the absolute CN state of genomic segment \u003cem\u003ej\u003c/em\u003e is defined as the state corresponding to the maximum posterior probability:\u003cdiv id=\"Equn\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equn\" name=\"EquationSource\"\u003e\n$$\\:P\\left({\\widehat{q}}_{j}\\right)=p\\left({q}_{j}|{\\widehat{\\theta\\:}}_{q}\\right)\\bullet\\:\\frac{N\\left({R}_{j}|{\\widehat{\\mu\\:}}_{q},{\\widehat{\\sigma\\:}}^{2}\\right)}{L\\left({R}_{j}|\\widehat{\\mu\\:},{\\widehat{\\sigma\\:}}^{2},\\widehat{\\theta\\:},{w}_{j}\\right)}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cdiv id=\"Sec23\" class=\"Section3\"\u003e \u003ch2\u003eOptimizing copy number estimation\u003c/h2\u003e \u003cp\u003eEach subgroup corresponds to an estimated expected copy-ratio \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\mu\\:={\\{\\mu\\:}_{q\\in\\:Q}\\}\\)\u003c/span\u003e\u003c/span\u003e associated with the absolute CN state set \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:Q\\)\u003c/span\u003e\u003c/span\u003e. Ideally, the interval between two consecutive CN states, defined as \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\Delta\\:}={\\mu\\:}_{q}-{\\mu\\:}_{q-1}\\:(\\forall\\:q\\in\\:Q)\\)\u003c/span\u003e\u003c/span\u003e, should remain consistent across the values derived from the same subgroup. However, bias in the estimation of \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\mu\\:\\)\u003c/span\u003e\u003c/span\u003e may arise due to inhomogeneous observed values. To address this potential bias, we iterate through all possible values of \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\Delta\\:}\\)\u003c/span\u003e\u003c/span\u003e for each subgroup, and include \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{k}_{1}\\bullet\\:{\\Delta\\:}\\:({where\\:k}_{1}=2)\\)\u003c/span\u003e\u003c/span\u003e for the \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\Delta\\:}\u0026lt;0.3\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{k}_{2}\\bullet\\:{\\Delta\\:}\\:({k}_{2}=\\frac{1}{2})\\)\u003c/span\u003e\u003c/span\u003e for \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\Delta\\:}\u0026gt;0.6\\)\u003c/span\u003e\u003c/span\u003e to update expected copy-ratio \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\mu\\:={\\{\\mu\\:}_{q\\in\\:Q}\\}\\)\u003c/span\u003e\u003c/span\u003e. We calculate the Akaike Information Criterion (AIC) value for all candidate \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\mu\\:\\)\u003c/span\u003e\u003c/span\u003e and determine the optimal CN estimation for each corresponding subgroup based on the minimum AIC value.\u003c/p\u003e \u003cp\u003eSubsequently, the overall ploidy of the subgroup is then calculated based on the estimated integer CN states of segments weighted by their genomic fraction:\u003cdiv id=\"Equo\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equo\" name=\"EquationSource\"\u003e\n$$\\:ploidy=\\sum\\:_{j=1}^{l}{w}_{j}\\bullet\\:{\\widehat{q}}_{j}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003c/div\u003e \u003c/div\u003e \u003cdiv id=\"Sec24\" class=\"Section2\"\u003e \u003ch2\u003eScoring the confidence of cell subgroup with estimated copy number profile\u003c/h2\u003e \u003cp\u003eTo assess the confidence of identified homogeneous subgroups, we score each subgroup using a combination of mean squared error (MSE) between observed and expected copy ratio, and the proportion of genomic segments explained by the expected copy ratios.\u003c/p\u003e \u003cp\u003eFirst, the MSE between the observed and expected copy-ratios is calculated as:\u003cdiv id=\"Equp\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equp\" name=\"EquationSource\"\u003e\n$$\\:MSE=\\frac{1}{l}\\sum\\:_{i=1}^{l}\\left({R}_{i}-{\\mu\\:}_{i}\\right)$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cem\u003el\u003c/em\u003e is the total number of segments and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{\\mu\\:}_{i}\\)\u003c/span\u003e\u003c/span\u003e corresponds to the expected ratio of CN state for segment \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:i\\)\u003c/span\u003e\u003c/span\u003e.\u003c/p\u003e \u003cp\u003eNext, we calculate the proportion of genomic segments explained by the expected copy ratios:\u003cdiv id=\"Equq\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equq\" name=\"EquationSource\"\u003e\n$$\\:F=\\sum\\:_{q\\in\\:Q}\\left(\\sum\\:_{i\\in\\:\\left\\{l\\right|\\left|{R}_{l}-{\\mu\\:}_{q}\\right|\u0026lt;d\\}}\\frac{{length}_{i}}{L}\\right)\\:$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:L\\)\u003c/span\u003e\u003c/span\u003e is the total length of the genome. \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:d\\)\u003c/span\u003e\u003c/span\u003e is the allowed maximal difference between observed and expected ratio (default is 0.1), and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{length}_{i}\\)\u003c/span\u003e\u003c/span\u003e is the length of segment \u003cem\u003ei\u003c/em\u003e.\u003c/p\u003e \u003cp\u003eThe overall score, weighted by a parameter \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\beta\\:\\)\u003c/span\u003e\u003c/span\u003e, is defined as:\u003cdiv id=\"Equr\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equr\" name=\"EquationSource\"\u003e\n$$\\:score=F\\bullet\\:\\left(-logMSE\\right)\\bullet\\:\\beta\\:$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ewhere:\u003cdiv id=\"Equs\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equs\" name=\"EquationSource\"\u003e\n$$\\:\\beta\\:=\\frac{{\\delta\\:}_{{\\Delta\\:}}}{log\\left(1+l\\right)\\bullet\\:{e}^{ploidy-2}}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eand:\u003cdiv id=\"Equt\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equt\" name=\"EquationSource\"\u003e\n$$\\:{\\delta\\:}_{{\\Delta\\:}}=\\left\\{\\begin{array}{c}\\varDelta\\:,\\:if\\:\\varDelta\\:\u0026lt;0.3\\\\\\:1,\\:otherwise\\end{array}\\right.$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eThe weighted parameter \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\beta\\:\\)\u003c/span\u003e\u003c/span\u003e decreases with smaller interval between consecutive CN states, a larger number of segments, and greater difference in ploidy between estimated and diploidy.\u003c/p\u003e \u003cp\u003eTherefore, a higher score indicates a more reliable and highly homogeneous subgroup, reflecting that more genomic regions can be explained by the integer CN states with smaller MSE, appropriate interval between consecutive CN states and consistent genomic segmentation.\u003c/p\u003e \u003cdiv id=\"Sec25\" class=\"Section3\"\u003e \u003ch2\u003eSubclone partitioning\u003c/h2\u003e \u003cp\u003eWe define an initial adjacency matrix \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:A\\)\u003c/span\u003e\u003c/span\u003e to represent whether two subgroups can be merged into one clone.\u003cdiv id=\"Equu\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equu\" name=\"EquationSource\"\u003e\n$$\\:{A}_{i,j}=\\left\\{\\begin{array}{c}1,\\:i=j\\\\\\:0,i\\ne\\:j\\end{array}\\right.$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003ehere, a value of 1 indicates that subgroup \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:i\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:j\\)\u003c/span\u003e\u003c/span\u003e can be merged, while a value of 0 indicates they cannot.\u003c/p\u003e \u003cp\u003eTo update the adjacency matrix \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:A\\)\u003c/span\u003e\u003c/span\u003e, we compare the peak signals of subgroup in a pair-wise manner using the peak value at the group level, defined as follow:\u003cdiv id=\"Equv\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equv\" name=\"EquationSource\"\u003e\n$$\\:{{Y}^{{\\prime\\:}}}_{{s}_{i}}={\\left[\\stackrel{-}{Y}{{\\prime\\:}}_{1},\\dots\\:,\\stackrel{-}{Y}{{\\prime\\:}}_{n}\\right]}^{{s}_{i}}$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equw\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equw\" name=\"EquationSource\"\u003e\n$$\\:{{Y}^{{\\prime\\:}}}_{{s}_{j}}={\\left[\\stackrel{-}{Y}{{\\prime\\:}}_{1},\\dots\\:,\\stackrel{-}{Y}{{\\prime\\:}}_{n}\\right]}^{{s}_{j}}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eHere,\u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:\\:{{Y}^{{\\prime\\:}}}_{{s}_{i}}\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{{Y}^{{\\prime\\:}}}_{{s}_{j}}\\)\u003c/span\u003e\u003c/span\u003e represents the average values of peaks across cells from subgroup \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{i}\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{j}\\)\u003c/span\u003e\u003c/span\u003e, respectively. We then calculate the odds ratios between subgroup \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{i}\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{j}\\)\u003c/span\u003e\u003c/span\u003e:\u003cdiv id=\"Equx\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equx\" name=\"EquationSource\"\u003e\n$$\\:{R}_{{s}_{i},{s}_{j}}=\\frac{{{Y}^{{\\prime\\:}}}_{{s}_{i}}}{{{Y}^{{\\prime\\:}}}_{{s}_{j}}}\\:$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eNext, genome segmentation is performed based on \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{R}_{{s}_{i},{s}_{j}}\\)\u003c/span\u003e\u003c/span\u003e using the same approach described previously. For each segment, we compare the distribution of peak values (equivalent to copy-ratio) between subgroup \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{i}\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{j}\\)\u003c/span\u003e\u003c/span\u003e using a two-sided t-test, and adjust the resulting P values for multiple hypotheses through the Benjamini-Hochberg method. If the difference in any one segment is significant (adjusted p-value\u0026thinsp;\u0026lt;\u0026thinsp;0.05) and the estimated integer CN is different between \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{i}\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{s}_{j}\\)\u003c/span\u003e\u003c/span\u003e, we set:\u003cdiv id=\"Equy\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equy\" name=\"EquationSource\"\u003e\n$$\\:{A}_{i,j}=\\:{A}_{j,i}=0$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eOtherwise, we update the matrix as:\u003cdiv id=\"Equz\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equz\" name=\"EquationSource\"\u003e\n$$\\:{A}_{i,j}=\\:{A}_{j,i}=1$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eBased on the updated adjacent matrix, we obtain the final subclonal partitioning. For the final subclone substructure, we update the corresponding integer CN profile using the same method as previously described.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec26\" class=\"Section3\"\u003e \u003ch2\u003eApplication of other CNV inference approaches\u003c/h2\u003e \u003cp\u003eWe detected CNVs using epiAneufinder\u003csup\u003e\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e\u003c/sup\u003e and Copy-scAT from scATAC-seq of three ccRCC samples which have sample-matched bulk WGS data. The analysis was performed with each method\u0026rsquo;s default parameters to identify CNVs at single cell resolution.\u003c/p\u003e \u003cp\u003eIn addition, for two ccRCC samples with scATAC\u0026amp;RNA-seq, we employed inferCNV\u003csup\u003e\u003cspan citationid=\"CR13\" class=\"CitationRef\"\u003e13\u003c/span\u003e\u003c/sup\u003e on the gene expression data, using the same reference as TeaCNV. We followed the recommended parameters for 10X (denoise\u0026thinsp;=\u0026thinsp;TRUE, cutoff\u0026thinsp;=\u0026thinsp;0.1) and performed CNV calling using the \u0026lsquo;consensus\u0026rsquo; i6 HMM mode.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec27\" class=\"Section3\"\u003e \u003ch2\u003eEvaluating performance of copy number estimation\u003c/h2\u003e \u003cp\u003eTo evaluate the performance of CNV detection, we used CNVs identified by bulk WGS as the ground truth and calculated precision, recall, ACC and F1 score based on the overlap between the predicted and true CNVs. The genome was divided into 100-kb bins, excluding those overlapping breakpoints in either ground truth or inferred results. We defined true positive (TP) as bins with CNVs present in both the group truth and the inferred results, false positive (FP) as bins with CNVs present in the inferred results but absent in the ground truth, and false negative (FN) as bins with CNVs in the ground truth but not detected in the inferred results. Precision, recall, ACC and F1 score were then calculated using the following formulas:\u003cdiv id=\"Equaa\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equaa\" name=\"EquationSource\"\u003e\n$$\\:precision=\\:\\frac{TP}{TP+FP}$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equab\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equab\" name=\"EquationSource\"\u003e\n$$\\:recall=\\:\\frac{TP}{TP+FN}$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equac\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equac\" name=\"EquationSource\"\u003e\n$$\\:ACC=\\:\\frac{TP+TN}{TP+TN+FP+FN}$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equad\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equad\" name=\"EquationSource\"\u003e\n$$\\:F1=\\:\\frac{2\\bullet\\:precision\\bullet\\:recall}{precision+recall}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eTeaCNV reports CNVs at clonal level, we averaged the precision, recall and F1 score across the subclones derived from each sample. For epiAneufinder and Copy-scAT, CNV events were considered detected at the bulk level if they were identified in a specific proportion of cells, ranging from 10\u0026ndash;100% in 10% increments. The average precision, recall, and F1 score were then calculated across these varying cell proportion thresholds.\u003c/p\u003e \u003cp\u003eTo evaluate the accuracy of inferred CN profiles, we calculated the deviation of estimations from WGS using the root mean square error (RMSE) metric:\u003cdiv id=\"Equae\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equae\" name=\"EquationSource\"\u003e\n$$\\:RMSD=\\:\\sqrt{\\frac{1}{n}\\bullet\\:\\sum\\:_{i=1}^{n}{\\left({CN}_{i}-{CN}_{i}^{true}\\right)}^{2}}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eHere, \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{CN}_{i}\\)\u003c/span\u003e\u003c/span\u003e represents the inferred integer CN for bin \u003cem\u003ei\u003c/em\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{CN}_{i}^{true}\\)\u003c/span\u003e\u003c/span\u003e is the true integer CN derived from WGS. Because epiAneuFinder and Copy-scAT do not report integer CNs directly, we centered the inferred CNV scores to align with integer CNs for bins with distinct absolute CN states in WGS data. The definition of \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{CN}_{i}\\)\u003c/span\u003e\u003c/span\u003e was as follow:\u003cdiv id=\"Equaf\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equaf\" name=\"EquationSource\"\u003e\n$$\\:{CN}_{i}=\\left\\{\\begin{array}{c}{CN}_{i}\\:\\:\\:for\\:TeaCNV\\\\\\:{CNV\\:score}_{i}-\\underset{j\\in\\:\\left\\{j|{{CN}}_{j}^{{true}}={{CN}}_{i}^{{true}}\\right\\}}{\\text{average}}\\left({CNV\\:score}_{j}\\right)+{CN}_{i}^{true}\\:\\:\\:for\\:other\\:methods\\end{array}\\right.$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eFor TeaCNV, \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{CN}_{i}\\)\u003c/span\u003e\u003c/span\u003e is taken directly from the inferred integer CNs. For other methods, \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{CN}_{i}\\)\u003c/span\u003e\u003c/span\u003e is adjusted by subtracting the average of the CNV scores from bins with the same true CN state, followed by adding the true absolute CN. This adjustment allows for a comparison of inferred profiles against the ground truth.\u003c/p\u003e \u003cp\u003eIn theory, the inferred values should be confined to discrete ranges defined by distinct absolute CN states. Ideally, the density peaks in the distribution of inferred CNs for genomic regions with the same CN state should not overlap with those for regions with different CN states. The dispersion in the distribution of inferred value from the genomic regions with distinct CN states reflects the confidence of estimation. To quantify the dispersion of CN estimations, we use the following formula:\u003cdiv id=\"Equag\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equag\" name=\"EquationSource\"\u003e\n$$\\:dispersion=\\:\\sum\\:_{s}\\left[{\\int\\:}_{A}^{}max\\left({p}_{s}\\left(x\\right),{p}_{s+1}\\left(x\\right)\\right)dx\\right]$$\u003c/div\u003e\u003c/div\u003e\u003cdiv id=\"Equah\" class=\"Equation\"\u003e\u003cdiv format=\"TEX\" class=\"mathdisplay\" id=\"FileID_Equah\" name=\"EquationSource\"\u003e\n$$\\:A=\\left\\{x\\right|{p}_{s}\\left(x\\right)\\ne\\:0\\:\\\u0026amp;\\:{p}_{s+1}\\left(x\\right)\\ne\\:0\\}$$\u003c/div\u003e\u003c/div\u003e\u003c/p\u003e \u003cp\u003eWhere \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:s\\)\u003c/span\u003e\u003c/span\u003e represents unique absolute CN states from WGS data. \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{p}_{s}\\left(x\\right)\\)\u003c/span\u003e\u003c/span\u003e is the probability density of the inferred values for the genomic regions with integer CN \u003cem\u003es\u003c/em\u003e in WGS data. For TeaCNV, the inferred value is output integer CN profiles. For other methods, the inferred value is the estimated CNV score. The interval \u003cem\u003eA\u003c/em\u003e denotes the range of integration that includes all \u003cem\u003ex\u003c/em\u003e values where both \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{p}_{s}\\left(x\\right)\\)\u003c/span\u003e\u003c/span\u003e and \u003cspan class=\"InlineEquation\"\u003e\u003cspan class=\"mathinline\"\u003e\\(\\:{p}_{s+1}\\left(x\\right)\\)\u003c/span\u003e\u003c/span\u003e are non-zero. A larger dispersion indicates a greater difference in inferred values between distinct absolute CN states, suggesting stronger confidence in the separation of these states.\u003c/p\u003e \u003c/div\u003e \u003c/div\u003e "},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eData availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe scATAC-seq, scATAC\u0026amp;RNA-seq aco-assays and bulk WGS validation data of ccRCC can be accessed through the link with token (https://zenodo.org/records/14190637?preview=1\u0026amp;token=eyJhbGciOiJIUzUxMiJ9.eyJpZCI6Ijg1YTczZDU4LTlkNjQtNGZkYS05ZjY3LTA3M2JhMjcwOTY2MiIsImRhdGEiOnt9LCJyYW5kb20iOiJlYWRjZmU3ZDU1M2I1ZjMxYTNlNGQ2MThiYzFiNDhhNCJ9.i4vMudpQCvXLnj3kkcBhJpLoS9GD_6rb-tQW0MjAq7AG1H2RHG84ZWWm5jhc2cmd5fMl-QOcmCXUk_7MVjib_A). scATAC-seq data of BRCA, PDAC, HNSCC, CRC and OV were downloaded from Gene Expression Omnibus (GEO) with accession number GSE240822.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCode availability\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe TeaCNV algorithm is available at https://github.com/ShaojunLab/TeaCNV. The analysis scripts used to reproduce results included in the paper are available at https://github.com/ShaojunLab/TeaCNVanalysis.\u003c/p\u003e\u003cp\u003e \u003ch2\u003eCompeting interests\u003c/h2\u003e \u003cp\u003eThe authors declare no competing interests.\u003c/p\u003e \u003c/p\u003e\u003ch2\u003eAuthor contributions\u003c/h2\u003e \u003cp\u003eY. W., S. Z. and F. W. formulated the study and the overall approach. Y. W. developed and implemented the computational algorithm with contribution form Y. Z.. X. Y. provide clinical samples. Y. D., Y. C. and Z. C. performed single-cell sequencing experiments and bulk WGS. Y.W. performed the analysis with the help from M. Z., X. W. and H. L.. Y. W., F. W. and S. Z. drafted the manuscript. All authors provided suggestions and corrections on the manuscript text.\u003c/p\u003e\u003ch2\u003eAcknowledgements\u003c/h2\u003e \u003cp\u003eThis work was supported by the National Natural Science Foundation of China (grants 32170666 to S. Zhang, grants 32400528 to Y. Wang; grant 32270699 to F. Wang), Guangdong Pearl River Program (grant 2021QN02Y180 to F. Wang).\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eEsteller M, Dawson MA, Kadoch C, Rassool FV, Jones PA, Baylin SB (2024) The Epigenetic Hallmarks of Cancer. Cancer Discov 14:1783\u0026ndash;1809\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChakravarthi BV, Nepal S, Varambally S (2016) Genomic and Epigenomic Alterations in Cancer. Am J Pathol 186:1724\u0026ndash;1735\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFriedlander ML, Hedley DW, Taylor IW (1984) Clinical and biological significance of aneuploidy in human tumours. J Clin Pathol 37:961\u0026ndash;974\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eClark GM, Dressler LG, Owens MA, Pounds G, Oldaker T, McGuire WL (1989) Prediction of relapse or survival in patients with node-negative breast cancer by DNA flow cytometry. N Engl J Med 320:627\u0026ndash;633\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKallioniemi OP, Punnonen R, Mattila J, Lehtinen M, Koivula T (1988) Prognostic significance of DNA index, multiploidy, and S-phase fraction in ovarian cancer. Cancer 61:334\u0026ndash;339\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChoma D, Daures JP, Quantin X, Pujol JL (2001) Aneuploidy and prognosis of non-small-cell lung cancer: a meta-analysis of published data. Br J Cancer 85:14\u0026ndash;22\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBielski CM et al (2018) Genome doubling shapes the evolution and prognosis of advanced cancers. Nat Genet 50:1189\u0026ndash;1195\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCarter SL et al (2012) Absolute quantification of somatic DNA alterations in human cancer. Nat Biotechnol 30:413\u0026ndash;421\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSalehi S et al (2021) Clonal fitness inferred from time-series modelling of single-cell cancer genomes. Nature 595:585\u0026ndash;590\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMinussi DC et al (2021) Breast tumours maintain a reservoir of subclonal diversity during expansion. Nature 592:302\u0026ndash;308\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCasasent AK et al (2018) Multiclonal Invasion in Breast Tumors Identified by Topographic Single Cell Sequencing. Cell 172:205\u0026ndash;217e212\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLei Y et al (2021) Applications of single-cell sequencing in cancer research: progress and perspectives. J Hematol Oncol 14:91\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePatel AP et al (2014) Single-cell RNA-seq highlights intratumoral heterogeneity in primary glioblastoma. Science 344:1396\u0026ndash;1401\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTirosh I et al (2016) Dissecting the multicellular ecosystem of metastatic melanoma by single-cell RNA-seq. Science 352:189\u0026ndash;196\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTirosh I et al (2016) Single-cell RNA-seq supports a developmental hierarchy in human oligodendroglioma. Nature 539:309\u0026ndash;313\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eVenteicher AS et al (2017) Decoupling genetics, lineages, and microenvironment in IDH-mutant gliomas by single-cell RNA-seq. Science 355\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePuram SV et al (2017) Single-Cell Transcriptomic Analysis of Primary and Metastatic Tumor Ecosystems in Head and Neck Cancer. Cell 171:1611\u0026ndash;1624e1624\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGao R et al (2021) Delineating copy number and clonal substructure in human tumors from single-cell transcriptomes. Nat Biotechnol 39:599\u0026ndash;608\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGao T et al (2023) Haplotype-aware analysis of somatic copy number variations from single-cell transcriptomes. Nat Biotechnol 41:417\u0026ndash;426\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFan J et al (2018) Linking transcriptional and genetic tumor heterogeneity through allele analysis of single-cell RNA-seq data. Genome Res 28:1217\u0026ndash;1227\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHou Y et al (2016) Single-cell triple omics sequencing reveals genetic, epigenetic, and transcriptomic heterogeneity in hepatocellular carcinomas. Cell Res 26:304\u0026ndash;319\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBuenrostro JD et al (2015) Single-cell chromatin accessibility reveals principles of regulatory variation. Nature 523:486\u0026ndash;490\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGrandi FC, Modi H, Kampman L, Corces MR (2022) Chromatin accessibility profiling by ATAC-seq. Nat Protoc 17:1518\u0026ndash;1552\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRamakrishnan A et al (2023) epiAneufinder identifies copy number alterations from single-cell ATAC-seq data. Nat Commun 14:5846\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNikolic A et al (2021) Copy-scAT: Deconvoluting single-cell chromatin accessibility of genetic subclones in cancer. Sci Adv 7:eabg6045\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWu CY et al (2021) Integrative single-cell analysis of allele-specific copy number alterations and chromatin accessibility in cancer. Nat Biotechnol 39:1259\u0026ndash;1269\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKillick R, Fearnhead P, Eckley IA (2012) Optimal Detection of Changepoints With a Linear Computational Cost. J Am Stat Assoc 107:1590\u0026ndash;1598\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRicketts CJ et al (2018) The Cancer Genome Atlas Comprehensive Molecular Characterization of Renal Cell Carcinoma. Cell Rep 23:313\u0026ndash;326e315\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTerekhanova NV et al (2023) Epigenetic regulation during cancer transitions across 11 tumour types. Nature 623:432\u0026ndash;441\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCancer Genome Atlas Research Network. Electronic address aadhe, Cancer Genome Atlas Research N. Integrated Genomic Characterization of Pancreatic Ductal Adenocarcinoma. Cancer Cell 32, 185\u0026ndash;203 e113 (2017)\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCancer Genome Atlas N (2012) Comprehensive molecular portraits of human breast tumours. Nature 490:61\u0026ndash;70\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCancer Genome Atlas N (2015) Comprehensive genomic characterization of head and neck squamous cell carcinomas. Nature 517:576\u0026ndash;582\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCancer Genome Atlas N (2012) Comprehensive molecular characterization of human colon and rectal cancer. Nature 487:330\u0026ndash;337\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCancer Genome Atlas Research N (2011) Integrated genomic analyses of ovarian carcinoma. Nature 474:609\u0026ndash;615\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRice AM, McLysaght A (2017) Dosage sensitivity is a major determinant of human copy number variant pathogenicity. Nat Commun 8\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang F et al (2019) Integrated transcriptomic-genomic tool Texomer profiles cancer tissues. Nat Methods 16:401\u0026ndash;404\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShi H, Williams MJ, Satas G, Weiner AC, McPherson A, Shah SP (2024) Allele-specific transcriptional effects of subclonal copy number alterations enable genotype-phenotype mapping in cancer cells. Nat Commun 15:2482\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShlyakhtina Y, Moran KL, Portal MM (2021) Genetic and Non-Genetic Mechanisms Underlying Cancer Evolution. Cancers (Basel) 13\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMarine JC, Dawson SJ, Dawson MA (2020) Non-genetic mechanisms of therapeutic resistance in cancer. Nat Rev Cancer 20:743\u0026ndash;756\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTellez-Gabriel M, Ory B, Lamoureux F, Heymann MF, Heymann D (2016) Tumour Heterogeneity: The Key Advantages of Single-Cell Analysis. Int J Mol Sci 17\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLee J, Hyeon DY, Hwang D (2020) Single-cell multiomics: technologies and data analysis methods. Exp Mol Med 52:1428\u0026ndash;1442\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBaysoy A, Bai Z, Satija R, Fan R (2023) The technological landscape and applications of single-cell multi-omics. Nat Rev Mol Cell Biol 24:695\u0026ndash;713\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang Y et al (2008) Model-based analysis of ChIP-Seq (MACS). Genome Biol 9:R137\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStuart T, Srivastava A, Madad S, Lareau CA, Satija R (2021) Single-cell chromatin state analysis with Signac. Nat Methods 18:1333\u0026ndash;1341\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHao Y et al (2021) Integrated analysis of multimodal single-cell data. Cell 184:3573\u0026ndash;3587e3529\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi H, Durbin R (2009) Fast and accurate short read alignment with Burrows-Wheeler transform. Bioinformatics 25:1754\u0026ndash;1760\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDanecek P et al (2021) Twelve years of SAMtools and BCFtools. \u003cem\u003eGigaScience\u003c/em\u003e 10\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMinussi DC et al (2022) Resolving clonal substructure from single cell genomic data using CopyKit. \u003cem\u003ebioRxiv\u003c/em\u003e, 2022.2003.2009.483497\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKillick R, Eckley IA (2014) changepoint: AnRPackage for Changepoint Analysis. J Stat Softw 58\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFearnhead P, Rigaill G (2018) Changepoint Detection in the Presence of Outliers. J Am Stat Assoc 114:169\u0026ndash;183\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":true,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-6609843/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-6609843/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eAccurate inference of absolute copy numbers beyond simple gains and losses from single-cell chromatin accessibility (scATAC-seq) data remains challenging, thereby obscuring the distinction between genetic and epigenetically driven oncogenic dependencies. Here, we present TeaCNV, a computational framework that reconstructs clonal absolute copy number profiles and tumor clonal architectures from scATAC-seq data without matched DNA baselines. Validated against bulk whole-genome sequencing in renal cell carcinomas, TeaCNV resolved subclonal absolute copy number profiles with less than 10% error and detected copy number variations (CNVs) with 98.6% accuracy, outperforming existing methods. Applied to six cancer types including renal, breast, pancreatic, head and neck, colorectal, and ovarian cancers, TeaCNV delineated polyclonal architectures and revealed distinct chromatin accessibility patterns driven by CNVs in key driver genes, including \u003cem\u003eAKT2, ZNF217\u003c/em\u003e and \u003cem\u003eSOX2.\u003c/em\u003e By enabling absolute copy number profiling and clonal deconvolution from epigenomic assays, TeaCNV bridges critical gaps in studying oncogenic dependencies and genotype-phenotype relationships at single-cell resolution.\u003c/p\u003e","manuscriptTitle":"TeaCNV: decoding tumor somatic absolute copy number and clonal architecture from single cell chromatin accessibility data","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-05-21 10:18:21","doi":"10.21203/rs.3.rs-6609843/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"cbe8a240-291e-410f-a58a-d82871719014","owner":[],"postedDate":"May 21st, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[{"id":48757831,"name":"Biological sciences/Computational biology and bioinformatics/Computational models"},{"id":48757832,"name":"Biological sciences/Computational biology and bioinformatics/Genome informatics"}],"tags":[],"updatedAt":"2025-10-14T09:26:08+00:00","versionOfRecord":[],"versionCreatedAt":"2025-05-21 10:18:21","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-6609843","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-6609843","identity":"rs-6609843","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}