生物信息学工具--bowtie&bowtie2

时间:2022-09-26 09:08:54

Bowtie和Bowtie2使用

【怪毛匠子整理】

Source URL: http://www.bbioo.com/lifesciences/40-112837-1.html

Bowtie和Bowtie2使用

碱基 序列 种子 前导链 错配 基因组 末端

标题: Bowtie和Bowtie2使用

摘要: [Bowtie和Bowtie2使用]bowtie 比对http: bowtie-bio sourceforge net index shtmlhttp: www ncrna net bowtie-short-sequence-alignment-tool-detailed-solution http: www plob org 2011 12 13 932 html 常见的短序列比对工具有很多,如fasta、blast、bowtie、shrimp、s…… [关键词:碱基 序列 种子 前导链 错配 基因组 末端]……

关键词: 碱基 序列 种子 前导链 错配 基因组 末端

bowtie 比对

http://bowtie-bio.sourceforge.net/index.shtml

http://www.ncrna.net/bowtie-short-sequence-alignment-tool-detailed-solution/

http://www.plob.org/2011/12/13/932.html

常见的短序列比对工具有很多,如fasta、blast、bowtie、shrimp、soap等。每个工具都有其自身的优点,但同时也具备了一些缺点。权衡利弊,我选择bowtie作为主要的短序列比对工具。它速度很快,比对结果也容易理解。

现在举个例子来探讨bowtie的使用方法:现在有GENOME.fa、高通量测序数据Reads.fa,我们希望将Reads.fa比对到基因组GENOME.fa上。

(一)、对Reference文件(GENOME.fa)建库

1

bowtie-build GENOME.fa GENOME.fa

$ bowtie2-build genome.fasta index

建库步骤可能需要1h甚至几个小时,建议在后台执行:

nohup bowtie-build GENOME.fa GENOME.fa &

(二)、将Reads.fa比对到GENOME.fa上,只能比对到正链,且匹配到基因组不多于20个不同位置,允许有1个错配

1

bowtie -f -a -m 20 -v 1 --al Reads_aligned --un Reads_unaligned --norc GENOME.fa Reads.fa Reads.bwt 2> log

$ bowtie2 -u 10000 -p 8 -x index -1 reads1.fq -2 reads2.fq -S out.sam

  • 注:
  • -f 指定query文件为fasta格式
  • -a 保留所有比对结果
  • -m 指定最大比对到基因组的次数
  • -v 允许最大错配数,为[0-2]
  • --al 能map到GENOME的reads,fasta格式
  • --un 不能map到GENOME的reads,fasta格式
  • --norc 不输出匹配到负链的结果;如果不想输出比对到正链的结果,则用"--nofw"。不指定该选项则正负链结果都输出
  • 后面依次写上GENOME索引文件,Reads文件,输出结果文件Reads.bwt,日志文件log。

(三)、bowtie输出结果的说明

1

2

sample001_x75 + Chr1 12453 ATCGGCCAATTACGGACTTAA IIIIIIIIIIIIIIIIIIIII 4 9:G>T

1 2 3 4 5 6 7 8

  • 1. query id
  • 2. "+"表示正向match;"-"表示对query作反向互补后match
  • 3. reference id
  • 4. 第2列为"+"时,表示query 第一个碱基map到reference(5"->3")上的位置,0-based(以0开始);第2列为"-"时,表示query的反向互补序列第 一 个碱基map到reference(5"->3")上的位置,0-based(以0开始)
  • 5. 如果第2列为"+",则和query序列一致;否则,和query序列反向互补
  • 6. 质量文件,如果query文件为fasta格式,则无法获取质量文件,用I代替,I的数量与query序列长度一致
  • 7. 当前query能map到GENOME的4个不同位置
  • 8. 如果存在第8列,表示有mismatch。第8列可以分为三个部分,最左端的数字,中间的碱基为reference碱基,最右端的碱基为query碱基,下面分情况讨论:
  • 第2列为"+"时:最左端的数字9表示query从5"端数起,第10个碱基为"T",而对应的reference为"G";

第2列为"-"时:最左端的数字9表示query先作反向互补,然后从3"端数起,第10个碱基为"T",而对应的reference为"G";

Bowtie2使用方法与参数详细介绍

http://bowtie-bio.sourceforge.net/bowtie2/manual.shtml

http://www.hzaumycology.com/chenlianfu_blog/?p=178

http://www.plob.org/2012/09/02/4540.html

懒人必看

Bowtie2 -q --phred33 --sensitive --end-to-end -I 0 -X 500 --fr --un unpaired --al aligned \ --un-conc unconc --al-conc alconc -p 6 --reorder -x{-1-2| -U} -S []

用法:

bowtie2 [options]* -x <bt2-idx> {-1 <m1> -2 <m2> | -U <r>} -S [<hit>]

必须参数:

