GATK4流程學(xué)習(xí)之背景知識與前期準(zhǔn)備 - 簡書
GATK4流程學(xué)習(xí)之DNA-Seq variant calling(Germline：SNP+INDEL) - 簡書
GATK4流程學(xué)習(xí)之RNA-Seq variant calling(SNP+INDEL) - 簡書
補(bǔ)：Mutect2+scRNAseq+cancer cell - 簡書

說明：由于一些原因，中途在一個新服務(wù)器賬號創(chuàng)建了GATK分析環(huán)境，故后面系列分析的路徑可能與在下文的路徑不一致，但數(shù)據(jù)與軟件都是一致的。

要點(diǎn)一、GATK學(xué)習(xí)

1、GATK簡介

• The GATK is the industry standard for identifying SNPs and indels in germline DNA and RNAseq data.
• Its scope is now expanding to include somatic short variant calling, and to tackle copy number (CNV) and structural variation (SV).

variant calling pipeline

• 簡單理解就是gatk4是根據(jù)bam文件，生成vcf文件的軟件；不僅如此，gatk開發(fā)團(tuán)隊(broad institute)對整個從fatsq→vcf分析流程都建立了標(biāo)準(zhǔn)的分析pipeline，即GATK Best Practices系列

關(guān)于SNP、INDEL等變異類型可參考之前的VCF格式詳解筆記
(插一句就是我登錄broad institute GATK頁面總是有問題，不知道其他小伙伴也遇到類似問題。)
生信格式之fasta、fastq - 簡書 http://www.itdecent.cn/p/5bd5848eb596
生信格式之sam、bam - 簡書 http://www.itdecent.cn/p/f0f1f293f0bd
生信格式之vcf格式 - 簡書 http://www.itdecent.cn/p/34c1e22c92c8

2、相關(guān)概念區(qū)別

2.1 DNA-seq與RNA-seq

https://sciberg.com/resources/bioinformatics-faqs/the-differences-between-dna-and-rna-sequencing.html
（1） DNA-seq

• 如下圖，DNA-seq包括三種測序手段，分別為Whole Genome Sequencing (WGS,全基因組測序), Whole Exome Sequencing (WES or WXS,全外顯子測序) and targeted sequencing(靶向測序).
• WGS是對樣本整個基因組的全部測序，而WES則僅對能攜帶遺傳信息，參與編碼mRNA的外顯子序列(僅占基因組大小的3%)進(jìn)行測序。
• 以WGS與WES為代表的DNA-seq，主要用于研究rare mutations and/or common variants associated with a disorder or phenotype.
  
  DNA-seq

（2）RNA-seq

• 如下圖，RNA-seq主要是捕捉DNA的轉(zhuǎn)錄產(chǎn)物mRNA以及非編碼RNA（lncRNA，circRNA和miRNA等），分為mRNA-seq、miRNA-seq、circRNA， Whole Transcriptome Sequencing (WTS,全轉(zhuǎn)錄組測序)。
• 相比DNA-seq的測序步驟，RNA-seq首先需要提取特定類型RNA，再反轉(zhuǎn)錄成cDNA(complementary DNA,互補(bǔ)DNA)，然后構(gòu)建文庫，進(jìn)行測序。
• 相比DNA-seq的測序分析，RNA-seq的研究包括the detection of changes in gene expression, alternative splicing, post-transcriptional modifications, gene fusions as well as detection of mutations and SNPs.

RNA-seq

2.2、germline mutation與somatic mutation

https://www.zhihu.com/question/38765318
（1）germline mutation 胚系突變

• germline mutation是指上一代的生殖細(xì)胞（germ cells）精子或卵子發(fā)生突變(如下圖左)，然后經(jīng)減數(shù)分裂，形成合子，在子代中不斷分化增殖(有絲分裂，直接復(fù)制)，從而在該個體的所有體細(xì)胞中都存在germline mutation。
• 即取正常組織測序，在某一特定位點(diǎn)，germline突變的頻率理論上只有2種：50%突變（精子或卵子一方突變），或100%突變（精子與卵子均突變）【該個體的生殖細(xì)胞也是帶有突變】
• 所以胚系突變的特點(diǎn)是可遺傳性。如下圖右是僅父代精子胚系突變，導(dǎo)致該個體產(chǎn)生的精子中會有50%的遺傳性

• germline mutation是遺傳性疾病的研究重點(diǎn)；只有一少部分癌癥，是與遺傳相關(guān)的（研究最廣泛的遺傳性癌癥就是乳腺癌，攜帶BRCA1/2基因的突變會導(dǎo)致患乳腺癌、卵巢癌的幾率增加）。

  
  
  image.png

（2）somatic mutation 體細(xì)胞突變

• 如上圖有，somatic mutation與精卵子配體是否發(fā)生突變無關(guān)，而是在胚胎后期發(fā)育過程中，體細(xì)胞分裂過程中發(fā)生的突變。由于體細(xì)胞已經(jīng)高度分化，僅影響該類體細(xì)胞(皮膚，肝臟，骨髓，眼睛等的細(xì)胞均為體細(xì)胞)相關(guān)區(qū)域。
• 由于大部分somatic mutation 不會影響生殖細(xì)胞，所以somatic mutation是不會遺傳的。
• 絕大多數(shù)癌癥，都是由于后天體細(xì)胞突變導(dǎo)致；研究時一般取癌組織與癌旁組織對比研究，即在Call Somatic mutations 的時候最好有同一個體的正常組織進(jìn)行參照。

