R/HellingerDist.R In distrEx: Extensions of Package 'distr'

```###############################################################################
## Method: HellingerDist
## Hellinger distance of two distributions
###############################################################################
setMethod("HellingerDist", signature(e1 = "AbscontDistribution",
e2 = "AbscontDistribution"),
function(e1, e2, rel.tol=.Machine\$double.eps^0.3,
TruncQuantile = getdistrOption("TruncQuantile"),
IQR.fac = 15, ..., diagnostic = FALSE){
## find sensible lower and upper bounds for integration
# (a) quantile based
low <- min(getLow(e1, eps = TruncQuantile), getLow(e2, eps = TruncQuantile))
up  <- max(getUp(e1, eps = TruncQuantile), getUp(e2, eps = TruncQuantile))
# (b) scale based
s0 <- min(IQR(e1),IQR(e2))*IQR.fac
low0 <- min(median(e1),median(e2))-s0
up0 <- max(median(e1),median(e2))+s0
# (a) & (b)
lower <- max(low,low0); upper <- min(up,up0)

o.warn <- getOption("warn"); options(warn = -1)
on.exit(options(warn=o.warn))
integrand <- function(x){ 0.5*(sqrt(d(e1)(x))-sqrt(d(e2)(x)))^2 }
res <- distrExIntegrate(integrand, lower = lower, upper = upper,
rel.tol = rel.tol, diagnostic = diagnostic)
names(res) <- "Hellinger distance"
if(diagnostic){
diagn <- attr(res,"diagnostic")
diagn[["call"]] <- match.call()
}
return(res)
})
setMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "DiscreteDistribution"),
function(e1, e2, ...){
o.warn <- getOption("warn"); options(warn = -1)
on.exit(options(warn=o.warn))
supp <- union(support(e1), support(e2))

res <- 0.5*sum((sqrt(d(e1)(supp))-sqrt(d(e2)(supp)))^2)
names(res) <- "Hellinger distance"

return(sqrt(res)) # ^.5 added P.R. 19-12-06
})

## new PR 08-09-16
setMethod("HellingerDist", signature(e1 = "DiscreteMVDistribution",
e2 = "DiscreteMVDistribution"),
function(e1, e2, ...){
o.warn <- getOption("warn"); options(warn = -1)
on.exit(options(warn=o.warn))
## replace univariate line  supp <- union(support(e1), support(e2))   by

supp <- unique(rbind(support(e1), support(e2)))

res <- 0.5*sum((sqrt(d(e1)(supp))-sqrt(d(e2)(supp)))^2)
names(res) <- "Hellinger distance"

return(sqrt(res))
})

setMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "AbscontDistribution"),
function(e1, e2, ...){
res <- 1
names(res) <- "Hellinger distance"

return(res)
})
setMethod("HellingerDist", signature(e1 = "AbscontDistribution",
e2 = "DiscreteDistribution"),
function(e1, e2, ...){
res <- 1
names(res) <- "Hellinger distance"

return(res)
})
## Hellinger distance
setMethod("HellingerDist", signature(e1 = "numeric",
e2 = "DiscreteDistribution"),
function(e1, e2, ... ){
d1 <- table(e1)/length(e1)
d2 <- d(e2)(sort(unique(e1)))
e21 <- setdiff(support(e2), unique(e1))
d21 <- d(e2)(e21)
res <- sqrt(1/2)*sqrt(sum((sqrt(d1)-sqrt(d2))^2) + sum(d21))
names(res) <- "Hellinger distance"

return(res)
})
setMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "numeric"),
function(e1, e2, ...){
return(HellingerDist(e2, e1))
})

## to avoid trivial distances (distance = 1)
## abs.cont. distributions may be discretized
## resp. empirical distributions may be smoothed
## (by convolution with a normal distribution)
setMethod("HellingerDist", signature(e1 = "numeric",
e2 = "AbscontDistribution"),
function(e1, e2, asis.smooth.discretize = "discretize", n.discr =
getdistrExOption("nDiscretize"), low.discr = getLow(e2),
up.discr = getUp(e2), h.smooth = getdistrExOption("hSmooth"),
rel.tol=.Machine\$double.eps^0.3,
TruncQuantile = getdistrOption("TruncQuantile"),
IQR.fac = 15, ..., diagnostic = FALSE){
res <- .asis.smooth.discretize.distance(e1, e2, asis.smooth.discretize, n.discr,
low.discr, up.discr, h.smooth, HellingerDist,
rel.tol = rel.tol, TruncQuantile = TruncQuantile,
IQR.fac = IQR.fac, ..., diagnostic = diagnostic)
if(diagnostic){
diagn <- attr(res,"diagnostic")
diagn[["call"]] <- match.call()
class(diagn)<- "DiagnosticClass"
attr(res,"diagnostic") <- diagn
}
return(res)
})
setMethod("HellingerDist", signature(e1 = "AbscontDistribution",
e2 = "numeric"),
function(e1, e2, asis.smooth.discretize = "discretize", n.discr =
getdistrExOption("nDiscretize"), low.discr = getLow(e1),
up.discr = getUp(e1), h.smooth = getdistrExOption("hSmooth"),
rel.tol=.Machine\$double.eps^0.3,
TruncQuantile = getdistrOption("TruncQuantile"),
IQR.fac = 15, ..., diagnostic = FALSE){
res <- HellingerDist(e2, e1, asis.smooth.discretize = asis.smooth.discretize,
low.discr = low.discr, up.discr = up.discr, h.smooth = h.smooth,
rel.tol = rel.tol, TruncQuantile = TruncQuantile,
IQR.fac = IQR.fac, ..., diagnostic = diagnostic)
if(diagnostic){
diagn <- attr(res,"diagnostic")
diagn[["call"]] <- match.call()
class(diagn)<- "DiagnosticClass"
attr(res,"diagnostic") <- diagn
}
return(res)
})

#### new from version 2.0 on: Distance for Mixing Distributions
setMethod("HellingerDist", signature(e1 = "AcDcLcDistribution",
e2 = "AcDcLcDistribution"),
function(e1, e2, rel.tol=.Machine\$double.eps^0.3,
TruncQuantile = getdistrOption("TruncQuantile"),
IQR.fac = 15, ..., diagnostic = FALSE){
if( is(e1,"AbscontDistribution"))
e1 <- as(as(e1,"AbscontDistribution"), "UnivarLebDecDistribution")
if( is(e2,"AbscontDistribution"))
e2 <- as(as(e2,"AbscontDistribution"), "UnivarLebDecDistribution")
if(is(e1,"DiscreteDistribution"))
e1 <- as(as(e1,"DiscreteDistribution"), "UnivarLebDecDistribution")
if(is(e2,"DiscreteDistribution"))
e2 <- as(as(e2,"DiscreteDistribution"), "UnivarLebDecDistribution")
ac1 <- acPart(e1); ac2 <- acPart(e2)
ac1d <- ac1@d; ac2d <- ac2@d
ac1@d <- function(x) ac1d(x)*acWeight(e1)
ac2@d <- function(x) ac2d(x)*acWeight(e2)
dc1 <- discretePart(e1); dc2 <- discretePart(e2)
dc1d <- dc1@d; dc2d <- dc2@d
dc1@d <- function(x) dc1d(x)*discreteWeight(e1)
dc2@d <- function(x) dc2d(x)*discreteWeight(e2)
da2 <- HellingerDist(ac1,ac2, rel.tol = rel.tol,
TruncQuantile = TruncQuantile, IQR.fac = IQR.fac, ..., diagnostic = diagnostic)
dd2 <- HellingerDist(dc1,dc2)
res <- (da2^2+dd2^2)^.5
names(res) <- "Hellinger distance"
if(diagnostic){
diagn <- attr(da2,"diagnostic")
diagn[["call"]] <- match.call()
class(diagn)<- "DiagnosticClass"
attr(res,"diagnostic") <- diagn
}
return(res)
})

setMethod("HellingerDist", signature(e1 = "LatticeDistribution",
e2 = "LatticeDistribution"),
getMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "DiscreteDistribution")))

setMethod("HellingerDist", signature(e1 = "LatticeDistribution",
e2 = "DiscreteDistribution"),
getMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "DiscreteDistribution")))

setMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "LatticeDistribution"),
getMethod("HellingerDist", signature(e1 = "DiscreteDistribution",
e2 = "DiscreteDistribution")))
```

Try the distrEx package in your browser

Any scripts or data that you put into this service are public.

distrEx documentation built on May 2, 2019, 5:16 p.m.