Plotting Bivariate Normal distribution pdf using random numbers generated from the distribution.

I need to plot a bi-variate Pdf by generating 1000 random variables from a given (mu and sigma) of bi-variate Normal distribution. I tried writing a code (Matlab), but there's some error in it. Can anyone please check it and help me correct it.

 function [ ] = bivariatedn( mu,sigma)
  % Generating 1000 random variables
  R=mvnrnd(mu,sigma,1000);
  X=R(:,1);
  Y=R(:,2);
  % Define the grid for visualization
  [X,Y]=meshgrid(R(:,1),R(:,2));

% Define the constant (To write pdf in a easy way)
const = (1/sqrt(2*pi))^2;
const = const/sqrt(det(sigma));
temp = [X(:)-mu(1) Y(:)-mu(2)];

pdf=zeros(length(X),1);

d=zeros(length(X),1);
d=diag(temp*inv(sigma)*temp');
for i=1:length(X)
  pdf(i)= const*exp(-0.5*d(i));
end
% plot the result
surfc(X, Y, pdf, 'LineStyle', 'none');
end

Error using surfc (line 39) The surface Z must contain more than one row or column

This is the error being shown. Can anyone look into this?

Answer 1

[John D'Errico]

To me, the simple way to plot an approximate pdf from a random sample would come from a 2-d histogram.

X = mvnrnd([0 0],eye(2),1000);
histogram2(X(:,1),X(:,2))

Of course, that does not get the z axis scaled properly, so use this version instead.

histogram2(X(:,1),X(:,2),'Normalization','pdf')

That gives you a PDF from a random sample. Are you looking to compute the actual PDF? If you know the mean and variance, then why do you need to generate a random sample? You could just use mvnpdf.

[sarath l]

This is for my assignment. I shouldn't use built-in functions. And I'm looking to compute and plot the actual pdf(contour and 3d plots).

[John D'Errico]

But why are you generating a random sample then? You DID indeed use a built-in function, mvnrnd. It looks like you are using stuff you are not allowed to use as it is. For that matter, you used inv. meshgrid? det? Even exp. So which built-in functions can you use?

[sarath l]

Sorry for the confusion. I shouldn't use the in built functions which directly give me the PDF of bi variate normal distribution. I should construct a function to calculate it(which I did by using "const" and "pdf" variables). I can use other in built functions such as inv, mvnrnd(to generate rv's). The question is, I should generate some (n) number of random variables from a given mu and sigma of bivariate normal distribution and plot the PDF using these rv's and distribution formula.

[John D'Errico]

But why are you using mvnrnd at all here? You just want to plot a pdf. What you did there makes no sense at all.

If you are generating n random variates, then create a histogram using a 2-d bar plot. That is all a histogram is, anyway. A properly scaled and created bar chart.

But if you actually want to create a surface, then using random variables makes little sense. Trying to use meshgrid on a random sample is rather silly, and meaningless. Sorry, but it is. (Yes, you COULD create a surface from randomly scattered points. But that is far beyond the level you are working at from what I can see, involving ideas and methods that you have not even touched upon.)

Use meshgrid. START WITH THAT. What arguments will we feed into meshgrid? It will look something like this:

[X,Y] = meshgrid(what?,what?)

The answer is NOT the output of mvnrnd. The idea is we might want x to vary uniformly over something like this interval:

[mu(1) - 3*sigma(1,1),mu(1) + 3*sigma(1,1)]

We could do the same for y. Generate a pair of uniformly spaced vectors, using linspace. So the code will look like this:

xvec = linspace(mu(1) - 5*sigma(1,1),mu(1) + 5*sigma(1,1),100);
yvec = linspace(mu(2) - 5*sigma(2,2),mu(2) + 5*sigma(2,2),100);
[X,Y] = meshgrid(xvec,yvec);

Now, combine them into one array.

XY = [X(:) - mu(1),Y(:) - mu(2)];

Here is a simple way to write the exponential in a vectorized form. It can be written in a less compact way, or even using a loop. But why?

Z = exp(-1/2*sum(XY*inv(sigma).*XY,2))/(2*pi*sqrt(det(sigma)));
Z = reshape(Z,size(X));

surf(X,Y,Z)

The point of all this is you are somehow conflating a randomly generated sample with the pdf itself. You need to do one or the other, but you should not try to mix things. That will just get you a bit confused. (Which is, I think, what you already are in this matter.)