如何使用bedtools根据染色体上的起止位置拿到基因symbol

如何使用bedtools根据染色体上的起止位置拿到基因symbol,相信很多没有经验的人对此束手无策,为此本文总结了问题出现的原因和解决方法,通过这篇文章希望你能解决这个问题。

成都创新互联服务项目包括灵山网站建设、灵山网站制作、灵山网页制作以及灵山网络营销策划等。多年来,我们专注于互联网行业,利用自身积累的技术优势、行业经验、深度合作伙伴关系等,向广大中小型企业、政府机构等提供互联网行业的解决方案,灵山网站推广取得了明显的社会效益与经济效益。目前,我们服务的客户以成都为中心已经辐射到灵山省份的部分城市,未来相信会继续扩大服务区域并继续获得客户的支持与信任!

第一步:将你的染色体位置坐标文件整理成bed格式。

bed格式文件至少包括前3列,分别是:染色体的名字、染色体上的起始位置、染色体上的终止位置。这一步无论用写字板、excel、R等进行处理都可以,文件的后缀名也不重要,因为强行将文件后缀改为bed时,在后面的Linux系统中进行bedtools处理时也会报错。所需的bed格式文件参见下图。

如何使用bedtools根据染色体上的起止位置拿到基因symbol


 
 

第二步:获得人类基因组的注释文件。

可从gencode中根据自己的需要下载hg38或者hg19版本的人类基因组注释文件(文章中以hg38为例)。这一步可以进gencode官网(https://www.gencodegenes.org/human/)进行本地下载,然后用filezilla等文件传输工具将下载的本地文件传输到服务器。也可以直接在服务器的Linux系统中进行ftp下载。

本地下载:

如何使用bedtools根据染色体上的起止位置拿到基因symbol

ftp下载:

如何使用bedtools根据染色体上的起止位置拿到基因symbol

获得下载链接后,在Linux系统中输入下面的代码进行ftp下载:

wget ftp://ftp.ebi.ac.uk/pub/databases/gencode/Gencode_human/release_34/gencode.v34.annotation.gtf.gz
   

第三步:在Linux系统中处理下载的基因组注释文件,得到人类的蛋白编码基因的位置坐标。

在Linux系统中输入下面的代码,得到hg38版本的人类蛋白编码基因的位置坐标:

zcat  gencode.v34.annotation.gtf.gz | grep   protein_coding |perl -alne '{next unless $F[2] eq "gene" ;/gene_name \"(.*?)\";/; print "$F[0]\t$F[3]\t$F[4]\t$1" }' >protein_coding.hg38.position
   

第四步:在Linux系统中将自己待处理的bed格式文件转换为Tab键分隔的文件。

先将待处理的坐标bed格式文件链接或复制到第三步得到的结果文件所在的目录下,然后修改这一文件的后缀名为bed,再将这一文件转化为Tab键分隔的后缀名为bed的文件,需输入下面的代码(motif1.bed是自己命名的待处理坐标文件):

mv motif1.tsv motif1.bed
perl -p -i -e 's/ /\t/g' motif1.bed
 

如果在第一步的时候已将待处理的bed格式文件保存为了Tab键分隔格式,但是在后面的处理中仍然报错,不妨再进行一次Tab键分隔处理。

 

第五步:在Linux系统中利用bedtools得到包含染色体位置坐标的蛋白编码基因。

首先需要启动自己安装了bedtools软件的conda小环境,然后输入下面的代码:

bedtools intersect -a motif1.bed  -b ~/dna/exercise/protein_coding.hg38.position  -wa -wb
 

如何使用bedtools根据染色体上的起止位置拿到基因symbol


 

也可以对结果进行汇总,将位于相同染色体坐标的基因symbol写在一块,此时只需要加上|后面的代码即可。| 之前的文件得到的结果有几列,-c后面的数字就写几。如我得到的有7列,-c后面就写7。

bedtools intersect -a motif1.bed  -b ~/dna/exercise/protein_coding.hg38.position  -wa -wb | bedtools groupby -i - -g 1-4 -c 7 -o collapse
 

如何使用bedtools根据染色体上的起止位置拿到基因symbol


 

也可以另存结果:

bedtools intersect -a motif1.bed  -b ~/dna/exercise/protein_coding.hg38.position  -wa -wb | bedtools groupby -i - -g 1-4 -c 7 -o collapse >gene.tsv
 

新保存的gene.tsv文件就是结果文件了,然后可以拿着结果进行后续处理啦~。

看完上述内容,你们掌握如何使用bedtools根据染色体上的起止位置拿到基因symbol的方法了吗?如果还想学到更多技能或想了解更多相关内容,欢迎关注创新互联行业资讯频道,感谢各位的阅读!


当前题目:如何使用bedtools根据染色体上的起止位置拿到基因symbol
当前网址:http://ybzwz.com/article/jgshhs.html