go.class的制作

方法:

1. 使用OBO-Edit来搜索GO二级水平term包含的terms,并生成相应的obo文件;
2. 使用编写的perl程序来提取obo文件中的内容。

过程:

1. 下载gene_ontology.1_2.obo文件。

下载网页:http://www.geneontology.org/GO.downloads.ontology.shtml

2. 安装OBO-Edit并运行。

在官网http://www.oboedit.org/上下载OBO-Edit并安装。
文章写的时间对应的版本是oboedit_2_3-b3_unix_install4j.sh (45.0 MB)
# sh oboedit_2_3-b3_unix_install4j.sh
按提示进行安装。
$ rm oboedit_config -rf
$ OBO-Edit
运行该程序需要java支持。

3. 使用gene_ontology.1_2.obo文件来search,得到二级水平的GO term的子terms。

1. File菜单下Load Ontologies,载入gene_ontology.1_2.obo文件;
2. File菜单下Save as,进入Write ontology画面;
3. 点击Advanced——点击Add——输入保存路径(biological adhesion为例)——打勾Filter terms——点击+号——Find terms that have a ID that equals the value GO:0022610 in Ancestor that can be reached via any tyep——点击OK
4. 等待一段时间,则会生成文件名为“biological ahension”的obo文件

4. 对每一个GO二级水平的term都进行这样的search和save,得到很多obo文件。使用perl程序make_go_class.pl来生成go.class文件。

$ perl make_go_class.pl *    #其中*代表生成的所有的obo文件
该程序以‘>>'的方式写入到文件句柄中,注意不要重复命令。

5. 还有GO二级terms本身没有加入到go.class中去。按以上的search和save方式进行其obo文件的制作,不同的是将Ancestor换成Parent。得到的obo文件使用make_go_class.pl来成才go.class文件。该go.class文件是错误的,需要手动修改后再整合到总的go.class文件中去。
6. 最后生成了go.class文件为最新的文件了。但是使用http://www.hzaumycology.com/chenlianfu_blog/?p=586中的方法做GO分类图的时候,依旧提示error的GO号,这些GO号都是废弃了的GO号。
7. go.class文件下载(2013.2.27)

Perl模块的使用

1. 随机打乱数组

use List::Util;
@array = List::Util::shuffle @array;

2. 对数组取平均值,标准差等

use Statistics::Descriptive;
my $stat = Statistics::Descriptive;
$stat->add_data(@insert_size);
my $mean = $stat->mean();
my $sd = $stat->standard_deviation();
my $sum = $stat->sum();
my $min = $stat->min();
my $max = $stat->max();

3.

GO富集分析

一、程序与文件准备
1. 下载并安装libgd,gd图形库。网址: https://bitbucket.org/libgd/gd-libgd。

# gdlib-config

使用上述命令能查看libgd的一些信息,比如版本和库文件和头文件所在的路径。
2. 在cpan中下载GD并安装。网址:http://search.cpan.org/。安装过程中输入命令如下:

# cd GD-2.XX
# perl Makefile.PL
# make
# make test (optional)
# make html (optional)
# make install

上面第2步中该perl程序会通过gdlib-config程序来检查libgd的版本,看是否满足安装条件。出了问题必须检查$PATH中是否有/usr/local/bin。
3. 安装GO-TermFinder。网址:http://search.cpan.org/。
4. 下载最新的Ontology文件gene_ontology_edit.obo或gene_ontology.obo。网址:http://obo.cvs.sourceforge.net/viewvc/obo/obo/ontology/genomic-proteomic/。

COG分类

本文讲述如何对基因组预测基因进行COG分类。

一 COG简介

COG,即Clusters of Orthologous Groups of proteins。

其网址主页为: http://www.ncbi.nlm.nih.gov/COG/。

网页版使用工具网址: http://www.ncbi.nlm.nih.gov/COG/old/xognitor.html。

使用说明文档网址: http://www.nlm.nih.gov/COG/old/COGhelp.html。

其FTP站点为: ftp://ftp.ncbi.nih.gov/pub/COG/。

通过观看其主页和说明文档,可以理解为COG是NCBI的数据库。COG的中文释义即“同源蛋白簇”。COG分为两类,一类是原核生物的,另一类是真核生物。原核生物的一般称为COG数据库;真核生物的一般称为KOG数据库。

COG注释作用:1. 通过已知蛋白对未知序列进行功能注释; 2. 通过查看指定的COG编号对应的protein数目,存在及缺失,从而能推导特定的代谢途径是否存在; 3. 每个COG编号是一类蛋白,将query序列和比对上的COG编号的proteins进行多序列比对,能确定保守位点,分析其进化关系。

二 将序列进行COG分类的步骤

1. COG的ftp里边提供了一个名为myva的文件,该文件里面为COG数据库的蛋白质序列,有192987条。将该序列文件使用使用ncbi-blast-2.2.26+中的blastdb程序制作出一个前缀为cog的蛋白质数据库。

2. 将需要进行COG注释并分类的DNA序列或protein序列分别使用blastx或blastp比对到上一步骤建好的cog数据库中。得出xml的比对结果。

3. 根据上一步骤的比对结果,得出与query序列相似的cog蛋白id。在COG的ftp里面有一个名为whog的文件,该文件中记录着COG数据中绝大部分的蛋白质id以及其所对应的以COG开头的protein编号,同时也记录这COG编号对应的功能分类编号。因此,可以得出query序列注释的COG编号以及其功能分类编号。

4. ftp中还有一个名为fun.txt的文件,该文件记录这COG的功能分类编号,及其对编号的功能描述。

5. 因此,我编写了一个脚本程序用于COG分类注释。该脚本程序名字为cog.pl。输入为所需要进行COG分类的fasta序列文件,得出序列的比对结果和分类统计。

三 将序列进行KOG分类

方法和COG的分类一致。值得注意的是KOG数据中protein的编号是以LSE或TWOG开头的。同样,我编写了一个kog.pl。

四 说明补充

值得说明的是,KOG数据库中蛋白质序列数目为112920条,但是其中有protein编号的只有27887条,占25%。而COG数据库中蛋白质序列条数为192987,其中有COG编号的有129326条,占67%。所以比对结果中,很多序列比对上了KOG数据库,但是没有protein编号;而在比对到COG数据库时会好很多。

所以,可以先将序列比对到COG数据库,得出分类数据;然后将没有分类编号的序列挑选出来,再比对到KOG数据库,得出分类数据;然后将两个分类数据进行整合,然后画出COG分类图。

安装blast2go数据库,所需知道的MySQL的编译和安装和简易使用

1 修改配置文件

# cp /usr/share/mysql/my-huge.cnf /etc/my.cnf
# vim /etc/my.cnf
在 [mysqld] 下加入一行:
datadir=/var/lib/mysql

2 创建系统自带的数据库和表

# mysql_install_db --user=mysql

3 启动mysql服务并将mysql的root用户密码修改为123

# /etc/init.d/mysqld restart
# mysqladmin -u root password '123'

4 删除测试数据库和匿名用户

# mysql_secure_installation

5 使用数据库,并创建b2g数据库

# mysql -u root -h 127.0.0.1 -p
mysql> SHOW DATABASES;
mysql> CREATE DATABASE b2g;
mysql> SHOW DATABASES;
mysql> quit;

primer3设计引物详解

1. Primer3 简介

参考自:http://primer3.wi.mit.edu/primer3web_help.htm。
Primer3的主页:http://primer3.sourceforge.net/。主页上能看到对于Primer3的简要介绍和下载通道。
Primer3的网页版:http://primer3.wi.mit.edu/
Primer3的源代码:http://sourceforge.net/projects/primer3/.
Primer3文献:Steve Rozen and Helen J. Skaletsky (2000) Primer3 on the WWW for general users and for biologist programmers. In: Krawetz S, Misener S (eds) Bioinformatics Methods and Protocols: Methods in Molecular Biology. Humana Press, Totowa, NJ, pp 365-386

