QTL-seq流程说明文档

1QTL-seq流程说明文档版本号v1.0撰写日期：2017.6.26撰写者：柯文斯1目录一、分析流程...............................................................................................................................2二、调参示例.........................................................................................................................72示例：/lustre/Work/project/genome/20170523_Oryza_sativa_BSA/01.QTL-seq/一、工作原理所需文件：1个亲本数据，2个混池数据1.将亲本数据比对到参考基因组，进行snp检测；2.将参考基因组的snp位点进行碱基替换，构建新的reference；3.将亲本数据比对到新的参考基因组，进行snp检测，用于后续混池的过滤；4.将混池数据比对到新的参考基因组进行snp检测，筛选出相对亲本特有的SNP位点；5.对两个混池特有的snp计算出snp-index值，利用窗口滑动的方法结合boost模拟曲线，定位性状关联区域。二、分析流程1.数据准备对一个亲本、2个混池的fastq文件进行数据链接。有多个lane数据的，先做数据合并，再链接合并的结果。合并的参考脚本：zcatL7_1.clean.fq.gzL7_2.clean.fq.gz|gzipm1_1.clean.fq.gz链接后的fastq文件命名方式：BA_1_1_sequence.txt.gz、BA_1_2_sequence.txt.gz为混池BA的fq1、fq2文件，其中BA_1的“1”是必须的，可以用其他数字代替。2.设置参数修改配置文件config.txt，根据需要设置相关的参数。exportPATH=${PATH}:/lustre/Work/software/common/fastx_toolkit/bin3运行命令：$./Bat_make_common.fnc.sh该脚本运行时间很短，直接在本地命令行运行即可。创建流程所需的目录，生成流程的参数配置文件0.common/common.fnc，用于后续调用参数。3.数据过滤$cd1.qualify_read/分别对亲本和两个混池进行数据过滤。运行命令：$./Run_all_Bats.sh0$./Run_all_Bats.sh1$./Run_all_Bats.sh9这一步运行时间较长，所以需要用qsub投递任务。任务脚本：run.0.pbs、run.1.pbs、run.9.pbs具体过滤条件：q30p90，即reads中不低于90%的碱基质量值大于30。选取能配对的reads。对于两个混池，选取同样大小的数据量，即从数据量较多的一个混池中随机提取与另一个混池相同的数据量。4.构建reference$cd2.make_consensus/运行命令：$./Run_all_Bats.sh这一步运行时间较长，所以需要用qsub投递任务。任务脚本：run.pbs。具体运行的步骤：利用bwaaln，将亲本过滤后的reads比对参考基因组。利用covalrefine，利用covalcall，检测亲本中的SNP、Indel。4对于亲本中检测到的SNP位点，替换参考基因组的碱基，从而得到一个新的reference。$cd90.align_to_this_fasta/运行命令：$./Run_all_Bats.sh这一步运行时间较长，所以需要用qsub投递任务。任务脚本：run.pbs。具体运行的步骤：利用bwaaln，将亲本过滤后的reads比对新reference。利用covalrefine，利用covalcall，检测亲本中的SNP、Indel。最后得到的结果文件：rice_q30p90_MSR_Cov_10_S-snp.pileup，用于过滤后面分析中的假阳性SNP。这个文件在起始的配置文件config.txt中的参数名称为PileupDB。5.检测SNP$cd3.alignment/分别对混池A和混池B进行SNP检测。运行命令：$./Run_all_Bats.sh0$./Run_all_Bats.sh1这一步运行时间较长，所以需要用qsub投递任务。任务脚本：run.0.pbs、run.1.pbs。具体运行的步骤：利用bwaaln，将某一混池过滤后的reads比对新reference。利用covalrefine，利用covalcall，检测混池中的SNP、Indel。提取混池和亲本（上一步亲本比对新reference）共同的SNP位点，540.exclude_common_snps/mybulk_BA_q30p90_MSR_Cov_2_S-snp-common-pos.pileup。提取混池中相对亲本特有的SNP位点，40.exclude_common_snps/mybulk_BA_q30p90_MSR_Cov_2_S-snp-rmc2snp.pileup。对混池特有的SNP位点进行过滤（最低深度、最小变异质量、SNP-index最小值），得到50.awk_custom/mut_index_2/mybulk_BA_q30p90_cov2_co3.txt。其中cov2的2是用于Covalmismatchfilter的阈值，co3的3是用于深度过滤的阈值。6.比较两个混池的SNP-index$cd4.search_for_pair/运行命令：$./Run_all_Bats.sh这一步运行时间较长，所以需要用qsub投递任务。任务脚本：run.pbs。具体运行的步骤：对于上一步最后得到的最后结果，将混池A和混池B的SNP位点进行分类，共有的：10.paired_or_unpaired/mut_index_2/paired_mybulk_BA_q30p90_cov2_co3.txt、特有的：10.paired_or_unpaired/mut_index_2/unpaired_mybulk_BA_q30p90_cov2_co3.txt。将混池A特有的SNP位点与混池B的bam文件进行对比，得到共有的以及特有的SNP位点。对共有的SNP位点，分别计算在2个混池中的SNP-index（2个混池中均满足最低深度要求）。将混池A和混池B共有的SNP位点的数据一起整合到上述结果中，得到所有相关SNP位点在两个混池中的信息，40.merge_paired/mut_index_2/merge_mybulk_BABm_q30p90_paired_cov2_co3.txt67.比较两个混池的SNP-index$cd5.compare/运行命令：./Run_all_Bats.sh这一步运行时间较长，所以需要用qsub投递任务。任务脚本：run.pbs。具体运行的步骤：模拟数据，计算置信区间的值。对所有相关的SNP位点，计算不同置信区间的值。过滤结果（最低深度、SNP-index的最大值及最小值）。滑窗口，计算每个窗口的SNP-index。8.文件格式文件：mybulk_BA_q30p90_MSR_Cov_2_S-snp-rmc2snp.pileupfield#Description1chromosome2coordinate3referencebase4consensusbase5consensusquality6SNPquality7mappingquality8thenumberofreadscoveringthesite9readbases10basequalities11SNP-index文件：filtered_merge_mybulk_BABm_q30p90_paired_pvalue_sldwnd2M50K_cov2_7co5.txtfield#Description1chromosome2coordinate=windowcenterposition(incrementedslidingshift)3actualwidthofwindow4slidingwindowaverageofthefield8thin(type4);bulk-A’sdepth5slidingwindowaverageofthefield9thin(type4);bulk-A’ssnp-index6slidingwindowaverageofthefield19thin(type4);bulk-B’sdepth7slidingwindowaverageofthefield20thin(type4);bulk-B’ssnp-index8slidingwindowaverageofthefield23rdin(type4);smallerdepth9slidingwindowaverageofthefield24ndin(type4);delta(SNP-index)10slidingwindowaverageofthefield30thin(type4);U9511slidingwindowaverageofthefield31stin(type4);L9512slidingwindowaverageofthefield32ndin(type4);L9513slidingwindowaverageofthefield33rdin(type4);L95~U9514slidingwindowaverageofthefield34thin(type4);U9515slidingwindowaverageofthefield35thin(type4);U9916slidingwindowaverageofthefield36thin(type4);L9917slidingwindowaverageofthefield37thin(type4);L9918slidingwindowaverageofthefield38thin(type4);L99~U9919slidingwindowaverageofthefield39thin(type4);U9920SNPscountsinthiswindow三、调参示例1.深度阈值调参第一步：修改config,txt文件8第二步：运行程序./Bat_make_common.fnc.sh这时确认配置文件0.common/common.fnc中深度参数是否调整过来：第三步：进入5.compare文件夹，运行程序：./Run_all_Bats.sh（运行时间较长）即可得到不同深度阈值对应的定位区域了。