-x <bt2-idx> 由bowtie2-build所生成的索引文件的前缀。首先 在当前目录搜寻,然后 在环境变量 BOWTIE2_INDEXES 中制定的文件夹中搜寻。 -1 <m1> 双末端测寻对应的文件1。可以为多个文件,并用逗号分开;多个文件必须和 -2 <m2> 中制定的文件一一对应。比如:"-1 flyA_1.fq,flyB_1.fq -2 flyA_2.fq,flyB _2.fq". 测序文件中的reads的长度可以不一样。 -2 <m2> 双末端测寻对应的文件2. -U <r> 非双末端测寻对应的文件。可以为多个文件,并用逗号分开。测序文件中的reads的 长度可以不一样。 -S <hit> 所生成的SAM格式的文件前缀。默认是输入到标准输出。

以下是可选项:

输入

-q 输入的文件为FASTQ格式文件,此项为默认值。 -qseq 输入的文件为QSEQ格式文件。 -f 输入的文件为FASTA格式文件。选择此项时,表示--ignore-quals也被选择了。 -r 输入的文件中,每一行代表一条序列,没有序列名和测序质量等。选择此项时,表示-- ignore-quals也被选择了。 -c 后直接为比对的reads序列,而不是包含序列的文件名。序列间用逗号隔开。选择此项时, 表示—ignore-quals也被选择了。 -s/--skip <int> input的reads中,跳过前<int>个reads或者pairs。 -u/--qupto <int> 只比对前<int>个reads或者pairs(在跳过前<int>个reads或者 pairs后)。Default: no limit. -5/--trim5 <int> 剪掉5"端<int>长度的碱基,再用于比对。(default: 0). -3/--trim3 <int> 剪掉3"端<int>长度的碱基,再用于比对。(default: 0). --phred33 输入的碱基质量等于ASCII码值加上33. 在最近的illumina pipiline中 得以运用。 --phred64 输入的碱基质量等于ASCII码值加上64. --solexa-quals 将Solexa的碱基质量转换为Phred。在老的GA Pipeline版本中得以 运用。Default: off. --int-quals 输入文件中的碱基质量为用“ ”分隔的数值,而不是ASCII码。比如 40 40 30 40...。Default: off.

–end-to-end模式

--very-fast Same as: -D 5 -R 1 -N 0 -L 22 -i S,0,2.50 --fast Same as: -D 10 -R 2 -N 0 -L 22 -i S,0,2.50 --sensitive Same as: -D 15 -R 2 -N 0 -L 22 -i S,1,1.15 (default in --end-to-end mode) --very-sensitive Same as: -D 20 -R 3 -N 0 -L 20 -i S,1,0.50

–loca模式下的预设参数

–loca模式下的预设参数 --very-fast-local Same as: -D 5 -R 1 -N 0 -L 25 -i S,1,2.00 --fast-local Same as: -D 10 -R 2 -N 0 -L 22 -i S,1,1.75 --sensitive-local Same as: -D 15 -R 2 -N 0 -L 20 -i S,1,0.75 (default in --local mode) --very-sensitive-local Same as: -D 20 -R 3 -N 0 -L 20 -i S,1,0.50

参数:

-N <int> 进行种子比对时允许的mismatch数. 可以设为0或者1. Default: 0. -L <int> 设定种子的长度. ************************************************************ 功能选项 给bowtie的一些参数设定值的时候,使用一个计算公式代替,于是值的大小与比对序列的长 度成一定关系。<func>有三部分组成: (a)计算方法, 包括常数(C),线性(L),平方根(S)和 自然对数(G); (b)一个常数; (c)一个系数. 例如: <func> 为 L,-0.4,-0.6 则计算公式为: f(x) = -0.4 + -0.6 * x <func> 为G,1,5.4 则计算公式为: f(x) = 1.0 + 5.4 * ln(x) ************************************************************ -i <func> 设定两个相邻种子间所间距的碱基数。 ************************************************************ 例如:如果read的长度为30, 种子的长度为10, 相邻种子的间距为6,则提取出的种子如下 所示: Read: TAGCTACGCTCTACGCTATCATGCATAAAC Seed 1 fw: TAGCTACGCT Seed 1 rc: AGCGTAGCTA Seed 2 fw: CGCTCTACGC Seed 2 rc: GCGTAGAGCG Seed 3 fw: ACGCTATCAT Seed 3 rc: ATGATAGCGT Seed 4 fw: TCATGCATAA Seed 4 rc: TTATGCATGA ************************************************************ 在--end-to-end模式中默认值为”-i S,1,1.15”.即表示f(x) = 1 + 1.15 * sqrt(x). 如果read长度为100, 则相邻种子的间距为12. --n-ceil <func> 设定read中允许含有不确定碱基(非GTAC,通常为N)的最大数目. Default: L,0,0.15. 计算公式为: f(x) = 0 + 0.15 * x, 表示长度为100的read 最多运行存在15个不确定碱基. 一旦不确定碱基数超过15, 则该条read会被过滤掉. --dpad <int> Default: 15. --gbar <int> 在read头尾<int>个碱基内不允许gap. Default: 4. --ignore-quals 计算错配罚分的时候不考虑碱基质量. 当输入序列的模式为-f, -r 或 者-c的时候, 该设置自动成为默认设置. --nofw/--norc –nofw设定read不和前导链(forward reference strand)进行比对; --norc设定不和后随链(reverse-complement reference strand)进行比对. Default: both strands enabled. --end-to-end 比对是将整个read和参考序列进行比对. 该模式--ma的值为0. 该模式为 默认模式, --local模式冲突. --local 该模式下对read进行局部比对, 从而, read两端的一些碱基不比对,从而使比 对得分满足要求. 该模式下 –ma默认为2.

