Normally, I use bwa-mem to map my fastq reads for a bunch of species using one protein-coding gene CDS from my reference species. This was fine when I had to try it for just a group of genes...I would obtain single CDS fastq files and use them as a reference. However, we now want to obtain protein-coding sequences for ALL GENES in my exome datasets. If I have one file containing all the CDS from my reference species, can that be used in bwa-mem? The end of my custom script produces one multiple sequence alignment file when I ran it using one gene only...my concern is that bwa-mem will somehow replace certain CDS or would I just end up with over 20000+ MSAs? Which is fine because I need them to be individual MSAs for each protein-coding gene.
I use a simple python script at the end to extract base information from the mpileup files generated to create my MSAs.