How do you run a BWA alignnment command in R?

Question

[carlabertulfo]

Hi, I'm trying to align my reads to a reference genome using the bwa alignment. I used this command wherein I already specified the path for my BWA folder ("path_bwa"):

system(paste(path_bwa, "/bwa aln -t ", proc, " ", genome, " FASTQ_FILTERED.fq > ", file_sai, sep=""))

However, it always shows this error: sh: /Applications/bwa/bwa: No such file or directory

Somehow, I think it assumes that the "/" in the "/bwa aln -t " is a directory and it is trying to look for it.

I tried deleting the "/" and just put "bwa aln -t " but it still doesn't work. Do you have any idea on what's wrong with the code?

Thank you so much for your help. It will be greatly appreciated.

Best, Carla

[Pierre Lindenbaum]

citing OP on twitter :

Hi, thanks for the response! I'm basically just running a published software/command line for my 4C seq data. The scripts are in R.

— Carla Bertulfo (@carlabertulfo) September 19, 2017

[Kevin Blighe]

I can easily invoke my BWA program by using full paths like this:

bwa.path <- "/Programs/bwa-0.7.12/"
system(paste(bwa.path, "bwa", " index", sep=""))

What is the value of path_bwa and proc ?

[h.mon]

In addition to all comments on why do you want to wrap this on R, are you sure you want to use bwa aln? bwa mem is the recommended algorithm for most uses, unless your reads are really short, bellow 70bp.

[Pierre Lindenbaum]

1) why R ? why ? why ? why not a good-old bash/shell script ?

2) sounds like a path problem. check paste(path_bwa, "/bwa . What is the output of

find  /Applications/ -name "bwa"

and

which bwa

[shoujun.gu]

just curious why you don't do it in shell? And whether it's working in shell?

[Jeremy Leipzig]

...or why not use one of the many fine pipeline frameworks actually designed for this?

Answer 1

[viniciushs.z]

Maybe the below wrapper can help

https://rdrr.io/github/flow-r/ngsflows/man/bwa.html