How to directly compare normalized read counts across subgroups within a metatranscriptome?

Question

[N15]

I would very much like feedback on this metatranscriptomic/DESEq2 issue from a biostat perspective.

I have 19 individual RNAseq samples corresponding to different locations in the environment. They were sequenced separately using Illumina Hiseq (barcoded), and then all reads were used to generate an assembled metatranscriptome.

Reads from each of the 19 samples were mapped to this metatranscriptome. Based on our annotations, we have many different microbial groups present.

What I would like to do is compare transcript abundance across locations for each major taxonomic group. I have subsetted out contigs and their raw counts corresponding to the groups A, B and C. I have normalized counts for each group in DESeq2 and exported the pseudo-counts.

I would like to generate transcriptome profiles (heatmaps) for each of these groups and look for changes in normalized transcript abundance across sites and among groups. The question is, is it incorrect to display these pseudo-counts side-by-side each other in the same heat map? (I think yes.) Is there a better way to do this if you wish to directly compare normalized transcript abundance across independently-normalized taxonomic groups?

I don't think the answer is to normalize groups A,B,C all together and then later subset them out, because there may be important differences in expression tendencies between groups.

I am not actually interested in a pairwise analysis here, mainly just the changes in transcript abundance with location and group.

Apologies if there is an obvious answer here. Thank you so much for your time and thoughts.

Answer 1

[theobroma22]

You can make a heatmap or plot however you want, it just has to be visually acceptable and fit your interpretation.

Also, if you want to compare changes in transcript abundance with location and group, you could use ANOVA, no?

[N15]

I don't think it's valid to compare normalized counts across data sets. Think about this: you have 1 data set consisting of 2 treatments and 3 reps each (EXPT1) and another data set with 2 different treatments and 3 reps (EXPT2). You normalize both independently. It wouldn't make any sense to then compare all pseudo counts between EXPT1 and EXPT2.

So how would I get around this problem in my metatranscriptome where I normalized each subset (taxa) separately?

[theobroma22]

Oh right, I overlooked you're using count data while writing my reply, so disregard my last. Well, in your case and your current data setup, you can use Venn Diagrams, bar plots and heatmaps or other graphical displays like radar plots. In other words, there is no more statistical analysis to do, only plot the data you have and write your results.

Answer 2

[pedrorodrigues]

In regard to group-level normalization rather than global, I think the method proposed by Klingenberg and Meinicke 2017 might be helpful. I have personally not been able to use their method yet, but it looks promising.