GLCM texture features

Calculates texture features from the input GLCMs

Description

The GLCMs are stored in a i x j x n matrix, where n is the number of GLCMs calculated usually due to the different orientation and displacements used in the algorithm. Usually the values i and j are equal to 'NumLevels' parameter of the GLCM computing function graycomatrix(). Note that matlab quantization values belong to the set {1,..., NumLevels} and not from {0,...,(NumLevels-1)} as provided in some references
http://www.mathworks.com/access/helpdesk/help/toolbox/images/graycomatrix.html

Although there is a function graycoprops() in Matlab Image Processing Toolbox that computes four parameters Contrast, Correlation, Energy, and Homogeneity. The paper by Haralick suggests a few more parameters that are also computed here. The code is not vectorized and hence is not an efficient implementation but it is easy to add new features based on the GLCM using this code. The code takes care of 3 dimensional glcms (multiple glcms in a single 3D array)

If you find that the values obtained are different from what you expect or if you think there is a different formula that needs to be used from the ones used in this code please let me know. A few questions which I have are listed in the link http://www.mathworks.com/matlabcentral/newsreader/view_thread/239608

I plan to submit a vectorized version of the code later and provide updates based on replies to the above link and this initial code.

% Features computed
% Autocorrelation: [2] (out.autoc)
% Contrast: matlab/[1,2] (out.contr)
% Correlation: matlab (out.corrm)
% Correlation: [1,2] (out.corrp)
% Cluster Prominence: [2] (out.cprom)
% Cluster Shade: [2] (out.cshad)
% Dissimilarity: [2] (out.dissi)
% Energy: matlab / [1,2] (out.energ)
% Entropy: [2] (out.entro)
% Homogeneity: matlab (out.homom)
% Homogeneity: [2] (out.homop)
% Maximum probability: [2] (out.maxpr)
% Sum of sqaures: Variance [1] (out.sosvh)
% Sum average [1] (out.savgh)
% Sum variance [1] (out.svarh)
% Sum entropy [1] (out.senth)
% Difference variance [1] (out.dvarh)
% Difference entropy [1] (out.denth)
% Information measure of correlation1 [1] (out.inf1h)
% Informaiton measure of correlation2 [1] (out.inf2h)
% Inverse difference (INV) is homom [3] (out.homom)
% Inverse difference normalized (INN) [3] (out.indnc)
% Inverse difference moment normalized [3](out.idmnc)

Haralick uses 'Symmetric' = true in computing the glcm. There is no Symmetric flag in the Matlab version I use hence I add the diagonally opposite pairs to obtain the Haralick glcm. Here it is assumed that the diagonally opposite orientations are paired one after the other in the matrix. If the above assumption is true with respect to the input glcm then setting the flag 'pairs' to 1 will compute the final glcms that would result by setting 'Symmetric' to true. If your glcm is computed using the
Matlab version with 'Symmetric' flag you can set the flag 'pairs' to 0

% References:
1. R. M. Haralick, K. Shanmugam, and I. Dinstein, Textural Features of Image Classification, IEEE Transactions on Systems, Man and Cybernetics, vol. SMC-3, no. 6, Nov. 1973
2. L. Soh and C. Tsatsoulis, Texture Analysis of SAR Sea Ice Imagery Using Gray Level Co-Occurrence Matrices, IEEE Transactions on Geoscience and Remote Sensing, vol. 37, no. 2, March 1999.
3. D A. Clausi, An analysis of co-occurrence texture statistics as a
function of grey level quantization, Can. J. Remote Sensing, vol. 28, no.1, pp. 45-62, 2002
4. http://murphylab.web.cmu.edu/publications/boland/boland_node26.html

% Example:
% Usage is similar to graycoprops() but needs extra parameter 'pairs' apart from the GLCM as input

>I = imread('circuit.tif');
>GLCM2 = graycomatrix(I,'Offset',[2 0;0 2]);
>stats = GLCM_features1(GLCM2,0)

The output is a structure containing all the features calculated from the different GLCMs

• Required Products: Image Processing Toolbox
• MATLAB release: MATLAB 7.0.1 (R14SP1)

[matlab res]

I am new to this environment.I dont know what parameters should be given & how to give.
I get the following error.Give me the step wise details of how to arrive at the output.
 
??? Error using ==> GLCM_Features1
Too many or too few input arguments. Enter GLCM and pairs.

[matlab res]

I get this error . What should be done ?
GLCM2 = graycomatrix(I,'Offset',[2 0;0 2]);
??? Undefined command/function 'graycomatrix'.

