是一个R包，它从scRNA-seq数据的原始计数矩阵开始，通过将肿瘤微环境的非恶性细胞与恶性细胞分离来自动对活检中存在的细胞进行分类，并表征这些恶性细胞的克隆结构。它识别具有不同拷贝数结构的细胞亚群，并报告每个亚群的特异性和共享改变。该工具的目的是通过允许以非常简单且完全无监督的方式执行整个分析来自动化整个分析。
1、简介 DOI: 10.1038/s41467-023-36790-9 链接：https://pubmed.ncbi.nlm.nih.gov/36841879/ 1.1 工作流程 SCEVAN算法是一个旨在从原始scRNA-seq数据中自动解析肿瘤克隆亚结构的“端到端”计算管道（图1）。整个流程始于输入的原始基因表达计数矩阵。
在对低质量细胞和基因进行初步过滤后，算法首先会自动识别出一小部分可信度高的非恶性细胞，并将它们的表达谱作为基准，计算出所有细胞相对于这个正常基线的基因表达差异矩阵。为了在降噪的同时保留真实的拷贝数变化边界，该差异矩阵会经过一个边缘保留非线性扩散滤波器（edge-preserving nonlinear diffusion filter）进行平滑处理。处理后的平滑矩阵是后续分析的核心。
SCEVAN会采用其关键创新(一个源于计算机视觉领域的联合分割算法(joint segmentation algorithm))对所有细胞的表达谱进行同步分割。该算法基于一个核心假设：同一克隆内的所有细胞共享相同的拷贝数断点，因此可以整合所有细胞的信号来更准确地识别这些断点。分割完成后，算法会自动将所有细胞分为恶性与非恶性两大类，其依据是包含初始可信正常细胞最多的细胞簇被定义为非恶性群体。
随后，仅针对恶性细胞群体，SCEVAN会进一步通过聚类来识别不同的亚克隆。每个亚克隆的拷贝数谱会再经过一次独立的精细分割，最终将各基因组区域划分为“缺失”、“丢失”、“中性”、“增加”或“扩增”五种拷贝数状态。流程的最后，SCEVAN会进行全面的下游分析，包括区分主干、共享和克隆特异性的拷贝数变异，并构建克隆演化树，从而完整地描绘出肿瘤的克隆结构和演化关系。
1.2 合成数据上进行恶性细胞分类 研究团队希望在可控的条件下定量评估SCEVAN算法的准确性，于是他们在模拟数据（synthetic data）上进行了一系列严格的基准测试。利用scDesign2工具，一个基于真实的多发性骨髓瘤数据集，生成了500个具有已知细胞身份的模拟scRNA-seq数据集。
这些数据集被设计用来模拟真实的肿瘤样本，包含了不同的肿瘤纯度和细胞组成，并且为了进一步增加挑战性，研究者还使用SPLATTER工具向数据中人为添加了三种不同水平的dropout噪声。在对这总计超过32万个模拟细胞的测试中，SCEVAN的表现全面优于当前先进的CopyKAT算法。研究人员使用F1分数作为准确性的衡量标准，结果显示SCEVAN在所有噪声水平下都取得了显著更高的分数。
例如，在无额外噪声的情况下，SCEVAN的F1分数为0.948，而CopyKAT为0.798；在最高噪声水平下，SCEVAN的F1分数仍有0.824，而CopyKAT则降至0.726。这一系列在模拟数据上的测试结果，有力地证明了SCEVAN在自动识别恶性细胞方面具有卓越的准确性和鲁棒性。1.3 真实数据中的恶性细胞分类性能验证 作者将SCEVAN和CopyKAT应用到了一个包含106个真实肿瘤样本的大型公开数据集中进行现实场景验证。
这些样本覆盖了胶质母细胞瘤（GBM）、头颈癌（HNSCC）和结直肠癌三种不同的癌症类型和多种测序技术，总细胞数超过93,000个。在此次测试中，作者以原始研究中经过专家手动注释的细胞身份作为“金标准”，来评判两种自动化算法的准确性。结果再次显示了SCEVAN的优越性（图2）：它在63%的样本中取得了比CopyKAT更高的分类分数，而CopyKAT仅在23%的样本中占优。
从总体性能上看，SCEVAN在所有样本中的总平均F1分数达到了0.90，远高于CopyKAT的0.63。作者也指出，在肿瘤细胞数量极少（1-15个）的样本中，SCEVAN的准确性会有所下降，但这一系列在多样化真实数据集上的测试结果，仍旧有力地证明了SCEVAN能够准确地从scRNA-seq数据中区分肿瘤细胞和正常细胞。
1.4 合成数据中的拷贝数断点分割性能验证 在证明其卓越的细胞分类能力后，作者进一步对SCEVAN的核心技术之变分联合分割算法（variational joint segmentation algorithm）进行验证。在识别拷贝数变异（CNV）断点上的准确性进行了深入验证。研究团队在模拟数据（synthetic data）上进行了测试以便于进行可控的定量评估。
他们构建了包含克隆性和亚克隆性拷贝数变异的复杂数据集，其中变异的位置和强度都是已知的。在与CopyKAT及其他分割算法的直接比较中，SCEVAN在所有测试场景下都展现出显著更高的准确性。无论是在F1分数还是在精确率-召回率曲线下面积（AUC）这两个关键指标上，SCEVAN的表现都全面胜出，这表明其联合分割策略能更有效地从嘈杂的单细胞表达数据中识别出真实的拷贝数边界。
1.5 参考数据的分割准确性验证 那么如何在真实世界数据中验证其核心分割算法的准确性？作者使用了一批同时拥有scRNA-seq数据和作为“金标准”的群体DNA测序（WGS）数据的肿瘤样本。通过比较各种算法从scRNA-seq中推断出的CNV谱与DNA测序得到的真实CNV谱之间的皮尔逊相关性（Pearson correlation），来评估算法的准确性。
结果再次证实了SCEVAN的优越性（图3）：在胶质母细胞瘤（GBM）和多发性骨髓瘤（MM）两个数据集中，SCEVAN推断出的CNV谱与“金标准”的相关性均显著高于inferCNV和CopyKAT。例如，在GBM数据集中，SCEVAN的平均相关性为0.57，而inferCNV为0.44，CopyKAT仅为-0.03。此外，作者还通过一系列对照实验证明，SCEVAN的这种优势主要来源于其更优秀的分割算法，而非仅仅是细胞分类的准确性。
并且，该算法对于参考正常细胞的轻微错误分类表现出高度的鲁棒性（robustness），这在处理复杂的真实世界样本时是一个巨大的优势。1.6. 计算效率评估 除了准确性，计算效率也是衡量一个算法实用性的关键指标，尤其是在处理日益庞大的单细胞数据集时。研究团队对SCEVAN算法的运行速度进行了评估，并与CopyKAT和inferCNV进行了比较。
他们指出，SCEVAN之所以高效，关键在于其核心的分割步骤采用了一种贪心区域增长算法（greedy region-growing algorithm）。在对算法不同阶段的计时测试中，SCEVAN在初始的恶性细胞分类步骤中，速度是竞争对手的2至7倍。在更为耗时的分割步骤中，SCEVAN的优势更为显著：在一个胶质母细胞瘤（GBM）数据集上，它比CopyKAT快2倍，比inferCNV快5倍。
特别是在处理细胞数量庞大的10x Genomics数据集（如多发性骨髓瘤样本）时，SCEVAN的效率优势被进一步放大，其速度分别比inferCNV和CopyKAT快了11倍和19倍。这些结果有力地证明，SCEVAN不仅在精度上超越了同类工具，在计算效率上也表现卓越，使其特别适合应用于当前的大规模单细胞数据分析。
1.7 胶质母细胞瘤中的应用：揭示肿瘤内部的克隆异质性 作者将其应用于分析胶质母细胞瘤（Glioblastoma, GBM）这一具有高度瘤内异质性（intratumoral heterogeneity）的恶性脑瘤，以便于展示SCEVAN在真实生物学问题中的应用价值。
在一个名为MGH105的样本中，SCEVAN成功地从scRNA-seq数据中识别出了四个亚克隆（four subclones），而这些亚克隆是标准scRNA-seq分析流程无法分辨的，其存在此前仅通过单细胞DNA甲基化这一更直接的基因组学技术才得以证实，这突显了SCEVAN算法的超高分辨率。在另一个名为BT1160的样本中，SCEVAN不仅识别出三个亚克隆并构建了它们的演化树（phylogenetic tree）（图4）。
更重要的是，它还将这种由拷贝数变异（CNA）定义的遗传异质性与功能异质性关联了起来：分析显示，这三个亚克隆分别对应了GBM的不同代谢亚型（神经元型、线粒体型和增殖/祖细胞型）。此外，通过对增殖型亚克隆特有的基因扩增区域进行分析，算法还锁定了一个潜在的驱动基因UBE2T。最后，分析也揭示了关键抑癌基因（如PTEN和CDKN2A）的缺失是亚克隆级别的，即只存在于部分而非全部的肿瘤细胞中（图5）。
这些应用实例充分证明，SCEVAN不仅能精确地识别肿瘤亚克隆，还能将遗传结构与细胞功能状态联系起来，从而从scRNA-seq数据中发掘出具有重要生物学和临床意义的见解。1.8 多区域胶质母细胞瘤的克隆演化 作者还将其应用于一个多区域胶质母细胞瘤（GBM）样本的分析中，这样可以进一步展示SCEVAN在解析复杂肿瘤生物学问题上的能力。
由于单个肿瘤活检可能无法完全代表整个肿瘤的异质性，因此对来自同一肿瘤不同空间位置的样本进行分析，这对于理解肿瘤演化至关重要。在本文中，作者分析了一位患者的七个活检样本，其中两个来自肿瘤外周（tumor periphery），五个来自肿瘤核心（tumor core）。通过对每个活检样本进行独立的克隆分析，SCEVAN成功地重构了该肿瘤内部复杂的空间克隆结构，并推断出了一个演化树（evolutionary tree）（图6）。
该演化树清晰地显示，来自肿瘤外周的克隆形成了一个独立的演化分支，它们缺少存在于核心区域克隆中的4号和8号染色体扩增。此外，一个位于2号染色体上的扩增事件在外周区域是克隆性的（所有肿瘤细胞都存在），但在部分核心区域则转变为亚克隆性的（仅部分肿瘤细胞存在），这揭示了肿瘤在空间扩张过程中的动态演化过程。这一应用实例表明，SCEVAN能够有效整合多区域scRNA-seq数据，从而描绘出肿瘤的地理克隆图谱并推断其演化路径。
1.9 原发灶与转移淋巴结的克隆结构 作者利用该算法比较了头颈癌（HNSCC）患者的原发肿瘤（primary tumor）与其对应的淋巴结转移灶（lymph node metastasis）之间的克隆结构关系。在分析的多个病例中，大部分情况下转移灶的克隆结构与原发灶高度相似，显示出很高的相关性，这表明转移细胞直接来源于原发灶的主要克隆。
然而，在一个名为HNSCC5的特定病例中，SCEVAN揭示了一个重要的演化差异：原发灶中存在一个7号染色体的扩增，而这个扩增在淋巴结转移灶中却消失了（图7）。进一步分析发现，这个区域包含一个名为GPNMB的关键基因，该基因在转移灶中也相应地表现出下调。作者指出，GPNMB在其他癌症中通常与促进肿瘤生长和转移有关，因此它在转移过程中的丢失可能是一个值得深入研究的演化事件。
这个案例充分说明，SCEVAN可以被用来精确比较原发灶和转移灶的克隆差异，从而为研究肿瘤转移过程中的克隆演化和筛选过程提供有力的工具。1.10 小结 综上所述，SCEVAN算法通过变分方法成功从scRNA-seq数据中解卷肿瘤克隆子结构，准确区分恶性细胞，并揭示了肿瘤异质性、克隆演化和微环境交互的关键机制。在脑肿瘤应用中，它识别了如UBE2T等潜在治疗靶点，并阐明了PTEN和CDKN2A等抑癌基因的克隆状态。
然而，该算法的局限在于依赖非整倍体假设，不适合如白血病或儿童癌症等基因组变异较少的肿瘤类型。未来研究可扩展到更多技术平台和数据集，结合多组学整合，进一步验证其在临床转化中的鲁棒性。总体而言，SCEVAN为肿瘤异质性研究注入了新活力，有望推动精准医学的进步。
2、安装R包 library (devtools) install_github ( "miccec/yaGST" ) install_github ( "AntonioDeFalco/SCEVAN" ) 3、加载R包 # 清空环境变量 rm ( list = ls ) gc # used (Mb) gc trigger (Mb) max used (Mb) # Ncells 9781130 522.4 14930285 797.4 14930285 797.4 # Vcells 17998057 137.4 644490092 4917.1 1006972481 7682.6 # 展示系统时间 Sys.time # [1] "2025-08-29 18:10:40 CST" library (SCEVAN) # 设置工作目录 setwd ( "/data/09_scRNA-seq_cnv/" ) 4、单样本分析 从胶质母细胞公共数据集 (GSE131928) 加载MGH106样本的scRNA-seq数据。
单个调用（pipelineCNA）能够执行整个关于克隆结构分类和表征的分析。+ count_mtx：具有基因行（允许使用基因符号或 Ensembl ID）和细胞列（列中包含细胞名称）的计数矩阵。
+ sample：保存结果的样本名称（可选） + par_cores：运行管道所需的内核数量（可选 - 默认值为 20） + norm_cells：用于作为可靠正常细胞的可能已知正常细胞向量（可选） + FIXED_NORMAL_CELLS：如果 norm_cell 向量用于作为固定参考，则为 TRUE，如果您仅对克隆结构感兴趣而不关心正常/肿瘤分类，则为 FALSE（默认值为 FALSE） + SUBCLONES：如果对克隆结构感兴趣，则为 TRUE，否则仅为恶性和非恶性细胞分类感兴趣，则为 FALSE（可选 - 默认值为 TRUE） + beta_vega：指定用于分割的 beta 参数，beta 值越高分割越粗粒化。
（可选 - 默认值为 0.5） + ClonalCN：从所有肿瘤细胞中获取克隆 CN 分析结果（可选） + plotTree：绘制系统发育树（可选 - 默认值为 FALSE） + AdditionalGeneSets：正常细胞类型附加签名的列表（可选） + SCEVANsignatures：如果只想使用 AdditionalGeneSets 中指定的签名，则为 FALSE（可选 - 默认值为 TRUE） + organism : 要进行分析的生物体（可选 - “mouse” 或 “human” - 默认为 “human”） 4.1 加载数据 # 示例：results <- pipelineCNA(count_mtx) # 加载数据 load ( "/data/09_scRNA-seq_cnv/data/MGH106_data.RData" ) 4.2 执行分析 single_result <- SCEVAN :: pipelineCNA ( count_mtx, # 单样本每个基因在每个细胞的表达矩阵 sample = "MGH106" , # 保存结果的样本名称 par_cores = 4 , # 分析所用核心数 SUBCLONES = TRUE , # 如果对克隆结构感兴趣，则为 TRUE，否则仅为恶性和非恶性细胞分类感兴趣，则为 FALSE（可选 - 默认值为 TRUE） plotTree = TRUE ) # 绘制系统发育树（可选 - 默认值为 FALSE） # 如果上一步出现报错提示：Error in f(...) : Graphics API version mismatch # 可以尝试在Rstudio 运行 install.packages("ragg") # 但是这时候可能还是会因为一些可能的权限原因提示无法安装 # 可以去linux端，执行sudo R ,启用管理员权限后再执行install.packages("ragg")，将ragg包的版本更新至‘1.4.0’ # 默认保存到当前工作目录下面的output文件夹，可以挪动结果目录到指定位置 unlink ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/" , recursive = TRUE ) dir.create ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample" , recursive = TRUE ) file.rename ( "output/" , "result_SCEVAN/single_sample/" ) 4.3 结果查看 # 返回结果包含每个细胞的分类（肿瘤/正常），它是否被用作确信的正常细胞，以及它属于哪个克隆亚群 head (single_result) # class confidentNormal subclone # MGH106-P7-E10 tumor <NA> 1 # MGH106-P2-B01 tumor <NA> 1 # MGH106-P2-B12 tumor <NA> 1 # MGH106-P2-H07 tumor <NA> 2 # MGH106-P2-E09 tumor <NA> 1 # MGH106-P2-E04 tumor <NA> 1 dim (single_result) # [1] 200 3 # 正常和肿瘤细胞数目 table (single_result $ class) # # normal tumor # 64 136 # 被用作确信的正常的正常细胞数目 table (single_result $ confidentNormal) # # yes # 30 # 肿瘤细胞每个克隆亚型的数目 table (single_result $ subclone) # # 1 2 3 # 58 57 21 拷贝数变异矩阵的热图，包含非恶性和恶性细胞的分类。
文件路径为：/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106heatmap.png 肿瘤细胞克隆亚群拷贝数变异矩阵的热图 (heatmap_subclones.png)，文件路径为：/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106heatmap_subclones.png 从亚克隆档案推断出的克隆树，文件路径为：/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106CloneTree.png 每个亚群体中存在的变化的紧凑共识图，文件路径为：/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106consensus.png OncoPrint样图突出显示特定改变、亚克隆之间的共享改变或克隆改变。
这些片段分为五种预定义的拷贝数状态之一：缺失、丢失、中性、增益或扩增。文件路径为：/data/09_scRNA-seq_cnv/img/MGH106OncoHeat2.png 特定改变的 DE 分析 (DEchr*_subclones.png)。火山图是从属于特定改变的基因的差异表达分析中获得的。
# 展示出指定目录下所有文件名称包含了MGH106-DEchr字符的文件 files <- list.files ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/" , pattern = "MGH106-DEchr" ) files # [1] "MGH106-DEchr1-155687960-160286348_subclones.png" # [2] "MGH106-DEchr1-28736621-36305357_subclones.png" # [3] "MGH106-DEchr1-74198212-92841924_subclones.png" # [4] "MGH106-DEchr11-87037844-118925608_subclones.png" # [5] "MGH106-DEchr11-95768942-117204782_subclones.png" # [6] "MGH106-DEchr12-109116595-120534434_subclones.png" # [7] "MGH106-DEchr12-280129-6961316_subclones.png" # [8] "MGH106-DEchr12-28133249-75511636_subclones.png" # [9] "MGH106-DEchr12-30709552-76853696_subclones.png" # [10] "MGH106-DEchr13-90060247-102876001_subclones.png" # [11] "MGH106-DEchr15-45587123-55508234_subclones.png" # [12] "MGH106-DEchr15-81159575-99135593_subclones.png" # [13] "MGH106-DEchr16-46955362-58283836_subclones.png" # [14] "MGH106-DEchr16-51135975-65126112_subclones.png" # [15] "MGH106-DEchr19-47713343-54102692_subclones.png" # [16] "MGH106-DEchr2-217730-27988087_subclones.png" # [17] "MGH106-DEchr2-46580937-74473322_subclones.png" # [18] "MGH106-DEchr21-17788967-46665124_subclones.png" # [19] "MGH106-DEchr21-44328944-46665124_subclones.png" # [20] "MGH106-DEchr22-29788608-49827512_subclones.png" # [21] "MGH106-DEchr3-113532296-139357223_subclones.png" # [22] "MGH106-DEchr3-183487531-198043720_subclones.png" # [23] "MGH106-DEchr4-108650584-127840733_subclones.png" # [24] "MGH106-DEchr4-153152342-183040119_subclones.png" # [25] "MGH106-DEchr4-68646630-77819615_subclones.png" # [26] "MGH106-DEchr8-73290267-93741017_subclones.png" # [27] "MGH106-DEchr9-88388419-98015943_subclones.png" # 提取第一个文件路径 first_file <- paste0 ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/" ,files[ 1 ]) 例如我这里展示第一个文件 亚克隆的富集分析 (pathway Analysis_subclones*.png)，利用GSEA获得的每个亚克隆的REACTOME途径活性。
# 展示出指定目录下所有文件名称包含了MGH106pathwayAnalysis字符的文件 files <- list.files ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/" , pattern = "MGH106pathwayAnalysis" ) files # [1] "MGH106pathwayAnalysis_subclones1.png" # [2] "MGH106pathwayAnalysis_subclones2.png" # [3] "MGH106pathwayAnalysis_subclones3.png" # 提取第一个文件路径 first_file <- paste0 ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/" ,files[ 1 ]) 这里同样以第一个文件为例 5、多样本分析 从胶质母细胞公共数据集 (GSE131928) 加载了三个样本 MGH102、MGH104 和 MGH105的scRNA-seq数据。
单个调用（multiSampleComparisonClonalCN）允许比较多个样品的克隆谱。
+ listCountMtx：要分析的样品的原始计数矩阵的命名列表 + listNormCells：要分析的每个样品的正常细胞列表 + analysisName：分析的名称（可选） + organism：可选“鼠”或“人”，默认为“人” + par_cores：内核数（默认为20） 5.1 加载数据 # 示例：multiSampleComparisonClonalCN(listCountMtx) # load(url("https://www.dropbox.com/s/esqvnltucdqajg1/listCountMtx.RData?raw=1")) # 此步骤可能出现读取失败的情况，因此我已经下载到了本地，可以直接从本地加载 load ( "/data/09_scRNA-seq_cnv/data/listCountMtx.RData" ) 5.2 执行分析 # 运行多样本分析步骤,约10min multi_result <- SCEVAN :: multiSampleComparisonClonalCN ( listCountMtx, # 多样本表达矩阵list analysisName = "all" , # 分析的名称 organism = "human" , # 分析物种 par_cores = 4 ) # 核心数 # [1] " raw data - genes: 23686 cells: 355" # [1] "1) Filter: cells > 200 genes" # [1] "2) Filter: genes > 10% of cells" # [1] "12813 genes past filtering" # [1] "3) Annotations gene coordinates" # [1] "found 0 confident non malignant cells" # [1] "11490 genes annotated" # [1] "4) Filter: genes involved in the cell cycle" # [1] "10928 genes past filtering " # [1] "5) Filter: cells > 5genes per chromosome " # [1] "6) Log Freeman Turkey transformation" # [1] "A total of 355 cells, 10928 genes after preprocessing" # [1] "7) Measuring baselines (pure tumor - synthetic normal cells)" # [1] "8) Smoothing data" # [1] "10) Adjust baseline" # [1] "11) plot heatmap" # [1] "found 355 tumor cells" # [1] "time classify tumor cells: 3.66527685721715" # [1] " raw data - genes: 23686 cells: 634" # [1] "1) Filter: cells > 200 genes" # [1] "2) Filter: genes > 10% of cells" # [1] "12189 genes past filtering" # [1] "3) Annotations gene coordinates" # [1] "found 30 confident non malignant cells" # [1] "10986 genes annotated" # [1] "4) Filter: genes involved in the cell cycle" # [1] "10481 genes past filtering " # [1] "5) Filter: cells > 5genes per chromosome " # [1] "6) Log Freeman Turkey transformation" # [1] "A total of 634 cells, 10481 genes after preprocessing" # [1] "7) Measuring baselines (confident normal cells)" # [1] "8) Smoothing data" # [1] "9) Segmentation (VegaMC)" # [1] "10) Adjust baseline" # [1] "11) plot heatmap" # [1] "found 427 tumor cells" # [1] "time classify tumor cells: 5.70456975301107" # [1] " raw data - genes: 23686 cells: 312" # [1] "1) Filter: cells > 200 genes" # [1] "2) Filter: genes > 10% of cells" # [1] "12359 genes past filtering" # [1] "3) Annotations gene coordinates" # [1] "found 30 confident non malignant cells" # [1] "11153 genes annotated" # [1] "4) Filter: genes involved in the cell cycle" # [1] "10641 genes past filtering " # [1] "5) Filter: cells > 5genes per chromosome " # [1] "6) Log Freeman Turkey transformation" # [1] "A total of 312 cells, 10641 genes after preprocessing" # [1] "7) Measuring baselines (confident normal cells)" # [1] "8) Smoothing data" # [1] "9) Segmentation (VegaMC)" # [1] "10) Adjust baseline" # [1] "11) plot heatmap" # [1] "found 216 tumor cells" # [1] "time classify tumor cells: 2.71259829600652" # [1] "Phylogenic tree cannot be plotted (ggtree requires dplyr version 1.0.5) try downgrading dprly to 1.0.5 or upgrading ggtree" # 同样如果上一步出现报错提示：Error in f(...) : Graphics API version mismatch # 可以尝试在Rstudio 运行 install.packages("ragg") # 但是这时候可能还是会因为一些可能的权限原因提示无法安装 # 可以去linux端，执行sudo R ,启用管理员权限后再执行install.packages("ragg")，将ragg包的版本更新至‘1.4.0’ # 默认保存到当前工作目录下面的output文件夹，可以挪动结果目录到指定位置 unlink ( "/data/09_scRNA-seq_cnv/result_SCEVAN/multi_sample" , recursive = TRUE ) dir.create ( "/data/09_scRNA-seq_cnv/result_SCEVAN/multi_sample" , recursive = TRUE ) file.rename ( "output" , "result_SCEVAN/multi_sample" ) # [1] TRUE 5.3 结果查看 class (single_result) # [1] "data.frame" # 与单样本不同的是，多样本分析生成的结果对象是一个list,总共包含了两个部分 class (multi_result) # [1] "list" # 包含了三个样本 names (multi_result[[ 1 ]]) # [1] "MGH104" "MGH105" "MGH102" # 每个样本包含了鉴定的正常/肿瘤细胞信息，以及当前样本细胞它是否被用作确信的正常细胞 head (multi_result[[ 1 ]][[ 1 ]]) # class confidentNormal # MGH104-P5-A02 tumor <NA> # MGH104-P5-A04 tumor <NA> # MGH104-P5-A05 tumor <NA> # MGH104-P5-A06 tumor <NA> # MGH104-P5-A07 tumor <NA> # MGH104-P5-A08 tumor <NA> # 包含了三个样本以及他们共享的位点 names (multi_result[[ 2 ]]) # [1] "all_MGH104" "all_MGH105" "all_MGH102" # [4] "all_shared" "all_shareSubshared" # 包含了染色体编号，染色体起点和终点以及拷贝数变异数值 head (multi_result[[ 2 ]][[ 1 ]]) # Chr Start End Alteration segm.mean # 5 3 47228026 49535618 2 0.174473 # 7 3 61561569 111665750 -2 -0.151862 # 11 4 152618632 168928457 -2 -0.225926 # 24 11 62832342 66958376 2 0.151061 # 27 12 53251251 57095476 2 0.262316 # 31 14 66507407 75282230 2 0.155383 # 所有样本都共有的变异位点 head (multi_result[[ 2 ]][[ 4 ]]) # Chr Start End Alteration segm.mean # 13 5 92151 147782848 2 0.040585 # 16 7 192969 159144957 2 0.088789 # 21 10 134465 133421307 -2 -0.114114 # 43 22 29767369 50801309 -2 -0.159207 # 部分样本共有的变异位点 head (multi_result[[ 2 ]][[ 5 ]]) # Chr Start End Alteration segm.mean sh_sub Mean # all_share2.6 2 1646536241817413 2 0.069152 2-3 0 # all_share2.11 5 148825245 181261139 -1 -0.098875 2-3 0 # all_share1.14 6 291630 28093664 -2 -0.099298 1-3 0 # all_share1.19 9 89318506 122913330 -2 -0.219364 1-3 0 # all_share2.23 11 78068304 134411918 1 0.098385 2-3 0 # all_share2.40 19 43646095 55442863 -2 -0.187000 2-3 0 三个样本拷贝数变异对比热图，文件路径为：/data/09_scRNA-seq_cnv/result_SCEVAN/multi_sample/all_compareClonalCN.png 三个样本之间的共享改变或不同克隆改变。
这些片段分为五种预定义的拷贝数状态之一：缺失、丢失、中性、增益或扩增。
文件路径为：/data/09_scRNA-seq_cnv/img/allOncoHeat2.png 6、联合Seurat对象 6.1 创建Seurat对象 library (Seurat) seurObj <- Seurat :: CreateSeuratObject (count_mtx, meta.data = single_result) # 或者可以使用下面方法增加metadata seurObj <- Seurat :: AddMetaData (seurObj, metadata = single_result) seurObj # An object of class Seurat # 23686 features across 200 samples within 1 assay # Active assay: RNA (23686 features, 0 variable features) # 2 layers present: counts, data 6.2 将拷贝数变异结果加入Seurat对象 # 加载前面的分析结果 # 基因注释结果 load ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106_count_mtx_annot.RData" ) # 拷贝数变异矩阵 load ( "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106_CNAmtx.RData" ) # 分别对应了基因名称、染色体编号、染色体对应起点和重点和symbol名称 count_mtx_annot[ 1 : 5 , 1 : 5 ] # seqnames start end gene_id gene_name # ENSG00000227232 1 14404 29570 ENSG00000227232 WASH7P # ENSG00000225880 1 826206 827522 ENSG00000225880 LINC00115 # ENSG00000230368 1 868071 876903 ENSG00000230368 FAM41C # ENSG00000188976 1 944204 959309 ENSG00000188976 NOC2L # ENSG00000187608 1 1001138 1014541 ENSG00000187608 ISG15 # 查看数据维度 dim (count_mtx_annot) # [1] 10566 5 # 列分别对应了细胞名称 CNA_mtx_relat[ 1 : 5 , 1 : 5 ] # MGH106-P7-E10 MGH106-P2-B01 MGH106-P2-B12 MGH106-P2-H07 MGH106-P2-E09 # [1,] 0.197378 -0.09338299 0.01316769 -0.07037455 -0.08075804 # [2,] 0.197378 -0.09338299 0.01316769 -0.07037455 -0.08075804 # [3,] 0.197378 -0.09338299 0.01316769 -0.07037455 -0.08075804 # [4,] 0.197378 -0.09338299 0.01316769 -0.07037455 -0.08075804 # [5,] 0.197378 -0.09338299 0.01316769 -0.07037455 -0.08075804 # 查看数据维度 dim (CNA_mtx_relat) # [1] 10566 200 # 从CNA_mtx_relat中筛选出染色体3上特定区域的数据，计算每列的均值，然后按照[email protected]的行名对结果进行重新排序和命名，最后将其转换为数据框 chr3 <- apply (CNA_mtx_relat[count_mtx_annot $ seqnames == 3 & count_mtx_annot $ start >= 158644278 & count_mtx_annot $ end <= 194498364 ,], 2 , mean) head (chr3) # MGH106-P7-E10 MGH106-P2-B01 MGH106-P2-B12 MGH106-P2-H07 MGH106-P2-E09 # -0.073484354 0.143933614 0.006989981 -0.029519201 -0.074089242 # MGH106-P2-E04 # 0.022922624 chr3 <- chr3[ rownames (seurObj @ meta.data)] names (chr3) <- rownames (seurObj @ meta.data) chr3 <- as.data.frame (chr3) head (chr3) # chr3 # MGH106-P7-E10 -0.073484354 # MGH106-P2-B01 0.143933614 # MGH106-P2-B12 0.006989981 # MGH106-P2-H07 -0.029519201 # MGH106-P2-E09 -0.074089242 # MGH106-P2-E04 0.022922624 # 将chr3结果加入Seurat对象中 seurObj <- AddMetaData (seurObj, metadata = chr3) 6.3 降维分群 seurObj <- NormalizeData (seurObj , normalization.method = "LogNormalize" , scale.factor = 10000 ) # 寻找高变基因 seurObj <- FindVariableFeatures (seurObj , selection.method = "vst" , nfeatures = 2000 ) # 归一化数据 seurObj <- ScaleData (seurObj , features = VariableFeatures (seurObj)) # 计算线粒体含量 seurObj [[ "percent.mt" ]] <- PercentageFeatureSet (seurObj , pattern = "^MT-" ) seurObj <- ScaleData (seurObj , vars.to.regress = "percent.mt" ) [["RNA"]]@scale.data seurObj <- RunPCA (seurObj , features = VariableFeatures ( object = seurObj )) ElbowPlot (seurObj ) # 寻找邻近值 seurObj <- FindNeighbors (seurObj , dims = 1 : 15 ) seurObj <- RunUMAP (seurObj , dims = 1 : 15 ) 6.4 可视化 # 展示正常/肿瘤的降维图 DimPlot (seurObj , reduction = "umap" , group.by = "class" , label = TRUE ) # 展示chr3染色体选取区域拷贝数变异VlnPlot # 将[email protected]$subclone的NA对应的正常细胞名称调整为normal seurObj @ meta.data $ subclone[ is.na (seurObj @ meta.data $ subclone)] <- "normal" seurObj @ meta.data $ subclone[seurObj @ meta.data $ subclone == 1 ] <- "clone1" seurObj @ meta.data $ subclone[seurObj @ meta.data $ subclone == 2 ] <- "clone2" seurObj @ meta.data $ subclone[seurObj @ meta.data $ subclone == 3 ] <- "clone3" table (seurObj @ meta.data $ subclone) # # clone1 clone2 clone3 normal # 58 57 21 64 VlnPlot (seurObj, "chr3" , group.by = "subclone" ) # 热图 plotCNA_withAnnotCells ( SampleName = "MGH106" , # 分析时用的样本名称 outputPATH= "/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/" , # 输出结果目录，和前面保持一致 metadata = seurObj @ meta.data, # 包含分类信息的列表 COLUMNS_TO_PLOT = c ( "class" , "subclone" )) # 需要展示分类变量名称 # [1] "class" # [1] "subclone" # png # 2 热图的文件路径为：/data/09_scRNA-seq_cnv/result_SCEVAN/single_sample/MGH106heatmap_with_annotation.png Sys.time # [1] "2025-08-29 18:38:18 CST" save.image ( "result_SCEVAN/SCEVAN.rdata" )