根据gff文件判断一段序列是否位于其内

最新推荐文章于 2020-12-30 08:11:17 发布
最新推荐文章于 2020-12-30 08:11:17 发布
阅读量826 收藏
点赞数
分类专栏： 生物信息 python 文章标签： 基因组-生物信息学 python gff3
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/msw521sg/article/details/52562229
版权

1.txt

UN227692 chr6B 558820383 558820604 . 100.0 99.1 216 2 0 4
UN113387 chr7A 683635472 683635624 . 100.0 100.0 153 0 0 0
UN128584 chr7D 27592786 27593326 . 100.0 100.0 541 0 0 0
UN170802 chr4B 505369881 505370146 . 98.9 99.6 265 1 0 0
UN181903 chr3A 106517703 106518022 . 100.0 100.0 320 0 0 0
UN076932 chr2B 452598011 452598795 . 99.1 99.2 781 4 2 0
UN067930 chr3D 23548729 23549162 . 100.0 100.0 434 0 0 0
UN066624 chr2B 514684148 514685994 + 100.0 100.0 697 0 0 0
UN025033 chr5A 66831930 66832944 . 100.0 99.9 1014 0 1 0
UN011472 chr2D 600422765 600422990 . 100.0 99.6 225 1 0 0
UN061273 chr6D 445614030 445614663 + 98.4 99.6 538 1 1 0
UN156200 chr2B 32058281 32059267 + 97.9 99.6 275 0 1 0
UN185796 chr5A 585001939 585002769 - 100.0 99.6 745 3 0 0
UN175374 chr3D 34709434 34709880 . 100.0 98.4 440 7 0 0
UN184680 chr6A 603337083 603337320 . 100.0 100.0 238 0 0 0
UN088028 chr1A 581028511 581028919 . 97.4 100.0 409 0 0 0

augustus.gff

