CRAM compression and TLEN SAM's field

[priesgo]

Hi there,

I'm evaluating the CRAM compression. So, I'm comparing the same BAM sample pre and post CRAM compression.

I found that in my data the TLEN (Template Length) field is sistematically modified summing 1 to this value.

For instance an original read like:

SRR107049.155163702     163     1       69171   44      76M     =       69404   [B][COLOR="Red"]308[/COLOR][/B]     CTATGGAGGAATCGTGTTTGGAAACCTTCTTATTGTCATAACAGTGGTATCTGACTCCCACCTTCACTCTCCCATG    S?;<CC>D=?@=C<>D?A@F?@CC
D3A?<C?>=>G?CA=:AE>E@HC><??CDAEABD??:9=<AAFCBCABE<>S    RG:Z:SRR107049  NM:i:0  OQ:Z:GGG@AFDE@E@BC@CBFGFFAFEDD4>B9DBD??EDBC?;AD<ADEC>?C=DBDE@AD>A;:;=>EGEBE>BE<A9

After compression/decompression looks like:

2       163     1       69171   44      76M     =       69404   [B][COLOR="Red"]309[/COLOR][/B]     CTATGGAGGAATCGTGTTTGGAAACCTTCTTATTGTCATAACAGTGGTATCTGACTCCCACCTTCACTCTCCCATG    S?;<CC>D=?@=C<>D?A@F?@CCD3A?<C?>=>G?CA=:AE>E@HC><??CDAEABD??:9=<AAFCBCABE<>S    RG:Z:SRR107049

Anybody knows something about this issue? I have gone through some posts that talk about differences among the aligners in the TLEN calculation.

Other changes that are OK for me are the read name which are replaced by numbers and tags which are removed, except the read group tag.


Thanks!
Pablo.

[maubp]

An issue like that was reported on the cram mailing list last year, actually a Picard bug:

https://github.com/vadimzalunin/crammer/issues/20

What version of the cram-tools do you have?

[BAMseek]

The magnitude of TLEN is "the number of bases from the leftmost mapped base to the rightmost mapped base".

So you should be able to check the other read and calculate TLEN to see which is right.

Assuming that the other read also has 76 matches, then I would think

TLEN = 69404 + 76 - 69171 = 309.

Which aligner did you use for the original BAM file?

Justin

[priesgo]

Thanks for your answers!
maubp, I am using version 1.00-b244 and it is linked to Picard 1.79.
Justin, this data is an Illumina exome from the 1000 Genomes Project. As it is reported in the project information, it was aligned with Mosaik (don't know the version).

By the way this is the mate-pair:

SRR107049.155163702     83      1       69404   39      76M     =       69171   -308    TGGTGACCCCCATAGCCATGGGCTGTGACAGATATAGAGCAATATGCAAGCCCCTACACTACACTACAATTATGTG    ##############################################EABFCC@G?>E=E;<H>F>>F@>B=>>B>S    RG:Z:SRR107049  NM:i:4  OQ:Z:##############################################?CBECCAEC=>;A<:FBDCBDEDDDDADDD

The value of 309 also makes sense for me. By the way, Justin, I saw you assumed that the mate-pair had a length of 76 just as his mate. And this is correct as you can see above, do mates usually have the same length?


Pablo.

[BAMseek]

Hi Pablo,
Yeah, it looks like Picard might be giving the right answer on this. Not sure which downstream tools rely on TLEN - the only thing I can think of is possibly the genomic viewers.

Paired-end reads won't necessarily be the same length (for example, if there was some trimming done to the reads to remove adapters). And even if they were the same length, they may not have the same number of aligned bases. To calculate the length of an alignment, I think it is #M (match/mismatch) + #D (deletions) + #N (skips). An alignment can have soft clips (bases at the ends of the read that aren't aligned), and those wouldn't be counted towards the length of the alignment. I just got lucky that it also happened to be 76M.

best,
Justin

[priesgo]

Thanks again Justin,

Yes, it seems to be an issue with the original file, don't know if coming from mosaik or somewhere else.

Pablo.