[Avinash Uppuluri]

>I = imread('circuit.tif');
>GLCM2 = graycomatrix(I,'Offset',[2 0;0 2]);
>stats = GLCM_features1(GLCM2,0)

The above code should work fine if you have Matlab Image Processing Toolbox with the function graycomatrix() to compute the Gray Level Cooccurrence Matrix (GLCM). [ Try > help graycomatrix and see what if the function is present in your version]

If you are using any other method to compute the GLCM (e.g., a self written function) you can use GLCM_Features1() with the first input as the GLCM (i x j x n matrix) and a parameter pairs which I use to indicate if a symmetric GLCM needs to be computed from opposite direction pairs. Go through the comments in the code to understand better. Please let me know if you still have a problem.

For the error

"??? Error using ==> GLCM_Features1
Too many or too few input arguments. Enter GLCM and pairs."

Let me know what kind of input parameters you gave the function and I will try to help.

Thanks

[matlab res]

>> help graycomatrix()

graycomatrix().m not found.

I get this error . What should i do?

[Avinash Uppuluri]

It looks like the version of Image Processing Toolbox of matlab that you are working with does not have graycomatrix().

You can write a function to compute the GLCM for a image using the details provided in http://www.mathworks.com/access/helpdesk/help/toolbox/images/graycomatrix.html and the references listed in the description.

[Bilwaj Gaonkar]

 

[Jonathan Pedron]

Good job ;)
It works perfectly well.
I am not an expert on Haralick coefficient but it can be done in different direction, am I right?
If it is right would be interesting to be able to set this parameter...

[Avinash Uppuluri]

Jonathan,
You are right! Haralick coefficients can be computed on GLCMs calculated in different directions. The matlab function graycomatrix() gives option to change this direction as required. http://www.mathworks.com/access/helpdesk/help/toolbox/images/graycomatrix.html

A good read for this topic: Image Processing: Dealing with texture by Maria Petrou and Pedro Garcia Sevilla. Let me know if you need any specifics about changing the orientation. Thanks.

[Avinash Uppuluri]

Hi Preeti,

Try something like this:
% Form GLCMs for each of the 354 images using the required offsets (and other parameters as required) using the MATLAB function graycomatrix()

e.g., GLCM2 = graycomatrix(Dicom_images(:,:,1),'Offset',[0 2]);

% Then compute the features for each GLCM/image. You can also have multiple GLCMs per image depending on the parameters you input to graycomatrix().

stats = GLCM_features1(GLCM2,0)

Hope this helps,
Avinash

Hi
 I have read your code GLCM texture features. I have a mat file of size 512*512*354. in which 354 DICOM images are saved.
I want to use your code to extract texture features of all images in .mat file. please help.
Preeti

[Preeti]

Thanks for reply. Actually i am not able to load the data from mat file. can you please help me in using the data in .mat file? what is the code to import data from .mat file one by one , then only i can form GLCM for every dicom image. Please also tell that how to save these features in a database

Looking forward to your reply.

[gasmi karim]

Please I want known is that in calculating the Haralick feature of such energy to the entire image or for each pixel because in your code I think you calculate for the whole image a single value for each feature so how it can be used for classification is for this I see the need to calculate for each pixel the feature

[gasmi karim]

I work on the MRI image but I do not know how I can use to calculate your classification is for the I told you that the feature is calculated for each pixel I await your response

[Avinash Uppuluri]

http://www.mathworks.com/matlabcentral/fileexchange/22354-glcmfeatures4-m-vectorized-version-of-glcmfeatures1-m-with-code-changes

The above link has a faster version of this code (named GLCM_features4.m). They should give the same outputs. Do let me know if you find any difference in the outputs.

Thanks,
Avinash

[Avinash Uppuluri]

There have been a few questions on how to compute features for each pixel in a image. The link below provides further explanation into such details. See "Creating a texture image"

http://www.fp.ucalgary.ca/mhallbey/texture_calculations.htm

Main points being:
+ Texture for a pixel needs the definition of a "window" around the pixel
+ Once such a window is defined; you can calculate the GLCM for each pixel (taking the respective window into consideration);
+ These GLCMs for each pixel should later be fed into the GLCM_features function to extract the feature vector for each pixel (and related window)
+ If needed, the values can later be put together to form a texture feature image

Hope this helps.

[kapil veera]

 

[Kavitha Ravi]

Hi
I tried the above code and got an error stating "Function undefined for input of type double".
What type should i convert my GLCM matrix for this function to work ?
Thanks for the help.

