Some elegant bash / awk here. +1
• 1 views
•
link
written
by
ThePresident
bwa for multiple fastq files
I would like to do bwa annotation for forward and reverse reads from different samples against a reference genome using for loap
I have different samples like:
A2_forward.fastq A2_reverse.fastq A3_forward.fastq A3_reverse.fastq A4_forward.fastq A4_reverse.fastq
Also, I need the output for annotation of each sample in separate, as I need to compare the effects of different treatments in different samples.
• 13,492 views
•
link
updated
by
bsp017 •
•
written
by
dmotelb87 •
•
written
by
dmotelb87 •
1 answer
Do you men to say, you want to do mapping of different sample against your reference genome ? You can use following script for mapping multiple sample using bwa. (Note: save the below script in .sh format (say "mapping.sh") and make sure that are no other fastq file rather then your sample should be there in your working directory)
script for mapping using bwa mem
total_files=`find -name '*.fastq' | wc -l`
arr=( $(ls *.fastq) )
echo "mapping started" >> map.log
echo "---------------" >> map.log
for ((i=0; i<$total_files; i+=2))
{
sample_name=`echo ${arr[$i]} | awk -F "_" '{print $1}'`
echo "[mapping running for] $sample_name"
printf "\n"
echo "bwa mem -t 12 A1_denovo_assembly.fasta ${arr[$i]} ${arr[$i+1]} > $sample_name.sam" >> map.log
bwa mem -t 12 A1_denovo_assembly.fasta ${arr[$i]} ${arr[$i+1]} > $sample_name.sam
}
Follow along using the numbers below as line numbers to the script above.
- Counting number of files that end in
*.fastq. If your files have different endings then substitute accordingly. - Assigning the names of the files to an array called
arr. - Just printing a message to screen and append to log file called
map.log. - Just printing a message to screen and append to log file called
map.log. - A
forloop to step through the array taking two file names at a time (this is a paired-end dataset soforwardandreversepair, will need to change depending on your use case) until you take all elements from array. Array elements start at position[0]. sample_nameextracts the actual name of the sample (stuff before first_e.g.Sample1fromSample1_fastq)- Printing the name you are working with to screen with
echo - Printing the actula command you are working with. Name of output file created using
$sample_nameas$sample_name.sam. - Run the actual
bwa memcommand.
I tried this
total_files=`find -wholename './sample_data/*.fastq' | wc -l`
arr=( $(ls ./sample_data/*.fastq) )
But then bwa mem runs for two pairs of input files in the folder but output only one file named sample.bam
However, I am expecting two output files, A1.bam and B1.bam
I would try replacing this line instead:
total_files=`ls /sample_data/*.fastq | wc -l`
Keep an eye on
echo "[mapping running for] $sample_name"
to make sure sample names are being printed properly. You should get as many BAM files as there are pairs of sample files.
So I tried
total_files=`ls ./sample_data/*.fastq | wc -l`
arr=( $(ls ./sample_data/*.fastq) )
for ((i=0; i<$total_files; i+=2))
{
sample_name=`echo ${arr[$i]} | awk -F "_" '{print $1}'`
echo "[mapping running for] $sample_name"
printf "\n"
bwa mem -t 12 ref.fa ${arr[$i]} ${arr[$i+1]} | samtools sort -o output/$sample_name.sorted.bam
Echo says
[mapping running for] ./sample
And got one sorted.bam file not with the name of input file but as sample.sorted.bam
It is reading sample as a file name which is a problem.
Log in to answer this question.
More posts like this
-
Harmony integration of integarted samples
written by fifty_fifty •
-
Question about gffcompare using multiple reference annotations from different databases
written by newbie
-
RNA-Seq: X samples and X batches
written by Behram Radmanesh •
-
Automating to convert multiple fastq files into one fastq file
written by zhou_1228 •
-
Can we concatenate two fastq files from same sample but different runs
written by Payal
-
should I merge fastq files for different lanes before do QC?
written by Lila M
-
How to merge paired-end reads from sam files?
written by aquaq •
-
How to Separate a SAM File Containing Forward and Reverse Reads into Two Separate Files
written by kdevakan •
-
Comparing assembled transcripts against a reference
written by Rob
-
Transcriptome Assembly-Gene Annotation
written by jackuser1979
Thanks Suraj. Yes. This what I mean. It is working now.
If an answer was helpful you should upvote it, if the answer resolved your question you should mark it as accepted.
Please use
ADD COMMENTorADD REPLYto answer to previous reactions, as such this thread remains logically structured and easy to follow. I have now moved your reaction but as you can see it's not optimal. Adding an answer should only be used for providing a solution to the question asked.