chr1A	AUGUSTUS	gene	4608	5033	0.95	+	.	g1
chr1A	AUGUSTUS	transcript	4608	5033	0.95	+	.	g1.t1
chr1A	AUGUSTUS	start_codon	4608	4610	.	+	0	transcript_id "g1.t1"; gene_id "g1";
chr1A	AUGUSTUS	CDS	4608	5033	0.95	+	0	transcript_id "g1.t1"; gene_id "g1";
chr1A	AUGUSTUS	exon	4608	5033	.	+	.	transcript_id "g1.t1"; gene_id "g1";
chr1A	AUGUSTUS	stop_codon	5031	5033	.	+	0	transcript_id "g1.t1"; gene_id "g1";
chr1A	AUGUSTUS	gene	7095	7403	0.97	+	.	g2
chr1A	AUGUSTUS	transcript	7095	7403	0.97	+	.	g2.t1
chr1A	AUGUSTUS	start_codon	7095	7097	.	+	0	transcript_id "g2.t1"; gene_id "g2";
chr1A	AUGUSTUS	CDS	7095	7403	0.97	+	0	transcript_id "g2.t1"; gene_id "g2";
chr1A	AUGUSTUS	exon	7095	7403	.	+	.	transcript_id "g2.t1"; gene_id "g2";
chr1A	AUGUSTUS	stop_codon	7401	7403	.	+	0	transcript_id "g2.t1"; gene_id "g2";
chr1A	AUGUSTUS	gene	7752	8069	0.59	-	.	g3
chr1A	AUGUSTUS	transcript	7752	8069	0.59	-	.	g3.t1
chr1A	AUGUSTUS	stop_codon	7752	7754	.	-	0	transcript_id "g3.t1"; gene_id "g3";
chr1A	AUGUSTUS	CDS	7752	8069	0.59	-	0	transcript_id "g3.t1"; gene_id "g3";
chr1A	AUGUSTUS	exon	7752	8069	.	-	.	transcript_id "g3.t1"; gene_id "g3";
chr1A	AUGUSTUS	start_codon	8067	8069	.	-	0	transcript_id "g3.t1"; gene_id "g3";
chr1A	AUGUSTUS	gene	19696	20019	0.75	+	.	g4
chr1A	AUGUSTUS	transcript	19696	20019	0.75	+	.	g4.t1
chr1A	AUGUSTUS	start_codon	19696	19698	.	+	0	transcript_id "g4.t1"; gene_id "g4";
chr1A	AUGUSTUS	CDS	19696	20019	0.75	+	0	transcript_id "g4.t1"; gene_id "g4";
chr1A	AUGUSTUS	exon	19696	20019	.	+	.	transcript_id "g4.t1"; gene_id "g4";
chr1A	AUGUSTUS	stop_codon	20017	20019	.	+	0	transcript_id "g4.t1"; gene_id "g4";
chr1A	AUGUSTUS	gene	22080	22388	1	+	.	g5
chr1A	AUGUSTUS	transcript	22080	22388	1	+	.	g5.t1
chr1A	AUGUSTUS	start_codon	22080	22082	.	+	0	transcript_id "g5.t1"; gene_id "g5";
chr1A	AUGUSTUS	CDS	22080	22388	1	+	0	transcript_id "g5.t1"; gene_id "g5";
chr1A	AUGUSTUS	exon	22080	22388	.	+	.	transcript_id "g5.t1"; gene_id "g5";
chr1A	AUGUSTUS	stop_codon	22386	22388	.	+	0	transcript_id "g5.t1"; gene_id "g5";
chr1A	AUGUSTUS	gene	22678	23055	0.95	-	.	g6
chr1A	AUGUSTUS	transcript	22678	23055	0.95	-	.	g6.t1
chr1A	AUGUSTUS	stop_codon	22678	22680	.	-	0	transcript_id "g6.t1"; gene_id "g6";
chr1A	AUGUSTUS	CDS	22678	23055	0.95	-	0	transcript_id "g6.t1"; gene_id "g6";
chr1A	AUGUSTUS	exon	22678	23055	.	-	.	transcript_id "g6.t1"; gene_id "g6";
chr1A	AUGUSTUS	start_codon	23053	23055	.	-	0	transcript_id "g6.t1"; gene_id "g6";
chr1A	AUGUSTUS	gene	24298	24594	0.58	+	.	g7
chr1A	AUGUSTUS	transcript	24298	24594	0.58	+	.	g7.t1
chr1A	AUGUSTUS	start_codon	24298	24300	.	+	0	transcript_id "g7.t1"; gene_id "g7";
chr1A	AUGUSTUS	CDS	24298	24594	0.58	+	0	transcript_id "g7.t1"; gene_id "g7";
chr1A	AUGUSTUS	exon	24298	24594	.	+	.	transcript_id "g7.t1"; gene_id "g7";
chr1A	AUGUSTUS	stop_codon	24592	24594	.	+	0	transcript_id "g7.t1"; gene_id "g7";

#!/usr/bin/env python
# encoding: utf-8

from collections import defaultdict
from multiprocessing import Pool, cpu_count
from functools import partial


def find_sth(f2, f1=None):
    start, end = int(f2[3]), int(f2[4])
    for uno1, start1, end1 in f1[f2[0]]:
        if (start1 <= start and start <= end1) or (start1 <= end and end <= end1) or (start1 >= start and end >= end1):
            with open("out.txt", "a") as fh:
                fh.write(uno1 + "\t" + f2[-1] + "\n")
                #print(uno1, f2[-1])


def main():
    with open('1.txt', 'r') as f1:
        infile1 = defaultdict(set)
        for uno1, chr1, start1, end1, *others in map(str.split, f1):
            infile1[chr1].add((uno1, int(start1), int(end1)))
    with open('genemark.gff3', 'r') as f2:
        infile2 = [x for x in map(str.split, f2) if x[2] == 'gene']
    pool = Pool(cpu_count())
    pool.map(partial(find_sth, f1=infile1), infile2)
    pool.close()
    pool.join()


if __name__ == "__main__":
    main()
from collections import defaultdict
from intervaltree import Interval, IntervalTree

with open('1.txt') as f:
    d1 = defaultdict(list)
    xs = map(lambda x: x.strip().split(), f)
    for x in xs:
        y = (x[0], int(x[2]), int(x[3]))
        d1[x[1]].append(y)

    for k, v in d1.items():
        d1[k] = IntervalTree(Interval(s, e, u) for u, s, e in v)

with open('genemark.gff3') as f:
    for line in f:
        line = line.strip().split()
        if line[2] == 'gene':
            chr, start, end = line[0], int(line[3]), int(line[4])
            for start1, end1, un1 in d1[chr][start-1:end+1]:
                print(un1, line[-1])
