生物信息学与算法生物信息学生物信息杂谈

GATK - Read groups

2018-11-13  本文已影响2人  JeremyL

运行GATK时,报错:java.lang.IllegalArgumentException: samples cannot be empty
这个问题在GATK Forum有讨论到;
IllegalArgumentException: samples cannot be empty


使用 ValidateSamFile确定bam文件是否有问题?

Picard ValidateSamFile

Validates a SAM or BAM file.

ValidateSamFile下有两种mode:
VERBOSE : 检查到100个错误之后退出,并且输出错误到终端;
SUMMARYL: 输出结果是一个表格,展示errors和warnings的数目;

$ java -jar picard.jar ValidateSamFile \
   I=input.bam \
   MODE=SUMMARY
   
Error Type  Count
ERROR:MISSING_READ_GROUP    1
WARNING:RECORD_MISSING_READ_GROUP   1287903

问题在bam文件 MISSING_READ_GROUP;请自动屏蔽WARNING结果;


获取GATK格式bam文件

  1. 方法法1

bwa 比对使用参数-R

-R STR: Complete read group header line. ’\t’ can be used in STR and will be converted to a TAB in the output SAM. The read group ID will be attached to every read in the output. An example is ’@RG\tID:foo\tSM:bar’. [null]
  1. 方法2
    Picard AddOrReplaceReadGroups

Replace read groups in a BAM file.This tool enables the user to replace all read groups in the INPUT file with a single new read group and assign all reads to this read group in the OUTPUT BAM file.

$ java -jar picard.jar AddOrReplaceReadGroups \
       I=input.bam \
       O=output.bam \
       RGID=4 \
       RGLB=library1 \
       RGPL=illumina \
       RGPU=unit1 \
       RGSM=sample1
Option  Description
INPUT (String)  Input file (BAM or SAM or a GA4GH url). Required.
OUTPUT (File)   Output file (BAM or SAM). Required.
SORT_ORDER (SortOrder)  Optional sort order to output in. If not supplied OUTPUT is in the same order as INPUT. Default value: null. Possible values: {unsorted, queryname, coordinate, duplicate, unknown}
RGID (String)   Read Group ID Default value: 1. This option can be set to 'null' to clear the default value.
RGLB (String)   Read Group library Required.
RGPL (String)   Read Group platform (e.g. illumina, solid) Required.
RGPU (String)   Read Group platform unit (eg. run barcode) Required.
RGSM (String)   Read Group sample name Required.

Read groups详细介绍

英文原版查看

There is no formal definition of what is a read group, but in practice, this term refers to a set of reads that were generated from a single run of a sequencing instrument.

测序时:
样本建一个库在同一条lane上完成测序产生reads sets可定义为一个Read group;
样本建成分开独立的库测序得到的reads sets也就分属于不同的Read groups;
存在于SAM/BAM /CRAM 文件中Read groups是由一系列标签组成;这些标签代表着测序中的一些技术特征;有了这些数据之后,大家就可以对bam文件进行重复序列标识和碱基质量重新矫正。
GATK要求输入的bam文件包含Read groups,如果没有就会报错。

例子:

samtools view -H sample.bam | grep '@RG'
@RG ID:H0164.2  PL:illumina PU:H0164ALXX140820.2    LB:Solexa-272222    PI:0    DT:2014-08-20T00:00:00-0400 SM:NA12878  CN:BI

GATk 要求read group的格式
  Read group是@RG开始

ID = Read group identifier
  每一个Read group独有的ID;
  Illumina 测序数据中,read group IDs由flowcell ,lane name 和number组成。
  在矫正碱基质量时,read group IDs对区分技术批次效应是必须的;在这过程中,同一read group的reads假 定为有一样的技术误差。

PU = Platform Unit
  Platform Unit由三部分组成: {FLOWCELL_BARCODE}.{LANE}.{SAMPLE_BARCODE}
  {FLOWCELL_BARCODE} refers to the unique identifier for a particular flow cell;
  The {LANE} indicates the lane of the flow cell ;
  The {SAMPLE_BARCODE} is a sample/library-specific identifier;
  GATK 使用时,PU不是必须要求的;但是PU与ID同时存在时,PU优先级高于ID。

SM = Sample
  reads属于的样品名;SM要设定正确,因为GATK产生的VCF文件也使用这个名字,
  PL = Platform/technology used to produce the read
  测序使用的平台: ILLUMINA, SOLID, LS454, HELICOS and PACBIO.
  LB = DNA preparation library identifier
  有时候,同一个库可能在不同的lane上完成测序;这种情况下进行重复序列标记时,需要使用LB。

从read names中提取ID与PU

H0164ALXX140820:2:1101:10003:23460
H0164ALXX140820:2:1101:15118:25288

H0164____________ #portion of @RG ID and PU fields indicating Illumina flow cell
_____ALXX140820__ #portion of @RG PU field indicating barcode or index in a multiplexed run
_______________:2 #portion of @RG ID and PU fields indicating flow cell lane

例子-Multi-sample and multiplexed example
三个样品:MOM, DAD, KID;
建库:每个样品分别建两个库,一个insert为200bp,一个insert为400bp;
上机:每个测序库使用Illumina HiSeq的 两条lane;
reads:来自 3 x 2 x 2 = 12条lane,可以产生12个bam文件,结果如下:

Dad's data:
@RG     ID:FLOWCELL1.LANE1      PL:ILLUMINA     LB:LIB-DAD-1 SM:DAD      PI:200
@RG     ID:FLOWCELL1.LANE2      PL:ILLUMINA     LB:LIB-DAD-1 SM:DAD      PI:200
@RG     ID:FLOWCELL1.LANE3      PL:ILLUMINA     LB:LIB-DAD-2 SM:DAD      PI:400
@RG     ID:FLOWCELL1.LANE4      PL:ILLUMINA     LB:LIB-DAD-2 SM:DAD      PI:400

#解释
样本DAD建了两个库LIB-DAD-1(200bp)和LIB-DAD-2(400bp),使用了FLOWCELL1上面的4条lane;

Mom's data:
@RG     ID:FLOWCELL1.LANE5      PL:ILLUMINA     LB:LIB-MOM-1 SM:MOM      PI:200
@RG     ID:FLOWCELL1.LANE6      PL:ILLUMINA     LB:LIB-MOM-1 SM:MOM      PI:200
@RG     ID:FLOWCELL1.LANE7      PL:ILLUMINA     LB:LIB-MOM-2 SM:MOM      PI:400
@RG     ID:FLOWCELL1.LANE8      PL:ILLUMINA     LB:LIB-MOM-2 SM:MOM      PI:400

Kid's data:
@RG     ID:FLOWCELL2.LANE1      PL:ILLUMINA     LB:LIB-KID-1 SM:KID      PI:200
@RG     ID:FLOWCELL2.LANE2      PL:ILLUMINA     LB:LIB-KID-1 SM:KID      PI:200
@RG     ID:FLOWCELL2.LANE3      PL:ILLUMINA     LB:LIB-KID-2 SM:KID      PI:400
@RG     ID:FLOWCELL2.LANE4      PL:ILLUMINA     LB:LIB-KID-2 SM:KID      PI:400

参考:
Read groups
Picard
从零开始完整学习全基因组测序(WGS)数据分析:第4节 构建WGS主流程

上一篇下一篇

猜你喜欢

热点阅读