On Mon, 16 Apr 2018, RITSUKO OIKAWA wrote:

Samtools does have an option to *filter* the reads according to regions
specified in the BED format.

$ samtools view -hL chr22_21804028_ref.Gene.bed PLCPRF5bis44.sam >

[E::sam_parse1] missing SAM header
[W::sam_read1] Parse error at line 2
[main_samview] truncated file.

But as above
I get an error message.
Is the command incorrect ?
Is sam format different?
We look forward to your reply.

It looks like your input sam file doesn't have any `@SQ` lines to describe the reference sequences that you're using. They should be at the start of the file, before any alignment records. What did you use to make the sam file?

If you have a fasta file with your reference sequences in it, you can use `samtools dict` to make suitable `@SQ` records. You can prepend the resulting headers to your original file with a bit of shell scripting. For example this should add in the @SQ lines while preserving any other header records that you already have:

head -n 1 PLCPRF5bis44.sam | grep '^@HD' > new.sam
samtools dict -H ref.fa >> new.sam
grep -Ev '^@(HD|SQ)' PLCPRF5bis44.sam >> new.sam

Rob Davies              r...@sanger.ac.uk
The Sanger Institute    http://www.sanger.ac.uk/
Hinxton, Cambs.,        Tel. +44 (1223) 834244
CB10 1SA, U.K.          Fax. +44 (1223) 494919

The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE.
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
Samtools-help mailing list

Reply via email to