I have raw paired-end RNA-seq reads for two novel eukaryotic species. Some background: the reads represent a copepod (arthropod) species each. The mRNA for each read set was obtained by extracting mRNA from multiple individuals. I am not exactly sure what sequencer was used (some Illumina short read machine as far as I can tell; I am not directly involved in this project).
I'd like to assemble these reads into transcriptomes de novo as no reference genomes are available as far as I know; my hope is to assemble with Trinity and then annotate with Dammit using its default databases (and check for completeness w/ BUSCO).
However, it looks like something went wrong during the sequencing runs (or during sample preparation) because the sequence GC content plots from FastQC look like this:
Here the curve with the tall peak is one sample, and the more flatter curve is the other. There are two closely related curves per sample because FastQC apparently seems to have recognized each read (forward and reverse) as a separate sample.
I first presumed that this is just rRNA contamination (I believe the mRNA was prepared via poly-A selection), and I filtered the reads against all of sortmerna's default databases (which I believe include sequences from SILVA and Rfam among others). Then I passed the sortmerna'd read through FastP to filter out adapters, paired read correction, and low complexity filtering (alongside whatever else FastP does by default). I had ran FastQC before FastP to check if the rRNA filtering did anything and it hadn't, so I also asked FastP to filter by length (--length_required 151) since there were small amounts of reads with lengths shorter than this value which I thought could be the source of the additional peaks in the data.
However, even after sortmerna and FastP together, nothing has really changed:
(There are only two traces/curves here instead of four because the paired reads were apparently recognized properly by FastQC now.)
My plan was to run Khmer to clean the reads a little bit more and then assemble them. I did do a quick assembly with Trinity (with the raw reads) and the BUSCO scores were pretty bad (only ~ 25% of the arthropod BUSCO genes present in the assemblies). I was hoping I'd be able to get Dammit to annotate whatever was left after the cleaning, despite the poor BUSCO prognosis.
Now I'm not entirely sure what to do if I can't clean the reads to begin with, since preventing nonsensical contigs would be absolutely important for a decent annotation.
So this leads me to my questions:
- Are these data even salvageable?
- Is there anything else I can do to clean the data?
- Is this even contamination to begin with that I'm seeing here? (I am presuming some kind of contamination because both samples share that peak at ca. 18% GC content, among other shared features.) Would mapping all the reads against each other and discarding the shared sequences help? If yes, is there a tool I can use for this?
- Any other thoughts in this regard?
I searched through the posts here and on Biostars as well as I could but I wasn't really able to find anything helpful in this regard (hence my own post).
I'm very new to all of this, so any inputs, tips, advice, and/or guidance would be most appreciated.
