WES分析流程-从fp.gz开始（一）
目录
基础背景
服务器：
linux
环境配置：
anaconda3
python3.8
openjdk version "1.8.0_282"
R version 3.2.3
软件：
fastp 0.20.1
bwa 0.7.17
samtools 1.7
gatk4
数据及存储位置：
wes文件夹下 - raw(双端测序数据) - SRR3399199.1_1.fastq.gz - SRR3399199.1_2.fastq.gz - truseq-dna-exome-targeted-regions-manifest-v1-2.bed - clean - align - genome - h19_VCF
一、质量控制
用到的软件：fastp
输入：这里用到的是双端测序，如果是单端测序只需要 -I -O，将双端测序的原始数据输入，得到两个质控后的文件S.1_1_clean.fq.gz S.1_2_clean.fq.gz，放在clean文件夹下。
fastp -i SRR3399199.1_1.fastq.gz -I SRR3399199.1_2.fastq.gz -o /nfs-test/hanwq/wes/clean/S.1_1_clean.fq.gz -O /nfs-test/hanwq/wes/clean/S.1_2_clean.fq.gz
参考：https://github.com/OpenGene/fastp#simple-usage
目前文件夹下内容：
wes文件夹下
- raw(双端测序数据)
- SRR3399199.1_1.fastq.gz
- SRR3399199.1_2.fastq.gz
- truseq-dna-exome-targeted-regions-manifest-v1-2.bed
- clean
- S.1_1_clean.fq.gz
- S.1_2_clean.fq.gz
- align
- genome
- h19_VCF
二、比对
1. 人类基因组参考文件下载
下载地址https://hgdownload.soe.ucsc.edu/goldenPath/hg19/bigZips/hg19.fa.gz 这里我选择的是hg19
下载下来的文件为fa的压缩格式：fa.gz，为了下面步骤的方便，在这里先对文件进行解压。
wget https://hgdownload.soe.ucsc.edu/goldenPath/hg19/bigZips/hg19.fa.gz #下载hg19
mkdir backup #新建文件夹backup
cp hg19.fa.gz backup/ #保存下载的原文件备份到backup
gunzip hg19.fa.gz #解压成fa文件2. 创建三种索引
第一种：FM-index用bwa创建
bwa index -a bwtsw hg19.fa
# gz-a bwtsw 这一参数适用于全基因组。其余的基因组直接用index就可以这一步耗时比较久，最后结束会生成5个文件：hg19.fa.ann 、hg19.fa.amb 、hg19.fa.bwt 、hg19.fa.pac 、hg19.fa.sa
第二种：samtools faidx
按照基因组位置创建索引文件，为了方便快速访问基因。
samtools faidx hg19.fa
#这一步会生成一个文件 h19.fa.fai第三种：生成dict文件
GATK CreateSequenceDictionary -R hg19.fa -O hg19.dict确保以上生成的文件都在同一文件夹genome下面，现在wes文件夹内容如下：
wes文件夹下
- raw(双端测序数据)
- SRR3399199.1_1.fastq.gz
- SRR3399199.1_2.fastq.gz
- truseq-dna-exome-targeted-regions-manifest-v1-2.bed
- clean
- S.1_1_clean.fq.gz
- S.1_2_clean.fq.gz
- align
- genome
- backup
- hg19.fa.gz
- hg19.dict
- hg19.fa.amb
- hg19.fa.ann
- hg19.fa.bwt
- hg19.fa.fai
- hg19.fa.pac
- hg19.fa.sa
- h19_VCF
3. 比对
准备工作完成之后就可以开始比对了，这里采用了官方最推荐的mem算法，输出文件为S-pe.sam放到文件align下面
注意：在这一步比对的时候不能使用nohup，否则在后面格式转换的时候，sam文件识别不了，原因是因为nohup的log信息会夹杂到S-pe.sam文件中，导致失败。
bwa mem ../genome/hg19.fa S.1_1_clean.fq.gz S.1_2_clean.fq.gz > ../align/S-pe.sam三、排序
1. sam->bam格式转换
samtools view -bS aln-pe.sam > S-pe.bam
2.排序
samtools sort -@ 4 -m 4G -O bam -o S-pe-sorted.bam ..S-pe.bam
四、标记PCR重复并建立索引
这里用到的gatk4的集成tools picard
# 去重，生成S-pe-sorted.bam 和 mark-dup.txt
picard markduplicationgatk MarkDuplicates -I S-pe-sorted.bam -O S-pe-markdup.bam -M mark-dup.txt
#建立新索引,生成S-pe-markdup.bam.bai
samtools index S-pe-markdup.bam完成这些步骤之后，对已经生成的文件进行梳理，现在wes文件夹下情况：
wes文件夹下
- raw(双端测序数据)
- SRR3399199.1_1.fastq.gz
- SRR3399199.1_2.fastq.gz
- truseq-dna-exome-targeted-regions-manifest-v1-2.bed
- clean
- S.1_1_clean.fq.gz
- S.1_2_clean.fq.gz
- align
- S-pe.sam
- S-pe.bam
- sort
- S-pe-sorted.bam
- markdup
- S-pe-markup.bam
- S-pe-markup.bam
- genome
- backup
- hg19.fa.gz
- hg19.dict
- hg19.fa.amb
- hg19.fa.ann
- hg19.fa.bwt
- hg19.fa.fai
- hg19.fa.pac
- hg19.fa.sa
- h19_VCF
未完待续。。。