ITK  5.3.0 Insight Toolkit
Examples/Statistics/ImageEntropy1.cxx
/*=========================================================================
*
*
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
*
* Unless required by applicable law or agreed to in writing, software
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
*
*=========================================================================*/
// Software Guide : BeginLatex
//
//
// This example shows how to compute the entropy of an image. More formally
// this should be said : The reduction in uncertainty gained when we measure
// the intensity of \emph{one} randomly selected pixel in this image, given
// that we already know the statistical distribution of the image intensity
// values.
//
// In practice it is almost never possible to know the real statistical
// distribution of intensities and we are forced to estimate it from the
// evaluation of the histogram from one or several images of similar nature.
// We can use the counts in histogram bins in order to compute frequencies and
// then consider those frequencies to be estimations of the probablility of a
// new value to belong to the intensity range of that bin.
//
// \index{Entropy!Images}
// \index{Image!Entropy}
// \index{Image!Amount of information}
// \index{Amount of information!Image}
//
// Software Guide : EndLatex
// Software Guide : BeginLatex
//
// Since the first stage in estimating the entropy of an image is to compute
// its histogram, we must start by including the headers of the classes that
// will perform such a computation. In this case, we are going to use a scalar
// image as input, therefore we need the
// \subdoxygen{Statistics}{ScalarImageToHistogramGenerator} class, as well as
// the image class.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
#include "itkImage.h"
// Software Guide : EndCodeSnippet
int
main(int argc, char * argv[])
{
if (argc < 3)
{
std::cerr << "Missing command line arguments" << std::endl;
std::cerr << "Usage : ImageEntropy1 inputImageFileName ";
std::cerr << "numberOfHistogramBins" << std::endl;
return EXIT_FAILURE;
}
// Software Guide : BeginLatex
//
// The pixel type and dimension of the image are explicitly declared and
// then used for instantiating the image type.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using PixelType = unsigned char;
constexpr unsigned int Dimension = 3;
// Software Guide : EndCodeSnippet
try
{
}
catch (const itk::ExceptionObject & excp)
{
std::cerr << "Problem encoutered while reading image file : " << argv[1]
<< std::endl;
std::cerr << excp << std::endl;
return EXIT_FAILURE;
}
// Software Guide : BeginLatex
//
// The image type is used as template parameter for instantiating the
// histogram generator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using HistogramGeneratorType =
HistogramGeneratorType::Pointer histogramGenerator =
HistogramGeneratorType::New();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The parameters of the desired histogram are defined, including the
// number of bins and the marginal scale. For convenience in this example,
// we read the number of bins from the command line arguments. In this way
// we can easily experiment with different values for the number of bins and
// see how that choice affects the computation of the entropy.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
const unsigned int numberOfHistogramBins = std::stoi(argv[2]);
histogramGenerator->SetNumberOfBins(numberOfHistogramBins);
histogramGenerator->SetMarginalScale(10.0);
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We can then connect as input the output image from a reader and trigger
// the histogram computation by invoking the \code{Compute()} method in the
// generator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
histogramGenerator->Compute();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The resulting histogram can be recovered from the generator by using the
// \code{GetOutput()} method. A histogram class can be declared using the
// \code{HistogramType} trait from the generator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
using HistogramType = HistogramGeneratorType::HistogramType;
const HistogramType * histogram = histogramGenerator->GetOutput();
// Software Guide : EndCodeSnippet
const unsigned int histogramSize = histogram->Size();
std::cout << "Histogram size " << histogramSize << std::endl;
for (unsigned int bin = 0; bin < histogramSize; ++bin)
{
std::cout << "bin = " << bin << " frequency = ";
std::cout << histogram->GetFrequency(bin, 0) << std::endl;
}
// Software Guide : BeginLatex
//
// We proceed now to compute the \emph{estimation} of entropy given the
// histogram. The first conceptual jump to be done here is to assume that
// the histogram, which is the simple count of frequency of occurrence for
// the gray scale values of the image pixels, can be normalized in order to
// estimate the probability density function \textbf{PDF} of the actual
// statistical distribution of pixel values.
//
// First we declare an iterator that will visit all the bins in the
// histogram. Then we obtain the total number of counts using the
// \code{GetTotalFrequency()} method, and we initialize the entropy
// variable to zero.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
HistogramType::ConstIterator itr = histogram->Begin();
HistogramType::ConstIterator end = histogram->End();
double Sum = histogram->GetTotalFrequency();
double Entropy = 0.0;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We start now visiting every bin and estimating the probability of a pixel
// to have a value in the range of that bin. The base 2 logarithm of that
// probability is computed, and then weighted by the probability in order to
// compute the expected amount of information for any given pixel. Note that
// a minimum value is imposed for the probability in order to avoid
// computing logarithms of zeros.
//
// Note that the $\log{(2)}$ factor is used to convert the natural
// logarithm in to a logarithm of base 2, and makes it possible to report
// the entropy in its natural unit: the bit.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
while (itr != end)
{
const double probability = itr.GetFrequency() / Sum;
if (probability > 0.99 / Sum)
{
Entropy += -probability * std::log(probability) / std::log(2.0);
}
++itr;
}
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The result of this sum is considered to be our estimation of the image
// entropy. Note that the Entropy value will change depending on the number
// of histogram bins that we use for computing the histogram. This is
// particularly important when dealing with images whose pixel values have
// dynamic ranges so large that our number of bins will always underestimate
// the variability of the data.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
std::cout << "Image entropy = " << Entropy << " bits " << std::endl;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// As an illustration, the application of this program to the image
//
// \begin{itemize}
// \item \code{Examples/Data/BrainProtonDensitySlice.png}
// \end{itemize}
//
// results in the following values of entropy for different values of number
// of histogram bins.
//
// \begin{center}
// \begin{tabular}{|l|r|r|r|r|r|}
// \hline
// Number of Histogram Bins & 16 & 32 & 64 & 128 & 255 \cr
// \hline
// Estimated Entropy (bits) & 3.02 & 3.98 & 4.92 & 5.89 & 6.88 \cr
// \hline
// \end{tabular}
// \end{center}
//
//
// This table highlights the importance of carefully considering the
// characteristics of the histograms used for estimating Information Theory
// measures such as the entropy.
//
// Software Guide : EndLatex
return EXIT_SUCCESS;
}