diff --git a/.Rbuildignore b/.Rbuildignore index 95b0e70..1324998 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -6,3 +6,9 @@ ^docs$ ^pkgdown$ ^CRAN-SUBMISSION$ +^doc$ +^Meta$ +^examples/temp/ +^tests/temp/ + + diff --git a/.gitignore b/.gitignore index 8c83bbc..ae997c6 100644 --- a/.gitignore +++ b/.gitignore @@ -7,3 +7,5 @@ .o .so .dll +/doc/ +/Meta/ diff --git a/CRAN-SUBMISSION b/CRAN-SUBMISSION index 40039da..efae06a 100644 --- a/CRAN-SUBMISSION +++ b/CRAN-SUBMISSION @@ -1,3 +1,3 @@ -Version: 2.1 -Date: 2024-09-14 10:57:51 UTC -SHA: eb5656ea33fc519c7f43ebe74b0963825b5943c0 +Version: 2.2 +Date: 2024-11-13 18:00:38 UTC +SHA: 323e03f56d02eec77b7c9d98aa46ed597b3f3b74 diff --git a/DESCRIPTION b/DESCRIPTION index b23cbc2..f12afd9 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: biplotEZ Title: EZ-to-Use Biplots -Version: 2.1 +Version: 2.3 Authors@R: c(person("Sugnet", "Lubbe", , "muvisu@sun.ac.za", role = c("aut", "cre","cph"), comment = c(ORCID = "0000-0003-2762-9944")), @@ -27,7 +27,7 @@ Description: Provides users with an EZ-to-use platform for representing data License: MIT + file LICENSE Encoding: UTF-8 Roxygen: list(markdown = TRUE) -RoxygenNote: 7.3.2 +RoxygenNote: 7.3.3 VignetteBuilder: knitr Config/testthat/edition: 3 Suggests: @@ -40,14 +40,17 @@ Suggests: knitr, MASS, R.devices, + RColorBrewer, rgl, - rmarkdown, - testthat (>= 3.0.0) + rmarkdown Imports: graphics, grDevices, plotrix, + Rcpp, splines, stats, withr NeedsCompilation: yes +LinkingTo: + Rcpp diff --git a/NAMESPACE b/NAMESPACE index 32ac46e..6697b7c 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -3,6 +3,7 @@ S3method(AoD,biplot) S3method(CA,biplot) S3method(CVA,biplot) +S3method(CatPCA,biplot) S3method(PCA,biplot) S3method(PCO,biplot) S3method(classification,biplot) @@ -12,15 +13,14 @@ S3method(regress,biplot) S3method(summary,biplot) export(AoD) export(CA) -export(CATPCA) export(CLPs) export(CLRs) export(CVA) +export(CatPCA) export(PCA) export(PCO) export(alpha.bags) export(axes) -export(axes_coordinates) export(biplot) export(classification) export(classify) @@ -34,6 +34,8 @@ export(legend.type) export(means) export(newaxes) export(newsamples) +export(nom.axes) +export(ord.axes) export(prediction) export(reflect) export(regress) @@ -41,6 +43,7 @@ export(rotate) export(samples) export(sqrtManhattan) export(translate_axes) +importFrom(Rcpp,sourceCpp) importFrom(grDevices,chull) importFrom(grDevices,colorRampPalette) importFrom(grDevices,grey) diff --git a/R/CA.R b/R/CA.R index 84aa44d..bab9577 100644 --- a/R/CA.R +++ b/R/CA.R @@ -18,13 +18,19 @@ #' \item{Z}{Combined data frame of the row and column coordinates.} #' \item{r}{Numer of levels in the row factor.} #' \item{c}{Numer of levels in the column factor.} +#' \item{Dr}{Diagonal matrix of row profiles.} +#' \item{Dc}{Diagonal matrix of column profiles.} +#' \item{Drh}{Weighted row profiles.} +#' \item{Dch}{Weighted column profiles.} #' \item{rowcoor}{Row coordinates based on the selected \code{variant}.} #' \item{colcoor}{Column coordinates based on the selected \code{variant}.} -#' \item{P}{Correspondence Matrix} +#' \item{P}{Correspondence Matrix.} #' \item{Smat}{Standardised Pearson residuals.} #' \item{SVD}{Singular value decomposition solution: \code{d, u, v}.} -#' \item{qual}{Quality of the approximation.} +#' \item{e.vects}{Depending on what was specified in \code{CA} argument.} +#' \item{dim.biplot}{The dimension of the biplot.} #' \item{lambda.val}{The computed lambda value if lambda-scaling is requested.} +#' \item{gamma}{Contribution of the singular values, based on the CA variant.} #' #' @seealso [biplot()] #' @@ -130,10 +136,10 @@ CA.biplot <- function(bp, dim.biplot = c(2,1,3), e.vects = 1:ncol(bp$X), variant bp$Z <- Z bp$r <- r bp$c <- c + bp$Dr <- Dr bp$Dc <- Dc - bp$Rc <- Dr - bp$Dch <- Dch bp$Drh <- Drh + bp$Dch <- Dch bp$rowcoor <- rowcoor bp$colcoor <- colcoor bp$P <- P diff --git a/R/CATPCA.R b/R/CATPCA.R index 57fd524..aae96a0 100644 --- a/R/CATPCA.R +++ b/R/CATPCA.R @@ -1,40 +1,253 @@ -#' Categorical Principal Component Analysis +# ------------------------------------------------------------------------------------------- +#' Perform Categorical Principal Components Analysis (CatPCA) +#' +#' @description +#' This function appends the \code{biplot} object with elements resulting from performing CatPCA. #' #' @param bp an object of class \code{biplot} obtained from preceding function \code{biplot()}. -#' @param dim.biplot dimension of the biplot. Only values 1, 2 and 3 are accepted, with default \code{2}. -#' @param e.vects which eigenvectors (principal components) to extract, with default \code{1:dim.biplot}. -#' @param group.aes vector of the same length as the number of rows in the data matrix +#' @param cat.type a vector indicating whether each categorical variables is \code{nominal} or \code{ordinal}. +#' @param dim.biplot the dimension of the biplot. Only values \code{1}, \code{2} and \code{3} are accepted, with default \code{2}. +#' @param e.vects the vector indicating which eigenvectors (principal components) should be plotted in the biplot, with default \code{1:dim.biplot}. +#' @param group.aes a vector of the same length as the number of rows in the data matrix #' for differentiated aesthetics for samples. -#' @param show.class.means logical, indicating whether group means should be plotted in the biplot. +#' @param show.class.means a logical value indicating whether group means should be plotted in the biplot. +#' @param epsilon convergence criteria for iterative algorithm with defaults 1e-6. +#' +#' @return An object of class CatPCA with the following elements: +#' \item{X}{the matrix of the centered and scaled numeric variables.} +#' \item{Xcat}{the data frame of the categorical variables.} +#' \item{raw.X}{the original data.} +#' \item{classes}{the vector of category levels for the class variable. This is to be used for \code{colour}, \code{pch} and \code{cex} specifications.} +#' \item{na.action}{the vector of observations that have been removed.} +#' \item{center}{a logical value indicating whether \eqn{\mathbf{X}} is centered.} +#' \item{scaled}{a logical value indicating whether \eqn{\mathbf{X}} is scaled.} +#' \item{means}{the vector of means for each numerical variable.} +#' \item{sd}{the vector of standard deviations for each numerical variable.} +#' \item{n}{the number of observations.} +#' \item{p}{the number of variables.} +#' \item{group.aes}{the vector of category levels for the grouping variable. This is to be used for \code{colour}, \code{pch} and \code{cex} specification.} +#' \item{g.names}{the descriptive names to be used for group labels.} +#' \item{g}{the number of groups.} +#' \item{Title}{the title of the biplot rendered.} +#' \item{Z}{the matrix with each row containing the details of the points that are plotted (i.e. coordinates).} +#' \item{Lmat}{the matrix for transformation to the principal components.} +#' \item{Linv}{the inverse of \eqn{\mathbf{L}}.} +#' \item{eigenvalues}{the vector of eigenvalues of the covariance matrix of \eqn{\mathbf{X}}.} +#' \item{all.ax.one.unit}{one unit in the positive direction of each biplot axis.} +#' \item{ax.type}{one of numeric, ordinal or nominal.} +#' \item{ax.one.unit}{one unit in the positive direction of each numeric biplot axis.} +#' \item{nom.levels}{a list with the levels of each nominal variable} +#' \item{ord.levels}{a list with the levels of each ordinal variable} +#' \item{e.vects}{the vector indicating which principal components are plotted in the biplot.} +#' \item{dim.biplot}{the dimension of the biplot.} +#' \item{class.means}{a logical value indicating whether group means are plotted in the biplot.} +#' \item{Zmeans}{the matrix of class mean coordinates that are plotted in the biplot.} +#' +#' @usage CatPCA(bp, cat.type = NULL, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), +#' group.aes = NULL, show.class.means = FALSE, epsilon = 1e-6) +#' @aliases CatPCA #' -#' @return -#' an object of class \code{biplot} -#' #' @export -#' @usage CATPCA(bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), -#' group.aes = NULL, show.class.means = FALSE) -#' @aliases CATPCA +#' +#' @examples +#' biplot(iris) |> CatPCA() +#' # create a CatPCA biplot +#' biplot(data = iris) |> CatPCA() |> plot() #' +CatPCA <- function (bp, cat.type = NULL, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), group.aes=NULL, + show.class.means = FALSE, epsilon = 1e-6) +{ + UseMethod("CatPCA") +} + +#' Calculate elements for the CatPCA biplot +#' +#' @description This function performs calculations for the construction of a CatPCA biplot. +#' +#' @inheritParams CatPCA +#' +#' @return an object of class \code{CatPCA}, inherits from classes \code{PCA} and \code{biplot}. +#' @export +#' #' @examples -#' CATPCA (iris) -CATPCA<-function(bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), - group.aes = NULL, show.class.means = FALSE){ - indmat <- function (groep.vec) - { elements <- levels(factor(groep.vec)) - Y <- matrix(0, nrow = length(groep.vec), ncol = length(elements)) - dimnames(Y) <- list(NULL, paste(elements)) - for (i in 1:length(elements)) Y[groep.vec == elements[i], i] <- 1 - return(Y) - } - - biplot.check.G <- function (G, n) - { if (is.null(G)) { G <- matrix(indmat(rep(1, n)), ncol = 1) - dimnames(G) <- list(1:n, "AllData") } - if (nrow(G) != n) stop("number of rows of X and G differ") - if (is.null(dimnames(G))) dimnames(G) <- list(NULL, paste("class", 1:ncol(G), sep = "")) - if (length(dimnames(G)[[2]]) == 0) dimnames(G)[[2]] <- paste("class", 1:ncol(G), sep = "") - if (ncol(G) == 1) class.vec <- rep(dimnames(G)[[2]], n) else class.vec <- apply(t(apply(G, 1, function(x) x == max(x))), 1, function(s, G) dimnames(G)[[2]][s], G = G) - G - } - ##Possibly combine these two functions to obtain indicator matrices and overall indicator matrix (could be useful for MCA as well) +#' biplot(iris) |> CatPCA() +#' # create a CatPCA biplot +#' biplot(data = iris) |> CatPCA() |> plot() +#' +CatPCA.biplot <- function (bp, cat.type = NULL, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), group.aes=NULL, + show.class.means = FALSE, epsilon = 1e-6) +{ + dim.biplot <- dim.biplot[1] + if (dim.biplot != 1 & dim.biplot != 2 & dim.biplot != 3) stop("Only 1D, 2D and 3D biplots") + e.vects <- e.vects[1:dim.biplot] + if (!is.null(group.aes)) + { bp$group.aes <- factor(group.aes) + bp$g.names <-levels(factor(group.aes)) + bp$g <- length(bp$g.names) + } + + if (!bp$center) + { warning("Any form of PCA requires a centred datamatrix. Your numeric data was centred before computation. Use center = TRUE in the call to biplot()") + bp <- biplot (bp$X, center = TRUE, scaled=bp$scaled, bp$classes, bp$group.aes) + } + X <- bp$X + n <- bp$n + p <- bp$p + p2 <- bp$p2 + Xcat <- bp$Xcat + if (is.null(cat.type)) cat.type <- rep("nom", p2) + if (length(cat.type) != p2) stop ("The length of cat.type needs to match the number of columns in Xcat") + cat.type <- c("nominal","ordinal")[pmatch (cat.type, c("nominal", "ordinal"), duplicates.ok = TRUE)] + if (any(is.na(cat.type))) stop ("Only 'nominal' or 'ordinal' allowed for cat.type.") + + Gk <- vector("list", p2) + Lk.mat <- vector("list", p2) + Lk <- rep(NA, p2) + zlist <- vector("list", p2) + Hmat <- matrix(nrow = n, ncol = p2) + CatLabels <- NULL + + for (k in 1:p2) + { Gk[[k]] <- indmat(Xcat[, k]) + CatLabels <- append(CatLabels, dimnames(Gk[[k]])[[2]]) + Lk[k] <- ncol(Gk[[k]]) + Lk.mat[[k]] <- t(Gk[[k]]) %*% Gk[[k]] + } + + Gmat <- matrix(unlist(Gk), nrow = n) + L.mat <- diag(diag(t(Gmat) %*% Gmat)) + L.mat.min.half <- diag(diag(L.mat)^-0.5) + N.mat <- diag(n) - matrix(1, nrow = n, ncol = n)/n + svd.mca <- svd(L.mat.min.half %*% t(Gmat) %*% N.mat %*% Gmat %*% L.mat.min.half/p) + z.ini <- svd.mca$v[, 1] + names(z.ini) <- CatLabels + count <- 0 + for (k in 1:p2) + { zlist[[k]] <- as.vector(z.ini[(count + 1):(count + ncol(Gk[[k]]))]) + names(zlist[[k]]) <- CatLabels[(count + 1):(count + ncol(Gk[[k]]))] + count <- count + ncol(Gk[[k]]) + if (identical(cat.type[k], "ordinal")) + { + gradient <- sign(stats::lm(zlist[[k]] ~ I(1:length(zlist[[k]])))$coef[2]) + if (gradient>0) zlist[[k]] <- ties (d = 1:length(zlist[[k]]), delta=zlist[[k]]) + if (gradient<0) zlist[[k]] <- rev(ties (d = 1:length(zlist[[k]]), delta=rev(zlist[[k]]))) + if (gradient==0) zlist[[k]] <- mean(zlist[[k]]) + } + } + for (k in 1:p2) + { labs <- names(zlist[[k]]) + scaling.factor <- as.numeric(sqrt(t(zlist[[k]]) %*% t(Gk[[k]]) %*% N.mat %*% Gk[[k]] %*% zlist[[k]])) + if (is.na(1/scaling.factor)) scaling.factor <- 1 + scaling.vec <- rep(scaling.factor, length(zlist[[k]])) + if (identical(cat.type[k], "nominal")) zlist[[k]] <- as.vector(zlist[[k]]/scaling.vec) + if (identical(cat.type[k], "ordinal")) + { u.temp <- sum(diag(Lk.mat[[k]])) + v.temp <- sum(zlist[[k]] * diag(Lk.mat[[k]])) + w.temp <- c(t(zlist[[k]]) %*% Lk.mat[[k]] %*% zlist[[k]]) + b.temp <- sqrt(1/(w.temp - v.temp^2/u.temp)) + zlist[[k]] <- b.temp * zlist[[k]] - b.temp * v.temp/u.temp + } + names(zlist[[k]]) <- labs + Hmat[, k] <- Gk[[k]] %*% zlist[[k]] + } + Hmat <- as.matrix(Hmat) + dimnames(Hmat) <- dimnames(Xcat) + again <- TRUE + RSS.old <- NULL + Hmat.cent <- cbind(N.mat %*% Hmat, X) + iter <- 0 + while (again) + { iter <- iter + 1 + svd.out <- svd(Hmat.cent) + if (dim.biplot>1) + Yr <- scale(svd.out$u[, e.vects, drop=F] %*% diag(svd.out$d[e.vects]) %*% t(svd.out$v[, e.vects, drop=F]), scale = FALSE, center = TRUE) + else + Yr <- scale((svd.out$u[, e.vects, drop=F] * svd.out$d[e.vects]) %*% t(svd.out$v[, e.vects, drop=F]), scale = FALSE, center = TRUE) + RSS <- sum((Hmat.cent - Yr)^2) + if (!is.null(RSS.old)) again <- ((RSS.old - RSS) > epsilon) + RSS.old <- RSS + if (again) { + for (k in 1:p2) { + labs <- names(zlist[[k]]) + if (identical(cat.type[k], "nominal")) { + scale.factor <- as.numeric(sqrt(t(Yr[, k]) %*% Gk[[k]] %*% solve(Lk.mat[[k]]) %*% t(Gk[[k]]) %*% Yr[, k])) + zlist[[k]] <- as.vector((solve(Lk.mat[[k]]) %*% t(Gk[[k]]) %*% Yr[, k])/scale.factor) + } + if (identical(cat.type[k], "ordinal")) { + zlist[[k]] <- as.vector((solve(Lk.mat[[k]]) %*% t(Gk[[k]]) %*% Yr[, k])) + gradient <- sign(stats::lm(zlist[[k]] ~ I(1:length(zlist[[k]])))$coef[2]) + if (gradient>0) zlist[[k]] <- ties (d = 1:length(zlist[[k]]), delta=zlist[[k]]) + if (gradient<0) zlist[[k]] <- rev(ties (d = 1:length(zlist[[k]]), delta=rev(zlist[[k]]))) + if (gradient==0) zlist[[k]] <- mean(zlist[[k]]) + u.temp <- sum(diag(Lk.mat[[k]])) + v.temp <- sum(zlist[[k]] * diag(Lk.mat[[k]])) + w.temp <- c(t(zlist[[k]]) %*% Lk.mat[[k]] %*% zlist[[k]]) + b.temp <- sqrt(1/(w.temp - v.temp^2/u.temp)) + zlist[[k]] <- b.temp * zlist[[k]] - b.temp * v.temp/u.temp + } + names(zlist[[k]]) <- labs + Hmat[, k] <- Gk[[k]] %*% zlist[[k]] + } + } + Hmat.cent <- cbind(N.mat %*% Hmat, X) + } + Lz <- rep(NA, p2) + for (i in 1:p2) Lz[i] <- sum(Lk.mat[[i]] %*% zlist[[i]]) + zLz <- rep(NA, p2) + for (i in 1:p2) zLz[i] <- sum(t(zlist[[i]]) %*% Lk.mat[[i]] %*% zlist[[i]]) + Hmat.cent.svd <- svd(Hmat.cent) + + V.mat <- Hmat.cent.svd$v + U.mat <- Hmat.cent.svd$u + Sigma.mat <- diag(Hmat.cent.svd$d) + Lmat <- Hmat.cent.svd$v + Vr <- Hmat.cent.svd$v[, e.vects, drop = FALSE] + Z <- Hmat.cent %*% Hmat.cent.svd$v[, e.vects, drop=F] + rownames(Z) <- rownames(X) + + ax.one.unit <- 1/(diag(Vr %*% t(Vr))) * Vr + nom.levels <- ord.levels <- NULL + nom.count <- ord.count <- 0 + for (k in 1:p2) + { + if (cat.type[k] == "nominal") + { + nom.count <- nom.count + 1 + nom.levels[[nom.count]] <- zlist[[k]] + } + if (cat.type[k] == "ordinal") + { + ord.count <- ord.count + 1 + ord.levels[[ord.count]] <- zlist[[k]] + } + } + + ax.type <- c(cat.type, rep("numeric",p)) + num.levels <- apply(Xcat, 2, function(x) nlevels(factor(x))) + num.levels <- c(num.levels, rep(1,p)) + ax.type <- data.frame (ax.type, num.levels) + + bp$Z <- Z + bp$Lmat <- Lmat + bp$eigenvalues <- svd.out$d^2 + bp$all.ax.one.unit <- ax.one.unit + bp$ax.type <- ax.type + bp$ax.one.unit <- ax.one.unit[ax.type$ax.type=="numeric",] + bp$nom.levels <- nom.levels + bp$ord.levels <- ord.levels + bp$e.vects <- e.vects + bp$Vr <- Vr + bp$dim.biplot <- dim.biplot + if (bp$g == 1) bp$class.means <- FALSE else bp$class.means <- show.class.means + if (bp$class.means) + { + G <- indmat(bp$group.aes) + Xmeans <- solve(t(G)%*%G) %*% t(G) %*% X + Zmeans <- Xmeans %*% Lmat[,e.vects] + bp$Zmeans <- Zmeans + } + + class(bp) <- append(class(bp),"PCA") + class(bp) <- append(class(bp),"catPCA") + bp } diff --git a/R/PCA.R b/R/PCA.R index ae79c83..db3b562 100644 --- a/R/PCA.R +++ b/R/PCA.R @@ -14,6 +14,7 @@ #' optimally approximated in the biplot. If \code{TRUE}, the correlations between #' variables are optimally approximated by the cosine of the angles between #' axes. Default is \code{FALSE}. +#' @param ... additional arguments. #' #' @return An object of class PCA with the following elements: #' \item{X}{the matrix of the centered and scaled numeric variables.} @@ -44,8 +45,9 @@ #' #' @seealso [biplot()] #' -#' @usage PCA(bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), -#' group.aes = NULL, show.class.means = FALSE, correlation.biplot = FALSE) +#' @usage PCA(bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), +#' group.aes=NULL, show.class.means = FALSE, +#' correlation.biplot=FALSE, ...) #' @aliases PCA #' #' @export @@ -58,8 +60,9 @@ #' # create a PCA biplot #' biplot(data = iris) |> PCA() |> plot() #' -PCA <- function (bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), group.aes=NULL, show.class.means = FALSE, - correlation.biplot=FALSE) +PCA <- function (bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), + group.aes=NULL, show.class.means = FALSE, + correlation.biplot=FALSE, ...) { UseMethod("PCA") } @@ -79,7 +82,7 @@ PCA <- function (bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), group.aes= #' biplot(data = iris) |> PCA() |> plot() #' PCA.biplot <- function (bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), group.aes=NULL, - show.class.means = FALSE, correlation.biplot=FALSE) + show.class.means = FALSE, correlation.biplot=FALSE, ...) { dim.biplot <- dim.biplot[1] @@ -105,13 +108,15 @@ PCA.biplot <- function (bp, dim.biplot = c(2, 1, 3), e.vects = 1:ncol(bp$X), gro Lmat <- svd.out$v Vr <- svd.out$v[, e.vects, drop = FALSE] - if (correlation.biplot) - { - - if (dim.biplot>1) lambda.r <- diag(svd(t(X) %*% X)$d[1:dim.biplot]) - else lambda.r <- matrix(svd(t(X) %*% X)$d, nrow=1, ncol=1) + if (correlation.biplot) { + svals <- svd(t(X) %*% X)$d + if (dim.biplot > 1) + lambda.r <- diag(svals[1:dim.biplot]) + else lambda.r <- matrix(svals, nrow = 1, ncol = 1) + if (any(diag(lambda.r) < 1e-14)) stop (paste("Your data is of rank <", dim.biplot)) Z <- sqrt(n - 1) * X %*% Vr %*% sqrt(solve(lambda.r)) - Lmat <- sqrt(n-1) * Lmat %*% sqrt(solve(diag(svd(t(X) %*% X)$d))) + if (n > p) Lmat <- sqrt(n - 1) * Lmat %*% diag(sqrt(ifelse(svals < 1e-14, 0, 1/svals))) + else Lmat <- sqrt(n - 1) * Lmat %*% diag(sqrt(ifelse(svals < 1e-14, 0, 1/svals)))[1:n,1:n] } else { Z <- X %*% Vr } rownames(Z) <- rownames(X) diff --git a/R/RcppExports.R b/R/RcppExports.R new file mode 100644 index 0000000..d943b0d --- /dev/null +++ b/R/RcppExports.R @@ -0,0 +1,11 @@ +# Generated by using Rcpp::compileAttributes() -> do not edit by hand +# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393 + +alfunc <- function(BVEC, X, Y, M, MU, LAMBDA, CONST1, CONST2, U, V) { + .Call(`_biplotEZ_alfunc`, BVEC, X, Y, M, MU, LAMBDA, CONST1, CONST2, U, V) +} + +optimize_spline <- function(BVEC, X, Y, M, MU, LAMBDA, CONST1, CONST2, U, V, TAU, FTOL, TINY, ITMAX) { + .Call(`_biplotEZ_optimize_spline`, BVEC, X, Y, M, MU, LAMBDA, CONST1, CONST2, U, V, TAU, FTOL, TINY, ITMAX) +} + diff --git a/R/aesthetics.R b/R/aesthetics.R index 639154c..b281a6b 100644 --- a/R/aesthetics.R +++ b/R/aesthetics.R @@ -78,7 +78,6 @@ samples <- function (bp, which = 1:bp$g, col = ez.col, pch = 16, g <- bp$g n <- bp$n p <- bp$p - if(is.null(which) & length(col)==0) col <- ez.col if(inherits(bp,"CA")) @@ -194,7 +193,8 @@ samples <- function (bp, which = 1:bp$g, col = ez.col, pch = 16, #' This function allows the user to format the aesthetics for the class means or group means. #' #' @param bp an object of class \code{biplot}. -#' @param which a vector containing the groups or classes for which the means should be displayed, with default \code{bp$g}. +#' @param which a vector containing the groups or classes for which the means should be displayed, with default +#' the \code{which} argument of samples. #' @param col the colour(s) for the means, with default as the colour of the samples. #' @param pch the plotting character(s) for the means, with default \code{15}. #' @param cex the character expansion(s) for the means, with default \code{1}. @@ -225,7 +225,7 @@ samples <- function (bp, which = 1:bp$g, col = ez.col, pch = 16, #' @seealso [biplot()] #' #' @usage -#' means (bp, which = NULL, col = NULL, pch = 15, cex = 1, label = FALSE, +#' means (bp, which = bp$samples$which, col = NULL, pch = 15, cex = 1, label = FALSE, #' label.col = NULL,label.cex = 0.75, label.side = "bottom", label.offset = 0.5, #' opacity = 1, shade.darker = TRUE) #' @aliases means @@ -234,11 +234,10 @@ samples <- function (bp, which = 1:bp$g, col = ez.col, pch = 16, #' #' @examples biplot(iris[,1:4]) |> PCA() |> #' means(col = "purple", pch = 15, cex = 2) |> plot() -means <- function (bp, which = NULL, col = NULL, +means <- function (bp, which = bp$samples$which, col = NULL, pch = 15, cex = 1, label = FALSE, label.col=NULL, label.cex = 0.75, label.side = "bottom", label.offset = 0.5, opacity = 1, shade.darker = TRUE) { - bp$show.class.means <- TRUE if (is.null(bp$samples)) bp <- samples(bp) g <- bp$g @@ -250,10 +249,9 @@ means <- function (bp, which = NULL, col = NULL, which <- which[which <= g] which <- which[which > 0] } - else if (!is.null(bp$samples$which)) which <- bp$samples$which else which <- 1:g class.num <- length(which) - + if (!is.null(which)) if (is.null(col)) { if (max(which)<=length(bp$samples$col)) col <- bp$samples$col[which] else col <- bp$samples$col @@ -487,6 +485,419 @@ axes <- function (bp, X.names=colnames(bp$X), which = 1:bp$p, col = grey(0.7), l bp } +# ---------------------------------------------------------------------------------------------- +#' Format aesthetics for the nominal biplot axes +#' +#' @description +#' This function allows the user to format the aesthetics for the nominal biplot axes. +#' +#' +#' @param bp an object of class \code{biplot}. +#' @param X.names a vector of column names of \code{bp} to specify which axes should be labelled. +#' @param which a vector containing the columns or variables for which the axes should be displayed, with default \code{1:p}. +#' @param col a list with each component the colour(s) for the levels of that axis. +#' @param lwd the line width(s) for the axes, with default \code{1}. +#' @param lty the line type(s) for the axes, with default \code{1}. +#' @param label.dir a character string indicating the placement of the axis titles to the side of the figure. One of "\code{Orthog}" for axis titles to appear orthogonal to the side of the figure (default) , "\code{Hor}" for axis titles to appear horizontally or "\code{Paral}" for axis titles to appear parallel to the side of the figure. +#' @param label.col the colour(s) for the axis labels, with default, \code{col}. +#' @param label.cex the label expansion for the axis labels, with default \code{0.75}. +#' @param label.line the distance of the axis title from the side of the figure, with default \code{0.1}. +#' @param label.offset a four-component numeric vector controlling the distances axis titles are displayed from the side of the figure, with default \code{rep(0,4)}. Sides are numbered \code{1} to \code{4} according to \code{R} conventions. +#' @param ticks an integer-valued vector indicating the number of tickmarks for each axis, with default \code{5} for each axis. +#' @param tick.col the colour(s) for the tick marks, with default \code{col}. +#' @param tick.size a vector specifying the sizes of tick marks for each axis, with default \code{1} for each . +#' @param tick.label a logical value indicating whether the axes should be labelled, with default \code{TRUE}. +#' @param tick.label.side a character string indicating the position of the tick label. One of "\code{below}" for the label to appear below the tick mark (default) or "\code{above}" for the label to appear above the tick mark. +#' @param tick.label.col the colour(s) for the tick mark labels, with default \code{tick.col}. +#' @param tick.label.cex the label expansion for the tick mark labels, with default \code{0.6}. +#' @param predict.col the colour(s) for the predicted samples, with default \code{col}. +#' @param predict.lwd the line width(s) for the predicted samples, with default \code{lwd}. +#' @param predict.lty the line type(s) for the predicted samples, with default \code{lty}. +#' @param ax.names a vector of size \code{p} containing user defined titles for the axes. +#' @param orthogx a numeric vector of size \code{p} specifying the x-coordinate of the parallel transformation of each axis, with default \code{0} for each axis. This is only used when \code{dim.biplot = 2}. +#' @param orthogy a numeric vector of size \code{p} specifying the y-coordinate of the parallel transformation of each axis, with default \code{0} for each axis. This is only used when \code{dim.biplot = 2}. +#' +#' @return The object of class \code{biplot} will be appended with a list called \code{axes} containing the following elements: +#' \item{which}{a vector containing the columns for which the axes are displayed.} +#' \item{col}{the colour(s) of the axes.} +#' \item{lwd}{the line width(s) of the axes.} +#' \item{lty}{the line type(s) of the axes.} +#' \item{label.dir}{the placement of the axis titles to the side of the figure.} +#' \item{label.col}{the colour(s) of the axis titles.} +#' \item{label.cex}{the label expansion(s) of the axis titles.} +#' \item{label.line}{the distance(s) of the axis titles from the side of the figure.} +#' \item{ticks}{the number of tick marks per axis.} +#' \item{tick.col}{the colour(s) of the tick marks.} +#' \item{tick.size}{the size(s) of the tick marks.} +#' \item{tick.label}{logical value(s) indicating whether axes are labelled.} +#' \item{tick.label.side}{the position of the tick mark labels.} +#' \item{tick.label.col}{the colour(s) of the tick mark labels.} +#' \item{tick.label.cex}{the expansion(s) of the tick mark labels.} +#' \item{predict.col}{the colour(s) of the predicted samples.} +#' \item{predict.lty}{the line type(s) of the predicted samples.} +#' \item{predict.lwd}{the line width(s) of the predicted samples.} +#' \item{names}{the user defined axis titles.} +#' \item{orthogx}{the horizontal translations for each axis.} +#' \item{orthogy}{the vertical translations for each axis.} +#' \item{vectors}{a logical value indicating whether calibrated axes are plotted.} +#' +#' @seealso [biplot()] +#' +#' @usage +#' nom.axes(bp, X.names=colnames(bp$Xcat)[bp$ax.type$ax.type=="nominal"], +#' which =1:sum(bp$ax.type$ax.type=="nominal"), col = NULL, lwd = 1, +#' lty = 1, label.dir = "Orthog", label.col = "black", label.cex = 0.75, +#' label.line = 0.1, label.offset=rep(0,4), ticks = 5, tick.col = "black", +#' tick.size = 1, tick.label = TRUE, tick.label.side = "below", +#' tick.label.col = tick.col, tick.label.cex = 0.6, predict.col = col, +#' predict.lwd = lwd, predict.lty = lty, ax.names = X.names, +#' orthogx = 0, orthogy = 0) +#' @aliases nom.axes +#' +#' @export +#' +#' @examples +#' biplot(iris) |> CatPCA() |> nom.axes(col=c("blue","magenta","gold")) |> plot() +#' +nom.axes <- function (bp, X.names=colnames(bp$Xcat)[bp$ax.type$ax.type=="nominal"], + which =1:sum(bp$ax.type$ax.type=="nominal"), col = NULL, lwd = 1, + lty = 1, label.dir = "Orthog", label.col = "black", + label.cex = 0.75, + label.line = 0.1, label.offset=rep(0,4), ticks = 5, tick.col = "black", + tick.size = 1, tick.label = TRUE, tick.label.side = "below", + tick.label.col = tick.col, tick.label.cex = 0.6, + predict.col = col, predict.lwd = lwd, predict.lty = lty, ax.names = X.names, + orthogx = 0, orthogy = 0) +{ + if (!requireNamespace("RColorBrewer", quietly = TRUE)) { + stop("Package 'RColorBrewer' is required for this function. Please install it.", call. = FALSE) + } + + p <- sum(bp$ax.type$ax.type=="nominal") + num.levels.vec <- bp$ax.type$num.levels[bp$ax.type$ax.type == "nominal"] + if (!is.null(p)) + { + if (!all(is.numeric(which))) which <- match(which, X.names, nomatch = 0) + which <- which[which <= p] + which <- which[which > 0] + ax.num <- length(which) + + if (is.null(col)) + col <- c(RColorBrewer::brewer.pal(12,"Set3")[-2], RColorBrewer::brewer.pal(7,"Dark2")) + if (!is.list(col)) + { + temp <- vector ("list", ax.num) + for (j in 1:ax.num) temp[[j]] <- col + col <- temp + } + if (length(col) < ax.num) + { + temp <- vector ("list", ax.num) + col.length <- length(col) + j <- 0 + while (j < ax.num) + { for (k in 1:col.length) if (k+j ax.num) + { + temp <- vector ("list", ax.num) + for (j in 1:ax.num) temp[[j]] <- col[[j]] + col <- temp + } + num.levels <- num.levels.vec[which] + + for (j in 1:ax.num) + col[[j]] <- col[[j]][ifelse (1:num.levels[j] %% length(col[[j]]) == 0, length(col[[j]]), + 1:num.levels[j] %% length(col[[j]]))] + + lwd.len <- length(lwd) + lwd <- lwd[ifelse(1:ax.num%%lwd.len==0,lwd.len,1:ax.num%%lwd.len)] + if(is.null(lwd)){lwd <- rep(0, ax.num)} + + lty.len <- length(lty) + lty <- lty[ifelse(1:ax.num%%lty.len==0,lty.len,1:ax.num%%lty.len)] + if(is.null(lty)){lty <- rep(0, ax.num)} + if (label.dir != "Orthog" & label.dir != "Hor" & label.dir != "Paral") + stop("Incorrect specification of axis label direction") + + label.col.len <- length(label.col) + label.col <- label.col[ifelse(1:ax.num%%label.col.len==0,label.col.len,1:ax.num%%label.col.len)] + if(is.null(label.col)){label.col <- rep(0, ax.num)} + + label.cex.len <- length(label.cex) + label.cex <- label.cex[ifelse(1:ax.num%%label.cex.len==0,label.cex.len,1:ax.num%%label.cex.len)] + if(is.null(label.cex)){label.cex <- rep(0, ax.num)} + + label.line.len <- length(label.line) + label.line <- label.line[ifelse(1:ax.num%%label.line.len==0,label.line.len,1:ax.num%%label.line.len)] + if(is.null(label.line)){label.line <- rep(0, ax.num)} + + ticks.len <- length(ticks) + ticks <- ticks[ifelse(1:ax.num%%ticks.len==0,ticks.len,1:ax.num%%ticks.len)] + if(is.null(ticks)){ticks <- rep(0, ax.num)} + + tick.col.len <- length(tick.col) + tick.col <- tick.col[ifelse(1:ax.num%%tick.col.len==0,tick.col.len,1:ax.num%%tick.col.len)] + if(is.null(tick.col)){tick.col <- rep(0, ax.num)} + + tick.size.len <- length(tick.size) + tick.size <- tick.size[ifelse(1:ax.num%%tick.size.len==0,tick.size.len,1:ax.num%%tick.size.len)] + if(is.null(tick.size)){tick.size <- rep(0, ax.num)} + + tick.label.len <- length(tick.label) + tick.label <- tick.label[ifelse(1:ax.num%%tick.label.len==0,tick.label.len,1:ax.num%%tick.label.len)] + if(is.null(tick.label)){tick.label <- rep("", ax.num)} + + tick.label.col.len <- length(tick.label.col) + tick.label.col <- tick.label.col[ifelse(1:ax.num%%tick.label.col.len==0,tick.label.col.len,1:ax.num%%tick.label.col.len)] + if(is.null(tick.label.col)){tick.label.col <- rep(0, ax.num)} + + tick.label.cex.len <- length(tick.label.cex) + tick.label.cex <- tick.label.cex[ifelse(1:ax.num%%tick.label.cex.len==0,tick.label.cex.len,1:ax.num%%tick.label.cex.len)] + if(is.null(tick.label.cex)){tick.label.cex <- rep(0, ax.num)} + + tick.label.side.len <- length(tick.label.side) + tick.label.side <- tick.label.side[ifelse(1:ax.num%%tick.label.side.len==0,tick.label.side.len,1:ax.num%%tick.label.side.len)] + if(is.null(tick.label.side)){tick.label.side <- rep(0, ax.num)} + + predict.col.len <- length(predict.col) + predict.col <- predict.col[ifelse(1:ax.num%%predict.col.len==0,predict.col.len,1:ax.num%%predict.col.len)] + if(is.null(predict.col)){predict.col <- rep(0, ax.num)} + + predict.lwd.len <- length(predict.lwd) + predict.lwd <- predict.lwd[ifelse(1:ax.num%%predict.lwd.len==0,predict.lwd.len,1:ax.num%%predict.lwd.len)] + if(is.null(predict.lwd)){predict.lwd <- rep(0, ax.num)} + + predict.lty.len <- length(predict.lty) + predict.lty <- predict.lty[ifelse(1:ax.num%%predict.lty.len==0,predict.lty.len,1:ax.num%%predict.lty.len)] + if(is.null(predict.lty)){predict.lty <- rep(0, ax.num)} + + ax.names.len <- length(ax.names) + ax.names <- ax.names[ifelse(1:p%%ax.names.len==0,ax.names.len,1:p%%ax.names.len)] + if(is.null(ax.names)){ax.names <- rep("", p)} + + orthogx.len <- length(orthogx) + orthogx <- orthogx[ifelse(1:p%%orthogx.len==0,orthogx.len,1:p%%orthogx.len)] + if(is.null(orthogx)){orthogx <- rep("", p)} + + orthogy.len <- length(orthogy) + orthogy <- orthogx[ifelse(1:p%%orthogy.len==0,orthogy.len,1:p%%orthogy.len)] + if(is.null(orthogy)){orthogy <- rep("", p)} + + bp$nom.axes = list(which = which, col = col, lwd = lwd, lty = lty, label.dir = label.dir, + label.col = label.col, label.cex = label.cex, label.line = label.line, + ticks = ticks, tick.col = tick.col, tick.size = tick.size, + tick.label = tick.label, tick.label.col = tick.label.col, + tick.label.cex = tick.label.cex, tick.label.side=tick.label.side, + predict.col = predict.col, predict.lty = predict.lty, + predict.lwd = predict.lwd, names = ax.names, orthogx = orthogx, + orthogy = orthogy) + + } + else bp$nom.axes <- list(which = NULL) + bp +} + +# ---------------------------------------------------------------------------------------------- +#' Format aesthetics for the ordinal biplot axes +#' +#' @description +#' This function allows the user to format the aesthetics for the ordinal biplot axes. +#' +#' +#' @param bp an object of class \code{biplot}. +#' @param X.names a vector of column names of \code{bp} to specify which axes should be labelled. +#' @param which a vector containing the columns or variables for which the axes should be displayed, with default \code{1:p}. +#' @param col the colour(s) for the axes, with default \code{grey(0.7)}. Alternatively, provide a vector of colours corresponding to \code{X.names}. +#' @param lwd.factor a numeric value for increasing line width for each consecutive level +#' @param reverse logical indicator to switch the narrow to wide calibration on the axis +#' @param lwd the line width(s) for the axes, with default \code{1}. +#' @param lty the line type(s) for the axes, with default \code{1}. +#' @param label.dir a character string indicating the placement of the axis titles to the side of the figure. One of "\code{Orthog}" for axis titles to appear orthogonal to the side of the figure (default) , "\code{Hor}" for axis titles to appear horizontally or "\code{Paral}" for axis titles to appear parallel to the side of the figure. +#' @param label.col the colour(s) for the axis labels, with default, \code{col}. +#' @param label.cex the label expansion for the axis labels, with default \code{0.75}. +#' @param label.line the distance of the axis title from the side of the figure, with default \code{0.1}. +#' @param label.offset a four-component numeric vector controlling the distances axis titles are displayed from the side of the figure, with default \code{rep(0,4)}. Sides are numbered \code{1} to \code{4} according to \code{R} conventions. +#' @param ticks an integer-valued vector indicating the number of tickmarks for each axis, with default \code{5} for each axis. +#' @param tick.col the colour(s) for the tick marks, with default \code{col}. +#' @param tick.size a vector specifying the sizes of tick marks for each axis, with default \code{1} for each . +#' @param tick.label a logical value indicating whether the axes should be labelled, with default \code{TRUE}. +#' @param tick.label.side a character string indicating the position of the tick label. One of "\code{below}" for the label to appear below the tick mark (default) or "\code{above}" for the label to appear above the tick mark. +#' @param tick.label.col the colour(s) for the tick mark labels, with default \code{tick.col}. +#' @param tick.label.cex the label expansion for the tick mark labels, with default \code{0.6}. +#' @param predict.col the colour(s) for the predicted samples, with default \code{col}. +#' @param predict.lwd the line width(s) for the predicted samples, with default \code{lwd}. +#' @param predict.lty the line type(s) for the predicted samples, with default \code{lty}. +#' @param ax.names a vector of size \code{p} containing user defined titles for the axes. +#' @param orthogx a numeric vector of size \code{p} specifying the x-coordinate of the parallel transformation of each axis, with default \code{0} for each axis. This is only used when \code{dim.biplot = 2}. +#' @param orthogy a numeric vector of size \code{p} specifying the y-coordinate of the parallel transformation of each axis, with default \code{0} for each axis. This is only used when \code{dim.biplot = 2}. +#' +#' @return The object of class \code{biplot} will be appended with a list called \code{axes} containing the following elements: +#' \item{which}{a vector containing the columns for which the axes are displayed.} +#' \item{col}{the colour(s) of the axes.} +#' \item{lwd.factor}{the factor with which to increase the line width for each level.} +#' \item{reverse}{whether the axis goes narrow to wide or vice versa.} +#' \item{lwd}{the line width(s) of the axes.} +#' \item{lty}{the line type(s) of the axes.} +#' \item{label.dir}{the placement of the axis titles to the side of the figure.} +#' \item{label.col}{the colour(s) of the axis titles.} +#' \item{label.cex}{the label expansion(s) of the axis titles.} +#' \item{label.line}{the distance(s) of the axis titles from the side of the figure.} +#' \item{ticks}{the number of tick marks per axis.} +#' \item{tick.col}{the colour(s) of the tick marks.} +#' \item{tick.size}{the size(s) of the tick marks.} +#' \item{tick.label}{logical value(s) indicating whether axes are labelled.} +#' \item{tick.label.side}{the position of the tick mark labels.} +#' \item{tick.label.col}{the colour(s) of the tick mark labels.} +#' \item{tick.label.cex}{the expansion(s) of the tick mark labels.} +#' \item{predict.col}{the colour(s) of the predicted samples.} +#' \item{predict.lty}{the line type(s) of the predicted samples.} +#' \item{predict.lwd}{the line width(s) of the predicted samples.} +#' \item{names}{the user defined axis titles.} +#' \item{orthogx}{the horizontal translations for each axis.} +#' \item{orthogy}{the vertical translations for each axis.} +#' \item{vectors}{a logical value indicating whether calibrated axes are plotted.} +#' +#' @usage +#' ord.axes(bp, X.names=colnames(bp$Xcat)[bp$ax.type$ax.type=="ordinal"], +#' which =1:sum(bp$ax.type$ax.type=="ordinal"), col = grey(0.7), +#' reverse = rep(FALSE, sum(bp$ax.type$ax.type=="ordinal")), lwd.factor = 1.5, +#' lwd = 1, lty = 1, label.dir = "Orthog", label.col = col, label.cex = 0.75, +#' label.line = 0.1, label.offset=rep(0,4), ticks = 5, tick.col = col, +#' tick.size = 1, tick.label = TRUE, tick.label.side = "below", +#' tick.label.col = tick.col, tick.label.cex = 0.6, +#' predict.col = col, predict.lwd = lwd, predict.lty = lty, ax.names = X.names, +#' orthogx = 0, orthogy = 0) +#' @aliases ord.axes +#' +#' @export +#' +#' @examples +#' biplot(iris) |> CatPCA(cat.type = "ord") |> ord.axes(col="purple") |> plot() +#' +ord.axes <- function (bp, X.names=colnames(bp$Xcat)[bp$ax.type$ax.type=="ordinal"], + which =1:sum(bp$ax.type$ax.type=="ordinal"), col = grey(0.7), + reverse = rep(FALSE, sum(bp$ax.type$ax.type=="ordinal")), lwd.factor = 1.5, + lwd = 1, lty = 1, label.dir = "Orthog", label.col = col, label.cex = 0.75, + label.line = 0.1, label.offset=rep(0,4), ticks = 5, tick.col = col, + tick.size = 1, tick.label = TRUE, tick.label.side = "below", + tick.label.col = tick.col, tick.label.cex = 0.6, + predict.col = col, predict.lwd = lwd, predict.lty = lty, ax.names = X.names, + orthogx = 0, orthogy = 0) +{ + p <- sum(bp$ax.type$ax.type=="ordinal") + num.levels.vec <- bp$ax.type$num.levels[bp$ax.type$ax.type == "ordinal"] + if (!is.null(p)) + { + if (!all(is.numeric(which))) which <- match(which, X.names, nomatch = 0) + which <- which[which <= p] + which <- which[which > 0] + ax.num <- length(which) + + col.len <- length(col) + col <- col[ifelse(1:ax.num%%col.len==0,col.len,1:ax.num%%col.len)] + if(is.null(col)){col <- rep(0, ax.num)} + + lwdf.len <- length(lwd.factor) + lwd.factor <- lwd.factor[ifelse(1:ax.num%%lwdf.len==0,lwdf.len,1:ax.num%%lwdf.len)] + if(is.null(lwd.factor)){lwd.factor <- rep(0, ax.num)} + + lwd.len <- length(lwd) + lwd <- lwd[ifelse(1:ax.num%%lwd.len==0,lwd.len,1:ax.num%%lwd.len)] + if(is.null(lwd)){lwd <- rep(0, ax.num)} + + rev.len <- length(reverse) + reverse <- reverse[ifelse(1:ax.num%%rev.len==0,rev.len,1:ax.num%%rev.len)] + if(is.null(reverse)){reverse <- rep(FALSE, ax.num)} + + lty.len <- length(lty) + lty <- lty[ifelse(1:ax.num%%lty.len==0,lty.len,1:ax.num%%lty.len)] + if(is.null(lty)){lty <- rep(0, ax.num)} + if (label.dir != "Orthog" & label.dir != "Hor" & label.dir != "Paral") + stop("Incorrect specification of axis label direction") + + label.col.len <- length(label.col) + label.col <- label.col[ifelse(1:ax.num%%label.col.len==0,label.col.len,1:ax.num%%label.col.len)] + if(is.null(label.col)){label.col <- rep(0, ax.num)} + + label.cex.len <- length(label.cex) + label.cex <- label.cex[ifelse(1:ax.num%%label.cex.len==0,label.cex.len,1:ax.num%%label.cex.len)] + if(is.null(label.cex)){label.cex <- rep(0, ax.num)} + + label.line.len <- length(label.line) + label.line <- label.line[ifelse(1:ax.num%%label.line.len==0,label.line.len,1:ax.num%%label.line.len)] + if(is.null(label.line)){label.line <- rep(0, ax.num)} + + ticks.len <- length(ticks) + ticks <- ticks[ifelse(1:ax.num%%ticks.len==0,ticks.len,1:ax.num%%ticks.len)] + if(is.null(ticks)){ticks <- rep(0, ax.num)} + + tick.col.len <- length(tick.col) + tick.col <- tick.col[ifelse(1:ax.num%%tick.col.len==0,tick.col.len,1:ax.num%%tick.col.len)] + if(is.null(tick.col)){tick.col <- rep(0, ax.num)} + + tick.size.len <- length(tick.size) + tick.size <- tick.size[ifelse(1:ax.num%%tick.size.len==0,tick.size.len,1:ax.num%%tick.size.len)] + if(is.null(tick.size)){tick.size <- rep(0, ax.num)} + + tick.label.len <- length(tick.label) + tick.label <- tick.label[ifelse(1:ax.num%%tick.label.len==0,tick.label.len,1:ax.num%%tick.label.len)] + if(is.null(tick.label)){tick.label <- rep("", ax.num)} + + tick.label.col.len <- length(tick.label.col) + tick.label.col <- tick.label.col[ifelse(1:ax.num%%tick.label.col.len==0,tick.label.col.len,1:ax.num%%tick.label.col.len)] + if(is.null(tick.label.col)){tick.label.col <- rep(0, ax.num)} + + tick.label.cex.len <- length(tick.label.cex) + tick.label.cex <- tick.label.cex[ifelse(1:ax.num%%tick.label.cex.len==0,tick.label.cex.len,1:ax.num%%tick.label.cex.len)] + if(is.null(tick.label.cex)){tick.label.cex <- rep(0, ax.num)} + + tick.label.side.len <- length(tick.label.side) + tick.label.side <- tick.label.side[ifelse(1:ax.num%%tick.label.side.len==0,tick.label.side.len,1:ax.num%%tick.label.side.len)] + if(is.null(tick.label.side)){tick.label.side <- rep(0, ax.num)} + + predict.col.len <- length(predict.col) + predict.col <- predict.col[ifelse(1:ax.num%%predict.col.len==0,predict.col.len,1:ax.num%%predict.col.len)] + if(is.null(predict.col)){predict.col <- rep(0, ax.num)} + + predict.lwd.len <- length(predict.lwd) + predict.lwd <- predict.lwd[ifelse(1:ax.num%%predict.lwd.len==0,predict.lwd.len,1:ax.num%%predict.lwd.len)] + if(is.null(predict.lwd)){predict.lwd <- rep(0, ax.num)} + + predict.lty.len <- length(predict.lty) + predict.lty <- predict.lty[ifelse(1:ax.num%%predict.lty.len==0,predict.lty.len,1:ax.num%%predict.lty.len)] + if(is.null(predict.lty)){predict.lty <- rep(0, ax.num)} + + ax.names.len <- length(ax.names) + ax.names <- ax.names[ifelse(1:p%%ax.names.len==0,ax.names.len,1:p%%ax.names.len)] + if(is.null(ax.names)){ax.names <- rep("", p)} + + orthogx.len <- length(orthogx) + orthogx <- orthogx[ifelse(1:p%%orthogx.len==0,orthogx.len,1:p%%orthogx.len)] + if(is.null(orthogx)){orthogx <- rep("", p)} + + orthogy.len <- length(orthogy) + orthogy <- orthogx[ifelse(1:p%%orthogy.len==0,orthogy.len,1:p%%orthogy.len)] + if(is.null(orthogy)){orthogy <- rep("", p)} + + bp$ord.axes = list(which = which, col = col, lwd.factor = lwd.factor, reverse = reverse, + lwd = lwd, lty = lty, label.dir = label.dir, + label.col = label.col, label.cex = label.cex, label.line = label.line, + ticks = ticks, tick.col = tick.col, tick.size = tick.size, + tick.label = tick.label, tick.label.col = tick.label.col, + tick.label.cex = tick.label.cex, tick.label.side=tick.label.side, + predict.col = predict.col, predict.lty = predict.lty, + predict.lwd = predict.lwd, names = ax.names, orthogx = orthogx, + orthogy = orthogy) + + } + else bp$ord.axes <- list(which = NULL) + bp +} + # ---------------------------------------------------------------------------------------------- #' Aesthetics for alpha-bags #' @@ -771,12 +1182,47 @@ newsamples <- function (bp, col = "darkorange1", pch = 1, cex = 1, any(label.side!="bottom") | any(label.offset !=0.5) | any(label.cex!=0.75)) label<-TRUE nn <- nrow(bp$Xnew) - while (length(col) < nn) col <- c(col, col) - col <- as.vector(col[1:nn]) - while (length(pch) < nn) pch <- c(pch, pch) - pch <- as.vector(pch[1:nn]) - while (length(cex) < nn) cex <- c(cex, cex) - cex <- as.vector(cex[1:nn]) + gg <- bp$new.g + + if (length(col)==1) + col <- rep(col, nn) + else + { + col.len <- length(col) + col <- col[ifelse(1:gg%%col.len==0,col.len,1:gg%%col.len)] + if(is.null(col)){col <- rep(NA, gg)} + new.col <- rep(NA, nn) + for (j in 1:bp$new.g) + new.col[bp$new.group.aes == bp$new.g.names[j]] <- col[j] + col <- new.col + } + + if (length(pch)==1) + pch <- rep(pch, nn) + else + { + pch.len <- length(pch) + pch <- pch[ifelse(1:gg%%pch.len==0,pch.len,1:gg%%pch.len)] + if(is.null(pch)){pch <- rep(NA, gg)} + new.pch <- rep(NA, nn) + for (j in 1:bp$new.g) + new.pch[bp$new.group.aes == bp$new.g.names[j]] <- pch[j] + pch <- new.pch + } + + if (length(cex)==1) + cex <- rep(cex, nn) + else + { + cex.len <- length(cex) + cex <- cex[ifelse(1:gg%%cex.len==0,cex.len,1:gg%%cex.len)] + if(is.null(cex)){cex <- rep(NA, gg)} + new.cex <- rep(NA, nn) + for (j in 1:bp$new.g) + new.cex[bp$new.group.aes == bp$new.g.names[j]] <- cex[j] + cex <- new.cex + } + if (label[1] == "ggrepel") { label <- label[1] diff --git a/R/biplot.R b/R/biplot.R index 158f93e..6f8c829 100644 --- a/R/biplot.R +++ b/R/biplot.R @@ -85,8 +85,10 @@ #' * [rotate()] #' #' @references -#' Gabriel, K.R. (1971) The biplot graphic display of matrices with application to principal component analysis. \emph{Biometrika.} 58(3):453–467.

-#' Gower, J., Gardner-Lubbe, S. & Le Roux, N. (2011, ISBN: 978-0-470-01255-0) \emph{Understanding Biplots.} Chichester, England: John Wiley & Sons Ltd.

+#' Gabriel, K.R. (1971) The biplot graphic display of matrices with application to principal component analysis. \emph{Biometrika.} 58(3):453–467. +#' +#' Gower, J., Gardner-Lubbe, S. & Le Roux, N. (2011, ISBN: 978-0-470-01255-0) \emph{Understanding Biplots.} Chichester, England: John Wiley & Sons Ltd. +#' #' Gower, J.C. & Hand, D.J.(1996, ISBN: 0-412-71630-5) \emph{Biplots.} London: Chapman & Hall. #' #' @usage biplot(data, classes = NULL, group.aes = NULL, center = TRUE, scaled = FALSE, @@ -740,12 +742,14 @@ summary.biplot <- function (object, adequacy = TRUE, axis.predictivity = TRUE, #' containing supplementary data points to be added onto the biplot. #' @param newvariable a new data set, similar in structure to the data set supplied to \code{biplot()} #' containing supplementary variables to be added onto the biplot. +#' @param new.group.aes aesthetics for the new samples #' #' @return The object of class \code{biplot} will be appended with the following elements: #' \item{Xnew.raw}{the new data.} #' \item{Xnew}{the matrix of the centered and scaled new numeric variables of new data.} #' \item{Xnew.cat}{the matrix of the categorical variables of new data.} #' \item{Znew}{the matrix of the coordinates of the new data in the biplot.} +#' \item{new.group.aes}{aeshetics for the new samples} #' #' For an object of class \code{CA} the following additional elements will be appended: #' \item{newrowcoor}{the matrix of row coordinates of the new data in the biplot.} @@ -758,7 +762,8 @@ summary.biplot <- function (object, adequacy = TRUE, axis.predictivity = TRUE, #' biplot(HairEyeColor[,,2], center = FALSE) |> CA(variant = "Symmetric") |> #' interpolate(newdata = HairEyeColor[,,1]) |> plot() #' -interpolate <- function (bp, newdata=NULL, newvariable=NULL) +interpolate <- function (bp, newdata=NULL, newvariable=NULL, + new.group.aes = NULL) { # New samples if(!is.null(newdata)) @@ -811,7 +816,19 @@ interpolate <- function (bp, newdata=NULL, newvariable=NULL) Znew <- Ynew[,bp$e.vects] } - + if (!is.null(new.group.aes)) + { + bp$new.group.aes <- factor(new.group.aes) + bp$new.g.names <-levels(factor(new.group.aes)) + bp$new.g <- length(bp$new.g.names) + } + else + { + bp$new.group.aes <- factor(rep(1, nrow(Xnew))) + bp$new.g.names <-levels(bp$new.group.aes) + bp$new.g <- length(bp$new.g.names) + } + bp$Xnew.raw <- Xnew.raw bp$Xnew <- Xnew bp$Xnew.cat <- Xnew.cat @@ -1013,6 +1030,9 @@ prediction <- function (bp, predict.samples=NULL,predict.means=NULL,which=1:bp$p #' classify(col=c("red","blue","orange"),opacity=0.1) |> plot() classify <- function(bp, classify.regions = TRUE, col=ez.col, opacity=0.4, borders = FALSE) { + if (!requireNamespace("caret", quietly = TRUE)) { + stop("Package 'caret' is required for this function. Please install it.", call. = FALSE) + } if (inherits(bp, "CVA")) { while (length(col) < bp$g) col <- c(col, col) diff --git a/R/biplotEZ-package.R b/R/biplotEZ-package.R index d95718a..0b92bb4 100644 --- a/R/biplotEZ-package.R +++ b/R/biplotEZ-package.R @@ -54,5 +54,7 @@ #' @importFrom graphics points #' @importFrom graphics text #' @importFrom grDevices colorRampPalette +#' @importFrom Rcpp sourceCpp +#' @useDynLib biplotEZ, .registration = TRUE ## usethis namespace: end NULL diff --git a/R/calibrate_axes.R b/R/calibrate_axes.R deleted file mode 100644 index 5aa838f..0000000 --- a/R/calibrate_axes.R +++ /dev/null @@ -1,59 +0,0 @@ -#' Calibrate Biplot Axes -#' -#' Convenience function to obtain the coordinates of the calibrated ticks marks on the biplot -#' -#' @param x an object of class \code{biplot} -#' -#' @return An ordered list containing the coordinates the of tick marks to plotted on the biplot -#' @export -#' -#' @examples -#' x<-biplot(iris) |> PCA() -#' coordinates<-axes_coordinates(x) -axes_coordinates<-function(x){ - Z <- x$Z - if(!inherits(x,"CA")) { - if (is.null(x$axes)) x <- axes(x)} - ax.aes <- x$axes - if (length(ax.aes$which) > 0) - { - if (!is.null(x$Lmat)) - if (nrow(x$Lmat) == ncol(x$Lmat)) - Xhat <- x$Z %*% solve(x$Lmat)[x$e.vects,] - else Xhat <- x$X - else - Xhat <- x$X - if (x$scaled) Xhat <- scale(Xhat, center=FALSE, scale=1/x$sd) - if (x$center) Xhat <- scale(Xhat, center=-1*x$means, scale=FALSE) - - if(!is.null(x$PCOaxes)) - { if (x$PCOaxes == "splines") # Only for PCO - if axes (type) is set to splines. - { - z.axes <- lapply(1:length(ax.aes$which), biplot.spline.axis, Z, x$raw.X, - means=x$means, sd=x$sd, n.int=ax.aes$ticks, - spline.control=x$spline.control) - - - } - else if(x$PCOaxes == "regression") # Only for PCO - if axes (type) is set to regression. - { - z.axes <- lapply(1:length(ax.aes$which), .calibrate.axis, Xhat, x$means, x$sd, x$ax.one.unit, ax.aes$which, - ax.aes$ticks, ax.aes$orthogx, ax.aes$orthogy) - - for(i in 1:length(z.axes)){ - z.axes[[i]]<-z.axes[[i]][[1]] - } - } - } - else - { # Otherwise calibrate linear axes - z.axes <- lapply(1:length(ax.aes$which), .calibrate.axis, Xhat, x$means, x$sd, x$ax.one.unit, ax.aes$which, - ax.aes$ticks, ax.aes$orthogx, ax.aes$orthogy) - for(i in 1:length(z.axes)){ - z.axes[[i]]<-z.axes[[i]][[1]] - } - } - } - - return(z.axes) -} \ No newline at end of file diff --git a/R/clouds.R b/R/clouds.R index 70a5873..af25f24 100644 --- a/R/clouds.R +++ b/R/clouds.R @@ -153,7 +153,7 @@ density1D <- function(bp,which = NULL, h = "nrd0", kernel="gaussian", #' \item{lwd}{vector of line widths for the \eqn{\alpha}-bags.} #' #' @references -#' Gower, J., Gardner-Lubbe, S. & Le Roux, N. (2011, ISBN: 978-0-470-01255-0) \emph{Understanding Biplots.} Chichester, England: John Wiley & Sons Ltd.

+#' Gower, J., Gardner-Lubbe, S. & Le Roux, N. (2011, ISBN: 978-0-470-01255-0) \emph{Understanding Biplots.} Chichester, England: John Wiley & Sons Ltd. #' #' @export #' @usage alpha.bags(bp, alpha = 0.95, which = NULL, col = ez.col[which], lty = 1, @@ -168,6 +168,9 @@ alpha.bags <- function(bp, alpha=0.95, which = NULL, col = ez.col[which], lty = lwd = 1, max = 2500, trace = TRUE, opacity = 0.25, outlying=FALSE) { + if (!requireNamespace("geometry", quietly = TRUE)) { + stop("Package 'geometry' is required for this function. Please install it.", call. = FALSE) + } g <- bp$g g.names <- bp$g.names if (is.null(which)) which <- 1:g @@ -564,6 +567,14 @@ calc.alpha.bags <- function (x, y, aa=0.95, na.rm = TRUE, approx.limit = 2500, p } xydata <- if (missing(y)) x else cbind(x, y) + if (all(!is.na(apply(xydata,2,stats::sd)))) + if (any(apply(xydata,2,stats::sd)<1e-15)) + { + xy <- apply(xydata, 2, mean) + xy <- jitter(rbind (xy, xy, xy)) + return (list(xy=xy)) + } + if (is.data.frame(xydata)) xydata <- as.matrix(xydata) if (any(is.na(xydata))) { if (na.rm) { diff --git a/R/plot2D.R b/R/plot2D.R index 1842a16..9a7d28c 100644 --- a/R/plot2D.R +++ b/R/plot2D.R @@ -74,7 +74,6 @@ group.aes = group.aes, pch = rep(NA,n), col = rep(NA,n), cex = rep(NA,n), cex.vec, Z) - for(j in 1:length(sample.aes$which)) { ZZ$pch[group.aes==g.names[sample.aes$which[j]]] = sample.aes$pch[j] @@ -176,41 +175,45 @@ if(!any(invals)) return() ZZ <- Z[sample.aes$which,] - - toetsers<-invals[sample.aes$which] - if (sample.aes$label[1]=="ggrepel") + + if (nrow(ZZ) > 0) { - for (j in 1:nrow(ggrepel.labs$coords)){ - graphics::text(ggrepel.labs$coords[j, 1], ggrepel.labs$coords[j, 2], - labels = ggrepel.labs$coords[j,3], - cex = sample.aes$label.cex[ggrepel.labs$visible[j]], - col = sample.aes$label.col[ggrepel.labs$visible[j]]) + toetsers<-invals[sample.aes$which] + + if (sample.aes$label[1]=="ggrepel") + { + for (j in 1:nrow(ggrepel.labs$coords)){ + graphics::text(ggrepel.labs$coords[j, 1], ggrepel.labs$coords[j, 2], + labels = ggrepel.labs$coords[j,3], + cex = sample.aes$label.cex[ggrepel.labs$visible[j]], + col = sample.aes$label.col[ggrepel.labs$visible[j]]) + } + for (j in ggrepel.labs$textlines) + { + label.val <- rownames(ZZ)[j] + label.xy <- ggrepel.labs$coords[match(label.val, ggrepel.labs$coords[,3]),1:2] + graphics::lines (x=c(label.xy[1],ZZ[j,1]), y=c(label.xy[2],ZZ[j,2]), col=sample.aes$label.col[j]) + } } - for (j in ggrepel.labs$textlines) - { - label.val <- rownames(ZZ)[j] - label.xy <- ggrepel.labs$coords[match(label.val, ggrepel.labs$coords[,3]),1:2] - graphics::lines (x=c(label.xy[1],ZZ[j,1]), y=c(label.xy[2],ZZ[j,2]), col=sample.aes$label.col[j]) + else + { + Z.labels <- rownames(Z)[sample.aes$which] + for (j in 1:length(sample.aes$label.side)) + { if(!toetsers[j]) + next + text.pos <- match(sample.aes$label.side[j], c("bottom", "left", "top", "right")) + if (sample.aes$label[j]) graphics::text(ZZ[j, 1], ZZ[j, 2], labels = Z.labels[j], + cex = sample.aes$label.cex[j], col = sample.aes$label.col[j], + pos = text.pos, offset = sample.aes$label.offset[j]) + } } - } - else - { - Z.labels <- rownames(Z)[sample.aes$which] - for (j in 1:length(sample.aes$label.side)) - { if(!toetsers[j]) + for (j in 1:length(sample.aes$which)){ + if(!toetsers[j]) next - text.pos <- match(sample.aes$label.side[j], c("bottom", "left", "top", "right")) - if (sample.aes$label[j]) graphics::text(ZZ[j, 1], ZZ[j, 2], labels = Z.labels[j], - cex = sample.aes$label.cex[j], col = sample.aes$label.col[j], - pos = text.pos, offset = sample.aes$label.offset[j]) + graphics::points(x = Z[sample.aes$which[j], 1], y = Z[sample.aes$which[j], 2], + pch = sample.aes$pch[j], col = sample.aes$col[j], cex = sample.aes$cex[j]) } } - for (j in 1:length(sample.aes$which)){ - if(!toetsers[j]) - next - graphics::points(x = Z[sample.aes$which[j], 1], y = Z[sample.aes$which[j], 2], - pch = sample.aes$pch[j], col = sample.aes$col[j], cex = sample.aes$cex[j]) - } } #' Title @@ -352,6 +355,61 @@ + +#---------- +#' .calibrate.cat.axis +#' +#' @param j Index of axis to calibrate +#' @param axes.rows +#' @param ax.which +#' @param ax.orthogxvec +#' @param ax.orthogyvec +#' @param markers +#' +#' @noRd +.calibrate.cat.axis <- function (j, axes.rows, ax.which, ax.orthogxvec, ax.orthogyvec, + markers) +{ + ax.num <- ax.which[j] + ax.direction <- axes.rows[ax.num,] + r <- ncol(axes.rows) + ax.orthog <- rbind(ax.orthogxvec, ax.orthogyvec) + if (nrow(ax.orthog) < r) ax.orthog <- rbind(ax.orthog, 0) + if (nrow(axes.rows) > 1) phi.vec <- diag(1 / diag(axes.rows %*% t(axes.rows))) %*% axes.rows %*% ax.orthog[, ax.num] + else phi.vec <- (1 / (axes.rows %*% t(axes.rows))) %*% axes.rows %*% ax.orthog[, ax.num] + + std.range <- range(markers[[ax.num]]) + std.ax.tick.label.min <- min(markers[[ax.num]]) - (std.range[2] - std.range[1]) + std.ax.tick.label.max <- max(markers[[ax.num]]) + (std.range[2] - std.range[1]) + axis.vals <- c(std.ax.tick.label.min, markers[[ax.num]], std.ax.tick.label.max) + + number.points <- length(axis.vals) + axis.points <- matrix(0, nrow = number.points, ncol = r) + for (i in 1:r) + axis.points[, i] <- ax.orthog[i, ax.num] + (axis.vals - phi.vec[ax.num]) * ax.direction[i] + rownames(axis.points) <- c("", names(markers[[ax.num]]), "") + + #slope = delta y / delta x of two datapoints + slope <- (axis.points[1, 2] - axis.points[2, 2]) / (axis.points[1, 1] - axis.points[2, 1]) + #if slope is infinite then all x-values are same + v <- NULL + if (is.na(slope)) + { v <- axis.points[1, 1] + slope = NULL + } + else if (abs(slope) == Inf) + { v <- axis.points[1, 1] + slope = NULL + } + + #y=mx+c... c=y-mx + intercept <- axis.points[1, 2] - slope * axis.points[1, 1] + + retvals <- list(coords = axis.points, a = intercept, b = slope, v = v) + return(retvals) +} + + #' Insert tick marks and label - called by .marker.func #' #' @param x x-coo of tick @@ -366,14 +424,14 @@ #' @param offset REMOVE #' @param label.col color of label #' @param cex Bigness of label -#' @param label_pos 1 for above -1 for below +#' @param label.pos 1 for above -1 for below #' #' @return inserts tick marks and label to plot #' @noRd .marker.label <- function(x, y, grad, marker.val, expand = 1, col, label.on.off, label.col, cex, - label_pos=-1,usr) + label.pos=-1,usr) { #get fraction of plot width and expand @@ -382,14 +440,14 @@ if (grad == "v"){ graphics::lines(rep(x, 2), c(y - d, y + d), col = col) if (label.on.off == 1) - graphics::text(x, y - label_pos(d+2.5*mm), label = marker.val, + graphics::text(x, y - label.pos(d+2.5*mm), label = marker.val, col = label.col, cex = cex) } if (grad == "h"){ graphics::lines(c(x - d, x + d), rep(y, 2), col = col) if (label.on.off == 1) - graphics::text(x - label_pos*(d+2.5*mm), y, label = marker.val, + graphics::text(x - label.pos*(d+2.5*mm), y, label = marker.val, col = label.col, cex = cex, srt=90) } @@ -406,14 +464,14 @@ if (label.on.off == 1){ if(grad>0) - graphics::text(x+label_pos*(2.5*x_shift+b), - y+label_pos*(2.5*y_shift+a), + graphics::text(x+label.pos*(2.5*x_shift+b), + y+label.pos*(2.5*y_shift+a), label = marker.val, col = label.col, cex = cex,srt=180+(90+atan(grad)*180/pi)) if(grad<0) - graphics::text(x-label_pos*(2.5*x_shift+b), y-label_pos*(2.5*y_shift+a), + graphics::text(x-label.pos*(2.5*x_shift+b), y-label.pos*(2.5*y_shift+a), label = marker.val, col = label.col, cex = cex,srt=(90+atan(grad)*180/pi)) @@ -439,20 +497,17 @@ #' @return prints markers to plot #' @noRd .marker.func <- function(vec, coef, col, tick.size, - label_pos, label.col, cex, + label.pos, label.col, cex, usr) { - if(label_pos=="below") - label_pos<--1 - if(label_pos=="above") - label_pos<- 1 + if(label.pos=="below") label.pos<--1 + if(label.pos=="above") label.pos<- 1 - ## vec = [x_coo, y_coo, marker val, show] ## coef = [intercept, slope] - x <- vec[1] - y <- vec[2] + x <- as.numeric(vec[1]) + y <- as.numeric(vec[2]) marker.val <- vec[3] label.on.off <- vec[4] if (is.na(coef[2])) @@ -460,20 +515,20 @@ expand = tick.size, col = col, label.on.off = label.on.off, label.col = label.col, - cex = cex,usr=usr, label_pos=label_pos) + cex = cex,usr=usr, label.pos=label.pos) else if (coef[2] == 0) # line is verticle .marker.label(x, y, grad = "v", marker.val, expand = tick.size, col = col, label.on.off = label.on.off, label.col = label.col, - cex = cex, usr=usr,label_pos=label_pos) + cex = cex, usr=usr,label.pos=label.pos) else .marker.label(x, y, grad = -1/coef[2], marker.val, expand = tick.size, col = col, label.on.off = label.on.off, label.col = label.col, - cex = cex, usr=usr,label_pos=label_pos) + cex = cex, usr=usr,label.pos=label.pos) } @@ -503,60 +558,90 @@ marker.mat <- marker.mat[rev(order(marker.mat[, 3])), ] x.vals <- marker.mat[, 1] y.vals <- marker.mat[, 2] + k <- nrow(marker.mat) - lin.coef<-c(a=this.axis$a,b=this.axis$b) if (is.null(this.axis$b)) graphics::abline(v = this.axis$v, col = ax.aes$col[i], lwd = ax.aes$lwd[i], lty = ax.aes$lty[i]) else - graphics::abline(coef=lin.coef, col = ax.aes$col[i], lwd = ax.aes$lwd[i], lty = ax.aes$lty[i]) - if (ax.aes$label.dir == "Hor") { graphics::par(las = 1) - adjust <- c(0.5, 1, 0.5, 0) } - if (ax.aes$label.dir == "Orthog") { graphics::par(las = 2) - adjust <- c(1, 1, 0, 0) } - if (ax.aes$label.dir == "Paral") { graphics::par(las = 0) - adjust <- c(0.5, 0.5, 0.5, 0.5) } - h <- nrow(marker.mat) - if (is.null(this.axis$b)) - { if (y.vals[1] < y.vals[h]) - graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], adj = adjust[1], at = x.vals[1], col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) + graphics::abline(coef=c(this.axis$a, this.axis$b), col = ax.aes$col[i], lwd = ax.aes$lwd[i], lty = ax.aes$lty[i]) + if (ax.aes$label.dir == "Hor") + { graphics::par(las = 1) + adjust <- c(0.5, 1, 0.5, 0) + } + if (ax.aes$label.dir == "Orthog") + { graphics::par(las = 2) + adjust <- c(1, 1, 0, 0) + } + if (ax.aes$label.dir == "Paral") + { graphics::par(las = 0) + adjust <- c(0.5, 0.5, 0.5, 0.5) + } + + if (is.null(this.axis$b)) + { if (y.vals[1] < y.vals[k]) + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) else - graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], adj = adjust[3], at = y.vals[1], col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) } - else - { y1.ster <- lin.coef[2] * usr[1] + lin.coef[1] - y2.ster <- lin.coef[2] * usr[2] + lin.coef[1] - x1.ster <- (usr[3] - lin.coef[1])/lin.coef[2] - x2.ster <- (usr[4] - lin.coef[1])/lin.coef[2] - if (lin.coef[2] == 0) - { if (x.vals[1] < x.vals[h]) - graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], adj = adjust[2], at = y.vals[1], col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else - graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], adj = adjust[4], at = y.vals[1], col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - } - if (lin.coef[2] > 0) - { if (x.vals[1] < x.vals[h]) - if (y1.ster <= usr[4] & y1.ster >= usr[3]) - graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else - graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else if (y2.ster <= usr[4] & y2.ster >= usr[3]) - graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else - graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - } - if (lin.coef[2] < 0) - { if (x.vals[1] < x.vals[h]) - if (y1.ster <= usr[4] & y1.ster >= usr[3]) - graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else - graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else if (y2.ster <= usr[4] & y2.ster >= usr[3]) - graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) - else - graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], cex = ax.aes$label.cex[i]) + else + { y1.ster <- this.axis$b * usr[1] + this.axis$a + y2.ster <- this.axis$b * usr[2] + this.axis$a + x1.ster <- (usr[3] - this.axis$a)/this.axis$b + x2.ster <- (usr[4] - this.axis$a)/this.axis$b + if (this.axis$b == 0) + { if (x.vals[1] < x.vals[k]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + if (this.axis$b > 0) + { if (x.vals[1] < x.vals[k]) + if (y1.ster <= usr[4] & y1.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else if (y2.ster <= usr[4] & y2.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) } - } - + if (this.axis$b < 0) + { if (x.vals[1] < x.vals[k]) + if (y1.ster <= usr[4] & y1.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else if (y2.ster <= usr[4] & y2.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + } + invals <- x.vals < usr[2] & x.vals > usr[1] & y.vals < usr[4] & y.vals > usr[3] std.markers <- marker.mat[invals, 3] if (is.numeric(std.markers)) std.markers <- zapsmall(std.markers) @@ -565,18 +650,308 @@ if (ax.aes$tick.label[i]) label.on.off <- rep(1, sum(invals)) else label.on.off <- rep(0, sum(invals)) if (!ax.aes$tick.label[i]) label.on.off[c(1, length(label.on.off))] <- 1 if(sum(invals)>0) apply(data.frame(x.vals, y.vals, std.markers, label.on.off), 1, .marker.func, - coef = lin.coef, col = ax.aes$tick.col[i], tick.size = ax.aes$tick.size[i], - label.col = ax.aes$tick.label.col[i], label_pos=ax.aes$tick.label.side[i], + coef = c(this.axis$a, this.axis$b), col = ax.aes$tick.col[i], tick.size = ax.aes$tick.size[i], + label.col = ax.aes$tick.label.col[i], label.pos=ax.aes$tick.label.side[i], cex = ax.aes$tick.label.cex[i], usr=usr) if(ax.num %in% predict_which) { - if (!is.null(predict.mat)) apply(cbind(predict.mat,y.vals[1]), 1, .predict.func, coef = lin.coef,col=ax.aes$predict.col,lwd=ax.aes$predict.lwd,lty=ax.aes$predict.lty) + if (!is.null(predict.mat)) apply(cbind(predict.mat,y.vals[1]), 1, .predict.func, coef = c(this.axis$a, this.axis$b), + col=ax.aes$predict.col,lwd=ax.aes$predict.lwd,lty=ax.aes$predict.lty) + } + } +} + +#' Plot linear nominal axes on biplots +#' +#' @param z.axes list containing all the info to draw axis. see below +#' @param ax.aes Axis aestetics +#' @param predict.mat Xhat +#' @param too.small cutoff: predictivity smaller than cutoff not plotted +#' @param usr plot dim +#' +#' z.axes: List containing following components for each axis: +#' $coords: Tick mark coordinates +#' $a $b $v: Intercept, slope, vertical +#' +#' @noRd +.nom.axes.plot <- function(z.axes, ax.aes, predict.mat, too.small, usr, predict_which) +{ + for (i in 1:length(ax.aes$which)) + { ax.num <- ax.aes$which[i] + if (!is.null(too.small)) if (ax.num %in% too.small) next + this.axis <- z.axes[[i]] + marker.mat <- this.axis$coords + marker.mat <- marker.mat[order(marker.mat[, 1]), ] + + ends.mat <- c(usr[1], NA) + k <- nrow(marker.mat) + for (j in 2:nrow(marker.mat)) + ends.mat <- rbind (ends.mat, (marker.mat[j,1:2] + marker.mat[j-1,1:2])/2) + ends.mat <- rbind (ends.mat, c(usr[2],NA)) + x.vals <- marker.mat[, 1] + y.vals <- marker.mat[, 2] + + if (ends.mat[1,1] < ends.mat[2,1]) ends.mat[2,] <- ends.mat[1,] + ends.mat <- ends.mat[-1,] + if (ends.mat[k,1] > ends.mat[k-1,1]) ends.mat[k-1,] <- ends.mat[k,] + ends.mat <- ends.mat[-k,] + + if (is.null(this.axis$b)) + { + ends.mat[1,1] <- ends.mat[nrow(ends.mat),1] <- 0 + ends.mat[1,2] <- usr[3] + ends.mat[nrow(ends.mat),2] <- usr[4] + # graphics::abline(v = this.axis$v, col = ax.aes$col[i], lwd = ax.aes$lwd[i], lty = ax.aes$lty[i]) + } + else + { + ends.mat[1,2] <- this.axis$a + ends.mat[1,1]*this.axis$b + ends.mat[nrow(ends.mat),2] <- this.axis$a + ends.mat[nrow(ends.mat),1]*this.axis$b + } + + if (ax.aes$label.dir == "Hor") + { graphics::par(las = 1) + adjust <- c(0.5, 1, 0.5, 0) + } + if (ax.aes$label.dir == "Orthog") + { graphics::par(las = 2) + adjust <- c(1, 1, 0, 0) + } + if (ax.aes$label.dir == "Paral") + { graphics::par(las = 0) + adjust <- c(0.5, 0.5, 0.5, 0.5) + } + + if (is.null(this.axis$b)) + { if (y.vals[1] < y.vals[k]) + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + else + { y1.ster <- this.axis$b * usr[1] + this.axis$a + y2.ster <- this.axis$b * usr[2] + this.axis$a + x1.ster <- (usr[3] - this.axis$a)/this.axis$b + x2.ster <- (usr[4] - this.axis$a)/this.axis$b + if (this.axis$b == 0) + { if (x.vals[1] < x.vals[k]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + if (this.axis$b > 0) + { if (x.vals[1] < x.vals[k]) + if (y1.ster <= usr[4] & y1.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else if (y2.ster <= usr[4] & y2.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + if (this.axis$b < 0) + { if (x.vals[1] < x.vals[k]) + if (y1.ster <= usr[4] & y1.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else if (y2.ster <= usr[4] & y2.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + } + for (j in 2:nrow(ends.mat)) + { + graphics::lines (x = ends.mat[(-1:0)+j, 1], y = ends.mat[(-1:0)+j, 2], + col = ax.aes$col[[i]][j-1], lwd = 2*ax.aes$lwd[i], lty=ax.aes$lty[i]) + .marker.func (vec = c(marker.mat[j,], rownames(marker.mat)[j], 1), + coef = c(this.axis$a, this.axis$b), col = ax.aes$tick.col, + tick.size = 0, label.pos = ax.aes$tick.label.side[i], label.col = ax.aes$tick.label.col, + cex = ax.aes$tick.label.cex, usr = usr) + } + + if(ax.num %in% predict_which) + { + if (!is.null(predict.mat)) apply(cbind(predict.mat,y.vals[1]), 1, .predict.func, + coef = c(this.axis$a, this.axis$b), col=ax.aes$predict.col, + lwd=ax.aes$predict.lwd,lty=ax.aes$predict.lty) } } } +#' Plot linear ordinal axes on biplots +#' +#' @param z.axes list containing all the info to draw axis. see below +#' @param ax.aes Axis aestetics +#' @param predict.mat Xhat +#' @param too.small cutoff: predictivity smaller than cutoff not plotted +#' @param usr plot dim +#' +#' z.axes: List containing following components for each axis: +#' $coords: Tick mark coordinates +#' $a $b $v: Intercept, slope, vertical +#' +#' @noRd +.ord.axes.plot <- function(z.axes, ax.aes, predict.mat, too.small, usr, predict_which) +{ + for (i in 1:length(ax.aes$which)) + { ax.num <- ax.aes$which[i] + if (!is.null(too.small)) if (ax.num %in% too.small) next + this.axis <- z.axes[[i]] + marker.mat <- this.axis$coords + marker.mat <- marker.mat[order(marker.mat[, 1]), ] + + ends.mat <- c(usr[1], NA) + k <- nrow(marker.mat) + for (j in 2:nrow(marker.mat)) + ends.mat <- rbind (ends.mat, (marker.mat[j,1:2] + marker.mat[j-1,1:2])/2) + ends.mat <- rbind (ends.mat, c(usr[2],NA)) + x.vals <- marker.mat[, 1] + y.vals <- marker.mat[, 2] + + if (ends.mat[1,1] < ends.mat[2,1]) ends.mat[2,] <- ends.mat[1,] + ends.mat <- ends.mat[-1,] + if (ends.mat[k,1] > ends.mat[k-1,1]) ends.mat[k-1,] <- ends.mat[k,] + ends.mat <- ends.mat[-k,] + + if (is.null(this.axis$b)) + { + ends.mat[1,1] <- ends.mat[nrow(ends.mat),1] <- 0 + ends.mat[1,2] <- usr[3] + ends.mat[nrow(ends.mat),2] <- usr[4] + } + else + { + ends.mat[1,2] <- this.axis$a + ends.mat[1,1]*this.axis$b + ends.mat[nrow(ends.mat),2] <- this.axis$a + ends.mat[nrow(ends.mat),1]*this.axis$b + } + + if (ax.aes$label.dir == "Hor") + { graphics::par(las = 1) + adjust <- c(0.5, 1, 0.5, 0) + } + if (ax.aes$label.dir == "Orthog") + { graphics::par(las = 2) + adjust <- c(1, 1, 0, 0) + } + if (ax.aes$label.dir == "Paral") + { graphics::par(las = 0) + adjust <- c(0.5, 0.5, 0.5, 0.5) + } + + if (is.null(this.axis$b)) + { if (y.vals[1] < y.vals[k]) + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + else + { y1.ster <- this.axis$b * usr[1] + this.axis$a + y2.ster <- this.axis$b * usr[2] + this.axis$a + x1.ster <- (usr[3] - this.axis$a)/this.axis$b + x2.ster <- (usr[4] - this.axis$a)/this.axis$b + if (this.axis$b == 0) + { if (x.vals[1] < x.vals[k]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y.vals[1], col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + if (this.axis$b > 0) + { if (x.vals[1] < x.vals[k]) + if (y1.ster <= usr[4] & y1.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else if (y2.ster <= usr[4] & y2.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + if (this.axis$b < 0) + { if (x.vals[1] < x.vals[k]) + if (y1.ster <= usr[4] & y1.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 2, line = ax.aes$label.line[i], + adj = adjust[2], at = y1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 3, line = ax.aes$label.line[i], + adj = adjust[3], at = x2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else if (y2.ster <= usr[4] & y2.ster >= usr[3]) + graphics::mtext(text = ax.aes$names[i], side = 4, line = ax.aes$label.line[i], + adj = adjust[4], at = y2.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + else + graphics::mtext(text = ax.aes$names[i], side = 1, line = ax.aes$label.line[i], + adj = adjust[1], at = x1.ster, col = ax.aes$label.col[i], + cex = ax.aes$label.cex[i]) + } + } + if (ax.aes$reverse[i]) lwd.vec <- (nrow(ends.mat)-1):1 else lwd.vec <- 1:(nrow(ends.mat)-1) + lwd.vec <- lwd.vec * ax.aes$lwd.factor[i] + for (j in 2:nrow(ends.mat)) + { + graphics::lines (x = ends.mat[(-1:0)+j, 1], y = ends.mat[(-1:0)+j, 2], + col = ax.aes$col[i], lwd = lwd.vec[j-1], lty=ax.aes$lty[i]) + .marker.func (vec = c(marker.mat[j,], rownames(marker.mat)[j], 1), + coef = c(this.axis$a, this.axis$b), col = ax.aes$tick.col, + tick.size = 0, label.pos = ax.aes$tick.label.side[i], label.col = ax.aes$tick.label.col, + cex = ax.aes$tick.label.cex, usr = usr) + } + + if(ax.num %in% predict_which) + { + if (!is.null(predict.mat)) apply(cbind(predict.mat,y.vals[1]), 1, .predict.func, + coef = c(this.axis$a, this.axis$b), col=ax.aes$predict.col, + lwd=ax.aes$predict.lwd,lty=ax.aes$predict.lty) + } + } +} + .CLPs.plot <- function (coords, colvec, cexvec) { text (x=coords[,1], y=coords[,2], labels=rownames(coords), col=colvec, cex=cexvec) @@ -638,6 +1013,15 @@ #' @noRd .get.ggrepel.coords <- function(df) { + if (!requireNamespace("ggplot2", quietly = TRUE)) { + stop("Package 'ggplot2' is required for this function. Please install it.", call. = FALSE) + } + if (!requireNamespace("ggrepel", quietly = TRUE)) { + stop("Package 'ggrepel' is required for this function. Please install it.", call. = FALSE) + } + if (!requireNamespace("grid", quietly = TRUE)) { + stop("Package 'grid' is required for this function. Please install it.", call. = FALSE) + } pp <- ggplot2::ggplot (df, ggplot2::aes(df$x,df$y,label=df$z)) + ggplot2::geom_point() + ggrepel::geom_text_repel() @@ -686,12 +1070,11 @@ } - #' calibrate spline based axes on biplots #' #' @param j Index of the axis to be calibrated in the data #' @param X Coordinates of the samples on the biplot space -#' @param Ytilde Raw data used to construct the biplot +#' @param Y Raw data used to construct the biplot #' @param means column means of the raw data #' @param sd Column standard deviations of the data #' @param n.int the number of tick marks per axis @@ -699,16 +1082,15 @@ #' @param dmeth Argument unused #' @param ... additional arguments #' -#' @useDynLib biplotEZ, .registration = TRUE -#' #' @noRd + + + biplot.spline.axis <- function(j, X, Y, means, sd, - n.int, spline.control, dmeth=0, ... ) -{ + n.int, spline.control, dmeth=0, ... ){ n <- nrow(X) p <- ncol(X) - if (n > 103) - { + if (n > 103){ my.sample <- sample (1:n, size=103, replace=F) X <- X[my.sample,] Y <- Y[my.sample,] @@ -743,8 +1125,7 @@ biplot.spline.axis <- function(j, X, Y, means, sd, mu <- unique(mu) nmu <- length(mu) - if (v>0) - { + if (v>0){ knots <- seq.int(from=0,to=1,length.out=v+2)[-c(1,v+2)] knots <- stats::quantile(y,knots) M <- splines::bs(mu,knots=knots,degree=u,intercept=FALSE) @@ -757,56 +1138,44 @@ biplot.spline.axis <- function(j, X, Y, means, sd, const1 <- sum(y^2) const2 <- sum(X^2)/(n*p) TotalNumberOfLossFunctionCalls <- 0 - - optimtouse <- function(Bvec) - { + optimtouse <- function(Bvec) { timetemp <- proc.time()[3] - LOSS <- 1.0 - LOSS1 <- 1.0 - Ind <- rep(1,n) - pred <- rep(0,nmu) - deltmp <- 0 - tau <- tau - #.5 # the choice of tau seems to affect perfomance quite substantially. - # tau is used to specify the points on the inital simplex. - Ay <- rep(0,(u+v)*p+1) - TEMPVK <- rep(0,(u+v)*p) - iter1 <- 0 - iter <- 0 - ERRO <- 0 - - # Prepare for Fortran subroutine - storage.mode(X) <- "double" - storage.mode(Ind) <- "integer" - storage.mode(mu) <- "double" - storage.mode(pred) <- "double" - storage.mode(y) <- "double" - storage.mode(M) <- "double" - storage.mode(Bvec) <- "double" - storage.mode(Ay) <- "double" - storage.mode(TEMPVK) <- "double" + # Use Rcpp optimize_spline function + returned_data <- optimize_spline( + BVEC = Bvec, + X = X, + Y = y, + M = M, + MU = mu, + LAMBDA = lambda, + CONST1 = const1, + CONST2 = const2, + U = u, + V = v, + TAU = tau, + FTOL = ftol, + TINY = tiny, + ITMAX = itmax + ) - returned_data <-.Fortran('L',LOSS=as.double(LOSS),X=X,n=as.integer(n),p=as.integer(p),nmu=as.integer(nmu),Ind=Ind, - mu=mu,pred=pred,lambda=as.double(lambda),y=y,const1=as.double(const1),const2=as.double(const2),u=as.integer(u), - v=as.integer(v),M=M,Bvec=Bvec,tau=as.double(tau),Ay=Ay,TEMPVEK=TEMPVK,iter=as.integer(iter), - ftol=as.double(ftol),LOSS1=as.double(LOSS1),iter1=as.integer(iter1),fout = as.integer(ERRO), - const3=as.double(tiny), itmax=as.integer(itmax)) - if(returned_data$fout > 0) - { - cat("Fout is: ", returned_data$fout, "\n") - warning("Increase itmax for Fortran \n") + if(returned_data$ERRO > 0) { + cat(" Warning: Error code =", returned_data$ERRO, "\n") } - B <- matrix(returned_data$Bvec,ncol=p) - Z <- M%*%B + aa <- list( + BestValue = returned_data$LOSS, + BestSolution = returned_data$BVEC, + ConvergenceCode = returned_data$ERRO, + iter1 = returned_data$ITER1, + iter = returned_data$ITER, + TimeTaken = proc.time()[3] - timetemp + ) - aa <- list(BestValue=returned_data$LOSS,BestSolution=returned_data$Bvec,ConvergenceCode=returned_data$fout, iter1=returned_data$iter1, - iter=returned_data$iter,TimeTaken=proc.time()[3]-timetemp) aa } - EuclidDist2 <- function (X, Y) - { + + EuclidDist2 <- function (X, Y) { n <- nrow(X) m <- nrow(Y) bx <- rowSums(X^2) @@ -838,29 +1207,24 @@ biplot.spline.axis <- function(j, X, Y, means, sd, test.iter1 <- temp$iter1 ### Last best coefficients perturbed - for (gammacounter in 2:(gamma+1)) - { - if (nSameSolutionConsecutively>=bigsigmaactivate) - { + for (gammacounter in 2:(gamma+1)){ + if (nSameSolutionConsecutively>=bigsigmaactivate){ temp <- optimtouse(outBestSolutions[,which.min(outBestValues)]+stats::rnorm((u+v)*2,mean=0,sd=bigsigma)) BigSigmaActivations <- c(BigSigmaActivations,gammacounter) } else temp <- optimtouse(outBestSolutions[,which.min(outBestValues)]+stats::rnorm((u+v)*2,mean=0,sd=smallsigma)) outTimeTaken[gammacounter] <- temp$TimeTaken tempSquaredDistances <- EuclidDist2(matrix(temp$BestSolution,nrow=1),t(outBestSolutions[,1:DistinctSolutions])) - if (any(tempSquaredDistances PCA() |> plot() plot.biplot <- function(x, exp.factor=1.2, axis.predictivity=NULL, sample.predictivity=NULL, - zoom=FALSE, xlim = NULL, ylim = NULL, ...) + zoom = FALSE, add = FALSE, xlim = NULL, ylim = NULL, ...) { + #----- See all the internal functions in utility_2D.R if (is.null(x$Z)) stop ("Add a biplot method before generating a plot") else Z <- x$Z #aesthetics for samples if (is.null(x$samples)) x <- samples(x) - + + if (add) zoom <- FALSE if(zoom) grDevices::dev.new() @@ -63,6 +66,9 @@ plot.biplot <- function(x, exp.factor=1.2, axis.predictivity=NULL, sample.predic ggrepel.new <- ggrepel.means <- ggrepel.samples <- NULL if (do.ggrepel) { + if (!requireNamespace("R.devices", quietly = TRUE)) { + stop("Package 'R.devices' is required for this function. Please install it.", call. = FALSE) + } out <- R.devices::suppressGraphics(.get.ggrepel.coords(df)) if (n.newsamples>n.means) ggrepel.new <-list(coords = out$coords[out$visible>n.newsamples & out$visiblen.newsamples & out$visible 0) { - if (!is.null(x$Lmat)) - if (nrow(x$Lmat) == ncol(x$Lmat)) - Xhat <- x$Z %*% solve(x$Lmat)[x$e.vects,] - else Xhat <- x$X - else - Xhat <- x$X - if (x$scaled) Xhat <- scale(Xhat, center=FALSE, scale=1/x$sd) - if (x$center) Xhat <- scale(Xhat, center=-1*x$means, scale=FALSE) - if(!is.null(x$PCOaxes)) { if (x$PCOaxes == "splines") # Only for PCO - if axes (type) is set to splines. { @@ -201,15 +198,17 @@ plot.biplot <- function(x, exp.factor=1.2, axis.predictivity=NULL, sample.predic } else if(x$PCOaxes == "regression") # Only for PCO - if axes (type) is set to regression. { - z.axes <- lapply(1:length(ax.aes$which), .calibrate.axis, Xhat, x$means, x$sd, x$ax.one.unit, ax.aes$which, - ax.aes$ticks, ax.aes$orthogx, ax.aes$orthogy) + z.axes <- lapply(1:length(ax.aes$which), .calibrate.axis, x$X, x$means, + x$sd, x$ax.one.unit, ax.aes$which, ax.aes$ticks, + ax.aes$orthogx, ax.aes$orthogy) .lin.axes.plot(z.axes, ax.aes, predict.mat, too.small,usr=usr,predict_which=x$predict$which) } } else { # Otherwise calibrate linear axes - z.axes <- lapply(1:length(ax.aes$which), .calibrate.axis, Xhat, x$means, x$sd, x$ax.one.unit, ax.aes$which, - ax.aes$ticks, ax.aes$orthogx, ax.aes$orthogy) + z.axes <- lapply(1:length(ax.aes$which), .calibrate.axis, x$X, x$means, x$sd, + x$ax.one.unit, ax.aes$which, ax.aes$ticks, + ax.aes$orthogx, ax.aes$orthogy) .lin.axes.plot(z.axes, ax.aes, predict.mat, too.small,usr=usr,predict_which=x$predict$which) } } @@ -226,13 +225,56 @@ plot.biplot <- function(x, exp.factor=1.2, axis.predictivity=NULL, sample.predic if (length(new.ax.aes$which) > 0) { - z.axes.new <- lapply(1:length(new.ax.aes$which), .calibrate.axis, - x$newvariable, x$new.means, x$new.sd, x$new.ax.one.unit, new.ax.aes$which, + z.axes.new <- lapply(1:length(new.ax.aes$which), .calibrate.axis, x$newvariable, + x$new.means, x$new.sd, x$new.ax.one.unit, new.ax.aes$which, new.ax.aes$ticks, new.ax.aes$orthogx, new.ax.aes$orthogy) - .lin.axes.plot(z.axes.new, new.ax.aes, predict.mat, too.small, usr=usr, predict_which=x$predict$which) + .lin.axes.plot(z.axes.new, new.ax.aes, predict.mat, too.small, usr=usr, + predict_which=x$predict$which) } } - + + # Nominal and ordinal axes + if(inherits(x,"catPCA")) + { + if (sum(x$ax.type == "nominal") > 0) + { + nom.one.unit <- x$all.ax.one.unit[x$ax.type$ax.type == "nominal",,drop = FALSE] + if (is.null(x$nom.axes)) x <- nom.axes(x) + + z.axes <- lapply(1:length(x$nom.axes$which), .calibrate.cat.axis, + nom.one.unit, x$nom.axes$which, x$nom.axes$orthogx, + x$nom.axes$orthogy, x$nom.levels) + .nom.axes.plot(z.axes, x$nom.axes, predict.mat, too.small, + usr = usr, predict_which = x$predict$which) + } + + if (sum(x$ax.type == "ordinal") > 0) + { + ord.one.unit <- x$all.ax.one.unit[x$ax.type$ax.type == "ordinal",,drop = FALSE] + if (is.null(x$ord.axes)) x <- ord.axes(x) + for (j in 1:length(x$ord.levels)) + { + current.markers <- x$ord.levels[[j]] + i <- 1 + while (i < length(current.markers) - 1) + { + if (abs(current.markers[i]-current.markers[i+1]) < 1.5*.Machine$double.eps) + { + names(current.markers[i]) <- paste(names(current.markers[i]), + names(current.markers[i+1]), sep="*") + current.markers <- current.markers[-(i+1)] + } + else i <- i + 1 + } + x$ord.levels[[j]] <- current.markers + } + z.axes <- lapply(1:length(x$ord.axes$which), .calibrate.cat.axis, + ord.one.unit, x$ord.axes$which, x$ord.axes$orthogx, + x$ord.axes$orthogy, x$ord.levels) + .ord.axes.plot(z.axes, x$ord.axes, predict.mat, too.small, + usr = usr,predict_which = x$predict$which) + } + } # Fit measures too.small <- NULL @@ -329,6 +371,9 @@ plot.biplot <- function(x, exp.factor=1.2, axis.predictivity=NULL, sample.predic plot3D <- function(bp, exp.factor = 1.2,...) { + if (!requireNamespace("rgl", quietly = TRUE)) { + stop("Package 'rgl' is required for this function. Please install it.", call. = FALSE) + } if (is.null(bp$Z)) stop ("Add a biplot method before generating a plot") else Z <- bp$Z @@ -408,7 +453,9 @@ if(!inherits(bp,"CA")){ } - + # Nominal and ordinal axes + if(inherits(bp,"catPCA")) stop ("catPCA not implemented for 3D biplots") + # Bags @@ -520,6 +567,9 @@ plot1D <- function(bp, exp.factor = 1.2,...) } } + # Nominal and ordinal axes + if(inherits(bp,"catPCA")) stop ("catPCA not implemented for 1D biplots") + # Predictions of points with lines down to the axes. if (!is.null(predict.mat)){ diff --git a/R/translate_axes.R b/R/translate_axes.R index c732687..67495a6 100644 --- a/R/translate_axes.R +++ b/R/translate_axes.R @@ -35,6 +35,9 @@ #' dist[7] <- 0.4 #' bp |> translate_axes(delta = 0.2, distances=dist) |> plot() translate_axes<-function(bp, delta=0,swop=FALSE,distances=NULL){ + if (!requireNamespace("cluster", quietly = TRUE)) { + stop("Package 'cluster' is required for this function. Please install it.", call. = FALSE) + } p<-bp$p m<-bp$ax.one.unit[,2]/bp$ax.one.unit[,1] #slopes of the axes thetas<-atan(m) #angles of the axes diff --git a/R/utility.R b/R/utility.R index 2752c5f..9e9e86e 100644 --- a/R/utility.R +++ b/R/utility.R @@ -413,4 +413,34 @@ extended.matching.coefficient <- function (X) Dmat <- Dmat + matrix (1, nrow=n, ncol=n) - Gk %*% t(Gk) } stats::as.dist(sqrt(Dmat)) -} \ No newline at end of file +} + +#' Applies primary or secondary approach to ties +#' +#' @param d observed dissimilarities +#' @param delta given dissimilarities to match (order) +#' @param method either \code{primary} or \code{secondary} approach to ties +#' +#' @returns a vector of ordered values +#' +ties <- function (d, delta, method = c("primary", "secondary")) +{ + method <- method[1] + if (!(identical(method, "primary") | identical(method, "secondary"))) stop("method must be either 'primary' or 'secondary' ") + if (is.null(names(delta))) labels <- 1:length(delta) else labels <- names(delta) + pairs <- data.frame(d, delta, labels) + if (identical(method, "primary")) pairs <- pairs[order(d, delta), ] + if (identical(method, "secondary")) pairs <- pairs[order(d), ] + z <- as.vector(pairs[, 2]) + names(z) <- pairs[, 3] + while (!all(order(z) == (1:length(z)))) + { block.vec <- rep(1, length(z)) + block.num <- 1 + for (i in 2:length(z)) + { if (z[i] > z[i - 1]) block.num <- block.num + 1 + block.vec[i] <- block.num + } + for (i in 1:max(block.vec)) z[block.vec == i] <- mean(z[block.vec == i]) + } + z +} diff --git a/README.Rmd b/README.Rmd index 55cde2b..ba90985 100644 --- a/README.Rmd +++ b/README.Rmd @@ -15,9 +15,6 @@ knitr::opts_chunk$set( # biplotEZ - - - The goal of biplotEZ is to provide users an EZ-to-use platform for visually representing their data with biplots. Currently, this package includes principal component analysis (PCA) and canonical variate diff --git a/README.md b/README.md index 362eca6..65d3167 100644 --- a/README.md +++ b/README.md @@ -13,6 +13,7 @@ analysis (CVA) biplots. This is accompanied by various formatting options for the samples and axes. Alpha-bags and concentration ellipses are included for visual enhancements and interpretation. + ## Installation You can install the development version of biplotEZ like this: diff --git a/_pkgdown.yml b/_pkgdown.yml new file mode 100644 index 0000000..5b1d790 --- /dev/null +++ b/_pkgdown.yml @@ -0,0 +1,21 @@ +template: + bootstrap: 5 +navbar: + components: + articles: + text: Articles + menu: + - text: biplotEZ + href: articles/biplotEZ.html + - text: biplotEZ enhancements + href: articles/biplotEZ_enhancements.html + - text: Multi-dimensional scaling + href: articles/MDS.html + - text: Biplots for separating classes + href: articles/Class_separation.html + - text: Categorical biplots + href: articles/CA_in_biplotEZ.html + - text: One dimensional biplots + href: articles/Biplots_in_1D.html + - text: Three dimensional biplots + href: articles/Biplots_in_3D.html \ No newline at end of file diff --git a/biplotEZ.Rproj b/biplotEZ.Rproj index 21a4da0..4b726bb 100644 --- a/biplotEZ.Rproj +++ b/biplotEZ.Rproj @@ -1,4 +1,5 @@ Version: 1.0 +ProjectId: 19a4037d-9088-49bb-bddb-bbc5a80cabdc RestoreWorkspace: Default SaveWorkspace: Default diff --git a/docs/404.html b/docs/404.html index 54624f2..547fff5 100644 --- a/docs/404.html +++ b/docs/404.html @@ -4,127 +4,88 @@ - + Page not found (404) • biplotEZ - - - - - - - - - - - - + + + + + + + + + + - - - - -
-
-
- +
+
+
Content not found. Please use links in the navbar. -
- - - +
- -
- - + diff --git a/docs/CNAME b/docs/CNAME new file mode 100644 index 0000000..09e131f --- /dev/null +++ b/docs/CNAME @@ -0,0 +1 @@ +www.biplotez.org \ No newline at end of file diff --git a/docs/LICENSE-text.html b/docs/LICENSE-text.html index 7479159..d776643 100644 --- a/docs/LICENSE-text.html +++ b/docs/LICENSE-text.html @@ -1,98 +1,65 @@ -License • biplotEZ - - -
-
-
- +
+
+
YEAR: 2024
 COPYRIGHT HOLDER: biplotEZ authors
 
-
+
- +
- - -
-
- - diff --git a/docs/LICENSE.html b/docs/LICENSE.html index ddc7ee5..817c01d 100644 --- a/docs/LICENSE.html +++ b/docs/LICENSE.html @@ -1,68 +1,43 @@ -MIT License • biplotEZ - - -
-
-
- +
+
+
@@ -73,30 +48,22 @@

MIT License

THE SOFTWARE IS PROVIDED “AS IS”, WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

-
+
- +
- - -
-
- - diff --git a/docs/articles/Biplots_in_1D.html b/docs/articles/Biplots_in_1D.html index 3f117d7..0a41d86 100644 --- a/docs/articles/Biplots_in_1D.html +++ b/docs/articles/Biplots_in_1D.html @@ -4,101 +4,72 @@ - + Biplots in 1D • biplotEZ - - - - - - - - - - - - - - - + + + + + + + + + + - - + + Skip to contents -
-
-
-

The print() function provides the misclassification rate, as well as the confusion matrix, of the CVA classifier applied to @@ -160,8 +135,28 @@

Classification using classify()

-

Adding new samples +

One-dimensional CA biplot +

+

Similarly for a CA 1d biplot, dim.biplot = 1 is +specified in the call to CA(). The one-dimensional biplot +is constructed from the first columns of +𝐔𝚲𝛄\mathbf{U\Lambda^\gamma} +and $&\mathbf{V\Lambda^{1-\gamma}}$.

+

Consider the HairEyeColor example again as discussed in +CA in biplotEZ:

+
+biplot(HairEyeColor[,,2], center = FALSE) |> CA(variant = "Princ", 
+  dim.biplot=1, lambda.scal = T) |> plot()
+

+
+
+

The function interpolate()

+

interpolate() allows for new observations or axes to be +added to the biplot.

+
+

The function newsamples() +

The process of adding new samples to the biplot, called interpolation utilises the functions interpolate() and newsamples(). These functions work in the same way as in @@ -171,14 +166,14 @@

Adding new samplesnewsamples() operates the same way as samples() in that the user can specify the aesthetics of the interpolated samples.

-
+
 biplot(iris[c(1:50,101:150),1:4])|> PCA(dim.biplot = 1) |> axes(col="black") |> 
   interpolate(newdata = iris[51:100,1:4]) |> newsamples(col="purple") |> plot()
-

+

-
-

Adding new variables -

+
+

The function newaxes() +

To interpolate new variables to the biplot, the function interpolate() and newaxes() are called. The function interpolate() accepts the argument @@ -187,15 +182,19 @@

Adding new variablesnewaxes() allows the user to specify the aesthetics of the interpolated variables.

-
+
 biplot(iris[,1:3])|> PCA(dim.biplot = 1) |> axes(col="black") |> 
   interpolate(newvariable = iris[,4]) |> 
   newaxes(col="darkred",X.new.names = "Petal.Width") |> plot()
-

+

+

-

Predicting samples +

The function prediction()

+
+

Predicting Samples +

To add the prediction of samples on the biplot, the prediction() function is used. The predict.samples argument takes in a vector indicating @@ -204,16 +203,16 @@

Predicting samplesaxes() function: predict.col and predict.lwd.

-
+
 biplot(iris) |> 
   PCA(group.aes = iris$Species,dim.biplot = 1,show.class.means = TRUE) |> 
   axes(col="black",predict.col = "darkred") |> 
   prediction(predict.samples=100:150) |> plot()
-

+

-
-

Predicting group means -

+
+

Predicting Group Means +

Similarly, to add the prediction of group means, the function prediction() is used. The argument predict.means takes in a vector specifying which group @@ -221,12 +220,13 @@

Predicting group meansshow.class.means must be set to TRUE in the PCA() function.

-
+
 biplot(iris) |> 
   PCA(group.aes = iris$Species,dim.biplot = 1,show.class.means = TRUE) |> 
-  axes(col="black",predict.col = "darkred") |> 
+  axes(col="black",predict.col = "darkred") |> means(label=TRUE,which=1:3)|>
   prediction(predict.means = 1) |> plot()
-

+

+

Ellipses and Alpha bags @@ -236,49 +236,117 @@

Ellipses and Alpha bags -
+
 biplot(iris) |> PCA(group.aes = iris[,5],dim.biplot = 1) |> 
   axes(col="black") |> 
   ellipses() |> plot()
 #> Computing 1.96 -ellipse for setosa 
 #> Computing 1.96 -ellipse for versicolor 
 #> Computing 1.96 -ellipse for virginica
-

+

The one-dimensional representation of an Alpha bag will simply be an empirical interval. The empirical interval is indicated using rectangles spanning the range of the interval.

-
+
 biplot(iris) |> PCA(group.aes = iris[,5],dim.biplot = 1) |> 
   axes(col="black") |> 
   alpha.bags(alpha = 0.7) |> plot()
 #> Computing 0.7 -bag for setosa 
 #> Computing 0.7 -bag for versicolor 
 #> Computing 0.7 -bag for virginica
-

+

-

Group Densities +

The function density1D()

-
+

Overlapping points make the distribution of points on the scatter +line difficult to identify. density1D() uses kernel density +estimation (KDE), which adds a density plot to the one-dimensional +biplot.

+
+biplot(iris) |> PCA(dim.biplot = 1) |> 
+  axes(col="black") |>
+  density1D() |> plot()
+

+

This KDE may be too smooth to display the distribution of the data. +By changing the parameters of the KDE, we are able to address this +issue. The bandwidth and kernel used in density1D() are +controlled by the arguments h= and kernel=, +respectively. The bandwidth h can take any positive value, +see ?stats::density for more detail. kernel +can take on any kernel supported by stats::density().

+
+biplot(iris) |> PCA(dim.biplot = 1) |> axes(col='black') |>
+  density1D(h = 0.5 ,kernel = "triangular") |> plot()
+

+

The high concentration of observations in the right of the plot now +becomes evident. To further explore the distributions of the +observations, we may want to explore the density of groupings in the +data. To do this, simply specify the group.aes= argument in +PCA(). Here the density of the three species of iris is +displayed.

+
 biplot(iris) |> PCA(group.aes = iris[,5],dim.biplot = 1) |> 
   axes(col="black") |>
-  density1D(which = c(1,2,3)) |> plot()
-

+ density1D() |> plot()
+

+

To only display the density of certain groups, use the +which= argument in density1D().

+
+biplot(iris) |> PCA(group.aes = iris[,5],dim.biplot = 1) |> 
+  axes(col="black") |>
+  density1D(which = c(2,3)) |> plot()
+

-

Legend +

The function legend.type()

-
-biplot(iris) |> CVA(classes  = iris[,5],dim.biplot = 1, show.class.means = TRUE) |> 
-  axes(col="black") |> density1D() |> samples(opacity=0.5)|> classify(borders = TRUE)|> alpha.bags()|>
-  legend.type(samples = TRUE,means = TRUE, bags = TRUE, ellipses = TRUE,new = TRUE) |> plot()
+

legend.type adds a legend to the plot. A separate legend +is created for each of the elements by setting each of +samples, means, bags and +ellipses equal to TRUE. Here, we add a legend +for the samples.

+
+biplot(iris) |> PCA(group.aes  = iris[,5],dim.biplot = 1, show.class.means = TRUE) |>
+  axes(col="black") |> density1D() |> samples(opacity=0.5)|> alpha.bags()|>
+  legend.type(samples = TRUE) |> plot()
+#> Computing 0.95 -bag for setosa 
+#> Computing 0.95 -bag for versicolor 
+#> Computing 0.95 -bag for virginica
+

+

If other legends are added, they will overlap with the elements of +the plot as displayed below.

+
+biplot(iris) |> PCA(group.aes  = iris[,5],dim.biplot = 1, show.class.means = TRUE) |>
+  axes(col="black") |> density1D() |> samples(opacity=0.5)|> alpha.bags()|>
+  legend.type(samples = TRUE,means = TRUE, bags = TRUE) |> plot()
+#> Computing 0.95 -bag for setosa 
+#> Computing 0.95 -bag for versicolor 
+#> Computing 0.95 -bag for virginica
+

+

By specifying new=TRUE, the legends will be displayed on +a new plot.

+
+biplot(iris) |> PCA(group.aes  = iris[,5],dim.biplot = 1, show.class.means = TRUE) |>
+  axes(col="black") |> density1D() |> samples(opacity=0.5)|> alpha.bags()|>
+  legend.type(samples = TRUE,means = TRUE, bags = TRUE, new=TRUE) |> plot()
 #> Computing 0.95 -bag for setosa 
 #> Computing 0.95 -bag for versicolor 
 #> Computing 0.95 -bag for virginica
-

+

+

For the CVA biplot legend.type also displays a legend +for classification regions if `

+
+bp <- biplot(iris) |> CVA(classes  = iris[,5],dim.biplot = 1, show.class.means = TRUE) |>
+  axes(col="black") |> classify() |> density1D() |> samples(opacity=0.5)|> alpha.bags()|>
+  legend.type(samples = TRUE,means = TRUE, bags = TRUE, regions = TRUE, new=TRUE) |> plot() #
+#> Computing 0.95 -bag for setosa 
+#> Computing 0.95 -bag for versicolor 
+#> Computing 0.95 -bag for virginica
+

-

Fit measures +

The function fit.measures()

fit.measures() calculates measures of fit for the biplot. Passing a biplot object which has been piped to @@ -289,21 +357,21 @@

Fit measures
+
 a <- biplot(iris) |> PCA(group.aes = iris[,5],dim.biplot = 1) |> fit.measures()
 summary(a)
 #> Object of class biplot, based on 150 samples and 5 variables.
 #> 4 numeric variables.
 #> 1 categorical variable.
 #> 
-#> Quality of fit = 92.5% 
-#> Adequacy of variables:
+#> Quality of fit in 1 dimension(s) = 92.5% 
+#> Adequacy of variables in 1 dimension(s):
 #> Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
 #>  0.130600269  0.007144055  0.733884527  0.128371149 
-#> Axis predictivity:
+#> Axis predictivity in 1 dimension(s):
 #> Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
 #>    0.8053299    0.1590003    0.9957524    0.9342141 
-#> Sample predictivity:
+#> Sample predictivity in 1 dimension(s):
 #>            1            2            3            4            5            6 
 #> 0.9859320755 0.9885358364 0.9974032622 0.9858701194 0.9843168833 0.8999107868 
 #>            7            8            9           10           11           12 
@@ -359,35 +427,26 @@ 

Fit measuresReferences

-

- - - +
-
- - + diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-10-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-10-1.png index d8c8b67..9e855a8 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-10-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-10-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-11-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-11-1.png index 3da834c..a51863e 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-11-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-11-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-12-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-12-1.png index 41ad959..dfabd38 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-12-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-12-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-13-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-13-1.png index decc530..1e34061 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-13-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-13-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-14-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-14-1.png new file mode 100644 index 0000000..74f4aba Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-14-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-15-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-15-1.png new file mode 100644 index 0000000..eaa98e7 Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-15-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-16-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-16-1.png new file mode 100644 index 0000000..c8fe91f Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-16-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-17-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-17-1.png new file mode 100644 index 0000000..4edacb2 Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-17-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-18-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-18-1.png new file mode 100644 index 0000000..3461014 Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-18-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-19-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-19-1.png new file mode 100644 index 0000000..d0f5bed Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-19-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-19-2.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-19-2.png new file mode 100644 index 0000000..1c6e9be Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-19-2.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-2-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-2-1.png index 75a41d9..80cec2c 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-2-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-2-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-1.png new file mode 100644 index 0000000..5430e1a Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-2.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-2.png new file mode 100644 index 0000000..1c6e9be Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-2.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-3.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-3.png new file mode 100644 index 0000000..142c83d Binary files /dev/null and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-20-3.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-3-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-3-1.png index 0e68c16..67b4e74 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-3-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-3-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-4-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-4-1.png index ff84862..4821d86 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-4-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-4-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-6-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-6-1.png index e9dc6b8..3e7bb92 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-6-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-6-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-7-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-7-1.png index 438571b..880f082 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-7-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-7-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-8-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-8-1.png index 3887827..d88e20a 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-8-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-8-1.png differ diff --git a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-9-1.png b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-9-1.png index 37fee28..06da112 100644 Binary files a/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-9-1.png and b/docs/articles/Biplots_in_1D_files/figure-html/unnamed-chunk-9-1.png differ diff --git a/docs/articles/Biplots_in_3D.html b/docs/articles/Biplots_in_3D.html index 5cba260..22f1c5d 100644 --- a/docs/articles/Biplots_in_3D.html +++ b/docs/articles/Biplots_in_3D.html @@ -4,108 +4,78 @@ - + Biplots in 3D • biplotEZ - - - - - - - - - - - - - - - + + + + + + + + + + - - + + Skip to contents -
-
-
-
-