Hi all,

I am developing an analysis of protein microarray data and have found that limma performs very well in terms of identifying positive controls in the assay. Unfortunately, there are a large number of known, non-specific interactors that previous people have been identified in earlier experiments. For the pipeline I would like to down-weight, but not exclude, these proteins from the analysis. I have an idea, but I'm not sure if it's a good one:

The thought is to add a value to the intensity of a probe that depends on its average across past experiments. This would essentially maintain the variance while reducing the fold changes between groups for that probe. It kind of places a higher burden of proof on the that interaction for it to be real.

The problem is that I won't be seeing true fold changes anymore when I look at my output and that seems pretty bad. Does anyone out there have any other solutions to this sort of problem or can think of a way to moderate the t values by incorporating old data? Or is there some feature of limma/eBayes that could handle probe weights when calculating p-values?

Any thoughts would be appreciated!

Thanks, Andrew

The main problem here is that it's hard to convert a vague sense of "I don't like these proteins" into a concrete, quantitative weight. You say that you want to downweight the problematic interactors, but by how much? What does it mean to be twice as problematic? How many angels dance on the head of a pin?

Now let's consider your proposed approach. It seems that you're proposing to add a constant value - dependent on the average intensity of the probe from past experiments - to all observations for a probe. If this is done on the log-scale, then it will have no effect on anything (trend-based shrinkage aside), as the added value will cancel out when computing variances or log-fold changes. If this is done on the raw scale, then you're adding an arbitrary pseudo-count that will affect both the log-fold changes and variances of the log-values - not good for empirical Bayes shrinkage. More generally, there is no sound rationale for adding the average intensity. Why not half the average? Or the square root of the average? Or any number I might pull out of thin air?

The better approach is to make use of the definition of these problematic interactors. You say that these were identified as being non-specific in previous experiments. If this is non-specific DE, presumably you have log-fold changes for these proteins from previous experiments as well. This means you can test for changes in expression beyond this non-specific log-fold change using the treat function. Some testing suggests that lfc can support numeric vector inputs, so you can give a different lfc per gene. If the log-fold changes are larger than what would be expected due to non-specific activity, then there's probably a genuine change in expression.

Or you could take the safe route and just discard the problematic interactors. This would be my preferred approach if there's not too many of them. Why take the risk of reporting them when they're known to be non-specific? I know my collaborators would be pretty skeptical.

From your description, it sounds like you are simply worried about probes with large variances. limma already downweights probes with large variances. Belinda Phipson and I have even written a paper specifically about how to deal with probes that have particularly large variances (which we call "hypervariable" genes):

Robust hyperparameter estimation protects against hypervariable genes and improves power to detect differential expression

and have implemented a method for dealing with them efficiently through the limma empirical Bayes procedure. To activate this strategy, specify robust=TRUE when running eBayes() or treat().

So limma already does what you want. It will generally deal with hypervariable probes very sensibly, requiring them to show a much larger fold change than other genes before they are counted as DE.

If you have prior data that truly represents the variability of the different probes, then you could consider adding that to your experimental design and analyse it at the same time as you data of interest. That does assume however that the probes will continue to have the same underlying variances in the new experiment as in the old, and I think that it is a very strong assumption. I think this will usually be unnecessary and a straightforward analysis will be enough. Generally speaking, your attempts to micro-manage the process will just get in the way of limma doing its job.