如何根据class_code筛选转录本?
问题描述
链特异文库鉴定长链非编码RNA(lncRNA)的基本步骤是
- hisat2将原始测序数据比对到参考基因组
- samtools获得排序的bam文件
- stringtie每个样本分别组装得到转录本,获得是一个gtf文件
stringtie -p 12 --rf -o transcripts.gtf sorted.bam
- gffcompare合并所有样本的gtf文件
ls *.gtf > list.merged.txt
~/Biotools/gffcompare/gffcompare -r reference.gtf -i list.merged.txt -o merged
得到一个 merged.combined.gtf这个文件里给每一个转录本分配了一个class_code用来表示转录本相对于参考基因组的位置
以上图片来源于论文 GFF Utilities: GffRead and GffCompare
长链非编码RNA通常是选择class_code为u/x/i,比如论文 Global identification of Arabidopsis lncRNAs revealsthe regulation of MAF4 by a natural antisense RNA 中就提到 only transcripts with TAIR10 annotation [Cufflinks class codes ‘u’ (intergenic transcripts),’x’ (Exonic overlap with reference on the opposite strand),’i’ (transcripts entirely within intron) were retained.
那么问题就来了,如何利用 merged.combined.gtf 这个文件获得 class_code 为 u、x和i的转录本的gtf文件呢
找到了一个办法,python中有一个模块 pyGTF,github链接是https://github.com/chengcz/pyGTF
直接使用pip安装
pip install pyGTF
可以解析gft格式的注释文件
利用这个模块来写一个简单的脚本
import sys
from pyGTF import Transcript
from pyGTF import GTFReader
in_gft = sys.argv[1]
class_code = sys.argv[2]
out_gtf = sys.argv[3]
fw = open(out_gtf,'w')
with GTFReader(in_gtf,flag_stream=True) as fi:
for i in fi:
if i._attri['class_code'] == class_code:
i.to_gtf(fw)
fw.close()
使用方法是
python 01.py in.gtf i out.gtf
####今天学到的另外一个知识点: samtools统计fasta文件序列长度,根据序列名提取序列
参考
https://www.cnblogs.com/xudongliang/p/5200655.html
使用命令
samtools faidx input.fasta
会生成一个input.fasta.fai的文件,文件的内容总共有5列 第一列是序列名,第二列是序列长度,第四列是每行多少个碱基
根据序列名提取序列 这里好像只能提取单条序列
samtools faidx input.fasta TCONS_00000018 > TCONS_00000018.fa
还可以加上指定的位置
samtools faidx input.fasta TCONS_00000018:1-10
>TCONS_00000018:1-10
TGGGCGAACG
- JavaScript 教程
- JavaScript 编辑工具
- JavaScript 与HTML
- JavaScript 与Java
- JavaScript 数据结构
- JavaScript 基本数据类型
- JavaScript 特殊数据类型
- JavaScript 运算符
- JavaScript typeof 运算符
- JavaScript 表达式
- JavaScript 类型转换
- JavaScript 基本语法
- JavaScript 注释
- Javascript 基本处理流程
- Javascript 选择结构
- Javascript if 语句
- Javascript if 语句的嵌套
- Javascript switch 语句
- Javascript 循环结构
- Javascript 循环结构实例
- Javascript 跳转语句
- Javascript 控制语句总结
- Javascript 函数介绍
- Javascript 函数的定义
- Javascript 函数调用
- Javascript 几种特殊的函数
- JavaScript 内置函数简介
- Javascript eval() 函数
- Javascript isFinite() 函数
- Javascript isNaN() 函数
- parseInt() 与 parseFloat()
- escape() 与 unescape()
- Javascript 字符串介绍
- Javascript length属性
- javascript 字符串函数
- Javascript 日期对象简介
- Javascript 日期对象用途
- Date 对象属性和方法
- Javascript 数组是什么
- Javascript 创建数组
- Javascript 数组赋值与取值
- Javascript 数组属性和方法