#' Orthonormal Contrast Matrices for Bayesian Estimation
#'
#' Returns a design or model matrix of orthonormal contrasts such that the
#' marginal prior on all effects is identical. Implementation from Singmann &
#' Gronau's [`bfrms`](https://github.com/bayesstuff/bfrms/), following
#' the description in Rouder, Morey, Speckman, & Province (2012, p. 363).
#' \cr\cr
#' Though using this factor coding scheme might obscure the interpretation of
#' parameters, it is essential for correct estimation of Bayes factors for
#' contrasts and order restrictions of multi-level factors (where `k>2`). See
#' info on specifying correct priors for factors with more than 2 levels in
#' [the
#' Bayes factors vignette](https://easystats.github.io/bayestestR/articles/bayes_factors.html).
#'
#' @inheritParams stats::contr.treatment
#'
#' @details
#' When `contrasts = FALSE`, the returned contrasts are equivalent to
#' `contr.treatment(, contrasts = FALSE)`, as suggested by McElreath (also known
#' as one-hot encoding).
#'
#' @references
#' - McElreath, R. (2020). Statistical rethinking: A Bayesian course with
#'   examples in R and Stan. CRC press.
#'
#' - Rouder, J. N., Morey, R. D., Speckman, P. L., & Province, J. M. (2012).
#'   Default Bayes factors for ANOVA designs. *Journal of Mathematical
#'   Psychology*, 56(5), 356-374. https://doi.org/10.1016/j.jmp.2012.08.001
#'
#' @return A `matrix` with n rows and k columns, with k=n-1 if contrasts is
#'   `TRUE` and k=n if contrasts is `FALSE`.
#'
#' @aliases contr.bayes
#'
#' @examples
#' contr.orthonorm(2) # Q_2 in Rouder et al. (2012, p. 363)
#'
#' contr.orthonorm(5) # equivalent to Q_5 in Rouder et al. (2012, p. 363)
#'
#' ## check decomposition
#' Q3 <- contr.orthonorm(3)
#' Q3 %*% t(Q3) ## 2/3 on diagonal and -1/3 on off-diagonal elements
#' @export
contr.orthonorm <- function(n, contrasts = TRUE, sparse = FALSE) {
  contr <- stats::contr.treatment(n,
    contrasts = FALSE, base = 1,
    sparse = sparse & !contrasts
  )

  if (contrasts) {
    n <- ncol(contr)
    I_a <- diag(n)
    J_a <- matrix(1, nrow = n, ncol = n)
    Sigma_a <- I_a - J_a / n
    contr <- eigen(Sigma_a)$vectors[, seq_len(n - 1), drop = FALSE]
  }

  contr
}

# ----------

#' @export
contr.bayes <- function(n, contrasts = TRUE) {
  .Deprecated(new = "contr.orthonorm", old = "contr.bayes")
  contr.orthonorm(n, contrasts = contrasts)
}