Adding Read Group info to a set of Bam files

[sarvidsson]

thanks for you quick reply,

my rg.txt is ID:FP SM:F11-12 LB:FP PL:Illumina
ID:FP SM:L6-6 LB:FP PL:Illumina

with your code "\@RG... looks like. should it be like this?

@RG ID:FP SM:F11-12 LB:FP PL:Illumina
@RG ID:FP SM:L6-6 LB:FP PL:Illumina

with try leaving -r in samtools,

appreciate your quick reply

Yes, that looks OK. Varsågod!

[kapr0007]

thanks for you quick reply,

my rg.txt is ID:FP SM:F11-12 LB:FP PL:Illumina
ID:FP SM:L6-6 LB:FP PL:Illumina

with your code "\@RG... looks like. should it be like this?

@RG ID:FP SM:F11-12 LB:FP PL:Illumina
@RG ID:FP SM:L6-6 LB:FP PL:Illumina

with try leaving -r in samtools,

appreciate your quick reply

[sarvidsson]

Hi, this is the command i used

perl -e ’print "@RG\tID:FP\tSM:F11-12\tLB:FP\tPL:Illumina\[email protected]\tID:FP\tSM:L6-6\tLB:FP\tPL:Illumina\n"’ > rg1.txt
samtools merge -rh rg1.txt merged_f11-l6 l6-6.bam f11-12.bam

Check your rg1.txt, you need to escape the '@'s:

perl -e ’print "\@RG\tID:FP\tSM:F11-12\tLB:FP\tPL:Illumina\n\@RG\tID:FP\tSM:L6-6\tLB:FP\tPL:Illumina\n"’ > rg1.txt

gives a correct header file. I guess you can leave out "-r" in the samtools call, didn't test that however (given that the single files have correctly set read groups already, if not you can skip the -h and the rg1.txt file and only use the -r feature).

[kapr0007]

Hi, this is the command i used

perl -e ’print "@RG\tID:FP\tSM:F11-12\tLB:FP\tPL:Illumina\[email protected]\tID:FP\tSM:L6-6\tLB:FP\tPL:Illumina\n"’ > rg1.txt
samtools merge -rh rg1.txt merged_f11-l6 l6-6.bam f11-12.bam

[sarvidsson]

Thank you for your reply,

I merged into one bam file before i suppled to freebayes, the above grep for my bam file looks empty, but when it view it has rg tags RG:Z:f11-12

thanks

Then something went wrong in merging those BAMs (what did you use for this step?); your header got truncated. If the individual BAM files have intact headers - with the read groups in there, you can simply feed Freebayes with several BAM files.

[dpryan]

What tool did you use to merge the BAM files and what do you get if you perform the grep command that sarvidsson posted on one of the unmerged BAM files? The problem you're running into is due to not having read group information in the BAM header. Most tools will read this to create a dictionary of valid read groups that alignments can have.

[kapr0007]

Thank you for your reply,

I merged into one bam file before i suppled to freebayes, the above grep for my bam file looks empty, but when it view it has rg tags RG:Z:f11-12

thanks

[sarvidsson]

Did you merge into one bam file (and if so, with which tool) or supply freebayes with several BAM files?

Could you post the output from

samtools view -H file.bam | grep "^@RG"

for either your merged file (if you have one) or a few of your single sample files and the exact freebayes call you used?

[kapr0007]

Hello all,

Thank you all for your posts,
first I would like to tell you all, I am trying to add RG tags and sample ID to a list of sorted bam files, i have tried all the above mentioned scripts. i was able to add the RG tags to all my files, later calling with freebayes for SNP, my ouput vcf file has only the header part
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT L6-6

the rest is empty. could some one help me with this, all i need is after i call for SNP, i would like to have all the genotypes labeled with the specific sample IDs.

Thank you all

[blancha]

I'm not an expert on the SAM format, and it's not a particularly interesting subject, to be honest.
Some of the RG tags only appear in the header, I believe.
In my BAM files, only the RGID appears on each alignment line, but you set it to null.
RGLB, RGPL, RGSM and RGPU only appear in the header.

If you have trouble falling asleep, you can read the SAM format specification.

http://samtools.github.io/hts-specs/SAMv1.pdf

Only ID is specified as being used in the RG tags of alignment records.

If I'm wrong, feel free to correct me.

[pliang]

Thanks for your quick reply, Blancha! The use 12 threads was unintended and was really unnecessary, and I know the use -H/h in samtools view, which just displays the header. But I was talking about the actual RG tags added to each alignment line, which was missing from mine bam file generated by picard AddOrReplaceReadGroups. I am looking into the perl script posted by Brugge.

[blancha]

Why in the world would you use 12 threads to add or replace read groups?

-XX:ParallelGCThreads=12

More importantly, you need to specify the parameter -H (or -h) to view the header.

samtools view -H A1_M1.bam
or
samtools view -h A1_M1.bam | more

[pliang]

Hello, Wjeck and others on this blog:

Picard AddOrReplaceReadGroups does not seem to work for me. Below is the command I used and the first part of bam file it generated. No @RG was added.

Am I missing anything? Your help is greatly appreciated.


