##### Contents

back to iatexture

# Function histstat

## Synopse

The histstat function computes the statistics from the image histogram.

• h,stats = histstat(f, mask=[])
• Output
• h: histogram
• stats: array with histogram' statistics
• Input
• f: input image

## Description

The first block of histstat function computes the histogram from a specific image. So, the second block

extracts some relevant statistics about the image from the histogram that was found before. These statistics are

called Mazda's attributes and they are:

• Mean;
• Variance;
• Skewness;
• Kurtosis;
• Perc.01%;
• Perc.10%;
• Perc.50%;
• Perc.90%;
• Perc.99%;
• Entropy;
• Mode;
• Median;

## Function Code

``` 1 import ia636 as ia
2 import numpy as np
3
5
6
7     ### compute histogram ###
8
9     if mask !=[]:
10         h = ia.iahistogram(f[mask])
11     else:
12         h = ia.iahistogram(f)
13     hn = 1.0*h/h.sum()                                                                   # compute the normalized image histogram
14     stats = np.zeros(12)                                                                 # number of statistics
15
16
17     ### compute statistics ###
18
19     n = len(h) # number of gray values
20     stats[0]  = np.sum((np.arange(n)*hn))                                                # mean
21     stats[1]  = np.sum(np.power((np.arange(n)-stats[0]),2)*hn)                           # variance
22     stats[2]  = np.sum(np.power((np.arange(n)-stats[0]),3)*hn)/(np.power(stats[1],1.5))  # skewness
23     stats[3]  = np.sum(np.power((np.arange(n)-stats[0]),4)*hn)/(np.power(stats[1],2))-3  # kurtosis
24     stats[4:9]  = np.round(ia.iah2percentile(h,np.array([1,10,50,90,99])))               # 1,10,50,90,99% percentile
25
26
27     ## extra attributes ##
28
29     stats[9] =  -1*(hn*np.log10(hn+np.spacing(1))).sum()                                 # entropy
30     stats[10] = np.argmax(hn)                                                            # mode
31     stats[11] = np.where(np.cumsum(hn) >= 0.5)[0][0]                                     # median
32
33     if  stats[1] == 0:                                                                   # if std equal to zero, both skewness and kurtosis are zero
34         stats[2] = 0
35         stats[3] = 0
36
37     return h,stats
```

## Examples

### Numeric Example

``` 1 import ia636 as ia
2 import numpy as np
3 import scipy
4
5 from iatexture.histstat import histstat
6
7 f = np.array([1,1,1,0,1,2,2,2,1])
8 h,stats = histstat(f)
9 print 'histogram',h
10 print
11 print 'statistics =',stats
12 print scipy.stats.entropy(f)
```
```histogram [1 5 3]

statistics = [ 1.22222222  0.39506173 -0.20992233 -0.62109375  0.          1.          1.
2.          2.          0.40688542  1.          1.        ]
2.01981499249
```

### Numeric Example using a mask

``` 1 import numpy as np
2 import ia636 as ia
3
4 from iatexture.histstat import histstat
5
6 f = np.array([1,1,1,0,1,2,2,2,1])
7 mask = f>0
8 h,stats = histstat(f,mask)
9 print 'histogram',h
10 print
11 print 'statistics =',stats)
```
```------------------------------------------------------------
*** Exception while evaluating code:
File "<string>", line 11
print 'statistics =',stats)
^
SyntaxError: invalid syntax

------------------------------------------------------------

```

### Numeric Example 3D

``` 1 import numpy as np
2 import ia636 as ia
3
4 from iatexture.histstat import histstat
5
6 f = np.array([1,1,1,0,1,2,2,2,1,0,0,0,1,1,1,1,2,0,1,1,2,1,0,0]).reshape(2,3,4)
7 h,stats = histstat(f)
8 print 'histogram',h
9 print
10 print 'statistics =', stats
```
```histogram [ 7 12  5]

statistics = [ 0.91666667  0.49305556  0.11700664 -0.97242611  0.          0.          1.
2.          2.          0.44851494  1.          1.        ]
```

### Numeric Example 3D

``` 1 import numpy as np
2 import ia636 as ia
3
4 from iatexture.histstat import histstat
5
6 f = array([[[1,2,2,0,0,1],[0,0,1,2,2,1],[1,1,0,0,0,2],[1,1,1,2,2,2],[1,1,2,2,0,0]],
7     [[1,2,2,0,0,1],[0,0,1,2,2,1],[1,1,0,0,0,2],[1,1,1,2,2,2],[1,1,2,2,0,0]],
8     [[1,2,2,0,0,1],[0,0,1,2,2,1],[1,1,0,0,0,2],[1,1,1,2,2,2],[1,1,2,2,0,0]]], dtype=uint8)
9
10 h,stats = histstat(f)
11 print 'histogram',h
12 print
13 print 'statistics =', stats
```
```histogram [27 33 30]

statistics = [ 1.03333333  0.63222222 -0.05953097 -1.41606308  0.          0.          1.
2.          2.          0.47567118  1.          1.        ]
```

