linux提取fasta文件的id,FASTA序列文件处理一网打尽
推荐两个地方:
地方一都是小脚本,但实用,大伙也可以自己练习写。
地方二成熟软件SeqKit,也很实用。
一、小脚本
大家可以在这里下载以下脚本:
https://github.com/jorvis/biocode/tree/master/fasta
各脚本作用信息如下:
|-- append_to_fasta_header.py
每个序列ID添加后缀
|-- check_for_embedded_fasta_headers.py
检查此种错误类型fasta序列
>gi1006569 DnaA [Vibrio harveyi]
MSSSLWLQCLQQLQEELPATEFSMWVRPLQAEVLHAC>gi409247 DNA repair protein
MVSLTFKNFKKEKVPLDLEPSNTILETKTKLAQSISCEESQIKLIYSGKVLQDSKTVSECGLKDGDQVVF
|-- compare_two_fastas.pl
比较两个文件相同序列与不同序列数目
|-- convert_fasta_contigs_to_gff3.py
将contig序列转换为GFF格式
|-- convert_fastq_to_fasta.py
将FASTQ转换为FASTA
|-- create_fasta_pseudomolecules.pl
从遗传图或者Hi-C确定的map信息文件中获得假染色体序列(将contig连成染色体序列)
|-- extract_fasta_regions.py
提取特定区域的序列
|-- fasta_base_content.py
碱基含量统计
|-- fasta_size_distribution.pl
GC含量和长度分布统计
|-- fasta_size_distribution_plot.py
序列长度分布作图
|-- fasta_size_report.pl
序列长度信息统计
|-- filter_fasta_by_header_regex.py
使用正则表达式获取FASTA序列
|-- filter_fasta_by_ids.pl
根据ID提取序列
|-- filter_fasta_by_size.pl
根据长度过滤序列
|-- merge_fasta_files_and_uniquify_ids.py
你有多个序列ID相同的文件,想合并一起,保证序列ID唯一性
>CL1Contig1
>CL1Contig2
|-- merge_masked_fasta_files.py
将不同软件冰壁重复序列结果合并
|-- prepend_to_fasta_header.py
你有多个序列ID相同的文件,想合并一起,在序列前面添加前缀,保证序列ID唯一性
|-- reformat_fasta_residue_lengths.py
按照一行固定字符数目重新显示FASTA序列
|-- remove_duplicate_sequences.py
移除掉重复的序列
|-- remove_empty_sequences.pl
移除掉文件中空的序列
|-- reorient_fasta.pl
从某一位置截取序列到另一个文件中
|-- reorient_sequences_by_id.py
将指定序列反向或者反向互补
|-- split_fasta_into_even_files.py
分割成指定数目的文件
|-- split_interleaved_sequence_file.py
将FASTQ/FASTA序列中成对的分成R1,R2两个文件,单独非成对的生成一个文件
|-- subsample_fasta.py
随机抽取指定数目的FASTA序列
|-- validate_fasta.py
主要检查以下内容:
- ERROR: Entry with 0-length sequence
- ERROR: No > symbols embedded in sequence residues (suggests merged records)
- ERROR: Nonzero-length identifier after > symbol until first whitespace
- ERROR: Multiple records within an individual file with same identifier
- ERROR: Incorrect record count (if --expected_record_count is passed)
- ERROR: Make sure the file ends with a blank line (else commonly used tools like cat will cause errors)
- WARNING: Residue lines > 60bp
- WARNING: Comment lines with # (most parsers will fail to handle this)
- WARNING: Homopolymers of length > --homopolymer_limit which are not Ns
- WARNING: Internal stops (*) detected in protein FASTA (if --check_internal_stops is passed)
`-- write_fasta_from_gff.pl
从GFF提取CDS和蛋白序列
二、SeqKit
处理FASTA和FASTQ的神器,window\linux系统版本都有.对于没有编程基础的小伙伴们,我们照样可以轻松操作序列文件。
该软件功能强大,小编只罗列部分模块功能,更详细功能参见软件网站:
http://bioinf.shenwei.me/seqkit/usage/。
一、序列操作。
seqkit seq [flags] file
参数:
-p, --complement 取互补序列
--dna2rna DNA to RNA
-l, --lower-case 将序列以小写字母形式输出
-g, --remove-gaps 移除组装序列中的gap
-r, --reverse 取反向序列
--rna2dna RNA to DNA
-u, --upper-case 将序列以大写字母形式输出
-w, --line-width int 以每行指定长度输出序列 (0 for no wrap) (default 60)
举例:
seqkit seq test.fa -w 0#将此文件fasta序列转换成一行输出
seqkit seq -w 100 test.fa#将此文件fasta序列转换成100个碱基一行输出
seqkit seq --dna2rna test.fa#将此文件fasta序列dna转换成rna
seqkit seq -w 100 -p -r test.fa#将此文件fasta序列反向互补输出,每行100碱基
二、Fasta/q之间及与tab格式互换
1、FASTQ转换成FASTA: seqkit fq2fa
举例:
seqkit fq2fa reads_1.fq -o reads_1.fa
2、FASTA/FASTQ转换成tab格式。seqkit fx2tab
举例:
seqkit fx2tab test.fa>test.fa.tab.fa
seqkit fx2tab test.fq>test.fq.tab.fq
tab格式:ID sequence
三、序列信息统计
1、序列碱基含量及序列长度信息统计
seqkit fx2tab [flags]
参数:
-B, --base-content value 要输出的碱基含量e.g. -B AT -B N
-g, --gc print GC content
-l, --length print sequence length
-n, --name only print names
-i, --only-id print ID instead of full head
举例:
seqkit fx2tab -l -g -n -i -H test.fa
输出结果:
#name seq qual length GC
gene1 30 40.00
2、序列长度分布统计
Usage:
seqkit stat [flags]
举例:
seqkit stat test.fa
输出结果
file format type num_seqs sum_len min_len avg_len max_len
test.fa FASTA DNA 1 30 30 30 30
四、根据ID或特定的motif筛选提取序列
seqkit grep [flags]
参数:
-n, --by-name 匹配整个序列的名字,包含description部分,而不是序列id。
-s, --by-seq 匹配序列
-d, --degenerate pattern/motif 包含简并碱基
-i, --ignore-case 忽略大小写
-v, --invert-match 输出不匹配此模式的内容
-p, 匹配模式,支持连续写多个模式,匹配任一模式即输出。如-p ^ATG -p TAA$。注意该功能仅能正向匹配,不能实现对互补链匹配。
-f, --pattern-file string 支持匹配模式写到一个文件中,如要提取的序列ID。
-R, --region string 匹配位置选择。e.g 1:12 for first 12 bases, -12:-1 for last 12 bases
-r, --use-regexp 使用正则表达式,必须加入此参数,如^匹配首端。同-p联合使用。
示例:
seqkit grep -s -r -i -p ^atg cds.fa#选取有起始密码子的序列
seqkit grep -f list test.fa > new.fa#根据ID提取序列
seqkit grep -s -d -i -p TTSAA#简并碱基使用。S 代表C or G.
seqkit grep -s -R 1:30 -i -r -p GCTGG##匹配限定到某区域
五、motif定位
对grep的拓展,可以正反链同时匹配,输出匹配的位置。
seqkit locate [flags]
参数
-d, --degenerate pattern/motif contains degenerate base
-i, --ignore-case ignore case
-P, --only-positive-strand only search at positive strand
-p, --pattern value search pattern/motif
-f, --pattern-file string pattern/motif file (FASTA format)
举例
seqkit locate -i -d -p AUGGACUN test.fa
输出结果
seqID patternName pattern strand start end matched
cel-mir-58a AUGGACUN AUGGACUN + 81 88 AUGGACUG
ath-MIR163 AUGGACUN AUGGACUN - 122 129 AUGGACUC
六、多个序列文件比较寻找相同的序列或者ID相同的序列
seqkit common [flags]
参数:
-n, --by-name 匹配整个序列的名字,包含description部分,而不是序列id
-s, --by-seq match by sequence
-i, --ignore-case ignore case
-m, --md5 use MD5 reduce memory usage
举例:
1、By ID (default,>后面,空格之前的名字)输出ID名字相同的。
seqkit common test1.fa test2.fa -o common.fasta
2、By full name(整个序列的名字,包含description部分)。输出序列名字相同的。
seqkit common test1.fa test2.fa -n -o common.fasta
3、输出要比较的文件中序列相同的序列
seqkit common test1.fa test2.fa -s -i -o common.fasta
4、输出要比较的文件中序列相同的序列 (for large sequences)
seqkit common test1.fa test2.fa -s -i -o common.fasta --md5
七、提取部分序列
如随机抽取10000条FASTQ序列做NT污染评估。同时他也可以对FASTA序列提取
seqkit sample [flags]
参数:
-n, --number int sample by number (result may not exactly match)
-p, --proportion float sample by proportion
-s, --rand-seed int rand seed for shuffle (default 11)
-2, --two-pass 2-pass mode lower memory
举例:随机抽取序列
seqkit sample -n 10000 -s 11 test1_1.fq -o sample.fq
seqkit sample -p 0.1 -s 11 test1_1.fq -o sample.fq
八、排序输出命令
seqkit sort [flags]
参数:
-l, --by-length 按照序列长度排序
-n, --by-name by full name
-s, --by-seq 按照序列排序
-i, --ignore-case 按序列排序时忽略大小写
-r, --reverse 反向排序
-2, --two-pass 对于FASTA序列排序可以减少内存
举例:
seqkit sort -l test.fa
九、文件切割
seqkit split [flags]
参数:
-i, --by-id split squences according to sequence ID
-p, --by-part int 将一个文件分割成N 份
-s, --by-size int 将一个文件按照N 条序列一个文件进行分割
-O, --out-dir string output directory (default value is infile.split)
-2, --two-pass two-pass mode to lower memory usage(only FAST)
举例:
seqkit split hairpin.fa.gz -p 4
linux提取fasta文件的id,FASTA序列文件处理一网打尽相关推荐
- php 蓝奏网盘上传文件,蓝奏云_文件上传_API
~~~[api] post:https://api-lanzou.anyhulian.ml/v2/upload.php *content-type=multipart/form-data; bound ...
- linux提取fasta文件的id,从大的fasta文件中提取特定的fasta序列
我想使用以下脚本从大的fasta文件中提取特定的fasta序列,但输出为空.从大的fasta文件中提取特定的fasta序列 transcripts.txt文件包含我想从assembly.fasta到s ...
- linux中fq格式转fa,fasta/fq文件处理万能工具——Seqkit学习记录
shenwei爪哥开发的处理Fasta/Fastq文件的万能工具.之前处理fq/fa文件时花时间写的一些脚本发现在seqkit里直接能一行命令就解决.实在是提升效率,整合流程中十分好的工具.本文是对S ...
- linux系统fasta程序,快速计算fasta序列长度的方法
最近看了一下进入PLoB的网页来路分析,看到有同学搜索计算fasta序列长度.其实自己在之前的数据分析中也遇到过相关的问题,这里给大家分享两种我常用的方法. 方法一:linux下用awk计算fasta ...
- linux提取基因名称和序列,一种批量提取基因组基因信息并翻译比对分析序列的方法与流程...
技术特征: 1.一种批量提取基因组基因信息并翻译比对分析序列的方法,其特征在于,将某一物种的转录本id或者基因id,依据供试基因组cds文件.蛋白质文件.gff文件和染色体fasta文件信息,通过6个 ...
- 新手写python脚本提取最长转录本ID以及最长转录本序列(此代码已过时,最新代码见主页)
原始蛋白序列长这个样子: 写脚本的具体思路: 1.创建转录本ID和序列的字典(键值对) 2.提取每个转录本的序列长度,形成三列,第一列是转录本ID,第二列是转录本长度,第三列是gene ID 3.这一 ...
- linux 提取cpio_15. Linux提取RPM包文件(cpio命令)详解
在讲解如何从 RPM 包中提取文件之前,先来系统学习一下 cpio 命令. cpio 命令用于从归档包中存入和读取文件,换句话说,cpio 命令可以从归档包中提取文件(或目录),也可以将文件(或目录) ...
- Linux提取文件夹中文件名,linux命令提取文件夹内特定文件的路径与文件名
最近需要实现自动化搜寻特定文件夹下的特定文件,并且需要分别保存文件路径与文件名.算然使用python的walk能够实现,但是感觉复杂了些.于是想看看linux自带的命令是否能完成这项工作. 环境 需要 ...
- linux 提取cpio_【rpm】从rpm包中提取文件:rpm2cpio和cpio的使用
rpm2cpio命令可以用于将rpm格式的文件转为cpio格式的文件. rpm是Linux中常用的文件格式,方便了用户的安装,但没有cpio格式灵活. cpio是用来建立,还原备份档的工具程序,它可以 ...
最新文章
- 构建企业级业务高可用的延时消息中台
- 中国工程院《全球工程前沿2020》报告在京发布
- 按15分钟取数据_步行15分钟能获得什么?这组数据能告诉你……
- 阿云搭建php博客,阿云的应聘故事
- picpick尺子像素大小精度不够准确_矿用电子皮带秤该如何维护,以提高使用精度?...
- 服务器CPU X86 ARM PowerPC RISC介绍
- 甘特图:有效的项目管理计划怎么制定?
- multimedia教学设计计算机英语,高中英语教学设计案例
- 基带、频带、宽带、带宽
- 人工神经网络的结构基本上分为两类
- 【Python百日基础系列】Day03 - Python 数据类型
- 华为手机左侧快捷方式_华为手机的这六个快捷键,让使用更简便!
- 怎么做一份漂亮的地质图
- Java 从lambda 表达式引用的本地变量必须是最终变量或实际上的最终变量问题解决
- mysql的467_数据库报错排解
- linux mt命令,Cheat—— 给Linux初学者和管理员一个终极命令行备忘单
- 微信小程序调用阿里OCR识别
- 绿盾启用限制终端无法上网策略,关闭策略后终端不能正常恢复上网功能
- secKill项目 --- 总结 + 推荐阅读顺序 + 源码地址
- EGE基础入门篇(二):开始使用EGE