如何统计id很复杂的fasta文件的长度？

最新推荐文章于 2022-03-18 17:33:09 发布
原创 最新推荐文章于 2022-03-18 17:33:09 发布 · 1.9k 阅读 · 0 ·
CC 4.0 BY-SA版权
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
文章标签：

#linux #perl #生物信息 #统计

本文探讨了在处理含有复杂ID的FASTA文件时遇到的问题，原本针对简单ID的脚本无法正确统计复杂的基因组数据，导致染色体数量统计出现偏差。为解决此问题，需要对现有脚本进行调整。

摘要生成于 C知道 ，由 DeepSeek-R1 满血版支持， 前往体验 >

对于一般的fasta文件的格式是：

>chr1

AAGCCATCCGG

但是最近两天遇到id很复杂的fasta，对于这样的fa文件用现有的脚本却统计错误


原本人的hg38染色体有23条染色体，但是统计出来却只有一条染色体的长度：


因此，需要重新修改length.pl

用python编写统计fasta格式的序列的长度脚本
tangxc10的专栏
09-30 9435
如果用perl来编写统计fasta序列的长度脚本，很简单的几行代码就可以搞定，但是想了想，觉得用python写更时候处理大的文件，尤其是想用python实现多线程处理。因此，就有了用python来编写最初版的统计fasta序列长度的脚本的想法。  运行方法：nohup python stat_length.py input.fasta > input.len & 运行结果： >Aquca_0
生信脚本练习（12）求fasta文件各序列长度并统计作图
49
08-20 6636
题目要求是要从一个fasta文件中统计出每条序列的长度分布，并作图。代码如下：import os import getpass import matplotlib.pyplot as plt usr = getpass.getuser() os.chdir('c:/Users/' + usr + '/Desktop') seq_len = {} # 把fasta文件全部读取做成字典，键是带‘>’的那
使用Python计算fasta文件的序列长度
zhangyikai2017的博客
11-30 4234
使用Python计算fasta文件的序列长度 在这里插入代码片使用Python计算fasta文件的序列长度 #！/usr/bin/python #-- coding:utf-8 -- import sys f = open(sys.argv[1],‘r’) out = open(sys.argv[2],‘w’) def chr_length(infile,outfile): f = open(sy...
linux统计染色体长度,利用bedtools统计VCF文件中各染色体不同区域中的变异数量...
weixin_31834275的博客
05-16 2912
1.获得基因组的各染色体的长度信息$ samtools faidx genome.fa#生成gemome.fa.fai文件第一列为染色体，第二列为对应染色体长度chr01 44488843 7 44488843 44488844chr02 38522657 44488858 38522657 38522658...
linux提取fasta文件的id,FASTA序列文件处理一网打尽
weixin_33093403的博客
05-14 5809
推荐两个地方：地方一都是小脚本，但实用，大伙也可以自己练习写。地方二成熟软件SeqKit，也很实用。一、小脚本大家可以在这里下载以下脚本：https://github.com/jorvis/biocode/tree/master/fasta各脚本作用信息如下：|--append_to_fasta_header.py每个序列ID添加后缀|--check_for_embedded_fasta_he...
用Python语言编写程序，统计Botrytis_cinerea_protein.fa文件中蛋白序列的长度，并将长度小于300的蛋白序列以fasta的格式输出文件
03-13
好的，我现在需要帮用户解决如何使用Python处理FASTA文件，统计蛋白序列长度并过滤保存短于300个氨基酸的序列的问题。首先，我得回忆一下FASTA文件的格式。FASTA文件通常以“>”开头作为序列标识符的行，接着是换行...
使用seqkit处理fasta文件，来进行EBI远程interproscan
最新发布
08-02
如果FASTA文件中的序列很长或者数量较多，建议分批提交或者使用本地InterProScan（但本地安装复杂，资源消耗大）。 以上步骤完成了使用seqkit处理FASTA文件并调用EBI远程服务进行InterProScan分析。 相关问题：...
import os from Bio import SeqIO from Bio.Seq import Seq from Bio.SeqRecord import SeqRecord # 定义文件路径 genome_file = r"C:\Users\Lenovo\Desktop\Botrytis_cinerea_genome.fa" gff_file = r"C:\Users\Lenovo\Desktop\Botrytis_cinerea.gff3" output_file = r"C:\Users\Lenovo\Desktop\Botrytis_cinerea_genes.fasta" # 读取基因组文件 genome_dict = SeqIO.to_dict(SeqIO.parse(genome_file, "fasta")) # 解析GFF3文件，提取基因位置信息 def parse_gff(gff_file): genes = {} with open(gff_file, "r") as gff: for line in gff: if line.startswith("#"): continue # 跳过注释行 fields = line.strip().split("\t") if len(fields) < 9: continue # 跳过格式不正确的行 seqid, source, feature, start, end, score, strand, phase, attributes = fields if feature == "gene": gene_id = None for attr in attributes.split(";"): if attr.startswith("ID="): gene_id = attr.split("=")[1] break if gene_id: genes[gene_id] = { "seqid": seqid, "start": int(start), "end": int(end), "strand": strand } return genes # 提取基因序列 def extract_genes(genome_dict, genes): gene_sequences = [] for gene_id, gene_info in genes.items(): seqid = gene_info["seqid"] start = gene_info["start"] end = gene_info["end"] strand = gene_info["strand"] if seqid not in genome_dict: print(f"Warning: Sequence ID {seqid} not found in genome file.") continue gene_seq = genome_dict[seqid].seq[start - 1:end] # GFF3位置是1-based if strand == "-": gene_seq = gene_seq.reverse_complement() gene_record = SeqRecord(gene_seq, id=gene_id, description="") gene_sequences.append(gene_record) return gene_sequences 逐行修改以上文件
03-28
比如检查FASTA记录的ID和序列长度，GFF3记录的属性和字段是否正确解析。可以用assert语句或单元测试框架。 最后，整合FASTA和GFF3数据，比如根据GFF3中的位置信息提取对应的序列区域。这需要将序列数据与注释结合，...
bam获取序列_如何从BAM文件中提取fastq
weixin_39722196的博客
12-20 2443
虽然高通量测序分析最常用的操作是将fastq比对到参考基因组得到BAM文件，但偶尔我们也需要提取BAM文件中特定区域中fastq。最开始我认为这是一个非常简单的操作，因为samtools其实已经提供了相应的工具samtools fastq.以biostar handbook的Ebola病毒数据为例，首先获取比对得到的BAM文件。# 建立文件夹mkdir -p refs# 根据Accession下载...
【生物序列数据的处理】FASTA与GENBANK格式的定义与特点
[【生物序列数据的处理】FASTA与GENBANK格式的定义与特点](https://www.gensas.org/sites/default/files/images/Fig6.jpg) # 1. 生物序列数据处理概述 生物序列数据处理是生物信息学领域的一项基础而关键的技术，它...
使用perl计算fasta序列长度
weixin_43840576的博客
12-03 3509
Fasta序列长度的计算的原理就是将不同染色体长度存入哈希，最后用length去计算哈希元素的长度。 perl脚本如下： #!/usr/bin/env perl use strict; my($name, %seq ); open IN, $ARGV[0]; die “Need the inut fasta.” if( ! defined $ARGV[0] ); while(){ chomp; i...
python 从fastq文件中挑选出序列长度在规定范围的序列
niuhuihui_fei的博客
05-23 7562
python 从fastq文件中挑选出序列长度在规定范围的序列 FASTQ文件长度过滤
根据序列ID提取fasta序列
热门推荐
周欣的博客
06-15 2万+
根据序列ID快速抽提fasta序列在进行完序列比对以及，在微生物多样性分析中经常要根据物种信息抽提特定ID的fasta格式的序列文件。传统方法大家肯定是一个个ID去查找，复制，粘贴。这种方法虽然很精确但是效率实在是太低了。如果能用几行命令解决不仅简单高效，而且还能一劳永逸，妈妈再也不用担心我花几天时间在查找序列上了。下面就正式教大家，怎么把两个文件名字不一样的文件中特定想要的序列文件一一匹配并提取...
根据ID从FASTA文件中批量提取序列【Python】
duoluka的博客
03-18 4874
当你有一个ID文件，该如何批量提取所对应的fasta序列？
根据ID列表提取fasta序列
薛猫_柳叶刀
08-08 1万+
根据基因列表，从总的fasta文件中提取
人类基因组测序文件chromFa.tar （hg19），统计每条染色体的长度，GC含量，N比值，N、A、T、C、G各碱基数量。
qq_38790244的博客
08-16 4242
import os import time begin = time.time() os.chdir(r'G:\tem\fa') def count_n_and_gc(file): content = [] chromsome = [] g = 0; c = 0; n = 0; a = 0; t = 0 with open(file) as f: r...
C语言读取fasta中核酸序列
in_the_shadow的博客
05-30 3345
C语言是面向过程的一门语言，本文程序用C语言编写，实现了生物信息学中fasta文件的读取，并读取其中核酸序列以及一小段描述符，并将其地址存取于指针变量中，方便对其中每个字符进行操作。FASTA格式（又称为Pearson格式），是一种基于文本用于表示核苷酸序列或氨基酸序列的格式。在这种格式中碱基对或氨基酸用单个字母来编码，且允许在序列前添加序列名及注释。
fasta多行文件处理
weixin_30562507的博客
10-23 665
1.创建fa文件，如下，命令为1.fa>SOX2 ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGAC>POU5F1 CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGTCGGAAGGTAGTCGTCAGTGCAGCGAGT...
Python处理多行文本问题--一个简单方法读取多行fasta文件
weixin_33984032的博客
03-10 1383
在处理fasta序列时，常常会遇到一条序列多行排列的现象，如下所示： $cat test.fasta >test_1 TGGGGAATCTTGGACAATGGGGGCAACCCTGATCCAGCCATGCCGCGTGAGCGATGAAGGCCTTAGGGTTGTAAAGCTCT TTCAGCTGGGAAGATAATGACGGTACCAGCAGAAGAAGCCCCGGCTAACTC...
利用hifiasm组装基因组，怎么bam文件转fasta文件？
12-18
对于`.cns`文件，你可以找到对应的`.fa`或`.fasta`文件，如果没有，可以直接通过`cns`文件生成。例如，命令如下： ``` samtools faidx your_contigs.cns > your_contigs.fasta ``` 3. 确保你的`samtools`已经...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

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

抵扣说明：

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

余额充值