Error remove 0x1 from alternate contigs tutorial
Answered
Hi!
Hi! I am trying to follow this tutorial: https://gatk.broadinstitute.org/hc/en-us/articles/360037498992--How-to-Map-reads-to-a-reference-with-alternate-contigs-like-GRCH38 and am stuck at the step to Remove the 0x1 bitwise flag from alignments.
a) GATK version used: gatk-4.1.9.0
b) Exact command used:
samtools view -h /families/bam/CS7.bam | gawk '{printf "%s\t", $1; if(and($2,0x1)) {t=$2-0x1}else{t=$2}; printf "%s\t" , t; for (i=3; iprintf "%s\n",$NF}'| samtools view -Sb - > /families/bam/CS7.rm0x1.sam
The error:
gawk: {printf "%s\t", $1; if(and($2,0x1)) {t=$2-0x1}else{t=$2}; printf "%s\t" , t; for (i=3; iprintf "%s\n",$NF}
gawk: ^ syntax error
[main_samview] fail to read the header from "-".
I also found the tutorial from the GATK3 version on Github (https://github.com/broadinstitute/gatk-docs/blob/master/gatk3-tutorials/(How_to)_Map_reads_to_a_reference_with_alternate_contigs_like_GRCh38.md) which had a slightly different command:
samtools view -h /families/bam/CS7.bam | gawk '{printf "%s\t", $1; if(and($2,0x1)){t=$2-0x1}else{t=$2}; printf "%s\t" , t; for (i=3; i<NF; i++){printf "%s\t", $i} ; printf "%s\n",$NF}'| samtools view -Sb - > /families/bam/CS7.rm0x1.sam
But still got an error:
[E::sam_hrecs_error] Malformed key:value pair at line 2: "@SQ SN:chr1 LN:248956422 AS:38 M5:6aef897c3d6ff0c78aff06ac189178dd UR:/seq/references/Homo_sapiens_assembly38/v0/Homo_sapiens_assembly38.fasta SP:Homo sapiens"
samtools view: failed to add PG line to the header
I am not sure why I got these errors and am not sure how to get past these errors to move onto the HaplotypeCaller step.
Thank you in advance!
-
Hi Linda Do, which samtools version are you using? It looks like the version used for this tutorial is v1.3.1.
-
I am currently using v1.11. Do I need to downgrade to v1.3.1 for this tutorial?
-
Hi Linda Do,
Since this is a documentation question, I am going to move it to the Documentation Questions topic.
We don't support samtools so I don't have any insight into getting around the error, and also that tutorial has not been updated so the commands are out of date.
To get it to work though, it probably would be good to try the older version of samtools. Other community members might have insight as well.
Best,
Genevieve
-
Hi Genevieve Brandt. Unfortunately, after downgrading samtools to v1.3.1, the code using the tutorial ( https://gatk.broadinstitute.org/hc/en-us/articles/360037498992--How-to-Map-reads-to-a-reference-with-alternate-contigs-like-GRCH38 ) still did not work.
samtools view -h /families/bam/CS6.bam | gawk '{printf "%s\t", $1; if(and($2,0x1)){t=$2-0x1}else{t=$2}; printf "%s\t" , t; for (i=3; iprintf "%s\n",$NF}'| samtools view -Sb - > /vipbg-data/home/shianglab/20180828_NG/VariantFiltration/families/bam/CS6.rm0x1.bamThe error I recieved:
gawk: {printf "%s\t", $1; if(and($2,0x1)){t=$2-0x1}else{t=$2}; printf "%s\t" , t; for (i=3; iprintf "%s\n",$NF}
gawk: ^ syntax errorBut when I run (from the GitHub GATK3 tutorial)
samtools view -h /families/bam/CS6.bam | gawk '{printf "%s\t", $1; if(and($2,0x1)){t=$2-0x1}else{t=$2}; printf "%s\t" , t; for (i=3; i<NF; i++){printf "%s\t", $i} ; printf "%s\n",$NF}'| samtools view -Sb - > /vipbg-data/home/shianglab/20180828_NG/VariantFiltration/families/bam/CS6.rm0x1.bamIt does output a file but I am still unsure if that file will work with the rest of the workflow (I will update this post later if it works or not) since it is from the GitHub GATK3 tutorial and not the same as in the tutorial on the this (GATK4) website. The only difference between the two tutorials in the gawk part so does this mean that the code in the tutorial is wrong and that the gawk code from the GATK3 GitHub is correct?
Thank you in advance!
-
Linda Do It looks like that might be a copy-paste typo, thank you for finding that! I'll try to get the GATK4 website version fixed.
Please sign in to leave a comment.
5 comments