Mapping by bwa error

[asmmahmoud]

Hi Lijia Yu,
When I was running sveditor.py script with an inversion variant as an input, the job was killed when the step of mapping by bwa starts

File "/home/asmmahmoud/ctg/varben/VarBen-master/bin/../varben/common/bamconvert.py", line 14, in _call
raise Exception("Cmd Error: %s" % cmd)
Exception: Cmd Error: bwa mem -R "@rg\tID:A4-NA24143-WEST_L1\tLB:RGLB\tPU:1\tSM:A4-NA24143-WEST" -t 1 /home/asmmahmoud/ctg/ref/GRCh38_full_analysis_set_plus_decoy_hla.fa /home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.to_1.fq /home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.to_2.fq >/home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.sam

Could you help me fix it ?

[yulijia]

Hi,
Can you show me if you run this by command line

bwa mem -R "@rg\tID:A4-NA24143-WEST_L1\tLB:RGLB\tPU:1\tSM:A4-NA24143-WEST" -t 1 /home/asmmahmoud/ctg/ref/GRCh38_full_analysis_set_plus_decoy_hla.fa /home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.to_1.fq /home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.to_2.fq >/home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.sam

What's error it will return?

[asmmahmoud]

Hi,
this is error msg: [E::bwa_set_rg] the read group line is not started with @rg

[yulijia]

I think the problem is the read group tag should be @RG but not @rg.

can you run the line below to see if it works?

bwa mem -R "@RG\tID:A4-NA24143-WEST_L1\tLB:RGLB\tPU:1\tSM:A4-NA24143-WEST" -t 1 /home/asmmahmoud/ctg/ref/GRCh38_full_analysis_set_plus_decoy_hla.fa /home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.to_1.fq /home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.to_2.fq >/home/asmmahmoud/ctg/varben/inv_test/tempDir/edit.remap.sam

[asmmahmoud]

It gave this output:
[M::bwa_idx_load_from_disk] read 0 ALT contigs
[M::process] read 2 sequences (0 bp)...
[M::mem_pestat] # candidate unique pairs for (FF, FR, RF, RR): (0, 0, 0, 0)
[M::mem_pestat] skip orientation FF as there are not enough pairs
[M::mem_pestat] skip orientation FR as there are not enough pairs
[M::mem_pestat] skip orientation RF as there are not enough pairs
[M::mem_pestat] skip orientation RR as there are not enough pairs
[mem_sam_pe] paired reads have different names: "/1", "/2"
Could you help me to run the script without changing the read group tag ?
I used same bam file with another variant list, it was completed successfully without this error at mapping by bwa step
bwa mem -R "@rg\tID:A4-NA24143-WEST_L1\tLB:RGLB\tPU:1\tSM:A4-NA24143-WEST" -t 1 /home/asmmahmoud/ctg/ref/GRCh38_full_analysis_set_plus_decoy_hla.fa /home/asmmahmoud/ctg/varben/A4/Hom1_mut/tempDir/edit.remap.to_1.fq /home/asmmahmoud/ctg/varben/A4/Hom1_mut/tempDir/edit.remap.to_2.fq >/home/asmmahmoud/ctg/varben/A4/Hom1_mut/tempDir/edit.remap.sam

[yulijia]

I think the problem is similar to #23 now.
I will check with the demo in a few weeks.
But I cannot help you generate the edited data in a short time.
If you still need to get the edited data recently, can you try using bamsurgeon?
https://github.com/adamewing/bamsurgeon

[asmmahmoud]

Hi Lijia Yu,
Thank you for your reply. I can't use bamsurgeon because I need to be consistent with other edited bam files generated by VarBen. I will wait for your update :)

[yulijia]

I have tried to edit a CNV in the example data EGFR_sort.bam, and I didn't encounter any error. My bwa version is 0.7.17-r1198-dirty.

chr7    55248067        55249038        cnv     2.5     gain

I am not able to help you check what's the problem in your bam file.