组装参数

参数文件中参数设置对组装过程中每一步都有影响。设定的默认值适用于大部分的组装过程,仅需个别进行调整(比如,文件路径),但应用到个例时需要根据实际情况进行调整。ipyrad的一个主要特点就是参数调整简单易控。

以下是参数设置的具体解释、影响的步骤以及示例文件内容。

0. Assembly name

所有输出文件的前缀,必须是唯一的,是当前计算数据集的代表。多种数据时,一般可使用连续的命名法(e.g. data1, data2)或加上一些关键参数(e.g. data_clust90, data_clust85)。当使用“-n”产生组装名称后,名称不可变,但新的组装可以通过分支命名(见分支模块)。

影响步骤: 1-7; 示例:新的组装通过-n 或-b 产生:

>>> ipyrad -n data1                       ## 产生命名为data1的新的组装任务
>>> ipyrad -p params-data1.txt -b data2   ## 产生命名为data2的分支组装任务

1. Project dir

项目目录;一个项目目录可用来组合不同的相关组装。若目录原先不存在,通过该步骤会产生相应的工作目录。一个好的工作目录名称能够更好地开展工作。在运行过程中,路径一般不能修改,除非整个目录移到了不同的位置/设备上。

影响步骤: 1-7; 示例:

/home/deren/ipyrad/tests/finches   ## [1] 产生/使用名称为finches的项目目录
finches                            ## [1] 产生/使用名称为finches的项目目录

2. Raw fastq path

raw data(未解编数据)的路径(若数据已分好个体,则留空)。输入文件可以使gzip压缩格式(即,具文件末尾有.gz)。若设置了该路径,则必须设置barcode路径。若选择多个或所有文件,使用*;

影响步骤:1; 示例:

/home/deren/ipyrad/tests/data/*.fastq.gz     ## [2] 选择所有的gzip文件
~/ipyrad/tests/data/*.fastq.gz               ## [2] 选择所有的gzip文件
./ipsimdata/rad_example*.fastq.gz            ## [2] 选择特定文件 `rad_example`

3. Barcodes path

barcodes的路径;step1中用来解编的,也能在step2中移除和识别adapter/primer的序列。若数据已解编,留空。

影响步骤: 1-2; 示例:

/home/deren/ipsimdata/rad_example_barcodes.txt    ## [3] 选择barcode文件
./ipsimdata/rad_example_barcodes.txt              ## [3] 选择barcode文件

4. Sorted fastq path

已解编数据路径;若数据已解编,则在第1步直接运行该命令,进行已解编数据的分类,可通过通配符来选择多个文件(*)

影响步骤: 1; 示例:

/home/deren/ipyrad/tests/ipsimdata/*.fastq.gz    ## [4] 选择所有gzip文件
~/ipyrad/tests/ipsimdata/*.fastq                 ## [4] 选择所有fastq文件
./ipsimdata/rad_example*.fastq.gz                ## [4] 选择特定文件 `rad_example`

5. Assembly method

选择组装方法;4种:从头组装(denovo),参考序列(reference),从头组装+参考序列(denovo+reference),从头组装-参考序列(denovo-reference)。无参考基因组,选denovo。后三者都需要提供参考序列(param #6),案例可见运行ipyrad的软件运行指导。

影响步骤: 3, 6; 示例:

denovo                            ## [5] denovo组装
reference                         ## [5] reference组装
denovo+reference                  ## [5] 加reference组装
denovo-reference                  ## [5] 减reference组装

6. Reference sequence

参考基因组路径;参考基因组必须是fasta格式,但无需是全基因组,可以使用任何一种map RAD数据的类型,例如转录组或质体基因组;若无参考基因组,留空。

~/ipyrad/tests/ipsimdata/rad_example_genome.fa   ## [6] 选择fasta文件
./data/finch_full_genome.fasta                   ## [6] 选择fasta文件

7. Datatype

选择数据类型;现在有多种数据类型,主要有普通rad(单酶切),ddrad(双酶切),gbs(GBS)等方法;(详见数据类型)

rad                       ## [7] rad数据类型(1个酶切位点,声波降解法)
pairddrad                 ## [7] 双端ddrad类型(2个不同的酶切位点)
pairgbs                   ## [7] 双端gbs类型(1个酶切位点,双端)

8. Restriction_overhang

酶切序列;若filter_adapter这个参数保持开启状态,则参数可用于第一步数据分类和第二步过滤接头和酶切位点,确定正确的限制酶切位点序列是比较困难的,无需键入酶切识别位点序列,只需键入酶切后序列上的一部分序列即可。例如,酶PstI具有以下会问序列,^表示酶切位点:

5'...C TGCA^G...'3
3'...G^ACGT C...'5

酶切后,DNA片段在酶切位点附近具有CTGCA序列,反向序列的开头则具有TGCAG。最简单的寻找酶切序列的方法:查看解编数据时酶切序列一般在各片段中开头的一部分并保持不变;在未解编数据中,则在barcode文件的最后一部分。使用一下命令来查看fastQ文件的不变的序列:

## 解压缩文件,
## -c标记表示预览,
##  `less` 表示只预览前100行
zless 100 my_R1_input_file.fastq.gz

如此,会得到以下文件。通过查看每一行的开头就能发现限制性酶切序列了。在双酶切中,则不会得到如图所示的第二行;若是双末端序列,_R2_文件将以第二个限制性酶切序列开头。第二个酶切序列只在filter_adapters参数设置为1时用来检测adapters/primers。第二个酶切序列可留空。

@HWI-ST609:152:D0RDLACXX:2:2202:18249:93964 1:N:0:
TGCAGCAGCAAGTGCTATTCGTACAGTCATCGATCAGGGTATGCAACGAGCAGAAGTCATGATAAAGGGTCCCGGTCTAGGAAGAGACGCAGCATTA
+
BDFFHHHHHJJJHIJJJJJIJIGJJJJJJJJJJJJJJJJDGGIJJJJJHIIIJJJJHIIHIGIGHHHHFFFFEDDDDDACCDCDDDDDDDDDBBDC:
@HWI-ST609:152:D0RDLACXX:2:2202:18428:93944 1:N:0:
TGCAGGATATATAAAGAATATACCAATCCTAAGGATCCATAGATTTAATTGTGGATCCAACAATAGAAACATCGGCTCAACCCTTTTAGTAAAAGAT
+
ADEFGHFHGIJJJJJJIJJIIJJJIJJIJGIJJJJJJJJIJJIJJJJIIIGGHIEGHJJJJJJG@@CG@@DDHHEFF>?A@;>CCCDC?CDDCCDDC
@HWI-ST609:152:D0RDLACXX:2:2202:18489:93970 1:N:0:
TGCAGCCATTATGTGGCATAGGGGTTACATCCCGTACAAAAGTTAATAGTATACCACTCCTACGAATAGCTCGTAATGCTGCGTCTCTTCCTAGACC
+
BDFFHHHHHJJJJIJJIJJJJJJJHIJJJJJJJJHIIJJJJIFHIJJJJFGIIJFHIJJIJJIFBHHFFDFEBACEDCDDDDBBBDCCCDDDCDDC:
@HWI-ST609:152:D0RDLACXX:2:2202:18714:93960 1:N:0:
TGCAGCATCTGGAAATTATGGGGTTATTTCACAGAAGCTGGAATCTCTTGGGCAATTTCACAGAATCTGGGAATATCTGGGGTAAATCTGCAAGATC
+
BDEFHHHHHJJJIJJJJJJJJJJCGIIJJGHJJJJJJJJJJIJJIJJJIHIJJJJJJJHHIIJJJJJJJHGHHHFEFFFDEDABDDFDDEDDDDDDA
@HWI-ST609:152:D0RDLACXX:2:2202:20484:93905 1:N:0:
TGCAGAGGGGGATTTTCTGGAGTTCTGAGCATGGACTCGTCCCGTTTGTGCTGCTCGAACACTGACGTTACTTCGTTGATCCCTATGGACTTGGTCA
+
ADEEHHHHGJJHIIJJJJJJIJFHJIIIJJJJIIIJJJJEHHIIGHIGGHGHHHHFFEEEDDDD;CBBDDDACC@DC?<CDDCCCCCCA@CCD>>@:
@HWI-ST609:152:D0RDLACXX:2:2202:20938:93852 1:N:0:
TGCAGAAGCTGGAGATTCTGGGGCAGCTTTGCAGCAAGCTGAAAATTCTGGGGGTCGATCTGCAAGATCGGAAGAGCGGTTCAGCAGGAATGCCGAG

在某些情况下,限制性内切酶具有不止一个序列,例如Apol:AATTY(即 AATTC & AATTT)。此时,可通过模糊代码处理。

影响步骤:1,2; 示例:

TGCAG                     ## [8] 单酶切位点(如 rad, gbs)
TGCAG, AATT               ## [8] 双酶切(如 ddrad, pairddrad)
CWGC                      ## [8] 单酶切 w/ degenerate base

NB:3RAD和SeqCap数据有4个限制性内切酶,若使用这种数据,将所有序列罗列即可:

CTAGA, CTAGC, AATTC               ## [8] 3rad 数据(多个酶切位点)

9. max_low_qual_bases

设置允许一个序列片段里最多有几个低质量碱基;在step2中,当质量低于20时(可在phred_Qscore_offset中设置),3’端的碱基将被剔除。然而,序列片段仍可能包含一些未被剔除的模糊碱基,这会影响下游组装的效率和准确性。该参数限制了允许序列片段中最大的Ns数。默认参数为5,即若一个序列片段中有超过5个低质量碱基的话,则过滤该序列片段。强烈建议不要增加该数值。

影响步骤:2; 示例:

0                      ## [9] 在一个片段中不允许包含任何低质量碱基
5                      ## [9] 在一个片段中允许包含5个低质量碱基

10. Phred_Qscore_offset

测序质量低于20的片段从3’端剔除的碱基个数。默认值为33。一些早先的数据默认值为64,但是这增加了稀有性(数据利用不充分)。读者可通过调控offset的收来改变修剪的阈值。例如,若将offset从33调至23,就相当于将序列片段最低质量分数从20调至10,碱基正确率大概有95%。

影响步骤:2; 示例:

33                 ## [10] 默认值为33,相当于最小分数20
43                 ## [10] offset增加10,相当于最小分数为30
64                 ## [10] 此处的offset适用于以前的数据,相当于最小分数为20

11. mindepth_statistical

第5步生成一致性序列时统计碱基的最小搜索深度。默认参数是6,该值为鉴别测序错误和杂合的最小值;即如果某一位点的深度(覆盖次数)小于6,就可能是测序错误。

影响步骤: 4,5; 示例:

6                 ## [11] 设置最小深度为6
10                ## [11] 设置为10

12. mindepth_majrule

第5步获得一致性序列时,执行碱基多数规则的最小搜索深度。默认参数与mindepth_statistical保持一致,如此,只保留统计碱基。这个值必须<= mindepth_statistical。若低于,则覆盖度>=mindepth_majrule并且<mindepth_statistical的位点就会使用多数规则。若数据集具有非常低的覆盖度以至于许多聚类都由于低的测序深度被排除了,那么降低mindepth_majrule是一个增加大量可用碱基的高效方法,然要注意的是,多数规则下得到的一致性碱基将会降低杂合性。即,若一个数据集的覆盖度大多在5,参数11设置为6,参数12若设置为6,则大部分数据被剔除,可用性很低;若参数12设置为3,则可以保留大部分的数据;但是,若参数12设的非常低,为1,则基本都保留,保留的结果就是降低了杂合性。

影响步骤: 4, 5; 示例:

6                 ## [12] 与mindepth_stat相同的较高值
2                 ## [12] 低于最高值

13. Maxdepth

最高覆盖度。由于文库准备时片段间扩增的差异性,覆盖度间常常是不均匀的;不同文库间差异各异。此外,基因组的重复区因高度相似聚类成高深度的cluster。设置最大深度可避免上述后者的问题,但同时也具有潜在的排除好的聚类的代价。默认值为设的非常高(10,000),可自行调整。

影响步骤:4,5; 示例:

10000             ## [13] 高于聚类去除的最大深度

14. clust_threshold

聚类阈值。判断两序列为同源序列的最小相似度,若高于,则归为一个cluster。该值必须为小数(e.g. 0.90)。不建议高于0.95,因同源序列可能存在的Ns,indels, 测序错误,多态性等使得很难高于这个阈值。

影响步骤: 3, 6; 示例:

0.90              ## [14] 聚类阈值为90%
0.85              ## [14] 聚类阈值为85%

15. max_barcodes_mismatch

设置barcodes的错配率。默认参数为0。barcodes之间通常会有至少2个碱基的不同,故不建议该值设置大于2;(无barcodes文件的不用调)

影响步骤: 1; 示例:

0              ## [15] 错配率设置为0
1              ## [15] 允许1个碱基错配

16. filter_adapters

去除Illumina adapters是非常重要的。这取决于文库准备时选择的大小方案的保真度,至少有一定比例的序列长度要比真实的DNA片段长,这些往往具有primer/adapter在read的末端。这在双酶切(GBS, ddRAD)中十分常见,这类数据用相同的酶切,尤其在GBS数据集里存在很大问题,GBS中序列片段很小,却需要从两段进行测序。 该参数具有三个设置值(0, 1, 2)。 若为0,则仅当序列片段质量低于参数【9】max_low_qual_bases时才被删除; 若为1,则每一个Qscore<20的序列片段都会修剪开头的碱基并剔除含有过多Ns的序列片段。 若为2,则所有序列片段根据Illumina测序接头、第二位的反向剪接位点(若提供)和barcodes(若提供),并依其去除。此过滤器使用cutadapt程序,有一定的容错率。

影响步骤: 2; 示例:

0                ## [16] 不需要过滤adapter
1                ## [16] 过滤条件基于碱基质量
2                ## [16] 严格过滤adapters

17. filter_min_trim_len

定过滤后序列片段的最短长度。在step2中,若filter_adapter>0, 序列片段就会被剪切成更短的长度,该参数即为剪切后最小的长度。默认35,即,若过滤后序列片段短于35bp,则剔除该序列片段。可调至更高。

影响步骤: 2; 示例:

50                ## [17] 最小剪切的序列长度为50
75                ## [17] 最小剪切的序列长度为75

18. max_alleles_consens

测序错误评估后在每个一致性序列片段(个体)中允许存在的唯一性等位基因最大值。默认值为2,适用于二倍体。在该设置下,任何一个基因座中大于2的样本都会被过滤或剔除。若设置为1(单倍体),则在step4中的H-错误率和杂合性估计就固定为0.0,此时碱基是用来估计错误率的。若设置>2,允许存在更多的等位基因,然而,估计杂合性的碱基仍然是以二倍体为假设前提的,即,等位基因杂合频率=50%

影响步骤: 4, 7; 示例:

2                ## [18] 二倍体,去除 >2 等位基因的碱基
1                ## [18] 单倍体,去除 >1 等位基因的碱基
4                ## [18] 二倍体-碱基,去除 >4 等位基因的碱基

19. max_Ns_consens

联合序列(R1, R2)生成时,允许存在的模糊碱基的最大数量。若一个碱基不能够被确定地(统计地)定义,则成为模糊碱基(N)。一致性序列中过多的Ns会影响聚类效果,也会导致排列错误。默认参数为5,5。单末端只用前面的一个数据,双末端用两者。即俩序列片段聚类发生时,若有5个以上的模糊碱基,该组聚类剔除。

影响步骤: 5; 示例:

2                ## [19] 允许一致性序列中最多存在2 Ns
5, 5             ## [19] 允许一致性序列中最多存在8 Ns (R1, R2)

20. max_Hs_consens

联合序列(R1, R2)生成时,允许存在的杂合位点最大值。该过滤帮助去除一些低质量联合过程中的过多的杂合性。默认参数是8,8,单末端只用前面的一个数据,双末端用两者。即,若联合序列后有8个以上的杂合位点,剔除该聚类。

影响步骤: 5; 示例:

2                ## [20] 允许一致性序列中最多存在2 Hs
8, 8             ## [20] 允许一致性序列中最多存在8 Hs (R1, R2)

21. min_samples_locus

个体间比较时,同一个位点上出现的样本数量。默认参数为4,即至少需要4个样本中含有该位点,否则剔除位点。若设成总样本数,即只保留在所有样本中都出现的位点。可调,不要太小。若popfile中出现min_samples值,该参数被覆盖。

影响步骤: 7; 示例:

4                ## [21] 个体间组装最小个体数为4
12               ## [21] 个体间组装最小个体数为12

22. max_SNPs_locus

个体间比较时,允许同一个序列片段上出现的最多SNP数量。可剔除由于多重复区域剔除了过多N snps导致的低质量联合。默认参数是20,20,即同一个序列片段上最多出现20个SNP。降低该值只能增加一些混乱数据集。单末端只用前面的一个数据,双末端用两者。

影响步骤: 7; 示例:

20                  ## [22] 允许在单端的位点上最大的SNP数为20
20, 30              ## [22] 允许在双端的位点上最大的SNP数为20及30

23. max_Indels_locus

个体间比较时,允许同一个序列片段上出现的最多插入缺失(indels)数量。该参数有助于过滤低质量联合,尤其是双末端数据。默认参数是8,8。单末端只用前面的一个数据,双末端用两者。

影响步骤: 7; 示例:

5                ## [23] 允许在单端的位点上最大的indels数为5
5, 10            ## [23] 允许在双端的位点上最大的indels数为5及10

24. max_shared_Hs_locus

个体间比较时,同一个位点上,允许多态性位点出现在多少个体上的个数/比率。该参数用来检测潜在的旁系同源基因,多个样本间共享的杂合位点可能代表了具有固定差异的旁系同源基因的聚类而非真正的杂合位点。默认参数是0.5,即在该位点上杂合子最多出现在50%的个体中,否则就不判定为杂合(可能是旁系同源)。单末端只用前面的一个数据,双末端用两者。

影响步骤: 7; 示例:

0.25             ## [24] 允许异源位点最多出现在25%的样本中
10               ## [24] 允许异源位点最多出现在10个样本中

25. trim_序列片段

设定数据过滤的时候把开始几个碱基删除。有时fastq文件有一些测序错误,offset并不是预期的那样,此时,可通过设定该参数来认为修剪数据集。

影响步骤: 2; 示例:

0, 0, 0, 0       ## [25] 不做任何处理
5, 0, 0, 0       ## [25] 在R1端减去前5个碱基
5, -5, 0, 0      ## [25] 在R1端减去前5个和末端5个碱基
5, 80, 0, 0      ## [25] 在R1端减去前5个碱基,并且剪切的最大长度为80
5, 75, 5, 75     ## [25] 在R1端和R2端减去前5个碱基,并且剪切的最大长度为75

26. trim_loci

比对(alignment)好的数据,在末端出现低质量的可能性比较高,所以可以设置把边缘的碱基删掉几个。例如设成0,5,0,0,即把R1 中右边缘的5个碱基删除。

影响步骤: 7; 示例:

0, 0, 0, 0     ## [26] 没有位点边缘被剪切
5, 0, 0, 0     ## [26] 剪去R1端比对位点的前5个碱基
0, 5, 5, 0     ## [26] 剪去R1端比对位点的后5个碱基和R2端的前5个碱基

27. output_formats

设定输出文件的格式。*是输出所有格式。例如:(详见output formats章节)

p: PHYLIP (Full dataset)
s: PHYLIP (SNPs only)
u: PHYLIP (One SNP per locus)
n: NEXUS
k: STRUCTURE
g: EIGENSTRAT .geno
G: G-PhoCS
v: VCF (SNPs only)

影响步骤: 7; 示例:

*                     ## [27] 生成所有格式的输出文件
n, v, g               ## [27] 只生成 nexus, vcf 和 geno 格式的输出文件
u,k                   ## [27] 只生成 unlinked snps in phylip 和 structure 格式的文件

28. pop_assign_file

根据给定的群体分组文件路径,把样本分组。

/home/user/ipyrad/popfile.txt        ## [28] example...

种群分类信息文件应先在空白txt文件中编辑,用空格区分个体及群体数据。注意拼写和大小写。每一行应包含一个样本名字以及种群名字。其余的解释信息应用“#”做注释。这些信息会在最后的组装过程中覆盖min_samples_locus值。

示例:

Sample1 pop1
Sample2 pop1
Sample3 pop1
Sample4 pop2
Sample5 pop2

# pop1:2 pop2:2