Holy freakin' smokes! What? It's time for StatQuest!
Hello, and welcome to StatQuest. StatQuest is brought to you by the friendly folks in the Genetics Department at the University of North Carolina at Chapel Hill. Today, we're going to be talking about false discovery rates, or FDR.
If you've ever seen or done anything with high-throughput sequencing, chances are you've heard of false discovery rates, FDR, before. You may have even used them. But where do false discovery rates come from? And how do they work?
Before we get down to the nitty-gritty, let me blurt out the main idea of this whole StatQuest. False discovery rates are a tool to weed out bad data that looks good.
Now let's get down to the nitty-gritty.
Let's start with an example of measuring gene expression using RNA sequencing. Here, we're going to plot the measurements, or the read counts for a gene called gene X, which is an imaginary gene, on a graph with the Y-axis being gene counts and the X-axis being the samples. For this example, imagine that we are looking at normal, wild-type mice. Later on, we'll be comparing them to mice that have been treated with a drug.
Isn't it funny that normal mice are called wild-type? If someone said I was a wild type, I don't think they would also think I was normal.
Anyways, here's our measurement for the first mouse that we do this RNA sequencing on. And here's the measurement for the second mouse that we do the RNA sequencing on. RNA sequencing isn't perfect, and different samples are always a little different. So each time we measure expression, we'll get slightly different values.
Here's the measurement for the third mouse. And here are the measurements for a bunch of mice. If we measured all normal mice, we'd be able to calculate the average for all normal mice.
Most of the values are going to be close to the mean. Rarely, we'll get a value that is much larger than the mean. And rarely, we'll get a value that is much smaller than the mean.
We can summarize the distribution of the measurements using this bell-shaped curve. Most of the measurements, which are close to the mean, will come from the middle of this curve. The rare measurement that is significantly larger than the average would come from the right side of the bell-shaped curve. And the rare measurement that is significantly less than the average would come from the left side of the bell-shaped curve.
Now, imagine that we do RNA-Seq on three mice. Collectively, we'll call these three measurements sample number one. Because these measurements are close to the mean, they come from the middle of the distribution. Now, imagine we compare sample number one to another three measurements taken from normal, wild-type mice. We'll call these new measurements sample number two.
Again, these three measurements come from the middle of the distribution. If we did a statistical test to compare sample number one to sample number two, the P-value would be large, greater than 0.05, because the two samples overlap. Very rarely, we'll get two samples that do not overlap. When this happens, the P-value will be less than 0.05. This is called a false positive, because the small P-value suggests that the samples are from two types of mice, or two separate distributions, and this is false.
Normally, false positives are rare. Unless you're a P-hacker, but that's another StatQuest already on YouTube. Anyways, normally, false positives are rare. 95% of the time, the two samples will overlap. This will result in a P-value greater than 0.05. 5% of the time, they don't. This will result in a false positive with a P-value less than 0.05.
But human and mouse cells have at least 10,000 transcribed genes. If we took two samples from the same type of mice and compared all 10,000 genes, well, 5% of 10,000 equals 500 false positives. That means there will be 500 genes that appear to be interesting, even when they are not. 500 false positives is a lot. Can we do something about them?
The false discovery rate can control the number of false positives. Technically, the false discovery rate is not a method to limit false positives, but the term is used interchangeably with the methods. In particular, it is used for the Benjamini-Hochberg method. Now there's a high probability that I just mispronounced Benjamini or Hochberg, and if I did, I apologize.
Before we talk about the details of the Benjamini-Hochberg method, let's review the concepts that it's based on. We'll start by generating 10,000 P-values from samples taken from the same distribution. That is to say, we'll start with test number one, and we'll use wild-type mice, and we'll take two samples from them. We'll then compare the two samples with a statistical test and calculate the P-value. In this case, the P-value is large, it's 0.83. This is exactly what we expect, because both samples are taken from the same type of mice. And then we repeat the procedure for test number two and calculate another P-value. This time it's 0.98. Again, this is as expected. To make a long story short, we just repeat this procedure 10,000 times. No big deal.
Here, I've drawn a histogram of the 10,000 P-values generated by testing samples taken from the same distribution. On the X-axis, we have possible values for P-values. On the Y-axis, we have the number of P-values in each bin. 510 P-values, or 5.1%, are less than 0.05. Close to 5% of the P-values are between 0.05 and 0.1. Actually, each bin contains about 5% of the P-values, about 500 P-values per bin.
Since the P-values are uniformly distributed, there is an equal probability that a test P-value falls into any one of these bins.
Now, let's look at how P-values are distributed when they come from two different distributions. And by two different distributions, I mean two different types of mice, or we have wild-type versus knockout, or control versus drug, we're just comparing two different situations.
Like before, we start off with test number one, but now we have two different distributions. The black distribution is for our control mice. The red distribution is for mice that have been treated with a drug. In this example, the drug increases this gene's transcription.
Like before, we take two samples. Since the samples are now coming from two separate distributions, there's a higher likelihood that the two samples will be separated and not overlap. When we do the statistical test, in this case, we get a P-value that equals 0.03. And then we do the exact same thing for test number two. Notice that both of the P-values are less than 0.05, so they're statistically significant. Since the samples were taken from two separate distributions, this is what we'd expect.
Like before, we repeat this process 10,000 times.
Here, I've drawn a histogram of the 10,000 P-values generated by testing samples taken from two different distributions. Most of the P-values are less than 0.05. This is what we'd expect. The P-values greater than 0.05 are false negatives from where the samples overlapped. Pssst. You can reduce the number of false negatives by increasing the sample size.
To summarize what we know so far, when the samples come from the same distribution, the P-values are uniformly distributed. But when the samples come from different distributions, the P-values are heavily skewed and closer to zero.
Now, imagine we're doing an experiment where we are testing all of the active genes in neuronal cells. One set of neuronal cells is treated with a drug, the other is not. The drug might affect 1,000 genes. The measurements for these genes will come from two different distributions. The black sample is from the control cells, and the red sample is from the cells treated with the drug. Since the samples come from different distributions, the P-values are skewed. The remaining 9,000 active genes might not be affected by the drug. This means the measurements for most of the genes will come from the same distribution. The P-values for these genes should be uniformly distributed.
The histogram of P-values we obtain from all 10,000 genes is the sum of the two separate histograms. The uniformly distributed P-values come from the genes unaffected by the drug. The P-values on the left side are a mixture from genes affected by the drug and genes unaffected by the drug. By eye, we can see where the P-values are uniformly distributed and determine how many tests are in each bin. Here, I've drawn a line indicating that about 450 P-values are in each bin in the uniformly distributed part of the histogram. We can extend this line and use it as a cut-off to identify the true positives.
Since we usually use a cut-off of 0.05, we're going to focus on these P-values. Roughly 450 P-values less than 0.05 are above the dotted line. And roughly 450 P-values less than 0.05 are below the dotted line. One way to isolate the true positives (genes affected by the drug) from the false positives would be to only consider the smallest 450 P-values. This procedure works fairly well because the P-values within the bins are skewed for the genes affected by the drug (Note: This histogram is for P-values between 0 and 0.05) and spread evenly for the genes not affected by the drug. Bam!
If you can understand these concepts, then you understand more about false discovery rates and the Benjamini-Hochberg method than most people who use it. All Benjamini and Hochberg did is convert this procedure that we just did by eye into a mathematical formula. So now let's talk about the details of the Benjamini-Hochberg method. Like I just said, it's based on the eyeball method we just saw.
The Benjamini-Hochberg method adjusts P-values in a way that limits the number of false positives that are reported as significant. Adjusted P-values means that it makes them larger. For example, before the false discovery rate correction, your P-value might be 0.04, i.e., significant. After the FDR correction, your P-value might be 0.06, no longer significant.
If your cut-off for significance is FDR values less than 0.05, then less than 5% of the significant results will be false positives. In other words, these are the genes with P-values less than 0.05. The black box shows the genes with FDR modified P-values less than 0.05. Notice that not all of the true positive genes are inside the box. However, only 5% of the modified P-values in the box are false positives. The remaining 95% are true positives.
Why don't all of the true positive genes have adjusted FDR P-values less than 0.05? Because not all true positive genes will have super small P-values. Here's the histogram of true positive P-values less than 0.05. These genes on the right side of the histogram probably won't remain significant after adjustment. Surprisingly, the math behind the Benjamini-Hochberg method is simple. Let's take a look.
Let's start with another simple example. We'll take 10 pairs of samples taken from the same distribution, i.e., 10 genes that were not affected by the drug. And here are the P-values from those 10 tests. The first thing we do is we order the P-values from smallest to largest. Notice that one of the P-values is a false positive. That is to say, it's less than 0.05.
Let's see what the Benjamini-Hochberg method does to it. The second step is to rank the P-values. And let's make spaces for the FDR adjusted P-values that we're going to create. The largest FDR adjusted P-value and the largest P-value are the same. The next largest adjusted P-value is the smaller of two options: either the previous adjusted P-value, which in this case equals 0.91, or the current P-value times the total number of P-values divided by the P-value rank. In this case, the current P-value is 0.81, the total number of P-values is 10, and the P-value rank is 9. Plugging these numbers in, we get 0.90. Since we select the smaller option, we're going to go with 0.90.
For the next largest adjusted P-value, we just repeat step four. That is to say, we select the smaller of these two options. Plugging in the numbers gives us the choice of 0.90 or 0.89. Since we use the smaller value, we go with 0.89. And for the next largest adjusted P-value. Okay, you get the idea. We just repeat until we've adjusted the remaining P-values.
And here I've plugged in the remaining adjusted P-values. That false positive P-value is no longer significant. Hooray!
Now let's look at a huge example. The blue boxes represent the P-values from when the samples came from two separate distributions. That is to say, these P-values are for genes that were affected by the drug. I've made these P-values relatively small to reflect their normal skew. The P-values in red boxes came from samples taken from the same distribution. Note, we've got some false positives.
The eyeball method suggests we draw a line at the top of the uniformly distributed P-values, and extend it to separate the false positives from the true positives. These are the P-values that the eyeball method suggests are true positives. Now let's see what the Benjamini-Hochberg method does.
Here, I've shown you the adjusted P-values. The false positives are now all greater than 0.05. But these true positives remain less than 0.05. Double bam!
Hooray, we made it to the end. Tune in next time for another exciting StatQuest.
