bivariate_normal.py

# -*- coding: UTF-8 -*-
"""
https://blog.ittraining.com.tw/2024/10/normal-distribution.html

"""


import numpy as np
import matplotlib.pyplot as plt
#from matplotlib.mlab import bivariate_normal


def bivariate_normal(X, Y, sigmax=1.0, sigmay=1.0,
                 mux=0.0, muy=0.0, sigmaxy=0.0):
    """
    Bivariate Gaussian distribution for equal shape *X*, *Y*.
    See `bivariate normal
    <http://mathworld.wolfram.com/BivariateNormalDistribution.html>`_
    at mathworld.
    """
    Xmu = X-mux
    Ymu = Y-muy

    rho = sigmaxy/(sigmax*sigmay)
    z = Xmu**2/sigmax**2 + Ymu**2/sigmay**2 - 2*rho*Xmu*Ymu/(sigmax*sigmay)
    denom = 2*np.pi*sigmax*sigmay*np.sqrt(1-rho**2)
    return np.exp(-z/(2*(1-rho**2))) / denom

def generateData(n, mean, cov):
    """
    generate normal distibution data
    """
    np.random.seed(2033)
    data = np.random.multivariate_normal(mean, cov, size=n)
    return data


def drawData(ax, mu, cov):

    data = generateData(150, mu, cov)
    ax.scatter(data[:, 0], data[:, 1])
    x = np.arange(-10, 10, 0.1)
    
    X, Y = np.meshgrid(x, x) # 依據 meshgird(x,x)產生的點, 共有x*x個點,
                             # 傳回這些點的X及Y, 故 X's shape (x,x), Y's shape(x,x) 
    Z = bivariate_normal(X, Y, cov[0, 0], cov[1, 1], mu[0], mu[1], cov[0, 1])
    ax.contour(X, Y, Z)
    ax.set_xlim([-10, 10])
    ax.set_ylim([-10, 10])
    ax.get_yaxis().set_visible(False)
    ax.get_xaxis().set_visible(False)


def visualize():
    """
    visualize 
    """

    fig = plt.figure(figsize=(10, 10), dpi=80)

    ax = fig.add_subplot(2, 2, 1)
    cov = np.array([[1., 0.], 
                    [0., 1.]])
    mu = np.array([0., 0.])
    
    drawData(ax, mu, cov)
    ax = fig.add_subplot(2, 2, 2)
    cov = np.array([[4., 0.], 
                    [0., 4.]])
    mu = np.array([0., 0.])
    
    drawData(ax, mu, cov)
    ax = fig.add_subplot(2, 2, 3)
    cov = np.array([[4., 3.], 
                    [3., 4.]])
    mu = np.array([0., 0.])
    
    drawData(ax, mu, cov)
    ax = fig.add_subplot(2, 2, 4)
    cov = np.array([[4., -3.], 
                    [-3., 4.]])
    mu = np.array([0., 0.])
   
    drawData(ax, mu, cov)
    plt.show()


if __name__ == "__main__":
    visualize()