使用cnvkit来对大批量wes样本找cnv

楼主  收藏   举报   帖子创建时间:  2018-12-16 00:00 回复:1 关注量:130

cnvkit被设计来处理同一个批次的多个肿瘤配对样本测序情况,首先对所有的normal数据进行bin处理拿到背景值,然后就这个背景值来处理所有的tumor测序数据计算拷贝数变异情况。

该软件使用比较复杂,建议读一读官网教程。所有的命令都被包装到一个python脚本里面,使用该脚本调用一系列字命令,如下:

每个命令都有自己的特殊功能,需要仔细阅读。

流程图:

流程代码如下:

可以看到软件提供的命令基本上的都用到了,coverage—>fix—>segment—>segment

其实上面这么一大串的命令已经被包装成了一句命令,就是:

这个一句话命令与上面的多行代码是等效的,默认的segment算法是 circular binary segmentation algorithm (CBS),也可以用-m切换使用其它算法,比如: faster HaarSeg (haar) or Fused Lasso (flasso)

上面得到的只是segment的结果,还可以call一下:

这个时候需要考虑到已有的vcf变异文件,或者计算好的tumor纯度,或者倍性等等。把segment计算得到的log2 ratio值还原成 0,1,2,3,4这样的拷贝数。

但是,事实上上面的代码一般来说是不能直接使用的,因为我们的测序数据通常是WES数据,需要加上很多参数。

实践运行cnvkit

上面流程很复杂,命令也很多,但是不知道也没关系,用起来其实就一个batch命令即可,当然这个batch命令本身参数也不少,而且被设计用来处理不同的数据情况。

值得注意的就是,如果是全基因组测序数据,用batch --method wgs ,如果是捕获基因组测序,包括全外显子,就用 batch --method amplicon ,然后一定要提供捕获区域的bed文件,一般是外显子加上其侧翼上下游的50bp长度。

人类外显子长度平均是200bp,所以默认的bin是267bp,这样可以把比较长的exon给拆分开来。

还有就是 access 参数需要的文件,

至于最后得到cnv片段该如何注释到对应区域的基因这种小事,就不在本文讨论范围啦。

输入输出文件

其中 coverage 命令会对每一个 normal 样本都计算 *.targetcoverage.cnn and *.antitargetcoverage.cnn files , 说明是: target and antitarget coverage tables (.cnn)

这些文件需要合并起来:

然后再校正区域测序深度及GC含量,之后变成 copy number ratios (.cnr) 文件。

最后把 copy number ratios (.cnr) 文件用segment算法跑一下即可,输出cns后缀文件的 segment信息。

对于normal样本只需要输出cnn即可,合并成Reference.cnn,然后把一个个tumor样品,根据这个Reference.cnn来计算 cnr,进而计算 cns 。

最后对于cns可以进行call找到真正的拷贝数。

可以看到 cns文件内容如下,其中第4列是注释的基因,因为太多,看不清楚,我就没有秀给大家。

上面的segment结果还可以call一下,如果有需要的话。

可视化结果:

很明显可以看到有拷贝数变异的区域了。

  • ***2XP 2018-04-24 02:10
    #1

    谢谢提醒。问题已经修复。