|
|
| Line 1: |
Line 1: |
| The clusters are not computed correctly.
| | {{warning|1=The media wiki content on this page is no longer maintained. The examples presented on the https://itk.org/Wiki/* pages likely require ITK version 4.13 or earlier releases. In many cases, the examples on this page no longer conform to the best practices for modern ITK versions.}} |
| | |
| ==ExpectationMaximizationMixtureModelEstimator_Image.cxx==
| |
| <source lang="cpp">
| |
| #include "itkVector.h"
| |
| #include "itkListSample.h"
| |
| #include "itkGaussianMixtureModelComponent.h"
| |
| #include "itkExpectationMaximizationMixtureModelEstimator.h"
| |
| #include "itkSampleClassifierFilter.h"
| |
| #include "itkMaximumDecisionRule2.h"
| |
| #include "itkImageToListSampleFilter.h"
| |
| #include "itkCovariantVector.h"
| |
| #include "itkImageRegionIterator.h"
| |
| #include "itkImageFileReader.h"
| |
| #include "itkSimpleFilterWatcher.h"
| |
| | |
| typedef itk::CovariantVector<unsigned char, 3> PixelType;
| |
| typedef itk::Image<PixelType, 2> ImageType;
| |
| | |
| static void ControlledImage(ImageType::Pointer image);
| |
| static void RandomImage(ImageType::Pointer image);
| |
| | |
| int main(int argc, char*argv[])
| |
| { | |
| | |
| ImageType::Pointer image = ImageType::New();
| |
| | |
| RandomImage(image);
| |
| //ControlledImage(image);
| |
| | |
| typedef itk::Statistics::ImageToListSampleFilter<ImageType> ImageToListSampleFilterType;
| |
| ImageToListSampleFilterType::Pointer imageToListSampleFilter = ImageToListSampleFilterType::New();
| |
| imageToListSampleFilter->SetInput(image);
| |
| imageToListSampleFilter->Update();
| |
| | |
| unsigned int numberOfClasses = 3;
| |
| | |
| typedef itk::Array< double > ParametersType;
| |
| ParametersType params( numberOfClasses + numberOfClasses*numberOfClasses ); // 3 for means and 9 for 3x3 covariance
| |
| | |
| // Create the first set (for the first cluster/model) of initial parameters
| |
| std::vector< ParametersType > initialParameters( numberOfClasses );
| |
| for(unsigned int i = 0; i < 3; i++)
| |
| {
| |
| params[i] = 5.0; // mean of dimension i
| |
| }
| |
| unsigned int counter = 0;
| |
| for(unsigned int i = 0; i < 3; i++)
| |
| {
| |
| for(unsigned int j = 0; j < 3; j++)
| |
| {
| |
| if(i == j)
| |
| {
| |
| params[3+counter] = 5; // diagonal
| |
| }
| |
| else
| |
| {
| |
| params[3+counter] = 0; // off-diagonal
| |
| }
| |
| counter++;
| |
| }
| |
| }
| |
| | |
| | |
| initialParameters[0] = params;
| |
| | |
| // Create the second set (for the second cluster/model) of initial parameters
| |
| params[0] = 210.0;
| |
| params[1] = 5.0;
| |
| params[2] = 5.0;
| |
| counter = 0;
| |
| for(unsigned int i = 0; i < 3; i++)
| |
| {
| |
| for(unsigned int j = 0; j < 3; j++)
| |
| {
| |
| if(i == j)
| |
| {
| |
| params[3+counter] = 5; // diagonal
| |
| }
| |
| else
| |
| {
| |
| params[3+counter] = 0; // off-diagonal
| |
| }
| |
| counter++;
| |
| }
| |
| }
| |
| initialParameters[1] = params;
| |
| | |
| // Create the third set (for the third cluster/model) of initial parameters
| |
| params[0] = 5.0;
| |
| params[1] = 210.0;
| |
| params[2] = 5.0;
| |
| counter = 0;
| |
| for(unsigned int i = 0; i < 3; i++)
| |
| {
| |
| for(unsigned int j = 0; j < 3; j++)
| |
| {
| |
| if(i == j)
| |
| {
| |
| params[3+counter] = 5; // diagonal
| |
| }
| |
| else
| |
| {
| |
| params[3+counter] = 0; // off-diagonal
| |
| }
| |
| counter++;
| |
| }
| |
| }
| |
| initialParameters[2] = params;
| |
| | |
| std::cout << "Initial parameters: " << std::endl;
| |
| for ( unsigned int i = 0 ; i < numberOfClasses ; i++ )
| |
| {
| |
| std::cout << initialParameters[i] << std::endl;
| |
| }
| |
| | |
| typedef itk::Statistics::GaussianMixtureModelComponent< ImageToListSampleFilterType::ListSampleType >
| |
| ComponentType;
| |
| | |
| std::cout << "Number of samples: " << imageToListSampleFilter->GetOutput()->GetTotalFrequency() << std::endl;
| |
| | |
| // Create the components
| |
| std::vector< ComponentType::Pointer > components;
| |
| for ( unsigned int i = 0 ; i < numberOfClasses ; i++ )
| |
| {
| |
| components.push_back( ComponentType::New() );
| |
| (components[i])->SetSample( imageToListSampleFilter->GetOutput() );
| |
| (components[i])->SetParameters( initialParameters[i] );
| |
| }
| |
|
| |
| | |
| typedef itk::Statistics::ExpectationMaximizationMixtureModelEstimator<
| |
| ImageToListSampleFilterType::ListSampleType > EstimatorType;
| |
| EstimatorType::Pointer estimator = EstimatorType::New();
| |
| | |
| estimator->SetSample( imageToListSampleFilter->GetOutput() );
| |
| estimator->SetMaximumIteration( 200 );
| |
| | |
| itk::Array< double > initialProportions(numberOfClasses);
| |
| initialProportions[0] = 0.33;
| |
| initialProportions[1] = 0.33;
| |
| initialProportions[2] = 0.33;
| |
| | |
| std::cout << "Initial proportions: " << initialProportions << std::endl;
| |
| | |
| estimator->SetInitialProportions( initialProportions );
| |
| | |
| for ( unsigned int i = 0 ; i < numberOfClasses ; i++)
| |
| {
| |
| estimator->AddComponent( components[i]);
| |
| }
| |
| //itk::SimpleFilterWatcher watcher(estimator);
| |
| estimator->Update();
| |
| | |
| // Output the results
| |
| 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: ";
| |
| // Outputs: // mean of dimension 1, mean of dimension 2, covariance(0,0), covariance(0,1), covariance(1,0), covariance(1,1)
| |
| std::cout << " " << estimator->GetProportions()[i] << std::endl;
| |
| }
| |
| | |
| // Display the membership of each sample
| |
| typedef itk::Statistics::SampleClassifierFilter< ImageToListSampleFilterType::ListSampleType > FilterType;
| |
| | |
| typedef itk::Statistics::MaximumDecisionRule2 DecisionRuleType;
| |
| DecisionRuleType::Pointer decisionRule = DecisionRuleType::New();
| |
| | |
| typedef FilterType::ClassLabelVectorObjectType ClassLabelVectorObjectType;
| |
| typedef FilterType::ClassLabelVectorType ClassLabelVectorType;
| |
| | |
| ClassLabelVectorObjectType::Pointer classLabelsObject = ClassLabelVectorObjectType::New();
| |
| ClassLabelVectorType & classLabelVector = classLabelsObject->Get();
| |
| | |
| typedef FilterType::ClassLabelType ClassLabelType;
| |
| | |
| ClassLabelType class0 = 0;
| |
| classLabelVector.push_back( class0 );
| |
| | |
| ClassLabelType class1 = 1;
| |
| classLabelVector.push_back( class1 );
| |
| | |
| ClassLabelType class2 = 2;
| |
| classLabelVector.push_back( class2 );
| |
| | |
| FilterType::Pointer sampleClassifierFilter = FilterType::New();
| |
| sampleClassifierFilter->SetInput( imageToListSampleFilter->GetOutput() );
| |
| sampleClassifierFilter->SetNumberOfClasses( numberOfClasses );
| |
| sampleClassifierFilter->SetClassLabels( classLabelsObject );
| |
| sampleClassifierFilter->SetDecisionRule( decisionRule );
| |
| sampleClassifierFilter->SetMembershipFunctions( estimator->GetOutput() );
| |
| sampleClassifierFilter->Update(); | |
| | |
| const FilterType::MembershipSampleType* membershipSample = sampleClassifierFilter->GetOutput();
| |
| FilterType::MembershipSampleType::ConstIterator iter = membershipSample->Begin();
| |
| | |
| while ( iter != membershipSample->End() )
| |
| {
| |
| std::cout << (int)iter.GetMeasurementVector()[0] << " " << (int)iter.GetMeasurementVector()[1] << " " << (int)iter.GetMeasurementVector()[2]
| |
| << " : " << iter.GetClassLabel() << std::endl;
| |
| ++iter;
| |
| }
| |
| | |
| return EXIT_SUCCESS;
| |
| }
| |
| | |
| void ControlledImage(ImageType::Pointer image)
| |
| {
| |
| // Create an image
| |
| ImageType::RegionType region;
| |
| ImageType::IndexType start;
| |
| start[0] = 0;
| |
| start[1] = 0;
| |
| | |
| ImageType::SizeType size;
| |
| size[0] = 10;
| |
| size[1] = 10;
| |
| | |
| region.SetSize(size);
| |
| region.SetIndex(start);
| |
| | |
| image->SetRegions(region);
| |
| image->Allocate();
| |
| | |
| // Make a red and a green square
| |
| itk::CovariantVector<unsigned char, 3> green;
| |
| green[0] = 0;
| |
| green[1] = 255;
| |
| green[2] = 0;
| |
| | |
| itk::CovariantVector<unsigned char, 3> red;
| |
| red[0] = 255;
| |
| red[1] = 0;
| |
| red[2] = 0;
| |
| | |
| itk::CovariantVector<unsigned char, 3> black;
| |
| black[0] = 0;
| |
| black[1] = 0;
| |
| black[2] = 0;
| |
| | |
| itk::ImageRegionIterator<ImageType> imageIterator(image,region);
| |
| imageIterator.GoToBegin();
| |
| | |
| while(!imageIterator.IsAtEnd())
| |
| {
| |
| if(imageIterator.GetIndex()[0] > 2 && imageIterator.GetIndex()[0] < 5 &&
| |
| imageIterator.GetIndex()[1] > 2 && imageIterator.GetIndex()[1] < 5)
| |
| {
| |
| imageIterator.Set(green);
| |
| }
| |
| else if(imageIterator.GetIndex()[0] > 6 && imageIterator.GetIndex()[0] < 9 &&
| |
| imageIterator.GetIndex()[1] > 6 && imageIterator.GetIndex()[1] < 9)
| |
| {
| |
| imageIterator.Set(red);
| |
| }
| |
| else
| |
| {
| |
| imageIterator.Set(black);
| |
| }
| |
| ++imageIterator;
| |
| }
| |
| | |
| }
| |
| | |
| void RandomImage(ImageType::Pointer image)
| |
| {
| |
| // Create an image
| |
| ImageType::RegionType region;
| |
| ImageType::IndexType start;
| |
| start[0] = 0;
| |
| start[1] = 0;
| |
| | |
| ImageType::SizeType size;
| |
| size[0] = 10;
| |
| size[1] = 10;
| |
| | |
| region.SetSize(size);
| |
| region.SetIndex(start);
| |
| | |
| image->SetRegions(region);
| |
| image->Allocate();
| |
| | |
| itk::ImageRegionIterator<ImageType> imageIterator(image,region);
| |
| imageIterator.GoToBegin();
| |
| | |
| while(!imageIterator.IsAtEnd())
| |
| {
| |
| // Get a random color
| |
| itk::CovariantVector<unsigned char, 3> pixel;
| |
| pixel[0] = rand() * 255;
| |
| pixel[1] = rand() * 255;
| |
| pixel[2] = rand() * 255;
| |
| imageIterator.Set(pixel);
| |
| ++imageIterator;
| |
| }
| |
| | |
| }
| |
| </source>
| |
| | |
| {{ITKCMakeLists|{{SUBPAGENAME}}}}
| |