2. 安装 Primer3

$ wget http://downloads.sourceforge.net/project/primer3/primer3/2.3.6/primer3-src-2.3.6.tar.gz
$ tar zxf primer3-src-2.3.6.tar.gz -C /opt/biosoft
$ cd primer3-2.3.6/src/
$ make all
$ make test   此项用于检测 primer3 是否正确安装,确保输出中不能有 FAILED。
$ echo 'PATH=$PATH:/opt/biosoft/primer3-2.3.6/src/' >> ~/.bashrc
$ source ~/.bashrc

3. 使用 Primer3

Primer3 的主程序是 primer3_core .

3.1 primer3_core 的常用例子和参数

常用例子:

$ perl -n -e 's/\s*#.*//;print' /opt/biosoft/primer3-2.3.6/my_default_settings.txt > p3_settings_file 
$ primer3_core -p3_settings_file p3_settings_file -strict_tags input_file > result.p3.out
如果没有input file,则 primer3_core 从标准输入读取数据。

常用参数:

-about
显示primer3版本号并退出

-default_version=n
默认设置为 n=2,让primer3使用最新的默认设置;n=1 则让程序使用 2.23 或之前版本的设置

-format_output
让 primer3_core 产生人类易读的结果,否则产生机器易读的结果(Boulder-IO 格式结果)。

-strict_tags
要求 input file 中的标签要全部正确。设置该参数后,如果有标签不能被程序识别,则报错并停止执行。不设置此参数,则软件忽略不识别的参数。

-p3_settings_file=file_path
指定 primer_core 的配置文件,该配置文件的设定会取代默认设置。当然,input file 的设置也能取代这个文件的设置。此文件的格式为:第 1 行固定为 "Primer3 File - http://primer3.sourceforge.net"; 第 2 行固定为 "P3_FILE_TYPE=settings"; 第 3 行是个空行;从第 4 行开始则是标准的 Boulder-IO 格式内容。

-echo_settings_file
打印出 p3_settings_file 中的设置信息。如果没有指定设置文件,或含有 -format_output,则该参数失效。

-io_version=n 
此参数现在仅能设置为默认值 4 。表示此参数算是取消了。

-output=file_path
指定输出文件路径,如果不指定,则输出到标准输出。

-error=file_path
指定错误信息输出路径,如果不指定,则输出到stderr中。

3.2 Primer3 的输入和输出

默认设置下,primer3 的输入文件为 Boulder-IO 格式的,该格式有利于程序间的相互交流。同时输出文件的格式也为 Boulder-IO。该文件格式是文本形式的,每一个引物设计的记录的结尾以“=\n”进行分隔。每个记录由多个标签和对应的值构成。例如:

SEQUENCE_ID=example
SEQUENCE_TEMPLATE=GTAGTCAGTAGACNATGACNACTGACGATGCAGACNACACACACACACACAGCACACAGGTATTAGTGGGCCATTCGATCCCGACCCAAATCGATAGCTACGATGACG
SEQUENCE_TARGET=37,21
PRIMER_TASK=generic
PRIMER_PICK_LEFT_PRIMER=1
PRIMER_PICK_INTERNAL_OLIGO=1
PRIMER_PICK_RIGHT_PRIMER=1
PRIMER_OPT_SIZE=18
PRIMER_MIN_SIZE=15
PRIMER_MAX_SIZE=21
PRIMER_MAX_NS_ACCEPTED=1
PRIMER_PRODUCT_SIZE_RANGE=75-100
P3_FILE_FLAG=1
SEQUENCE_INTERNAL_EXCLUDED_REGION=37,21
PRIMER_EXPLAIN_FLAG=1
=

3.3 常用 primer3 标签

primer3发布2.0后,许多的Boulder-IO的标签有了变动,而且新增了许多标签。标签分三类:分别以SEQUENCE_, PRIMER_ 和 P3_开头。关于标签的变动,具体参考本文开头的参考网址。

SEQUENCE 标签 (以 SEQUENCE_ 开头):

SEQUENCE_ID (string; default empty)
序列的名称,为序列的简要描述。

SEQUENCE_TEMPLATE (nucleotide sequence; default empty)
进行引物设计的模板琏序列,方向必须是从5'到3'。碱基序列最好用大写字母,并且序列中的所有碱基都在一行上书写。如果用小写字母,则能被标签 PRIMER_LOWERCASE_MASKING 进行屏蔽。

SEQUENCE_INCLUDED_REGION (interval list; default empty)
该标签的值为:<start>,<length>。即在从序列的第<start>位点的碱基开始,到往后<length>长度的子序列中设计引物。

SEQUENCE_TARGET (interval list; default empty)
该标签的值为:<start>,<length>。即从序列的第<start>位点的碱基开始,到往后<length>长度的区域是要扩增的目的片段,设计的引物在此区域的侧翼才能去得较好的分值。可以在多行上设置多次该标签,从而设定多个目标区间。

SEQUENCE_EXCLUDED_REGION (interval list; default empty)
该标签的值为:<start><length>。设计引物时要避开此段子序列区域。比如此区域含有诸如ALUs或者LINEs等重复序列时。

SEQUENCE_PRIMER_PAIR_OK_REGION_LIST (semicolon separated list of integer "quadruples"; default empty)
该标签的值为多个,值和值之间用“;“隔开,单个值为:<left_start>,<left_length>,<right_start>,<right_length>。
比如:
SEQUENCE_PRIMER_PAIR_OK_REGION_LIST=100,50,300,50 ; 900,60,, ; ,,930,100
表明有引物设计有三种选择:
左引物在100~150bp区间进行设计,右引物在300~350bp的区间进行设计;
左引物在900~960bp区间进行设计,右引物随意;
右引物在930~1030bp区间进行设计,左引物随意。

SEQUENCE_OVERLAP_JUNCTION_LIST (space separated integers; default empty)
指定引物在此位置要有重叠。与标签 PRIMER_MIN_3_PRIME_OVERLAP_OF_JUNCTION 或 PRIMER_MIN_5_PRIME_OVERLAP_OF_JUNCTION 搭配使用。例如:此值设为 20, 默认下,设计的引物的 3' 端必须要包含有 20~23 区段的碱基,或引物的 5' 端必须包含有 20~26 区段的碱基。

SEQUENCE_INTERNAL_EXCLUDED_REGION (interval list; default empty)
该标签的值为:<start>,<length>。中间oligos的设计避开此段子序列。

SEQUENCE_PRIMER (nucleotide sequence; default empty)
该标签的值为左引物序列,该序列必须是 SEQUENCE_TEMPLATE 的一段子序列;程序运行时会检查该子序列并在其附近设计右引物或中间oligos.

SEQUENCE_INTERNAL_OLIGO (nucleotide sequence; default empty)
该标签的值为中间oligo序列,该序列必须是SEQUENCE_TEMPLATE的一段子序列;程序运行时会检查该子序列并在其附近设计左右引物。

SEQUENCE_PRIMER_REVCOMP (nucleotide sequence; default empty)
该标签的值为右引物序列,该序列必须是SEQUENCE_TEMPLATE的一段反向互补子序列;程序运行时会检查该子序列并在其附近设计左引物。

SEQUENCE_START_CODON_POSITION (int; default -2000000)
标记起始密码子的第一个碱基位点,从而让primer3在起始密码子的左边设计左引物,只有在必要时才在起始密码子的右边设计左引物。可以设置为负数,表示起始密码子在整条序列不可见的左边。如果该值为非负数且该位点处的碱基不是ATG,则primer3会发出一个错误信号。如果值小于或等于-10^6,则忽略该参数。primer3从起始密码子开始扫描,直到获得一个终止密码子,理想状况下设计的右引物结束于该终止密码子或其后的碱基。

SEQUENCE_QUALITY (space separated integers; default empty)
该标签的值是一系列的以空格隔开的碱基质量值,每一个碱基质量值对应着模板琏上的碱基。

