histmed22Y {T4transport}

Wasserstein Median of Histograms by You et al. (2022)

Description

Given multiple histograms represented as "histogram" S3 objects, compute the Wasserstein median of order 2. We need one requirement that all histograms in an input list hists must have same breaks. See the example on how to construct a histogram on predefined breaks/bins.

Usage

histmed22Y(hists, weights = NULL, lambda = NULL, ...)

Arguments

hists

a length-N list of histograms ("histogram" object) of same breaks.

weights

a weight of each image; if NULL (default), uniform weight is set. Otherwise, it should be a length-N vector of nonnegative weights.

lambda

a regularization parameter; if NULL (default), a paper's suggestion would be taken, or it should be a nonnegative real number.

...

extra parameters including

abstol

stopping criterion for iterations (default: 1e-8).

init.vec

an initial weight vector (default: uniform weight).

maxiter

maximum number of iterations (default: 496).

nthread

number of threads for OpenMP run (default: 1).

print.progress

a logical to show current iteration (default: FALSE).

Value

a "histogram" object of the Wasserstein median histogram.

Examples

#----------------------------------------------------------------------
#                      Binned from Two Gaussians
#
# EXAMPLE : small example for CRAN for visualization purpose.
#----------------------------------------------------------------------
# GENERATE FROM TWO GAUSSIANS WITH DIFFERENT MEANS
set.seed(100)
x  = stats::rnorm(1000, mean=-4, sd=0.5)
y  = stats::rnorm(1000, mean=+4, sd=0.5)
bk = seq(from=-10, to=10, length.out=20)

# HISTOGRAMS WITH COMMON BREAKS
histxy = list()
histxy[[1]] = hist(x, breaks=bk, plot=FALSE)
histxy[[2]] = hist(y, breaks=bk, plot=FALSE)

# COMPUTE
hmean = histbary15B(histxy)
hmeds = histmed22Y(histxy)

# VISUALIZE
opar <- par(no.readonly=TRUE)
par(mfrow=c(1,3))
barplot(histxy[[1]]$density, col=rgb(0,0,1,1/4), 
        ylim=c(0, 1.05), main="Two Histograms")
barplot(histxy[[2]]$density, col=rgb(1,0,0,1/4), 
        ylim=c(0, 1.05), add=TRUE)
barplot(hmean$density, main="Barycenter",
        ylim=c(0, 1.05))
barplot(hmeds$density, main="Wasserstein Median",
        ylim=c(0, 1.05))
par(opar)