Use a variable instead of a file for bedtools getfasta BED input?

338 views Asked by At

I would like to input a variable in the place of a file in the command below, is this possible?

I would like to generate some FASTA files using the programme bedtools getfasta. Typically I would order the input files first using the awk to bedtools sort commands below. The output is then piped into an output file. This output file would then follow the -bed flag in bedtools getfasta so I can create the FASTA files I require, e.g.

awk '{if ($2>$3)print $1,$3,$2,".",".","-";else print $1,$2,$3,".",".","+";}' OFS='\t' Infile.bed | 
awk '{a=$2-1;print $1,a,$3,$4,$5,$6;}' OFS='\t' | 
bedtools sort > OrderedFile.bed


bedtools getfasta -s -fi Infile.fasta -bed OrderedFile.bed -fo Outfile.fasta

However, I have a lot of files I would like to use bedtools getfasta for. I was hoping to avoid creating the additional OrderedFile.bed files by setting the output of the initial awk to bedtools sort commands as a variable (see below)

swapped=$(awk '{if ($2>$3)print $1,$3,$2,".",".","-";else print $1,$2,$3,".",".","+";}' OFS='\t' Infile.bed | awk '{a=$2-1;print $1,a,$3,$4,$5,$6;}' OFS='\t' | bedtools sort) 

This works quite nicely:

echo "${swapped}"
HEADING_1   4   12  .   .   +
HEADING_2   4   12  .   .   -

When I use the variable in the bedtools getfasta command no output is generated. Is there a way to for a variable to be read like a file? I have tried the following, but it is still not working:

  1. bedtools getfasta -s -fi Infile.fasta -bed "${swapped}" -fo Outfile.fasta
  2. bedtools getfasta -s -fi Infile.fasta -bed <(echo "${swapped}") -fo Outfile.fasta
  3. bedtools getfasta -s -fi Infile.fasta -bed <(<<< "${swapped}") -fo Outfile.fasta

Basically, can I use a variable in place of a file as an argument for a command?

I hope that this makes sense

Thanks,

Jamie

1

There are 1 answers

0
Charles Duffy On BEST ANSWER

If you look in the test suite for the bedtools getfasta command, you'll see that it passes the word stdin as a filename when it wants the BED input to be read from stdin. For example:

LINES=$(echo $'chr1\t1\t10' | $BT getfasta -fi t.fa -bed stdin -fo - | awk 'END{ print NR }')

So, we just need to do that same thing in your script:

bedtools getfasta -s -fi Infile.fasta -bed stdin -fo Outfile.fasta <<<"$swapped"

By the way -- in most cases, your second attempt would have worked:

bedtools getfasta -s -fi Infile.fasta -bed <(echo "${swapped}") -fo Outfile.fasta

...insofar as the <(...) expression is replaced by a filename from which the output at hand can be read. (There are some caveats: It's typically passed through a /dev/fd link, so any program that closes file descriptors other than the default stdin, stdout and stderr won't be able to read from content given that way; also, insofar as that filename is an end of a FIFO, anything that needs to be able to seek around in input, read it more than once, check its size before reading, etc. won't work).