# hdrlevels: Highest Density Region (HDR) Levels In mclust: Gaussian Mixture Modelling for Model-Based Clustering, Classification, and Density Estimation

## Description

Compute the levels of Highest Density Regions (HDRs) for any density and probability levels.

## Usage

 1 hdrlevels(density, prob) 

## Arguments

 density A vector of density values computed on a set of (observed) evaluation points. prob A vector of probability levels in the range [0,1].

## Details

From Hyndman (1996), let f(x) be the density function of a random variable X. Then the 100(1-α)\% HDR is the subset R(f_α) of the sample space of X such that

R(f_α) = {x : f(x) ≥ f_α }

where f_α is the largest constant such that Pr( X \in R(f_α)) ≥ 1-α

## Value

The function returns a vector of density values corresponding to HDRs at given probability levels.

L. Scrucca

## References

Rob J. Hyndman (1996) Computing and Graphing Highest Density Regions. The American Statistician, 50(2):120-126.

plot.densityMclust

## Examples

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 # Example: univariate Gaussian x <- rnorm(1000) f <- dnorm(x) a <- c(0.5, 0.25, 0.1) (f_a <- hdrlevels(f, prob = 1-a)) plot(x, f) abline(h = f_a, lty = 2) text(max(x), f_a, labels = paste0("f_", a), pos = 3) mean(f > f_a[1]) range(x[which(f > f_a[1])]) qnorm(1-a[1]/2) mean(f > f_a[2]) range(x[which(f > f_a[2])]) qnorm(1-a[2]/2) mean(f > f_a[3]) range(x[which(f > f_a[3])]) qnorm(1-a[3]/2) # Example 2: univariate Gaussian mixture set.seed(1) cl <- sample(1:2, size = 1000, prob = c(0.7, 0.3), replace = TRUE) x <- ifelse(cl == 1, rnorm(1000, mean = 0, sd = 1), rnorm(1000, mean = 4, sd = 1)) f <- 0.7*dnorm(x, mean = 0, sd = 1) + 0.3*dnorm(x, mean = 4, sd = 1) a <- 0.25 (f_a <- hdrlevels(f, prob = 1-a)) plot(x, f) abline(h = f_a, lty = 2) text(max(x), f_a, labels = paste0("f_", a), pos = 3) mean(f > f_a) # find the regions of HDR ord <- order(x) f <- f[ord] x <- x[ord] x_a <- x[f > f_a] j <- which.max(diff(x_a)) region1 <- x_a[c(1,j)] region2 <- x_a[c(j+1,length(x_a))] plot(x, f, type = "l") abline(h = f_a, lty = 2) abline(v = region1, lty = 3, col = 2) abline(v = region2, lty = 3, col = 3) 

### Example output

Package 'mclust' version 5.4.3
Type 'citation("mclust")' for citing this R package in publications.
50%        75%        90%
0.31878769 0.19683029 0.08894239
[1] 0.5
[1] -0.6671619  0.6681078
[1] 0.6744898
[1] 0.75
[1] -1.184999  1.187529
[1] 1.150349
[1] 0.9
[1] -1.726935  1.731267
[1] 1.644854
75%
0.09291342
[1] 0.75


mclust documentation built on Nov. 20, 2020, 5:09 p.m.