推荐两个地方: 地方一都是小脚本,但实用,大伙也可以自己练习写。 地方二成熟软件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./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碱基 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 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##匹配限定到某区域 对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 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 |
|
来自: 海陵岛的鱼 > 《Bioinformation》