pysam alignment error

Question

[rthapa]

Hi,

I am trying to read bam file using a command: samfile = pysam.AlignmentFile(os.path.join(path,bam),'rb') but I keep on getting an error:

Traceback (most recent call last):
  File "<stdin>", line 3, in <module>
  File "pysam/libcalignmentfile.pyx", line 741, in pysam.libcalignmentfile.AlignmentFile.__cinit__
  File "pysam/libcalignmentfile.pyx", line 990, in pysam.libcalignmentfile.AlignmentFile._open
ValueError: file has no sequences defined (mode='rb') - is it SAM/BAM format? Consider opening with check_sq=False

I checked my bam files and they look fine. All of my bam files have headers too. Could anyone provide any suggestion?

Thanks

[WouterDeCoster]

I checked my bam files and they look fine.

How did you check and how did you conclude the files are fine?

[rthapa]

I used samtools flagstat to see if the alignment is okay and I got the output. It didn't throw me any error,

25906192 + 0 in total (QC-passed reads + QC-failed reads)
0 + 0 secondary
0 + 0 supplementary
0 + 0 duplicates
20865552 + 0 mapped (80.54% : N/A)
25906192 + 0 paired in sequencing
12953096 + 0 read1
12953096 + 0 read2
20348218 + 0 properly paired (78.55% : N/A)
20501230 + 0 with itself and mate mapped
364322 + 0 singletons (1.41% : N/A)
152966 + 0 with mate mapped to a different chr
152966 + 0 with mate mapped to a different chr (mapQ>=5)

[finswimmer]

Hello,

can you use samtools to have a look into the the bam file?

$ samtools view -h input.bam

Or does this give you an error as well?

fin swimmer

[rthapa]

When I used samtools view -h input.bam - I get the following output,

@SQ SN:super_3314   LN:1110
@SQ SN:super_3316   LN:1100
@SQ SN:super_3319   LN:1072
@SQ SN:super_3320   LN:1071
@SQ SN:super_3322   LN:1061
@SQ SN:super_3323   LN:1039
@SQ SN:super_3326   LN:1005
@SQ SN:Chr01    LN:80884392
@SQ SN:Chr02    LN:77742459
@SQ SN:Chr03    LN:74386277
@SQ SN:Chr04    LN:68658214
@SQ SN:Chr05    LN:71854669
@SQ SN:Chr06    LN:61277060
@SQ SN:Chr07    LN:65505356
@SQ SN:Chr08    LN:62686529
@SQ SN:Chr09    LN:59416394

[finswimmer]

Ok,

just another guess: Are you sure you have a bam and not a sam file? What's the output of head input.bam?

fin swimmmer

[rthapa]

The output of head .bam looks like this,

$ head SRR5748777.bam 
?BC?З/hBAM?X@HD VN:1.0  SO:unsorted
@PG ID:GSNAP    PN:gsnap    VN:2018-03-25   CL:gsnap.sse42 -t 4 -m 0.02 -B 5 -n 1 -Q --nofails -d data -D /work/rthapa/GO/gsnap_indexes -A sam SRR5748777_1.fastq SRR5748777_2.fastq
@SQ SN:Chr10    LN:61233695
@SQ SN:super_16 LN:1391575
@SQ SN:super_18 LN:1172180
@SQ SN:super_20 LN:888508
@SQ SN:super_22 LN:639901
@SQ SN:super_25 LN:641756
@SQ SN:super_26 LN:650379
@SQ SN:super_27 LN:526259

[finswimmer]

This is not a bam file. It is a sam file.

pysam.AlignmentFile(os.path.join(path,bam),'rb')

Just remove the 'rb' part. Because this opens the file in binary mode, which you don't have.

[rthapa]

I used samtools view -h input.sam > output.bam to convert the sam file to bam file. I do not know why it didn't convert to bam file. As you suggested I opened the file removing 'rb',

pysam.AlignmentFile(os.path.join(path,bam))

But still getting the following error for

samfile.count(mychr, start, stop)

Traceback (most recent call last):
  File "<stdin>", line 39, in <module>
  File "pysam/libcalignmentfile.pyx", line 1431, in pysam.libcalignmentfile.AlignmentFile.count
  File "pysam/libcalignmentfile.pyx", line 1112, in pysam.libcalignmentfile.AlignmentFile.fetch
ValueError: fetching by region is not available for SAM files

[finswimmer]

The fetch() method doesn't allow random access. Convert to bam and coordinate sort your sam file:

$ samtools sort -o out.bam input.bam

---

samtools view -h input.sam > output.bam

The standard output by samtools is sam. So this command doesn't convert to bam. You have to use the parameter -Ob to get this. Or instead of redirecting to a file (> output.bam), use the -o parameter to define the output file. samtools recognizes the file extension, and will convert.

[rthapa]

thanks a lot, it finally worked.

[rthapa]

I am not sure if updating the newest version of pysam solved this issue. After updating the version, it is not throwing me the error. But I am getting other error while counting in the sam file, samfile.count(mychr, start, stop).

Traceback (most recent call last):
  File "<stdin>", line 39, in <module>
  File "pysam/libcalignmentfile.pyx", line 1431, in pysam.libcalignmentfile.AlignmentFile.count
  File "pysam/libcalignmentfile.pyx", line 1112, in pysam.libcalignmentfile.AlignmentFile.fetch
ValueError: fetching by region is not available for SAM files

Answer 1

[skbrimer]

I have not used pysam before but have you tried the suggestion the code gives samfile = pysam.AlignmentFile(os.path.join(path,bam), check_sq=False, 'rb') ? also are you using Windows to process the file the rb call is for windows, linux/mac would be just r

Answer 2

[Jianyu]

A similar error on github: https://github.com/pysam-developers/pysam/issues/442

Maybe updating the version of your pysam can fix the problem?