# How to find Sigma Matrix for Bivariate Distribution Given Data

7 views (last 30 days)
Thomas Rodriguez on 10 Apr 2022
Commented: Thomas Rodriguez on 12 Apr 2022
Hi,
I'm trying to construct a Bivariate Normal Distribution with Given Data the Location (Lattiude, Longitude).
So far, I have the following code:
xc = 32.7419; % Center of Latitude Coordinate
yc = -117.0904; % Center of Longtitude Coordinate
pop = 72994; % Population at Specified Location
% Creating a 2 dimensional space
x1 = linspace(xc-0.05,xc+0.05,31);
x2 = linspace(yc-0.05,yc+0.05,31);
[X1,X2] = meshgrid(x1,x2);
X = [X1(:) X2(:)];
% Estabilishing Mean Vector,mu & Covariance Matrix,Sigma
mu = [xc yc]
[muHat, Sigma] = normfit(X);
S1 = Sigma(1);
S2 = Sigma(2);
Sigma = [S1 0; 0 S2]
% Sigma = [0.25 0.3; 0.3 1];
% Sigma = [0.85 0.1; 0.1 0.85];
% Evaluating pdf of the normal distribution:
Y = mvnpdf(X,mu,Sigma)*pop;
Y = reshape(Y,length(x2),length(x1));
% Plotting:
figure(1);
surf(x1,x2,Y)
xlabel('x1')
ylabel('x2')
zlabel('Probability Density')
title("Example")
The code produces a distribution, but I would like to create a distribution with minimal standard deviation so that it resembles a "typical" Bivariate Normal Distribution. But I'm having a hard finding a Sigma Matrix that does such.

Jeff Miller on 11 Apr 2022
You can't really use normfit to estimate sigma in the manner you are doing it, especially because x1 & x2 have such different numerical values (i.e., the standard deviation of your x = [X1(:) X2(:)] vector is far more than the sd of either vector on its own).
I suggest getting rid of normfit and doing something like this:
x1range = 0.1;
x2range = 0.1;
x1 = linspace(xc-x1range/2,xc+x1range/2,31);
x2 = linspace(yc-x2range/2,yc+x2range/2,31);
[X1,X2] = meshgrid(x1,x2);
X = [X1(:) X2(:)];
mu = [xc yc]
% 2 notes about the following line:
% (a) the constant 8 reflects the fact that the range of a normal
% distribution is about 8 sd's, from Z=-4 to Z=+4.
% Admittedly, 6 might be a better choice Z=-3 to Z=+3.
% (b) mvnpdf wants variances, not sd's
Sigma = [(x1range/8)^2 0; 0 (x2range/8)^2];
Y = mvnpdf(X,mu,Sigma)*pop;
Thomas Rodriguez on 12 Apr 2022
Thank you so much. It worked!
xmin = min(Lat);
xmax = max(Lat);
ymin = min(Lon);
ymax = max(Lon);
x = linspace(xmin,xmax,31);
y = linspace(ymin,ymax,31);
[X1,X2] = meshgrid(x,y);
X = [X1(:) X2(:)];
% Mu:
mu = [xc yc];
% Sigma:
Sigma = [(x1range/8)^2 0; 0 (x2range/8)^2];
% Recompute PDF:
Y = mvnpdf(X,mu,Sigma)*pop;
Y = reshape(Y,length(x),length(y));
figure(2)
surf(x,y,Y)
xlabel('x1')
ylabel('x2')
zlabel('Probability Density')
title("Example")

R2020b

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!