Hi

I've been using ROAST and Camera, but I am bit unsure of how much to trust each. Let's say that I want to test if a just a a couple of gene sets are enriched in an RNAseq dataset. roast/mroast should do the trick, right? Well, I am under the impression that it's very easy for roast to produce interesting p-values, even if it's just random data.

As an example, let's randomly generate gene expression for 10000 genes, and make 1000 of these DE expressed between 2 groups. Then generate two lists with 100 genes: one where 50% of the genes are DE, and the other where genes are picked entirely at random.

set.seed(0)

# generate random expression 10000 genes
y <- matrix(rnorm(10000 * 4), 10000,4)
design <- cbind(Intercept = 1, Group = c(0, 0, 1, 1))

# make 1000 genes differentially expressed
y[1:1000, 3:4] <- y[1:1000, 3:4] + 3

# index 1 - 50% genes are DE
DE_50 <- 951:1050

# index 2 - random genes
DE_random <- sample(1:10000, 100)

mroast(y, list(DE_50 = DE_50, DE_random = DE_random),
       design = design, contrast = 2, set.statistic = 'mean50')

This gives:

          NGenes PropDown PropUp Direction PValue   FDR PValue.Mixed FDR.Mixed
DE_50        100     0.02   0.55        Up  0.002 0.003        0.002     0.003
DE_random    100     0.11   0.17        Up  0.030 0.030        0.009     0.009

Quite interesting p-values for a random list! ROAST appears to have a huge false positive rate!

In contrast, Camera behaves as expected:

> camera(y, DE_50, design = design, contrast = 2)
     NGenes Direction       PValue
set1    100        Up 5.534736e-14
> camera(y, DE_random, design = design, contrast = 2)
     NGenes Direction    PValue
set1    100      Down 0.3403619

This makes me quite unsure about trusting ROAST. Camera is/was said to be quite conservative, but partly this could be just because ROAST gives many false positives. Are there any real advantages of using ROAST over Camera?

Thanks, Nuno

ROAST is doing its job here. In ROAST, the null hypothesis is that none of the genes in the set are DE. You will find that the intersection of DE_random and the 1:1000 true DE genes is non-empty, which means that the null hypothesis is not true. So ROAST is correct in returning a low p-value. You can argue about the sensitivies of the different set statistics, given that the documentation suggests that using "mean50" would only detect DE proportions above 25%, but there is nothing inherently wrong with ROAST itself.

If you read the documentation in ?roast, you will find that ROAST and CAMERA test different hypotheses. ROAST performs a self-contained gene set test, where it looks for any DE within the set of genes. CAMERA performs a competitive gene set test, where it compares the DE within the gene set to the DE outside of the gene set. Thus, it is not surprising that the results of these two procedures will differ.

Personally, I find the results of ROAST to be easier to interpret than the results of CAMERA. In many cases, I just want to know whether the genes in my gene set or pathway are systematically DE. For this, a self-contained gene set test is more suitable. In contrast, the output of a competitive gene set test depends not only on the DE of genes within the set, but also on DE in genes outside the set. This can lead to counterintuitive behaviour where the significance of one gene set depends on the significance (or lack thereof) of other sets.

Do you care about whether a gene set is "more DE" than other sets? If so, then you should use a competitive gene set test like CAMERA. If you just want to know whether the gene set is DE in one (or any) direction, use a self-contained gene set test like ROAST.

Hi Nuno,

Aaron has explained why your simulation is not showing any false discoveries. But I do wonder what documentation you are reading. The help pages for roast and camera are only a few paragraphs long and they explain how the two functions test different hypotheses. Do you understand the difference between a self-contained gene set test and a competitive test? The help pages try to explain this in simple terms. If you understand this, then the results from your simulation shouldn't be at all surprising. The version of camera() that you are using is BTW not conservative in practical situations. Again, the documentation would tell you that.

Your simulation is actually a good illustration of the (intended!) differences between two functions. Your "random" set is DE in absolute terms, with 11 out of 100 genes being strongly DE all in the same direction. In many circumstances we would want roast() to detect this set as DE.

Your question fails to mention that, by default, roast does not detect the random set as DE. Instead it gives P=0.33 for the random set:

> mroast(y, GeneList, design)
          NGenes PropDown PropUp Direction PValue   FDR PValue.Mixed FDR.Mixed
DE_50        100     0.05   0.53        Up  0.003 0.005        0.007     0.013
DE_random    100     0.09   0.13        Up  0.330 0.330        0.022     0.022

This is because we have tuned roast() by default so that it does not give small p-values when the proportion of DE genes in a set is small. So, rather than illustrating a "problem", your simulated example actually shows that roast() does exactly what you seem to want it to do. Indeed roast gives almost the same p-value for the random set that camera does.

But, rather than relying on default settings, you have specially asked roast to try to detect a smaller proportion of DE genes, by setting set.statistic="mean50". So it is difficult to understand why you should be at all surprised that roast then does exactly what you asked it to do, and does then succeed in detecting the 11% DE genes in the set.

I think the real message here is: don't ask roast() to be sensitive to very small proportions of DE genes unless that is what you actually want.

Personally, I like both roast() and camera() and use them both regularly. But to see a situation when camera() would give you unintuitive results, let's make a small change to your simulation. In this simulation, all the genes are DE:

> set.seed(0)
> y <- matrix(rnorm(10000 * 4), 10000,4)
> design <- cbind(Intercept = 1, Group = c(0,0,1,1))
> y[, 3:4] <- y[, 3:4] + 3
> GeneList <- list(DE_50 = 951:1050, DE_random = sample(1:10000, 100))
> mroast(y, GeneList, design)
          NGenes PropDown PropUp Direction PValue   FDR PValue.Mixed FDR.Mixed
DE_random    100        0   0.91        Up  0.001 0.001        0.001     0.001
DE_50        100        0   0.96        Up  0.002 0.002        0.002     0.002
> camera(y, GeneList, design)
          NGenes Direction PValue   FDR
DE_50        100        Up  0.116 0.132
DE_random    100      Down  0.132 0.132

Now, roast finds both sets to be significant but camera finds neither. Now tell me, how you are going to explain to your lab head or biological collaborator why camera finds neither pathway to be DE, even though it is obvious that every single gene in both pathways is strongly DE, with every gene DE in the same direction by the same amount? Note that camera's failure to detect DE here has nothing to do with conservativeness - camera is not conservative. It has everything to do with the meaning of a self-contained test vs a competitive gene set test.

Best

Gordon