I've got a question about the data from a biological experiment.

Three times the same 1024 different proteins are spotted on one testing chip. Target of the experiment is to see whether certain proteins in the lysate bind to particular proteins on the chip. They are spotted three times to increase reliability. The chip is graphically analysed after radioactive ATP made a print by a third party.

After analysis I get a matrix of 3 columns and 1024 rows, where the rows are the values of the same protein, but they differ a bit. Depicted below is are the highest data points in a boxplot, plotted on their lowest measured value. The first 600 are similar to the ones on the left on the chart; they run almost in a flat line.

enter image description here

In brown there is the median, in blue the median + 3*MAD, in yellow the Inter Quartile Range * 2,2 (note, there is 1,5 in the title, it should be 2,2)

To tell if a protein is bound (gets and "ON-call"). I use the Median Absolute Deviation-method on the lowest values of the measurements. I use this method because the number of on-calls might be significant (200 of 1024 proteins), and greatly influence the mean in other outlier detection methods — like you can see it does with the IQR. This gives me a criterion: median + 3*(1,4826*MAD). I have taken the 1,4826 factor from literature. If a value is higher, it is an outlier; in this case that means the protein is bound.

The question: Since I only use the lowest values of the measurements, I basically ignore the other values. In particular I am unsure of the measurements of a protein, where the lowest point is under the criterion, but the average is not. After removing outliers the data is normalish distributed (a tat skewed).

enter image description here

Or maybe I am completely wrong and there is an entirely other way to detect values that correspond to bound proteins. Please tell me what you think!

In general: how to detect outliers between points after three measurements?

  • 2
    $\begingroup$ MAD is ambiguous: median absolute deviation from the median is one interpretation but as m stands for mean too, I can imagine three other interpretations at least. You write of the MAD method as if it were standard and known. I doubt that there is any such standard method and even if there were it would be helpful to explain it. Note further that many people here prefer to see data examples upfront, not on password-requiring sites or in proprietary format datasets. Otherwise put, no-one flagged this as unclear, but lack of clarity is probably why it has received no response to date. $\endgroup$ – Nick Cox Apr 4 '16 at 12:55
  • 1
    $\begingroup$ In addition to clarifying about what you mean by "MAD," "flagged outlier," and "on-call, "and providing perhaps a small sample of data directly, it would help a lot if you could give a bit more background about the nature of the measurements (are these 3 replicates under the same conditions, or evaluations under 3 different experimental conditions?) and the scientific question you are trying to answer. Answers might also be different, for example, if data were ELISA results versus microarray or RNAseq data. $\endgroup$ – EdM Apr 4 '16 at 16:28
  • $\begingroup$ In general, the first measurement is the lowest and the third (or last) measurement is the highest. Did you rearrange the order within the rows to reflect this or is this an artifact of the process used? Why is outlier detection important? Why is the method you've described, and that most of us find ambiguous, not sufficient? $\endgroup$ – Mike Hunter Apr 4 '16 at 17:32
  • 1
    $\begingroup$ Are the median and the MAD values based on all 3072 data points? What are those values? Am I correct that by an "outlier" you mean a protein that is substantially bound so that you will give it an "ON-call"? Are you primarily interested in discovering all the proteins that might reasonably be thought to be bound, or rather in using this procedure to identify a few candidates for more detailed analysis? $\endgroup$ – EdM Apr 4 '16 at 19:47
  • 1
    $\begingroup$ Please edit the question itself to clarify, as there Is some non-standard terminology. It seems that you want to know which of these proteins have high signals in this assay system, and what you are calling "outliers" are really the high-signal proteins you care about. The word "outliers" more typically means values that are suspiciously out of line with some type of expectation. This MAD method seems also to be non-standard, perhaps from some type of quality-control-chart analysis that might not be useful on this type of data. Details on these "chips" and what the values represent would help. $\endgroup$ – EdM Apr 5 '16 at 9:00

There are three important points here. First, the procedure you are using to identify bound proteins is a heuristic approach rather than a formal statistical test. Second, since your data are fundamentally based on rates of radioactive decay, more formal statistical tests might be possible. Third, regardless of whether you stick with heuristics or move to statistical tests, the answer depends on your tradeoffs between false-positive and false-negative results.

Heuristics. The "median + 3* MAD" approach originates from control-chart analyses formulated nearly a century ago at Bell Labs, in a different context. The original context was to determine when something might be going wrong with a manufacturing process. If values measured on a process were drawn from a normal distribution, then fewer than 3 out of 1000 values are expected to be more than 3 standard deviations (SD) away from the mean value. Even if the underlying distribution is not normal, this "3-sigma" criterion provides a useful guide for prompting an evaluation of whether the "outlier" values represents a process that is out of control. But it's a heuristic guide to further action, not a formal statistical test.

