数据质量控制
把ubuntu同步到D盘
数据质量对后续分析至关重要,不可忽视。
数据质量主要包括以下方面:
- read各个位置的碱基质量值分布
- 碱基的总体质量值分布
- read各个位置上碱基分布比例,目的是为了分析碱基的分离程度
- GC含量分布
- read各位置的N含量
- read是否还包含测序的接头序列
- read重复率,这个是实验的扩增过程所引入的
好的数据应该是达到Q20的碱基要在95%以上(最差不低于90%),Q30要求大于85%(最差也不要低于80%)
质量查看-fastp
对输入和输出数据进行详尽的质量剖析,生成人性化的报告。(quality curves, base contents, KMER, Q20/Q30, GC Ratio, duplication, adapter contents...)
过滤掉低质量,太短和太多N的"bad reads"。
从头部(5')或尾部(3')计算滑动窗内的碱基质量均值,并切除低质量碱基。(类似Trimmomatic,但是非常快)。
头尾剪裁。(上下机序列可能不稳定)
去除接头。(可以不用输入接头序列,算法可以自动识别接头序列并进行剪裁)
不匹配碱基对矫正。(对于双端测序(PE)的数据,软件能够矫正重叠区域的不匹配碱基对。如:重叠区域内,一个碱基质量很高,另一个非常低,fastp会依据质量值进行校正。
修剪尾部(3')的polyG。(修剪ployX尾,可以去掉不想要的ploy。如:mRNA-Seq 中的ployA。)
处理使用了唯一分子标识符(UMI)的数据,并将UMI转换为序列名称。
不仅生成质控和过滤结果的HTML报告(存储在fastp.html,可用-h可指定),还生成了程序可读性非常强的JSON报告(fastp.json,可用-j指定)。
为了适合并行处理,fastp将输出进行分拆。有两种模式可选,a. 指定拆分的个数; b. 指定拆分后每个文件的行数。
fastp支持gzip的输入和输出,同时支持SE和PE数据处理;支持PacBio/Nanopore的长reads数据;也支持交错输入。
用bioconda安装fastp
conda install -c bioconda fastp
fastp -h #查看是否安装成功,运行帮助命令
fastp全指令说明
-尝试对本地下载数据进行质控并生成报告
1.设定工作文件夹
sudo mount -t drvfs D: /mnt/d #把d盘挂到wsl下
cd /mnt/d/seq #打开进入seq这个文件夹 我的数据放在seq文件夹下
ls -l #显示文件夹下文件,看一下是否运行成功
2.按默认条件设置过滤流程
fastp -i SRR31045179_1.fastq -o out.R1.fq -I SRR31045179_2.fastq -O out.R2.fastq #默认会生成fastp.html的质控报告
1.基于短reads组装
利用KmerGenie计算最适合组装的kmer值
conda install -c bioconda kmergenie #安装 KmerGenie
kmergenie --version #查看安装是否成功
mkdir result_haploid #创建结果输出文件
kmergenie fastq_list.txt -o ./result_haploid/haploid -l 15 -k 105 -s 10 -t 4 > haploid.log1.txt 2> haploid.log2.txt
# fastq_list.txt 写有绝对路径的txt文件
#kmergenie: 这是 kmergenie 工具的执行命令。
#fastq_list.txt: 这是包含 fastq 文件列表的文件名。该文件包含了要进行组装的测序reads的路径信息,每行对应一个 fastq 文件。
#-o ./result_haploid/haploid: -o 是输出目录选项,指定输出文件和结果的存放位置。这里指定输出目录为 ./result_haploid/haploid,其中 . 表示当前目录。
#-l 15: -l 是指定 k-mer 长度范围的选项。这里指定 k-mer 长度范围为 15,表示 kmergenie 将尝试使用长度从 15 到 105 的 k-mer 进行组装。
#-k 105: -k 是指定最大 k-mer 长度的选项。这里指定最大 k-mer 长度为 105。
#-s 10: -s 是指定步长大小的选项。这里指定步长大小为 10,表示 kmergenie 将以 10 为步长尝试不同的 k-mer 长度。
#-t 4: -t 是指定线程数的选项。这里指定使用 4 个线程进行组装,可以加快计算速度。
#> haploid.log1.txt 2> haploid.log2.txt: 这两个命令将标准输出(stdout)和错误输出(stderr)分别重定向到 haploid.log1.txt 和 haploid.log2.txt 文件中。> 表示将标准输出重定向,2> 表示将错误输出重定向。
a利用spades进行基因组组装
conda install -c bioconda spades#安装spades
spades.py --pe1-1 Bacillus_subtilis.filt_R1.fastq.gz --pe1-2 Bacillus_subtilis.filt_R2.fastq.gz -t 200 -k 71 -m 600 --careful --phred-offset 33 -o 4-assembly1/Bacillus_subtilis
#--pe1-1、--pe1-2:pair-end测序得到的两个文件名,这里我用的是经过trim后的数据。
#-t:线程数(对于CPU密集型的程序,应该使用少一点的线程,对于IO密集型任务,应该使用多一点的线程)。
#-k:kmer选择的值,一般不需要设置,默认的话21,33,55,77,但是对于我的data(水稻)组装结果明显不好,所以我决定以前面kmergenie预测出来的结果为参考设置k值。
#-m:内存。---跑不下去加到跑的下去,不能超。单位为Gb,注:free -m可以查看最大的内存。
#--careful:运行BayesHammer(read纠错)、SPAdes、MismatchCorrector。
#--phred-offset:碱基质量体系,33或64。
#-o:输出文件目录。
#note:spades还有一个优势,如果你有其他组装工具产生的contigs的话可以,可以使用--trusted-contigs或 --untrusted-contigs选项。前者是当获得高质量组装中使用的,这些contigs将用于graph的构建、间隙填补、处理重复序列。第二个选项用于相对不太可靠的contigs,这些contigs用于间隙填补、处理重复序列。
b评价组装结果
conda install -c bioconda quast#安装软件
quast.py -o ./quast_results/PB_501_trim_quast -R ./0-refer/Bacillus_subtilis.str168.fasta -t 70 ./4-assembly1/Bacillus_subtilis/scaffolds.fasta
#-o输出目录
#-R参考基因组序列
#-t线程数
#后面为要比较的contig/scaffold所在目录。
简书参考