## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(echo = TRUE, warning = FALSE, message = FALSE)

## ----eval = FALSE-------------------------------------------------------------
# install.packages("DEmixR")

## -----------------------------------------------------------------------------
# Load the package
library(DEmixR)

# Simulation data
set.seed(123)
x <- c(rnorm(3000, mean = 0, sd = 1),
       rnorm(2000, mean = 3, sd = 1.2))

# Diagnostic plots
prelim_plots(x, which = c("hist", "qq"))

## -----------------------------------------------------------------------------
# Automatically choose the best mixture family
best_model <- select_best_mixture(x, n_runs = 3, NP = 50, itermax = 2000)
best_model$best$family
best_model$BICs

## -----------------------------------------------------------------------------
# Fit a two-component normal mixture
fit <- fit_norm2(x, n_runs = 5, quiet = 2, parallelType = 0)

fit$par      # estimated parameters
fit$logLik   # log-likelihood
fit$AIC      # Akaike Information Criterion
fit$BIC      # Bayesian Information Criterion

## -----------------------------------------------------------------------------
# Simulation data
set.seed(123)
y <- c(rlnorm(3000, meanlog = 0, sdlog = 0.5),
       rlnorm(2000, meanlog = 1, sdlog = 0.7))

# Fit a two-component lognormal mixture
fit_ln <- fit_lognorm2(y, n_runs = 5, parallelType = 0)
fit_ln$par

## -----------------------------------------------------------------------------
# Bootstrap confidence intervals
boot_res <- bootstrap_mix2(fit_ln, B = 100, parametric = TRUE, parallelType = 0)
boot_res$central
boot_res$ci

## -----------------------------------------------------------------------------
evaluate_init(par_init = c(0.5, 0, 0.6, 1.1, 0.8), y, family = "lognormal")

## -----------------------------------------------------------------------------
# Simulated biomarker data with two populations
set.seed(456)
biomarker <- c(rlnorm(4000, meanlog = 2.5, sdlog = 0.4),  # Healthy group
               rlnorm(1000, meanlog = 3.8, sdlog = 0.6))  # Disease group

# Initial exploration
prelim_plots(biomarker, which = c("hist", "logqq"))

# Fit lognormal mixture
fit_bio <- try(fit_lognorm2(biomarker, n_runs = 5, parallelType = 0, quiet = 2))

if (!inherits(fit_bio, "try-error")) {
  cat("Biomarker analysis results:\n")
  print(fit_bio$par)
  
  # Estimate proportion in each group
  cat("\nEstimated proportion in group 1 (healthy):", 
      round(fit_bio$par[1], 3), "\n")
  cat("Estimated proportion in group 2 (disease):", 
      round(1 - fit_bio$par[1], 3), "\n")
  
  # Bootstrap for confidence intervals
  boot_bio <- try(bootstrap_mix2(fit_bio, B = 100, quiet = 2))
  if (!inherits(boot_bio, "try-error")) {
    cat("\n95% CI for proportion in group 1:\n")
    print(boot_bio$ci[, "p"])
  }
}