得分罚分参数

--ma <int> 设定匹配得分. --local模式下每个read上碱基和参考序列上碱基匹配, 则 加<int>分. 在—end-to-end模式中无效. Default: 2. --mp MX,MN 设定错配罚分. 其中MX为所罚最高分, MN为所罚最低分. 默认设置下罚分与 碱基质量相关. 罚分遵循的公式为: MN + floor( (MX-MN)(MIN(Q, 40.0)/40.0) ). 其中Q为碱基的质量值. 如果设置了—ignore-qual参数, 则错配总是罚最高分. Default: MX = 6, MN = 2. --np <int> 当匹配位点中read, reference上有不确定碱基(比如N)时所设定的罚分值. Default: 1. --rdg <int1>,<int2> 设置在read上打开gap 罚分<int1>, 延长gap罚分<int2>. Default: 5, 3. --rfg <int1>,<int2> 设置在reference上打开gap 罚分<int1>, 延长gap罚分 <int2>. Default: 5, 3. --score-min <func> 设定成为有效比对的最小分值. 在—end-to-end模式下默认值为: L,-0.6,-0.6; 在--local模式下默认值为: G,20,8.

报告参数

-k <int> 默认设置下, bowtie2搜索出了一个read不同的比对结果, 并报告其中最好的 比对结果(如果好几个最好的比对结果得分一致, 则随机挑选出其中一个). 而在该模式下, bowtie2最多搜索出一个read <int>个比对结果, 并将这些结果按得分降序报告出来. -a 和-k参数一样, 不过不限制搜索的结果数目. 并将所有的比对结果都按降序报告出来. 此参数和-k参数冲突. 值得注意的是: 如果基因组含有很多重复序列时, 该参数会导致程序 运行极其缓慢.

Effort 参数

-D <int> 比对时, 将一个种子延长后得到比对结果, 如果不产生更好的或次好的比对结果, 则该次比对失败. 当失败次数连续达到<int>次后, 则该条read比对结束. Bowtie2才会 继续进行下去. Default: 15. 当具有-k或-a参数, 则该参数所产生的限制会自动调整. -R <int> 如果一个read所生成的种子在参考序列上匹配位点过多. 当每个种子平均匹配超 过300个位置, 则通过一个不同的偏移来重新生成种子进行比对. <int>则是重新生成种子 的次数. Default: 2.

Paired-end 参数

-I/--minins <int> 设定最小的插入片段长度. Default: 0. -X/--maxins <int> 设定最长的插入片段长度. Default: 500. --fr/--rf/--ff 设定上下游reads和前导链paired-end比对的方向. --fr: 匹配时, read1在5"端上游, 和前导链一致, read2在3"下游, 和前导链反向互补. 或者read2在 上游, read1在下游反向互补; --rf: read1在5"端上游, 和前导链反向互补, read2在 3"端下游, 和前导链一致; --fr: 两条reads都和前导链一致. Default: --fr. 默认 设置适合于Illumina的paired-end测序数据; 若是mate-paired, 则要选择—rf参数. --no-mixed 默认设置下, 一对reads不能成对比对到参考序列上, 则单独对每个read进 行比对. 该选项则阻止此行为. --no-discordant 默认设置下, 一对reads不能和谐比对(concordant alignment, 即满足-I, -X, --fr/--rf/--ff的条件)到参考序列上, 则搜寻其不和谐比对(discon cordant alignment, 即两条reads都能独一无二地比对到参考序列上, 但是不满足-I, -X,--fr/--rf/--ff的条件). 该选项阻止此行为. --dovetail read1和read2的关系为dovetail的时候,该状况算为和谐比对. 默认情况 下dovetail不算和谐比对. --no-contain read1和read2的关系为包含的时候, 该状况不算为和谐比对. 默认情况 下包含关系算为和谐比对. --no-overlap read1和read2的关系为有重叠的时候, 该状况不算为和谐比对. 默认情 况下两个reads重叠算为和谐比对.

输出参数

