Distribute Sampling Using GMM EM#

Warning

Fix Problem Contains problems not fixed from original wiki.

Synopsis#

Compute distributions of samples using GMM EM.

Results#

Note

Help Wanted Implementation of Results for sphinx examples containing this message. Reconfiguration of CMakeList.txt may be necessary. Write An Example <https://itk.org/ITKExamples/Documentation/Contribute/WriteANewExample.html>

Code#

C++#

#include "itkVector.h"
#include "itkListSample.h"
#include "itkGaussianMixtureModelComponent.h"
#include "itkExpectationMaximizationMixtureModelEstimator.h"
#include "itkNormalVariateGenerator.h"

int
main()
{
  unsigned int numberOfClasses = 2;
  using MeasurementVectorType = itk::Vector<double, 1>;
  using SampleType = itk::Statistics::ListSample<MeasurementVectorType>;
  auto sample = SampleType::New();

  using NormalGeneratorType = itk::Statistics::NormalVariateGenerator;
  auto normalGenerator = NormalGeneratorType::New();

  normalGenerator->Initialize(101);

  MeasurementVectorType mv;
  double                mean = 100;
  double                standardDeviation = 30;
  for (unsigned int i = 0; i < 10; ++i)
  {
    mv[0] = (normalGenerator->GetVariate() * standardDeviation) + mean;
    std::cout << "m[" << i << "] = " << mv[0] << std::endl;
    sample->PushBack(mv);
  }

  normalGenerator->Initialize(3024);
  mean = 200;
  standardDeviation = 30;
  for (unsigned int i = 0; i < 10; ++i)
  {
    mv[0] = (normalGenerator->GetVariate() * standardDeviation) + mean;
    std::cout << "m[" << i << "] = " << mv[0] << std::endl;
    sample->PushBack(mv);
  }

  using ParametersType = itk::Array<double>;
  ParametersType params1(2);

  std::vector<ParametersType> initialParameters(numberOfClasses);
  params1[0] = 110.0;
  params1[1] = 50.0;
  initialParameters[0] = params1;

  ParametersType params2(2);
  params2[0] = 210.0;
  params2[1] = 50.0;
  initialParameters[1] = params2;

  using ComponentType = itk::Statistics::GaussianMixtureModelComponent<SampleType>;

  std::vector<ComponentType::Pointer> components;
  for (unsigned int i = 0; i < numberOfClasses; ++i)
  {
    components.push_back(ComponentType::New());
    components[i]->SetSample(sample);
    components[i]->SetParameters(initialParameters[i]);
  }

  using EstimatorType = itk::Statistics::ExpectationMaximizationMixtureModelEstimator<SampleType>;
  auto estimator = EstimatorType::New();

  estimator->SetSample(sample);
  estimator->SetMaximumIteration(500);

  itk::Array<double> initialProportions(numberOfClasses);
  initialProportions[0] = 0.5;
  initialProportions[1] = 0.5;

  estimator->SetInitialProportions(initialProportions);

  for (unsigned int i = 0; i < numberOfClasses; ++i)
  {
    estimator->AddComponent((ComponentType::Superclass *)components[i].GetPointer());
  }

  estimator->Update();

  for (unsigned int i = 0; i < numberOfClasses; ++i)
  {
    std::cout << "Cluster[" << i << "]" << std::endl;
    std::cout << "    Parameters:" << std::endl;
    std::cout << "         " << components[i]->GetFullParameters() << std::endl;
    std::cout << "    Proportion: ";
    std::cout << "         " << estimator->GetProportions()[i] << std::endl;
  }

  return EXIT_SUCCESS;
}

Classes demonstrated#

template<typename TSample> class ExpectationMaximizationMixtureModelEstimator : public itk::Object

This class generates the parameter estimates for a mixture model using expectation maximization strategy.

The first template argument is the type of the target sample data. This estimator expects one or more mixture model component objects of the classes derived from the MixtureModelComponentBase. The actual component (or module) parameters are updated by each component. Users can think this class as a strategy or a integration point for the EM procedure. The initial proportion (SetInitialProportions), the input sample (SetSample), the mixture model components (AddComponent), and the maximum iteration (SetMaximumIteration) are required. The EM procedure terminates when the current iteration reaches the maximum iteration or the model parameters converge.

Recent API changes: The static const macro to get the length of a measurement vector, MeasurementVectorSize has been removed to allow the length of a measurement vector to be specified at run time. It is now obtained at run time from the sample set as input. Please use the function GetMeasurementVectorSize() to get the length.

See

MixtureModelComponentBase, GaussianMixtureModelComponent

ITK Sphinx Examples:

All ITK Sphinx Examples

2D Gaussian Mixture Model Expectation Maximum

Distribution Of Pixels Using GMM EM

Distribute Sampling Using GMM EM

See itk::Statistics::ExpectationMaximizationMixtureModelEstimator for additional documentation.