克隆排序和进化可视化R包：ClonEvol
导语
GUIDE ╲
ClonEvol是一个用于克隆排序和克隆进化可视化的R包。
背景介绍
进化树在生物学中，用来表示物种之间的进化关系。生物分类学家和进化论者根据各类生物间的亲缘关系的远近，把各类生物安置在有分枝的树状的图表上，简明地表示生物的进化历程和亲缘关系。在进化树上每个叶子结点代表一个物种，如果每一条边都被赋予一个适当的权值，那么两个叶子结点之间的最短距离就可以表示相应的两个物种之间的差异程度。
在肿瘤细胞里染色体存在着克隆演化(clone evolution)，形成现异质性(heterogeneity)， 演变为多克隆性，其中占主导数目的克隆构成肿瘤干系(stemline)，系肿瘤细胞的染色体数目称为众数，除干系外还有一些非主导的克隆，称为旁系。
R包简介
R包ClonEvol利用其他方法预先聚类的变异来推断和可视化克隆进化树。它还可以可视化由其他方法识别的树。它输入的数据是其他工具识别出的杂合变异的聚类，从而推断一致性的克隆进化树，并估计个体样本克隆中的癌细胞比例（也称为克隆频率）。ClonEvol可以处理测序和数据分析中可能会歪曲个体变异的细胞流行率评估的统计不确定性和错误。
细胞流行率（cellular prevalence）可以从样本和肿瘤两个角度来定义：
①样本的细胞流行率=包含事件的整个样本中细胞的比例
②肿瘤细胞流行率=所有肿瘤细胞中包含该事件的细胞比例
ClonEvol使用bootstrap重采样方法来估计克隆的癌细胞部分(CCF)，通过下面的公式(sum rule)，给定变异及clusters的CCF。
根据这个公式可以得到克隆Y(给定其直接子克隆Xi)的CCF抽样分布的bootstrap估计，可用于估计：①克隆Y的CCF的置信区间，②克隆Y的CCF为负值(或非负值)的概率。
克隆Y的负CCF评估表明克隆进化模型中所有克隆{Xi}直接来自克隆Y是一个违反了总和规则的顺序。一个负CCF评估也可能来自统计上的不确定性和错误（上面提到的），在ClonEvol输入的细胞流行率和变异clusters中存在。因此，在ClonEvol中，只有当克隆的CCF的置信区间严重向负值偏移时，才认为违反了和规则，意味着负CCF的概率较高(或CCF非负的概率较低)。
ClonEvol评估所有潜在树中所有亲本无性系的克隆排序。如果一个树的任何克隆都没有违反上面描述的总和规则，那么它会被报告输出。ClonEvol可以产生多种可视化效果，包括：
①Bell绘图来呈现随时间推移的克隆动态（基于Fishplot建立）
②使用细胞球来表示样本的克隆亚群
③对以节点为基础和分枝为基础的树进行注释，以表示样本间的克隆关系和种子模式
使用方法
01
从Github安装
install.packages('devtools')
library(devtools)
install_github('hdng/clonevol')
install.packages('gridBase')
install.packages('gridExtra')
install.packages('ggplot2')
install.packages('igraph')
install.packages('packcircles')
install_github('hdng/trees')02
克隆进化分析流程
使用测序数据进行克隆进化推断的典型分析工作流包括以下步骤:
①准备一套全面和可靠的体细胞突变数据
②基于细胞流行率对突变进行聚类
③评估聚类结果
④推断无性系进化树(如克隆排序)
⑤可视化数据和克隆进化树
⑥解释结果
ClonEvol为步骤③-⑥提供了工具。应该在运行之前使用其他工具(简要描述如下)完成步骤①-②。
下面对以上步骤进行详细介绍
（1）Step 1: 为克隆进化推断准备variants
测序的深度、样本的数量和质量以及体细胞变异的数量和质量可能对由此产生的克隆进化模型产生深远的影响。在理想情况下，对于工作最好需要：
①大规模的样本量
②大规模的变异量（外显子组测序是可以的，但是全基因组测序提供了更好的passenger体细胞突变覆盖率）
③突变时间点
④多重区域样品（由于瘤内异质性）
⑤深度测序，如有针对性的验证
但是并不意味着在没有这样一个理想的数据集的情况下不能研究癌症的克隆进化。但是，在高度异质性的患者/肿瘤中，你的数据可能产生低估真实的模型。
（2）Step 2: 变异聚类
基于样本中细胞流行率的变异聚类是一个关键步骤。变异聚类的目的是识别克隆。由于肿瘤的异质性，不同克隆的细胞流行率在样本之间可能存在不同的频率（如样本A有90%的克隆X和10%的克隆Y，而样本B有50%的克隆X和50%的克隆）。因此，一个克隆的克隆marker变异将在不同样本中显示相似的细胞流行率，并可能聚集在一起。变异聚类工具的最终产品是聚类，每个聚类由一个不同克隆的克隆marker变异组成。
在聚类算法中使用的变异细胞流行率通常由变异等位基因频率（VAF）来衡量，由携带变异基因的读数与位点总读数的比率来计算。聚类算法工作的假设是，VAF提供了很好的变异细胞分数评估，即携带变异的细胞比例。在二倍体杂合变异（拷贝数中性）的情况下，一个变异的癌细胞分数（CCF）可以计算为其VAF的两倍。然而，拷贝数变异在癌症中是很常见的，使它们的VAF进一步偏离，如果不加以纠正，会导致他们的CCF与实际的细胞流行率不一致。二倍体杂合子和拷贝变异都可以用于聚类算法。如果只使用二倍体杂合子变异，可以使用sciClone算法进行聚类，在ClonEvol中也可以使用VAF。如果使用拷贝变异variants，则应该使用拷贝数识别工具（如Pyclone）来进行聚类。在后一种情况下，ClonEvol中更倾向于使用聚类工具提供的拷贝数校正细胞流行率估计，特别是在增扩增和缺失事件之间存在偏差时。
（3）Step 3: 评估变异聚类结果
由于每个类代表一个克隆，一个类缺失或被错误推断可能会阻碍我们成功地构建进化模型。因此，获得良好的聚类结果是非常重要的。ClonEvol提供了跨多个样本的变异聚类的方便可视化，以帮助评估聚类结果，特别是在没有推断的树的情况下。
假设我们已经有一个聚类结果，包括聚类的识别和个体变异的细胞流行率估计。它可以存储为 tabular text文件，并使用read.table读取到数据框中。
下面介绍的例子使用了ClonEvol的复发性髓细胞白血病（AML1）的数据(Ding et al. 2012)。包含两个样本（一个原发样本和一个复发样本）的全基因组测序和深度靶向验证。
①加载ClonEvol和AML1聚类数据
library(clonevol)
data(aml1)
#数据aml1包含以下部分
x <- aml1$variants
#fix(x)
②准备数据聚类
其他工具输出的聚类结果可能不利于ClonEvol分析和可视化。ClonEvol需要输入数据框，该数据框至少包含一个聚类列和一个或多个变异细胞流行率列，每个列对应于一个样本。聚类应该用从1开始的连续整数命名。为了更好地显示，细胞流行率列的名称应该简短。
# 尽量缩短vaf列名
vaf.col.names <- grep('.vaf', colnames(x), value=T)
sample.names <- gsub('.vaf', '', vaf.col.names)
x[, sample.names] <- x[, vaf.col.names]
vaf.col.names <- sample.names
# 准备样本分组
sample.groups <- c('P', 'R')
names(sample.groups) <- vaf.col.names
# 设置要在不同的plot中显示的聚类的顺序
x <- x[order(x$cluster),]③选择克隆的颜色
ClonEvol有内置的颜色，用来区分大约20个不同的克隆。用户也可以指定自己的颜色。为了设置将在整个可视化过程中使用的聚类/克隆的颜色，创建一个颜色矢量，如下所示。在这种情况下，选择了与 Ding et al (2012)原始图形相匹配的颜色。
clone.colors <- c('#999793', '#8d4891', '#f8e356', '#fe9536', '#d7352e')
#如果自己设置颜色，只需将其设置为NULL
clone.colors <- NULL④可视化变异类
可以绘制跨聚类和样本的变异的细胞流行率(CCF或V AF)，使用jitter、box和violin plots来对聚类类进行的密切调查。它可以同时可视化许多样本和clusters
pdf('C:/Users/DELL/Desktop/box.pdf', width = 4, height = 4, useDingbats = FALSE, title='')
pp <- plot.variant.clusters(x,
cluster.col.name = 'cluster',
show.cluster.size = FALSE,
cluster.size.text.color = 'blue',
vaf.col.names = vaf.col.names,
vaf.limits = 70,
sample.title.size = 20,
violin = FALSE,
box = FALSE,
jitter = TRUE,
jitter.shape = 1,
jitter.color = clone.colors,
jitter.size = 3,
jitter.alpha = 1,
jitter.center.method = 'median',
jitter.center.size = 1,
jitter.center.color = 'darkgray',
jitter.center.display.value = 'none',
highlight = 'is.driver',
highlight.shape = 21,
highlight.color = 'blue',
highlight.fill.color = 'green',
highlight.note.col.name = 'gene',
highlight.note.size = 2,
order.by.total.vaf = FALSE)
dev.off()
在评估clusters时，寻找潜在的离群值clusters（如具有少量变异的clusters）、潜在的合并clusters(如具有在多个样本中变异的VAF从零延伸到非零值，以进一步分裂成多个clusters)和噪声clusters（如在样本之间显示非常相似和低的V AF，说明错误的变异调用）。随着更深层次的测序，这些条件可以更宽松，因为聚类可以更准确。对于AML1的情况，存在几个小的clusters，但是超深测序使它们具有可解释性。
⑤跨样本绘制成对的VAFs或CCFs图
如果需要对配对样本检查变异clusters，下面的语句可以对VAF或CCF的成对绘图实现。
plot.pairwise(x, col.names = vaf.col.names,
out.prefix = 'variants.pairwise.plot',
colors = clone.colors)
⑥绘制样本间clusters的平均值/中位数(cluster流)
pdf('C:/Users/DELL/Desktop/flow.pdf', width=5, height=5, useDingbats=FALSE, title='')
plot.cluster.flow(x, vaf.col.names = vaf.col.names,
sample.names = c('Primary', 'Relapse'),
colors = clone.colors)
dev.off()
（4）Step 4: ClonEvol进行克隆排序
①推断克隆进化树
这部分是ClonEvol执行克隆排序并构建一致树。在AML1的例子中，使用变异的VAF。如果你的数据包含由聚类工具（如Pyclone）估计的copy-altered变异和拷贝数校正的变异，可以通过infer.clonal.models中的ccf.col.names参数向ClonEvol提供正确的CCF ，或计算等量拷贝数校正VAF作为CCF的一半。
#infer.clonal.models功能提取聚类结果，并评价各克隆序次以重构克隆进化树，并评估克隆在个体样本中的CCF。
y = infer.clonal.models(variants = x, #variants,变量聚类数据框
cluster.col.name = 'cluster',
#cluster.col.name,聚类列名
vaf.col.names = vaf.col.names,
#vaf.col.names,VAF列名
sample.groups = sample.groups,
cancer.initiation.model='monoclonal',
subclonal.test = 'bootstrap',
subclonal.test.model = 'non-parametric',
num.boots = 1000,
founding.cluster = 1,
cluster.center = 'mean',
ignore.clusters = NULL,
clone.colors = clone.colors,
min.cluster.vaf = 0.01,
sum.p = 0.05,
#sum.p,样本中一个克隆的CCF估计在可接受的克隆排序中是非负的最小概率
alpha = 0.05
# alpha,alpha水平，或克隆的CCF估计的[1 -(置信水平)]
)②将driver事件映射到树中
如果前面的步骤成功并提供了一棵树或几棵树，接下来我们可以将一些driver事件映射到树中，以确保稍后它们可以被可视化。对于AML1样本，is.driver指示该变量是否为（潜在的）driver事件。将使用gene列中的基因名称来注释树中的变异。
y <- transfer.events.to.consensus.trees(y,
x[x$is.driver,],
cluster.col.name = 'cluster',
event.col.name = 'gene')③将基于节点的树转换为基于分支的树
ClonEvol既可以绘制基于节点的树（每个克隆是一个节点），也可以绘制基于分支的树（每个分支代表一个克隆从其亲本克隆进化而来，每个节点代表克隆established/founded的一个点）。在绘制基于分支的树之前，需要准备它。
y <- convert.consensus.tree.clone.to.branch(y, branch.scale = 'sqrt')④树储存在什么地方?
一致树存储在ymatchedmerge .trees中。这是一个数据框列表，每个数据框描述一个带有各种注释的推断树，包括CCF估计和图形化参数。根据不违反sum rule的概率对树进行排序，因此是ymatchedmerge .trees[[1]]是最好的得分树。几个重要的列包括:
• lab: 克隆的标签，与cluster标签匹配。
• parent: 树中克隆的亲本
• sample.with.nonzero.cell.frac.ci: 估计的克隆有阳性CCF的样本
（5）Step 5: 可视化结果
①将多个图形和树一起展示
用以下命令用来绘制变异clusters、bell plots和克隆进化树。
plot.clonal.models(y,
#箱线图参数
box.plot = TRUE,
fancy.boxplot = TRUE,
fancy.variant.boxplot.highlight = 'is.driver',
fancy.variant.boxplot.highlight.shape = 21,
fancy.variant.boxplot.highlight.fill.color = 'red',
fancy.variant.boxplot.highlight.color = 'black',
fancy.variant.boxplot.highlight.note.col.name = 'gene',
fancy.variant.boxplot.highlight.note.color = 'blue',
fancy.variant.boxplot.highlight.note.size = 2,
fancy.variant.boxplot.jitter.alpha = 1,
fancy.variant.boxplot.jitter.center.color = 'grey50',
fancy.variant.boxplot.base_size = 12,
fancy.variant.boxplot.plot.margin = 1,
fancy.variant.boxplot.vaf.suffix = '.VAF',
# bell图参数
clone.shape = 'bell',
bell.event = TRUE,
bell.event.label.color = 'blue',
bell.event.label.angle = 60,
clone.time.step.scale = 1,
bell.curve.step = 2,
#基于节点的一致树参数
merged.tree.plot = TRUE,
tree.node.label.split.character = NULL,
tree.node.shape = 'circle',
tree.node.size = 30,
tree.node.text.size = 0.5,
merged.tree.node.size.scale = 1.25,
merged.tree.node.text.size.scale = 2.5,
merged.tree.cell.frac.ci = FALSE,
#基于分支的一致树参数
merged.tree.clone.as.branch = TRUE,
mtcab.event.sep.char = ',',
mtcab.branch.text.size = 1,
mtcab.branch.width = 0.75,
mtcab.node.size = 3,
mtcab.node.label.size = 1,
mtcab.node.text.size = 1.5,
#细胞总体参数
cell.plot = TRUE,
num.cells = 100,
cell.border.size = 0.25,
cell.border.color = 'black',
clone.grouping = 'horizontal',
#meta参数
scale.monoclonal.cell.frac = TRUE,
show.score = FALSE,
cell.frac.ci = TRUE,
disable.cell.frac = FALSE,
#输出图参数
out.dir = 'C:/Users/DELL/Desktop/output',
out.format = 'pdf',
#out.prefix = 'C:/Users/DELL/Desktop/output.plot',
overwrite.output = TRUE,
width = 8,
height = 4,
#每个面板从左到右缩放的宽度向量
panel.widths = c(3,4,2,4,2)
)
②绘制树
如果你只想绘制树：
pdf('C:/Users/DELL/Desktop/trees.pdf', width = 3, height = 5, useDingbats = FALSE)
plot.all.trees.clone.as.branch(y, branch.width = 0.5,
node.size = 1, node.label.size = 0.5)
dev.off()
③可视化其他工具预测的树
为了将其他工具预测的树可视化，需要准备两个文件:
A. variants.tsv: 类似于ClonEvol输入，使用其他工具分配的clusters
B. tree.tsv: 预测树，包括:
a.至少3列：克隆，亲本，sample.with.nonzero.cell.frac.ci
b.其他列：颜色，事件
#首先从文件读取树和变异list
y = import.tree('tree.tsv', 'variants.tsv')
#然后，可以根据扩展到克隆marker variants的分支长度，准备一个注释的基于分枝的树
y = convert.consensus.tree.clone.to.branch(y, branch.scale = 'sqrt')
#然后，还可以将driver events映射到树中(如果准备的变异文件有“cluster”、“is.driver”和“gene”列):
y <- transfer.events.to.consensus.trees(y,
y$variants[y$variants$is.driver,],
cluster.col.name = 'cluster',
event.col.name = 'gene')
#绘制树
pdf('C:/Users/DELL/Desktop/imported-tree.pdf', width=3, height=5, useDingbats=F)
plot.all.trees.clone.as.branch(y, branch.width = 0.5,
node.size = 1, node.label.size = 0.5)
dev.off()小编总结
R包ClonEvol利用其他方法预先聚类的变异来推断和可视化克隆进化树。它还可以可视化由其他方法识别的树。它输入的数据是其他工具识别出的杂合变异的聚类，从而推断一致性的克隆进化树，并估计个体样本克隆中的癌细胞比例（也称为克隆频率）。ClonEvol可以处理测序和数据分析中可能会歪曲个体变异的细胞流行率评估的统计不确定性和错误。
腾讯云开发者
