<jupyter_start><jupyter_text>--- Example of classification using two gaussian distributions.---In this example, we generate a global probability distribution function (PDF) that is the sum of two gaussian PDF. Each one is weighted with its a priori class probability $P(C_{i})$:  $P(\bf{x}) = P(\bf{x}|C_{1}) P(C_{1}) + P(\bf{x}|C_{2}) P(C_{2})$ Our goal is to locate the influence zone of each class i. This corresponds to the zone where  $P(\bf{x}|C_{i}) P(C_{i}) > P(\bf{x}|C_{j}) P(C_{j})$  The 2-D gausian PDFIn its matrix-form, the equation for the 2-D gausian PDF reads like this:   $P(\bf{x}) = \frac{1}{2\pi |\Sigma|^{0.5}} \exp{[-\frac{1}{2}(\bf{x}-\bf{\mu})^\top \Sigma^{-1} (\bf{x}-\bf{\mu})]}$  where  $\begin{align}\bf{x} &= [x_{1} x_{2}]^\top \\\bf{\mu} &= [\mu_{1} \mu_{2}]^\top  \\\Sigma &= \begin{pmatrix} \sigma_{x_{1}}^2 & \rho\sigma_{x_{1}}\sigma_{x_{2}} \\ \rho\sigma_{x_{1}}\sigma_{x_{2}} & \sigma_{x_{2}}^2 \end{pmatrix}\end{align}$                where $\rho$ is the correlation factor between the $x_{1}$ and $x_{2}$ data.  A useful trick to generate a covariance matrix $\Sigma$ with desired features Instead of guessing the values of $\sigma_{x_{1}}$, $\sigma_{x_{2}}$ and $\rho$ to create a covariance matrix $\Sigma$ for visualization purpose, we can design it with some desired characteristics. Those are the principal axis variances $\sigma_{1}^2$ and $\sigma_{2}^2$, and the rotation angle $\theta$.First, we generate the covariance matrix of a correlated PDF with its principal axes oriented along the $x_{1}$ and $x_{2}$ axes.  $\Sigma_{PA} = \begin{pmatrix} \sigma_{1}^2 & 0 \\ 0 & \sigma_{2}^2 \end{pmatrix}$ Next, we generate the rotation matrix for the angle $\theta$:  $R = \begin{pmatrix} \cos{\theta}  & -\sin{\theta} \\ \sin{\theta} & \cos{\theta} \end{pmatrix}$ The covariance matrix we are looking for is   $\Sigma = R \Sigma_{PA} R^\top $ Hence, we only have to specify the values of $\sigma_{1}$ and $\sigma_{2}$ along principal axes and the rotation angle $\theta$. The correlation coefficient $\rho$ depends on the value of $\theta$:  $\rho>0$ when $\theta>0$   $\rhoN.B. For ease of reading, we use below the variables x and y instead of x1 and x2. The final results are shown with x1 and x2.<jupyter_code>print(__doc__)

# Author: <NAME> <<EMAIL>>
# License: BSD

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import colors
from matplotlib import cm
from matplotlib.ticker import NullFormatter
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from scipy.stats import multivariate_normal

from skimage.transform import resize

import seaborn as sns
sns.set(color_codes=True)

# Used for reproductibility of the results
np.random.seed(43)<jupyter_output>Automatically created module for IPython interactive environment<jupyter_text>Here are the characteristics of the 2-D PDF we want to generate:<jupyter_code># Origins of the PDF
Mu = np.zeros((2,2))
Mu[0,:] = [6., 3.]
Mu[1,:] = [6., 7.]

# Standard deviations along the x1 and x2 directions (or x and y)
sigma = np.zeros((2,2))
sigma[0,:] = [2., 0.7]
sigma[1,:] = [1.2, 0.7]

# Rotation angles for the principal axis
theta = np.array([-30., 0.]) 

# A priori class probabilities
prob_C = np.array([0.5, 0.5])<jupyter_output><empty_output><jupyter_text>Compute the covariance matrix $\Sigma$<jupyter_code>Sigma = np.zeros((2,2,2))
for i in range(2):
    angle = np.radians(theta[i])
    c, s = np.cos(angle), np.sin(angle)

    # Rotation matrix
    R = np.array(((c, -s), (s, c)))        

    # Covariance matrix for a PDF with its principal axes oriented along the x and y directions
    sig = np.array([[sigma[i,0]**2, 0.],[0., sigma[i,1]**2]])

    # Covariance matrix after rotation
    Sigma[:,:,i] = R.dot( sig.dot(R.T) )<jupyter_output><empty_output><jupyter_text>Global PDF computationGenerate a spatial grid where the PDF will be evaluated locally.<jupyter_code>x_min, x_max = 0., 10.