3、筆記內(nèi)容

• 基于GATK Best Practices的identifying SNPs and indels in germline DNA and RNAseq data的流程學(xué)習(xí)；
• 主要以用為主，通過示例數(shù)據(jù)操作為主，同時再盡量解釋清楚每一步的含義，但背后深入算法還是并不太明白，例如pairHMM算法。

1、下載相關(guān)軟件

• sra-tools、aspera 是兩個常用的下載公共數(shù)據(jù)庫測序數(shù)據(jù)的軟件；
• fastqc、trimmomatic是對fastq測序文件質(zhì)控的兩個軟件；
• bwa、star是兩個常用的比對軟件，各有所長；
• GATK4是variant calling的常規(guī)軟件，目前已發(fā)布第4版本；
• 其它seqtk、tree......

軟件安裝一般到官網(wǎng)或者github主頁，根據(jù)提示下載安裝即可；有的是解壓即用，有的需要make之類的操作(編譯)一下。建議選擇合適的文件路徑，方便以后管理方便。

https://github.com/lh3/bwa

2、conda創(chuàng)建GATK分析環(huán)境

• 區(qū)別上述的方法，conda環(huán)境下可以軟件命令的操作更加方便，不需要考慮環(huán)境變量因素。
• conda的基礎(chǔ)學(xué)習(xí)可參考前面的筆記--Linux的conda軟件管家http://www.itdecent.cn/p/84a0d5c407aa

conda create -n GATK python=3
conda activate GATK
conda install -c bioconda -y sra-tools  seqtk
conda install -c bioconda -y fastqc trimmomatic samtools
conda install -c bioconda -y bwa gatk4
# aspera比較特殊，需從hcc channel源下載
conda install -c hcc aspera-cli
conda list

• 根據(jù)后面的踩坑教訓(xùn)，有兩個軟件需要安裝指定版本才可以

conda install -c bioconda -y star=2.7.1a
conda install -c bioconda  -y sambamba=0.6.6
conda list

但是還是建議手動安裝下上述所有軟件，我是分別建立了一個GATKconda環(huán)境與biosoft文件加下安裝了上述軟件。

3、下載參考數(shù)據(jù)庫

#部分?jǐn)?shù)據(jù)集特別大，耗時，建議后臺運(yùn)行
mkdir -p ~/path/to/GATK/bundle/hg38
cd ~/path/to/GATK/bundle/hg38

（1）下載參考基因組

nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/Homo_sapiens_assembly38.fasta.gz >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/Homo_sapiens_assembly38.fasta.fai >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/Homo_sapiens_assembly38.dict >/dev/null 2>&1 &

如下，bwa與star是兩個測序數(shù)據(jù)比對軟件，比對時需要建立索引文件。根據(jù)GATK流程推薦，bwa適合DNA-seq數(shù)據(jù)找變異；star適合RNA-seq數(shù)據(jù)找變異

（2）bwa建立參考基因組(human)索引

#比較耗時，1-2h
mkdir bwa_index 
cd bwa_index
nohup ~/biosoft/bwa/bwa-0.7.15/bwa index -a bwtsw -p gatk_hg38 ../Homo_sapiens_assembly38.fasta >/dev/null 2>&1 &

（3）下載star的參考基因組(human)索引

• 由于STAR建立索引十分耗資源，因此這里下載搭建好的STAR軟件比對人類參考基因組數(shù)據(jù)的全套數(shù)據(jù)(31G)。因?yàn)檫@套數(shù)據(jù)里的比對索引是star 2.7.1a建立的，故后面比對時需要使用對應(yīng)版本的star，以及找變異時使用版本一致的基因組文件。

mkdir  /home/shensuo/biosoft/star/STAR-2.7.7a/db/
cd /home/shensuo/biosoft/star/STAR-2.7.7a/db/
#網(wǎng)速好的話，一晚上可以下載好
wget -c https://data.broadinstitute.org/Trinity/CTAT_RESOURCE_LIB/GRCh38_gencode_v22_CTAT_lib_Apr032020.plug-n-play.tar.gz
# -c參數(shù)表示斷點(diǎn)續(xù)傳，下載大文件時建議使用
tar -zcvf GRCh38_gencode_v22_CTAT_lib_Apr032020.plug-n-play.tar.gz
cd GRCh38_gencode_v22_CTAT_lib_Apr032020.plug-n-play/ctat_genome_lib_build_dir/
gatk CreateSequenceDictionary -R ref_genome.fa
ls

image.png

（4）下載人類基因組參考變異注釋數(shù)據(jù)

nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/dbsnp_146.hg38.vcf.gz >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/dbsnp_146.hg38.vcf.gz.tbi >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz.tbi >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/1000G_phase1.snps.high_confsampleence.hg38.vcf.gz >/dev/null 2>&1 &
nohup wget ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/hg38/1000G_phase1.snps.high_confsampleence.hg38.vcf.gz.tbi >/dev/null 2>&1 &

nohup搭配&是后臺不斷線的下載。因?yàn)橛械臄?shù)據(jù)比較大，以及建立索引都比較耗時。
此外都是人類測序的相關(guān)分析數(shù)據(jù)。