Muszę obliczyć obszar, na którym nakładają się dwie funkcje. Używam normalnych dystrybucji w tym konkretnym uproszczonym przykładzie, ale potrzebuję bardziej ogólnej procedury, która również dostosowuje się do innych funkcji.

Zobacz obraz poniżej, aby dowiedzieć się, co mam na myśli, gdzie czerwony obszar jest tym, za czym się podążam:

To jest jak dotąd MWE:

import matplotlib.pyplot as plt
import numpy as np
from scipy import stats

# Generate random data uniformly distributed.
a = np.random.normal(1., 0.1, 1000)
b = np.random.normal(1., 0.1, 1000)

# Obtain KDE estimates foe each set of data.
xmin, xmax = -1., 2.
x_pts = np.mgrid[xmin:xmax:1000j]
# Kernels.
ker_a = stats.gaussian_kde(a)
ker_b = stats.gaussian_kde(b)
# KDEs for plotting.
kde_a = np.reshape(ker_a(x_pts).T, x_pts.shape)
kde_b = np.reshape(ker_b(x_pts).T, x_pts.shape)


# Random sample from a KDE distribution.
sample = ker_a.resample(size=1000)

# Compute the points below which to integrate.
iso = ker_b(sample)

# Filter the sample.
insample = ker_a(sample) < iso

# As per Monte Carlo, the integral is equivalent to the
# probability of drawing a point that gets through the
# filter.
integral = insample.sum() / float(insample.shape[0])

print integral

plt.xlim(0.4,1.9)
plt.plot(x_pts, kde_a)
plt.plot(x_pts, kde_b)

plt.show()

gdzie się ubiegamMonte Carlo aby uzyskać całkę.

Problem z tą metodą polega na tym, że gdy oceniam próbki w obu dystrybucjach za pomocąker_b(sample) (lubker_a(sample)), Dostaję wartościbezpośrednio nad linia KDE. Z tego powodu nawet wyraźnie pokrywające się rozkłady, które powinny zwracać wspólną / nakładającą się wartość powierzchni bardzo zbliżoną do 1. zwracają zamiast niej małe wartości (całkowity obszar każdej krzywej wynosi 1., ponieważ są to oszacowania gęstości prawdopodobieństwa).

Jak mogę naprawić ten kod, aby uzyskać oczekiwane wyniki?

W ten sposób zastosowałem odpowiedź Zhenyi

# Calculate overlap between the two KDEs.
def y_pts(pt):
    y_pt = min(ker_a(pt), ker_b(pt))
    return y_pt
# Store overlap value.
overlap = quad(y_pts, -1., 2.)