y_min, y_max = 0., 10.
nx, ny = 300, 300
x = np.linspace(x_min, x_max, nx)
y = np.linspace(y_min, y_max, ny)
xx, yy = np.meshgrid(x,y) 
pos = np.dstack((xx, yy))<jupyter_output><empty_output><jupyter_text>Compute the global PDF<jupyter_code># Weighted PDF
model = multivariate_normal(Mu[0,:], Sigma[:,:,0]) 
pdf0 = prob_C[0]*model.pdf(pos)

model = multivariate_normal(Mu[1,:], Sigma[:,:,1]) 
pdf1 = prob_C[1]*model.pdf(pos)

# Global PDF 
pdf = pdf0 + pdf1

# Define a mask to identify the influence zone of each class: Mask=False for class 0 and Mask=True for class 1
mask = pdf1 > pdf0<jupyter_output><empty_output><jupyter_text>Display the global PDF<jupyter_code>z_offset = -0.1
view=[22., -30.]


fig = plt.figure(figsize = (15,10))
ax = fig.gca(projection='3d')


# ---------Generate overlying 3-D surface with colors indicating the influence zone of each class
rstride, cstride = 5, 5

s = ax.plot_surface(xx, yy, pdf, rstride=rstride, cstride=cstride, linewidth=.5, antialiased=True, color='gray', edgecolors='k')       
a1 = s.__dict__['_original_facecolor']
b1 = s.__dict__['_facecolors']
c1 = s.__dict__['_facecolors3d']

s = ax.plot_surface(xx, yy, pdf, rstride=rstride, cstride=cstride, linewidth=.5, antialiased=True, color='w', edgecolors='k')
a2 = s.__dict__['_original_facecolor']
b2 = s.__dict__['_facecolors']
c2 = s.__dict__['_facecolors3d']

Lx = int(nx/rstride)
Ly = int(ny/cstride)

mask2 = resize(mask, (Lx,Ly), order=0)
indx = np.argwhere(mask2)
idx = indx[:,0]*Lx + indx[:,1]

a = a1
b = b1
c = c1
for i in idx:
    a[i,:] = a2[i,:]
    b[i,:] = b2[i,:]
    c[i,:] = c2[i,:]
    
s.__dict__['_original_facecolor'] = a
s.__dict__['_facecolors'] = b
s.__dict__['_facecolors3d'] = c


# --------- Generate underlying contours with borders between each class
cset = ax.contourf(xx, yy, pdf, zdir='z', offset=z_offset, cmap='viridis')
ax.contour(xx, yy, mask, [0.5], offset=z_offset, linewidths=2., colors='white') 

# N.B. For convenience, we use classes 1 and 2 instead of classes 0 and 1 in the figure
ax.text(2, 1, z_offset,'$\mathcal{R}_{1}$', zdir='y', color='w',fontsize=18)        
ax.text(9, 8, z_offset,'$\mathcal{R}_{2}$', zdir='y', color='w',fontsize=18)    
ax.text(6, 2.5, 0.05, '$p(\mathbf{X}|C_{1})p(C_{1})$', zdir='y', horizontalalignment='right', verticalalignment='bottom', 
        color='k',fontsize=18)    
ax.text(6, 6., 0.08, '$p(\mathbf{X}|C_{2})p(C_{2})$', zdir='y', horizontalalignment='right', verticalalignment='bottom', 
        color='gray',fontsize=18)    

ax.set_zlim(z_offset, np.max(pdf))
ax.view_init(view[0], view[1])

ax.set_ylabel('$x_{2}$', fontsize=18)
ax.xaxis.set_rotate_label(False)  
ax.set_xlabel('$x_{1}$', rotation=10, fontsize=18)
ax.set_zlabel('$Z$', rotation=0, fontsize=18)
ax.set_title('$Z = PDF(x_{1},x_{2})$', fontsize=20)
    
fig.tight_layout()

plt.savefig('Influence_zones_in_gaussian_PDF_classification.png')
plt.savefig('Influence_zones_in_gaussian_PDF_classification.pdf')

plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1734752589182969, 'fraction_numerical': 0.02977560414269275, 'mean_word_length': 4.694512694512695, 'pattern_counts': {'":': 0, '<': 27, '<?xml version=': 0, '>': 31, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '646459', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 2680, 'n_tokens_neox': 2550, 'n_words': 712, 'path': 'generate_explanatory_diagram_for_influence_zones_in_two_class_segmentation.ipynb', 'repo_name': 'AstroPierre/Scripts-for-figures-courses-GIF-4101-GIF-7005', 'size': '6979'}