Running pe-utils to find insert length for miso failing

Question

[spacer]

I am running pe-utils to determine the insert length for Miso. However, firstly I am getting messages like the following:

WARNING: NB501137:15:HHYT3BGXY:1:23310:15215:17414 left mate starts later than right mate WARNING: NB501137:15:HHYT3BGXY:4:11608:2223:10635 left mate starts later than right mate WARNING: NB501137:15:HHYT3BGXY:1:21306:5350:10090 left mate starts later than right mate WARNING: NB501137:15:HHYT3BGXY:1:12103:10191:8692 left mate starts later than right mate

The run then ends with the following error and the insert_len file is empty: Error: Could not find any properly mated pairs to compute insert length with. Are you sure your BAM reads are properly paired and map the chromosome headers in the constitutive exon file?

I don't really know what all this mean. Why would the left mate start later than the right mate in the first place?

The fastq files were aligned with star and sorted by coordinate. I could not find anything in the Miso manual about how the files bam files should be sorted. I might have missed it.

If someone could perhaps help me solve the problem it would be greatly appreciated.

[rjoliver]

Does anyone have a fix for this? My chromosome headers were different before (chr1 vs 1) which I changed all to chr1 but then I get this same error message with the left mate/right mate issue. I have tried to run the miso protocol up to this step with a sorted and an unsorted file bam file. Any suggestions? I will gladly provide more information if required.

Answer 1

[Brian Bushnell]

You can calculate the insert size distribution from a sam or bam file with the BBMap package like this:

reformat.sh in=mapped.sam ihist=ihist.txt

This requires the latest version of BBMap (36.84 or higher). However, calculating the insert size distribution requires that your reads were mapped paired.

[spacer]

Hi Brian

Thanks for that. I'll give it a go. But I still would like to know why pe-utils won't work!?!

My reads were mapped as follows:

STAR --genomeDir ${GENDIR} --readFilesIn ${FASTQ_R1} ${FASTQ_R2} --readFilesCommand gunzip -c --runThreadN ${NSLOTS} --outSAMtype BAM Unsorted --twopassMode Basic

I then sorted by coordinate using samtools:

samtools sort ${FILENAME}.sorted.coord.bam -O bam ${FILENAME}.bam

Is there anything else I should have done for miso/pe-utils to work?

[spacer]

reformat with the ihist parameter does not work for bam files and looking at the help I can't find an alternative switch for use with bam files. Does someone perhaps have any idea why pe-utils fails or have a suggestion for calculating the insert size from bam files with something else?

[Brian Bushnell]

For the BBMap package to handle bam files, you need to have samtools installed and in your path.

I don't know why pe-utils doesn't work, though; I have never used that program. Perhaps it wants an unsorted bam file with the paired reads interleaved?

[GenoMax]

Perhaps PE utils expect your BAM file to be sorted (not an unreasonable expectation). Have you tried to sort your files with samtools and then used PE utils?

[spacer]

Like I mentioned and showed in my comment above to Brian, the files are sorted with samtools using:

samtools sort ${FILENAME}.sorted.coord.bam -O bam ${FILENAME}.bam

[GenoMax]

Sorry I missed that line in your post above. I have highlighted the code in that post.

Are you using the latest samtools (v.1.3.1)? Your command looks a bit strange since it would be samtools sort -o file_sorted.bam unsorted_file.bam with new samtools.