java -Xmx4g -XX:ParallelGCThreads=12 -jar /work/nrap1100/bin/picard-tools-1.78/AddOrReplaceReadGroups.jar I=Sample_RS-01812720/merged.bam O=A1_M1.bam RGID=null RGLB=$d RGPL=Illumina RGSM=A1_M1 RGPU=TAAGGCG

samtools view A1_M1.bam |head
8LSZMS1:104:C4KJNACXX:1:1203:8209:39039 145 chr10 3100000 37 101M chr12 59597178 0 AGAATTCTCACCTGAGAAATACCGAATGGCAGAGAAACACCTGAATAAAATGTTCAACATCCTTAATCATCAGGGAAATGCAAATCAAAACAACACTGAGA EDDEECEEEDFEFFFHHHHHHIJJJJJJJJJJIJIIJIHJJIJJIGIJJIIJIIGIJJIJJJJJJJIJJIJJJJJJJJJJJJJJIJJJHHHHHFEFDF@BB X0:i:1 X1:i:0 MD:Z:101 XG:i:0 AM:i:0 NM:i:0 SM:i:37 XM:i:0 XN:i:1 XO:i:0 XT:A:U
8LSZMS1:104:C4KJNACXX:1:2105:3168:69401 99 chr10 3100097 57 101M = 3100236 240 GAGATTCCACTTCACTCCAGTTAGAATGGCTAAGATCAAAAACTCAGGTGACAACAGATGTTGGCGAGGATGTGGAGAAAGGGGAACACTCCTCCATTGTT CC@FFFFFHHHHHJJJJJJJIIJGIIJJJJJJJJIIIIJJJJJGGIGHBGHIJJJJJJJIIGHIJJJJHFFFDEFEDECCCBDDDDDDDACDDDDDDDDED X0:i:1 X1:i:2 XA:Z:chr4,+72093156,101M,1;chr7,+144726933,101M,1; MD:Z:101 XG:i:0 AM:i:20 NM:i:0 SM:i:20 XM:i:0 XO:i:0 XT:A:U
8LSZMS1:104:C4KJNACXX:1:1212:1749:88430 99 chr10 3100168 29 101M = 3100238 171 GTGGAGAAAGGGGAACACTCCTCCATTGTTGGTGGGATTGCAAGCTTGTACAACCACTCTGGAAATCAGTCTGGCGGTTCCTCAGAAAATTGGACATAGTA CCCFFFFFHGGHHJJIJJJJIJJJJJJJIJJJGGIIGIIJIIGIJJJJGHIJJJJJJJJJJHHHHHFFFFFFEEEBD@BDDDDDDDDDDDDDDDDDDDDEE X0:i:452 MD:Z:101 XG:i:0 AM:i:0 NM:i:0 SM:i:0 XM:i:0 XO:i:0 XT:A:R

[jester112358]

I run a bash script with samtools merge -h rg.txt -r

below is the list of my rg.txt

@RG ID:Kapo93+94+95_1_09-1107_091207_HWI-EAS418_9_s_7 DS:091207_HWI-EAS418_9 SM:09-1107
@RG ID:Kapo93+94+95_1_09-1107_100618_ILLUMINA-8C38E9_0112_s_3 DS:100618_ILLUMINA-8C38E9_0112 SM:09-1107
@RG ID:Kapo93+94+95_1_09-1107_110103_HWUSI-EAS1785_0223_s_2 DS:110103_HWUSI-EAS1785_0223 SM:09-1107
@RG ID:Kapo93+94+95_1_09-1107_110616_SN588_0054_AC00T5ABXX_s_4 DS:110616_SN588_0054_AC00T5ABXX SM:09-1107
@RG ID:Kapo93+94+95_1_09-1107_111028_SN653_0108_BB0418ABXX_s_3 DS:111028_SN653_0108_BB0418ABXX SM:09-1107
@RG ID:Kapo93+94+95_1_09-1107_120216_SN670_0098_AC04HRABXX_s_1 DS:120216_SN670_0098_AC04HRABXX SM:09-1107
@RG ID:Kapo93-95-1_09-1107_110714_SN588_0055_AB0A5UABXX_s_1 DS:110714_SN588_0055_AB0A5UABXX SM:09-1107

I get this error and can't merge the bams:

##### ERROR MESSAGE: SAM/BAM file SAMFileReader{/xxxxx/xxxx/xxx/xxxxx/Kaposi/Kapo93+94+95/09-1107/Kapo_93-95_09-1107_merged_s_99.nodup.bam} is malformed: Read ILLUMINA-8C38E9_0112:3:84:1445:15512#0 is either missing the read group or its read group is not defined in the BAM header, both of which are required by the GATK. Please use http://gatkforums.broadinstitute.org...lacereadgroups to fix this problem
##### ERROR

I have tried pretty much everything that I can at this point and it would be really nice to get this work. Any help is appreciated.

[Clare S]

I think, tough this is not with any authority, that

ID = id name for the readgroup
SM = sample name
LB = label? dunno about this one
PL = platform

These are not currently standardized (I think) but ARE used by the Broad GATK, which means getting them right may be important for your pipeline

LB is library.

For many tools the really important field is ID, which must be unique to the read group. Reads are considered to be from different experimental conditions if they have different read groups. So we usually form our readgroup IDs by concatenating all the relevant information that uniquely identifies the experimental conditions (run/flowcell, lane, etc).