SEQUENCE_FORCE_LEFT_START (int; default -1000000)
强制左引物的5'端从指定的位点开始,不管是否违背了一些其它的条件。默认值表示左引物的起始位点是任意位置。

SEQUENCE_FORCE_LEFT_END (int; default -1000000)
强制左引物的3'端到指定位点结束,不管是否违背了一些其它的条件。默认值表示左引物的结束位点是任意位置。

SEQUENCE_FORCE_RIGHT_START (int; default -1000000)
强制右引物的5'端从指定的位点开始,不管是否违背了一些其它的条件。默认值表示左引物的起始位点是任意位置。

SEQUENCE_FORCE_RIGHT_END (int; default -1000000)
强制右引物的3'端到指定位点结束,不管是否违背了一些其它的条件。默认值表示左引物的结束位点是任意位置。

全局输入标签 (以 PRIMER_ 开头):

PRIMER_TASK (string; default generic)
该标签确定 primer3 设计什么引物,其有效的值有:
1. generic
最常用的值,即设计引物对。
2. pick_detection_primers
该值已取消。和 generic 等同。
3. check_primers
对引物进行检测,必须要在 SEQUENCE_PRIMER, SEQUENCE_INTERNAL_OLIGO 和 SEQUENCE_PRIMER_REVCOMP 标签中提供引物序列。
4. pick_primer_list
仅仅独立地进行 left primers, right primers 和 internal oligos 的提取。
5. pick_sequencing_primers
设计引物适合于对目标区域进行测序。可以使用 SEQUENCE_TARGET 来制定多个目标区间。为了获得最佳的测序结果,来计算 primer 的位置。此时会忽略 PRIMER_PRODUCT_SIZE_RANGE 所设定的值。
6. pick_cloning_primers
设计引物用于基因的克隆。必须使用 SEQUENCE_INCLUDED_REGION 来指定模板序列中起始和终止的核苷酸位置(ORF的起始与终止)。从而,程序仅能调节引物序列的长度。同时再进行 PRIMER_PICK_ANYWAY=1 而不考虑和其它引物设计条例相冲突,来强制设计引物。
7. pick_discriminative_primers
设计多态性扩增引物,用于将一个引物的序列强制设计到一个多态性位点上。和上一个方法原理一致,需要设置 SEQUENCE_INCLUDED_REGION 和 PRIMER_PICK_ANYWAY=1 标签。
8. pick_pcr_primers
该值已取消。等同于:
PRIMER_TASK=generic
PRIMER_PICK_LEFT_PRIMER=1
PRIMER_PICK_INTERNAL_OLIGO=1
PRIMER_PICK_RIGHT_PRIMER=1
9. pick_pcr_primers_and_hyb_probe
该值已取消。 和 pick_pcr_primers 完全一致。
10. pick_left_only
该值已取消。仅设计 left primer.
11. pick_right_only
该值已取消。仅设计 right primer.
12. pick_hyb_probe_only
该值已取消。仅设计 interanl oligo.

PRIMER_PICK_LEFT_PRIMER (boolean; default 1)
是否设计 left primer.

PRIMER_PICK_INTERNAL_OLIGO (boolean; default 0)
是否设计 internal oligo.

PRIMER_PICK_RIGHT_PRIMER (boolean; default 1)
是否设计 right primer.

PRIMER_NUM_RETURN (int; default 5)
最大返回的引物(对)数目。这些引物(对)按其质量进行排序。增大该值将增加运算时间。

PRIMER_MIN_3_PRIME_OVERLAP_OF_JUNCTION (int; default 4)
默认设置下,引物的 3' 端序列和模板指定位置必须要有 4 个碱基的重叠。

PRIMER_MIN_5_PRIME_OVERLAP_OF_JUNCTION (int; default 7)
默认设置下,引物的 5' 端序列和模板指定位置必须要有 4 个碱基的重叠。

PRIMER_MUST_MATCH_FIVE_PRIME (ambiguous nucleotide sequence; default empty)
该值是个 5bp 长度的碱基序列(可以使用各种模糊碱基字符),用于指定引物 5' 端的序列。

PRIMER_INTERNAL_MUST_MATCH_FIVE_PRIME (ambiguous nucleotide sequence; default empty)
和上一个标签一致,用于指定 internal oligo 的 5' 端序列。

PRIMER_MUST_MATCH_THREE_PRIME (ambiguous nucleotide sequence; default empty)
该值是个 5bp 长度的碱基序列(可以使用各种模糊碱基字符),用于指定引物 3' 端的序列。

PRIMER_INTERNAL_MUST_MATCH_THREE_PRIME (ambiguous nucleotide sequence; default empty)
和上一个标签一致,用于指定 internal oligo 的 3' 端的序列。

PRIMER_PRODUCT_SIZE_RANGE (size range list; default 100-300)
指定引物产物的长度范围。

PRIMER_PRODUCT_OPT_SIZE (int; default 0)
指定引物产物最佳的长度值。

PRIMER_PAIR_WT_PRODUCT_SIZE_LT (float; default 0.0)
如果产物长度低于 PRIMER_PRODUCT_OPT_SIZE 的值,所得到的罚分。

PRIMER_PAIR_WT_PRODUCT_SIZE_GT (float; default 0.0)
如果产物长度高于 PRIMER_PRODUCT_OPT_SIZE 的值,所得到的罚分。

PRIMER_MIN_SIZE (int; default 18)
所允许的最小的 primer 长度。

PRIMER_INTERNAL_MIN_SIZE (int; default 18)
所允许的最小的 internal oligo 长度。

PRIMER_OPT_SIZE (int; default 20)
最佳的 primer 长度。

PRIMER_INTERNAL_OPT_SIZE (int; default 20)
最佳的 internal oligo 长度。

PRIMER_MAX_SIZE (int; default 27)
所允许的最大的 primer 长度。

PRIMER_INTERNAL_MAX_SIZE (int; default 27)
所允许的最大的 internal oligo 长度。

PRIMER_WT_SIZE_LT (float; default 1.0)
primer 比最佳长度小所获得的罚分。

PRIMER_INTERNAL_WT_SIZE_LT (float; default 1.0)
internal oligo 比最佳长度小所获得的罚分。

PRIMER_WT_SIZE_GT (float; default 1.0)
primer 比最佳长度大所获得的罚分。

PRIMER_INTERNAL_WT_SIZE_GT (float; default 1.0)
internal oligo 比最佳长度大所获得的罚分。

PRIMER_MIN_GC (float; default 20.0)
primer 中所允许的最小的 GC 含量的百分比。

PRIMER_INTERNAL_MIN_GC (float; default 20.0)
internal oligo 中所允许的最小的 GC 含量的百分比。

PRIMER_OPT_GC_PERCENT (float; default 50.0)
primer 中最佳的 GC 含量值。

PRIMER_INTERNAL_OPT_GC_PERCENT (float; default 50.0)
internal oligo 中最佳的 GC 含量值。

PRIMER_MAX_GC (float; default 80.0)
primer 中所允许的最大的 GC 含量的百分比。

PRIMER_INTERNAL_MAX_GC (float; default 80.0)
internal oligo 中所允许的最大的 GC 含量的百分比。

PRIMER_WT_GC_PERCENT_LT (float; default 0.0)
primer 比最佳 GC 含量值低所获得的罚分。

PRIMER_INTERNAL_WT_GC_PERCENT_LT (float; default 0.0)
internal oligo 比最佳 GC 含量值低所获得的罚分。

PRIMER_WT_GC_PERCENT_GT (float; default 0.0)
primer 比最佳 GC 含量值高所获得的罚分。

PRIMER_INTERNAL_WT_GC_PERCENT_GT (float; default 0.0)
internal oligo 比最佳 GC 含量值高所获得的罚分。

