StatQuest: edgeR, part 1, Library Normalization

Key Takeaways

that quest is super cool that class won't make you truth that quest stand quest guaranteed not to make you drool hello and welcome to stack quest stack quest is brought to you by the friendly folks in the genetics department at the University of North Carolina at Chapel Hill today we're starting part 1 of our exploration of how edge our works we're going to talk about library normalization just like Dec to our does not use rpkm TPM or any of the other standard normalization techniques this is because it needs to adjust for two things sequencing depth that's what rpkm and all those other methods deal with it also needs to take library composition into account this means different samples can contain different active genes and that can change things I cover this concept in depth in the stat quest on Dec 2 part 1 library normalization so check that out if you have any more questions without further ado let's talk about how edge are normalizes libraries step 1 remove all untranscribed genes that is to say remove all genes with 0 read counts in all samples here's an example we've got 3 samples and 5 genes the last gene gene 5 has 0 read counts in all 3 samples and so we will remove gene 5 from our list of genes step two pick one sample to be the reference sample this is the sample that we will use to normalize all of the other samples against imagine these blue dots are samples or libraries a jar chooses one sample to be the reference sample and uses it to normalize all of the remaining samples here's a question for you what's a good reference sample alternatively what's a bad reference sample here's an example of an extremely bad reference sample sample three would be a terrible reference sample scaling would be based on a single potentially very noisy measurement to avoid choosing extreme samples edge our attempts to identify the most average sample let's see how it does this in order to pick a reference sample the first thing it does is it scales each sample by its total read counts so here if we imagine each sample has four genes we just sum the read counts for each sample and we divide the original read counts for each gene in each sample by each samples total and here's what the scaled read counts look like the second part of picking a reference sample is for each sample determine the value such that 75% of the scaled data or equal to or smaller than it in sample number one three of the four values 75% are less than or equal to zero point two six in sample number two three of the four values 75% are less than or equal to 0.36 lastly in sample number three three of the four values or 75% are less than or equal to zero point one three now we've calculated the 75th quantiles for each sample the third part of picking a reference sample is to average the 75th quantiles in this case the average is 0.25 the reference sample is the one whose 75th quantile is closest to the average in this case sample number 1 will be the reference sample that's because its 75th quantile 0.26 is closest to the average 0.25 now that we've picked a reference sample we now need to select the jeans for calculating the scaling factors this is done separately for each sample relative to the reference sample here we see the reference sample we selected in step number two and now we will select a set of genes to create a scaling factor for sample number two and we will select another set of genes to create a scaling factor for sample number three this is one of the ramifications of edge ARS approach different samples use different genes to derive their scaling factors let's see how we select genes for sample number two we'll start by looking at the different types of genes to choose from he's an XY plot that will demonstrate the different types of genes we have to choose from on the left side we see a gene primarily transcribed in the reference sample on the right side we see a gene primarily transcribed in sample number two way up high we see a gene with a ton of reads in both samples and way down low we see a gene with hardly any reads in both samples genes in the middle don't have much bias and have a moderate number of reads mapped to them in both samples edge our selects the genes in the middle with more effort put into excluding biased genes let's see how it does this we'll start with the scaled read counts these are the read counts for each gene divided by the samples total the next few steps only makes sense if there are a lot of genes too we need to put on the screen so that's what the dot-dot-dot and gene n are all about they represent lots more genes edge our filters out biased genes by looking at log fold differences remember with logs if the reference is way high relative to sample number two you'll get a value way up here and if sample number two is way high relative to the reference we'll get a value way down here ultimately we'll pick a threshold like plus or minus six and remove all genes with more extreme biases for more information about logs check out the stat quest on logs we'll start by calculating the log to ratio for gene number one plugging in the scaled read counts we get the log two of zero divided by zero point nine the log 2 of 0 is defined as negative infinity by are the programming language that edge R was written in so we put negative infinity as the value for gene 1 in our table of log 2 ratios now let's calculate the log 2 ratio for gene number 2 plugging in the scaled read counts gives us the log base 2 of 0.8 which equals negative 0.32 and we just do that for every single gene in our list of genes now we remove genes with infinite values ie genes without any reads mapping to them in one or both samples in this example we'll remove gene number 1 now that we have a table of log ratios to identify biased genes let's make a table to identify genes that are highly and lowly transcribed in both samples to identify genes that are high and low in both samples first calculate the geometric mean for each gene remember the geometric mean is cool since it is not easily influenced by outliers more details on the geometric mean are in the stat quest for logs anyway here's how we calculate the geometric mean for gene number one we take the average of the log two of the scaled read counts for each sample plugging in the numbers gives us negative infinity that's because there's zero read counts in the reference for gene number one and we put that value in our table of the mean of the logs technically we are not calculating the geometric mean since we are not converting the average back into normal number land but the effect is the same outliers are less influential on the value that we are keeping track of here's how we calculate the mean of the logs for genome or two and then we calculate the mean of the logs for all of the other genes in the list now we remove all genes with infinite values ie genes without any reads mapping to one or more samples in this case that means we will be removing gene number one now we have two tables one to identify biased genes that's the log of the ratio between the two samples and one to identify genes that are highly and lowly transcribed in both samples and that's the mean of the logs now sort both tables from low to high filter out the top 30% and the bottom 30% of the biased genes and filter out the top 5% and the bottom 5% of the highly and lowly transcribed genes genes that are still in both lists are used to calculate the scaling factor unfortunately the only genes in this example that are in both lists are the dot-dot-dot genes well you get the idea hooray we figured out which genes to use to scale sample number two now we need to figure out which genes to use for sample number three I'll leave this as an exercise for the reader now that we figured out which genes we're going to use to create a scaling factor we can move on a step four calculate the weighted average of the remaining log two ratios for your information edge our calls this the weighted trimmed mean of the log to ratios because we trimmed off the most extreme genes by excluding extreme genes we avoid the effect of outliers sort of like using the geometric mean okay so we're back to talking about sample number two and imagine genes a through Z are the dot-dot-dot genes these are the genes that made it through the filters in step 3 once you have selected which genes will be used to calculate the scaling factor just calculate the weighted average of their log2 ratios genes with more reads map to them get more weight this is because log ratios have more variance with low read counts let's look at an illustration here we have a table of read counts for sample number one and sample number two we've got six genes and on the right side we show the log to ratio between the two samples the first three genes have relatively high read counts however the difference between gene number one and gene number three is for the bottom three genes have relatively low read counts however the difference between gene four and gene six is also four now in the right hand column we see the log two ratios and when there are lots of reads small differences in recount do not make big log fold differences when there are only a few reads small differences and read counts can make big log fold differences so we're back to calculating the weighted average of the log two ratios genes with more reads mapped to them get more weight because they're less noisy and then we do the same thing for sample number three we calculate the weighted average of the log two ratios step 5 convert the weighted average of the log two ratios to normal numbers that is to say we raise to two the weighted average of the sample two log two ratios that gives us a scaling factor however this is not quite the scaling factor that edge our reports we also calculate a scaling factor for sample number three and we continue to do that for all the samples in our data set here are actual scaling factors that I calculated from an rna-seq experiment we see that WT 2 was used as the reference sample the other samples were scaled to it 1 step 2 step 6 Center the scaling factors around 1 here our raw scaling factors the values are roughly centered around 0.95 and here are the final centered scaling factors they've been shifted over so that they are centered on one the sintering is done by dividing each raw scaling factor by their geometric mean although sintering does not change the results mark robinson the guy who wrote this method says it gives the scaling factor some nice mathematical properties so I guess it's a sort of artistic signature on a mathematical formula that's it now you know how scaling factors are calculated in a jar tune in next time for another exciting stat quest