%% Variable Aleatoria Bidimensional Discreta
% Copyright (c) 2023 G.T.A.S. (Grupo Tratamiento Avanzado de Señal). Universidad de Cantabria.
% Métodos Matemáticos para Telecomunicaciones.
%
% En este ejemplo se representan gráficamente la función densidad de probabilidad conjunta y la función de distribución conjunta de diversas variables aleatorias bidimensionales discretas. También se representan las distintas marginales y su relación con la conjunta.
%% Generación de las realizaciones de la v.a. bidimensional
% Puede modificar el tipo de variable aleatoria bidimensional que se generará sin más que modificar el código de la función matlab denominada genera_va_bidimensional_discreta que se encuentra al final del Live Script.

N = 100000;     % Número de realizaciones
[x,y] = genera_va_bidimensional_discreta(N);   % Vectores con N realizaciones de la v.a. bidimensional

%% Estimación de la Función Densidad de Probabilidad
% Se estima la fdp a partir del histograma bidimensional. Recuerde que, las verdaderas fdp's estarán compuestas por Deltas de Dirac que aquí se representan y visualizan mediante barras estrechas de altura igual a la constante que las multiplica (es decir, igual a su área y, por lo tanto, a la probabilidad).

histogram2(x,y,[min(x)-1:0.2:max(x)+1],[min(y)-1:0.2:max(y)+1],'Normalization','probability')
xlabel('x');ylabel('y');zlabel('estima f_X_Y(x,y)');title('Función Densidad de Probabilidad Conjunta')

% También es posible representar los histogramas de las variables individuales, es decir, las fdp's marginales.

subplot(211);histogram(x,'Normalization','probability','BinWidth',0.01);
xlabel('x');ylabel('estima f_X(x)');title('Función Densidad de Probabilidad Marginal de X')
subplot(212);histogram(y,'Normalization','probability','BinWidth',0.01);
xlabel('y');ylabel('estima f_Y(y)');title('Función Densidad de Probabilidad Marginal de Y')

% Nótese que la marginal de X se puede calcular como la integral de la fdp conjunta desde y=-\infty hasta y=+\infty para cada valor x, es decir:
% P(X=x_i)=\sum_{y_k \in \Omega_Y} P(X=x_i,Y=y_k)=\sum_{k}p_{ik}
% Observe en la figura inferior (fdp conjunta representada como mapa de calor) cómo la marginal deX para un valor concretox se puede obtener integrando dicha fdp conjunta a lo largo de la línea vertical que pasa por dicho valor de x.

subplot(111);histogram2(x,y,'DisplayStyle','tile','Normalization','probability','BinWidth',[0.1 0.1]);
xlabel('x');ylabel('y');title('Función Densidad de Probabilidad Conjunta');colorbar;axis equal

%% Estimación de la Función de Distribución
% Se estima la Función de Distribución conjunta a partir del histograma acumulativo bidimensional

subplot(111);histogram2(x,y,[min(x)-1:max(x)+1],[min(y)-1:max(y)+1],'Normalization','cdf');
xlabel('x');ylabel('y');zlabel('estima F_X_Y(x,y)');title('Función de Distribución Conjunta')

% También es posible representar los histogramas acumulativos de las variables individuales, es decir, las funciones de distribución marginales.

subplot(211);histogram(x,[min(x)-1:max(x)+1],'Normalization','cdf','DisplayStyle',"stairs");
xlabel('x');ylabel('estima F_X(x)');title('Función de Distribución Marginal de X');grid
subplot(212);histogram(y,[min(y)-1:max(y)+1],'Normalization','cdf','DisplayStyle',"stairs");
xlabel('y');ylabel('estima F_Y(y)');title('Función de Distribución Marginal de Y');grid

% Nótese que la función distribución marginal de X se puede calcular como la función de distribución conjunta particularizada en y=+\infty, es decir, F_X(x)=F_{XY}(x,+\infty)
% Puede comprobarlo girando la figura de la función de distribución conjunta hasta que desdaparezca el eje y y se observe de frente el eje x. En esa posición se observará la la función distribución marginal de X.

return

%
%% Función matlab genera_va_bidimensional_discreta
% Modifique el código de esta función para generar la variable aleatoria bidimensional discreta deseada.

function [x,y] = genera_va_bidimensional_discreta(N)
% Función de generación de variables aleatorias bidimensionales discretas
% Modifique el código para generar la variable aleatoria bidimensional deseada
    
% Ejemplo lanzamiento de dos monedas
    z = randi([0,1],1,N); % lanzamiento moneda 1
    w = randi([0,1],1,N); % lanzamiento moneda 2
    x = z;    % moneda 1
    y = w;    % moneda 2
    %y = z+w;  % moneda 1 + moneda 2

% Ejemplo lanzamiento de dos dados
%     z = randi([1,6],1,N); % lanzamiento dado 1
%     w = randi([1,6],1,N); % lanzamiento dado 2
%     x = abs(z-w);   % Diferencia de puntuación de los dados
%     y = z+w;        % Suma de puntuaciones de dados
        
end