diff --git a/DESCRIPTION b/DESCRIPTION
index d7b4663..0984cc4 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -30,10 +30,14 @@ URL: https://github.com/stmcg/learner
BugReports: https://github.com/stmcg/learner/issues
Imports: ScreeNOT,
Rcpp (>= 1.0.11),
- RcppEigen
+ RcppEigen,
+ graphics,
+ grDevices,
+ stats
LinkingTo: Rcpp, RcppEigen
Suggests:
- testthat (>= 3.0.0)
+ testthat (>= 3.0.0),
+ ggplot2 (>= 3.4.0)
Config/testthat/edition: 3
Depends:
R (>= 2.10)
diff --git a/NAMESPACE b/NAMESPACE
index 1ab3843..02814b0 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -3,6 +3,10 @@
export(cv.learner)
export(dlearner)
export(learner)
+export(plot_contributions)
+export(plot_cv)
+export(plot_objective)
+export(plot_source_spaces)
import(RcppEigen)
importFrom(Rcpp,sourceCpp)
useDynLib(learner, .registration = TRUE)
diff --git a/NEWS.md b/NEWS.md
index e9ac59c..3631f30 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -6,6 +6,8 @@
populations with independently generated noise.
* Add separate row- and column-space penalties in `learner()` and
cross-validation over these penalties in `cv.learner()`.
+* Add plotting functions for objective function values, cross-validation
+ results, population latent spaces, and variable contributions.
### learner version 1.0.0 (2025-03-02)
diff --git a/R/plot-spaces.R b/R/plot-spaces.R
new file mode 100644
index 0000000..d4287bc
--- /dev/null
+++ b/R/plot-spaces.R
@@ -0,0 +1,189 @@
+# Internal helpers shared by the space and contribution plots.
+check_plot_matrix <- function(x, name) {
+ if (!is.matrix(x) || !is.numeric(x) || is.complex(x) || any(dim(x)==0) || any(!is.finite(x))) {
+ stop(name, ' must be a complete finite real numeric matrix.')
+ }
+}
+
+plot_indices <- function(index, n, max_display, name) {
+ if(is.null(index)) return(unique(as.integer(round(seq(1,n,length.out=min(n,max_display))))))
+ if(!is.numeric(index)||is.complex(index)||length(index)==0||any(!is.finite(index))||
+ any(index!=floor(index))||any(index<1|index>n)||anyDuplicated(index)) {
+ stop(name, ' must contain distinct valid integer indices.')
+ }
+ as.integer(index)
+}
+
+plot_svd <- function(x, r) {
+ s<-svd(x,nu=r,nv=r)
+ if(s$d[r] <= max(dim(x))*s$d[1]*.Machine$double.eps) {
+ stop('r includes numerically zero singular directions; use a smaller rank.')
+ }
+ s
+}
+
+plot_limits_pair <- function(limits, automatic) {
+ if(is.null(limits))return(automatic)
+ if(!is.numeric(limits)||is.complex(limits)||length(limits)!=2L||
+ any(!is.finite(limits))||any(limits<=0))stop('limits must contain two finite positive values: row then column.')
+ stats::setNames(as.numeric(limits),c('row','column'))
+}
+
+# Use a raster to avoid drawing hundreds of thousands of individual rectangles.
+draw_matrix_panel <- function(values, x_index, y_index, limits, colors, title, xlab, ylab) {
+ nx<-nrow(values);ny<-ncol(values)
+ clipped<-pmax(limits[1],pmin(limits[2],values))
+ bins<-1L+floor((clipped-limits[1])/diff(limits)*(length(colors)-1L))
+ pixels<-matrix(colors[pmax(1L,pmin(length(colors),bins))],nx,ny)
+ graphics::plot.new()
+ # Reserve a constant fraction of the panel width for the color key.
+ graphics::plot.window(xlim=c(0.5,nx+0.5+nx*0.15),ylim=c(0.5,ny+0.5),xaxs='i',yaxs='i')
+ graphics::rasterImage(grDevices::as.raster(t(pixels[,ny:1,drop=FALSE])),
+ 0.5,0.5,nx+0.5,ny+0.5,interpolate=FALSE)
+ graphics::rect(0.5,0.5,nx+0.5,ny+0.5,border='gray40',lwd=0.7)
+ ix<-unique(as.integer(round(seq(1,nx,length.out=min(5,nx)))))
+ iy<-unique(as.integer(round(seq(1,ny,length.out=min(8,ny)))))
+ graphics::axis(1,at=ix,labels=x_index[ix],cex.axis=0.85)
+ graphics::axis(2,at=iy,labels=y_index[iy],las=1,cex.axis=0.85)
+ graphics::title(main=title,xlab=xlab,ylab=ylab,cex.main=1.05,cex.lab=0.9)
+ draw_color_key(nx+0.5+nx*0.05,nx+0.5+nx*0.075,0.5,ny+0.5,limits,colors)
+}
+
+#' Compare target and source latent spaces
+#'
+#' Plot column-space and row-space projection submatrices for the target and
+#' each source, to compare their retained latent spaces.
+#'
+#' @param Y_target Complete finite target matrix, or a completed estimate.
+#' Missing values are rejected; the function does not impute them.
+#' @param Y_source Aligned complete source matrix or nonempty list of matrices.
+#' @param r Positive integer giving the retained rank for every matrix.
+#' @param row_index,column_index Optional distinct integer indices selecting
+#' variables to display, in display order. The same indices are used in all
+#' populations. SVD always uses the full matrices before subsetting.
+#' @param max_display Maximum number of displayed variables on each side when
+#' indices are omitted; default 200. Deterministic evenly spaced indices are
+#' used, without changing the random-number state.
+#' @param limits Optional two positive color limits, in row-side then column-side
+#' order. Each side uses a symmetric scale about zero shared across populations.
+#' Values beyond supplied limits are clipped for display and counts are returned.
+#' @param row_label,column_label Axis labels. Use Variant index and Phenotype
+#' index only when these accurately describe the input rows and columns.
+#' @param target_name Panel label describing the supplied target matrix.
+#' @return Invisibly, a list with \code{row_projectors}, \code{column_projectors},
+#' \code{row_index}, \code{column_index}, \code{limits}, and \code{clipped_count}.
+#' Each projector is the displayed submatrix, computed without allocating the
+#' full row-by-row projector. Its population name identifies the input.
+#' @details
+#' For an SVD basis B, the displayed matrix is a submatrix of B B-transpose.
+#' Column-side panels appear on top, row-side panels below. Red denotes negative
+#' values and blue positive values. Color scales are shared within each side,
+#' including across pages, but can differ between the two sides. At most three
+#' populations appear per page. Source names are taken from the input list when
+#' supplied. No name-based alignment or matching of factors is performed.
+#' These are individual population projectors, not the weighted source operator.
+#' A submatrix need not itself be an idempotent projector.
+#' Numerically zero retained singular directions are rejected. With repeated
+#' singular values at the truncation boundary the retained space is not unique.
+#' The input order must already encode any desired phenotype or genomic ordering.
+#' @examples
+#' plot_source_spaces(dat_multisource$Y_target, dat_multisource$Y_source,
+#' r = 3, max_display = 20)
+#' @export
+plot_source_spaces <- function(Y_target, Y_source, r, row_index=NULL, column_index=NULL,
+ max_display=200L, limits=NULL, row_label='Row index',
+ column_label='Column index', target_name='Target population') {
+ check_plot_matrix(Y_target,'Y_target')
+ sources<-prepare_sources(Y_source,Y_target,allow_missing_target=FALSE)$matrices
+ for(s in sources)check_plot_matrix(s,'Y_source')
+ r<-resolve_rank(r,Y_target);validate_count(max_display,'max_display')
+ ri<-plot_indices(row_index,nrow(Y_target),max_display,'row_index')
+ ci<-plot_indices(column_index,ncol(Y_target),max_display,'column_index')
+ labels<-names(sources)
+ if(is.null(labels))labels<-rep('',length(sources))
+ labels[is.na(labels)]<-''
+ labels<-ifelse(nzchar(labels),paste('Source',labels),paste('Source',seq_along(sources)))
+ populations<-c(list(Y_target),sources);names(populations)<-make.unique(c(target_name,labels))
+ bases<-lapply(populations,plot_svd,r=r)
+ up<-lapply(bases,function(s)tcrossprod(s$u[ri,,drop=FALSE]))
+ vp<-lapply(bases,function(s)tcrossprod(s$v[ci,,drop=FALSE]))
+ auto<-c(row=max(vapply(up,function(p)max(abs(p)),numeric(1))),
+ column=max(vapply(vp,function(p)max(abs(p)),numeric(1))))
+ auto[auto==0]<-1e-8
+ cap<-plot_limits_pair(limits,auto)
+ counts<-list(row=vapply(up,function(p)sum(abs(p)>cap[1]),integer(1)),
+ column=vapply(vp,function(p)sum(abs(p)>cap[2]),integer(1)))
+ colors<-grDevices::colorRampPalette(c('firebrick3','white','dodgerblue3'))(101)
+ old<-graphics::par(no.readonly=TRUE);on.exit(graphics::par(old),add=TRUE)
+ for(page in split(seq_along(populations),ceiling(seq_along(populations)/3))) {
+ graphics::par(mfrow=c(2,length(page)),mar=c(3.5,3.8,2.5,4),oma=c(2,0,0,0),mgp=c(2.2,0.6,0))
+ for(i in page)draw_matrix_panel(vp[[i]],ci,ci,c(-cap[2],cap[2]),colors,names(populations)[i],column_label,column_label)
+ for(i in page)draw_matrix_panel(up[[i]],ri,ri,c(-cap[1],cap[1]),colors,names(populations)[i],row_label,row_label)
+ caption<-paste('Red: negative; blue: positive. Shared scales within each row.',
+ if(sum(unlist(counts))>0)paste(sum(unlist(counts)),'displayed entries clipped across all panels.')else '')
+ graphics::mtext(caption,outer=TRUE,side=1,line=0.5,cex=0.8)
+ }
+ invisible(list(row_projectors=up,column_projectors=vp,row_index=ri,column_index=ci,
+ limits=cap,clipped_count=counts))
+}
+
+#' Plot variable contributions to latent components
+#'
+#' Show which rows and columns contribute to each latent component of a fitted
+#' target matrix, using squared singular-vector entries. Higher scores are darker red.
+#'
+#' @param x A complete estimated matrix, a \code{learner()} result, or a
+#' \code{dlearner()} result. No rotation of the singular vectors is performed.
+#' @param r Number of components. Required for a matrix; defaults to the stored
+#' rank for a fitting result. Numerically zero directions are rejected.
+#' @inheritParams plot_source_spaces
+#' @param main Title describing the input estimate.
+#' @return Invisibly, a list with full \code{row_scores} and \code{column_scores},
+#' displayed \code{row_index} and \code{column_index}, \code{limits}, and
+#' \code{clipped_count}. Score matrices have variables in rows and components
+#' in columns; all values are returned without display clipping.
+#' @details
+#' An SVD of the estimated matrix supplies U and V. Row scores are U squared
+#' elementwise; column scores are V squared elementwise. Scores sum to one per
+#' component over all variables on the corresponding side. They are not squared
+#' entries of the raw optimization factors, effect sizes, or correlations.
+#' Subsetting is applied only for display and scores are not renormalized.
+#' The left panel displays columns and the right panel rows. Color limits are
+#' separate for the two sides and start at zero. Optional limits specify row
+#' then column upper limits, with clipped-entry counts reported.
+#' Signs do not affect squared scores, but rotations within repeated singular
+#' subspaces can change component-specific scores. Components from different
+#' populations are not automatically aligned or given biological labels.
+#' @examples
+#' fit <- dlearner(dat_multisource$Y_source, dat_multisource$Y_target, r = 3)
+#' plot_contributions(fit,max_display=50)
+#' @export
+plot_contributions <- function(x, r=NULL, row_index=NULL, column_index=NULL,
+ max_display=200L, limits=NULL, row_label='Row index',
+ column_label='Column index', main='Estimated target contributions') {
+ if(is.list(x)) {
+ estimate<-if(!is.null(x$learner_estimate))x$learner_estimate else x$dlearner_estimate
+ if(is.null(r))r<-x$r
+ } else estimate<-x
+ check_plot_matrix(estimate,'x')
+ if(is.null(r))stop('Supply r for a matrix without a stored rank.')
+ r<-resolve_rank(r,estimate);validate_count(max_display,'max_display')
+ ri<-plot_indices(row_index,nrow(estimate),max_display,'row_index')
+ ci<-plot_indices(column_index,ncol(estimate),max_display,'column_index')
+ s<-plot_svd(estimate,r)
+ row_scores<-s$u^2;col_scores<-s$v^2
+ cap<-plot_limits_pair(limits,c(row=max(row_scores[ri,,drop=FALSE]),column=max(col_scores[ci,,drop=FALSE])))
+ cap[cap==0]<-1e-8
+ counts<-c(row=sum(row_scores[ri,,drop=FALSE]>cap[1]),column=sum(col_scores[ci,,drop=FALSE]>cap[2]))
+ colors<-grDevices::colorRampPalette(c('white','firebrick3'))(100)
+ old<-graphics::par(no.readonly=TRUE);on.exit(graphics::par(old),add=TRUE)
+ graphics::par(mfrow=c(1,2),mar=c(3.8,3.8,2.2,4),oma=c(2,0,2,0),mgp=c(2.2,0.6,0))
+ draw_matrix_panel(col_scores[ci,,drop=FALSE],ci,seq_len(r),c(0,cap[2]),colors,'Column contributions',column_label,'Component')
+ draw_matrix_panel(row_scores[ri,,drop=FALSE],ri,seq_len(r),c(0,cap[1]),colors,'Row contributions',row_label,'Component')
+ graphics::mtext(main,outer=TRUE,side=3,line=0.4,font=2,cex=1.1)
+ caption<-paste('Squared singular-vector entries; separate color scales.',
+ if(sum(counts)>0)paste(sum(counts),'displayed entries clipped.')else '')
+ graphics::mtext(caption,outer=TRUE,side=1,line=0.4,cex=0.8)
+ invisible(list(row_scores=row_scores,column_scores=col_scores,row_index=ri,column_index=ci,
+ limits=cap,clipped_count=counts))
+}
diff --git a/R/plotting.R b/R/plotting.R
new file mode 100644
index 0000000..2c3d185
--- /dev/null
+++ b/R/plotting.R
@@ -0,0 +1,431 @@
+#' Plot the objective function values of a LEARNER fit
+#'
+#' Plots the values of the objective function at each iteration of the
+#' numerical optimization algorithm in \code{learner}.
+#'
+#' @param x A result returned by \code{learner}.
+#' @param iterations Iterations to show: \code{"all"} (the default),
+#' \code{"near_optimum"}, or a vector of iteration numbers such as
+#' \code{1:1000}. \code{"near_optimum"} shows a window centered on the first
+#' smallest finite recorded objective value. The window contains 10\% of the
+#' recorded iterations (rounded up), with a minimum of 10, or all iterations
+#' for shorter histories. At either boundary the window shifts to stay within
+#' the recorded history. Both axes adjust to the selected window.
+#' @param main Title of the plot. By default, \code{"Objective function"},
+#' with \code{"(near optimum)"} added when \code{iterations = "near_optimum"}.
+#'
+#' @return Invisibly, a data frame of the plotted iterations and objective
+#' values.
+#'
+#' @details The red point marks the smallest objective value, which is also
+#' given above the plot. The objective value at the initialization is not
+#' included. With \code{iterations = "near_optimum"}, the plot shows a shaded
+#' band spanning each pair of consecutive objective values, so that values
+#' that alternate between iterations remain visible, and a line at their
+#' average.
+#'
+#' @examples
+#' fit <- learner(dat_highsim$Y_source, dat_highsim$Y_target,
+#' r = 3, lambda_1 = 1, lambda_2 = 1, step_size = 0.003,
+#' control = list(max_iter = 1000))
+#' plot_objective(fit)
+#' plot_objective(fit, iterations = "near_optimum")
+#' plot_objective(fit, iterations = 1:100)
+#'
+#' @export
+plot_objective <- function(x, iterations = "all", main = NULL) {
+ objective <- x$objective_values
+ if (!is.numeric(objective) || length(objective) == 0L ||
+ !any(is.finite(objective))) {
+ stop('x must be a result returned by learner().')
+ }
+ n <- length(objective)
+ finite <- is.finite(objective)
+ minimum <- which.min(replace(objective, !finite, Inf))
+ near <- identical(iterations, 'near_optimum')
+ if (identical(iterations, 'all')) {
+ shown <- seq_len(n)
+ } else if (near) {
+ width <- min(n, max(10L, ceiling(0.1 * n)))
+ start <- max(1L, min(minimum - floor((width - 1L) / 2), n - width + 1L))
+ shown <- seq.int(start, length.out = width)
+ } else if (is.numeric(iterations) && length(iterations) > 0L &&
+ all(iterations %in% seq_len(n))) {
+ shown <- sort(unique(as.integer(iterations)))
+ } else {
+ stop("iterations must be 'all', 'near_optimum', or iteration numbers between 1 and ",
+ n, '.')
+ }
+ if (!any(finite[shown])) {
+ stop('The selected iterations contain no finite objective values; choose a different range.')
+ }
+ if (is.null(main)) {
+ main <- if (near) 'Objective function (near optimum)' else 'Objective function'
+ }
+ values <- replace(objective[shown], !finite[shown], NA_real_)
+
+ old <- graphics::par(mar = c(5, 4.5, 4, 1.5))
+ on.exit(graphics::par(old), add = TRUE)
+ xlim <- range(shown)
+ graphics::plot(shown, values, type = 'n', xlim = xlim, xlab = 'Iteration',
+ ylab = 'Objective value', bty = 'l')
+ graphics::title(main, line = 2)
+ graphics::mtext(paste0('Minimum: ', format(objective[minimum], digits = 6),
+ ' at iteration ', minimum),
+ side = 3, line = 0.6, cex = 0.8, col = 'gray30')
+ if (near && length(shown) > 1L) {
+ # The band spans each pair of consecutive values, so alternating values
+ # appear as a band rather than as dense vertical lines.
+ nxt <- c(values[-1], values[length(values)]) # the last value pairs with itself
+ keep <- !is.na(values) & !is.na(nxt)
+ graphics::polygon(c(shown[keep], rev(shown[keep])),
+ c(pmin(values, nxt)[keep], rev(pmax(values, nxt)[keep])),
+ col = grDevices::adjustcolor('#2166AC', 0.4), border = NA)
+ graphics::lines(shown[keep], ((values + nxt) / 2)[keep], col = '#2166AC', lwd = 1.5)
+ } else {
+ graphics::lines(shown, values, col = '#2166AC', type = if (length(shown) == 1L) 'p' else 'l')
+ }
+ if (minimum %in% shown) {
+ graphics::points(minimum, objective[minimum], pch = 16, cex = 0.9, col = '#B2182B')
+ }
+ # Symbol keys distinguish the band, mean line, and minimum without colored text.
+ band_color <- grDevices::adjustcolor('#2166AC', 0.4)
+ labels <- if (near && length(shown) > 1L)
+ c('Consecutive-iteration range', 'Pairwise average') else 'Objective value'
+ colors <- if (near && length(shown) > 1L) c(band_color, '#2166AC') else '#2166AC'
+ symbols <- if (near && length(shown) > 1L) c(22, NA) else NA
+ line_types <- if (near && length(shown) > 1L) c(0, 1) else 1
+ if (minimum %in% shown) {
+ labels <- c(labels, 'Minimum recorded objective')
+ colors <- c(colors, '#B2182B')
+ symbols <- c(symbols, 16)
+ line_types <- c(line_types, 0)
+ }
+ graphics::legend('topright', inset = 0.02, legend = labels,
+ col = colors, text.col = 'gray30', pch = symbols,
+ pt.bg = band_color, pt.cex = 1.1, lty = line_types,
+ lwd = 1.5, cex = 0.8, y.intersp = 1.3, bty = 'n')
+ invisible(data.frame(iteration = shown, objective = objective[shown]))
+}
+
+#' Plot the LEARNER cross-validation scores
+#'
+#' Draw a common-penalty heatmap or separate-penalty heatmaps, with a shared
+#' color scale and a black box at the first minimum in R array order.
+#'
+#' @param x A result returned by \code{cv.learner()}.
+#' @param lambda_1_all Common-penalty candidates in the exact order used for CV.
+#' Required for a two-dimensional score array; omit for separate penalties.
+#' @param lambda_2_all Balance-penalty candidates in the exact order used for CV.
+#' @param lambda_1_row_all,lambda_1_col_all Separate-penalty candidates in the
+#' exact order used for CV. Both are required for a three-dimensional array.
+#' @param n_observed Optional positive integer giving the number of observed
+#' target entries used in CV. If supplied, scores are divided by this count and
+#' labeled Validation MSE. Otherwise the original sums are labeled Validation
+#' SSE. Missing target entries must not be included in this count.
+#' @param lambda_2_index Optional vector of distinct integer indices into
+#' \code{lambda_2_all}, selecting panels to draw in separate-penalty mode.
+#' \code{NULL} draws all panels. Not used in common-penalty mode.
+#'
+#' @param log_labels Logical scalar. If TRUE, display base-10 logarithms of
+#' the candidates on the heatmap axes.
+#' All candidates on those axes must be strictly positive. Log labels are
+#' rounded to two decimal places. Defaults to FALSE.
+#'
+#' @param relative_limit Optional positive relative increase cutoff, default 0.3
+#' (30 percent). Larger finite errors are gray and labeled above the cutoff,
+#' never as diverged. NULL keeps all finite values in the color scale.
+#' @param show_values Whether to print scores inside cells. Defaults to TRUE.
+#' @param main Plot title.
+#'
+#' @return Invisibly, a list with \code{scores}, \code{scale}, \code{grids},
+#' \code{best_index}, \code{shown_lambda_2_index}, \code{color_limits}, and
+#' \code{unavailable_count}, \code{log_labels}, \code{score_range}, \code{minimum},
+#' \code{exact_minima}, \code{boundary}, and \code{lambda_2_slice_range}.
+#' Indices refer to the original candidate order. Boundary flags indicate a
+#' selected numerical endpoint on a grid with at least two distinct values.
+#' Also returns \code{display_data}, \code{display_limits}, \code{relative},
+#' \code{above_limit_count}, \code{plot} (the first ggplot page), and \code{plots}
+#' (all ggplot pages). The original absolute score diagnostics are unchanged.
+#'
+#' @details
+#' Existing CV results do not store the full candidate grids, so callers must
+#' provide them. Lengths and the selected penalty values are checked, but these
+#' checks cannot detect every incorrectly supplied grid. Use the original grids,
+#' without sorting or removing duplicates. Axes are equally spaced candidate
+#' positions labeled with actual values, or their base-10 logarithms when
+#' \code{log_labels=TRUE}. This changes labels only, not scores, candidate order,
+#' or cell spacing. The reversed magma palette uses pale colors for lower error
+#' and dark colors for higher error. Cell labels show increases over the minimum.
+#'
+#' Common mode shows \code{lambda_1_all} on the horizontal axis and
+#' \code{lambda_2_all} on the vertical axis. Separate mode shows the U-space
+#' penalty horizontally and the V-space penalty vertically, with one panel for
+#' each balance penalty. All panels share limits computed from the full score
+#' array, even when only some panels are selected. Up to six panels are drawn
+#' per page. Gray cells have non-finite scores and are excluded from selection.
+#' Only the first minimum is marked; a panel need not contain a black box.
+#' Captions report the score range, exact ties, and boundary selections. These are
+#' descriptive diagnostics, not significance or convergence tests.
+#'
+#' The stored \code{mse_all} values are sums of held-out squared errors. MSE
+#' conversion requires the correct \code{n_observed}; the function cannot infer
+#' this count from the scores. Neither this plot nor CV certifies convergence
+#' of the individual fits. No fitting is performed and \code{x} is unchanged.
+#' Requires optional ggplot2 >= 3.4.0; install it with
+#' \code{install.packages("ggplot2")}. Estimation and other plots do not require it.
+#' Draws with ggplot2 without modifying base graphics parameters. Each page shares
+#' one legend. A sufficiently large device is recommended for multiple panels.
+#' Relative increases are undefined when the minimum is zero; the display then
+#' uses absolute scores and ignores relative_limit. High relative error alone
+#' never establishes optimizer divergence. Non-finite scores are unavailable.
+#' Use ggplot2::ggsave(plot = result$plot, ...) to export the first page, or
+#' export each member of result$plots for multiple pages.
+#'
+#' @examples
+#' if (requireNamespace("ggplot2", quietly = TRUE) &&
+#' utils::packageVersion("ggplot2") >= "3.4.0") {
+#' row_grid <- c(0.1, 1)
+#' col_grid <- c(0.1, 2)
+#' balance_grid <- c(0.1, 1)
+#' target <- dat_highsim$Y_target[1:12, 1:8]
+#' set.seed(42)
+#' cv <- cv.learner(dat_highsim$Y_source[1:12, 1:8], target,
+#' r = 2, lambda_1_row_all = row_grid,
+#' lambda_1_col_all = col_grid, lambda_2_all = balance_grid,
+#' step_size = 0.0003, n_folds = 2,
+#' control = list(max_iter = 20))
+#' plot_cv(cv, lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+#' lambda_2_all = balance_grid, n_observed = sum(!is.na(target)))
+#'
+#' }
+#'
+#' @export
+plot_cv <- function(x, lambda_1_all = NULL, lambda_2_all = NULL,
+ lambda_1_row_all = NULL, lambda_1_col_all = NULL,
+ n_observed = NULL, lambda_2_index = NULL, log_labels = FALSE,
+ relative_limit = 0.3, show_values = TRUE,
+ main = 'Cross-validation error') {
+ d <- prepare_cv_plot(x, lambda_1_all, lambda_2_all, lambda_1_row_all,
+ lambda_1_col_all, n_observed, lambda_2_index, log_labels)
+ if (!is.null(relative_limit) && (!is.numeric(relative_limit) ||
+ is.complex(relative_limit) || length(relative_limit) != 1L ||
+ !is.finite(relative_limit) || relative_limit <= 0)) {
+ stop('relative_limit must be NULL or a finite positive number.')
+ }
+ if (!is.logical(show_values) || length(show_values) != 1L || is.na(show_values)) {
+ stop('show_values must be TRUE or FALSE.')
+ }
+ if (!is.character(main) || length(main) != 1L || is.na(main)) {
+ stop('main must be one nonmissing character string.')
+ }
+ if (!requireNamespace("ggplot2", quietly = TRUE)) {
+ stop('plot_cv() requires the optional ggplot2 package. Install it with install.packages("ggplot2").',
+ call. = FALSE)
+ }
+ if (utils::packageVersion("ggplot2") < "3.4.0") {
+ stop('plot_cv() requires ggplot2 >= 3.4.0. Update it with install.packages("ggplot2").',
+ call. = FALSE)
+ }
+ separate <- length(dim(d$scores)) == 3L
+ indices <- as.data.frame(arrayInd(seq_along(d$scores), dim(d$scores)))
+ names(indices) <- if (separate) c('x_index', 'y_index', 'slice') else c('x_index', 'y_index')
+ if (!separate) indices$slice <- 1L
+ indices$score <- as.vector(d$scores)
+ indices$available <- is.finite(indices$score)
+ # Relative increases are undefined when the best score is zero.
+ relative <- d$minimum > 0
+ indices$display <- if (relative) indices$score / d$minimum - 1 else indices$score
+ indices$above_limit <- indices$available & relative & !is.null(relative_limit)
+ if (!is.null(relative_limit)) indices$above_limit <- indices$above_limit & indices$display > relative_limit + 8 * .Machine$double.eps * max(1, relative_limit)
+ indices$fill_value <- indices$display
+ indices$fill_value[!indices$available | indices$above_limit] <- NA_real_
+ finite_fill <- indices$fill_value[is.finite(indices$fill_value)]
+ upper <- max(finite_fill)
+ limits <- c(0, if (upper > 0) upper else if (relative) 0.01 else 1e-8)
+ # One global scale across all slices, even if only a subset is displayed.
+ label_digits <- if (limits[2] <= 0.05) 2L else 1L
+ indices$label <- if (relative) sprintf(paste0('%.', label_digits, 'f%%'), 100 * indices$display) else
+ formatC(indices$display, digits = 4, format = 'g')
+ indices$label[!indices$available] <- 'NA'
+ if (any(indices$above_limit)) indices$label[indices$above_limit] <-
+ paste0('>', format(100 * relative_limit, trim = TRUE), '%')
+ indices$text_color <- ifelse(is.finite(indices$fill_value) &
+ indices$fill_value > limits[2] * 0.55, 'white', '#202020')
+ best_slice <- if (separate) d$best_index[3] else 1L
+ indices$selected <- indices$x_index == d$best_index[1] &
+ indices$y_index == d$best_index[2] & indices$slice == best_slice
+ panels <- if (separate) d$shown_lambda_2_index else 1L
+ labels_x <- if (separate) d$grids$row else d$grids$common
+ labels_y <- if (separate) d$grids$col else d$grids$balance
+ xlab <- if (separate) expression(lambda[1*','*plain(row)]) else expression(lambda[1])
+ ylab <- if (separate) expression(lambda[1*','*plain(col)]) else expression(lambda[2])
+ if (log_labels) {
+ labels_x <- round(log10(labels_x), 2)
+ labels_y <- round(log10(labels_y), 2)
+ xlab <- bquote(log[10](.(xlab[[1]])))
+ ylab <- bquote(log[10](.(ylab[[1]])))
+ }
+ metric <- if (is.null(n_observed)) 'SSE' else 'MSE'
+ selected_values <- vapply(seq_along(d$grids), function(i) d$grids[[i]][d$best_index[i]], numeric(1))
+ selected_names <- if (separate) c('row', 'col', 'balance') else c('lambda1', 'lambda2')
+ subtitle <- paste0('Selected ', paste(paste0(selected_names, ' = ',
+ trimws(formatC(selected_values, digits = 4, format = 'g'))), collapse = ', '),
+ ' | Minimum ', metric, ' = ', formatC(d$minimum, digits = 6, format = 'g'))
+ caption <- paste0('Black box: selected parameters. Cell labels rounded.\n', metric, ' range: ',
+ formatC(d$score_range[1], digits = 6, format = 'g'), ' to ',
+ formatC(d$score_range[2], digits = 6, format = 'g'), '.')
+ if (any(d$boundary)) caption <- paste0(caption, ' Selection on grid boundary.')
+ if (d$exact_minima > 1) caption <- paste0(caption, '\n', d$exact_minima,
+ ' exact minima; first in CV order selected.')
+ if (any(indices$above_limit)) caption <- paste0(caption, '\nGray >',
+ format(100 * relative_limit, trim = TRUE), '%: high error, not a divergence diagnosis.')
+ if (d$unavailable_count > 0) caption <- paste0(caption, '\nGray NA: unavailable score.')
+ if (!relative) caption <- paste0(caption, '\nMinimum is zero; colors and labels show absolute ', metric, '.')
+ breaks <- seq(0, limits[2], length.out = 5L)
+ key_labels <- if (relative) sprintf(paste0('%.', label_digits, 'f%%'), 100 * breaks) else
+ formatC(breaks, digits = 4, format = 'g')
+ # Declare aesthetic bindings for package checks; ggplot evaluates them in data.
+ x_index <- y_index <- fill_value <- label <- text_color <- NULL
+ plots <- list()
+ for (page in split(panels, ceiling(seq_along(panels) / 6))) {
+ data <- indices[indices$slice %in% page, , drop = FALSE]
+ facet_labels <- paste0('lambda[2] == ', format(d$grids$balance[page], digits = 4, trim = TRUE))
+ if (anyDuplicated(facet_labels)) {
+ facet_labels <- paste0('lambda[2] == ', format(d$grids$balance[page], digits = 17, trim = TRUE))
+ }
+ # Grid indices identify panels; rounded display labels must not merge them.
+ data$panel <- factor(data$slice, levels = page)
+ panel_labels <- stats::setNames(facet_labels, as.character(page))
+ plot <- ggplot2::ggplot(data, ggplot2::aes(x = x_index, y = y_index, fill = fill_value)) +
+ ggplot2::geom_tile(color = 'white', linewidth = 0.35, width = 1, height = 1) +
+ ggplot2::scale_fill_viridis_c(option = 'magma', direction = -1,
+ na.value = 'grey85', limits = limits, breaks = breaks, labels = key_labels,
+ name = if (relative) paste0('CV ', metric, '\nabove minimum') else paste('CV', metric)) +
+ ggplot2::scale_x_continuous(breaks = seq_along(labels_x), labels = labels_x,
+ expand = ggplot2::expansion(mult = 0)) +
+ ggplot2::scale_y_continuous(breaks = seq_along(labels_y), labels = labels_y,
+ expand = ggplot2::expansion(mult = 0)) +
+ ggplot2::labs(title = main, subtitle = subtitle, x = xlab, y = ylab, caption = caption) +
+ ggplot2::theme_minimal(base_size = 11) +
+ ggplot2::theme(panel.grid = ggplot2::element_blank(),
+ plot.title = ggplot2::element_text(face = 'bold', size = 14),
+ plot.subtitle = ggplot2::element_text(size = 10, margin = ggplot2::margin(b = 10)),
+ plot.caption = ggplot2::element_text(size = 9, hjust = 0, margin = ggplot2::margin(t = 10)),
+ legend.title = ggplot2::element_text(size = 10),
+ legend.text = ggplot2::element_text(size = 9),
+ axis.text = ggplot2::element_text(color = '#303030', size = 10),
+ strip.text = ggplot2::element_text(size = 11),
+ plot.margin = ggplot2::margin(12, 18, 10, 12))
+ if (show_values) plot <- plot +
+ ggplot2::geom_text(ggplot2::aes(label = label, color = text_color), size = 2.8) +
+ ggplot2::scale_color_identity()
+ plot <- plot + ggplot2::geom_tile(data = data[data$selected, , drop = FALSE],
+ fill = NA, color = 'black', linewidth = 0.9, width = 0.94, height = 0.94)
+ if (separate) plot <- plot + ggplot2::facet_wrap(~panel, ncol = min(3L, length(page)),
+ labeller = ggplot2::as_labeller(panel_labels, default = ggplot2::label_parsed))
+ print(plot)
+ plots[[length(plots) + 1L]] <- plot
+ }
+ d$display_data <- indices
+ d$display_limits <- limits
+ d$relative <- relative
+ d$above_limit_count <- sum(indices$above_limit)
+ d$plots <- plots
+ d$plot <- plots[[1L]]
+ invisible(d)
+}
+
+prepare_cv_plot <- function(x, lambda_1_all = NULL, lambda_2_all = NULL,
+ lambda_1_row_all = NULL, lambda_1_col_all = NULL,
+ n_observed = NULL, lambda_2_index = NULL, log_labels = FALSE) {
+ if (!is.list(x) || !is.numeric(x$mse_all) || is.complex(x$mse_all) ||
+ !(length(dim(x$mse_all)) %in% c(2L, 3L)) || any(dim(x$mse_all) < 1L)) {
+ stop('x must contain a nonempty two- or three-dimensional mse_all array.')
+ }
+ if (!is.logical(log_labels) || length(log_labels) != 1L || is.na(log_labels)) {
+ stop('log_labels must be TRUE or FALSE.')
+ }
+ scores <- x$mse_all
+ finite <- is.finite(scores)
+ if (!any(finite)) stop('No finite validation scores are available to plot.')
+ if (any(scores[finite] < 0)) stop('Validation squared-error scores cannot be negative.')
+ separate <- length(dim(scores)) == 3L
+ check_grid <- function(values, size, name) {
+ validate_penalty_grid(values, name)
+ if (!is.null(dim(values)) || length(values) != size) {
+ stop(name, ' must match its score-array dimension in the original CV order.')
+ }
+ }
+ if (separate) {
+ if (!is.null(lambda_1_all)) stop('Do not supply lambda_1_all for separate penalties.')
+ check_grid(lambda_1_row_all, dim(scores)[1], 'lambda_1_row_all')
+ check_grid(lambda_1_col_all, dim(scores)[2], 'lambda_1_col_all')
+ check_grid(lambda_2_all, dim(scores)[3], 'lambda_2_all')
+ grids <- list(row = lambda_1_row_all, col = lambda_1_col_all,
+ balance = lambda_2_all)
+ panels <- if (is.null(lambda_2_index)) seq_along(lambda_2_all) else lambda_2_index
+ if (!is.numeric(panels) || length(panels) == 0L || any(!is.finite(panels)) ||
+ any(panels != floor(panels)) || any(panels < 1 | panels > length(lambda_2_all)) ||
+ anyDuplicated(panels)) {
+ stop('lambda_2_index must contain distinct valid balance-grid indices.')
+ }
+ } else {
+ if (!is.null(lambda_1_row_all) || !is.null(lambda_1_col_all) ||
+ !is.null(lambda_2_index)) stop('Separate-penalty arguments require a 3D score array.')
+ check_grid(lambda_1_all, dim(scores)[1], 'lambda_1_all')
+ check_grid(lambda_2_all, dim(scores)[2], 'lambda_2_all')
+ grids <- list(common = lambda_1_all, balance = lambda_2_all)
+ panels <- 1L
+ }
+ axis_grids <- if (separate) c(grids$row, grids$col) else c(grids$common, grids$balance)
+ if (log_labels && any(axis_grids <= 0)) {
+ stop('log_labels requires strictly positive candidates on both heatmap axes.')
+ }
+ selectable <- scores
+ selectable[!finite] <- Inf
+ best <- as.integer(arrayInd(which.min(selectable), dim(scores)))
+ expected <- if (separate) {
+ c(lambda_1_row_min = grids$row[best[1]], lambda_1_col_min = grids$col[best[2]],
+ lambda_2_min = grids$balance[best[3]])
+ } else c(lambda_1_min = grids$common[best[1]], lambda_2_min = grids$balance[best[2]])
+ for (name in names(expected)) {
+ if (!is.numeric(x[[name]]) || length(x[[name]]) != 1L ||
+ !isTRUE(all.equal(unname(x[[name]]), unname(expected[[name]]), tolerance = 1e-12))) {
+ stop('Supplied grids do not match the selected penalties stored in x: ', name)
+ }
+ }
+ scale <- 'Validation SSE'
+ if (!is.null(n_observed)) {
+ validate_count(n_observed, 'n_observed', upper = Inf)
+ scores <- scores / n_observed
+ scale <- 'Validation MSE'
+ }
+ limits <- range(scores[finite])
+ if (diff(limits) == 0) {
+ pad <- max(abs(limits[1]) * 0.01, 1e-8)
+ limits <- c(max(0, limits[1] - pad), limits[2] + pad)
+ }
+ score_range <- range(scores[finite])
+ minimum <- min(scores[finite])
+ boundary <- vapply(seq_along(grids), function(i) {
+ g <- grids[[i]]
+ length(unique(g)) > 1L && g[best[i]] %in% range(g)
+ }, logical(1))
+ names(boundary) <- names(grids)
+ slice <- if (separate) scores[best[1], best[2], ] else scores[best[1], ]
+ list(scores = scores, scale = scale, grids = grids, best_index = best,
+ shown_lambda_2_index = if (separate) panels else NULL,
+ color_limits = limits, unavailable_count = sum(!finite), log_labels = log_labels,
+ score_range = score_range, minimum = minimum, exact_minima = sum(scores[finite] == minimum),
+ boundary = boundary, lambda_2_slice_range = range(slice[is.finite(slice)]))
+}
+
+# Five ticks include the numerical endpoints even for a narrow score range.
+draw_color_key <- function(left, right, bottom, top, limits, colors) {
+ edges <- seq(bottom, top, length.out=length(colors)+1L)
+ graphics::rect(left,edges[-length(edges)],right,edges[-1],col=colors,border=NA)
+ ticks <- seq(limits[1],limits[2],length.out=5)
+ graphics::axis(4,pos=right,at=seq(bottom,top,length.out=5),
+ labels=formatC(ticks,format='g',digits=5),las=1,cex.axis=0.9,tck=-0.015)
+}
diff --git a/README.Rmd b/README.Rmd
index 595e3ff..a145bf8 100644
--- a/README.Rmd
+++ b/README.Rmd
@@ -4,12 +4,10 @@ output: github_document
-```{r, include = FALSE}
-knitr::opts_chunk$set(
- collapse = TRUE,
- comment = "#>",
- out.width = "100%"
-)
+```{r setup, include=FALSE}
+knitr::opts_chunk$set(collapse = TRUE, comment = "#>", out.width = "100%",
+ fig.width = 8, fig.height = 5,
+ fig.path = "README_files/figure-gfm/")
```
@@ -20,72 +18,171 @@ knitr::opts_chunk$set(
# learner
-The `learner` package implements transfer learning methods for low-rank matrix estimation. These methods leverage similarity in the latent row and column spaces between the source and target populations to improve estimation in the target population. The methods include the LatEnt spAce-based tRaNsfer lEaRning (LEARNER) method and the direct projection LEARNER (D-LEARNER) method described by [McGrath et al. (2024)](https://doi.org/10.48550/arXiv.2412.20605).
+The `learner` package implements transfer learning methods for low-rank matrix
+estimation. These methods use similarity in the latent row and column spaces
+of source and target populations to improve estimation in the target population.
+The methods include LatEnt spAce-based tRaNsfer lEaRning (LEARNER) and direct
+projection LEARNER (D-LEARNER), described by
+[McGrath et al. (2024)](https://doi.org/10.48550/arXiv.2412.20605).
## Installation
-You can install the released version of `learner` from [CRAN](https://CRAN.R-project.org/package=learner) with:
-``` r
+Install the released version from CRAN:
+
+```r
install.packages("learner")
```
-You can install the development version of `learner` from [GitHub](https://github.com/stmcg/learner) with:
+The examples below use the multiple-source and plotting extensions on the
+`plot-diagnostics` development branch:
-``` r
-# install.packages("devtools")
-devtools::install_github("stmcg/learner")
+```r
+# install.packages("remotes")
+remotes::install_github("stmcg/learner", ref = "plot-diagnostics")
+# Optional: required for plot_cv(), but not for fitting or the other plots.
+install.packages("ggplot2") # version 3.4.0 or later
```
-## Example
+## One source population
-We illustrate an example of how `learner` can be used. We first load the package.
-```{r, results="hide", message=FALSE}
+`dat_highsim` contains a simulated source and target with similar latent spaces.
+LEARNER fits the target while penalizing departures from the source spaces.
+Here `lambda_1` controls both space penalties, and `lambda_2` controls the
+balance between the two target factors.
+
+```{r single-source, message=FALSE}
library(learner)
+fit_single <- learner(
+ Y_source = dat_highsim$Y_source, Y_target = dat_highsim$Y_target,
+ r = 3, lambda_1 = 10, lambda_2 = 1, step_size = 0.003,
+ control = list(max_iter = 2000)
+)
+dim(fit_single$learner_estimate)
```
-In this illustration, we will use one of the toy data sets in the package (`dat_highsim`) that has a high degree of similarity between the latent spaces of the source and target populations. The object `dat_highsim` is a list which contains the observed source population data matrix `Y_source` and the target population data matrix `Y_target`. Since the data was simulated, the true values of the matrices are included in `dat_highsim` as `Theta_source` and `Theta_target`.
+D-LEARNER projects the target data directly onto the source spaces, without
+iterative optimization or penalty selection:
-#### LEARNER Method
+```{r single-direct}
+direct_single <- dlearner(dat_highsim$Y_source, dat_highsim$Y_target, r = 3)
+```
+
+## Multiple sources and separate penalties
+
+`dat_multisource` contains one target and two sources, each with 100 rows and
+50 columns. The true signals have rank 3 and share their latent spaces.
+Independent normal noise is added to each population: the noise variances are
+0.10 for the target, 0.01 for source A, and 0.50 for source B. The true signal
+matrices are also included, but are not needed for fitting.
+
+Pass the sources as a list. The weights below are illustrative choices; they
+are held fixed during cross-validation. Larger row or column penalties encourage
+the corresponding target factors to follow the combined source spaces more closely.
+
+```{r multiple-data}
+sources <- dat_multisource$Y_source
+target <- dat_multisource$Y_target
+weights <- c(0.7, 0.3)
+row_grid <- c(0, 1, 10)
+col_grid <- c(0, 1, 3)
+balance_grid <- c(0.1, 1)
+```
-We can apply the LEARNER method via the `learner` function.
+Select penalties by cross-validation. This small grid and iteration budget keep
+the example quick; use suitable grids and optimization settings for an analysis.
-This method allows for flexible patterns of heterogeneity across the source and target populations. It consequently requires specifying the tuning parameters `lambda_1` and `lambda_2` which control the degree of transfer learning between the populations. These values can be selected based on cross-validation via the `cv.learner` function. For example, we can specify candidate values of 1, 10, and 100 for `lambda_1` and `lambda_2` and select the optimal values based on cross-validation as follows:
-```{r}
-res_cv.learner <- cv.learner(Y_source = dat_highsim$Y_source,
- Y_target = dat_highsim$Y_target,
- lambda_1_all = c(1, 10, 100),
- lambda_2_all = c(1, 10, 100),
- step_size = 0.003)
-res_cv.learner$lambda_1_min
-res_cv.learner$lambda_2_min
+```{r multiple-cv}
+set.seed(815)
+cv <- cv.learner(
+ Y_source = sources, Y_target = target, r = 3, source_weights = weights,
+ lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, n_folds = 3, step_size = 0.0003,
+ control = list(max_iter = 2000)
+)
+c(row = cv$lambda_1_row_min, column = cv$lambda_1_col_min,
+ balance = cv$lambda_2_min)
```
-Next, we apply the `learner` function with these values of `lambda_1` and `lambda_2`:
-```{r}
-res_learner <- learner(Y_source = dat_highsim$Y_source,
- Y_target = dat_highsim$Y_target,
- lambda_1 = 100, lambda_2 = 1,
- step_size = 0.003)
+For a common space penalty, supply `lambda_1_all` instead of the separate row
+and column grids. Cross-validation selects penalties, not the rank or source weights.
+
+### Cross-validation results
+
+The heatmap shows how much the CV error increases above its minimum, with a
+black box marking the selected combination. Supplying `n_observed` converts
+the stored sums of squared errors to mean squared errors (MSE).
+
+```{r cv-heatmap, fig.height=5.5, fig.cap="Cross-validation error across the candidate penalties."}
+plot_cv(
+ cv, lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, n_observed = sum(!is.na(target))
+)
```
-The LEARNER estimate is given by the `learner_estimate` component in the output of the `learner` function, e.g.,
-```{r}
-res_learner$learner_estimate[1:5, 1:5]
+
+A selected value on a grid boundary can be a reason to expand that grid.
+CV scores can also come from fits that reached the iteration limit, so inspect
+the optimization settings before interpreting the selection.
+
+### Refit and inspect optimization
+
+Refit using the selected penalties and all observed target entries:
+
+```{r multiple-refit}
+fit <- learner(
+ Y_source = sources, Y_target = target, r = cv$r, source_weights = weights,
+ lambda_1_row = cv$lambda_1_row_min,
+ lambda_1_col = cv$lambda_1_col_min, lambda_2 = cv$lambda_2_min,
+ step_size = 0.0003, control = list(max_iter = 5000)
+)
+fit$convergence_criterion
```
-We can check the convergence of the numerical optimization algorithm by plotting the values of the objective function at each iteration:
-```{r, out.width="75%"}
-plot(res_learner$objective_values,
- xlab = 'Iteration', ylab = 'Objective function',
- type = 'l', col = 'blue', lwd = 2)
+
+```{r objective-history, fig.cap="Objective function values near the optimum for the refitted target model."}
+plot_objective(fit, iterations = "near_optimum")
```
+Use `iterations = "all"` (the default) for the full history, or a vector such
+as `1:500` for a specific range.
-#### D-LEARNER Method
+D-LEARNER accepts the same source list and weights:
-We can apply the D-LEARNER method via the `dlearner` function. This method makes stronger assumptions on the heterogeneity across the source and target populations. It consequently does not rely on choosing tuning parameters. The `dlearner` function can be applied as follows:
-```{r}
-res_dlearner <- dlearner(Y_source = dat_highsim$Y_source,
- Y_target = dat_highsim$Y_target)
-res_dlearner$dlearner_estimate[1:5, 1:5]
+```{r multiple-direct}
+direct <- dlearner(sources, target, r = 3, source_weights = weights)
+dim(direct$dlearner_estimate)
```
+### Compare source choices
+
+For this simulated dataset, the true target signal is available. Compare
+D-LEARNER fits at the same rank to isolate the choice of source spaces:
+
+```{r source-comparison}
+estimates <- list(
+ A_only = dlearner(sources$A, target, r = 3)$dlearner_estimate,
+ B_only = dlearner(sources$B, target, r = 3)$dlearner_estimate,
+ Weighted_A_B = direct$dlearner_estimate
+)
+source_errors <- vapply(estimates, function(estimate) {
+ norm(estimate - dat_multisource$Theta_target, "F")
+}, numeric(1))
+round(source_errors, 3)
+```
+
+Lower Frobenius error is better. In this example, source A alone gives the
+smallest error, while the weighted combination improves on source B alone.
+The weights 0.7/0.3 are illustrative, not optimized. These results describe one
+simulated dataset; they are not CV errors or evidence that combining sources
+always improves estimation.
+
+### Latent spaces and contributions
+
+Compare the target and source spaces, or inspect which variables contribute
+to the components of a fitted target matrix:
+
+```{r additional-plots, eval=FALSE}
+plot_source_spaces(target, sources, r = 3)
+plot_contributions(fit)
+```
+See `?learner`, `?cv.learner`, `?dlearner`, and `?dat_multisource` for details.
+More runnable examples are included in `inst/examples/` in the source repository.
diff --git a/README.md b/README.md
index a059c2c..3258f07 100644
--- a/README.md
+++ b/README.md
@@ -1,4 +1,3 @@
-
@@ -12,121 +11,202 @@ coverage](https://codecov.io/gh/stmcg/learner/graph/badge.svg)](https://app.code
# learner
The `learner` package implements transfer learning methods for low-rank
-matrix estimation. These methods leverage similarity in the latent row
-and column spaces between the source and target populations to improve
-estimation in the target population. The methods include the LatEnt
-spAce-based tRaNsfer lEaRning (LEARNER) method and the direct projection
-LEARNER (D-LEARNER) method described by [McGrath et
-al. (2024)](https://doi.org/10.48550/arXiv.2412.20605).
+matrix estimation. These methods use similarity in the latent row and
+column spaces of source and target populations to improve estimation in
+the target population. The methods include LatEnt spAce-based tRaNsfer
+lEaRning (LEARNER) and direct projection LEARNER (D-LEARNER), described
+by [McGrath et al. (2024)](https://doi.org/10.48550/arXiv.2412.20605).
## Installation
-You can install the released version of `learner` from
-[CRAN](https://CRAN.R-project.org/package=learner) with:
+Install the released version from CRAN:
``` r
install.packages("learner")
```
-You can install the development version of `learner` from
-[GitHub](https://github.com/stmcg/learner) with:
+The examples below use the multiple-source and plotting extensions on
+the `plot-diagnostics` development branch:
``` r
-# install.packages("devtools")
-devtools::install_github("stmcg/learner")
+# install.packages("remotes")
+remotes::install_github("stmcg/learner", ref = "plot-diagnostics")
+# Optional: required for plot_cv(), but not for fitting or the other plots.
+install.packages("ggplot2") # version 3.4.0 or later
```
-## Example
+## One source population
-We illustrate an example of how `learner` can be used. We first load the
-package.
+`dat_highsim` contains a simulated source and target with similar latent
+spaces. LEARNER fits the target while penalizing departures from the
+source spaces. Here `lambda_1` controls both space penalties, and
+`lambda_2` controls the balance between the two target factors.
``` r
library(learner)
+fit_single <- learner(
+ Y_source = dat_highsim$Y_source, Y_target = dat_highsim$Y_target,
+ r = 3, lambda_1 = 10, lambda_2 = 1, step_size = 0.003,
+ control = list(max_iter = 2000)
+)
+dim(fit_single$learner_estimate)
+#> [1] 100 50
```
-In this illustration, we will use one of the toy data sets in the
-package (`dat_highsim`) that has a high degree of similarity between the
-latent spaces of the source and target populations. The object
-`dat_highsim` is a list which contains the observed source population
-data matrix `Y_source` and the target population data matrix `Y_target`.
-Since the data was simulated, the true values of the matrices are
-included in `dat_highsim` as `Theta_source` and `Theta_target`.
+D-LEARNER projects the target data directly onto the source spaces,
+without iterative optimization or penalty selection:
+
+``` r
+direct_single <- dlearner(dat_highsim$Y_source, dat_highsim$Y_target, r = 3)
+```
-#### LEARNER Method
+## Multiple sources and separate penalties
-We can apply the LEARNER method via the `learner` function.
+`dat_multisource` contains one target and two sources, each with 100
+rows and 50 columns. The true signals have rank 3 and share their latent
+spaces. Independent normal noise is added to each population: the noise
+variances are 0.10 for the target, 0.01 for source A, and 0.50 for
+source B. The true signal matrices are also included, but are not needed
+for fitting.
-This method allows for flexible patterns of heterogeneity across the
-source and target populations. It consequently requires specifying the
-tuning parameters `lambda_1` and `lambda_2` which control the degree of
-transfer learning between the populations. These values can be selected
-based on cross-validation via the `cv.learner` function. For example, we
-can specify candidate values of 1, 10, and 100 for `lambda_1` and
-`lambda_2` and select the optimal values based on cross-validation as
-follows:
+Pass the sources as a list. The weights below are illustrative choices;
+they are held fixed during cross-validation. Larger row or column
+penalties encourage the corresponding target factors to follow the
+combined source spaces more closely.
``` r
-res_cv.learner <- cv.learner(Y_source = dat_highsim$Y_source,
- Y_target = dat_highsim$Y_target,
- lambda_1_all = c(1, 10, 100),
- lambda_2_all = c(1, 10, 100),
- step_size = 0.003)
-res_cv.learner$lambda_1_min
-#> [1] 100
-res_cv.learner$lambda_2_min
-#> [1] 1
+sources <- dat_multisource$Y_source
+target <- dat_multisource$Y_target
+weights <- c(0.7, 0.3)
+row_grid <- c(0, 1, 10)
+col_grid <- c(0, 1, 3)
+balance_grid <- c(0.1, 1)
+```
+
+Select penalties by cross-validation. This small grid and iteration
+budget keep the example quick; use suitable grids and optimization
+settings for an analysis.
+
+``` r
+set.seed(815)
+cv <- cv.learner(
+ Y_source = sources, Y_target = target, r = 3, source_weights = weights,
+ lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, n_folds = 3, step_size = 0.0003,
+ control = list(max_iter = 2000)
+)
+c(row = cv$lambda_1_row_min, column = cv$lambda_1_col_min,
+ balance = cv$lambda_2_min)
+#> row column balance
+#> 10 3 1
+```
+
+For a common space penalty, supply `lambda_1_all` instead of the
+separate row and column grids. Cross-validation selects penalties, not
+the rank or source weights.
+
+### Cross-validation results
+
+The heatmap shows how much the CV error increases above its minimum,
+with a black box marking the selected combination. Supplying
+`n_observed` converts the stored sums of squared errors to mean squared
+errors (MSE).
+
+``` r
+plot_cv(
+ cv, lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, n_observed = sum(!is.na(target))
+)
```
-Next, we apply the `learner` function with these values of `lambda_1`
-and `lambda_2`:
+
+
+A selected value on a grid boundary can be a reason to expand that grid.
+CV scores can also come from fits that reached the iteration limit, so
+inspect the optimization settings before interpreting the selection.
+
+### Refit and inspect optimization
+
+Refit using the selected penalties and all observed target entries:
``` r
-res_learner <- learner(Y_source = dat_highsim$Y_source,
- Y_target = dat_highsim$Y_target,
- lambda_1 = 100, lambda_2 = 1,
- step_size = 0.003)
+fit <- learner(
+ Y_source = sources, Y_target = target, r = cv$r, source_weights = weights,
+ lambda_1_row = cv$lambda_1_row_min,
+ lambda_1_col = cv$lambda_1_col_min, lambda_2 = cv$lambda_2_min,
+ step_size = 0.0003, control = list(max_iter = 5000)
+)
+fit$convergence_criterion
+#> [1] 1
```
-The LEARNER estimate is given by the `learner_estimate` component in the
-output of the `learner` function, e.g.,
+``` r
+plot_objective(fit, iterations = "near_optimum")
+```
+
+
+
+Use `iterations = "all"` (the default) for the full history, or a vector
+such as `1:500` for a specific range.
+
+D-LEARNER accepts the same source list and weights:
``` r
-res_learner$learner_estimate[1:5, 1:5]
-#> [,1] [,2] [,3] [,4] [,5]
-#> [1,] 0.1468889 -1.63414688 1.1857308 0.008093962 0.2196070
-#> [2,] 0.1087766 -0.24392770 -1.0061737 -0.058119721 -0.4012974
-#> [3,] -0.6213028 1.04206753 -0.6907776 0.905709331 -0.5324437
-#> [4,] 1.5811588 -3.00001490 0.3978080 -2.124285120 0.6909386
-#> [5,] 0.5046133 0.01845138 -1.4798209 -0.765601712 -0.2230227
+direct <- dlearner(sources, target, r = 3, source_weights = weights)
+dim(direct$dlearner_estimate)
+#> [1] 100 50
```
-We can check the convergence of the numerical optimization algorithm by
-plotting the values of the objective function at each iteration:
+### Compare source choices
+
+For this simulated dataset, the true target signal is available. Compare
+D-LEARNER fits at the same rank to isolate the choice of source spaces:
``` r
-plot(res_learner$objective_values,
- xlab = 'Iteration', ylab = 'Objective function',
- type = 'l', col = 'blue', lwd = 2)
+estimates <- list(
+ A_only = dlearner(sources$A, target, r = 3)$dlearner_estimate,
+ B_only = dlearner(sources$B, target, r = 3)$dlearner_estimate,
+ Weighted_A_B = direct$dlearner_estimate
+)
+source_errors <- vapply(estimates, function(estimate) {
+ norm(estimate - dat_multisource$Theta_target, "F")
+}, numeric(1))
+round(source_errors, 3)
+#> A_only B_only Weighted_A_B
+#> 2.323 17.214 5.324
```
-
+Lower Frobenius error is better. In this example, source A alone gives
+the smallest error, while the weighted combination improves on source B
+alone. The weights 0.7/0.3 are illustrative, not optimized. These
+results describe one simulated dataset; they are not CV errors or
+evidence that combining sources always improves estimation.
-#### D-LEARNER Method
+### Latent spaces and contributions
-We can apply the D-LEARNER method via the `dlearner` function. This
-method makes stronger assumptions on the heterogeneity across the source
-and target populations. It consequently does not rely on choosing tuning
-parameters. The `dlearner` function can be applied as follows:
+Compare the target and source spaces, or inspect which variables
+contribute to the components of a fitted target matrix:
``` r
-res_dlearner <- dlearner(Y_source = dat_highsim$Y_source,
- Y_target = dat_highsim$Y_target)
-res_dlearner$dlearner_estimate[1:5, 1:5]
-#> [,1] [,2] [,3] [,4] [,5]
-#> [1,] 0.0959171 -1.72143637 1.15500149 0.005478273 0.3258085
-#> [2,] 0.1077137 -0.32243480 -1.03557177 -0.108195636 -0.4093049
-#> [3,] -0.7237275 0.86634922 -0.43528206 1.105802194 -0.5961351
-#> [4,] 1.6347696 -2.91985925 0.04876288 -2.357504089 0.8127403
-#> [5,] 0.4757450 0.05665543 -1.50373470 -0.744076401 -0.2876033
+plot_source_spaces(target, sources, r = 3)
+plot_contributions(fit)
```
+
+See `?learner`, `?cv.learner`, `?dlearner`, and `?dat_multisource` for
+details. More runnable examples are included in `inst/examples/` in the
+source repository.
diff --git a/README_files/figure-gfm/cv-heatmap-1.png b/README_files/figure-gfm/cv-heatmap-1.png
new file mode 100644
index 0000000..13c0365
Binary files /dev/null and b/README_files/figure-gfm/cv-heatmap-1.png differ
diff --git a/README_files/figure-gfm/objective-history-1.png b/README_files/figure-gfm/objective-history-1.png
new file mode 100644
index 0000000..5489c95
Binary files /dev/null and b/README_files/figure-gfm/objective-history-1.png differ
diff --git a/README_files/figure-gfm/unnamed-chunk-6-1.png b/README_files/figure-gfm/unnamed-chunk-6-1.png
deleted file mode 100644
index 360e568..0000000
Binary files a/README_files/figure-gfm/unnamed-chunk-6-1.png and /dev/null differ
diff --git a/README_files/figure-gfm/unnamed-chunk-7-1.png b/README_files/figure-gfm/unnamed-chunk-7-1.png
deleted file mode 100644
index 360e568..0000000
Binary files a/README_files/figure-gfm/unnamed-chunk-7-1.png and /dev/null differ
diff --git a/inst/examples/plot-diagnostics.R b/inst/examples/plot-diagnostics.R
new file mode 100644
index 0000000..fb19984
--- /dev/null
+++ b/inst/examples/plot-diagnostics.R
@@ -0,0 +1,44 @@
+# This heatmap demonstration requires the optional ggplot2 package.
+if (!requireNamespace("ggplot2", quietly = TRUE) ||
+ utils::packageVersion("ggplot2") < "3.4.0") {
+ stop('Install ggplot2 >= 3.4.0 to run this example: install.packages("ggplot2")')
+}
+# Run after installing the plotting development version of learner.
+# Synthetic data demonstrate the interface, not statistical performance.
+# Use a small aligned subset to keep this demonstration quick.
+sources <- lapply(learner::dat_multisource$Y_source,
+ function(x) x[1:20, 1:10])
+target <- learner::dat_multisource$Y_target[1:20, 1:10]
+weights <- c(0.7, 0.3)
+row_grid <- c(0, 1, 10)
+col_grid <- c(0, 1, 3)
+balance_grid <- c(0.1, 1)
+
+set.seed(815)
+cv <- learner::cv.learner(
+ sources, target, r = 3, source_weights = weights,
+ lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, step_size = 0.0003, n_folds = 3,
+ control = list(max_iter = 500))
+learner::plot_cv(
+ cv, lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, n_observed = sum(!is.na(target)))
+
+fit <- learner::learner(
+ sources, target, r = cv$r, source_weights = weights,
+ lambda_1_row = cv$lambda_1_row_min, lambda_1_col = cv$lambda_1_col_min,
+ lambda_2 = cv$lambda_2_min, step_size = 0.0003,
+ control = list(max_iter = 2000))
+learner::plot_objective(fit, iterations = "near_optimum")
+
+# Common-penalty results are supported as well.
+set.seed(815)
+common_grid <- c(0, 1, 10)
+common_cv <- learner::cv.learner(
+ sources, target, r = 3, source_weights = weights,
+ lambda_1_all = common_grid, lambda_2_all = balance_grid,
+ step_size = 0.0003, n_folds = 3, control = list(max_iter = 500))
+learner::plot_cv(common_cv, lambda_1_all = common_grid,
+ lambda_2_all = balance_grid, n_observed = sum(!is.na(target)))
+
+# Inspect stopping codes as well as curves. Plotting does not prove convergence.
diff --git a/inst/examples/plot-paper-displays.R b/inst/examples/plot-paper-displays.R
new file mode 100644
index 0000000..00539df
--- /dev/null
+++ b/inst/examples/plot-paper-displays.R
@@ -0,0 +1,8 @@
+# Population spaces and contributions using simulated multiple-source data.
+# Run after loading the plotting development version of learner.
+space_target <- learner::dat_multisource$Y_target
+space_sources <- learner::dat_multisource$Y_source
+learner::plot_source_spaces(space_target, space_sources, r=3,
+ target_name='Target (simulated)')
+space_fit <- learner::dlearner(space_sources,space_target,r=3,source_weights=c(0.7,0.3))
+learner::plot_contributions(space_fit, main='D-LEARNER contributions (simulated data)')
diff --git a/inst/examples/plotting-guide.md b/inst/examples/plotting-guide.md
new file mode 100644
index 0000000..1e80e7c
--- /dev/null
+++ b/inst/examples/plotting-guide.md
@@ -0,0 +1,91 @@
+# Plotting guide
+
+These functions display existing fits or measured results. Plot appearance alone
+cannot establish convergence or statistical superiority.
+
+## Objective history
+
+```r
+plot_objective(fit) # All iterations
+plot_objective(fit, iterations = "near_optimum") # Around the smallest recorded value
+plot_objective(fit, iterations = 1:500) # A specific range
+```
+
+`"near_optimum"` centers a window on the first smallest finite recorded
+objective value. It shows 10% of the recorded iterations (rounded up), at least
+10, or the full history if fewer than 10 were recorded. At either end the window
+shifts to stay within the history. Both axes zoom to the selected window.
+The shaded band spans consecutive objective values and the line shows their
+pairwise average; the band is not an uncertainty interval.
+
+## CV heatmaps
+
+Only the heatmap requires the optional package ggplot2 >= 3.4.0. Install it
+with `install.packages("ggplot2")`. Estimation and all other plotting functions
+work without ggplot2.
+
+```r
+plot_cv(cv, lambda_1_all = common_grid, lambda_2_all = balance_grid,
+ n_observed = sum(!is.na(target)), log_labels = TRUE)
+```
+
+Use the exact grids supplied to `cv.learner()`. With independent penalties,
+replace `lambda_1_all` with `lambda_1_row_all` and `lambda_1_col_all`.
+Heatmaps then have one panel per balance penalty. Zero candidates require
+`log_labels = FALSE`. Heatmap positions are equally spaced grid slots.
+
+The CV object's `mse_all` contains pooled held-out squared-error sums. Supplying
+the total number of observed target entries converts these to MSE. Without that
+number the plots explicitly show SSE. Source entries are not in the denominator.
+Heatmaps use a reversed magma palette: pale cells have lower error. Labels and
+colors show percentage increases over the minimum, with a black selection box.
+Finite errors more than 30% above the minimum are gray and labeled >30%, not
+as optimizer divergence. Use `relative_limit = NULL` to color all finite errors
+or `show_values = FALSE` to hide cell text. If the minimum is zero, the heatmap
+uses absolute scores instead. Cell labels are rounded, so matching printed
+percentages do not establish an exact tie.
+
+The plots report the score range, minimum, exact ties, and boundary selections.
+A boundary minimum deserves investigation; a small change across a parameter
+is not an exact tie or evidence that the parameter is irrelevant.
+
+For a customized export, capture the returned object:
+
+```r
+result <- plot_cv(cv, lambda_1_all = common_grid, lambda_2_all = balance_grid)
+ggplot2::ggsave("cv-heatmap.png", plot = result$plot, width = 8, height = 5.5, dpi = 200)
+```
+
+For multiple pages, `result$plots` contains every ggplot page. Inspect exported
+images at the intended size to check text and legend placement.
+
+## Population latent spaces
+
+```r
+plot_source_spaces(target, list(A = source_A, B = source_B), r = 3)
+```
+
+For each population, compute an SVD and retain orthonormal bases B_U and B_V.
+Plot B_U B_U^T for the row side and B_V B_V^T for the column side. These are
+individual population projectors, not raw data matrices or the weighted average
+used by the optimizer. Shared color limits within each side support comparisons.
+Negative entries are red and positive entries blue.
+
+SVD uses full matrices before display subsetting. `row_index`, `column_index`,
+or `max_display` limit plotting cost without changing the fitted spaces. All
+populations must have aligned variables in the same order. Labels default to
+Row index and Column index; use biological labels only when appropriate.
+
+## Variable contributions
+
+```r
+plot_contributions(fit, r = 3)
+```
+
+Recompute the SVD of the estimated signal. For component j, a row's score is
+B_U[i,j]^2 and a column's score is B_V[i,j]^2. The scores sum to one per component
+on each side over the full set of variables. They are not squared optimization
+factors, which include scaling. The plot uses separate color scales on the two
+sides and reports display clipping if explicit limits are supplied.
+These component-specific scores depend on basis choice in repeated singular
+subspaces; the function does not match components or assign biological meaning.
diff --git a/man/plot_contributions.Rd b/man/plot_contributions.Rd
new file mode 100644
index 0000000..6f5637f
--- /dev/null
+++ b/man/plot_contributions.Rd
@@ -0,0 +1,69 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/plot-spaces.R
+\name{plot_contributions}
+\alias{plot_contributions}
+\title{Plot variable contributions to latent components}
+\usage{
+plot_contributions(
+ x,
+ r = NULL,
+ row_index = NULL,
+ column_index = NULL,
+ max_display = 200L,
+ limits = NULL,
+ row_label = "Row index",
+ column_label = "Column index",
+ main = "Estimated target contributions"
+)
+}
+\arguments{
+\item{x}{A complete estimated matrix, a \code{learner()} result, or a
+\code{dlearner()} result. No rotation of the singular vectors is performed.}
+
+\item{r}{Number of components. Required for a matrix; defaults to the stored
+rank for a fitting result. Numerically zero directions are rejected.}
+
+\item{row_index, column_index}{Optional distinct integer indices selecting
+variables to display, in display order. The same indices are used in all
+populations. SVD always uses the full matrices before subsetting.}
+
+\item{max_display}{Maximum number of displayed variables on each side when
+indices are omitted; default 200. Deterministic evenly spaced indices are
+used, without changing the random-number state.}
+
+\item{limits}{Optional two positive color limits, in row-side then column-side
+order. Each side uses a symmetric scale about zero shared across populations.
+Values beyond supplied limits are clipped for display and counts are returned.}
+
+\item{row_label, column_label}{Axis labels. Use Variant index and Phenotype
+index only when these accurately describe the input rows and columns.}
+
+\item{main}{Title describing the input estimate.}
+}
+\value{
+Invisibly, a list with full \code{row_scores} and \code{column_scores},
+displayed \code{row_index} and \code{column_index}, \code{limits}, and
+\code{clipped_count}. Score matrices have variables in rows and components
+in columns; all values are returned without display clipping.
+}
+\description{
+Show which rows and columns contribute to each latent component of a fitted
+target matrix, using squared singular-vector entries. Higher scores are darker red.
+}
+\details{
+An SVD of the estimated matrix supplies U and V. Row scores are U squared
+elementwise; column scores are V squared elementwise. Scores sum to one per
+component over all variables on the corresponding side. They are not squared
+entries of the raw optimization factors, effect sizes, or correlations.
+Subsetting is applied only for display and scores are not renormalized.
+The left panel displays columns and the right panel rows. Color limits are
+separate for the two sides and start at zero. Optional limits specify row
+then column upper limits, with clipped-entry counts reported.
+Signs do not affect squared scores, but rotations within repeated singular
+subspaces can change component-specific scores. Components from different
+populations are not automatically aligned or given biological labels.
+}
+\examples{
+fit <- dlearner(dat_multisource$Y_source, dat_multisource$Y_target, r = 3)
+plot_contributions(fit,max_display=50)
+}
diff --git a/man/plot_cv.Rd b/man/plot_cv.Rd
new file mode 100644
index 0000000..f665d1b
--- /dev/null
+++ b/man/plot_cv.Rd
@@ -0,0 +1,121 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/plotting.R
+\name{plot_cv}
+\alias{plot_cv}
+\title{Plot the LEARNER cross-validation scores}
+\usage{
+plot_cv(
+ x,
+ lambda_1_all = NULL,
+ lambda_2_all = NULL,
+ lambda_1_row_all = NULL,
+ lambda_1_col_all = NULL,
+ n_observed = NULL,
+ lambda_2_index = NULL,
+ log_labels = FALSE,
+ relative_limit = 0.3,
+ show_values = TRUE,
+ main = "Cross-validation error"
+)
+}
+\arguments{
+\item{x}{A result returned by \code{cv.learner()}.}
+
+\item{lambda_1_all}{Common-penalty candidates in the exact order used for CV.
+Required for a two-dimensional score array; omit for separate penalties.}
+
+\item{lambda_2_all}{Balance-penalty candidates in the exact order used for CV.}
+
+\item{lambda_1_row_all, lambda_1_col_all}{Separate-penalty candidates in the
+exact order used for CV. Both are required for a three-dimensional array.}
+
+\item{n_observed}{Optional positive integer giving the number of observed
+target entries used in CV. If supplied, scores are divided by this count and
+labeled Validation MSE. Otherwise the original sums are labeled Validation
+SSE. Missing target entries must not be included in this count.}
+
+\item{lambda_2_index}{Optional vector of distinct integer indices into
+\code{lambda_2_all}, selecting panels to draw in separate-penalty mode.
+\code{NULL} draws all panels. Not used in common-penalty mode.}
+
+\item{log_labels}{Logical scalar. If TRUE, display base-10 logarithms of
+the candidates on the heatmap axes.
+All candidates on those axes must be strictly positive. Log labels are
+rounded to two decimal places. Defaults to FALSE.}
+
+\item{relative_limit}{Optional positive relative increase cutoff, default 0.3
+(30 percent). Larger finite errors are gray and labeled above the cutoff,
+never as diverged. NULL keeps all finite values in the color scale.}
+
+\item{show_values}{Whether to print scores inside cells. Defaults to TRUE.}
+
+\item{main}{Plot title.}
+}
+\value{
+Invisibly, a list with \code{scores}, \code{scale}, \code{grids},
+\code{best_index}, \code{shown_lambda_2_index}, \code{color_limits}, and
+\code{unavailable_count}, \code{log_labels}, \code{score_range}, \code{minimum},
+\code{exact_minima}, \code{boundary}, and \code{lambda_2_slice_range}.
+Indices refer to the original candidate order. Boundary flags indicate a
+selected numerical endpoint on a grid with at least two distinct values.
+Also returns \code{display_data}, \code{display_limits}, \code{relative},
+\code{above_limit_count}, \code{plot} (the first ggplot page), and \code{plots}
+(all ggplot pages). The original absolute score diagnostics are unchanged.
+}
+\description{
+Draw a common-penalty heatmap or separate-penalty heatmaps, with a shared
+color scale and a black box at the first minimum in R array order.
+}
+\details{
+Existing CV results do not store the full candidate grids, so callers must
+provide them. Lengths and the selected penalty values are checked, but these
+checks cannot detect every incorrectly supplied grid. Use the original grids,
+without sorting or removing duplicates. Axes are equally spaced candidate
+positions labeled with actual values, or their base-10 logarithms when
+\code{log_labels=TRUE}. This changes labels only, not scores, candidate order,
+or cell spacing. The reversed magma palette uses pale colors for lower error
+and dark colors for higher error. Cell labels show increases over the minimum.
+
+Common mode shows \code{lambda_1_all} on the horizontal axis and
+\code{lambda_2_all} on the vertical axis. Separate mode shows the U-space
+penalty horizontally and the V-space penalty vertically, with one panel for
+each balance penalty. All panels share limits computed from the full score
+array, even when only some panels are selected. Up to six panels are drawn
+per page. Gray cells have non-finite scores and are excluded from selection.
+Only the first minimum is marked; a panel need not contain a black box.
+Captions report the score range, exact ties, and boundary selections. These are
+descriptive diagnostics, not significance or convergence tests.
+
+The stored \code{mse_all} values are sums of held-out squared errors. MSE
+conversion requires the correct \code{n_observed}; the function cannot infer
+this count from the scores. Neither this plot nor CV certifies convergence
+of the individual fits. No fitting is performed and \code{x} is unchanged.
+Requires optional ggplot2 >= 3.4.0; install it with
+\code{install.packages("ggplot2")}. Estimation and other plots do not require it.
+Draws with ggplot2 without modifying base graphics parameters. Each page shares
+one legend. A sufficiently large device is recommended for multiple panels.
+Relative increases are undefined when the minimum is zero; the display then
+uses absolute scores and ignores relative_limit. High relative error alone
+never establishes optimizer divergence. Non-finite scores are unavailable.
+Use ggplot2::ggsave(plot = result$plot, ...) to export the first page, or
+export each member of result$plots for multiple pages.
+}
+\examples{
+if (requireNamespace("ggplot2", quietly = TRUE) &&
+ utils::packageVersion("ggplot2") >= "3.4.0") {
+row_grid <- c(0.1, 1)
+col_grid <- c(0.1, 2)
+balance_grid <- c(0.1, 1)
+target <- dat_highsim$Y_target[1:12, 1:8]
+set.seed(42)
+cv <- cv.learner(dat_highsim$Y_source[1:12, 1:8], target,
+ r = 2, lambda_1_row_all = row_grid,
+ lambda_1_col_all = col_grid, lambda_2_all = balance_grid,
+ step_size = 0.0003, n_folds = 2,
+ control = list(max_iter = 20))
+plot_cv(cv, lambda_1_row_all = row_grid, lambda_1_col_all = col_grid,
+ lambda_2_all = balance_grid, n_observed = sum(!is.na(target)))
+
+}
+
+}
diff --git a/man/plot_objective.Rd b/man/plot_objective.Rd
new file mode 100644
index 0000000..7b7ce87
--- /dev/null
+++ b/man/plot_objective.Rd
@@ -0,0 +1,47 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/plotting.R
+\name{plot_objective}
+\alias{plot_objective}
+\title{Plot the objective function values of a LEARNER fit}
+\usage{
+plot_objective(x, iterations = "all", main = NULL)
+}
+\arguments{
+\item{x}{A result returned by \code{learner}.}
+
+\item{iterations}{Iterations to show: \code{"all"} (the default),
+\code{"near_optimum"}, or a vector of iteration numbers such as
+\code{1:1000}. \code{"near_optimum"} shows a window centered on the first
+smallest finite recorded objective value. The window contains 10\% of the
+recorded iterations (rounded up), with a minimum of 10, or all iterations
+for shorter histories. At either boundary the window shifts to stay within
+the recorded history. Both axes adjust to the selected window.}
+
+\item{main}{Title of the plot. By default, \code{"Objective function"},
+with \code{"(near optimum)"} added when \code{iterations = "near_optimum"}.}
+}
+\value{
+Invisibly, a data frame of the plotted iterations and objective
+values.
+}
+\description{
+Plots the values of the objective function at each iteration of the
+numerical optimization algorithm in \code{learner}.
+}
+\details{
+The red point marks the smallest objective value, which is also
+given above the plot. The objective value at the initialization is not
+included. With \code{iterations = "near_optimum"}, the plot shows a shaded
+band spanning each pair of consecutive objective values, so that values
+that alternate between iterations remain visible, and a line at their
+average.
+}
+\examples{
+fit <- learner(dat_highsim$Y_source, dat_highsim$Y_target,
+ r = 3, lambda_1 = 1, lambda_2 = 1, step_size = 0.003,
+ control = list(max_iter = 1000))
+plot_objective(fit)
+plot_objective(fit, iterations = "near_optimum")
+plot_objective(fit, iterations = 1:100)
+
+}
diff --git a/man/plot_source_spaces.Rd b/man/plot_source_spaces.Rd
new file mode 100644
index 0000000..7e341fc
--- /dev/null
+++ b/man/plot_source_spaces.Rd
@@ -0,0 +1,71 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/plot-spaces.R
+\name{plot_source_spaces}
+\alias{plot_source_spaces}
+\title{Compare target and source latent spaces}
+\usage{
+plot_source_spaces(
+ Y_target,
+ Y_source,
+ r,
+ row_index = NULL,
+ column_index = NULL,
+ max_display = 200L,
+ limits = NULL,
+ row_label = "Row index",
+ column_label = "Column index",
+ target_name = "Target population"
+)
+}
+\arguments{
+\item{Y_target}{Complete finite target matrix, or a completed estimate.
+Missing values are rejected; the function does not impute them.}
+
+\item{Y_source}{Aligned complete source matrix or nonempty list of matrices.}
+
+\item{r}{Positive integer giving the retained rank for every matrix.}
+
+\item{row_index, column_index}{Optional distinct integer indices selecting
+variables to display, in display order. The same indices are used in all
+populations. SVD always uses the full matrices before subsetting.}
+
+\item{max_display}{Maximum number of displayed variables on each side when
+indices are omitted; default 200. Deterministic evenly spaced indices are
+used, without changing the random-number state.}
+
+\item{limits}{Optional two positive color limits, in row-side then column-side
+order. Each side uses a symmetric scale about zero shared across populations.
+Values beyond supplied limits are clipped for display and counts are returned.}
+
+\item{row_label, column_label}{Axis labels. Use Variant index and Phenotype
+index only when these accurately describe the input rows and columns.}
+
+\item{target_name}{Panel label describing the supplied target matrix.}
+}
+\value{
+Invisibly, a list with \code{row_projectors}, \code{column_projectors},
+\code{row_index}, \code{column_index}, \code{limits}, and \code{clipped_count}.
+Each projector is the displayed submatrix, computed without allocating the
+full row-by-row projector. Its population name identifies the input.
+}
+\description{
+Plot column-space and row-space projection submatrices for the target and
+each source, to compare their retained latent spaces.
+}
+\details{
+For an SVD basis B, the displayed matrix is a submatrix of B B-transpose.
+Column-side panels appear on top, row-side panels below. Red denotes negative
+values and blue positive values. Color scales are shared within each side,
+including across pages, but can differ between the two sides. At most three
+populations appear per page. Source names are taken from the input list when
+supplied. No name-based alignment or matching of factors is performed.
+These are individual population projectors, not the weighted source operator.
+A submatrix need not itself be an idempotent projector.
+Numerically zero retained singular directions are rejected. With repeated
+singular values at the truncation boundary the retained space is not unique.
+The input order must already encode any desired phenotype or genomic ordering.
+}
+\examples{
+plot_source_spaces(dat_multisource$Y_target, dat_multisource$Y_source,
+ r = 3, max_display = 20)
+}
diff --git a/tests/testthat/helper-plotting.R b/tests/testthat/helper-plotting.R
new file mode 100644
index 0000000..2e9ca4f
--- /dev/null
+++ b/tests/testthat/helper-plotting.R
@@ -0,0 +1,15 @@
+# Devices are temporary and closed even when an assertion fails.
+plot_device <- function() {
+ path <- tempfile(fileext = '.pdf')
+ grDevices::pdf(path, width = 10, height = 5)
+ path
+}
+
+cv_plot_example <- function(separate = TRUE) {
+ if (separate) {
+ scores <- array(c(8, 7, 6, 5, 4, 3, 2, 1), c(2, 2, 2))
+ list(mse_all = scores, lambda_1_min = NA_real_, lambda_1_row_min = 0,
+ lambda_1_col_min = 3, lambda_2_min = 1, r = 2)
+ } else list(mse_all = matrix(c(8, 7, 6, 5), 2, 2),
+ lambda_1_min = 0, lambda_2_min = 1, r = 2)
+}
diff --git a/tests/testthat/test-cv-heatmap-style.R b/tests/testthat/test-cv-heatmap-style.R
new file mode 100644
index 0000000..1c9e0e8
--- /dev/null
+++ b/tests/testthat/test-cv-heatmap-style.R
@@ -0,0 +1,75 @@
+test_that('relative heatmaps preserve score alignment and distinguish high error from missing', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ x <- list(mse_all = matrix(c(100, 110, 140, NA_real_), 2), lambda_1_min = 1, lambda_2_min = 0.1)
+ args <- list(x = x, lambda_1_all = c(1, 2), lambda_2_all = c(.1, 1))
+ z <- do.call(plot_cv, args)
+ expect_s3_class(z$plot, 'ggplot')
+ expect_equal(z$display_data$score, c(100, 110, 140, NA_real_))
+ expect_equal(z$display_data$label, c('0.0%', '10.0%', '>30%', 'NA'))
+ expect_equal(z$display_data$selected, c(TRUE, FALSE, FALSE, FALSE))
+ expect_equal(z$above_limit_count, 1)
+ expect_equal(z$unavailable_count, 1)
+ expect_equal(z$display_limits, c(0, .1))
+ expect_match(z$plot$labels$caption, 'not a divergence diagnosis', fixed = TRUE)
+ mse <- do.call(plot_cv, c(args, list(n_observed = 10)))
+ expect_equal(mse$display_data$display, z$display_data$display)
+ expect_equal(mse$minimum, 10)
+ full <- do.call(plot_cv, c(args, list(relative_limit = NULL, show_values = FALSE)))
+ expect_equal(full$above_limit_count, 0)
+ expect_equal(full$display_limits, c(0, .4))
+ expect_equal(length(full$plot$layers), 2)
+ expect_error(do.call(plot_cv, c(args, list(relative_limit = 0))), 'relative_limit')
+ expect_error(do.call(plot_cv, c(args, list(show_values = NA))), 'show_values')
+})
+
+test_that('zero minima use absolute scores and exact ties remain first in CV order', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ x <- list(mse_all = matrix(c(0, 0, 2, 4), 2), lambda_1_min = 1, lambda_2_min = .1)
+ z <- plot_cv(x, lambda_1_all = c(1, 2), lambda_2_all = c(.1, 1))
+ expect_false(z$relative)
+ expect_equal(z$above_limit_count, 0)
+ expect_equal(z$display_data$display, c(0, 0, 2, 4))
+ expect_equal(z$exact_minima, 2)
+ expect_equal(which(z$display_data$selected), 1L)
+ expect_match(z$plot$labels$caption, 'absolute SSE', fixed = TRUE)
+ x$mse_all[] <- 0
+ expect_gt(diff(plot_cv(x, lambda_1_all = c(1, 2), lambda_2_all = c(.1, 1))$display_limits), 0)
+})
+
+test_that('separate heatmaps share scales and return printable pages', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ x <- list(mse_all = array(7:1, c(1, 1, 7)), lambda_1_row_min = 1,
+ lambda_1_col_min = 2, lambda_2_min = 7)
+ args <- list(x = x, lambda_1_row_all = 1, lambda_1_col_all = 2, lambda_2_all = 1:7,
+ relative_limit = NULL)
+ full <- do.call(plot_cv, args)
+ subset <- do.call(plot_cv, c(args, list(lambda_2_index = c(7, 2))))
+ expect_length(full$plots, 2)
+ expect_equal(full$display_limits, subset$display_limits)
+ expect_equal(as.character(subset$plot$data$panel), c('2', '7'))
+ expect_equal(levels(subset$plot$data$panel), c('7', '2'))
+ expect_equal(sum(subset$plot$data$selected), 1)
+})
+
+
+test_that('close or repeated balance candidates remain distinct heatmap panels', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ for (balance in list(c(1, 1.00001), c(1, 1))) {
+ x <- list(mse_all = array(1:8, c(2, 2, 2)), lambda_1_row_min = 1,
+ lambda_1_col_min = 1, lambda_2_min = 1)
+ result <- plot_cv(x, lambda_1_row_all = c(1, 2),
+ lambda_1_col_all = c(1, 2), lambda_2_all = balance)
+ built <- ggplot2::ggplot_build(result$plot)
+ expect_equal(nrow(built$layout$layout), 2L)
+ expect_equal(as.integer(table(built$data[[1]]$PANEL)), c(4L, 4L))
+ expect_equal(as.integer(built$data[[3]]$PANEL), 1L)
+ }
+})
diff --git a/tests/testthat/test-cv-plot-errors.R b/tests/testthat/test-cv-plot-errors.R
new file mode 100644
index 0000000..ec0239f
--- /dev/null
+++ b/tests/testthat/test-cv-plot-errors.R
@@ -0,0 +1,48 @@
+test_that('CV heatmaps reject invalid titles and separate arguments in common mode', {
+ args <- list(x = cv_plot_example(FALSE), lambda_1_all = c(10, 0),
+ lambda_2_all = c(0, 1))
+ for (title in list(NA_character_, character(), c('a', 'b'), 1)) {
+ expect_error(do.call(plot_cv, c(args, list(main = title))),
+ 'main must be one nonmissing character string', fixed = TRUE)
+ }
+ for (extra in list(list(lambda_1_row_all = c(1, 2)),
+ list(lambda_1_col_all = c(1, 2)),
+ list(lambda_2_index = 1))) {
+ expect_error(do.call(plot_cv, c(args, extra)),
+ 'Separate-penalty arguments require a 3D score array', fixed = TRUE)
+ }
+})
+
+test_that('CV heatmaps explain how to install the optional plotting dependency', {
+ # Older supported testthat versions do not provide local_mocked_bindings().
+ skip_if_not_installed('testthat', minimum_version = '3.2.0')
+ original_require <- base::requireNamespace
+ testthat::local_mocked_bindings(
+ requireNamespace = function(package, ...) {
+ if (identical(package, 'ggplot2')) return(FALSE)
+ original_require(package, ...)
+ }, .package = 'base')
+ expect_error(plot_cv(cv_plot_example(FALSE), lambda_1_all = c(10, 0),
+ lambda_2_all = c(0, 1)),
+ 'Install it with install.packages("ggplot2")', fixed = TRUE)
+})
+
+test_that('CV heatmaps explain when the plotting dependency needs updating', {
+ skip_if_not_installed('testthat', minimum_version = '3.2.0')
+ # Simulate availability/version checks without changing installed packages.
+ original_require <- base::requireNamespace
+ original_version <- utils::packageVersion
+ testthat::local_mocked_bindings(
+ requireNamespace = function(package, ...) {
+ if (identical(package, 'ggplot2')) return(TRUE)
+ original_require(package, ...)
+ }, .package = 'base')
+ testthat::local_mocked_bindings(
+ packageVersion = function(pkg, ...) {
+ if (identical(pkg, 'ggplot2')) return(package_version('3.3.6'))
+ original_version(pkg, ...)
+ }, .package = 'utils')
+ expect_error(plot_cv(cv_plot_example(FALSE), lambda_1_all = c(10, 0),
+ lambda_2_all = c(0, 1)),
+ 'requires ggplot2 >= 3.4.0. Update it', fixed = TRUE)
+})
diff --git a/tests/testthat/test-paper-plots.R b/tests/testthat/test-paper-plots.R
new file mode 100644
index 0000000..c4c6a53
--- /dev/null
+++ b/tests/testthat/test-paper-plots.R
@@ -0,0 +1,56 @@
+test_that('CV annotations distinguish small differences from exact ties and numeric boundaries', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path<-plot_device();on.exit({grDevices::dev.off();unlink(path)},add=TRUE)
+ x<-list(mse_all=matrix(c(1.06,1.051,1.07,1.06,1.052,1.07),3),
+ lambda_1_min=50,lambda_2_min=0.01)
+ d<-plot_cv(x,lambda_1_all=c(10,50,100),lambda_2_all=c(0.01,1))
+ expect_equal(d$exact_minima,1)
+ expect_equal(d$score_range,c(1.051,1.07))
+ expect_equal(d$lambda_2_slice_range,c(1.051,1.052))
+ expect_identical(d$boundary,c(common=FALSE,balance=TRUE))
+ x$mse_all[2,2]<-1.051
+ expect_equal(plot_cv(x,lambda_1_all=c(10,50,100),lambda_2_all=c(0.01,1))$exact_minima,2)
+ x$mse_all<-matrix(c(3,1,2),3);x$lambda_1_min<-1;x$lambda_2_min<-1
+ expect_identical(plot_cv(x,lambda_1_all=c(3,1,2),lambda_2_all=1)$boundary,
+ c(common=TRUE,balance=FALSE))
+})
+
+test_that('space plots compute full-matrix spaces before taking common submatrices', {
+ path<-plot_device();on.exit({grDevices::dev.off();unlink(path)},add=TRUE)
+ # This rank-one matrix has analytically known normalized left/right vectors.
+ u<-c(1,2,3);v<-c(2,1);y<-outer(u,v)
+ set.seed(1);before<-.Random.seed
+ old<-graphics::par(c('mfrow','mar'))
+ d<-plot_source_spaces(y,list(A=y,B=2*y),r=1,row_index=c(3,1),column_index=2:1)
+ expect_equal(d$row_projectors[[1]],outer(u,u)[c(3,1),c(3,1)]/sum(u^2))
+ expect_equal(d$column_projectors[[1]],outer(v,v)[2:1,2:1]/sum(v^2))
+ expect_equal(d$row_projectors[[2]],d$row_projectors[[3]])
+ expect_identical(.Random.seed,before)
+ expect_equal(graphics::par(c('mfrow','mar')),old)
+ expect_equal(d$limits,c(row=9/14,column=4/5))
+ clipped<-plot_source_spaces(y,y,r=1,limits=c(.1,.1))
+ expect_gt(sum(unlist(clipped$clipped_count)),0)
+ expect_equal(clipped$row_projectors[[1]],outer(u,u)/sum(u^2))
+ expect_error(plot_source_spaces(y,y,r=2),'zero singular')
+ expect_error(plot_source_spaces(y,y,r=1,row_index=c(1,1)),'distinct')
+ expect_error(plot_source_spaces(y,y,r=1,limits=c(0,.1)),'positive')
+ y[1]<-NA
+ expect_error(plot_source_spaces(y,y,r=1),'complete finite')
+})
+
+test_that('contributions use squared normalized singular vectors, not raw factors', {
+ path<-plot_device();on.exit({grDevices::dev.off();unlink(path)},add=TRUE)
+ u<-c(1,2,3);v<-c(2,1);y<-outer(u,v)
+ d<-plot_contributions(y,r=1,row_index=c(3,1),column_index=2)
+ expect_equal(drop(d$row_scores),u^2/sum(u^2))
+ expect_equal(drop(d$column_scores),v^2/sum(v^2))
+ expect_equal(colSums(d$row_scores),1)
+ expect_equal(colSums(d$column_scores),1)
+ expect_equal(d$row_index,c(3L,1L))
+ expect_equal(plot_contributions(-y,r=1)$row_scores,d$row_scores)
+ expect_equal(plot_contributions(list(learner_estimate=y,r=1))$column_scores,d$column_scores)
+ expect_equal(plot_contributions(list(dlearner_estimate=y,r=1))$row_scores,d$row_scores)
+ expect_error(plot_contributions(y),'Supply r')
+ expect_error(plot_contributions(y,r=2),'zero singular')
+ expect_error(plot_contributions(y,r=1,column_index=3),'valid integer')
+})
diff --git a/tests/testthat/test-plotting.R b/tests/testthat/test-plotting.R
new file mode 100644
index 0000000..7101c73
--- /dev/null
+++ b/tests/testthat/test-plotting.R
@@ -0,0 +1,204 @@
+test_that('objective plots show the requested iterations', {
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ # Decrease to the minimum at iteration 31, then values that stay above it.
+ fit <- list(objective_values = c(seq(100, 10, length.out = 30), rep(c(9, 9.5), 50)),
+ convergence_criterion = 2L)
+ original <- fit
+ result <- withVisible(plot_objective(fit))
+ expect_false(result$visible)
+ expect_equal(result$value$iteration, 1:130)
+ expect_equal(result$value$objective, fit$objective_values)
+ expect_equal(plot_objective(fit, iterations = 2:4)$iteration, 2:4)
+ expect_equal(plot_objective(fit, iterations = c(5, 2, 2))$iteration, c(2, 5))
+ # A 13-iteration window centered on the first minimum at iteration 31.
+ expect_equal(plot_objective(fit, iterations = 'near_optimum')$iteration, 25:37)
+ expect_identical(fit, original)
+})
+
+test_that('near-optimum windows include both sides and zoom the horizontal axis', {
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ # Interior minimum followed by a long increase.
+ fit <- list(objective_values = (seq_len(101) - 51)^2)
+ plot_objective(fit)
+ full <- graphics::par('usr')[1:2]
+ result <- plot_objective(fit, iterations = 'near_optimum')
+ expect_equal(result$iteration, 46:56)
+ expect_gt(graphics::par('usr')[1], full[1])
+ expect_lt(graphics::par('usr')[2], full[2])
+ # At the beginning, shift the window right; exact ties use the first minimum.
+ fit$objective_values <- seq_len(100)
+ expect_equal(plot_objective(fit, iterations = 'near_optimum')$iteration, 1:10)
+ fit$objective_values <- rep(1, 100)
+ expect_equal(plot_objective(fit, iterations = 'near_optimum')$iteration, 1:10)
+})
+
+test_that('near-optimum plots show the final stretch when the minimum is at the end', {
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ # A history that is still decreasing at the last iteration.
+ fit <- list(objective_values = seq(100, 10, length.out = 100), convergence_criterion = 1L)
+ expect_equal(plot_objective(fit, iterations = 'near_optimum')$iteration, 91:100)
+ # Short histories are shown in full.
+ fit <- list(objective_values = c(9, 5, 4), convergence_criterion = 1L)
+ expect_equal(plot_objective(fit, iterations = 'near_optimum')$iteration, 1:3)
+})
+
+test_that('near-optimum plots handle alternating and missing values and keep graphics settings', {
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ # Alternating values after a decrease; near_optimum draws them as a band.
+ values <- c(seq(50, 11, length.out = 40), 10 + rep(c(0, 0.5), 200))
+ values[100] <- NA
+ fit <- list(objective_values = values, convergence_criterion = 2L)
+ old <- graphics::par('mar')
+ result <- plot_objective(fit)
+ expect_equal(result$iteration, seq_along(values))
+ expect_equal(plot_objective(fit, iterations = 'near_optimum')$iteration, 20:63)
+ expect_equal(graphics::par('mar'), old)
+})
+
+test_that('objective plots handle flat, singleton, and non-finite histories', {
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ for (values in list(4, c(4, 4, 4), c(4, Inf, NaN), c(Inf, 3, 1, NA))) {
+ fit <- list(objective_values = values, convergence_criterion = 2L)
+ expect_equal(plot_objective(fit)$objective, values)
+ expect_true(length(plot_objective(fit, iterations = 'near_optimum')$iteration) >= 1L)
+ }
+ fit <- list(objective_values = c(2, 1, 3), convergence_criterion = 2L)
+ for (bad in list(0, 4, 1.5, NA, numeric(), 'near', TRUE)) {
+ expect_error(plot_objective(fit, iterations = bad), 'iterations')
+ }
+ expect_error(plot_objective(list(objective_values = c(Inf, NA))), 'learner')
+ expect_error(plot_objective(list()), 'learner')
+})
+
+test_that('CV plots retain candidate order, mark the global minimum, and scale SSE', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ x <- cv_plot_example()
+ original <- x
+ old <- graphics::par(c('mfrow','mar'))
+ args <- list(x = x, lambda_1_row_all = c(10, 0),
+ lambda_1_col_all = c(0.2, 3), lambda_2_all = c(0, 1))
+ sse <- do.call(plot_cv, args)
+ mse <- do.call(plot_cv, c(args, list(n_observed = 40)))
+ expect_identical(sse$grids$row, c(10, 0))
+ expect_identical(sse$best_index, c(2L, 2L, 2L))
+ expect_equal(sse$scores, x$mse_all)
+ expect_equal(mse$scores, x$mse_all / 40)
+ expect_equal(sse$scale, 'Validation SSE')
+ expect_equal(mse$scale, 'Validation MSE')
+ expect_equal(mse$color_limits, range(x$mse_all) / 40)
+ expect_identical(x, original)
+ expect_equal(graphics::par(c('mfrow','mar')), old)
+ subset <- do.call(plot_cv, c(args, list(lambda_2_index = 1)))
+ expect_equal(subset$color_limits, sse$color_limits)
+ expect_equal(subset$shown_lambda_2_index, 1)
+ expect_equal(subset$best_index, sse$best_index)
+ common <- plot_cv(cv_plot_example(FALSE), lambda_1_all = c(10, 0), lambda_2_all = c(0, 1))
+ expect_identical(common$best_index, c(2L, 2L))
+})
+
+test_that('CV plots handle singleton axes, duplicate candidates, ties and unavailable scores', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ x <- list(mse_all = array(c(5, 5, Inf, NA), c(1, 2, 2)),
+ lambda_1_row_min = 0, lambda_1_col_min = 2, lambda_2_min = 1)
+ z <- plot_cv(x, lambda_1_row_all = 0, lambda_1_col_all = c(2, 2),
+ lambda_2_all = c(1, 10))
+ expect_identical(z$best_index, c(1L, 1L, 1L))
+ expect_equal(z$unavailable_count, 2)
+ expect_gt(diff(z$color_limits), 0)
+ x <- list(mse_all = matrix(0, 1, 1), lambda_1_min = 0, lambda_2_min = 0)
+ z <- plot_cv(x, lambda_1_all = 0, lambda_2_all = 0)
+ expect_identical(z$best_index, c(1L, 1L))
+ expect_gt(diff(z$color_limits), 0)
+ # Seven balance slices exercise pagination without inventing new candidates.
+ x <- list(mse_all = array(7:1, c(1, 1, 7)), lambda_1_row_min = 0,
+ lambda_1_col_min = 2, lambda_2_min = 7)
+ z <- plot_cv(x, lambda_1_row_all = 0, lambda_1_col_all = 2, lambda_2_all = 1:7)
+ expect_equal(z$shown_lambda_2_index, 1:7)
+ expect_identical(z$best_index, c(1L, 1L, 7L))
+})
+
+test_that('CV plots reject missing, mismatched, and ambiguous grid information', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ x <- cv_plot_example()
+ args <- list(x = x, lambda_1_row_all = c(10, 0),
+ lambda_1_col_all = c(0.2, 3), lambda_2_all = c(0, 1))
+ expect_error(plot_cv(x), 'lambda_1_row_all')
+ bad <- args; bad$lambda_1_row_all <- c(0,10)
+ expect_error(do.call(plot_cv,bad), 'selected penalties')
+ bad <- args; bad$lambda_1_row_all <- 0
+ expect_error(do.call(plot_cv,bad), 'score-array dimension')
+ expect_error(do.call(plot_cv,c(args,list(lambda_1_all=1))), 'Do not supply')
+ expect_error(do.call(plot_cv,c(args,list(n_observed=0))), 'n_observed')
+ expect_error(do.call(plot_cv,c(args,list(lambda_2_index=c(1,1)))), 'distinct')
+ expect_error(do.call(plot_cv,c(args,list(lambda_2_index=3))), 'distinct')
+ expect_error(do.call(plot_cv,c(args,list(lambda_2_index=NA_real_))), 'distinct')
+ x$mse_all[] <- Inf
+ expect_error(plot_cv(x), 'No finite')
+ x$mse_all[] <- -1
+ expect_error(plot_cv(x), 'cannot be negative')
+ expect_error(plot_cv(list(mse_all=1:3)), 'two- or three-dimensional')
+})
+
+test_that('plot functions accept actual fits and both CV modes without changing outputs', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ set.seed(623)
+ s <- matrix(rnorm(48), 8, 6)
+ t <- matrix(rnorm(48), 8, 6); t[c(1,3)] <- NA_real_
+ fit <- learner(s,t,r=2,lambda_1=1,lambda_2=1,step_size=0.001,
+ control=list(max_iter=6))
+ expect_equal(plot_objective(fit)$objective,fit$objective_values)
+ for (separate in c(FALSE,TRUE)) {
+ args <- list(Y_source=list(s,s*0.9),Y_target=t,r=2,step_size=0.001,
+ lambda_2_all=c(0,1),n_folds=2,control=list(max_iter=6))
+ grid <- if(separate) list(lambda_1_row_all=c(0,1),lambda_1_col_all=c(0.5,2)) else list(lambda_1_all=c(0,1))
+ cv <- do.call(cv.learner,c(args,grid))
+ original <- cv
+ result <- do.call(plot_cv,c(list(x=cv,lambda_2_all=c(0,1),n_observed=sum(!is.na(t))),grid))
+ expect_equal(result$scores,cv$mse_all/46)
+ expect_identical(cv,original)
+ }
+})
+
+
+test_that('paper-style log labels preserve scores and candidate selection', {
+ testthat::skip_if_not_installed("ggplot2", minimum_version = "3.4.0")
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ x <- list(mse_all = matrix(c(4, 3, 2, 1), 2, 2),
+ lambda_1_min = 100, lambda_2_min = 1)
+ raw <- plot_cv(x, lambda_1_all = c(10, 100), lambda_2_all = c(0.01, 1))
+ log <- plot_cv(x, lambda_1_all = c(10, 100), lambda_2_all = c(0.01, 1), log_labels = TRUE)
+ expect_equal(log$scores, raw$scores)
+ expect_equal(log$best_index, raw$best_index)
+ expect_equal(log$grids, raw$grids)
+ expect_true(log$log_labels)
+ expect_error(plot_cv(x, lambda_1_all = c(0, 100), lambda_2_all = c(0.01, 1),
+ log_labels = TRUE), 'strictly positive')
+ expect_error(plot_cv(x, log_labels = NA), 'TRUE or FALSE')
+ expect_error(plot_cv(x, log_labels = 'yes'), 'TRUE or FALSE')
+})
+
+test_that('objective plots reject ranges with no finite values before drawing', {
+ path <- plot_device()
+ on.exit({grDevices::dev.off(); unlink(path)}, add = TRUE)
+ fit <- list(objective_values = c(3, 1, NA_real_, Inf, NaN))
+ old <- graphics::par(c('mar', 'usr'))
+ expect_error(
+ withCallingHandlers(plot_objective(fit, iterations = 3:5),
+ warning = function(w) stop('Unexpected warning: ', conditionMessage(w))),
+ 'selected iterations contain no finite objective values')
+ expect_equal(graphics::par(c('mar', 'usr')), old)
+ expect_equal(plot_objective(fit, iterations = 2:5)$objective,
+ fit$objective_values[2:5])
+})
+
+