# [R] Cross-checking a custom function for separability indices

Nikos Alexandris nikos.alexandris at felis.uni-freiburg.de
Wed Apr 21 19:42:37 CEST 2010

```Hi list!

I have prepared a custom function (below) in order to calculate separability
indices (Divergence, Bhattacharyya, Jeffries-Matusita, Transformed divergene)
between two samples of (spectral land cover) classes.

I need help to cross-compare results to verify that it works as expected
(since I don't know of any other foss-tool that will give me quickly some
results).

Does anybody use another mean/tool/function (within or out of R) to calculate
these distances?

For example,

# two samples
sample.1 <- c (1362, 1411, 1457, 1735, 1621, 1621, 1791, 1863, 1863, 1838)
sample.2 <- c (1354, 1458, 1458, 1458, 1550, 1145, 1428, 1573, 1573, 1657)

# running the custom function (below)
separability.measures ( sample.1 , sample.2 )

Divergence:             1.5658
Bhattacharryya:         0.1683
Jeffries-Matusita:      0.3098
Transformed divergence: 0.3555

# --- Code -------------------------------------------------------------------
# Custom function for various Separability Measures
# by Nikos Alexandris, Freiburg, 8.04.2010

# ( based on Divergence and Jeffries-Matusita, requires input variables
"as.matrices" )
separability.measures <- function ( Vector.1 , Vector.2 ) {

# convert vectors to matrices in case they are not
Matrix.1 <- as.matrix (Vector.1)
Matrix.2 <- as.matrix (Vector.2)

# define means
mean.Matrix.1 <- mean ( Matrix.1 )
mean.Matrix.2 <- mean ( Matrix.2 )

# define difference of means
mean.difference <- mean.Matrix.1 - mean.Matrix.2

# define covariances for supplied matrices
cv.Matrix.1 <- cov ( Matrix.1 )
cv.Matrix.2 <- cov ( Matrix.2 )

# define the halfsum of cv's as "p"
p <- ( cv.Matrix.1 + cv.Matrix.2 ) / 2

# --%<------------------------------------------------------------------------
# calculate the Bhattacharryya index
bh.distance <- 0.125 *
t ( mean.difference ) *
p^ ( -1 ) *
mean.difference +
0.5 * log10 (
det ( p ) / sqrt (
det ( cv.Matrix.1 ) *
det ( cv.Matrix.2 )
)
)

# --%<------------------------------------------------------------------------
# calculate Jeffries-Matusita
jm.distance <- 2 * ( 1 - exp ( - bh.distance ) )

# --%<------------------------------------------------------------------------
# calculate the divergence
# trace (is the sum of the diagonal elements) of a square matrix
trace.of.matrix <- function ( SquareMatrix ) {
sum ( diag ( SquareMatrix ) ) }

# term 1
divergence.term.1 <- 1/2 *
trace.of.matrix (
( cv.Matrix.1 - cv.Matrix.2 ) *
( cv.Matrix.2^ (-1) - cv.Matrix.1^ (-1) )
)

# term 2
divergence.term.2 <- 1/2 *
trace.of.matrix (
( cv.Matrix.1^ (-1) + cv.Matrix.2^ (-1) ) *
( mean.Matrix.1 - mean.Matrix.2 ) *
t ( mean.Matrix.1 - mean.Matrix.2 )
)

# divergence
divergence <- divergence.term.1 + divergence.term.2

# --%<------------------------------------------------------------------------
# and the transformed divergence
transformed.divergence  <- 2 * ( 1 - exp ( - ( divergence / 8 ) ) )

# --%<------------------------------------------------------------------------
# print results --- look at "separability.measures.pp"

# get "column" names (hopefully they exist) --- to use for/with... ??? ---
name.of.Matrix.1 <- colnames ( Matrix.1 )
name.of.Matrix.2 <- colnames ( Matrix.2 )

# create new objects to be used with... ?
assign ( "divergence.vector"             , divergence             ,
envir = .GlobalEnv)

assign ( "jm.distance.vector"            , jm.distance            ,
envir = .GlobalEnv)

assign ( "bh.distance.vector"            , bh.distance            ,
envir = .GlobalEnv )

assign ( "transformed.divergence.vector" , transformed.divergence ,
envir = .GlobalEnv )

# print some message --- ?
# a "just do it" solution using cat()
cat ( paste ( "Divergence:             " ,
round ( divergence , digits = 4 ) , sep = "" ) , "\n")

cat ( paste ( "Bhattacharryya:         " ,
round ( bh.distance , digits = 4 ) , sep = "" ) , "\n")

cat ( paste ( "Jeffries-Matusita:      " ,
round ( jm.distance , digits = 4 ) , sep = "" ) , "\n")

cat ( paste ( "Transformed divergence: " ,

round ( transformed.divergence , digits = 4 ) , sep = "" ) , "\n" )
cat ("\n")

}
# --- Code -------------------------------------------------------------------

# --%<------------------------------------------------------------------------
# Sources

# Richards, John A., Jia, X., "Remote Sensing Digital Image Analysis: ...",
Springer
# Wikipedia

```