2021-02-24--TCGA 拷贝数变异数据分析
GISTIC2 分析 TCGA拷贝数变异数据及可视化 全过程
一、数据下载
1、Segment_data数据下载
TCGA官网下载
下载地址
它有如下CNV类型，而在TCGA数据库里面我们通常关心的是somatic CNV，也就是那些肿瘤病人的拷贝数变异然后要剔除正常对照里面的CNV多态性信息，只有这些somatic CNV才更可能是跟肿瘤相关的。
Masked Copy Number Segment, 此表是在上面数据上过滤掉了一些与生殖和性染色体相关的数据。
image.png
①Copy Number Segment：
从筛选界面就能看出来这些信息：
Experimental Strategy：Genotyping Array；
Workflow Type：DNAcopy
Platform：affymetrix snp 6.0
下载一个文件查看一下发现它又一下六列信息；病历号，染色体位置，Num_probes，Segment_Mean；
Num_probes应该是指的探针号，暂时来说用不到这个信息。
Segment_Mean值--就是 log2(copy_number/ 2), 正常来说人是二倍体生物则此value值为0，如果拷贝数小于2（删除）则小于0，拷贝数大于2（扩增）则大于0.
image.png
用R语言下载和整理
官网下下来的文件整理比较麻烦，这里直接用R语言整理比较方便
BiocManager::install('TCGAbiolinks')
library(TCGAbiolinks)
query <- GDCquery(project = "TCGA-LUAD",
data.category = "Copy Number Variation",
data.type = "Masked Copy Number Segment")
GDCdownload(query, method = "api", files.per.chunk = 100)
segment_dat <- GDCprepare(query = query)
segment_dat$Sample <- substring(segment_dat$Sample,1,16)
segment_dat <- grep("01A$",segment_dat$Sample) %>%
segment_dat[.,]
segment_dat[,1] <- segment_dat$Sample
segment_dat <- segment_dat[,-7]
write.table(segment_dat,"MaskedCopyNumberSegment.txt",sep="\t",
quote = F,col.names = F,row.names = F)
segment_data长这样
image.png
2、Markers File
下载地址 ，选择文件 SNP6 GRCh38 Remapped Probeset File for Copy Number Variation Analysis，并注意提示 If you are using Masked Copy Number Segment for GISTIC analysis, please only keep probesets with freqcnv = FALSE ，所以只保留 freqcnv = FALSE 的行和前三列。
R语言整理Markers file 代码如下
snp<-read.table('snp6.na35.remap.hg38.subset.txt',header=T)
snp<-snp[snp$freqcnv=='FALSE',]
snp<-snp[,1:3]
colnames(snp)<-c("Marker name","chromosome","Marker position")
write.table(snp,'Marker.txt',sep="\t",
quote = F,col.names = F,row.names = F)
文件长这样
image.png
好了，我们所需要下载的文件就已经准备好了
二、GISTIC2的使用
1、在线分析
Genepattern 云分析平台中有GISTIC2模块，可以直接注册使用，缺点是网不好。
2、安装GISTIC2软件分析
在linux服务器上安装GISTIC2软件
wget -c ftp://ftp.broadinstitute.org/pub/GISTIC2.0/GISTIC_2_0_23.tar.gz
mkdir GISTIC2 && mv GISTIC_2_0_23.tar.gz GISTIC2/ && cd GISTIC2/
tar -xzvf GISTIC_2_0_23.tar.gz
cd MCR_Installer/ && unzip MCRInstaller.zip
./install -mode silent -agreeToLicense yes -destinationFolder /PATH/GISTIC2/MATLAB_Compiler_Runtime/
export XAPPLRESDIR=/PATH/GISTIC2/MATLAB_Compiler_Runtime/v83/X11/app-defaults:$XAPPLRESDIR
参数
软件的参数有很多，详细可以见软件包下的document，这里介绍几个常用的：
-b 输出目录，后面接一个输出文件前缀prefix
-seg 上面讲到的Segmentation文件
-refgene 上面讲到的参考基因组文件
-conf 置信水平，默认是0.75
TCGA官方参数
gistic2
-b <base_directory>
-seg <segmentation_file>
-mk <marker_file>
-refgene <reference_gene_file>
-conf 0.99
-td 0.1 -genegistic 1 -gcm extreme -js 4 -maxseg 2000 -ta 0.1 -armpeel 1 -brlen 0.7 -cap 1.5 -qvt 0.25 -rx 0 -savegene 1 (-broad 1)
输出文件
image.png
三、使用maftools读取GISTIC输出的CNV数据并统计
需要如下这四个文件输入到R包--maftools进行可视化分析
"all_lesions.conf_99.txt",
"amp_genes.conf_99.txt"
"del_genes.conf_99.txt"
"scores.gistic"
library(maftools)
luad.gistic <- readGistic(gisticAllLesionsFile="all_lesions.conf_99.txt", gisticAmpGenesFile="amp_genes.conf_99.txt", gisticDelGenesFile="del_genes.conf_99.txt", gisticScoresFile="scores.gistic", isTCGA=TRUE)
## -Processing Gistic files..
## --Processing amp_genes.conf_99.txt
## --Processing del_genes.conf_99.txt
## --Processing scores.gistic
## --Summarizing by samples
染色体突变图
gisticChromPlot(gistic=luad.gistic, markBands="all")
image.png
气泡图
gisticBubblePlot(gistic=luad.gistic)
image.png
也可以用oncoplot同时展示CNV和突变
luad.plus.gistic <- read.maf(maf="TCGA.LUAD.maf", gisticAllLesionsFile="all_lesions.conf_99.tx
t", gisticAmpGenesFile="amp_genes.conf_99.txt", gisticDelGenesFile="del_genes.conf_99.txt", gist
icScoresFile="scores.gistic")
# 因为样本数太多，使用borderCol=NULL不显示样本边界。可惜的是目前gisticOncoPlot还没有加上这个参数
oncoplot(maf=luad.plus.gistic, borderCol=NULL, top=10)
image.png