#!/usr/bin/env python
# -*- coding: utf-8 -*-
infile2 = open('augustus.gtf', 'r')
infile1 = set(line1.strip() for line1 in open('1.txt', 'r'))
for line in infile2:
    line = line.strip().split()
    if line[2] == 'transcript':
        chr, start, end = line[0], int(line[3]), int(line[4])
        for line1 in infile1:
            line1 = line1.split()
            chr1, start1, end1 = line1[1], int(line1[2]), int(line1[3])
            if chr1 == chr:
                if start1 < start < end1:
                    print line1[0], line[-1]
                    break
                if start1 < end < end1:
                    print line1[0], line[-1]
                    break
            else:
                continue
msw521sg
关注
  • 0
    点赞
  • 踩
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
提取基因序列文件
11-28
perl可编译代码，可用于提取目的序列文件中所需序列id
python解析gff文件中的转录本
Cassiel60的博客
05-27 6503
1.下载基因组注释文件，选择对应的版本：ftp://ftp.ncbi.nlm.nih.gov/genomes/Homo_sapiens/ARCHIVE/BUILD.37.3/GFF/ 2.GTF 为General Transfer Format ，熟悉格式http://www.huoyunjn.com/wuliuxinwen/2/33709819.htm。 第三列feature- 后面...
python提取GFF3文件信息
山海之间
02-22 8717
前言因为最近使用REPET这个程序对基因组进行重复序列的注释，但是最后输出的结果是GFF3格式的文件，缺少统计信息。因此用python写了个脚本，对GFF3的信息进行提取并统计。GFF3文件想要从GFF3文件中提取的信息为重复序列的种类，数量，以及bp数。 GFF3文件分为9列，这次用到是第2、4、5、9列，分别代表来源，序列起始位置，结束位置，以及属性。 重复序列的种类可以从第9列属性中提取，
根据GFF3文件统计外显子大小和数量以及内含子大小
xxxxx
01-17 1万+
根据GFF3文件统计外显子大小和数量以及内含子大小#!/usr/bin/env python # -*- coding: utf-8 -*- __author__ = "Sheng-Wei Ma" with open('TGACv1.cdna.gff3', 'r') as f: for line in f: lin = line.strip().split('\t')
Lecture-GFF:读取.gff文件，目前为NCBI
04-06
这段代码首先导入了`gffutils`库，然后创建了一个数据库连接，将GFF文件中的数据加载到内存中。`features_of_type`函数用于筛选特定类型的特征，这里我们选择了'gene'。最后，通过遍历所有基因特征，可以访问每个...
gff2gtf:用于将 GFF 文件转换为 GTF 的 AWK 脚本
06-10
awk -f gff2gtf.awk 文件.gff > 文件.gtf 适用于那些具有 AWK 解释器的操作系统。 该脚本依赖于记录的顺序。 规范不要求对 GFF 进行排序，但通常对记录进行排序。 对于 RSEM： awk -f gff2gtf-strict.awk FILE....
gffcompare：通过与参考注释GFF进行比较来对GFF文件进行分类，合并，跟踪和注释
02-11
对一个或多个GTF / GFF3文件中的笔录进行分类，因为它们与注释文件中提供的参考笔录有关（也为GTF / GFF3格式） 可以在以下文章中找到更多详细信息和用法示例： ，也可以引用该软件。 可在此处在线找到此实用程序...
gffutils:GFF和GTF文件操作和互转换
05-25
gffutils是一个Python软件包，用于处理和处理通常用于基因组注释的GFF和GTF格式文件。 文件被加载到sqlite3数据库中，与仅使用纯文本方法相比，可以对层次结构特征（例如，基因，转录本和外显子）进行更为复杂的操作...
nwn_langedit:一个小的gff文件语言编辑器（希望如此）
04-02
"nwn_langedit"就是这样一款专为编辑GFF文件而设计的小型工具，它由AutoHotkey编程语言编写，旨在简化游戏中的语言文件处理过程。 GFF，全称为“Generic Formatted File”，是一种广泛用于BioWare公司开发的游戏，...
基因组注释文件(GFF,GTF)下载的五种方法
热门推荐
Baimoc
04-17 5万+
文章目录NCBIEnsemblUCSCGeneCode NCBI Ncbi 里包含现在最全的参考基因组数据，可以进入FTP站点查看：ftp://ftp.ncbi.nlm.nih.gov/genomes/ 这里的文件夹名为物种的拉丁名，这里以 Human(Homo_sapiens) 为例，下载方法如下： wget ftp://ftp.ncbi.nlm.nih.gov/genomes/Homo_sa...
gff数据处理
ruby912的博客
06-04 3947
Perl脚本练习 通过数组和哈希引用读取存储复杂数据 要求 gff文件，记录每个基因的名称、起始终止位置、染色体、转录本信息 数据 #注释行 分割：\t 9列：gff格式 思路 构建一个哈希，结构 my%gene = ( geneID => { 'location' => [chrs, start, end] 'transcripts' => [ { ...
gff文件_pyhton读写数据文件
weixin_36464343的博客
12-30 600
1、打开文件获取文件对象2、操作文件3、关闭文件fr = open("test.txt","r") #打开文件ff=fr.read()#读取文件所有内容(不建议使用，如果文件内容巨大，内存会爆)print(ff)文件类型r#只读，默认模式 打开数据文件w#只写，不可读，若文件不存在则创建，若存在，则删除内容，写入新内容a #只追加，不可读，若文件不存在则创建，存在则追...
嵌合体序列
weixin_34195546的博客
03-03 949
嵌合体序列由来自两条或者多条模板链的序列组成，示意图如下： 在PCR反应中，在延伸阶段，由于不完全延伸，就会导致嵌合体序列的出现，以上图为例， 在扩增序列X的过程中，在序列延伸阶段，只产生了部分X序列延伸阶段就结束了，在下一轮的PCR反应中，这部分序列作为序列Y的引物接着延伸，扩增就会形成X和Y的嵌合体序列； 在放一张具体一点的示意图，不完全延伸产生的序列作为下一轮PCR反应的产物，进行延伸 ...
gff文件_GTF与GFF文件所有操作尽在这里
weixin_39640417的博客
12-18 8021
一、格式介绍(一)gtf文件。GTF 为General Transfer Format缩写，跟 GFF2格式类似。相信大家做转录组分析时候经常会看到Cufflinks或者Stringtie软件对转录组进行定量与组装会时产生一个gtf文件，里面包含的信息如下：每列信息的含义如下：seqname- 序列的ID，可以是染色体的ID也可以是Scaffold或者Contig的ID。source- 产生此文件...
从gff3中获取fasta序列
xxxxx
09-05 2267
#!/usr/bin/env python # -*- coding: utf-8 -*-from Bio import SeqIO record_dict = SeqIO.index(open('genome.fa', 'r'), "fasta") tmp = open('2.txt', 'w') gene = {} mRNA = {} exon = {} CDS = {} gene_sequen
两个文件的合并
xxxxx
02-18 2588
两个文件的合并分别有文件一和文件二，文件一和文件二都有一列内容是相同的。根据此列内容可将两者结合起来。文件一的格式如下 AEGTA04333 LOC_Os01g33520|LOC_Os01g33520.1 2.70E-25 LOC_Os01g33520|LOC_Os01g33520.1 PF02902 Ulp1 protease family, C-terminal cat
根据位置从gff中提取gene
GRC
06-29 1万+
小程序-根据位置从gff中提取gene Table of Contents 1 需要的文件格式2 程序 1 需要的文件格式 1: 染色体 位置 2 程序 1: use strict; 2: use warnings; 3: 4: my (@information,$chr,%hash,$key1,$key2,); 5
从gff3文件获取fasta序列(2)
xxxxx
09-06 4033
get fasta sequence from gff3#taken from https://www.biostars.org/p/46281/ and modified as needed 10/15/15 use strict; use warnings; use Bio::Seq; use Bio::SeqIO; use Bio::DB::Fasta;#add a help message
根据基因组文件和gff3文件提取启动子序列Python
最新发布
05-28
可以使用Biopython和pandas库来解析基因组文件和gff3文件，并提取启动子序列。 首先，需要从基因组文件中读取DNA序列。假设基因组文件是fasta格式，可以使用Biopython中的SeqIO模块读取序列： ```python from Bio ...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值