[KevinLam]

Thanks Li Heng. I just want to confirm that nothing is wrong with the sorted bam file.

Thanks a lot.

You may convert to sam and count the lines as a proxy for checking.
but yes as lh3 mentioned sorted files compress better due to the compression algorithm

[KevinLam]

There was a discussion on how the CS tag should be generated in sam files according to the specs. Is there a consensus on how it is to be done?

I have to write a script to append the CS tag info to the BWA alignment of SOLID reads. I am hoping to make it as painless as possible.

[bosTau2]

Is there sam2maq (sam to maq map format)? I am testing simulated data and like to convert ssaha2 alignment to maq output so that I can analyze them in maq. I cerated all the simu data in maq.
Thank you.

[wuhoucdc]

Dear All,

Do you know if samtools can perform multiple commands (>2) together? Here assuming that I have ten BAM files (result001.bam, result002.bam,......result010.bam) and want to merge them first and then sort and index them, the last step I hope is to extract the data for chromosome 1 (chr1), how can I edit the samtools command? I did it like this:
samtools merge result.bam result001.bam result002.bam ............result010.bam |\
samtools sort - result | \
samtools index result.bam | \
samtools view result.bam chr1 > resultchr1.bam

Is it right?

Thank you very much!

Wu

[sbaheti]

Hi lh3

Originally Posted by xguo
I got a list of candidate SNPs using BWA and samtools for RNASeq data, and am trying to weigh various filtering options. "samtools.pl varFilter" gives a list of filtering criterion with default setting. The maximum read depth is set at 100. Given that duplicate reads have been removed by "samtools rmdup", do I still need to limit the maximum read depth?

replied by lh3:
Yes, you need this unless you are doing target sequencing in which case the read depth is expected to vary a lot.

My question

If we are doing variant calling in exome capture analysis, do i need to limit the max read depth when using samtool varFilter tool?

[bosTau2]

It depends on regions. If a region is repetitive then you need to filter out possible duplications which can cause artificial hetero SNPs. "Repetitive" here means kmer uniqueness --- how many times a given kmer (30mer for example) can be found in a reference.

[aby]

I was trying to convert my single read illumina file to SAM format using the export2sam.pl script. It does not work and gives the following error:

Use of uninitialized value $t[21] in string ne at export2sam.pl line 67, <$fh1> line 14063.

Here is a sample data from my file.

HWUSI-EAS174_0025:2:1:5:488#0/1:CGGAGAATACGCTCCCATTCCCCCNGNANNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNTNTGATCTTAGATCGGA:aabbbbb_baaa_a]ab`_aBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB
HWUSI-EAS174_0025:2:1:5:1542#0/1:TGGATGCCTAGGCAATCAGAGGCGNANANNNNNNNNNNNNNNNNNNNNNNNNNNNNNNNANGTGATAAGCAGCGAA:abbbabbbbbbbbbbbbbbbBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBBB


How do I resolve it ? Is there any other way?

A.

[aby]

I converted the all.map file generated from MAQ tool to SAM format using maq2sam-long in samtools. The output file I named as Out. Now when I try to convert Out file (which is in SAM format) to BAM format using the command
samtools view -S Out

I get an error message that header @ not found. When I edit the Out file by manually adding @ at the top, more errors appear. What should I do? Is there a fault in the SAM file generated by the converter maq2sam ? Is there any other way to convert SAM file to BAM file? I tried 'samtool import' and that also does not work.

[NF_seq]

Question about the output of bwa?

I got the output, see below:
HWI-EAS307:1:54:758:902#0 20 19641_CLSZ1904.b1_P20.ab1_CLSZ_L._sativa_library_forward_335 301 20 36M * 0 0 CAAATCGGTGTGTTTTCACTGGTCGTGCTCGTTCCG aabaaaaaaaaababaa`aaaabaabaabbabaaaa XT:A:U NM:i:1 X0:i:1 X1:i:2 XM:i:1 XO:i:0 XG:i:0 MD:Z:35T0 XA:Z:13134_QGB27J17.yg.ab1_QGB_L._sativa_library_forward_448,-58,36M,2;7061_CLS_S3_Contig6993_CLS_S3_L._sativa_library_forward_968,-404,36M,2;

I can not understand the flag value 20. I used "samse" to process single reads.
"XT:A:U" indicates the read uniquely mapped to the reference, why i still got XA for alternative alignment inforamtion?
It is confused me. Someone could help me a bit for that? Thank you very much

It's also confused me that "XT:A:U" and "XA:..." information came up in one alignment. Could anyone please explain that? Thanks a lot! And if i only care about uniquely mapped reads, is this kind of reads what i want?

[aby]

sam to bam conversion not taking place

samtools import example.sam example.bam

This command does not work. Gives error:

Usage: bamtk import <in.ref_list> <in.sam> <out.bam>


What to do? What is in.ref_list ?

[AXW]

Does anybody know of a utility/script for converting .SAM/BAM files into .SOAP? I know that there are scripts out there to go from SOAP->SAM, but I can't find anything going the other way.

Cheers.

[swbarnes2]

samtools import example.sam example.bam

This command does not work. Gives error:

Usage: bamtk import <in.ref_list> <in.sam> <out.bam>


What to do? What is in.ref_list ?

Use faidx to make a .fai file, use that for the in.ref.list. It works for me.

[seq_GA]

Hi Heng,

Can you please explain about the new feature of samtools (ie) multisample pile up? Thanks.

[lh3]

http://samtools.sourceforge.net/mpileup.shtml

[aby]

Okay, I have solved my problems. Seems there is a script to add the header file, and different command options for conversion to Bam.