第一步：将你的染色体位置坐标文件整理成bed格式。bed格式文件至少包括前3列，分别是：染色体的名字、染色体上的起始位置、染色体上的终止位置。这一步无论用写字板、excel、R等进行处理都可以，文件的后缀名也不重要，因为强行将文件后缀改为bed时，在后面的Linux系统中进行bedtools处理时也会报错。所需的bed格式文件参见下图。 第二步：获得人类基因组的注释文件。可从gencode中根据自己的需要下载hg38或者hg19版本的人类基因组注释文件(文章中以hg38为例)。这一步可以进gencode官网(https://www./human/)进行本地下载，然后用filezilla等文件传输工具将下载的本地文件传输到服务器。也可以直接在服务器的Linux系统中进行ftp下载。 本地下载： ftp下载： 获得下载链接后，在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如果在第一步的时候已将待处理的bed格式文件保存为了Tab键分隔格式，但是在后面的处理中仍然报错，不妨再进行一次Tab键分隔处理。 第五步：在Linux系统中利用bedtools得到包含染色体位置坐标的蛋白编码基因。首先需要启动自己安装了bedtools软件的conda小环境，然后输入下面的代码： bedtools intersect -a motif1.bed -b ~/dna/exercise/protein_coding.hg38.position -wa -wb也可以对结果进行汇总，将位于相同染色体坐标的基因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 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也可以得到编码lncRNA的基因等，大家可以借鉴曾老师的教程进行探索。 如果你看这些软件用法以及Linux命令比较困难，而且时间很宝贵，建议参加生信技能树的学习班：
|
|
|