-t/--time --un <path> 将unpaired reads写入到<path>. --un-gz <path> 将unpaired reads写入到<path>, gzip压缩. --un-bz2 <path> 将unpaired reads写入到<path>, bz2压缩. --al <path> 将至少能比对1次以上的unpaired reads写入<path>. --al-gz <path> ... ,gzip压缩. --al-bz2 <path> ... ,bz2压缩. --un-conc <path> 将不能和谐比对的paired-end reads写入<path>. --un-conc-gz <path> ... ,gzip压缩. --un-conc-bz2 <path> ... ,bz2压缩. --al-conc <path> 将至少能和谐比对一次以上的paired-end reads写入<path>. --al-conc-gz <path> ... ,gzip压缩. --al-conc-bz2 <path>... ,bz2压缩. --quiet 安静模式,除了比对错误和一些严重的错误, 不在屏幕上输出任何东西. --met-file <path> 将bowtie2的检测信息(metrics)写入文件<path>. 用于debug. Default: metrics disabled. --met-stderr <path> 将bowtie2的检测信息(metrics)写入标准错误文件句柄. 和上 一个选项不冲突. Default: metrics disabled. --met <int> 每隔<int>秒写入一次metrics记录. Default: 1.

Sam 参数

--no-unal 不记录没比对上的reads. --no-hd 不记录SAM header lines (以@开头). --no-sq 不记录@SQ的SAM header lines. --rg-id <text> 设定read group Id到<text>. --rg <text> 增加<text>作为一行@RG.

性能参数

-o/--offrate <int> 无视index的offrate值, 以<int>取代之. Index默认的<int> 值为5. <int>值必须大于index的offrate值, 同时<int>越大, 耗时越长,耗内存越少. -p/--threads NTHREADS 设置线程数. Default: 1 --reorder 多线程运算时, 比对结果在顺序上会和文件中reads的顺序不一致, 使用该选 项, 则使其一致. --mm 使用内存定位的I/O来载入index, 而不是常规的文件I/O. 从而使多个bowtie程 序共用内存中同样的index, 节约内存消耗.

其它参数:

--qc-filter 滤除QSEQ fileter filed为非0的reads. 仅当有—qseq选项时有效.

Default: off.

--seed <int> 使用<int>作为随机数产生的种子. Default: 0.

--version 打印程序版本并退出

-h/--help 打印用法信息并推出

Bowtie2使用指南

常用步骤

  1. 对参考序列构建index

$ bowtie2-build genome.fasta index

  1. 尝试使用前10000个reads进行比对

$ bowtie2 -u 10000 -p 8 -x index -1 reads1.fq -2 reads2.fq -S out.sam

使用8个线程进行比对

$ bowtie2 -p 8 -x index -1 reads1.fq -2 reads2.fq -S out.sam

比对的sam结果中添加了read group信息

$ bowtie2 -p 8 --rg-id sample01 --rg "PL:ILLUMINA" --rg "SM:sample01" -x index -1 reads1.fq -2 reads2.fq -S out.sam

常用的参数进行比对,可以更改其中的参数获得更好的结果

$ bowtie2 -q --phred33 --sensitive --end-to-end -I 0 -X 500 --fr --un unpaired --al aligned --un-conc unconc --al-conc alconc -p 6 --reorder -x <bt2-idx> {-1 <m1gt; -2 <m2> | -U <r>} -S [<hit>]

用法:

bowtie2 [options]* -x <bt2-idx> {-1 <m1> -2 <m2> | -U <r>} -S [<hit>]

bowtie2-build用法

bowtie2-build默认情况下将fasta文件换成index的数据库。

$ bowtie2-build <fasta文件> <要生存的索引文件前缀名>

必须参数:

-x <bt2-idx>由bowtie2-build所生成的索引文件的前缀。首先 在当前目录搜寻,然后

在环境变量BOWTIE2_INDEXES中制定的文件夹中搜寻。

-1 <m1>双末端测寻对应的文件1。可以为多个文件,并用逗号分开;多个文件必须和-2

<m2>中制定的文件一一对应。比如:"-1 flyA_1.fq,flyB_1.fq -2 flyA_2.fq,flyB

_2.fq".测序文件中的reads的长度可以不一样。

-2 <m2>双末端测寻对应的文件2.

-U <r>非双末端测寻对应的文件。可以为多个文件,并用逗号分开。测序文件中的reads的

长度可以不一样。

-S <hit>所生成的SAM格式的文件前缀。默认是输入到标准输出。

以下是可选参数:

输入参数

-q输入的文件为FASTQ格式文件,此项为默认值。

-qseq输入的文件为QSEQ格式文件。

-f输入的文件为FASTA格式文件。选择此项时,表示--ignore-quals也被选择了。

-r输入的文件中,每一行代表一条序列,没有序列名和测序质量等。选择此项时,表示--

ignore-quals也被选择了。

-c后直接为比对的reads序列,而不是包含序列的文件名。序列间用逗号隔开。选择此项时,

表示—ignore-quals也被选择了。

-s/--skip <int> input的reads中,跳过前<int>个reads或者pairs。

-u/--qupto <int>只比对前<int>个reads或者pairs(在跳过前<int>个reads或者

pairs后)。Default: no limit.

-5/--trim5 <int>剪掉5'端<int>长度的碱基,再用于比对。(default: 0).

-3/--trim3 <int>剪掉3'端<int>长度的碱基,再用于比对。(default: 0).

