在处理前需要准备两个文件:
1.基因坐标的信息
需要下载gencode数据库或者其他数据库的人类基因组注释文件
下载 gencode.v25.annotation.gtf.gz的源代码
wget ftp://ftp.sanger.ac.uk/pub/gencode/Gencode_human/release_25/gencode.v25.annotation.gtf.gz
之后再将文件解压,再将其转换为基因组坐标文件
cat gencode.v25.annotation.gtf|perl -alne '{next unless $F[2] eq "gene";print}'|grep -w HAVANA |cut -f 1,4,5,9| cut -d";" -f 1,2,4|sed 's/gene_id//g'|sed 's/gene_type//g'|sed 's/gene_name//g'|sed 's/;//g'| sed 's/\"//g'|perl -alne '{/(ENSG\d+)/;print "$1\t$_"}' >human.gene.positions
得到文件如下:
ENSG00000223972 chr1 11869 14409 ENSG00000223972.5 transcribed_unprocessed_pseudogene DDX11L1
ENSG00000227232 chr1 14404 29570 ENSG00000227232.5 unprocessed_pseudogene WASH7P
ENSG00000243485 chr1 29554 31109 ENSG00000243485.4 lincRNA MIR1302-2
ENSG00000237613 chr1 34554 36081 ENSG00000237613.2 lincRNA FAM138A
ENSG00000268020 chr1 52473 53312 ENSG00000268020.3 unprocessed_pseudogene OR4G4P
ENSG00000240361 chr1 62948 63887 ENSG00000240361.1 unprocessed_pseudogene OR4G11P
ENSG00000186092 chr1 69091 70008 ENSG00000186092.4 protein_coding OR4F5
ENSG00000238009 chr1 89295 133723 ENSG00000238009.6 lincRNA RP11-34P13.7
ENSG00000239945 chr1 89551 91105 ENSG00000239945.1 lincRNA RP11-34P13.8
ENSG00000233750 chr1 131025 134836 ENSG00000233750.3 processed_pseudogene CICP27
2.准备单细胞测序转录组的表达矩阵,这个就不赘述了.
3.分析的方法和软件。
1.软件项目地址: https://github.com/broadinstitute/inferCNV.git
软件说明书:https://github.com/broadinstitute/inferCNV/wiki 是一个R包,安装即可,并且下载源代码。
2.另外除了软件分析,还可以自己写代码分析:方法是分染色体用100个步长移动计算CNV。所以每个染色体上面头尾50个基因是没办法计算CNV的,中间的基因的CNV就是它上下100个基因的表达量的平均值,就是这么简单。
之后就可以用CNV的矩阵归一化生成热图了。
网友评论