### Numeric Example 3D using a mask

``` 1 import numpy as np
2 import ia636 as ia
3
4 from iatexture.histstat import histstat
5
6 f = np.array([1,1,1,0,1,2,2,2,1,0,0,0,1,1,1,1,2,0,1,1,2,1,0,0]).reshape(2,3,4)
7 h,stats = histstat(f,mask = f>0)
8 print 'histogram',h
9 print
10 print 'statistics =', stats
```
```histogram [ 0 12  5]

statistics = [ 1.29411765  0.20761246  0.90369611 -1.18333333  1.          1.          1.
2.          2.          0.26309451  1.          1.        ]
```

### Image Example

``` 1 import numpy as np
2 import ia636 as ia
3
4 from iatexture.histstat import histstat
5
8 h,stats = histstat(f)
9
10 print 'histogram',h
11 print
12 print 'statistics =',stats
```
```histogram [   1  280 1229  952 1216  829  824  722  758  817  980  838 1172 1208 1353
1123  921  648  502  319  265  228  157  118  103   95   83   68   76   61
81   73   65   79   66   58   72   81   55   47   66   65   72   74   63
76   70   73   82   78   73   82   86   71   64   80   72   49   54   49
46   54   51   39   39   49   45   37   37   53   35   47   43   42   56
42   61   54   54   52   40   45   49   48   37   49   52   45   60   67
68   48   49   49   52   51   47   43   38   24   34   33   35   31   31
49   31   31   35   32   24   41   36   39   40   46   38   31   41   38
46   49   55   31   45   48   59   47   65   71   68   67   57   83   82
95   96   97   96   97   95  104  126  121  109  139  144  133  148  153
147  160  150  184  157  175  178  183  171  203  208  192  231  234  228
225  268  268  279  285  278  298  274  317  319  346  339  308  348  362
363  373  380  367  383  407  360  330  439  471  539  545  584  644  647
603  702  682  702  680  775  721  689  628  702  754  739  681  812  807
856  813  875  824  811  742  744  597  618  465  570  562  572  452  494
479  442  343  353  297  274  223  235  220  235  204  164  149  153  124
100   90   92   69   79   68   60   62   57   60   48  776]

statistics = [  1.37065933e+02   7.60425522e+03  -5.54976386e-01  -1.41520647e+00
2.00000000e+00   8.00000000e+00   1.80000000e+02   2.22000000e+02
2.51000000e+02   2.16964814e+00   1.40000000e+01   1.80000000e+02]
```

## Equation

The follow equations represents the statistics equations, where N denotes the number of intensity levels in the imagem ( usually between 0 and 255) and is the probability density function, i.e. the normalized histogram (h) of the image, and F is a random variable (any image) that represents the pixels intensities .

### Percentile

The k-th percentile is the value that corresponds to the cumulative frequency of , where N is the amostral size

### Mode

The value that appears most often.

### Median

It is the numerical value separating the higher half of a data sample or a probability distribution.

## Contributions

• Mariana Bento, August 2013: initial code
• Mariana Leite, September 2014: 3D examples

## Memory test

``` 1 import multiprocessing as mp
2 import numpy as np
3 import psutil
4 from ia636 import iacontour,iameshgrid,iaroi,iagshow
5 import glob, urllib, os
6
7 foldername = '/awmedia/www/media/p/LesionMRI/CAINdata/Original_Data/CAIN10320013/BASELINE - CAIN1/FLAIR_LongTR_601/*'
8 folder = glob.glob(foldername)
9 folder.sort()
10
11 img = zeros((len(folder),560,560),dtype=uint8)
12
13 for i in arange(len(folder)):
14      aux = dicom.read_file(find_attachment_file(folder[i]))
15      img[i] = aux.pixel_array.astype(float32)
16
17 # retalhar a imagem em retangulos
18 sizex = 30
19 sizey = 30
20 sizez = 30
21
22 z = arange(0,img.shape[0],sizez)
23 x = arange(0,img.shape[1],sizex)
24 y = arange(0,img.shape[2],sizey)
25 print 'initial memory use:  %9.3f MB' % (psutil.phymem_usage().used/1e06,)
26 for k in arange(z.shape[0]-1):
27     for i in arange(x.shape[0]-1):
28         for j in arange(y.shape[0]-1):
29             roi = np.zeros(img.shape, dtype = bool)
30             roi[z[k]:z[k+1],x[i]:x[i+1],y[j]:y[j+1]] =1
31             print 'memory used before:  %9.3f MB' % (psutil.phymem_usage().used/1e06,)
32             att_rec = histstat(img,mask = roi)
33             print 'memory used after:  %9.3f MB' % (psutil.phymem_usage().used/1e06,)
34
35
36
37 print 'final memory use:  %9.3f MB' % (psutil.phymem_usage().used/1e06,)
```
```------------------------------------------------------------
*** Exception while evaluating code:
File "<string>", line 14, in <module>
NameError: name 'dicom' is not defined

------------------------------------------------------------

```