scFAST-seq 突变分析：拟时序轨迹上的突变动态展示

作者: SeekGene

时长: 9 分钟

字数: 1.8k 字

更新: 2026-02-28

阅读: 0 次

全序列转录组 Notebooks 拟时序分析

右上角选择 monocle2(R)环境

确保已经完成了 monocle2 高级分析，这里直接读取 monocle2 的分析结果

library(monocle)

output

Loading required package: Matrix

Loading required package: Biobase

Loading required package: BiocGenerics


Attaching package: ‘BiocGenerics’


The following objects are masked from ‘package:stats’:

IQR, mad, sd, var, xtabs


The following objects are masked from ‘package:base’:

anyDuplicated, aperm, append, as.data.frame, basename, cbind,
colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
Position, rank, rbind, Reduce, rownames, sapply, setdiff, sort,
table, tapply, union, unique, unsplit, which.max, which.min


Welcome to Bioconductor

Vignettes contain introductory material; view with
'browseVignettes()'. To cite Bioconductor, see
'citation("Biobase")', and for packages 'citation("pkgname")'.


Loading required package: ggplot2

Loading required package: VGAM

Loading required package: stats4

Loading required package: splines

Loading required package: DDRTree

Loading required package: irlba

读取 monocle2 分析结果

修改下方的 rds 路径

路径为"data/流程 ID/advanced_analysis/output_高级分析任务 ID/monocle2/monocle_final.rds"

相对目录为../data/，绝对目录为/home/mambauser/data/

cds = readRDS("data/AY1732591902625/advanced_analysis/output_54833/monocle2/monocle_final.rds")

#查看meta信息
head(pData(cds))

A data.frame: 6 × 20

	orig.ident	nCount_RNA	nFeature_RNA	mito	resolution.0.2	resolution.0.5	resolution.0.8	resolution.1.1	resolution.1.4	Sample	raw_Sample	barcodes	CellAnnotation	singleR_blood	Neu_anno	singleR_human_all	Size_Factor	Pseudotime	State	barcode
	<chr>	<dbl>	<int>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<chr>	<dbl>	<dbl>	<fct>	<chr>
AAGTTCGTACTGGTTCT	PBMC	500	399	5.6	0	1	0	0	10	PBMC	PBMC	AAGTTCGTACTGGTTCT	T_cells	NK Cell	Unknown	NK_cell	0.4340618	16.367642	6	AAGTTCGTACTGGTTCT
ATCAGGTGCTTTCAGAC	PBMC	526	364	0	0	0	2	2	12	PBMC	PBMC	ATCAGGTGCTTTCAGAC	T_cells	NK Cell	Unknown	T_cells	0.4143317	9.843610	7	ATCAGGTGCTTTCAGAC
GACCCTTTACGTGTTCC	PBMC	528	180	3.2197	0	0	2	2	12	PBMC	PBMC	GACCCTTTACGTGTTCC	T_cells	NK Cell	Unknown	T_cells	0.1979931	13.107587	4	GACCCTTTACGTGTTCC
CTGTATTCGACAGCATT	PBMC	552	426	0.7246	0	0	2	2	12	PBMC	PBMC	CTGTATTCGACAGCATT	T_cells	NK Cell	Unknown	T_cells	0.4693683	11.759849	4	CTGTATTCGACAGCATT
GTGTAACCGATGCGGTC	PBMC	564	429	3.0142	1	4	5	10	11	PBMC	PBMC	GTGTAACCGATGCGGTC	T_cells	T Cell	Unknown	T_cells	0.4721374	0.646801	1	GTGTAACCGATGCGGTC
AGGACTTCGAATCTCTT	PBMC	565	414	3.1858	1	4	5	10	11	PBMC	PBMC	AGGACTTCGAATCTCTT	T_cells	T Cell	Unknown	T_cells	0.4631377	2.981440	1	AGGACTTCGAATCTCTT

#查看细胞的伪时间
plot_cell_trajectory(cds, color_by = "Pseudotime")

output

Warning message:
“\`select_()\` was deprecated in dplyr 0.7.0.
ℹ Please use \`select()\` instead.
ℹ The deprecated feature was likely used in the monocle package.
Please report the issue to the authors.”

#查看细胞注释结果，本数据仅作为示例
plot_cell_trajectory(cds, color_by = "singleR_blood")

读取突变数据

修改下方的突变矩阵文件名

对应该样本的两个突变矩阵（all/alt，可联系客户经理申请释放，仅限于全序列数据）

矩阵文件可以在左上角的 upload 处点击上传

snv_cover_mat = read.delim("PBMC.snp_indel.all_UMI.matrix", header = T, row.names = 1)
snv_cover_mat[1:3, 1:3]

snv_mut_mat = read.delim("PBMC.snp_indel.alt_UMI.matrix", header = T, row.names = 1)
snv_mut_mat[1:3, 1:3]

A data.frame: 3 × 3

	AAGTTCGTACTGGTTCT	CTGCAGGTACGGAGTAG	TAACGACCGACTGCGCA
	<int>	<int>	<int>
SDF4:chr1-1223263:T>G	0	0	0
SLC35E2B:chr1-1668373:C>T	0	0	0
CDK11A:chr1-1709071:C>T	0	0	0

A data.frame: 3 × 3

	AAGTTCGTACTGGTTCT	CTGCAGGTACGGAGTAG	TAACGACCGACTGCGCA
	<int>	<int>	<int>
SDF4:chr1-1223263:T>G	0	0	0
SLC35E2B:chr1-1668373:C>T	0	0	0
CDK11A:chr1-1709071:C>T	0	0	0

可以从突变矩阵行名来选择突变

all_mut = rownames(snv_mut_mat)
all_mut[1:3]

1. 'SDF4:chr1-1223263:T>G'
2. 'SLC35E2B:chr1-1668373:C>T'
3. 'CDK11A:chr1-1709071:C>T'

指定感兴趣的突变

例如 "NPM1:chr5-171405326:G>A", 将该突变的发生数作为新的标签添加到 cds 中

interest_mut = "NPM1:chr5-171405326:G>A"
mut_num = as.numeric(snv_mut_mat[interest_mut,])
names(mut_num) = colnames(snv_mut_mat)
mut_num = mut_num[colnames(cds)]

pData(cds)$mut_num = mut_num

在曲线上展示，mut_num 表示发生该突变的 UMI 数

plot_cell_trajectory(cds, color_by = "mut_num") +
scale_color_gradient(low = "grey", high = "red")

保存图片

p=plot_cell_trajectory(cds, color_by = "mut_num") +
scale_color_gradient(low = "grey", high = "red")
ggsave(p, file = "NPM1_G_to_A.png", width = 6, height = 6) #保存在当前目录下，可在左侧下载

END