--phred33输入的碱基质量等于ASCII码值加上33.在最近的illumina pipiline中

得以运用。最低碱基质量是“#”。

--phred64输入的碱基质量等于ASCII码值加上64.最低碱基质量是“B”。

--solexa-quals将Solexa的碱基质量转换为Phred。在老的GA Pipeline版本中得以

运用。Default: off.

--int-quals输入文件中的碱基质量为用“ ”分隔的数值,而不是ASCII码。比如40 40

30 40...。Default: off.

–end-to-end模式下的预设参数

--very-fast Same as: -D 5 -R 1 -N 0 -L 22 -i S,0,2.50

--fast Same as: -D 10 -R 2 -N 0 -L 22 -i S,0,2.50

--sensitive Same as: -D 15 -R 2 -N 0 -L 22 -i S,1,1.15 (default in --end-to-endmode)

--very-sensitive Same as: -D 20 -R 3 -N 0 -L 20 -i S,1,0.50

–loca模式下的预设参数

--very-fast-local Same as: -D 5 -R 1 -N 0 -L 25 -i S,1,2.00

--fast-local Same as: -D 10 -R 2 -N 0 -L 22 -i S,1,1.75

--sensitive-local Same as: -D 15 -R 2 -N 0 -L 20 -i S,1,0.75 (default in --local mode)

--very-sensitive-local Same as: -D 20 -R 3 -N 0 -L 20 -i S,1,0.50

比对参数:

-N <int>进行种子比对时允许的mismatch数.可以设为0或者1.Default: 0.

-L <int>设定种子的长度.

************************************************************

功能选项

给bowtie的一些参数设定值的时候,使用一个计算公式代替,于是值的大小与比对序列的长

度成一定关系。<func>有三部分组成: (a)计算方法,包括常数(C),线性(L),平方根(S)和

自然对数(G); (b)一个常数; (c)一个系数.

例如:

<func>为L,-0.4,-0.6则计算公式为: f(x) = -0.4 + -0.6 * x

<func> 为G,1,5.4则计算公式为: f(x) = 1.0 + 5.4 * ln(x)

************************************************************

-i <func>设定两个相邻种子间所间距的碱基数。

************************************************************

例如:如果read的长度为30,种子的长度为10,相邻种子的间距为6,则提取出的种子如下

所示:

Read: TAGCTACGCTCTACGCTATCATGCATAAAC

Seed 1 fw: TAGCTACGCT

Seed 1 rc: AGCGTAGCTA

Seed 2 fw: CGCTCTACGC

Seed 2 rc: GCGTAGAGCG

Seed 3 fw: ACGCTATCAT

Seed 3 rc: ATGATAGCGT

Seed 4 fw: TCATGCATAA

Seed 4 rc: TTATGCATGA

************************************************************

在--end-to-end模式中默认值为”-i S,1,1.15”.即表示f(x) = 1 + 1.15 *

sqrt(x).如果read长度为100,则相邻种子的间距为12.

--n-ceil <func>设定read中允许含有不确定碱基(非GTFAC,通常为N)的最大数目.

Default: L,0,0.15.计算公式为: f(x) = 0 + 0.15 * x,表示长度为100的read

最多运行存在15个不确定碱基.一旦不确定碱基数超过15,则该条read会被过滤掉.

--dpad <int> Default: 15.

--gbar <int>在read头尾<int>个碱基内不允许gap. Default: 4.

--ignore-quals计算错配罚分的时候不考虑碱基质量.当输入序列的模式为-f, -r或

者-c的时候,该设置自动成为默认设置.

--nofw/--norc –nofw设定read不和前导链(forward reference strand)进行比对;

--norc设定不和后随链(reverse-complement reference strand)进行比对.

Default: both strands enabled.

--end-to-end比对是将整个read和参考序列进行比对.该模式--ma的值为0.该模式为

默认模式, --local模式冲突.

--local该模式下对read进行局部比对,从而, read两端的一些碱基不比对,从而使比

对得分满足要求.该模式下 –ma默认为2.

得分罚分参数

--ma <int>设定匹配得分. --local模式下每个read上碱基和参考序列上碱基匹配,则

加<int>分.在—end-to-end模式中无效. Default: 2.

--mp MX,MN设定错配罚分.其中MX为所罚最高分, MN为所罚最低分.默认设置下罚分与

碱基质量相关.罚分遵循的公式为: MN + floor( (MX-MN)(MIN(Q, 40.0)/40.0) ).

其中Q为碱基的质量值.如果设置了—ignore-qual参数,则错配总是罚最高分. Default:

MX = 6, MN = 2.

--np <int>当匹配位点中read, reference上有不确定碱基(比如N)时所设定的罚分值.

Default: 1.

--rdg <int1>,<int2>设置在read上打开gap罚分<int1>,延长gap罚分<int2>.

Default: 5, 3.

--rfg <int1>,<int2>设置在reference上打开gap罚分<int1>,延长gap罚分

<int2>. Default: 5, 3.

