Another histogram modeling example

This notebook illustrates Gaussian mixture modeling for some 2d synthetic data, and assumes that the definitions in the notebook "statistical_modeling.nb" have already been evaluated.

Global user definitions

Colors representing two classes

Data colors

data1Color = Red;
data2Color = Blue;
centroidColor  =  Magenta;

style1  =  {data1Color, 
          PointSize[ .02]};
style2  =  {data2Color, 
          PointSize[ .02]};
centroidStyle  =  {centroidColor, PointSize[ .02]};

Classification region colors

class1Color  =  RGBColor[1, 1,  .5];
class2Color  =  RGBColor[ .5, 1, 1];

classColors  =  {class1Color,  class2Color};

Default plot range

pr  =  {{0, 1}, {0, 1}};

No frame ticks

nF = FrameTicks -> None;
FS = Sequence[yS,  nF];

Local definitions

Specialized data generation and visualization functions

AnnularRegion[center_,  r1_,  r2_,  n_] := 
        Module[{tst, x1, x2, count, xd, t1, t2, range1,  
            range2}, {x1, x2}  =  
            center; 
range1  =  {x1  -  r2,  
              x1 +  r2}; 
range2  =  {x2  -  r2,  
              x2  +  r2}; 
count  =  
            1; 
xd  =  
            Table[0,  {n}]; 
While[
            count  <=  
              n, {t1, t2}  =  {Random[Real,  range1],  
                Random[Real,  range2]}; 
            If[(tst  = 
                      Sqrt[(t1 - x1)^2  +  (t2 - x2)^2])  <  
                  r2  &&  tst  >  r1, 
              xd[[count++]]  =  {t1, 
                  t2}]
]; 
xd
];

AntiAnnularRegion[center_,  r1_,  r2_,  n_] := 
        Module[{tst, x1, x2, count, xd, t1, t2, range1,  
            range2}, {x1, x2}  =  
            center; 
range1  =  {0,  
              1}; 
range2  =  {0, 
              1}; 
count  =  
            1; 
xd  =  
            Table[0,  {n}]; 
While[
            count  <=  
              n, {t1, t2}  =  {Random[Real,  range1],  
                Random[Real,  range2]}; 
            If[(tst  = 
                      Sqrt[(t1 - x1)^2  +  (t2 - x2)^2])  >  
                  r2  ||  tst  <  r1, 
              xd[[count++]]  =  {t1, 
                  t2}]
]; 
xd
];

Region definitions

c1  =  { .5,  .5};  

    rr1  =   .12;  

    rr2  =   .35;

DrawAnnularRegion2[center_, r1_, r2_, pr_,  color_:  LightPink,  
          color2_:  White, 
          opt___] := 
Module[{xmin, xmax, ymin, 
            ymax}, 
{xmin, xmax, ymin, ymax} = 
            Flatten[pr] + {- .01,  .01, - .01,  .01}; 
Show[
Graphics[{{color, 
                  Disk[center,  r2]}, {color2, Disk[center, r1]}}],  
            AspectRatio -> Automatic,  
            PlotRange -> {{xmin, xmax}, {ymin, ymax}}, 
            Frame -> True, opt]
];

specialClassPlot  := 
        Show[Graphics[{class1Color, 
              Rectangle[{0, 0}, {1, 1}]}],   DrawAnnularRegion2[c1, rr1, 
            rr2,  pr, class2Color,  class1Color, nS] , 
          AspectRatio -> Automatic,  Frame -> True,  nF];

Visualize annular region in 2d

Draw annular region in 2d

DrawAnnularRegion[center_, r1_, r2_, pr_,  color_:  LightPink,  
          opt___] := 
Module[{xmin, xmax, ymin, 
            ymax}, 
{xmin, xmax, ymin, ymax} = 
            Flatten[pr] + {- .01,  .01, - .01,  .01}; 
Show[
Graphics[{{color, 
                  Disk[center,  r2]}, {GrayLevel[1. ], 
                  Disk[center, r1]}}],  AspectRatio -> Automatic,  
            PlotRange -> {{xmin, xmax}, {ymin, ymax}}, 
            Frame -> True, opt]
];

Plot multiple results for annular region

AnnularAllPlot[center_,  r1_,  
          r2_,  {μ_, Σ_, 
            priors_}] := 
Module[{pr, g1, g2, g3, g4, g5, 
            ga, gb}, 
pr  =  {{0, 1}, {0, 
                1}}; 
g1  =   DrawAnnularRegion[center, 
              rr1, rr2,  pr, LightPink,  nF, 
              nS]; 
g2  =  
            PlotMeans[{μ, Σ, priors},  pr,  Blue,  nF, 
              nS]; 
g3  =  
            PlotContourPDFOverlay[{μ, Σ, priors},  pr, 50, 
              nF, nS]; 
g4  =  
            PlotContourPDF[{μ, Σ, priors},  pr, 50, nF, 
              nS]; 
g5  =  
            PlotPDF[{μ, Σ, priors},  pr, 50, 
              Axes -> False, nS]; 
ga  =  
            Show[g1, g2]; 
gb  =  
            Show[g1, g3]; 
Show[
            GraphicsArray[{{ga, gb}, {g4, g5}}], yS,  
            imSize]
];

Data generation and visualization

Idealized data classes

idealPlot  =  specialClassPlot;

Generate data from above distributions

Number of points in each class

n1  =  1000;
n2  =  1000;

Generate data

data1  = AnnularRegion[c1,  rr1,  rr2,  n1];
data2  = AntiAnnularRegion[c1,  rr1,  rr2,  n2];

Visualize data

gd1  =  
        PlotData2D[data1, pr, style1, nS];
gd2  =  
        PlotData2D[data2, pr, style2, nS];
gd  =  Show[gd1, gd2,  FS];

Histogram modeling

Number of centroids for each class ((2^m))

m  =  5;

Compute data-dependent centroids

vqCentroids1 =  VQ[data1, m];
vqCentroids2  = VQ[data2, m];
vqCentroids  =  vqCentroids1~Join~vqCentroids2;

Visualize data

gCentroids =  
        PlotData2D[vqCentroids, pr,  centroidStyle,  
          FrameTicks -> None];

Visualize resultant partitioning of input space

gpart1  =  
        DrawPartitioning[vqCentroids,  pr, nS];
Show[gpart1, gCentroids, yS, FrameTicks -> None];

Histograms for two classes

Compute histograms based on data-dependent centroids

h1  = Histogram[data1, vqCentroids, pr, 150];
h2  = Histogram[data2, vqCentroids, pr, 150];

Visualize two histograms

gh1  = 
        DrawHistogram[h1, pr,  nS, 
          Axes -> False];
gh2  = 
        DrawHistogram[h2, pr,  nS, 
          Axes -> False];
Show[GraphicsArray[{gh1, gh2}], yS,  imSize];

New classification of input space

gdnoFrame  =  
        Show[gd,  Frame -> False,  nS];
gclass  =  
        DrawHistogramClass[{h1, h2},  {class1Color, class2Color}, 
          nS];
Show[GraphicsArray[{gdnoFrame, gclass}], yS,  imSize];