序列的提取是我们生物学分析的基础操作，在许多分析中都会使用到。在这里给大家提供几个提取序列的命令。
序列的提取是我们生物学分析的基础操作，在许多分析中都会使用到。在这里给大家提供几个提取序列的命令。
1. seqkit
使用 seqkit subseq 通过指定区域提取序列
-r 通过区域来截取序列，提取 in.fa 前12个字符，可以用--chr指定染色体
seqkit subseq -r 1:12 data/in.fa > out.fa
seqkit subseq -r 1:12 data/genome.fa --chr chrI
seqkit subseq --bed data/gene.bed data/genome.fa > out.fa
seqkit subseq --gtf data/gtf data/genome.fa > out.fa
使用 seqkit grep 命令，可以根据序列名称来提取特定的序列，也支持从文件列表中提取序列。
单个提取
seqkit grep --pattern YAL069W data/gene.fa > out.fa
批量提取 (list：ID列表 )
seqkit grep -f data/list data/gene.fa > out.fa
2. seqtk
使用 seqtk subseq 命令，可以根据序列名称来提取特定的序列，也支持从bed中提取序列。
位置信息文件 ( gene.bed )
seqtk subseq data/genome.fa data/gene.bed > out.fa
序列名称（ list：ID列表 ， -l可设定输出的每行长度）
seqtk subseq -l 80 data/gene.fa data/list > out.fa
使用 seqtk trimfq 截取序列
切除前5个和后10个碱基
seqtk trimfq -b 5 -e 10 data/in.fa > out.fa
3. Samtools
samtools faidx data/genome.fa chrI:1-20 > out.fa
4. Bedtools
awk ‘{ if($3==“gene”) print $1“\t”$4“\t”$5“\t”$NF}’ data/gff | cut -d “;” -f1 | sed “s/ID=//g” > data/gene.bed # 提取gene的bed文件
bedtools getfasta -fi data/genome.fa -bed data/gene.bed -nameOnly > gene.fa