--score-min <func>设定成为有效比对的最小分值.在—end-to-end模式下默认值为:

L,-0.6,-0.6;在--local模式下默认值为: G,20,8.

报告参数

-k <int>默认设置下, bowtie2搜索出了一个read不同的比对结果,并报告其中最好的

比对结果(如果好几个最好的比对结果得分一致,则随机挑选出其中一个).而在该模式下,

bowtie2最多搜索出一个read <int>个比对结果,并将这些结果按得分降序报告出来.

-a和-k参数一样,不过不限制搜索的结果数目.并将所有的比对结果都按降序报告出来.

此参数和-k参数冲突.值得注意的是:如果基因组含有很多重复序列时,该参数会导致程序

运行极其缓慢.

Effort 参数

-D <int>比对时,将一个种子延长后得到比对结果,如果不产生更好的或次好的比对结果,

则该次比对失败.当失败次数连续达到<int>次后,则该条read比对结束. Bowtie2才会

继续进行下去. Default: 15.当具有-k或-a参数,则该参数所产生的限制会自动调整.

-R <int>如果一个read所生成的种子在参考序列上匹配位点过多.当每个种子平均匹配超

过300个位置,则通过一个不同的偏移来重新生成种子进行比对. <int>则是重新生成种子

的次数. Default: 2.

Paired-end 参数

-I/--minins <int> 设定最小的插入片段长度. Default: 0.

-X/--maxins <int> 设定最长的插入片段长度. Default: 500.

--fr/--rf/--ff 设定上下游reads和前导链paired-end比对的方向. --fr: 匹配时,

read1在5'端上游, 和前导链一致, read2在3'下游, 和前导链反向互补. 或者read2在

上游, read1在下游反向互补; --rf: read1在5'端上游, 和前导链反向互补, read2在

3'端下游, 和前导链一致; --ff: 两条reads都和前导链一致. Default: --fr. 默认

设置适合于Illumina的paired-end测序数据; 若是mate-paired, 则要选择—rf参数.

--no-mixed 默认设置下, 一对reads不能成对比对到参考序列上, 则单独对每个read进

行比对. 该选项则阻止此行为.

--no-discordant 默认设置下, 一对reads不能和谐比对(concordant alignment,

即满足-I, -X, --fr/--rf/--ff的条件)到参考序列上, 则搜寻其不和谐比对(discon

cordant alignment, 即两条reads都能独一无二地比对到参考序列上, 但是不满足-I,

-X,--fr/--rf/--ff的条件). 该选项阻止此行为.

--dovetail read1和read2的关系为dovetail的时候,该状况算为和谐比对. 默认情况

下dovetail不算和谐比对.

--no-contain read1和read2的关系为包含的时候, 该状况不算为和谐比对. 默认情况

下包含关系算为和谐比对.

--no-overlap read1和read2的关系为有重叠的时候, 该状况不算为和谐比对. 默认情

况下两个reads重叠算为和谐比对.

输出参数

-t/--time --un <path> 将unpaired reads写入到<path>.

--un-gz <path> 将unpaired reads写入到<path>, gzip压缩.

--un-bz2 <path> 将unpaired reads写入到<path>, bz2压缩.

--al <path> 将至少能比对1次以上的unpaired reads写入<path>.

--al-gz <path> ... ,gzip压缩.

--al-bz2 <path> ... ,bz2压缩.

--un-conc <path> 将不能和谐比对的paired-end reads写入<path>.

--un-conc-gz <path> ... ,gzip压缩.

--un-conc-bz2 <path> ... ,bz2压缩.

--al-conc <path> 将至少能和谐比对一次以上的paired-end reads写入<path>.

--al-conc-gz <path> ... ,gzip压缩.

--al-conc-bz2 <path>... ,bz2压缩.

--quiet 安静模式,除了比对错误和一些严重的错误, 不在屏幕上输出任何东西.

--met-file <path> 将bowtie2的检测信息(metrics)写入文件<path>. 用于debug.

Default: metrics disabled.

--met-stderr <path> 将bowtie2的检测信息(metrics)写入标准错误文件句柄. 和上

一个选项不冲突. Default: metrics disabled.

--met <int> 每隔<int>秒写入一次metrics记录. Default: 1.

Sam 参数

--no-unal不记录没比对上的reads.

--no-hd不记录SAM header lines (以@开头).

--no-sq不记录@SQ的SAM header lines.

--rg-id <text>设定read group ID为text。在SAM文件的头中增加一行@RG,在输出的SAM

文件中添加Tag "RG:Z:text"。

--rg <text>使用text作为@RG的一列,比如"SM:Pool1"。在@RG中加入多列,则多次使用

该参数即可。在进行Variant calling的过程中需要@RG头,SM信息和Tag RG。

性能参数

-o/--offrate <int> 无视index的offrate值, 以<int>取代之. Index默认的<int>

值为5. <int>值必须大于index的offrate值, 同时<int>越大, 耗时越长,耗内存越少.

-p/--threads NTHREADS 设置线程数. Default: 1

