GATK - Read groups

运行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]

bwa -R 参数的使用

# 不加参数R
bwa mem ref.fa read1.fq read2.fq > aln-pe.sam
# 加参数R
bwa mem -R '@RG\tID:group_n\tLB:library_n\tPL:illumina\tPU:unit1\tSM:sample_n' ref.fa read1.fq read2.fq > aln-pe.R.sam
'@RG\tID:group_n\tLB:library_n\tPL:illumina\tPU:unit1\tSM:sample_n'
# ID:每一个Read group独有的ID;
# LB:DNA preparation library identifier
# PL:测序使用的平台: ILLUMINA, SOLID, LS454, HELICOS and PACBIO;
# PU:Platform Unit,由三部分组成: {FLOWCELL_BARCODE}.{LANE};
# SM:样品名
read group参数详细信息见本文后半部分

注:GATK 使用时,PU不是必须要求的;但是PU与ID同时存在时,PU优先级高于ID。

# 查看文件aln-pe.sam与aln-pe.R.sam中的group
samtools view -H aln-pe.sam  | grep '@RG'
# 不会有输出
samtools view -H aln-pe.R.sam  | grep '@RG'
@RG ID:group_n  LB:library_n    PL:illumina PU:unit1    SM:sample_n

# 多个比对的时候,使用循环
ls *.gRNA.fastq | while read id ;
do
bwa mem -R "@RG\tID:group_${id}\tLB:library_${id}\tPL:illumina\tSM:sample_${id}" ref.fa $id > ${id}.R.sam
done
  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
  对一个read group的reads进行重复序列标记时,需要使用LB来区分reads来自那条lane;有时候,同一个库可能在不同的lane上完成测序;为了加以区分,同一个或不同库只要是在不同的lane产生的reads都要单独给一个ID。

从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主流程

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 217,542评论 6 504
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 92,822评论 3 394
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 163,912评论 0 354
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 58,449评论 1 293
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 67,500评论 6 392
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 51,370评论 1 302
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 40,193评论 3 418
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 39,074评论 0 276
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 45,505评论 1 314
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 37,722评论 3 335
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 39,841评论 1 348
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 35,569评论 5 345
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 41,168评论 3 328
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 31,783评论 0 22
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,918评论 1 269
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 47,962评论 2 370
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 44,781评论 2 354

推荐阅读更多精彩内容