利用samtools截取基因组上指定位置的序列
利用samtools截取指定位置的序列
samtools faidx 能够对fasta 序列建立一个后缀为.fai 的文件，根据这个.fai 文件和原始的fastsa文件， 能够快速的提取任意区域的序列
用法：
samtools faidx input.fa
该命令对输入的fasta序列有一定要求：对于每条序列，除了最后一行外， 其他行的长度必须相同:
>chr1
ATGCATGCATGCATGCATGCATGCATGCAT
GCATGCATGCATGCATGCATGCATGCATGC
ATGCAT
>chr2 another chromosome
ATGCATGCATGCAT
GCATGCATGCATGC
可以使用seqtk转换：
seqtk seq -l 50 -L 50 input.fa >input1.fa
最后生成的.fai文件如下， 共5列，\t分隔；
chr1 66 5 30 31
chr2 28 98 14 15
第一列 NAME : 序列的名称，只保留“>”后，第一个空白之前的内容；
第二列 LENGTH: 序列的长度， 单位为bp；
第三列 OFFSET : 第一个碱基的偏移量， 从0开始计数，换行符也统计进行；
第四列 LINEBASES : 除了最后一行外， 其他代表序列的行的碱基数， 单位为bp；
第五列 LINEWIDTH : 行宽， 除了最后一行外， 其他代表序列的行的长度， 包括换行符， 在windows系统中换行符为\r\n, 要在序列长度的基础上加2；
提取序列：
samtools faidx input.fa chr1 > chr1.fa
samtools faidx input.fa chr1:10-20 > chr1.fa
更多生物信息课程：https://study.omicsclass.com/index
- 发表于 2020-02-21 19:36
- 阅读 ( 14828 )
- 分类：软件工具
