# Interpreting a kernel density function plot of gene expression

I need help interpreting a kernal density plot and a histogram of this data. I have provided some context below let me know if you need more information

A co-worker that I am collaborating with has a list of genes ranked by expression fold change. He says most of the genes in one dataset are noise with respect to his analysis so he is trying to filter this out. He has sent me a document of his methods for generating a new list of genes which he think are relevant. he has run a permutation test on the ranked genes I have provided a quote here:

"To determine which genes contributed more strongly in real data, as compared by permuted data, the rankings achieved for each gene in real data were compared to the shuffled rankings. For each gene, and each permutation type, a p-value was calculated by counting how many times that gene was ranked higher in shuffled data than in real data, and dividing it by the number of comparisons (16 * 256 = 4096). To achieve p<0.05 a gene’s real ranking would have to be higher than 95% of its shuffled rankings."

Following this he provides kernel density plots of the genes under different rankings.

"Below are kernel density estimation (kde) plots, overlaid on histograms, of the rankings of two genes randomly selected from the few genes with p<0.05 vs all three permutations." "Below are kde plots of the rankings of two genes randomly selected from the many genes with p>0.05 vs all three permutations." What am I looking at? What does this tell me about the genes ranking statistics.

• Are you asking what is a kernel density plot? The Matlab how-to is the most informative. Or are you asking why this enrichment analysis was favoured over others? – Matteo Ferla Nov 23 '19 at 22:52
• Both I guess. I don't understand what the x and y axes are in either of the graphs. Is the second chart showing probability density and the histogram the actual distribution? – Angus Campbell Nov 24 '19 at 5:37
• A kernel density plot is a like a histogram, but smoothed. Your coworker has given you rough data, e.g. unlabelled axes and little explanation. The x-axis is number of genes and the y-axis is the "density", which isn't "number of counts in a bin", but a number so that the area under the curve is one (it's continuous not discreet). So for now, I strongly suggest not worrying about that and concentrating on the gene enrichment analysis. The permutation test is to get statistical support of the enrichment set. Therefore, what you need is the enriched gene list not the metrics of support. – Matteo Ferla Nov 24 '19 at 10:11
• If x axis is the number of genes, how come there are negative values in the second plot? – haci Nov 24 '19 at 15:54
• No. The negative values are due to the kernel density function extrapolating beyond the 0 x-value in order to justify the high sample at zero. It is a common problem with KDE, or smoothing in general. Normally bounds are set beforehand to avoid this or, sloppily, the window set to exclude these. – Matteo Ferla Nov 25 '19 at 14:32

There are two topics here to be discussed. Gene Set Enrichment Analysis and a kernel Density plot.

## GSEA

The Wikipedia page about Gene Set Enrichment Analysis covers well the principle and the tools that can be used to calculate it. To cite it:

The general steps are summarized below:

Calculate the enrichment score (ES) that represents the amount to which the genes in the set are over-represented at either the top or bottom of the list. This score is a Kolmogorov–Smirnov-like statistic [a ranking based statistic].

Estimate the statistical significance of the ES. This calculation is done by a phenotypic-based permutation test in order to produce a null distribution for the ES. The P value is determined by comparison to the null distribution.

Calculating significance this way tests for the dependence of the gene set on the diagnostic/phenotypic labels

Basically, you want to test if the ranking is random or not. Therefore you make a random dataset and test whether your set is like that. One assumption underlying this is that the output of the gene whatever it may be is not normally distributed (parametric, follows a bell curve), hence the use of ranking, which removes this (you lose signal but you lose biases). I am assuming the non-parametric nature of your ouput values was tested or is a known non-parametric value, such as RPKM or log2 fold change (by virtue of being a ratio). In other words, the magnitude difference between your samples is lost, but their order is not.

There are several tools, some online, that do this enrichment testing for you depending on the task, some even work with gene ontologies.

## KDE

A kernel density plot is a like a histogram, but smoothed, albeit not in a moving average way. In a histogram you divide your x-axis into bins: it is discreet and you get an integer count per bin. In a kernel density plot the data is fitted to "probability density function", an equation which given x will give y, where the integral of the curve is one, hence why y is density not counts. The x axis will be your bins-like value, in this case rank. The binning window metric for a KS function is called bandwidth.

See Matlab documentation for a fairly non-opaque explanation.

However, if you want to compare histograms, I strongly suggest violin plots, not because they are strongly in vogue, but because you can compare histograms without worrying about overlaying them.

• Thank you for these resources. – Angus Campbell Nov 26 '19 at 3:30