Parzen窗估计

Octave/MATLAB代码

function [p] = parzen(data, range, h, f)
  
  if isempty(f)
    f = @(u)(1 / sqrt(2 * pi)) * exp(-0.5 * u.^2);
  end;
  
  N = size(data, 2);
  h = h / sqrt(N);
  [X Xi] = meshgrid(range, data);
  p = sum(f((X - Xi) / h)) / N / h;

解释

参数

预处理

结果

示例

Octave/MATLAB代码

S = 1000;
xi = [rand(1, S) * 0.5 - 2.5, rand(1, S) * 2 + 1];
x = linspace(-3, 4, 400);
p = parzen(xi, x, 22, @(u)(1 * (abs(u) < 0.5)));
plot(x, p);
figure();
p = parzen(xi, x, 5,[]);
plot(x, p);

解释

结果

正态核函数

f1

方核函数

f2

Kn近邻估计

Octave/MATLAB代码

function [p] = knn(data, range, kn)
  N = size(data, 2);
  [X Xi] = meshgrid(range, data);
  Dis = abs(X - Xi);
  for i = 1:kn-1
    [~, ind] = min(Dis);
    ind=sub2ind(size(Dis), ind, [1:size(Dis, 2)]);
    Dis(ind) = Inf;
  end
  Dis = min(Dis);
  V = 2 * Dis;
  p = kn * ones(1, size(range, 2)) ./ V / N;

解释

参数

预处理

寻找Kn近距离

通过min运算,且每次将得到的最小值改为正无穷,得到第近的距离。 最终Dis为各待估点的第近距离。

窗体积

其中

是待估点与近邻点的距离

所以仅考虑一维情况

结果

示例

Octave/MATLAB代码

S = 5000;
xi = [sqrt(rand(1, S)) + 0.5, -sqrt(rand(1, S)) + 2.5, rand(1, 2 * S) + 3];
x = linspace(0, 5, 100);
p = knn(xi, x, sqrt(S) * 20);
figure();
plot(x, p);

解释

结果

f3