WGS数据分析处理流程

1.拿到数据后先检查数据是否完整。用md5sum命令。

#生成md5文件
ls KPGP*| while read KPGP; do echo $KPGP;md5sum ${KPGP} >> ${KPGP}.md5; done

#检查完整性，全部显示OK即可
md5sum -c *.md5

2.对数据进行质检。

#质检
nohup fastqc -o /data/XXXX/WGS/01fastqc -t 10 *.fq.gz &

#multiqc合并质检报告查看
nohup multiqc * -o /data/XXXX/WGS/01fastqc/multiqc &
#若不合格还要进行用cutadaptor去接头等操作。此处合格，则不赘述。

3.与参考基因组进行比对。

为什么要比对？我们已经知道NGS测序下来的短序列（read）存储于FASTQ文件里面。虽然它们原本都来自于有序的基因组，但在经过DNA建库和测序之后，文件中不同read之间的前后顺序关系就已经全部丢失了。因此，FASTQ文件中紧挨着的两条read之间没有任何位置关系，它们都是随机来自于原本基因组中某个位置的短序列而已。
因此，我们需要先把这一大堆的短序列捋顺，一个个去跟该物种的 参考基因组比较，找到每一条read在参考基因组上的位置，然后按顺序排列好，这个过程就称为测序数据的比对。这 也是核心流程真正意义上的第一步，只有完成了这个序列比对我们才有下一步的数据分析。
我们这里将用于流程构建的BWA就是其中最优秀的一个，它将BW(Burrows-Wheeler)压缩算法和后缀树相结合，能够让我们以较小的时间和空间代价，获得准确的序列比对结果。

#bwa建立索引
bwa index -a bwtsw chrom.37.fa 

#bwa men
vi bwa.sh

for(( i=1 ; i<=6 ; i++ ))
do
  bwa mem -t 4 /data/XXXX/WGS/02hg19/chrom.37.fa \
   /data/XXXX/WGS/KPGP-00001_L${i}_R1.fq.gz \
   /data/XXXX/WGS/KPGP-00001_L${i}_R2.fq.gz > /data/XXXXX/WGS/03BWA/L${i}.sam 
done

nohup sh bwa.sh &

4.sam转为bam，并sort好。

#samtools sam2bam+sort
for(( i=1 ; i<=6 ; i++ ))
do
  samtools view L${i}.sam>L${i}.bam  
done

for(( i=1 ; i<=6 ; i++ ))
do
  samtools sort L${i}.sam L${i}.sort
done

为什么需要排序？FASTQ文件里面这些被测序下来的read是随机分布于基因组上面的，第一步的比对是按照FASTQ文件的顺序把read逐一定位到参考基因组上之后，随即就输出了，它不会也不可能在这一步里面能够自动识别比对位置的先后位置重排比对结果。因此࿰