Как нормализовать kde из scikit Learn?

Допустим, у меня есть массив формы (100000,1), представляющий выборки переменной X с равномерным распределением от 0 до 1. Я хочу аппроксимировать плотность вероятности этой переменной, и для этого я использую Scikit-Learn KernelDensity.

Проблема в том, что я получаю только результат, который не нормализован. Интеграл плотности вероятности не равен 1. Как мне сделать, чтобы нормализовать автоматически? Я делаю что-то неправильно ?

def kde_sklearn(data, grid, **kwargs):
    """
    Kernel Density Estimation with Scikit-learn

    Parameters
    ----------
    data : numpy.array
        Data points used to compute a density estimator. It
        has `n x p` dimensions, representing n points and p
        variables.
    grid : numpy.array
        Data points at which the desity will be estimated. It
        has `m x p` dimensions, representing m points and p
        variables.

    Returns
    -------
    out : numpy.array
        Density estimate. Has `m x 1` dimensions
    """
    kde_skl = KernelDensity(**kwargs)
    kde_skl.fit(data)
    # score_samples() returns the log-likelihood of the samples
    log_pdf = kde_skl.score_samples(grid)
    return np.exp(log_pdf) 

X = np.random.uniform(0,1,1000).reshape(-1,1)
X1 = np.linspace(0,1,100)[:,np.newaxis]

kde_sklearn(X,X1,kernel='tophat')

Out[43]: 
array([0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5,
       0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5])

Я ожидал, что вектор будет равен 1, так как сумма интеграла должна быть равна 1.


person Raphael Benezra    schedule 09.08.2019    source источник


Ответы (3)


Проблема не в нормализации, как я могу показать на примере. Предположим, что я запускаю следующий код, который подгоняет KDE к образцам из стандартного нормального распределения:

import numpy as np
import sklearn.neighbors as sn

# Sample from a standard normal distribution
XX = np.random.randn(1000).reshape(-1, 1)

# Fit a KDE
kde_sklg = sn.KernelDensity()
kde_sklg.fit(XX)

# Get estimated densities
XX1 = np.linspace(-4.0, 4.0, 100)[:, np.newaxis]
gdens = np.exp(kde_sklg.score_samples(XX1))

Затем я могу оценить площадь под PDF с помощью правила трапеций следующим образом:

my_area = 0.0
for i in range(1,gdens.shape[0]):
    my_area += 0.5*(gdens[i] + gdens[i-1])*(XX1[i,0] - XX1[i-1,0])

Расчетная площадь (my_area), которую я получаю, составляет около 0,996, что довольно близко к 1.

Проблема в том, что ваш KDE не обрабатывает скачки в вашем единообразном PDF, которые происходят на 0 и 1, поэтому он слишком их размывает. Примерно половина области под оценкой KDE вашего PDF-файла оказывается под этими размытыми областями. Если вы замените значение вашего X1, скажем, на X2 = np.linspace(-1,2,200)[:,np.newaxis], вы увидите значительную плотность в частях оценки KDE для PDF в интервалах [-1,0] и [1,2].

person jjramsey    schedule 09.08.2019
comment
Хороший ответ. Спасибо чувак :). Я постараюсь обучить свою модель с большим количеством примеров в моем образце, я считаю, что размытие должно исчезнуть. - person Raphael Benezra; 10.08.2019
comment
@RaphaelBenezra Я не уверен, но вам могут понадобиться образцы за пределами интервала [0,1], чтобы все работало. Вы также можете поэкспериментировать с различными ядрами, пропускной способностью и т. д. - person jjramsey; 12.08.2019

Опубликованные ответы, на мой взгляд, неясны, и поэтому я даю другой ответ.

Короче говоря, сумма integral равна 1, а не вероятности. Ниже я покажу 2 способа получить интеграл, который действительно равен 1.

import numpy as np
from sklearn.neighbors import KernelDensity

np.random.seed(1)

# some uniform data
X = np.random.uniform(-5,5,100).reshape(-1,1)

# grid to be used later0
grid = np.linspace(-5,5,1000)[:,np.newaxis]

# fit using the data
kde = KernelDensity(kernel = 'tophat', bandwidth= 0.5).fit(X)

# get log probailities of the grid
log_dens = kde.score_samples(grid)

# transform log prob to prob
probs = np.exp(log_dens)

# Integrate
print(np.trapz(probs.ravel(), grid.ravel()))
0.9732232232232225

plt.hist(X, density=True, bins=30)
plt.plot(grid.ravel(),probs.ravel())
plt.show()

Обратите внимание, что другой способ получить интеграл заключается в следующем, поскольку у нас есть тот же шаг в заданной сетке:

np.sum(probs*np.diff(grid.ravel())[0])
0.9732232232232225

введите здесь описание изображения

person seralouk    schedule 20.10.2020

это вероятности в каждой точке - что произойдет, если

X1 = np.linspace(0,1,10000000)[:,np.newaxis]

?

массив, который вы получаете, не является распределением/выборкой из случайной величины

person quester    schedule 09.08.2019
comment
Я получаю то же самое. Вектор 0,5. - person Raphael Benezra; 09.08.2019
comment
посмотрите, эти числа не являются вероятностями какого-то события, а являются p-значениями некоторых статистических данных, поэтому они не должны суммироваться с 1 - person quester; 09.08.2019
comment
@quester Это не p-значения, а плотности вероятности, и интеграл функции плотности вероятности по своей области должен быть равен 1. - person jjramsey; 09.08.2019
comment
@quester Кроме того, это не настоящий ответ. Это должен быть комментарий. - person jjramsey; 09.08.2019