This returns the result as a distance matrix such as produced by pdist2, except that all entries not corresponding to neighboring points are set to NaN. Neighbors are defined as points whose Voronoi polygons share one or more (finite) vertices.
[V,C]=voronoin(points);
F=all(isfinite(V),2);
D=pdist2(points,points);
D(~common_vertex(C,F))=nan,
function map = common_vertex(C,F)
n=numel(C);
map=false(n);
for i=1:n
for j=1:i-1
map(i,j)=any( F(intersect(C{i},C{j})) );
end
end
map=logical(map+map.');
map(1:n+1:end)=1;
end