mxPBF_cov: Change Point Detection in Covaraiance Structure using Maximum...

In hdbcp: Bayesian Change Point Detection for High-Dimensional Data

Change Point Detection in Covaraiance Structure using Maximum Pairwise Bayes Factor (mxPBF)

Description

This function detects change points in the covariance structure of multivariate Gaussian data using the Maximum Pairwise Bayes Factor (mxPBF). The function selects alpha that controls the empirical False Positive Rate (FPR), as suggested in the paper. One can conduct a multiscale approach using the function majority_rule_mxPBF().

Usage

mxPBF_cov(
  given_data,
  nws,
  alps,
  a0 = 0.01,
  b0 = 0.01,
  FPR_want = 0.05,
  n_sample = 300,
  n_cores = 1
)

Arguments

given_data	An (n \times p) data matrix representing n observations and p variables.
nws	A set of window sizes for change point detection.
alps	A grid of alpha values used in the empirical False Positive Rate (FPR) method.
a0	A hyperparameter a_0 used in the mxPBF (default: 0.01).
b0	A hyperparameter b_0 used in the mxPBF (default: 0.01).
FPR_want	Desired False Positive Rate for selecting alpha, used in the empirical FPR method (default: 0.05).
n_sample	Number of simulated samples to estimate the empirical FPR, used in the empirical FPR method (default: 300).
n_cores	Number of threads for parallel execution via OpenMP (default: 1).

Value

A list of length equal to the number of window sizes provided. Each element in the list contains:

Change_points

Locations of detected change points.

Bayes_Factors

Vector of calculated Bayes Factors for each middle points.

Selected_alpha

Optimal alpha value selected based on the method that controls the empirical FPR.

Window_size

Window size used for change point detection.

Examples

nws <- c(25, 60, 100)
alps <- seq(1,10,0.05)
## H0 data
mu <- rep(0,10)
sigma1 <- diag(10)
X <- mvrnorm_cpp(500, mu, sigma1)
res1 <- mxPBF_cov(X, nws, alps)

## H1 data
mu <- rep(0,10)
sigma2 <- diag(10)
for (i in 1:10) {
  for (j in i:10) {
    if (i == j) {
    next
    } else {
    cov_value <- rnorm(1, 1, 1)
    sigma2[i, j] <- cov_value
    sigma2[j, i] <- cov_value
    }
  }
}
sigma2 <- sigma2 + (abs(min(eigen(sigma2)$value))+0.1)*diag(10) # Make it nonsingular
Y1 <- mvrnorm_cpp(250, mu, sigma1)
Y2 <- mvrnorm_cpp(250, mu, sigma2)
Y <- rbind(Y1, Y2)
res2 <- mxPBF_cov(Y, nws, alps)

hdbcp documentation built on April 3, 2025, 7:39 p.m.