sid: A function for computing the spectral information divergence...

Description Usage Arguments Details Value Author(s) References See Also Examples

View source: R/sid.R

Description

\loadmathjax Experimental lifecycle

This function computes the spectral information divergence/dissimilarity between spectra based on the kullback-leibler divergence algorithm (see details).

Usage

1
2
3
4
5
6
7
8
sid(Xr, Xu = NULL,
    mode = "density",
    center = FALSE, scale = FALSE,
    kernel = "gaussian",
    n = if(mode == "density") round(0.5 * ncol(Xr)),
    bw = "nrd0",
    reg = 1e-04,
    ...)

Arguments

Xr

a matrix containing the spectral (reference) data.

Xu

an optional matrix containing the spectral data of a second set of observations.

mode

the method to be used for computing the spectral information divergence. Options are "density" (default) for computing the divergence values on the density distributions of the spectral observations, and "feature" for computing the divergence vales on the spectral variables. See details.

center

a logical indicating if the computations must be carried out on the centred X and Xu (if specified) matrices. If mode = "feature" centring is not carried out since this option does not accept negative values which are generated after centring the matrices. Default is FALSE. See details.

scale

a logical indicating if the computations must be carried out on the variance scaled X and Xu (if specified) matrices. Default is TRUE.

kernel

if mode = "density" a character string indicating the smoothing kernel to be used. It must be one of "gaussian" (default), "rectangular", "triangular", "epanechnikov", "biweight", "cosine" or "optcosine". See the density function of the stats package.

n

if mode = "density" a numerical value indicating the number of equally spaced points at which the density is to be estimated. See the density function of the stats package for further details. Default is round(0.5 * ncol(X)).

bw

if mode = "density" a numerical value indicating the smoothing kernel bandwidth to be used. Optionally the character string "nrd0" can be used, it computes the bandwidth using the bw.nrd0 function of the stats package (see bw.nrd0). See the density and the bw.nrd0 functions for more details. By default "nrd0" is used, in this case the bandwidth is computed as bw.nrd0(as.vector(X)), if Xu is specified the bandwidth is computed as bw.nrd0(as.vector(rbind(X, Xu))).

reg

a numerical value larger than 0 which indicates a regularization parameter. Values (probabilities) below this threshold are replaced by this value for numerical stability. Default is 1e-4.

...

additional arguments to be passed to the density function of the base package.

Details

This function computes the spectral information divergence (distance) between spectra. When mode = "density", the function first computes the probability distribution of each spectrum which result in a matrix of density distribution estimates. The density distributions of all the observations in the data sets are compared based on the kullback-leibler divergence algorithm. When mode = "feature", the kullback-leibler divergence between all the observations is computed directly on the spectral variables. The spectral information divergence (SID) algorithm (Chang, 2000) uses the Kullback-Leibler divergence (\mjeqnKLKL) or relative entropy (Kullback and Leibler, 1951) to account for the vis-NIR information provided by each spectrum. The SID between two spectra (\mjeqnx_ix_i and \mjeqnx_jx_j) is computed as follows:

\mjdeqn

sid(x_i,x_j) = KL(x_i \left |\right | x_j) + KL(x_j \left |\right | x_i)sid(x_i,x_j) = KL(x_i |\ | x_j) + KL(x_j |\ | x_i)

\mjdeqn

sid(x_i,x_j) = \sum_l=1^k p_l \ log(\fracp_lq_l) + \sum_l=1^k q_l \ log(\fracq_lp_l)sid(x_i,x_j) = \sum_l=1^k p_l \ log(p_l/q_l) + \sum_l=1^k q_l \ log(q_l/p_l)

where \mjeqnkk represents the number of variables or spectral features, \mjeqnpp and \mjeqnqq are the probability vectors of \mjeqnx_ix_i and \mjeqnx_ix_j respectively which are calculated as:

\mjdeqn

p = \fracx_i\sum_l=1^k x_i,lp = x_i/\sum_l=1^k x_i,l

\mjdeqn

q = \fracx_j\sum_l=1^k x_j,lq = x_j1/\sum_l=1^k x_j,l

From the above equations it can be seen that the original SID algorithm assumes that all the components in the data matrices are nonnegative. Therefore centering cannot be applied when mode = "feature". If a data matrix with negative values is provided and mode = "feature", the sid function automatically scales the matrix as follows:

\mjdeqn

X_s = \fracX-min(X)max(X)-min(X)X_s = X-min(X)/max(X)-min(X)

or

\mjdeqn

X_s = \fracX-min(X, Xu)max(X, Xu)-min(X, Xu)X_s = X-min(X, Xu)/max(X, Xu)-min(X, Xu)

\mjdeqn

Xu_s = \fracXu-min(X, Xu)max(X, Xu)-min(X, Xu)Xu_s = Xu-min(X, Xu)/max(X, Xu)-min(X, Xu)

if Xu is specified. The 0 values are replaced by a regularization parameter (reg argument) for numerical stability. The default of the sid function is to compute the SID based on the density distributions of the spectra (mode = "density"). For each spectrum in X the density distribution is computed using the density function of the stats package. The 0 values of the estimated density distributions of the spectra are replaced by a regularization parameter ("reg" argument) for numerical stability. Finally the divergence between the computed spectral histogramas is computed using the SID algorithm. Note that if mode = "density", the sid function will accept negative values and matrix centering will be possible.

Value

a list with the following components:

Author(s)

Leonardo Ramirez-Lopez

References

Chang, C.I. 2000. An information theoretic-based approach to spectral variability, similarity and discriminability for hyperspectral image analysis. IEEE Transactions on Information Theory 46, 1927-1932.

See Also

density

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
library(prospectr)

data(NIRsoil)

Xu <- NIRsoil$spc[!as.logical(NIRsoil$train), ]
Yu <- NIRsoil$CEC[!as.logical(NIRsoil$train)]
Yr <- NIRsoil$CEC[as.logical(NIRsoil$train)]
Xr <- NIRsoil$spc[as.logical(NIRsoil$train), ]

Xu <- Xu[!is.na(Yu), ]
Xr <- Xr[!is.na(Yr), ]

# Example 1
# Compute the SID distance between all the observations in Xr
xr_sid <- sid(Xr)
xr_sid

# Example 2
# Compute the SID distance between the observations in Xr and the observations
# in Xu
xr_xu_sid <- sid(Xr, Xu)
xr_xu_sid

resemble documentation built on Nov. 9, 2020, 5:08 p.m.