--reorder 多线程运算时, 比对结果在顺序上会和文件中reads的顺序不一致, 使用该选

项, 则使其一致.

--mm 使用内存定位的I/O来载入index, 而不是常规的文件I/O. 从而使多个bowtie程

序共用内存中同样的index, 节约内存消耗.

其它参数:

--qc-filter 滤除QSEQ fileter filed为非0的reads. 仅当有—qseq选项时有效.

Default: off.

--seed <int>使用<int>作为随机数产生的种子. Default: 0.

--version打印程序版本并退出

-h/--help 打印用法信息并推出

bowtie2是一款速度快,占用内存少的一款将段序列片段定位到参照序列上的一款软件。以下的内容基本上是对bowtie2的英文manual的翻译,并截取了其中比较中药的部分。

bowtie2与bowtie1的区别:

bowtie2并不是bowtie1的简单的升级版,与bowtie1有较大的区别:

· 由于近年来HTS的发展,产出的序列长度比之前更长一些。bowtie2在处理长度大于50bp的数据是速度更快;

· bowtie2支持gap形式的比对;

· bowtie2有local alignment和end-to-end alignment两个模式的比对,而bowtie1只有end-to-end 模式的比对;

· bowtie2对输入的序列的长度没有要求,而bowtie1最长要求为1000nt;

获得并安装bowtie2:

从download中下载,选择适合自己操作系统的版本。现有的bowtie2基本上支持所有的操作系统。在linux等操作系统下安装之后可以添加环境变量,一边方便之后的使用。

bowtie2对匹配质量的评价:

1.Socres:越高则相似度越大

score值是根据序列和参考序列的相似程度决定的,决定该值大小的因素有很多:选择的匹配的模式,序列的长度,错配以及是否存在gap,罚分设置等。罚分的设置可以通过--ma (匹配奖励), --mp (错配罚分), --np (出现N时罚分), --rdg (序列中出现gap的罚分) and --rfg (参照序列中出现gap的罚分) 。

· end-to-end alignment 模式:

出现一个错配默认罚6分,出现一个2nt的gap默认罚11分(默认gap open罚 5分,gap延长1nt罚3分)。

· local aligment模式:

出现一个错配默认罚6分,出现一个2nt的gap默认罚11分(默认gap open罚 5分,gap延长1nt罚3分)。假如配对成功的话,奖励默认2分。

2.Mapping qulity:越高则位点越独特

bowtie2有时并不能完全确定一个短的序列来自与参考序列的那个位置,特别是对于那些比较简单的序列。但是bowtie2会给出一个值来显示出这个段序列来自某个位点的概率值,这个值就是mapping qulity。Mapping qulity的计算方法是:Q=-10log10p,Q是一个非负值,p是这个序列不来自这个位点的估计值。

假如说一条序列在某个参考序列上找到了两个位点,但是其中一个位点的Q明显大于另一个位点的Q值,这条序列来源于前一个位点的可能性就比较大。Q值的差距越大,这独特性越高。

Q值的计算方法来自与SAM标准格式,请查看SAM总结。

3.Aliging pair

bowtie2支持paired-end或者是mate-pair形式的测序结果,具体在运行参数中设置 –1 指定第一个文件,-2 指定第二个文件就行。

报告

报告的模式决定了bowtie2寻找多少个位置的匹配,并报告这些匹配。Bowtie2有三种不同的报告模式。

1. 默认模式:寻找多个位点,并报告出最好的那一个

在默认情况下,bowtie2寻找明确的,序列的位点。当它找到一个确认的结果的时候,它将会继续寻找比这个结果更好的结果,或者至少一样好的结果。当bowtie2已经超过了寻找努力(通过-D 和 -R设置),或者它已经找到了所有的结果的时候,bowtie2将停止查找。

2. -k模式

在-k模式下,bowtie2最多寻找N(N为-k 的参数)个独特的位点。例如N为2,那么bowtie2就寻找最多2个位点,就结束查找,并报告出两个位点的结果。

3. -a模式

在-a模式下,bowtie2将会将查找所有可能的匹配位点,并报告出所有的结果。这样就可能会导致程序运行速度很慢。

命令行

1. 函数参数

在bowtie2中有些参数是通过函数来指定的,而非使用特定的数值或者设置来指定。这些函数需要用户指定三个参数来确定该函数:(1)函数类型 F,(2)一个常数 B,(3)一个系数 A。函数类型包括:常数型 C,线性型 L,平方根型 S,和自然对数型 G。这些参数被指定为 F,B,A,分别为函数类型,常量,和系数。例如:L,-0.4,-0.6就定义了一个函数:

f(x)=-0.4+-0.6x

2. 使用

bowtie2运行的参数一般为:

bowtie2 [options]* -x <bt2-idx> {-1 <m1> -2 <m2> | -U <r>} -S [<hit>]

3. 主要参数

· -x<bt2-idx> 为已经建立好的参考序列的索引文件,通过bowtie2_build建立;