Kavitha

[Avinash Uppuluri]

Hi Kavitha,

Can you please provide more information on how you generate your GLCM matrix (if possible also provide the image used) so that I can reproduce the error on my end.

Thanks,
Avinash

Hi
I tried the above code and got an error stating "Function undefined for input of type double".
What type should i convert my GLCM matrix for this function to work ?
Thanks for the help.

Kavitha

[Avinash Uppuluri]

Kapil,

I would appreciate if you can leave a short comment on how you think this work can be improved so that I can make it more useful for yourself and others.

Thanks,
Avinash

[bala k]

avinash uppuluri,

i'm doing image processing using GLCM function energy, contrast, entropy, inverse matrix, can you give how i generate for these function code.
my process is
query_image = imread('145.jpg');
grayimage = rgb2gray(query_image);
i need to calculate gray level co-occurrence matrix,
energy E,
entropy I,
contrast S,
inverse difference H,
F_query = [E I S H];

image database also used above process.

F_data = [E1 I1 S1 H1];

euclidean distance

d = norm(F_query - F_data);

i will get the best match images.
 
thanks,
bala

[Reem Md]

 

[Reem Md]

Hi Avinash
Good job ;)
Can you please help me, I'm trying to move a window, say 9X9 on a matrix, and create for every region (9X9) the glcm, and then move this 9X9 window until I reach the end of the matrix. How can I do that?

[Avinash Uppuluri]

Hi bala,

Please try to explain further. From what I understand it looks like you want to know how to generate the E I S H for each image in your database and use the same parameters of the input image to find the closest ( euclidean ) image to it in the database.

You can generate the GLCM for each image using the in built matlab function graycomatrix() and then pass the generated GLCM to the function you downloaded here GLCM_featuresX()

Hope that helps,
Avinash

[adzlan]

hi Avinash

im new in matlab. can you tell me how can i call only this 13 features from your code
• Angular Second Moment
• Contrast
• Correlation
• Variance
• Inverse Second Differential Moment
• Sum Average
• Sum Variance
• Sum Entropy
• Entropy
• Difference Variance
• Difference Entropy
• Measure of Correlation 1
• Measure of Correlation 2

thanks,
adzlan

[Hervé Chubaka Bagalwa]

Thank you. I have not checked if the output values are correct, but I'm just going to go and plug it in my classifier.

I will have more comments for you when i'm done with my project.

Thanks again

[leila]

 

[Will]

function [out] = cad_glcm_features(glcm)

% VECTORIZED CODE: FASTER

size_glcm_1 = size(glcm,1);
size_glcm_2 = size(glcm,2);
size_glcm_3 = size(glcm,3);

% checked
out.autoc = zeros(1,size_glcm_3); % Autocorrelation: [2]
out.contr = zeros(1,size_glcm_3); % Contrast: matlab/[1,2]
out.corrm = zeros(1,size_glcm_3); % Correlation: matlab
out.corrp = zeros(1,size_glcm_3); % Correlation: [1,2]
out.cprom = zeros(1,size_glcm_3); % Cluster Prominence: [2]
out.cshad = zeros(1,size_glcm_3); % Cluster Shade: [2]
out.dissi = zeros(1,size_glcm_3); % Dissimilarity: [2]
out.energ = zeros(1,size_glcm_3); % Energy: matlab / [1,2]
out.entro = zeros(1,size_glcm_3); % Entropy: [2]
out.homom = zeros(1,size_glcm_3); % Homogeneity: matlab
out.homop = zeros(1,size_glcm_3); % Homogeneity: [2]
out.maxpr = zeros(1,size_glcm_3); % Maximum probability: [2]
out.sosvh = zeros(1,size_glcm_3); % Sum of sqaures: Variance [1]
out.savgh = zeros(1,size_glcm_3); % Sum average [1]
out.svarh = zeros(1,size_glcm_3); % Sum variance [1]
out.senth = zeros(1,size_glcm_3); % Sum entropy [1]
out.dvarh = zeros(1,size_glcm_3); % Difference variance [4]
out.denth = zeros(1,size_glcm_3); % Difference entropy [1]
out.inf1h = zeros(1,size_glcm_3); % Information measure of correlation1 [1]
out.inf2h = zeros(1,size_glcm_3); % Informaiton measure of correlation2 [1]
out.indnc = zeros(1,size_glcm_3); % Inverse difference normalized (INN) [3]
out.idmnc = zeros(1,size_glcm_3); % Inverse difference moment normalized [3]

