Question: loop for samtools to deal with multiple files
0
gravatar for hfang4
19 months ago by
hfang40
hfang40 wrote:

Hello,

I have a script to run samtools pileup for multiple files. The script worked if I ran one file at a time without loop, but it gave me error saying "Illegal variable name" with the loop. So, I think it is something wrong with the loop. I have to use the old version of samtools with the pileup for the downstream analysis. Please anyone help me out with the issue.

 #! /bin/csh
cd /share/gwascotton/GBS_data/pileup_trial/
source /usr/local/apps/samtools/samtools017a.csh
for i in $(ls *.bam);do
samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/"$i" > /share/gwascotton/GBS_data/new/"$i".pileup
done ;

Thank you.

fhzh

snp software error • 1.7k views
ADD COMMENTlink modified 19 months ago • written 19 months ago by hfang40
1
gravatar for Kevin Blighe
19 months ago by
Kevin Blighe41k
Kevin Blighe41k wrote:

First get a listing of the BAM files and save in BAM.list, then:

#!/bin/bash
paste BAM.list | while read BAM ;
do
    samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/"${BAM}" > /share/gwascotton/GBS_data/new/"${BAM}".pileup
done ;
ADD COMMENTlink written 19 months ago by Kevin Blighe41k

Thank you for the code. I think I have to use csh because I have to submit my job to high performance computing (HPC) server and using the old version of samtools which the source code give the path to it. I modified your code a little bit as attached below, but I got error "while: Expression Syntax". Any ideas? Thanks

 #!/usr/bin/csh

cd /share/gwascotton/GBS_data/pileup_trial/
source /usr/local/apps/samtools/samtools017a.csh

ls *.bam > /share/gwascotton/GBS_data/pileup_trial/BAM.list ;
paste BAM.list | while read BAM ; 
do
    samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0. "${BAM}" > /share/gwascotton/GBS_data/new/"${BAM}".pileup
done;
ADD REPLYlink modified 19 months ago by genomax65k • written 19 months ago by hfang40

Every modern computing environment will support bash. So you should be able to use it. You will need to figure out what /usr/local/apps/samtools/samtools017a.csh is doing and replace that with equivalent bash commands/script.

ADD REPLYlink written 19 months ago by genomax65k

Oh dear, I noticed the old version of SAMtools as judged by the fact that you are using pileup, which has been replaced with mpileup.

If you really need to use the c-shell, then try:

#!/bin/csh
foreach BAM (`cat BAM.list`) ;
    samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/"${BAM}" > /share/gwascotton/GBS_data/new/"${BAM}".pileup
end ;

NB - cat BAM.list is surrounded by backticks (not single quotes).

ADD REPLYlink modified 19 months ago • written 19 months ago by Kevin Blighe41k
0
gravatar for dariober
19 months ago by
dariober10.0k
WCIP | Glasgow | UK
dariober10.0k wrote:

I think the csh shell doesn't support that for loop syntax (see for example https://www.cyberciti.biz/faq/csh-shell-scripting-loop-example/ ). I would use Bash unless you have to use csh.

ADD COMMENTlink modified 19 months ago by genomax65k • written 19 months ago by dariober10.0k

Thanks for the comments, but the link you provided can't be opened.

I tried bath script. It worked somehow because it produced .pileup files, but I also got error messages “/usr/local/apps/samtools/samtools017a.csh: line2: setenv: command not found” and “/usr/local/apps/samtools/samtools017a.csh: line8: syntax error: unexpected end of file”.

I think that because the bath script did not apply the samtools017a.csh which I wanted.

Anyway, thanks for your insight.

ADD REPLYlink modified 19 months ago • written 19 months ago by hfang40

Link in @Dario's post fixed.

ADD REPLYlink written 19 months ago by genomax65k

That's correct, but the syntax for looping in c-shell is:

foreach VAR (`cat file.list`) ;
    echo "do something" ;
end ;

See my solution above.

ADD REPLYlink modified 19 months ago • written 19 months ago by Kevin Blighe41k
0
gravatar for arfesta
19 months ago by
arfesta30
arfesta30 wrote:

In terminal why not just try:

  cd /share/gwascotton/GBS_data/pileup_trial/
  ls *.bam > ~/Desktop/list.of.bams
  while read file; do samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f 
      /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/$file >
      /share/gwascotton/GBS_data/new/$file.pileup; done < ~/Desktop/list.of.bams
ADD COMMENTlink written 19 months ago by arfesta30

Thanks for the idea. I run the following codes but got another error "Ambiguous output redirect". I think we are almost there. But any idea on the new error?

 #! /bin/csh

cd /share/gwascotton/GBS_data/pileup_trial/
source /usr/local/apps/samtools/samtools017a.csh
ls *.bam > /share/gwascotton/GBS_data/pileup_trial/list.of.bams |
while read file; do
samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/$file > /share/gwascotton/GBS_data/new/$file.pileup;
done < list.of.bams;
ADD REPLYlink modified 19 months ago by genomax65k • written 19 months ago by hfang40

See if putting this line in " helps in script above. "samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/$file > /share/gwascotton/GBS_data/new/$file.pileup"

ADD REPLYlink modified 19 months ago • written 19 months ago by genomax65k
0
gravatar for Kevin Blighe
19 months ago by
Kevin Blighe41k
Kevin Blighe41k wrote:

As per my comment from my thread above, here are the solutions for both bash and csh:

First get a listing of the BAM files and save in BAM.list, then:

BASH

#!/bin/bash
paste BAM.list | while read BAM ;
do
        samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/"${BAM}" > /share/gwascotton/GBS_data/new/"${BAM}".pileup
done ;

c-shell

#!/bin/csh
foreach BAM (`cat BAM.list`) ;
    samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/"${BAM}" > /share/gwascotton/GBS_data/new/"${BAM}".pileup
end ;

NB - cat BAM.list is surrounded by backticks (not single quotes).

ADD COMMENTlink written 19 months ago by Kevin Blighe41k

Please don't post duplicate answers. Instead edit your original answer with this material.

ADD REPLYlink written 19 months ago by genomax65k

Thanks for the scripts. I tried several files, but neither of them works. The bash one keeps running for ever and the csh one gave me the same error message "Ambiguous output redirect".

ADD REPLYlink written 19 months ago by hfang40

It works for me in both BASH and CSH, and even with SAMtools and a redirect. Just to be sure, cat BAM.list is enclosed by backticks, not quotes/single quotes.

This is looking like a local issue with your system architecture / program installation. What are the contents of /usr/local/apps/samtools/samtools017a.csh?

ADD REPLYlink written 19 months ago by Kevin Blighe41k
0
gravatar for hfang4
19 months ago by
hfang40
hfang40 wrote:

This final version of the script worked. Thank all for your suggestions, comments, codes and insight.

#! /bin/bash

cd /share/gwascotton/GBS_data/pileup_trial/
source /usr/local/apps/samtools/samtools017a.bash for i in $(ls *.bam);
do

samtools pileup -c -l /share/gwascotton/GBS_data/merged_files.sorted.bam.pileup -f /share/gwascotton/GBS_data/Ghirsutum_458_v1.0.fa /share/gwascotton/GBS_data/pileup_trial/${i} > /share/gwascotton/GBS_data/new/${i}.pileup

done ;

ADD COMMENTlink modified 19 months ago • written 19 months ago by hfang40
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: 1199 users visited in the last hour