Samtools gives "truncated file" when trying to view a "reheaded" file

Question

[olavur]

I have a BAM file that works fine, but after using samtools reheader, I cannot view it in samtools. The new file has the same size as the original.

I run the command: samtools reheader -P header.sam input.bam > output.bam

Where header.sam contains:

@HD VN:1.5 SO:coordinate
@SQ SN:10X_hg19_ucsc_v2.1.0 LN:3199910119
@CO Whole-exome sequencing of one individual.

input.bam is a BAM file with no existing header, which I can view as samtools view input.bam | less with no problem.

The reheader command finishes without giving any error. Then I try to view the file, using samtools view output.bam | less and

[main_samview] truncated file.

If I use the -h option (samtools view -h output.bam | less) I get:

[main_samview] truncated file.
@HD VN:1.5 SO:coordinate
@SQ SN:10X_hg19_ucsc_v2.1.0 LN:3199910119
@CO Whole-exome sequencing of one individual.

Answer 1

[Devon Ryan]

Reheader the file again, apparently there was an error the first time. The error means that the "magic number" at the end of all BAM files isn't there. That could indicate that there are alignments missing at the end of the reheadered file as well.

[olavur]

Have run it several times, no difference.

[Devon Ryan]

You must be running out of disk space then.

[olavur]

I have plenty of disk space, 100 TB available.

[Devon Ryan]

Then you have some sort of hardware issue. This isn't a reproducible problem by others.

[Panagiotis Moulos]

I can confirm this. I face a similar problem with samtools 1.5.

[Devon Ryan]

One of you will need to post a file that reproduces this issue somewhere.

[olavur]

Panagiotis, are you able to do this? The data I'm working on is confidential.

[Panagiotis Moulos]

Maybe this will help. https://github.com/pmoulos/samtools-reheader-fail

[Devon Ryan]

BTW, if you need to remove chromosomes that aren't at the end of the header then you need to either use cat (cat header <(samtools view file_with_alignment_removed.bam) | samtools view -bo reheadered.bam -) or a different tool (maybe something from Picard).

[Devon Ryan]

Your reheader command removed a chromosome that wasn't the last one, so you ended up with orphaned reads. That samtools sort then threw an error is more of a feature than a bug, since you wouldn't have known that you'd lost chrX (and reassigned its alignments in an illegal manner to chrY) otherwise.

[Devon Ryan]

I'm giving this a try now.

[Devon Ryan]

I can reproduce this example, I'll track down the issue.