% Indices
[i,j] = meshgrid(1:size_glcm_1,1:size_glcm_2);
idx1 = (i+j)-1;
idx2 = abs(i-j)+1;
ii = (1:(2*size_glcm_1-1))';
jj = (0:size_glcm_1-1)';

for k = 1:size_glcm_3 % number glcms
    % Normalize GLCM
    glcm_sum = sum(sum(glcm(:,:,k)));
    Pij = glcm(:,:,k)./glcm_sum; % Normalize each glcm
    glcm_mean = mean(Pij(:)); % compute mean after norm
    %
    p_x = squeeze(sum(Pij,2));
    p_y = squeeze(sum(Pij,1))';
    %
    u_x = sum(sum(i.*Pij));
    u_y = sum(sum(j.*Pij));
    %
    p_xplusy = zeros((2*size_glcm_1 - 1),1); %[1]
    p_xminusy = zeros((size_glcm_1),1); %[1]
    for aux = 1:max(idx1(:))
       p_xplusy(aux) = sum(Pij(idx1==aux));
    end
    for aux = 1:max(idx2(:))
       p_xminusy(aux) = sum(Pij(idx2==aux));
    end
    
    % Contrast
    out.contr(k) = sum(sum((abs(i-j).^2).*Pij));
    % Dissimilarity
    out.dissi(k) = sum(sum(abs(i-j).*Pij));
    % Energy
    out.energ(k) = sum(sum(Pij.^2));
    % Entropy
    out.entro(k) = -sum(sum(Pij.*log(Pij+eps)));
    % Homogeneity Matlab
    out.homom(k) = sum(sum(Pij./(1+abs(i-j))));
    % Homogeneity Paper
    out.homop(k) = sum(sum(Pij./(1+abs(i-j).^2)));
    % Sum of squares: Variance
    out.sosvh(k) = sum(sum(Pij.*((j-glcm_mean).^2)));
    % Inverse difference normalized
    out.indnc(k) = sum(sum(Pij./(1+(abs(i-j)./size_glcm_1))));
    % Inverse difference moment normalized
    out.idmnc(k) = sum(sum(Pij./(1+((i-j)./size_glcm_1).^2)));
    % Maximum probability
    out.maxpr(k) = max(Pij(:));
    % Sum average
    out.savgh(k) = sum((ii+1).*p_xplusy);
    % Sum entropy
    out.senth(k) = -sum(p_xplusy.*log(p_xplusy+eps));
    % Sum variance
    out.svarh(k) = sum((((ii+1) - out.senth(k)).^2).*p_xplusy);
    % Difference entropy
    out.denth(k) = -sum(p_xminusy.*log(p_xminusy+eps));
    % Difference variance
    out.dvarh(k) = sum((jj.^2).*p_xminusy);
    % Computes correlation
    hxy1 = -sum(sum(Pij.*log(p_x*p_y' + eps)));
    hxy2 = -sum(sum((p_x*p_y').*log(p_x*p_y' + eps)));
    hx = -sum(p_x.*log(p_x+eps));
    hy = -sum(p_y.*log(p_y+eps));
    hxy = out.entro(k);
    % Information measure of correlation 1
    out.inf1h(k) = (hxy-hxy1)/(max([hx,hy]));
    % Information measure of correlation 2
    out.inf2h(k) = (1-exp(-2*(hxy2-hxy)))^0.5;
    % Cluster Prominence
    out.cprom(k) = sum(sum(Pij.*((i+j-u_x-u_y).^4)));
    % Cluster Shade
    out.cshad(k) = sum(sum(Pij.*((i+j-u_x-u_y).^3)));
    %
    s_x = sum(sum(Pij.*((i-u_x).^2)))^0.5;
    s_y = sum(sum(Pij.*((j-u_y).^2)))^0.5;
    corp = sum(sum(Pij.*(i.*j)));
    corm = sum(sum(Pij.*(i-u_x).*(j-u_y)));
    % Autocorrelation
    out.autoc(k) = corp;
    % Correlation paper
    out.corrp(k) = (corp-u_x*u_y)/(s_x*s_y);
    % Correlation Matlab
    out.corrm(k) = corm/(s_x*s_y);
end

[Avinash Uppuluri]

Will, Thanks for the update.

[bhargav m]

Hi all

Can you explain the physical interpretation of the Inverse difference moment , Information measures of correlation , cluster prominence and shade in a image.
I was trying to get the physical meaning of them , but could not find properly.

Thanks

[giet]