个人基因组比对及其变异分析

这里用到生物学软件有

1.bwa（主要用途a.建立参考基因组索引文件的。b.进行比对）类似的建立索引的软件有bowtie，hisat

http://sourceforge.net/projects/bio-bwa/files/

2.samtools（主要用途a.排序b.合并c.建立基因组index）

http://sourceforge.net/projects/samtools/files/samtools/

该软件需要编译安装

3.bcftools（主要用途：view命令来进行SNP和Indel calling，它是samtools中附加的软件）

4.snpEff （主要用途：a.添加标签b.标记dbSNP的ID号）

5.freebayes（主要用途：添加一些标签）

(以上软件需要配置环境变量，添加方法见这里)

一、基因组下载

下载参考基因组hg19，所在的网站如下：

http://hgdownload.cse.ucsc.edu/goldenPath/hg19/bigZips/chromFa.tar.gz

下载至~/reference/genome/hg19/目录

下载过后解压，cat合并成hg19.fa一个文件

下载需要比对的个人基因组，这里我拿到的是韩国人个人基因组计划中编号为KPGP-00001的基因组，所在ftp站点：

ftp://biodisk.org/Release/KPGP/KPGP_Data_2016_Release_Candidate/WGS/KPGP-00001

下载至~/Data/KPGP-00001目录下

下载过后注意md5校验，md5sum + 文件名；

二、bwa软件建立索引

新建路径~/reference/index/bwa/，用以存放参考基因组的索引。

bwa index -a(STR BWT construction algorithm) bwtsw -p(STR prefix of the index [same as fasta name]) ~/reference/index/bwa/hg19 ~/reference/genome/hg19/hg19.fa 1>hg19.bwa.index.log 2>&1

因为运行时间较长，此命令最好在后台运行，(nohup time,加上time可以看见运行时间)

三、将fastq文件比对到hg19参考基因组

bwa mem(PE150测序比对模式) -t (INT number of threads) -M(mark shorter split hits as secondary 选择上一步建立好的参考基因组索引) ~/reference/index/bwa/hg19 KPGP-00001_R1.fq.gz KPGP-00001_R2.fq.gz 1>KPGP-00001.sam 2>KPGP-00001.bwa.log

若有多个fq文件，需要加for流程控制，得到多个sam文件。此命令在后台运行。

至此，我们用已经比对得到了sam格式比对结果文件，从sam文件中我们看到每条reads在基因组及染色体的位置

图一

若图片不清晰，可放大网页查看。