I am trying to find a way to calculate the overlapping (intersection) and non-overlapping areas of two distributions. Based on some posts (like this: FINDING AREA) I could figure out how to calculate the intersection of the two plots. But I had no luck finding the non-overlapping part.
Given this example distributions:
import pandas as pd
import numpy as np
import seaborn as sns
data1=np.random.normal(0, 0.1, 100)
data2=np.random.normal(0, 0.3, 100)
x0=pd.Series(data1)
x1=pd.Series(data2)
kwargs = dict(hist_kws={'alpha':.01}, kde_kws={'linewidth':2})
sns.distplot(x0,bins=3, color="dodgerblue", **kwargs)
sns.distplot(x1,bins=3, color="dodgerblue", **kwargs)
I wonder how I can calculate the areas between the two distributions at the edges?
I should mention that I have the data as well (if that is something that I could do with the data itself).
Advertisement
Answer
Taking inspiration from the answer you linked, from data you should compute kde curves through scipy.stats.gaussian_kde:
data1 = np.random.normal(0, 0.1, 100) data2 = np.random.normal(0, 0.3, 100) kde1 = gaussian_kde(data1, bw_method = 0.5) kde2 = gaussian_kde(data2, bw_method = 0.5) xmin = min(data1.min(), data2.min()) xmax = max(data1.max(), data2.max()) dx = 0.2*(xmax - xmin) xmin -= dx xmax += dx x = np.linspace(xmin, xmax, 1000) kde1_x = kde1(x) kde2_x = kde2(x)
Then you need to find intersection points between the kde curves:
idx = np.argwhere(np.diff(np.sign(kde1_x - kde2_x))).flatten()
where idx is the list of index of intersection points.
Finally, you can compute area through numpy.trapz, using index previously computed to slice x, kde1_x and kde2_x:
area1 = np.trapz(kde2_x[:idx[0]] - kde1_x[:idx[0]], x[:idx[0]]) # area under the lower kde, from the first leftmost point to the first intersection point area2 = np.trapz(kde2_x[idx[1]:] - kde1_x[idx[1]:], x[idx[1]:]) # area under the lower kde, from the second intersection point to the last rightmost point area3 = np.trapz(np.minimum(kde1_x, kde2_x), x) # intersection area between of the kde curves, between the two intersection points area4 = np.trapz(kde1_x[idx[0]:idx[1]] - kde2_x[idx[0]:idx[1]], x[idx[0]:idx[1]]) # area under the highest kde, excluding the lower one, between the two intersection points
Complete Code
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import gaussian_kde
data1 = np.random.normal(0, 0.1, 100)
data2 = np.random.normal(0, 0.3, 100)
kde1 = gaussian_kde(data1, bw_method = 0.5)
kde2 = gaussian_kde(data2, bw_method = 0.5)
xmin = min(data1.min(), data2.min())
xmax = max(data1.max(), data2.max())
dx = 0.2*(xmax - xmin)
xmin -= dx
xmax += dx
x = np.linspace(xmin, xmax, 1000)
kde1_x = kde1(x)
kde2_x = kde2(x)
idx = np.argwhere(np.diff(np.sign(kde1_x - kde2_x))).flatten()
area1 = np.trapz(kde2_x[:idx[0]] - kde1_x[:idx[0]], x[:idx[0]]) # area under the lower kde, from the first leftmost point to the first intersection point
area2 = np.trapz(kde2_x[idx[1]:] - kde1_x[idx[1]:], x[idx[1]:]) # area under the lower kde, from the second intersection point to the last rightmost point
area3 = np.trapz(np.minimum(kde1_x, kde2_x), x) # intersection area between of the kde curves, between the two intersection points
area4 = np.trapz(kde1_x[idx[0]:idx[1]] - kde2_x[idx[0]:idx[1]], x[idx[0]:idx[1]]) # area under the highest kde, excluding the lower one, between the two intersection points
fig, ax = plt.subplots()
ax.plot(x, kde1_x, color = 'dodgerblue', linewidth = 2)
ax.plot(x, kde2_x, color = 'orangered', linewidth = 2)
ax.fill_between(x[:idx[0]], kde2_x[:idx[0]], kde1_x[:idx[0]], color = 'dodgerblue', alpha = 0.3, label = 'area1')
ax.fill_between(x[idx[1]:], kde2_x[idx[1]:], kde1_x[idx[1]:], color = 'orangered', alpha = 0.3, label = 'area2')
ax.fill_between(x, np.minimum(kde1_x, kde2_x), 0, color = 'lime', alpha = 0.3, label = 'area3')
ax.fill_between(x[idx[0]:idx[1]], kde1_x[idx[0]:idx[1]], kde2_x[idx[0]:idx[1]], color = 'gold', alpha = 0.3, label = 'area4')
ax.plot(x[idx], kde2_x[idx], 'ko')
handles, labels = ax.get_legend_handles_labels()
labels[0] += f': {area1 * 100:.1f}%'
labels[1] += f': {area2 * 100:.1f}%'
labels[2] += f': {area3 * 100:.1f}%'
labels[3] += f': {area4 * 100:.1f}%'
ax.legend(handles, labels)
plt.show()
