使用Gatk Germline spns-indels Pipeline分析遗传病(耳聋)

这是GATK Best Practice系列学习文章中的一篇，本文尝试使用Gatk Germline spns-indels Pipeline来分析遗传病（耳聋）

数据

这次没有拿到遗传病的室间质评的数据，直接从NCBI上找一些数据来分析。NCBI上搜索deaf，点击第一条搜索结果，最后几经跳转找到数据下载页面：https://trace.ncbi.nlm.nih.gov/Traces/sra/?study=SRP218677

可以看到：Targeted next generation sequencing of 139 known deafness-related genes in 44 deaf patients with mono-allelic GJB2 mutations

这里一共44个样本，有兴趣的可以全部下载来分析。

奇怪的是，这里不是sra下载，是直接下载fastq文件，不管了。使用下载工具，最后得到的文件是C261 R1.fastq.1这种文件名，这里首先改名为C261_R1.fastq，然后使用gzip，bgzip压缩得到 C261_R1.fastq.gz；C261_R2.fastq.gz也是同样.

也可以直接从一下链接下载：

C261_R1.fastq.gz bc72a0b02876273e2457c753a87241aa 54M

C261_R2.fastq.gz a5f6fcac16cbbe2b2a0c98982243a949 54M

概览：

用到的分析系统及分析流程文件

流程用到的变量

分析流程

• 01 INPUT 输入文件

• 02 fastp,因为数据没有去接头，这里用fastp直接去接头，默认参数处理

• 03 Map Sort Dedup 一套操作

• 04-Recalibrator

• 05 Process Dup Bam

• 06-Call SNV INDEL

• 07-Hard Filter 因为call出来的突变数量不够运行VQSR，所以这里选择使用硬过滤（Hard Filter）， 过滤参数使用gatk推荐的数值， 见链接#https://gatk.broadinstitute.org/hc/en-us/articles/360035532412?id=11097

• 08-Annotation，使用Annovar注释结果

• 09-Annotation Filter对注释结果过滤，过滤工具使用python开发，见文章头部链接
  1. 初步过滤：过滤掉突变类型为NAN和synonymous SNV的结果
  2. 初步过滤：过滤掉突变位点深度低于10，或者reads数少于10的结果，因为遗传突变要么是50%要么100%，低于一定reads数该结果的可靠性不高。这里暂时定为10.
  3. 匹配blacklist.csv文件中记录的和耳聋相关的基因，过滤掉无关结果

  

• 10-QC 获取整体数据的QC状态，平均覆盖深度，总的reads数，等等。

• 11-Output 最后的输出文件：

• 12-分析结果匹配：

最后我们得到一个结果表格，根据Clinvar数据库的注释来判断的话，只有图中黄色一行的状态为：

Pathogenic/Likely_pathgenic(致病的/可能致病的)；

其余行记录均为：NAN或者Benign/Like_benign（良性的/可能良性的）

这与下载数据的内容相匹配：Targeted next generation sequencing of deaf patients with mono-allelic GJB2 mutations

说明我们的分析结果准确度还可以。

不足

• 因为未能下载到原数据用到150个和耳聋相关基因的bed文件，这里使用的是hg19全外的bed文件，所以QC工具最后计算的整体数据的QC数值不准确，平均测序深度、测序深度中位数值、1X，4X等深度下的覆盖度也很低，insert size更是不准确。这只有获得精确的bed文件之后可以得到改进。
• 因为流程分析中使用的是一个hg19全外的bed文件，未包含耳聋中部分线粒体突变的区域。这里最终结果不包含线粒体上的耳聋致病突变，
• 分析流程中使用的一个数据文件，blacklist.csv是搜集到的一些和耳聋相关的基因和突变位点，过滤的时候只匹配了基因。数据还相对比较粗糙，如果改进的话改进数据结构和提高数据的质量：如使用这个网站的数据https://hereditaryhearingloss.org/ 或者搜集更高质量的数据，注释时候关联更多的数据库如OMIM