PRIMER_GC_CLAMP (int; default 0)
要求 left primer 和 right primer 的 3' 末端序列中有连续指定数目的 Gs 或 Cs 碱基。

PRIMER_MAX_END_GC (int; default 5)
left primer 和 right primer 的 3' 末端的 5bp 碱基中最多允许有指定数目的 Gs 或 Cs 碱基。

PRIMER_MIN_TM (float; default 57.0)
primer 所允许的最小的 TM 值。

PRIMER_INTERNAL_MIN_TM (float; default 57.0)
internal oligo 所允许的最小的 TM 值。

PRIMER_OPT_TM (float; default 60.0)
primer 最优的 TM 值。

PRIMER_INTERNAL_OPT_TM (float; default 60.0)
internal oligo 最优的 TM 值。

PRIMER_MAX_TM (float; default 63.0)
primer 所允许的最大的 TM 值。

PRIMER_INTERNAL_MAX_TM (float; default 63.0)
internal oligo 所允许的最大的 TM 值。

PRIMER_WT_TM_LT (float; default 1.0)
primer 低于最优 TM 值的罚分。

PRIMER_INTERNAL_WT_TM_LT (float; default 1.0)
internal oligo 低于最优 TM 值的罚分。

PRIMER_WT_TM_GT (float; default 1.0)
primer 高于最优 TM 值的罚分。

PRIMER_INTERNAL_WT_TM_GT (float; default 1.0)
internal oligo 高于最优 TM 值的罚分。

PRIMER_PAIR_WT_DIFF_TM (float; default 0.0)
left primer 和 right primer 的 TM 值不等的罚分。

PRIMER_PRODUCT_MIN_TM (float; default -1000000.0)
扩增产物所允许最小的 TM 值。

PRIMER_PRODUCT_OPT_TM (float; default 0.0)
扩增产物最佳的 TM 值。 0 表示没有最佳的 TM 值。

PRIMER_PRODUCT_MAX_TM (float; default 1000000.0)
扩增产物所允许最大的 TM 值。计算其 TM 值的公式:Tm = 81.5 + 16.6(log10([Na+])) + 0.41*(%GC) - 600/length。其中[Na+]是钠盐的摩尔浓度;length 是序列的长度。此公式也用于计算 primer and internal oligo 的 Tm 值。

PRIMER_PAIR_WT_PRODUCT_TM_LT (float; default 0.0)
扩增产物低于最佳 TM 值的罚分。

PRIMER_PAIR_WT_PRODUCT_TM_GT (float; default 0.0)
扩增产物高于最佳 TM 值的罚分。

PRIMER_TM_FORMULA (int; default 1)
指定 TM 详细计算方法。 0 表示使用以前版本的兼容的算法。

PRIMER_SALT_MONOVALENT (float; default 50.0)
指定 PCR 反应体系中一价盐的浓度。默认是 50 mM。

PRIMER_INTERNAL_SALT_MONOVALENT (float; default 50.0)
和上一个标签一致。但是 for internal oligo.

PRIMER_SALT_DIVALENT (float; default 1.5)
指定 PCR 反应体系中二价盐(一般是 Mg2+)的浓度。默认是 1.6 mM。

PRIMER_INTERNAL_SALT_DIVALENT (float; default 0.0)
和上一个标签一致,但是 for internal oligo.

PRIMER_DNTP_CONC (float; default 0.6)
DNTP 的总浓度。

PRIMER_INTERNAL_DNTP_CONC (float; default 0.0)
和上一个标签一致,但是 for internal oligo.

PRIMER_SALT_CORRECTIONS (int; default 1)
对盐浓度的校正方法。有 0,1,2 三种方法。

PRIMER_DNA_CONC (float; default 50.0)
PCR 完毕后的 annealing oligo 的浓度。

PRIMER_INTERNAL_DNA_CONC (float; default 50.0)
和上一个标签一致,但是 for internal oligo.

PRIMER_THERMODYNAMIC_OLIGO_ALIGNMENT (boolean; default 1)
该值设为 1, 则程序会使用热力学模型来计算 oligos 形成发夹结构和二聚体的可能性。

PRIMER_THERMODYNAMIC_TEMPLATE_ALIGNMENT (boolean; default 0)
该值设为 1, 则程序会使用热力学模型来计算 oligos 退火时结合到模板链上非特异性位点的可能性。

PRIMER_THERMODYNAMIC_PARAMETERS_PATH (string; default ./primer3_config)
此标签指定进行热动力学方法进行分析所需要的相关参数文件的路径。在 Linux 下有 2 个默认的路径:./primer3_config/ and /opt/primer3_config/ ; 在 Windows 下,只有 1 个默认路径:.\primer3_config\ 。此参数在 Linux 下需要设置为 /opt/biosoft/primer3-2.3.6/src/primer3_config/ 。

PRIMER_MAX_SELF_ANY (decimal, 9999.99; default 8.00)
引物进行自我结合的可能性。自我匹配位点则记 1 分;匹配到 N 则记 -0.25 分;不匹配记 -1 分;打开一个 gap 记 -2 分。例如:
   5' ATCGNA 3'
      || | |
   3' TA-CGT 5'
这个产生的得分为 1.75 。得分越低,则表示越不可能发生引物内碱基匹配。

PRIMER_MAX_SELF_ANY_TH (decimal, 9999,99; default 47.00)
和 PRIMER_MAX_SELF_ANY 一样,但是使用热力学方法进行计算。此标签必须在 PRIMER_THERMODYNAMIC_OLIGO_ALIGNMENT=1 时生效。该标签的默认值比 PRIMER_MIN_TM 标签的值低 10 度。

PRIMER_INTERNAL_MAX_SELF_ANY (decimal, 9999.99; default 12.00)
和上上一个标签一致,但是 for internal oligo.

PRIMER_INTERNAL_MAX_SELF_ANY_TH (decimal, 9999.99; default 47.00)
和上上一个标签一致,但是 for internal oligo.

PRIMER_PAIR_MAX_COMPL_ANY (decimal, 9999.99; default 8.00)
和 PRIMER_MAX_SELF_ANY 一致,但是针对 left primer 结合到 right primer 上的可能性。

PRIMER_PAIR_MAX_COMPL_ANY_TH (decimal, 9999.99; default 47.00)
和上一个标签一样,但是使用热力学的方法进行计算。

PRIMER_WT_SELF_ANY (float; default 0.0)
罚分权重 for PRIMER_MAX_SELF_ANY 的值。

PRIMER_WT_SELF_ANY_TH (float; default 0.0)
罚分权重 for PRIMER_MAX_SELF_ANY_TH 的值。

PRIMER_INTERNAL_WT_SELF_ANY (float; default 0.0)
和上上一个标签一致,但是 for internal oligo.

PRIMER_INTERNAL_WT_SELF_ANY_TH (float; default 0.0)
和上上一个标签一致,但是 for internal oligo.

PRIMER_PAIR_WT_COMPL_ANY (float; default 0.0)
罚分权重 for PRIMER_MAX_SELF_ANY 的值。

PRIMER_PAIR_WT_COMPL_ANY_TH (float; default 0.0)
罚分权重 for PRIMER_MAX_SELF_ANY_TH 的值。

PRIMER_MAX_SELF_END (decimal, 9999.99; default 3.00)
引物的 3' 端可能结合到一起,从而形成一条序列,然后以之为模板进行扩增,形成引物二聚体,从而导致目标序列的扩增失败。该参数设置最多允许碱基在 3' 端进行互补配对的得分为 3。例如:
   5' ATGCCCTAGCTTCCGGATG 3'
                ||| |||||
             3' AAGTCCTACATTTAGCCTAGT 5'

