Оцените pdf вектора с помощью ядра Гаусса

Я использую ядро ​​Гаусса для оценки PDF-файла данных на основе уравнения введите здесь описание изображениягде K (. ) - гауссовское ядро, данные - заданный вектор. z - это ячейка от 1 до 256. Размер ячейки равен 1.

Я реализовал код Matlab. Однако результат показывает, что амплитуда моей оценки PDF (синий цвет) не похожа на реальный PDF-файл данных. Не могли бы вы увидеть мой код и прокомментировать мой код?

КОД MATLAB

function pdf_est=KDE()
close all;
%%Random values of 20 pixels, range=[1 256]
data=randi([1 256],1,20);

%% Estimate histogram%%%%% 
pdf_est=zeros(1,256);
z=256;

for i=1:z
    for j=1:length(data)
        pdf_est(i)=pdf_est(i)+Gaussian(i-data(j));
    end
end
%% Plot real histogram 1 to 256; binsize=1;
hold on
plot(imhist(uint8(data))./length(data),'r');
%% Plot histogram estimation
plot(pdf_est./length(data),'b');
hold off
function K=Gaussian(x)
   sigma=1;
   K=1./(sqrt(2*pi)*sigma)*exp(-x^2./(2*sigma^2));

РЕЗУЛЬТАТ СИНИЙ - мой результат, КРАСНЫЙ - настоящий pdf введите описание изображения здесь


matlab statistics probability probability-density kernel-density
person user3051460    schedule 16.02.2015    source источник
comment
Я предполагаю, что они нормализованы, верно? тем не менее, pdf выглядит довольно правильно стимулированным   -  person Ander Biguri    schedule 16.02.2015
comment
@AnderBiguri: Оба pdf разделены по длине (данные). Так что я думаю, что они уже нормализованы.   -  person user3051460    schedule 16.02.2015

Ответы (1)


arrow_upward
3
arrow_downward

У вас есть две проблемы:

  1. Смещение на 1 единицу между синим и красным графиками.
  2. Синие шипы шире и меньше красных.

Как решить каждую проблему:

  1. Это вызвано возможной путаницей между диапазоном данных 0, ..., 255 и интервалом индексации 1, ..., 256. Поскольку ваши данные представляют собой 8-битное изображение, значения должны быть 0, ..., 255 (не 1, ..., 256). Ваша горизонтальная ось на графике должна быть 0, ..., 255. То же самое и с переменной i в строке for. И затем, поскольку индексирование Matlab начинается с 1, вы должны использовать i+1 при индексировании pdf_est.

  2. Это нормальное поведение. Вы предполагаете дисперсию единиц в своем ядре. Чтобы увидеть более высокие синие шипы, вы можете уменьшить sigma, чтобы сделать ядро ​​уже и выше. Но вы никогда не получите ту же высоту, что и ваши данные (необходимое sigma будет зависеть от ваших данных).

    Фактически, у вас есть компромисс между высотой и шириной, которым управляет sigma. Но важно то, что область остается неизменной для любого sigma. Поэтому я предлагаю построить CDF (площадь) вместо pdf (плотность площади). Для этого постройте накопленную гистограмму и pdf (используя cumsum).

Код изменен согласно 1:

function pdf_est=KDE()
close all;
%%Random values of 20 pixels, range=[1 256]
data=randi([1 256],1,20)-1; %// changed: "-1"

%% Estimate histogram%%%%% 
pdf_est=zeros(1,256);
z=256;

for i=0:z-1 %// changed_ subtracted 1 
    for j=1:length(data)
        pdf_est(i+1)=pdf_est(i+1)+Gaussian(i-data(j)); %// changed: "+1" (twice)
    end
end
%% Plot real histogram 1 to 256; binsize=1;
hold on
plot(0:255, imhist(uint8(data))./length(data),'r'); %// changed: explicit x axis
%% Plot histogram estimation
plot(0:255, pdf_est./length(data),'b'); %// changed: explicit x axis
hold off
function K=Gaussian(x)
   sigma=1; %// change? Set as desired
   K=1./(sqrt(2*pi)*sigma)*exp(-x^2./(2*sigma^2));

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

Код изменен согласно пунктам 1 и 2:

function pdf_est=KDE()
close all;
%%Random values of 20 pixels, range=[1 256]
data=randi([1 256],1,20)-1; %// changed: "-1"

%% Estimate histogram%%%%% 
pdf_est=zeros(1,256);
z=256;

for i=0:z-1 %// changed: subtracted 1 
    for j=1:length(data)
        pdf_est(i+1)=pdf_est(i+1)+Gaussian(i-data(j)); %// changed: "+1" (twice)
    end
end
%% Plot real histogram 1 to 256; binsize=1;
hold on
plot(0:255, cumsum(imhist(uint8(data))./length(data)),'r'); %// changed: explicit x axis
                                                            %// changed: cumsum
%% Plot histogram estimation
plot(0:255, cumsum(pdf_est./length(data)),'b'); %// changed: explicit x axis
                                                %// changed: cumsum
hold off
function K=Gaussian(x)
   sigma=1; %// change? Set as desired
   K=1./(sqrt(2*pi)*sigma)*exp(-x^2./(2*sigma^2));

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

person Luis Mendo    schedule 16.02.2015
comment
@ user3051460 Рад, что смог с этим помочь! Это был интересный вопрос - person Luis Mendo; 16.02.2015
comment
да. Фактически, у нас есть функция, аналогичная идее в kr.mathworks.com/ помощь / статистика /. Однако в моем случае приведенная выше формула более проста и подходит для сегментации изображений. Но когда я запускаю вышеуказанный код, это занимает много времени (потому что у нас есть два цикла). Не могли бы вы увидеть мою формулу и попробовать оптимизировать, чтобы сократить время вычислений? - person user3051460; 16.02.2015
comment
@ user3051460 Я думаю, что некоторые петли (возможно, оба) можно устранить векторизацией. Но было бы сбивать с толку делать это здесь, так как это не то, что задавалось в вашем исходном вопросе. Предлагаю вам создать для этого отдельный вопрос. Выберите версию модифицированного кода, которую вы хотите использовать, опубликуйте ее и попросите векторизовать ее. Думаю, это легко сделать. Кто-то (возможно, я, если будет время) обязательно поможет с этим. - person Luis Mendo; 16.02.2015
comment
да. Я разместил его по адресу stackoverflow.com/questions/28545361/ - person user3051460; 16.02.2015