<html><head><style type="text/css"><!-- DIV {margin:0px;} --></style></head><body><div style="font-family:times new roman, new york, times, serif;font-size:12pt"><DIV>Dear All:</DIV>
<DIV> </DIV>
<DIV>I just figure out the problem. I need to use itk::GaussianMixtureModelComponent->Evaluate to make decision. The result is ok.</DIV>
<DIV> </DIV>
<DIV>but i am still not clear what happened if I used itk::GaussianDensityFunction->Evaluate() to get decision. Maybe I did not use it properly. Still interested to know what happened.</DIV>
<DIV> </DIV>
<DIV>Can somebody tell me whether my code is wrong or I did not used itk::GaussianDensityFunction->Evaluate() properly?</DIV>
<DIV> </DIV>
<DIV>Thanks</DIV>
<DIV> </DIV>
<DIV>Baoyun<BR></DIV>
<DIV style="FONT-SIZE: 12pt; FONT-FAMILY: times new roman, new york, times, serif"><BR>
<DIV style="FONT-SIZE: 12pt; FONT-FAMILY: times new roman, new york, times, serif"><FONT face=Tahoma size=2>
<HR SIZE=1>
<B><SPAN style="FONT-WEIGHT: bold">From:</SPAN></B> Baoyun Li <baoyun_li123@yahoo.com><BR><B><SPAN style="FONT-WEIGHT: bold">To:</SPAN></B> Luis Ibanez <luis.ibanez@kitware.com><BR><B><SPAN style="FONT-WEIGHT: bold">Cc:</SPAN></B> insight-users@itk.org<BR><B><SPAN style="FONT-WEIGHT: bold">Sent:</SPAN></B> Tuesday, March 17, 2009 12:24:04 PM<BR><B><SPAN style="FONT-WEIGHT: bold">Subject:</SPAN></B> cannot get proper classification using itk:GaussianDensityFunction after EM clustering<BR></FONT><BR>
<DIV style="FONT-SIZE: 12pt; FONT-FAMILY: times new roman, new york, times, serif">
<DIV>Dear ALL:</DIV>
<DIV> </DIV>
<DIV>I am using EM gassian mixutrue estimator to estimate the mean and covariance of vector image. I first adapt the image to sample, and the EM works well.</DIV>
<DIV>It produced mean value for Gray matter , whithe matter and CSF, the mean value looks ok to me.</DIV>
<DIV> </DIV>
<DIV>But when I try to assign the label to the image, I cannot get the proper labled image. I got the estimated all the voxle to the one class. And the pdf of GassianDensityFunction always have the biggest value for CSF, that should not happend.</DIV>
<DIV> </DIV>
<DIV>Below is part of my code.</DIV>
<DIV> </DIV>
<DIV>1, I first get the pararameter of each component by : finalparameters=(components[i])->GetFullParameters();<BR>2, Then I decompose the parameters to estimatedmean and estimatedcovariance</DIV>
<DIV>3, set the mean for GaussianDensityFunction by : membershipFunctions[i]->SetMean(&estimatedmean );<BR> membershipFunctions[i]->SetCovariance( &estimatedcovariance);</DIV>
<DIV> </DIV>
<DIV>4, In the lableing program, I used eveluateresult[i]=(membershipFunctions[i]->Evaluate(mv)) to get pdf for each class, mv is the measurement vector obtained of images by itk iterator</DIV>
<DIV>5, then I find the class with biggest pdf, assign the label .</DIV>
<DIV> </DIV>
<DIV>The code is listed below, can some body help me to find what is wrong.</DIV>
<DIV> </DIV>
<DIV> </DIV>
<DIV>Thanks</DIV>
<DIV> </DIV>
<DIV>Baoyun</DIV>
<DIV> </DIV>
<DIV> </DIV>
<DIV>///////give the mean value and covariance to GaussianDensityFunction </DIV>
<DIV> typedef itk::Statistics::GaussianDensityFunction<ComponentType::MeasurementVectorType > MembershipFunctionType;<BR> std::vector< MembershipFunctionType::Pointer > membershipFunctions(numberOfClasses);<BR> unsigned int paramIndex;<BR> ParametersType finalparameters((NumberOfComponents+1)*NumberOfComponents);</DIV>
<DIV> typedef MembershipFunctionType::MeanType MeanType;<BR> typedef MembershipFunctionType::CovarianceType CovarianceType;<BR> MeanType estimatedmean(NumberOfComponents);<BR> CovarianceType estimatedcovariance;</DIV>
<DIV> typedef ComponentType::MeasurementVectorSizeType MeasurementVectorSizeType;<BR> MeasurementVectorSizeType measurementVectorSize=NumberOfComponents;<BR> estimatedcovariance.SetSize( measurementVectorSize, measurementVectorSize );<BR> //estimatedmean.SetSize( measurementVectorSize, 1);<BR> // estimatedmean.Fill(0.0) ;<BR> for ( unsigned int i = 0 ; i < numberOfClasses ; i++ )<BR> {</DIV>
<DIV> // membershipFunctions.push_back(MembershipFunctionType::New());<BR> membershipFunctions[i]=MembershipFunctionType::New();<BR> std::cout << " Size of Gaussian : ";<BR> std::cout << " " << (membershipFunctions[i])->GetMeasurementVectorSize() << std::endl;<BR> <BR> finalparameters=(components[i])->GetFullParameters();<BR> paramIndex=0; <BR> for (unsigned int i1 = 0 ; i1 < measurementVectorSize ;
i1++)<BR> {<BR> <BR> estimatedmean[i1] = finalparameters[paramIndex] ;<BR> ++paramIndex ;<BR> }<BR> <BR> for (unsigned int i1 = 0 ; i1 < measurementVectorSize ; i1++ )<BR> {<BR> for (unsigned int j1 = 0 ; j1 < measurementVectorSize; j1++ )<BR> {<BR>
<BR> estimatedcovariance.GetVnlMatrix().put(i1, j1, finalparameters[paramIndex]) ;<BR> ++paramIndex ;<BR> }<BR> }<BR> <BR> // finalProportions[i]=(*estimator->GetProportions())[i];<BR> std::cout << "********************************" << std::endl;<BR> std::cout << "class[" << i << "]" << std::endl;<BR> std::cout << " Estimated Mean:" << std::endl;<BR> std::cout << " "
<< estimatedmean <BR> << std::endl;<BR> std::cout << " Estimated Covariance ";<BR> std::cout << " " << estimatedcovariance << std::endl;<BR> <BR> <BR> membershipFunctions[i]->SetMean(&estimatedmean );<BR> membershipFunctions[i]->SetCovariance( &estimatedcovariance);</DIV>
<DIV> std::cout << "***************************" << std::endl;<BR> std::cout << "class[" << i << "]" << std::endl;<BR> std::cout << " Gaussain Mean:" << std::endl;<BR> std::cout << " " <<*(membershipFunctions[i]->GetMean()+0)<< std::endl;<BR> // std::cout << " " <<*(membershipFunctions[i]->GetMean()+1)<< std::endl;<BR> std::cout << " Gaussian Covariance ";<BR> std::cout << " " << *(membershipFunctions[i]->GetCovariance()+0) <<
std::endl;<BR> // std::cout << " " << *(membershipFunctions[i]->GetCovariance()+1) << std::endl;<BR> // std::cout << " " << *(membershipFunctions[i]->GetCovariance()+3) << std::endl;</DIV>
<DIV> }</DIV>
<DIV> /////to lable the image </DIV>
<DIV> typedef unsigned char OutputPixelType;<BR> typedef itk::Image< OutputPixelType, 3 > OutputImageType;<BR> <BR> OutputImageType::Pointer outputPtr=reader1->GetOutput();<BR> typedef itk::ImageRegionIterator< OutputImageType > ImageIterator;</DIV>
<DIV> typedef VectorImageType::RegionType RegionType;<BR> RegionType region = outputPtr->GetBufferedRegion();<BR>// <BR>// <BR> ImageIterator pixel( outputPtr, region );</DIV>
<DIV> typedef itk::Vector< double, numberOfClasses > EvaluationVectorType;<BR> EvaluationVectorType eveluateresult;<BR> <BR> </DIV>
<DIV> //pixel.GoToBegin();</DIV>
<DIV>// itk::Array< double > finalProportions(numberOfClasses);<BR> for(unsigned char i=0;i<numberOfClasses;i++)<BR> {<BR> finalProportions[i]=(*estimator->GetProportions())[i];<BR> }<BR> <BR> double maxvalue;<BR> OutputPixelType classlabel;<BR> for ( pixel.GoToBegin(),constIterator.GoToBegin(); !pixel.IsAtEnd(); ++pixel,++constIterator )<BR> {<BR> mv=constIterator.Get();<BR> accumulate1=0;<BR>
maxvalue=0;<BR> for( unsigned int i=0;i<NumberOfComponents;i++)<BR> {<BR> accumulate1=accumulate1+mv[i];<BR> }</DIV>
<DIV> if(accumulate1>20)<BR> {<BR> </DIV>
<DIV> for(unsigned char i=0;i<numberOfClasses;i++)<BR> {<BR> eveluateresult[i]=(membershipFunctions[i]->Evaluate(mv))*finalProportions[i];<BR> }<BR> <BR> for(unsigned char i=0;i<numberOfClasses;i++)<BR> {<BR> if(eveluateresult[i]>=maxvalue)<BR>
{<BR> maxvalue=eveluateresult[i];<BR> classlabel=i+1;<BR> }<BR>// if(mv[0]==94&&mv[1]==77)<BR>// {<BR>// std::cout << "evaluation result" << eveluateresult << "]" << std::endl;<BR>// }</DIV>
<DIV> }<BR> pixel.Set(classlabel);<BR> }<BR> else<BR> {<BR> pixel.Set(0);<BR> }<BR> <BR> }<BR></DIV></DIV><BR></DIV></DIV></div><br>
</body></html>