RNAseq: Why would using a more complete scaffold reduce DEGs?

Question

I have a set of RNAseq data for a couple of relatively under-studied species.

They're fungi, of the order Hypocreales. The assemblies we're working with are de novo - reference seqs don't actually exist yet for most of our species. What the collaborator provided was updated genome assemblies in 60-90 fragments as opposed to the original 250-300, but obviously that's still not full resolution. As for the gene families, not a whole lot is known... in fact our species' genes are still not well described and we're mostly working off of predicted orthology.

During the analysis, a collaborator created new DNA scaffolds for these species which we'd now like to use.

However when I redid our RNAseq using the existing reads with the new scaffolds, one of our species had the number of differentially-expressed genes drop by anywhere from 20-50% per condition.

As the reads themselves haven't changed, we weren't expecting significant changes in differential gene expression, and indeed for our other species changes have been minimal.

Does anybody know what might cause this kind of change?

Thank you!

Answer 1

Fungi, ain't a kingdom most of us will know much about. The only thing I can guess at is that in bacteria there has resurgence in methods for de novo assembly specifically aimed at reducing the total number of contigs in an assembly. The 'new wave' is summarised in this response here, Improving prokaryotic assembly with other contig/scaffold-level data?. This is potentially very exciting, except the OP reported they didn't work ;-(

A reasonable guess is your colleague has deployed these approaches to reduce the number of contigs. HOWEVER, what is appropriate for bacteria, may not be appropriate for fungi, alternatively it may be an improvement. Both bateria and fungi are small genomes, but beyond that there could notable differences.

If fungi are dominated by multigene families these might be being removed in the new wave de novo assemblies (they are designed for bacteria). The collapse of multigene families absolutely occurs in other eukaryotes with small genomes, but they have to be really big on multigene families (>10% of genes in their genome). If that was true you would probably know about it.

Alternative, this could simply be poor gene identification in the first assembly resulting in an artificial jump in the number of genes.

Solution? Either way you should be able to figure which genes are different between the assemblies and precisely what type of genes are differing, e.g. are these large multi-gene families. Personally, I would pipe the predicted genes from both assemblies through eggNOG. This might be useful, but not really what it's designed for, but it will return the KEGG nos. More explicit methods based around comparing the RNAseq reads against the assembly via e.g. bedtools or these guys here recommend htseq. This might be a better way to rephrase the question.

Thus I recommend:

1. Finding out which de novo methods were used between first and second fungi genome assemblies. Do any of those methods match this post, Improving prokaryotic assembly with other contig/scaffold-level data?
2. Consider a more explicit genomics study to examine RNAseq reads against the genome, or comparison between predicted no. genes between genomes.
3. Having said all that ... once you've done your DEG you identify the KEGG no. and then simply compare the lists directly to find out which genes are different between the assemblies. This is an easy database operation to perform, e.g. via a merge function (separate question).

I'd check point 1 and then go for point 3 - because it's something I easily know how to do. Thats just my opinion, a purist would likely explore this via point 2.

Conclusion What I am certain about is that for a species where very little might be known it is better to perform an investigation rather than make an assumption that might prove later to be wrong.