or    5` AGGCTATGGGCCTCGCGA 3'
                  ||||||
               3' AGCGCTCCGGGTATCGGA 5'
前者得分为 7, 后者得分为 6 。

PRIMER_MAX_SELF_END_TH (decimal 9999.99; default 47.00)
和上一个标签一样,不过使用热力学方法进行计算。

PRIMER_INTERNAL_MAX_SELF_END (decimal 9999.99; default 12.00)
和上上个标签一样,但是 for internal oligo.

PRIMER_INTERNAL_MAX_SELF_END_TH (decimal 9999.99; default 47.00
和上一个标签一样,不过使用热力学方法进行计算。

PRIMER_PAIR_MAX_COMPL_END (decimal, 9999.99; default 3.00)
和 PRIMER_MAX_SELF_END 标签一样,但是针对 left primer 和 right primer 的 3' 端互补配对。

PRIMER_PAIR_MAX_COMPL_END_TH (decimal, 9999.99; default 47.00)
和上一个标签一样,不过使用热力学方法进行计算。

PRIMER_WT_SELF_END (float; default 0.0)
罚分权重 for PRIMER_MAX_SELF_END 的值。

PRIMER_WT_SELF_END_TH (float; default 0.0)
罚分权重 for PRIMER_MAX_SELF_END_TH 的值。

PRIMER_INTERNAL_WT_SELF_END (float; default 0.0)
和上上个标签一样,但是 for internal oligo.

PRIMER_INTERNAL_WT_SELF_END_TH (float; default 0.0)
和上上个标签一样,但是 for internal oligo.

PRIMER_PAIR_WT_COMPL_END (float; default 0.0)
罚分权重 for PRIMER_PAIR_MAX_COMPL_END 的值。

PRIMER_PAIR_WT_COMPL_END_TH (float; default 0.0)
罚分权重 for PRIMER_PAIR_MAX_COMPL_END_TH 的值。

PRIMER_MAX_HAIRPIN_TH (float; default 47.0)
和 PRIMER_MAX_SELF_ANY_TH 类似,该标签利用热力学方法计算形成发夹结构后所允许最大的 TM 值。该标签的默认值比 PRIMER_MIN_TM 标签的值低 10 度。

PRIMER_INTERNAL_MAX_HAIRPIN_TH (float; default 47.0)
和上个标签一样,但是 for internal oligo.

PRIMER_WT_HAIRPIN_TH (float; default 0.0)
罚分权重 for PRIMER_MAX_HAIRPIN_TH 的值。

PRIMER_INTERNAL_WT_HAIRPIN_TH (float; default 0.0)
罚分权重 for PRIMER_INTERNAL_MAX_HAIRPIN_TH 的值。

PRIMER_MAX_END_STABILITY (float, 999.9999; default 100.0)
left or right primer 的 5bp 3' 碱基最大的稳定性。该值选取 left or right primer 分别计算出的值中的较大值。计算方法为 PRIMER_TM_FORMULA 标签中指定的方法。默认方法下, GCGCG 碱基序列的值为6.86; TATAT 碱基序列的值为 0.86 。如果是使用旧版本的方法,此值则会变大,相应应该设置此值大些。

PRIMER_WT_END_STABILITY (float; default 0.0)
罚分因子 for PRIMER_MAX_END_STABILITY 的值。

PRIMER_MAX_NS_ACCEPTED (int; default 0)
所允许的 N 碱基的数目。

PRIMER_WT_NUM_NS (float; default 0.0)
罚分权重 for PRIMER_MAX_NS_ACCEPTED 的值。

PRIMER_MAX_POLY_X (int; default 5)
所允许的单核苷酸重复的次数,例如, 默认下 AAAAAA 是不允许的。

PRIMER_INTERNAL_MAX_POLY_X (int; default 5)
和上个标签一致,但是 for internal oligo.

PRIMER_MIN_LEFT_THREE_PRIME_DISTANCE (int; default -1)
在给出引物对的时候,任意两个 left primers 的 3' 末端在模板链上的位点相差的碱基数不超过此值。 -1 表示不作限制, 0 则表示两对引物的 left primers 的 3' 末端在模板链上的位点是相同的。

PRIMER_MIN_RIGHT_THREE_PRIME_DISTANCE (int; default -1)
和上一个标签一样,但是 for right primers。

PRIMER_MIN_THREE_PRIME_DISTANCE (int; default -1)
是上两个标签的综合体,用此标签相当于给上两个标签赋予同样的值。不能和上两个标签同时使用。

PRIMER_PICK_ANYWAY (boolean; default 0)
不管如何也要设计引物。

PRIMER_LOWERCASE_MASKING (int; default 0)
如果设置此标签值为 1, 则设计的引物的 3' 端将不会落在小写字符的碱基区间。

后面的还有一小部分,不想写了... 不是特别重要,而且用到的时候很少。

程序输入标签 (以 P3_ 开头):

P3_FILE_ID (string; default empty)
该标签用于指定 primer3 的配置文件路径。和主程序的 -p3_setting_file 的效果一致。
P3_FILE_FLAG (boolean; default 0)
如果设置该值为 1, 则程序会额外生成几个文件,文件名以 SEQUENCE_ID 标签的值为前缀,分别以 .for, .rev 或 .int 为后缀,对应着独立设计的 left primer, right primer 和 internal oligos 的序列。
P3_COMMENT (string; default empty)
该标签的值是无效的。可以用于输入信息的描述。

3.4 PRIMER3 打分

首先,通过一些限定阈值,程序找出所有的引物对。然后再对其进行打分。
默认情况下,引物的 TM 每偏离最优 TM 1 摄氏度,则罚 1 分; 引物的长度每偏离最优长度 1 bp, 则罚 1 分。而没有通过 GC 含量、发夹结构、3′ 端互补、自身互补配对、模板是重复序列 等方式进行罚分。要根据自身需求来进行罚分常数的设定。

4. /opt/biosoft/primer3-2.3.6/my_default_settings.txt

Primer3 File - http://primer3.sourceforge.net
P3_FILE_TYPE=settings

P3_FILE_ID=Settings for Batch task of searching SSR Primers
PRIMER_FIRST_BASE_INDEX=1
PRIMER_TASK=generic
PRIMER_NUM_RETURN=5
PRIMER_PICK_LEFT_PRIMER=1
PRIMER_PICK_INTERNAL_OLIGO=0
PRIMER_PICK_RIGHT_PRIMER=1
PRIMER_PICK_ANYWAY=1
PRIMER_THERMODYNAMIC_PARAMETERS_PATH=/opt/biosoft/primer3-2.3.6/src/primer3_config/


### 引物 TM 值设定: ###
PRIMER_TM_FORMULA=1                 # TM 的计算方法, 1 表示使用 the SantaLucia parameters (Proc Natl Acad Sci 95:1460-65)
PRIMER_MIN_TM=55.0
PRIMER_OPT_TM=60.0
PRIMER_MAX_TM=65.0
PRIMER_PAIR_MAX_DIFF_TM=5.0         # 两个引物之间的 TM 值最多相差 5 摄氏度
PRIMER_WT_TM_LT=1.0
PRIMER_WT_TM_GT=1.0                 # TM 每偏离最优 TM 值 1 摄氏度,则罚 1分。
PRIMER_PAIR_WT_DIFF_TM=0.0


### 引物长度设定: ###
PRIMER_MIN_SIZE=18
PRIMER_OPT_SIZE=20
PRIMER_MAX_SIZE=22
PRIMER_WT_SIZE_LT=1.0
PRIMER_WT_SIZE_GT=1.0                # 引物长度每偏离最优值 1 bp, 则罚 1 分。


### 引物 GC 含量设定:###
PRIMER_MIN_GC=30.0
PRIMER_MAX_GC=70.0
PRIMER_WT_GC_PERCENT_LT=0.0
PRIMER_WT_GC_PERCENT_GT=0.0


### 引物的热力学计算: ###
# 开启热力学计算
PRIMER_THERMODYNAMIC_OLIGO_ALIGNMENT=1

# 引物自身进行反向互补
PRIMER_MAX_SELF_ANY=8.00
PRIMER_WT_SELF_ANY=0.0
PRIMER_MAX_SELF_ANY_TH=45.00
PRIMER_WT_SELF_ANY_TH=0.0

# 引物自身进行 3' 端反向互补形成引物二聚体
PRIMER_MAX_SELF_END=3.00
PRIMER_WT_SELF_END=0.0
PRIMER_MAX_SELF_END_TH=35.00
PRIMER_WT_SELF_END_TH=0.0

# left primer 和 right primer 序列的反向互补
PRIMER_PAIR_MAX_COMPL_ANY=8.00
PRIMER_PAIR_WT_COMPL_ANY=0.0
PRIMER_PAIR_MAX_COMPL_ANY_TH=45.00
PRIMER_PAIR_WT_COMPL_ANY_TH=0.0

# left primer 和 right primer 进行 3' 端反向互补形成引物二聚体
PRIMER_PAIR_MAX_COMPL_END=3.00
PRIMER_PAIR_WT_COMPL_END=0.0
PRIMER_PAIR_MAX_COMPL_END_TH=35.00
PRIMER_PAIR_WT_COMPL_END_TH=0.0

# 发夹结构
PRIMER_MAX_HAIRPIN_TH=24.00
PRIMER_WT_HAIRPIN_TH=0.0

# 3' 端碱基的稳定性
PRIMER_MAX_END_STABILITY=9.0
PRIMER_WT_END_STABILITY=0.0


### 碱基序列设定: ###
PRIMER_LOWERCASE_MASKING=0              # 模板序列中包含小写字符不影响引物设计
PRIMER_MAX_POLY_X=4                     # 引物序列中不能包含单核苷酸连续长度超过 4 bp
PRIMER_MAX_NS_ACCEPTED=0                # 引物中允许的 N 的数目
PRIMER_WT_NUM_NS=0.0                    # 每个 N 的罚分
PRIMER_MAX_END_GC=5                     # 引物中 3' 端 5bp 碱基中允许的最大 Gs 或 Cs 的数目
PRIMER_GC_CLAMP=0                       # 引物中 3' 端碱基中不能出现连续的 Gs 和 Cs 序列
PRIMER_LIBERAL_BASE=1                   # 是否允许有诸如 N A R Y 等类型的碱基。必须在设定PRIMER_MAX_NS_ACCEPTED 不为 0 后方有效。
PRIMER_LIB_AMBIGUITY_CODES_CONSENSUS=0  # 如果设置为 1,则 C 能与 S 完美匹配,任意碱基能和 N 完美匹配。



### 碱基质量设定: ###
PRIMER_MIN_QUALITY=0                    # primer 序列允许最小的碱基质量
PRIMER_MIN_END_QUALITY=0
PRIMER_QUALITY_RANGE_MIN=0
PRIMER_QUALITY_RANGE_MAX=100
PRIMER_WT_SEQ_QUAL=0.0
PRIMER_WT_END_QUAL=0.0


## 引物位置设置: ###
PRIMER_SEQUENCING_LEAD=50               # 该参数仅在 PRIMER_TASK=pick_sequencing_primers 时有效. 表明引物的 3' 端距目标区域有 50bp。
PRIMER_SEQUENCING_SPACING=500           # 该参数仅在 PRIMER_TASK=pick_sequencing_primers 时有效. 该值决定了在同一条链上的两个引物的距离.
PRIMER_SEQUENCING_INTERVAL=250          # 该参数仅在 PRIMER_TASK=pick_sequencing_primers 时有效. 该值决定了在不同链上的两个引物的距离。
PRIMER_SEQUENCING_ACCURACY=20           # 该参数仅在 PRIMER_TASK=pick_sequencing_primers 时有效. 该值决定了引物设计的区间.
PRIMER_OUTSIDE_PENALTY=0
PRIMER_INSIDE_PENALTY=-1.0
PRIMER_WT_POS_PENALTY=0.0


### 非引物数据库设定:###
# 启用非引物数据库
#PRIMER_MISPRIMING_LIBRARY=

# 引物比对到非引物数据库
#PRIMER_MAX_LIBRARY_MISPRIMING=12.00
#PRIMER_WT_LIBRARY_MISPRIMING=0.0
#PRIMER_PAIR_MAX_LIBRARY_MISPRIMING=20.00
#PRIMER_PAIR_WT_LIBRARY_MISPRIMING=0.0

# 模板比对到非引物数据库
#PRIMER_MAX_TEMPLATE_MISPRIMING=12.00
#PRIMER_WT_TEMPLATE_MISPRIMING=0.0
#PRIMER_MAX_TEMPLATE_MISPRIMING_TH=40.00
#PRIMER_WT_TEMPLATE_MISPRIMING_TH=0.0
#PRIMER_PAIR_MAX_TEMPLATE_MISPRIMING=24.00
#PRIMER_PAIR_WT_TEMPLATE_MISPRIMING=0.0
#PRIMER_PAIR_MAX_TEMPLATE_MISPRIMING_TH=70.00
#PRIMER_PAIR_WT_TEMPLATE_MISPRIMING_TH=0.0

### PCR 反应体系设定: ###
PRIMER_SALT_MONOVALENT=50.0                 # 单价盐离子浓度(mM)
PRIMER_SALT_CORRECTIONS=1                   # PRIMER_SALT_CORRECTIONS=1 means use the salt correction in SantaLucia et al 1998
PRIMER_SALT_DIVALENT=1.5                    # 二价镁离子浓度(mM)
PRIMER_DNTP_CONC=0.6                        # 总dNTPs浓度(mM)
PRIMER_DNA_CONC=50.0                        # DNA产物浓度(mM)


### PCR 产物的设定: ###
PRIMER_PRODUCT_MIN_TM=-1000000.0
PRIMER_PRODUCT_OPT_TM=0.0
PRIMER_PRODUCT_MAX_TM=1000000.0
PRIMER_PAIR_WT_PRODUCT_TM_LT=0.0
PRIMER_PAIR_WT_PRODUCT_TM_GT=0.0
PRIMER_PRODUCT_OPT_SIZE=0
PRIMER_PAIR_WT_PRODUCT_SIZE_LT=0.0
PRIMER_PAIR_WT_PRODUCT_SIZE_GT=0.0


## 罚分因子: ###
PRIMER_PAIR_WT_PR_PENALTY=1.0                # left primer和right primer的罚分之和。此和乘以此系数,再加其它罚分作为最终罚分。
PRIMER_PAIR_WT_IO_PENALTY=0.0                # 将 internal oligo 的罚分乘以此系数,加入到引物的最终罚分中。


### internal oligo 的设定:###
# TM 值
PRIMER_INTERNAL_MIN_TM=57.0
PRIMER_INTERNAL_OPT_TM=60.0
PRIMER_INTERNAL_MAX_TM=63.0
PRIMER_INTERNAL_WT_TM_LT=1.0
PRIMER_INTERNAL_WT_TM_GT=1.0

# 长度
PRIMER_INTERNAL_MIN_SIZE=18
PRIMER_INTERNAL_OPT_SIZE=20
PRIMER_INTERNAL_MAX_SIZE=27
PRIMER_INTERNAL_WT_SIZE_LT=1.0
PRIMER_INTERNAL_WT_SIZE_GT=1.0

# GC 含量
PRIMER_INTERNAL_MIN_GC=20.0
PRIMER_INTERNAL_MAX_GC=80.0
PRIMER_INTERNAL_OPT_GC_PERCENT=50.0
PRIMER_INTERNAL_WT_GC_PERCENT_LT=0.0
PRIMER_INTERNAL_WT_GC_PERCENT_GT=0.0

# 热力学
PRIMER_INTERNAL_MAX_SELF_ANY=12.00
PRIMER_INTERNAL_WT_SELF_ANY=0.0
PRIMER_INTERNAL_MAX_SELF_END=12.00
PRIMER_INTERNAL_WT_SELF_END=0.0

# 碱基序列
PRIMER_INTERNAL_MAX_POLY_X=5
PRIMER_INTERNAL_MAX_NS_ACCEPTED=0

# 碱基质量
PRIMER_INTERNAL_MIN_QUALITY=0
PRIMER_INTERNAL_WT_END_QUAL=0.0
PRIMER_INTERNAL_WT_SEQ_QUAL=0.0

# 非引物数据库
#PRIMER_INTERNAL_MAX_LIBRARY_MISHYB=12.00
#PRIMER_INTERNAL_WT_LIBRARY_MISHYB=0.0

# PCR 反应体系
PRIMER_INTERNAL_SALT_MONOVALENT=50.0
PRIMER_INTERNAL_SALT_DIVALENT=1.5
PRIMER_INTERNAL_DNA_CONC=50.0
PRIMER_INTERNAL_DNTP_CONC=0.0
=

使用Misa结合Primer3来批量设计SSR引物

MISA,英文全称为MIcroSAtellite identification tool,即微卫星识别工具。

MISA是使用 perl 编写的一支程序,能识别出序列中的微卫星和复合微卫星(两个微卫星之间由由不多于100bp的碱基对隔开),并给出其所在位点。

MISA下载网址:http://pgrc.ipk-gatersleben.de/misa/misa.html

MISA用法:

$ misa.pl filename
misa.pl <FASTAfile>
其中,fastfile是序列文件,同时在运行程序的工作目录下必须有一个名称为“misa.ini”
的文件。该文件内容为:
definition(unit_size,min_repeats): 1-10 2-6 3-5 4-5 5-5 6-5
interruptions(max_difference_for_2_SSRs): 100 
该文件指定了misa的参数,即1个碱基重复10次及10次以上;2个碱基重复6次及6 次以上;
3个碱基重复5次及5次以上;4个碱基重复5次及5次以上;5个碱基重复5 次及5次以上;6碱
基重复5次及5次以上,这样的碱基重复序列才算是微卫星序列。 同时,两个微卫星之间的距
离小于100bp的时候,两个微卫星组成一个复合微卫星。

MISA的输出结果:

MISA会在 Fastafile 所在的文件夹下生成两个文件,分别是 “<FASTfile>.misa” 和 “<FASTfile>.statistics”

"<FASTfile>.misa" :以表格的形式列出微卫星的类型和位点;
"<FASTfile>.statistics" :统计微卫星的类型和频数。

在MISA的下载页面中,提供了3个附加的 perl 脚本,分别是:Get_est_trimmer.plp3_in.plp3_out.pl

由于MISA程序读取fasta文件中的序列ID,将序列ID中的空格用下划线 ”_” 填补了,所以在fasta文件中,其序列ID最好不要有空格。否则运行接下来的程序时,会出问题。

Get_est_trimmer.pl

针对EST序列,可以除去EST序列中短的序列和两端不明确的碱基。

p3_in.pl

输入 misa.pl 的输出结果,将引物设计的参数文件(模板,产物长度,目标区域等)导入到一个以“p3in”为后缀的文件中。

$ p3_in.pl filename.misa

 调用 primer3_core

该软件详细解说见:http://www.hzaumycology.com/chenlianfu_blog/?p=284,生成结果文件 filename.p3in。使用primer3-2.3.5版本的时候,MISA官网提供的p3_in.pl的结果不符合primer3-2.3.5的输入格式,故需要修改p3_in.pl和p3_out.pl文件。

 $ primer3_core -default_version=1 -output=filename.p3out filename.p3in

 p3_out.pl

对primer3产生的文件进行提取合,得到最后的结果文件 filename.result

$ p3_out.pl filename.p3out filename.misa

p3_in.pl 和 p3_out.pl 这两支程序需要修改才能正常使用。

修改过后的两支程序下载:p3_in.plp3_out.pl

结果文件示例

ID      SSR nr. SSR type        SSR     size    start   end     FORWARD PRIMER1 (5'-3') Tm(°C)  size    REVERSE PRIMER1 (5'-3') Tm(°C)  size    PRODUCT1 size (bp)      start (bp)      end (bp)        FORWARD PRIMER2 (5'-3') Tm(°C)  size    REVERSE PRIMER2 (5'-3') Tm(°C)  size    PRODUCT2 size (bp)      start (bp)      end (bp)        FORWARD PRIMER3 (5'-3') Tm(°C)  size    REVERSE PRIMER3 (5'-3') Tm(°C)  size    PRODUCT3 size (bp)      start (bp)      end (bp)
scaffold1_254817_bp     1       p3      (GCC)5  15      17114   17128   TGATGTCCTAGTGCGTCTCG    60.008  20      CATCCTGTCTTTGAACGGGT    59.966  20      226     17022   17247   TGATGTCCTAGTGCGTCTCG    60.008  20      ACATCCTGTCTTTGAACGGG    59.966  20      227     17022   17248   TGATGTCCTAGTGCGTCTCG    60.008  20      TGAGGGAGTTGTGGTGATGA    60.088  20      142     17022   17163
scaffold1_254817_bp     2       p1      (T)10   10      116694  116703  ATTGCAACCACCAAAGAAGG    59.971  20      CTCCAGGCGCTACGTTAATC    59.867  20      151     116600  116750  TCCCTACTGCATTGACCTCC    60.073  20      CTCCAGGCGCTACGTTAATC    59.867  20      223     116528  116750  ACACTGCCTTCGATTCATCC    60.081  20      CTCCAGGCGCTACGTTAATC    59.867  20      246     116505  116750
scaffold1_254817_bp     3       p3      (TTG)5  15      142162  142176  TCCACAACCCAATTTACGGT    60.088  20      CCAGAGTATGCCTGGTTCGT    60.134  20      212     142004  142215  GTCCACAACCCAATTTACGG    60.088  20      CCAGAGTATGCCTGGTTCGT    60.134  20      213     142003  142215  GCCAGTTTTGACAGGCGTAT    60.140  20      CCAGAGTATGCCTGGTTCGT    60.134  20      236     141980  142215
scaffold1_254817_bp     4       p3      (AGG)5  15      145861  145875  TTCGAGCTCGTCTGGTAGGT    60.012  20      ATTTATCGTCCAGTGCCCAG    59.955  20      212     145725  145936  GTTCGAGCTCGTCTGGTAGG    60.012  20      ATTTATCGTCCAGTGCCCAG    59.955  20      213     145724  145936  GGTTCGAGCTCGTCTGGTAG    60.012  20      ATTTATCGTCCAGTGCCCAG    59.955  20      214     145723  145936
scaffold1_254817_bp     5       p3      (CAT)6  18      182964  182981  TGTAGAGGGAGGCTGAGGAA    59.943  20      TTGCGAAAAGCAAGGAGAGT    60.132  20      270     182913  183182  TGTAGAGGGAGGCTGAGGAA    59.943  20      GCAAGGAGAGTCGGGTATGA    60.218  20      261     182913  183173  TGTAGAGGGAGGCTGAGGAA    59.943  20      AAAAGCAAGGAGAGTCGGGT    60.247  20      265     182913  183177
scaffold1_254817_bp     6       p3      (GAT)5  15      220964  220978  ATTGATACCGGTGGGTGAAA    60.051  20      TTGAAGGAACTTCGAATGGG    60.044  20      263     220929  221191  ATTGATACCGGTGGGTGAAA    60.051  20      TCGAATGGGATCAACTTTCC    59.871  20      252     220929  221180  ATTGATACCGGTGGGTGAAA    60.051  20      GGAACTTCGAATGGGATCAA    59.871  20      258     220929  221186
scaffold1_254817_bp     7       p2      (TA)7   14      249504  249517  TACCATGAGAAGGGGGAATG    59.744  20      TTTTCTCGACACGTCTGCAC    60.032  20      230     249469  249698  AGAAGGGGGAATGCAAAGTC    60.443  20      TTTTCTCGACACGTCTGCAC    60.032  20      223     249476  249698  GAGAAGGGGGAATGCAAAGT    60.443  20      TTTTCTCGACACGTCTGCAC    60.032  20      224     249475  249698
scaffold2_167145_bp     1       p3      (CTG)7  21      3622    3642    GGAGATATTTCCTCAGGGGC    59.866  20      AGGCAATGTCGATGCTATCC    60.066  20      240     3445    3684    AGAAGCAGAAGGAGGTGCAG    59.745  20      AGGCAATGTCGATGCTATCC    60.066  20      184     3501    3684    GGGAGATATTTCCTCAGGGG    59.722  20      AGGCAATGTCGATGCTATCC    60.066  20      241     3444    3684
scaffold2_167145_bp     2       p2      (AT)6   12      82759   82770   CCATCCCTCTTCCTCTTCCT    59.630  20      ACAAGGTGATGCACAATCCA    59.967  20      222     82648   82869   CCATCCCTCTTCCTCTTCCT    59.630  20      CACAAGGTGATGCACAATCC    59.967  20      223     82648   82870   CCATCCCTCTTCCTCTTCCT    59.630  20      CCACAAGGTGATGCACAATC    59.967  20      224     82648   82871
scaffold3_156598_bp     1       c       (GGA)5(GGT)5    30      70534   70563   GGTGGATGTATTGGATTGCC    60.021  20      GGCATGAACGACTTTTTGCT    60.257  20      183     70469   70651   GGGTGGATGTATTGGATTGC    60.021  20      GGCATGAACGACTTTTTGCT    60.257  20      184     70468   70651   TGGTCATCGAGCTGATGGTA    60.225  20      GGCATGAACGACTTTTTGCT    60.257  20      229     70423   70651
scaffold3_156598_bp     2       p3      (CGC)6  18      80301   80318   CCGAAAAGGCCATTAGTTCA    60.067  20      ACGACGAATGAAACCCTTTG    59.971  20      255     80206   80460   CCGAAAAGGCCATTAGTTCA    60.067  20      TTGGGGTGAGTTCCTTATCG    59.926  20      238     80206   80443   TCTTTTGACTTCGATGCCCT    59.813  20      ACGACGAATGAAACCCTTTG    59.971  20      216     80245   80460

SequenceServer的安装

出自http://www.sequenceserver.com/

SequenceServer这个软件开发出来不久。其作用是将 blast+ 整合到本地网络中。能自动识别出本地的数据库,界面简洁易用。和 wwwblast 功能类似。

1 安装需要Ruby (>= 1.8.7), RubyGems (>= 1.3.6), and NCBI BLAST+ (>= 2.2.25+).

# yum install *ruby*  #安装很多ruby相关软件,不然运行会出问题,特别是缺少
ruby-devel
NCBI BLAST+ ftp://ftp.ncbi.nih.gov/blast/executables/blast+/LATEST/

2 安装sequenceserver

# gem install sequenceserver

3  配置sequenceserver

# sequenceserver   运行软件生成带有注释的配置文件
# vim .sequenceserver.conf 将注释文件进行修改,加入blast路径
bin: ~/ncbi-blast-2.2.26+/bin/
database: /Users/me/blast_databases/

4 blast数据库的创建

<1> 使用 sequenceserver 来调用 makeblastdb 进行数据库创建

$ sequenceserver format-databases directory_with_fasta_files

<2> 使用 blast+ 本身所带 makeblastdb 来创建数据库

$ makeblastdb -dbtype <db type> -title <db title> -in <db> -parse_seqids

5 通过 passenger 在 Apache 或 Nginx 上运行 SequenceServer

# gem install passenger

for apache2 
# passenger-install-apache2-module   #按提示进行设置

for nginx 
# passenger-install-nginx-module

6 部署服务器设置

在 https://github.com/yannickwurm/sequenceserver 下载 sequencesercer 
将 sequenceserver-0.8.0.3.zip 解压到 /var/www/sequenceserver/. 使该文
件夹下存在 public 这个文件夹。
# vim /etc/http/conf/http.conf  并加入以下数行

<VirtualHost *:80>
    DocumentRoot /var/www/sequenceserver/public
    ServerName http://sequenceserver.hzaumycology.com
#前提条件是申请了sequenceserver.hzaumycology.com这个域名
    <Directory /var/www/sequenceserver/public>
        AllowOverride all
        Options -MultiViews
    </Directory>
</VirtualHost>

# /etc/init.d/httpd restart

7 运行 passenger 和 sequenceserver

# nohup passenger start &
# nohup sequenceserver &

若想要开机运行这两个命令,则
# vim /etc/rc.local  在末尾添加 passenger start &
sequenceserver &

 8 修改 .sequenceserver.conf 配置文件

可以选择修改端口,默认端口为4567,则要
# vim /etc/sysconfig/iptables   在相应位置加入一行
-A RH-Firewall-1-INPUT -m state --state NEW -m tcp -p tcp --dport 4567 -j ACCEPT
# /etc/init.d/iptables restart   重启iptables服务

根据服务器配置修改线程数,默认为1.

9  以上配置不出现错误,则可以在浏览器中进行 sequenceserver 访问了.

浏览器中输入: http://sequenceserver.hzaumycology.com:4567

则会出现一个简洁漂亮的界面了!

Bowtie2用法祥解

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

懒人必看

对参考序列构建index
$ bowtie2-build genome.fasta index

尝试使用前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> inputreads中,跳过前<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-end mode) 
--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中允许含有不确定碱基(GTAC,通常为N)的最大数目. 
Default: L,0,0.15. 计算公式为: f(x) = 0 + 0.15 * x, 表示长度为100read
最多运行存在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: 匹配时,
read15'端上游, 和前导链一致, read23'下游, 和前导链反向互补. 或者read2
上游, read1在下游反向互补; --rf: read15'端上游, 和前导链反向互补, read2
3'端下游, 和前导链一致; --ff: 两条reads都和前导链一致. Default: --fr. 默认
设置适合于Illuminapaired-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 read1read2的关系为dovetail的时候,该状况算为和谐比对. 默认情况
dovetail不算和谐比对.

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

--no-overlap read1read2的关系为有重叠的时候, 该状况不算为和谐比对. 默认情
况下两个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 不记录@SQSAM 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> 无视indexoffrate, <int>取代之. Index默认的<int> 
值为5. <int>值必须大于indexoffrate, 同时<int>越大, 耗时越长,耗内存越少.

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

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

--mm 使用内存定位的I/O来载入index, 而不是常规的文件I/O. 从而使多个bowtie
序共用内存中同样的index, 节约内存消耗.

 其它参数:

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

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

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

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


服务器行情

好久没关注服务器配置与价格了,今天通过淘宝查看了下行情,如下:

配件 名称 报价网址 价格 数量 价格
CPU Intel/英特尔XEON E5-2620 服务器CPU http://item.taobao.com/item.htm?spm=a1z10.3.17.21&id=14632262825& 2200 2 4400
主板 Intel/英特尔S2600CP2服务器主板 http://item.taobao.com/item.htm?spm=a1z10.3.17.21&id=14021314905& 2900 1 2900
内存 HY/现代16G DDR3 (RECC)REG ECC 1333服务器内存 http://item.taobao.com/item.htm?spm=a1z10.3.17.29&id=12489828234& 720 8 5760
硬盘 Seagate/希捷2TB SATA2 NS 企业级硬盘 http://item.taobao.com/item.htm?spm=a1z10.3.17.13&id=12389025877& 1620 2 3240
电源 航嘉HK600 额定功率500W 服务器电源 http://item.taobao.com/item.htm?spm=a1z10.3.17.21&id=15969008757& 480 1 480
机箱 道和—牧马人X800塔式服务器机箱/工作站机箱(标准版) http://item.taobao.com/item.htm?spm=a1z10.3.17.38&id=14798052919& 390 1 390
总计 17170