Conditioning bivariate gaussian distribution
显示 更早的评论
Hi,
I have a bivariate normal distribution as follows f(x):
m = 0;
c = [0.5 0.8; 0.8 2.0];
x1 = -4:0.2:4;
x2 = -4:0.2:4;
[X1, X2] = meshgrid(x1,x2);
X = X1(:)';
Y = X2(:)';
fun = @(X, Y) 1/(2*pi*(det(c))^(0.5))* exp(-(0.5)*sum(([X; Y]-m).*(inv(c)*([X; Y]-m))));
i = @(X)integral(@(Y)fun(X,Y),-inf,inf,'ArrayValued',true);
fplot(i)
And this gives output:
I want to find the f(x / y = 1.5)
I have tried to find f(x) and f(y) and then filter out y = 1.5 to get the x values from the distribution. but this method is not working and giving errors as follows:
mu_x = 0;
c = [0.5 0.8; 0.8 2.0];
x1 = -4:0.2:4;
x2 = -4:0.2:4;
[X1,X2] = meshgrid(x1,x2);
X = X1(:)';
Y = X2(:)';
fun = @(X, Y) 1/(2*pi*(det(c))^(0.5))* exp(-(0.5)*sum(([X; Y]-mu_x).*(inv(c)*([X; Y]-mu_x))));
px = @(X)integral(@(Y)fun(X,Y),-inf,inf,'ArrayValued',true);
py = @(Y)integral(@(X)fun(X,Y),-inf,inf,'ArrayValued',true);
%vq1 = interp1(-3:0.2:3,px,-3:0.2:3)
px = px([-3:0.2:3])
p = [px([-3:0.2:3]); py([-3:0.2:3])]
fplot(py)
The error :
Error using vertcat
Dimensions of arrays being concatenated are not consistent.
Error in pg2>@(X,Y)1/(2*pi*(det(c))^(0.5))*exp(-(1/2)*sum(([X;Y]-mu_x).*(inv(c)*([X;Y]-mu_x)))) (line 109)
fun = @(X, Y) 1/(2*pi*(det(c))^(0.5))* exp(-(1/2)*sum(([X; Y]-mu_x).*(inv(c)*([X; Y]-mu_x))));
Error in pg2>@(Y)fun(X,Y) (line 110)
px = @(X)integral(@(Y)fun(X,Y),-inf,inf,'ArrayValued',true);
Error in integralCalc/iterateArrayValued (line 156)
fxj = FUN(t(1)).*w(1);
Error in integralCalc/vadapt (line 130)
[q,errbnd] = iterateArrayValued(u,tinterval,pathlen);
Error in integralCalc (line 103)
[q,errbnd] = vadapt(@minusInfToInfInvTransform,interval);
Error in integral (line 88)
Q = integralCalc(fun,a,b,opstruct);
Error in pg2>@(X)integral(@(Y)fun(X,Y),-inf,inf,'ArrayValued',true) (line 110)
px = @(X)integral(@(Y)fun(X,Y),-inf,inf,'ArrayValued',true);
Error in pg2>dist (line 113)
px = px([-4:0.2:4])
How do i get the values f(x) and f(y) from px and py in range -4:0.2:4 so that i can find f(x / y = 1.5)?
采纳的回答
Bruno Luong
2019-9-15
编辑：Bruno Luong
2019-9-15
There a few issues with your code.
First INTEGRAL is integration of a scalar function, you cannot integrate a vector function. So in 1D you need to loop in the parameters (x) or (y).
Second parameters ans variable are no longer have the same length, you can not do vercat them.
Here is a modified code
mu_x = 0;
c = [0.5 0.8; 0.8 2.0];
% this block is irrelevant!!!
% x1 = -4:0.2:4;
% x2 = -4:0.2:4;
% [X1,X2] = meshgrid(x1,x2);
% X = X1(:)';
% Y = X2(:)';
fun = @(X, Y) gfun(X,Y,mu_x,c); % nested defined bellow
px = @(X) arrayfun(@(x) integral(@(y)fun(x,y),-inf,inf,'ArrayValued',true), X);
py = @(Y) arrayfun(@(y) integral(@(x)fun(x,y),-inf,inf,'ArrayValued',true), Y);
Pxy = [px([-3:0.2:3]); py([-3:0.2:3])]
close all
plot(Pxy')
%%
function f = gfun(X,Y,mu,c)
[X,Y] = ndgrid(X(:)',Y(:)'); % expanding the scalar to match the vector, regardless which is which
XY = [X; Y];
f = 1/(2*pi*(det(c))^(0.5))* exp(-(0.5)*sum((XY-mu).*(c\(XY-mu))));
end
6 个评论
Aishwarya Radhakrishnan
2019-9-15
编辑：Aishwarya Radhakrishnan
2019-9-15
Bruno Luong
2019-9-15
That's your pdf restricted on the line of equation
x = 1.5*y
Aishwarya Radhakrishnan
2019-9-16
Aishwarya Radhakrishnan
2019-9-16
Hi Bruno,
I have used the above methods and tried to find the conditional independence. Is this the correct procedure?
function [] = BivariateGaussianDistributions_ConditionalProbability()
fig = figure('Name','Marginalize','NumberTitle','off');
mu = 0;
c = [0.5 0.8; 0.8 2.0];
fun = @(X, Y) distfun(X,Y,mu,c); % nested defined bellow
py = @(Y) arrayfun(@(y) integral(@(x)fun(x,y),-inf,inf,'ArrayValued',true), Y);
Py2 = py([1.5]); % P(y=1.5)
[X1,X2] = meshgrid((1.5).*[-4:0.02:4],[-4:0.02:4]); % x = 1.5*y
X = [X1(:) X2(:)]';
PXY = 1/(2*pi*(det(c))^(0.5))* exp(-(1/2)*sum((X-mu).*(pinv(c)*(X-mu))));
PXY = reshape(PXY,size(X1)); % P(X,Y = 2)
PX_Y2 = PXY / Py2; % P(X | Y = 2) = P(X,Y = 2) / P(y=1.5)
plot(X1,PX_Y2)
title('Approximate Univariate Normal Distribution');
xlabel('x');
ylabel('Approximate p(x)');
end
function f = distfun(X,Y,mu,c)
[X,Y] = ndgrid(X(:)',Y(:)'); % expanding the scalar to match the vector, regardless which is which
XY = [X; Y];
f = 1/(2*pi*(det(c))^(0.5))* exp(-(0.5)*sum((XY-mu).*(c\(XY-mu))));
end
For this i'm getting the following curve:
Bruno Luong
2019-9-17
编辑：Bruno Luong
2019-9-17
No idea what you compute in your script and what want to plot.
The condition probability is simply the PDF projected on the line
m = 0;
c = [0.5 0.8; 0.8 2.0];
fun = @(X, Y) 1/(2*pi*(det(c))^(0.5))* exp(-(0.5)*sum(([X; Y]-m).*(c\([X; Y]-m))));
fun_condi_x = @(x1) fun(x1,1.5*x1);
x1_c = -4:0.05:4;
close all
plot(x1_c,fun_condi_x(x1_c)/integral(fun_condi_x,-Inf,Inf));
xlabel('x1');
ylabel('PDF of f(x1,x2) conditional x2/x1=1.5')
x1 = -4:0.05:4;
x2 = -4:0.05:4;
[X1, X2] = meshgrid(x1,x2);
Z = fun(X1(:)',X2(:)');
Z = reshape(Z,size(X1));
% First figure plot is the restriction of Z on the white line
% (normalized so that the integral == 1, since it's a PDF if we consider x1 as
% an independent parametrization variable.)
figure
imagesc(x1,x2,Z);
set(gca,'Ydir','normal')
hold on
x2_c = 1.5*x1_c;
plot(x1_c,x2_c,'w','linewidth',2);
Aishwarya Radhakrishnan
2019-9-18
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Uniform Distribution (Continuous) 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
