肿瘤拷贝数变异的差异分析应该以这个肺癌脑转移为标准
粉丝的问题很朴素，就是想把TCGA数据库里面的非小细胞肺癌里面的肺鳞癌区分成为是否有TP53这个基因的somatic突变的两个分组，然后去比较这两个组别里面的病人的肿瘤拷贝数变异，做一个差异分析。
关于这样的分析，我们推荐一个最近发表在《NATURE GENETICS》杂志，2020的文章：《Genomic characterization of human brain metastases identifies drivers of metastatic lung adenocarcinoma》，比较了BM-LUAD 和TCGA-LUAD 两个队列的拷贝数差异情况。
首先展现两个队列拷贝数变异全景图，而且是扩增情况和缺失情况分开展现 ：
拷贝数变异全景图
图例是：a, GISTIC amplification (top) and deletion (bottom) plots of BM-LUAD (n = 73) and covariate-matched samples in the TCGA-LUAD (n = 464) cohort
可以看到，这样的结果其实肉眼也能看到一些两个图的差异，但是毕竟没有统计学检验，这个时候两个拷贝数变异结果都可以导入GISTIC2软件，然后统计学检验后绘制散点图：
统计学检验后绘制散点图
图例是：b, Differentially amplified or deleted regions in BM-LUAD compared to TCGA-LUAD. Significant differential regions are labeled (FDR < 0.01; G-score difference > 0.5).
然后可以再具体看有差异的地方（染色体片段），整理展现4个基因组区域：
具体看有差异的地方
图例是：GISTIC plots of control region (NKX2-1) and candidate metastatic driver regions
可以看到 这个区域 region (NKX2-1) 在图上面是没有拷贝数差异的，但是另外的3个区域很明显在两个队列有肉眼可见的差异。有意思的是前面的差异分析散点图里面密码是4个区域，这个时候 candidate metastatic driver regions仅仅是展现了3个区域，缺了那个4q31区域，我猜想可能是这个4q31区域上面并没有很出名的基因，所以被研究者首先排除了。
最后箱线图去验证拷贝数差异区域上面的热点基因的拷贝数情况：
热点基因的拷贝数情况
图例是：d, Frequencies of amplifications or deletions of candidate metastatic drivers, adjusted by matching weights to control for confounding.
为了让在数据分析更可信，研究者们还加入了实验验证：
实验验证数据分析的拷贝数变异
图例是：e, Frequencies of amplifications of MYC and YAP1 in the BM-LUAD validation cohort (n = 105) as determined by FISH. TCGA-LUAD was reused as the control cohort.
最后作者的结论是：Amplified candidate drivers (MYC, MMP13, YAP1) tended to occur after the divergence of the metastatic and primary tumor lineages, in agreement with positive selection of these amplifica- tions contributing to a pro-metastatic phenotype
那么问题来了，这样的数据分析如何做呢？
这个还是有点难度的， 首先你需要拿到两个队列各自的cnv信息，如果是肿瘤外显子数据，首先需要对fastq数据进行质量控制，比对到参考基因组拿到bam文件，接着使用cnvkit处理这些bam文件进行cnv分析，CNVkit的标准流程，代码如下：
# From baits and tumor/normal BAMs
ref=/home/biosoft/GATK/resources/bundle/hg38/Homo_sapiens_assembly38.fasta
bed=/home/annotation/CCDS/human/hg38.baits.bed
source activate wes
## set --method amplicon if the data from PCR
cnvkit.py batch tumor_bam/*.bam \
--normal normal_bam/*.bam \
--targets $bed \
--fasta $ref \
--output-reference my_reference.cnn --output-dir T_results
运行起来还是非常快的，这个cnv分析对每个肿瘤外显子队列都是需要测肿瘤样品也同步测序其正常对照哦，这样分析才能是可靠的。
多样本的cns文件可以合并到IGV查看，很简单的一个awk脚本即可，如下：
awk '{print FILENAME"\t"$0}' *_recal.cns |grep -v chromosome |sed 's/_recal.cns//g' |awk '{print $1"\t"$2"\t"$3"\t"$4"\t"$8"\t"$6}' >final.seg
不同队列的 final.seg 文件，都需要去走前面提到的 GISTIC流程 ，一般来说每个 cnv文件（final.seg ）走了GISTIC流程 ，都会得到如下所示：
gistic_results
但是两个队列两个这样的图片，很难去肉眼比较，所以需要参考这个文章：《Genomic characterization of human brain metastases identifies drivers of metastatic lung adenocarcinoma》
如果你也恰好有两个队列的肿瘤外显子数据，或者有想挖掘的公共数据库里面的队列可以区分成为两个分组，我们提供这样的数据分析服务，绘制同样的4个图表：
- 肿瘤拷贝数扩增和缺失全景图
- 两个队列的cnv信息的的 GISTIC流程结果，并且差异分析散点图
- 局部差异区域可视化
- 局部差异区域基因可视化