## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(collapse = TRUE, comment = "#>", eval = FALSE)

## ----eval = TRUE--------------------------------------------------------------
library("msPCA")
data(snp500)

dim(snp500)
round(snp500[1:4, 1:4], 3)

## ----eval = TRUE--------------------------------------------------------------
ev <- eigen(snp500, symmetric = TRUE, only.values = TRUE)$values
sum(ev > 1e-8)          # 422: the market direction has been removed

## -----------------------------------------------------------------------------
# ks_grid <- seq(5, 35, by = 5)
# 
# results <- lapply(ks_grid, function(k) {
#   set.seed(42)
#   res <- mspca(snp500, r = 4, ks = rep(k, 4), verbose = FALSE,
#                maxIter = 100, feasibilityConstraintType = 0)
#   data.frame(
#     k      = k,
#     fve    = fraction_variance_explained(snp500, res$x_best),
#     orth   = feasibility_violation_off(snp500, res$x_best, 0),
#     pwcorr = feasibility_violation_off(snp500, res$x_best, 1)
#   )
# })
# results_df <- do.call(rbind, results)

## ----eval = TRUE--------------------------------------------------------------
res_grid <- read.csv(system.file("vignette-data", "snp_varyingk_results.csv",
                                 package = "msPCA"))

ks  <- sort(unique(res_grid$k))
by_constraint <- function(cn) {
  sub <- res_grid[res_grid$constraint == cn, ]
  sub[match(ks, sub$k), c("fve", "orth_violation")]
}

tab <- cbind(k = ks,
             by_constraint("orthogonality"),
             by_constraint("zero-correlation"),
             by_constraint("nsprcomp"))

knitr::kable(
  tab, digits = 4, row.names = FALSE,
  col.names = c("k", "FVE (msPCA - orth)", "orth viol (msPCA - orth)", "FVE (msPCA - zero-corr)",
                "orth viol (msPCA - zero-corr)", "FVE (nsprcomp)", "orth viol (nsprcomp)"),
  caption = paste("FVE and orthogonality violation across the sparsity grid.",
                  "Every violation column reports the orthogonality violation,",
                  "including for the fits run under the zero-correlation constraint.")
)

## ----echo = FALSE, eval = TRUE, out.width = "32%", fig.show = "hold"----------
knitr::include_graphics(c("figures/snp_fve.png",
                          "figures/snp_orth.png",
                          "figures/snp_pwcorr.png"))

## -----------------------------------------------------------------------------
# set.seed(42)
# res_orth <- mspca(snp500, r = 4, ks = rep(10, 4), verbose = FALSE,
#                   maxIter = 100, feasibilityConstraintType = 0)
# set.seed(42)
# res_corr <- mspca(snp500, r = 4, ks = rep(10, 4), verbose = FALSE,
#                   maxIter = 100, feasibilityConstraintType = 1)
# 
# print(res_orth)
# print(res_corr)

## -----------------------------------------------------------------------------
# summary(res_orth)
# summary(res_corr)

## ----echo = FALSE, eval = TRUE, out.width = "95%"-----------------------------
knitr::include_graphics("figures/snp_heatmap.png")