Question: mpileup for creating a consensus sequence
0
gravatar for phernet123
2.6 years ago by
phernet1230
phernet1230 wrote:

Dear colleagues

I desperately need your help with the samtools mpileup command. I am not variant calling, rather I want to create a consensus sequence of all the mapped reads in the Sorted.bam file. I am trying to run the following bash script:

#!/bin/bash
#PBS -l select=1:ncpus=24:mem=64GB
#PBS -P CBBI0911
#PBS -q smp
#PBS -l walltime=96:00:00
#PBS -o /mnt/lustre/users/felegbeleye/my_data/stdout.femi.txt
#PBS -e /mnt/lustre/users/felegbeleye/my_data/stderr.femi.txt
#PBS -N samtools
#PBS -M --removed--
module load chpc/BIOMODULES
module load samtools/0.1.19


EXE="samtools"

ARGS="mpileup -u -f 19Dec2016_S9zPn.fasta Sorted.bam |
bcftools call -c - |
vcfutils.pl vcf2fq -D 100 > femi.psmc.fq"

cd /mnt/lustre/users/felegbeleye/my_data/Ref_Genome

${EXE} ${ARGS}

I am getting the following error:

mpileup: invalid option -- 'c'
open: No such file or directory
[mpileup] failed to open |: No such file or directory

I do not have an option c in my mpileup command. My indexing of the reference was done in bwa and not in Samtools. Could that be it? I feel I am doing something basic wrong but I am new at this and would appreciate some guidance.

Thank you very much for your time and help.

Femi

ADD COMMENTlink modified 2.6 years ago by genomax83k • written 2.6 years ago by phernet1230
1

Try putting the entire command for ARGS on one line and see if that helps. (note: I redacted your email address from the script above).

ARGS="mpileup -u -f 19Dec2016_S9zPn.fasta Sorted.bam | bcftools call -c - | vcfutils.pl vcf2fq -D 100 > femi.psmc.fq"

ADD REPLYlink modified 2.6 years ago • written 2.6 years ago by genomax83k

For starters, you're using an old/deprecated version (v0.1.19) that did not include the BCFtools 'call' command (I believe it was 'view' back then). Update to the current version and try again.

ADD REPLYlink written 2.6 years ago by harold.smith.tarheel4.5k

Good catch. phernet123 : This may just be a matter of using right module load call. Check with module avail to see if newer versions of samtools (v.1.4+) are available.

ADD REPLYlink modified 2.6 years ago • written 2.6 years ago by genomax83k

My indexing of the reference was done in bwa and not in Samtools. Could that be it?

SAMtools and BWA will use different indices. If you haven't create the SAMtools indices, create it in addition to the ones for BWA with samtools faidx. Other programs use their own indices too.

I don't know if this is the exact problem. You should follow-up on what the other have already said.

ADD REPLYlink modified 2.6 years ago • written 2.6 years ago by Kevin Blighe60k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1420 users visited in the last hour