Your application of this heuristic approach is in a different context, as you expect deviations from the center of the distribution to be frequent, with many spots on your array bound by proteins in your lysate. For you, these "outliers" are desirable and are the point of your study. Using the median instead of the mean as the measure of the center of the distribution of values, and the median absolute deviation (MAD*) or the inter-quartile range instead of the SD as a measure of width, avoids some of the complications that thus arise, as these robust measures will be much less affected by the many high values that you expect. For example, in the context of screening all genes to see which ones have significant effects upon knockdown, median-MAD-based heuristics work better than mean-SD for identifying gene candidates for more detailed examination.

Formal statistical tests. In your application you have a pretty good idea of an important underlying statistical process, so you might be able to replace these heuristic assessments with formal statistical tests. Each of your values represents a result of the decay of radioactive phosphorus, a classic example of a Poisson process. A quick look at your data suggests that your values are not decay counts per se, as the within-protein variances seem to exceed the associated mean values as you would expect in a Poisson distribution (the variances do, however, tend to rise with mean values). There is presumably some post-decay detection/amplification process that leads to higher than Poisson variance, in addition to other sources of error such as random differences among the 3 spots for each protein.

A generalized linear model that treats your data as coming from an overdispersed Poisson distribution might lead to an analysis in which standard statistical significance tests are possible, providing within-protein errors (in a log scale for Poisson models) that are normally distributed and independent of the mean values. Modern statistical software provides such analyses readily. This approach has the advantage of removing reliance of the significance of the test result for each individual protein from the vagaries of its own particular 3 values, instead basing the estimate of the underlying measurement error on information pooled from measurements on all 3000+ spots.

False-positives versus false-negatives. Your initial question of whether to stick with the lowest of 3 values for each protein or to use the average in your heuristic assessment of protein binding points out a fundamental decision that you must make. What are your relative risks of making "ON-calls" that are erroneous, versus missing proteins that really are bound but don't pass your criterion for an "ON-call"?

The requirement that the lowest of the 3 values pass something similar to a 3-sigma test is quite stringent, and your top plot suggests that sticking with that requirement will miss a fair number of likely true binding events. Using the average instead, or changing the criterion to something like "median + 2*MAD" will provide more "ON-calls" but runs the risk of making "ON-calls" for more proteins that are not really bound. So there is no single correct answer to that part of your question; you have to decide on your tradeoffs.

If you move to formal statistical tests instead, you do not avoid making this type of decision but you may be able to present your results to others in a more defensible way. You then are faced with the classic multiple comparisons problem. If you perform 1000 significance tests with a criterion of p < 0.05 for significance, your will get about 50 positive results even if there are no true positive cases. There are two general ways of addressing this problem.

You can try to control the familywise error rate, trying to make sure than you will make a false-positive error in only 5% of similar experiments. A highly conservative example is the Bonferroni correction, which in this case would require, for 1000 tests, that you only accept significance of an individual test at p < (0.05/1000) = 0.00005. Less stringent corrections are possible, as described in the linked pages.

Alternatively, you can accept that you are going to make some false-positive calls and try to limit the fraction of such calls among all the calls that you make. This is called controlling the false discovery rate. This is a well-accepted approach for evaluating results of microarrays, of which your application is an example. You might, for example, choose to aim for a 10% false discovery rate, so that about 90% of your "ON-calls" would be true. That could be readily justified in a scientific publication.

* Note that there is an ambiguity in the use of "MAD" even in the context of "median absolute deviation." What is reported as the "MAD" is typically not the median of the absolute values of all deviations from the median, but that result scaled by 1.4826 so that the MAD and the SD should coincide for a normal distribution.

  • $\begingroup$ Greatly appreciated! Indeed i have used 1,4826 to scale the median absolute deviation. Forgot to mention it, so edited it. $\endgroup$ – Max Apr 6 '16 at 17:34
  • $\begingroup$ Would it be of any added value to this post to describe my final method, based on the advice you have given me? $\endgroup$ – Max Apr 7 '16 at 12:03
  • $\begingroup$ @Max, maybe you could leave a brief comment under this answer stating your choice of method. Thanks for working so hard at clarifying your question; that made it easier for me to answer and for others who see this question to understand the issues. $\endgroup$ – EdM Apr 7 '16 at 19:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.