· -1 <m1> 由逗号分隔的包含mate 1的多个文件,也可以输入一个文件;

· -2 <m2> 由逗号分隔的包含mate 2的多个文件,也可以只有一个文件;

· -U <r> 由逗号分隔的包含不是pair reads的文件;

· -S <hit> 结果输出的位置。

4. 设置

输入设置:

-q 输入的序列为FASTQ格式

–qseq 输入的序列为QSEQ格式

-f 输入的序列为FASTA格式

-r 输入的序列为每行一条序列,并不包含其他的信息

-c 从命令行输入序列

-s <int> 忽略输入序列中的前<int>条序列

-u <int> 只保留输入序列中的前<int>条序列

-5 <int> 剪切掉5`端的<int>个碱基

-3 <int> 剪切掉3`端的<int>个碱基

end-to-end模式下的预设值的选项

–very-fast 即: -D 5 -R 1 -N 0 -L 22 -i S,0,2.50

–fast 即: -D 10 -R 2 -N 0 -L 22 -i S,0,2.50

–sensitive 即: -D 15 -R 2 -L 22 -i S,1,1.1(为end-to-end模式下的默认设置)

--very-sensitive 即: -D 20 -R 3 -N 0 -L 20 -i S,1,0.50

local 模式下的预设选项

–very-fast-local 即: -D 5 -R 1 -N 0 -L 25 -i S,1,2.00

–fast-local 即: -D 10 -R 2 -N 0 -L 22 -i S,1,1.75

–sensitive-local 即: -D 15 -R 2 -N 0 -L 20 -i S,1,0.75 (为 --local 模式下的默认设置)

–very-sensitive-local 即: -D 20 -R 3 -N 0 -L 20 -i S,1,0.50

比对的设置

-N <int> 设置在种子配对是允许的错配数

-L <int> 设置种子的长度

--local bowite2进行local模式的匹配

--end-to-end bowtie2进行end-to-end模式匹配

还有几个参数这里就不介绍了

打分设置

--ma <int> 配对成功奖励

--mp MX,MN 设置当出现错配的时候最大的(MX)和最小的(MN)罚分,被扣的分数为MN+floor((mx-mn)(MIN)(Q,40)/40),Q为phred quality 值。

--np <int> 设置出现N的时候的罚分

--rfg <int1>,<int2> 设置gap open(<int1>)和延伸(<int2>)的罚分

报告设置

-k <int> 之前已经介绍过

-a 之前已经介绍过

输出设置

输出设置比较简单就不写了

性能设置

-p 设置线程数目

SAM 输出

一下是一个简单的对bowtie2输出的SAM格式文件的说明,更加详细的SAM格式文件说明请参考SAM format specification。

在SAM输出的结果中每一行都包括十二项通过Tab分隔,从左到右分别是:

1 序列的名字;

2 概括出一个合适的标记,各个数字分别代表

· 1 序列是一对序列中的一个

· 2 比对结果是一个pair-end比对的末端

· 4 没有找到位点

· 8 这个序列是pair中的一个但是没有找到位点

· 16 在这个比对上的位点,序列与参考序列反向互补

· 32 这个序列在pair-end中的的mate序列与参考序列反响互补

· 64 序列是 mate 1

· 128 序列是 mate 2

假如说标记为以上列举出的数目,就可以直接推断出匹配的情况。假如说标记不是以上列举出的数字,比如说83=(64+16+2+1),就是这几种情况值和。

3 参考序列的名字

4 在参考序列上的位置

5 mapping qulity

6 代表比对结果的CIGAR字符串,如37M1D2M1I,这段字符的意思是37个匹配,1个参考序列上的删除,2个匹配,1个参考序列上的插入。M代表的是alignment match(可以是错配)

7 mate 序列所在参考序列的名称

8 mate 序列在参考序列上的位置

9 估计出的片段的长度,当mate 序列位于本序列上游时该值为负值。

10 read的序列

11 ASCII码格式的序列质量

12 可选的区域

· AS:i 匹配的得分

· XS:i 第二好的匹配的得分

· YS:i mate 序列匹配的得分

· XN:i 在参考序列上模糊碱基的个数

· XM:i 错配的个数

· XO:i gap open的个数

· XG:i gap 延伸的个数

· NM:i 经过编辑的序列

· YF:i 说明为什么这个序列被过滤的字符串

· YT:Z

· MD:Z 代表序列和参考序列错配的字符串

bowtie2-build

在运行bowtie2之前,需要先建立参考序列的索引文件。bowtie2-build 输出六个文件分别是 .1.bt2,.2.bt2, .3.bt2, .4.bt2, .rev.1.bt2, 和.rev.2.bt2。

bowtie2-build运行的命令为:

bowtie2-build [options]* <reference_in> <bt2_base>

<reference_in> 参考序列的文件

<bt2_base> 要输出的索引文件的名称,运行bowtie2的时候使用这个名称

选项

-f 输入参考序列为fasta格式

-c 通过命令行输入参考序列

其余参数不再介绍