CODEWITHSUNDEEP
matlabplotgaussian
teaLeef
teaLeef

Reputation: 1989

Gaussian Probabilities plot around a trajectory

I am trying to write some code to generate a plot similar to the one below on matlab (taken from here):

I have a set of points on a curve (x_i,y_i,z_i). Each point generates a Gaussian distribution (of mean (x_i,y_i,z_i) and covariance matrix I_3). What I did is I meshed the space into npoint x npoints x npoints and computed the sum of the probability densities for each of the 'sources' (x_i,y_i,z_i) in each point (x,y,z). Then, if the value I get is big enough (say 95% of the maximum density), I keep the point. otherwise I discard it.

The problem with my code is that it is too slow (many for loops) and the graph I get doesn't look like the one below:

Does anyone know whether there is a package to get a similar plot as the one below?

enter image description here

Upvotes: 1

Views: 217

Answers (1)

Simon M
Simon M

Reputation: 266

Using isosurface we can do reasonably well. (Although I'm not honestly sure what you want, I think this is close:

% Create a path
points = zeros(10,3);
for ii = 2:10
    points(ii, :) = points(ii-1,:) + [0.8 0.04 0] + 0.5 * randn(1,3);
end

% Create the box we're interested in
x = linspace(-10,10);
y = x;
z = x;
[X,Y,Z] = meshgrid(x,y,z);

% Calculate the sum of the probability densities(ish)
V = zeros(size(X));
for ii = 1:10
    V = V + 1/(2*pi)^(3/2) * exp(-0.5 * (((X-points(ii,1)).^2 + (Y-points(ii,2)).^2 + (Z-points(ii,3)).^2)));
end

fv = isosurface(X,Y,Z,V, 1e-4 * 1/(2*pi)^(3/2), 'noshare');
fv2 = isosurface(X,Y,Z,V, 1e-5 * 1/(2*pi)^(3/2), 'noshare'); 
p = patch('vertices', fv.vertices, 'faces', fv.faces);
set(p,'facecolor', 'none', 'edgecolor', 'blue', 'FaceAlpha', 0.05)
hold on;
p2 = patch('vertices', fv2.vertices, 'faces', fv2.faces);
set(p2,'facecolor', 'none', 'edgecolor', 'red', 'FaceAlpha', 0.1)
scatter3(points(:,1), points(:,2), points(:,3));

Upvotes: 1

Related Questions