得到一个物种所有基因的TSS(转录起始位点)区域的bed文件。

首先在UCSC的table browser 里面下载下面这个文件:

可以看到我这里选择的mm10的refseq系统的所有基因,共有29037个不同的tss,36872个转录本,只有24540个基因,说明有部分基因有多个tss,这个其实挺麻烦的。

#bin    name    chrom   strand  txStart txEnd   cdsStart    cdsEnd  exonCount   exonStarts  exonEnds    score   name2   cdsStartStat    cdsEndStat  exonFrames0    NM_001282945    chr1    -   134199214   134235457   134202950   134234355   3   134199214,134234014,134235227,  134203590,134234446,134235457,  0   Adora1  cmpl    cmpl    2,0,-1,0    NM_001039510    chr1    -   134199214   134235457   134202950   134234355   3   134199214,134234014,134235227,  134203590,134234412,134235457,  0   Adora1  cmpl    cmpl    2,0,-1,0    NM_001291930    chr1    -   134199214   134235457   134202950   134203505   2   134199214,134235227,    134203590,134235457,    0   Adora1  cmpl    cmpl    0,-1,0    NM_001291928    chr1    -   134199214   134234856   134202950   134234733   2   134199214,134234662,    134203590,134234856,    0   Adora1  cmpl    cmpl    2,0,0    NM_001008533    chr1    -   134199214   134235457   134202950   134234355   2   134199214,134234014,    134203590,134235457,    0   Adora1  cmpl    cmpl    2,0,

其实里面可以设置直接下载所有基因的TSS区域的bed文件,可是我不会设置各种参数,也懒得去摸索,直接对上面的文件我可以写脚本处理得到需要的数据形式。 需要输出的是bed格式文件,如下: chrom / chromStart /chromEnd /name /score /strand 我这里定义的TSS(转录起始位点)区域上下游2.5kb,所以代码如下:

perl -alne '{next if /^#/;if($F[3] eq "+"){$start=$F[4]-2500;$end=$F[4]+2500}else{$start=$F[5]-2500;$end=$F[5]+2500}print join("\t",$F[2],$start,$end,$F[12],0,$F[3])}' ucsc.refseq.txt |sort -u >ucsc.refseq.tss.bed

最后得到的文件如下:

tail ucsc.refseq.tss.bed chrY    816212  821212  Uba1y   0   +chrY    81798997    81803997    Gm20747 0   -chrY    82222714    82227714    Gm20736 0   +chrY    83925411    83930411    Gm20854 0   -chrY    85527019    85532019    Gm20854 0   -chrY    8832669 8837669 Gm20815 0   -chrY    895287  900287  Kdm5d   0   +chrY    90752550    90757550    G530011O06Rik   0   -chrY    90782941    90787941    Erdr1   0   +chrY    90836906    90841906    G530011O06Rik   0   

这里面会有一个问题,对于部分基因在非正常染色体的,会出现如下诡异的结果,建议干脆删除掉。

chr4_GL456216_random    13380   18380   Dhrsx   0   +chr4_GL456350_random    -1369   3631    Ccl21c  0   -chr4_GL456350_random    -1369   3631    Gm10591 0   -chr4_GL456350_random    -1369   3631    Gm13304 0   -

记住,这个时候,部分基因还有多个tss哦,反正取决于你的下游分析流程啦。

原文发布于微信公众号 - 生信技能树(biotrainee)

原文发表时间:2017-07-26

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

发表于

我来说两句

0 条评论
登录 后参与评论

相关文章

来自专栏友弟技术工作室

solidity语言介绍以及开发环境准备

2223
来自专栏玉树芝兰

如何用Sikuli自动录入成绩?

手里明明有一份学生成绩Excel表格,却还得一一手动把它们输入到教务系统?类似这样的简单重复枯燥操作,其实你都可以一键让电脑自动替你完成。

672
来自专栏生信技能树

【直播】我的基因组 35:bam格式转化为bw格式看测序深度分布

我们在之前说到过bam文件还有55G大小,这样的文件,在很多时候都不方便查看和转移。而有些时候,我们只需要我们的测序数据在全基因组的测序深度这一个值,并不需要具...

3898
来自专栏魏艾斯博客www.vpsss.net

代码实现 WordPress 文章中英文数字间自动添加空格

1753
来自专栏写代码的海盗

初生牛犊不怕虎 golang入坑系列

读前必读,下面所有内容都是来自这里。 放到这里的目的,就是为了比对一下,哪里的读者多。。(X := '难' || X :='耐' || X := '好' ) 都...

2704
来自专栏数据分析

使用Visual Studio 2010 一步一步创建Powershell Module 和 Cmdlet

之前写了一个C# 调用PowerShell方法, 那么怎么反过来操作呢,也就是怎么样用C#写一个powershell命令呢? 现在就用C#写一个超级简单的Mod...

3849
来自专栏HaHack

Speed Up the Rendering Process of hexo 3

963
来自专栏数据魔术师

词云图:论一个精致猪猪男孩的数据修养

“词云”就是对网络文本中出现频率较高的“关键词”予以视觉上的突出,形成“关键词云层”或“关键词渲染”,从而过滤掉大量的文本信息,使浏览网页者只要一眼扫过...

703
来自专栏友弟技术工作室

solidity语言介绍以及开发环境准备

Solidity 是一门面向合约的、为实现智能合约而创建的高级编程语言。这门语言受到了 C++,Python 和 Javascript 语言的影响,设计的目的是...

3084
来自专栏coding

听说,撸代码,ide与vim更配哦vim折腾记vim常用命令

在选择编辑器上面,我是一个纠结的人,曾经年少的我执着地追求一款万能的编辑器,可以支持所有编辑语言,灵活可定制,可纯粹用键盘操作。符合这种条件的编辑器,非vim莫...

712

扫码关注云+社区