From 19296fdd7ef5262d939459de4aa9b4b794b20e48 Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Thu, 7 May 2026 11:09:35 +0200 Subject: [PATCH 1/9] Add splitPop and tests --- NAMESPACE | 1 + R/mergePops.R | 184 ++++++++++++++++++++++++++++++++ man/dot-splitAtLevels.Rd | 24 +++++ man/dot-splitPop.Rd | 21 ++++ man/splitPop.Rd | 67 ++++++++++++ tests/testthat/test-mergePops.R | 143 +++++++++++++++++++++++++ 6 files changed, 440 insertions(+) create mode 100644 man/dot-splitAtLevels.Rd create mode 100644 man/dot-splitPop.Rd create mode 100644 man/splitPop.Rd diff --git a/NAMESPACE b/NAMESPACE index 07534a87..5fd04848 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -112,6 +112,7 @@ export(solveRRBLUP_EM) export(solveRRBLUP_EM2) export(solveRRBLUP_EM3) export(solveUVM) +export(splitPop) export(unnameMultiPop) export(usefulness) export(varA) diff --git a/R/mergePops.R b/R/mergePops.R index c994fbd3..ecccad14 100644 --- a/R/mergePops.R +++ b/R/mergePops.R @@ -418,3 +418,187 @@ mergeMultiPops = function(..., level=0){ flatMultiPop = flattenMultiPop(multiPop) return(mergePops(flatMultiPop)) } + +#' @title Split Pop or MultiPop +#' +#' @description +#' Split a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object +#' into a \code{MultiPop} at specified nesting \code{level}(s) using one +#' or more grouping specifications (\code{by}). Grouping specs can be +#' atomic vectors (length \code{nInd(x)} or 1) or functions that return +#' such vectors. If \code{by} is a list, each element is applied +#' recursively to create nested \code{MultiPop} objects. +#' +#' @param x a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object +#' @param by a vector, function, or list of vectors/functions defining +#' groupings. Functions are called with a \code{Pop} and must return an +#' atomic vector of length \code{nInd(pop)} (or 1). +#' @param level A positive integer, a vector of positive integers, or +#' \code{Inf}. Only relevant when \code{x} is a \code{\link{MultiPop-class}} +#' object. If \code{level = Inf}, split any \code{\link{Pop-class}} within +#' \code{x} according to \code{by}. If \code{level = c(a, b)}, only +#' \code{Pop-class} objects at levels \code{a} and \code{b} (top \code{level=1}) +#' are split. An error is raised if any requested level is deeper than the +#' object's maximum nesting depth. +#' +#' @return Returns a \code{\link{MultiPop-class}} object. +#' +#' @examples +#' # Create founder haplotypes +#' founderPop = quickHaplo(nInd = 10, nChr = 1, segSites = 10) +#' +#' # Set simulation parameters +#' SP = SimParam$new(founderPop) +#' \dontshow{SP$nThreads = 1L} +#' SP$addTraitA(10) +#' +#' # Create population +#' pop = newPop(founderPop, simParam = SP) +#' +#' # Split pop into groups A/B deterministically +#' mp1 = splitPop(pop, by = rep(c("A", "B"), length.out = nInd(pop))) +#' mp1 +#' +#' # Nested split: First using a random grouping vector (three groups), then by family +#' mp2 = splitPop( +#' pop, +#' by = list( +#' sample(LETTERS[1:3], nInd(pop), replace = TRUE), +#' function(x) paste(x@mother, x@father, sep = "_") +#' ) +#' ) +#' mp2 +#' +#' # When x is a MultiPop, control which nesting levels are split +#' mp_nested = newMultiPop(pop[1:3], newMultiPop(pop[4:6], pop[7:9])) +#' splitPop(mp_nested, +#' by = function(p) rep(c("A","B"), length.out = nInd(p)), +#' level = 1) +#' +#' @export +splitPop = function(x, by, level = Inf) { + if (!isPop(x) && !isMultiPop(x)) { + stop("`x` must be a Pop or MultiPop object") + } + + # Validate by argument + if (!is.list(by)) { + by = list(by) + } + if (length(by) == 0) { + stop("`by` must have at least one grouping spec.") + } + + # Get max depth of nesting in MultiPop + md = ifelse(isMultiPop(x), .depthMultiPop(x), 1L) + + # Validate level argument + if (length(level) == 1L) { + if (is.infinite(level)) { + levels = Inf + } else { + if (!is.numeric(level) || is.na(level) || + level < 1 || level != as.integer(level)) { + stop("`level` must be a positive integer or Inf") + } + if (level > md) { + stop(sprintf("requested `level` exceeds max depth of `x` (%d)", md)) + } + levels = as.integer(level) + } + } else { + if (!is.numeric(level) || any(is.na(level)) || any(level <= 0)) { + stop("`level` must be a numeric vector of positive integers (no NA)") + } + if (any(is.infinite(level))) { + stop("cannot mix Inf with integer levels") + } + if (any(level != trunc(level))) { + stop("`level` must be a numeric vector of positive integers (no NA)") + } + if (any(level > md)) { + stop(sprintf("requested level(s) exceed max depth of `x` (%d)", md)) + } + levels = as.integer(unique(level)) + } + + # Recursively split at requested levels + res = .splitAtLevels(x, currentLevel = 1L, levels = levels, by = by) + return(res) +} + +#' Helper to apply splitting only at requested nesting levels +#' +#' @param obj A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object. +#' @param currentLevel Integer; current depth during recursion (top-level = 1). +#' @param levels Integer vector or \code{Inf}; levels at which to apply splits. +#' @param by A list of grouping specs (vectors or functions). +#' +#' @return The input object with splits applied at the requested levels. +#' +#' @keywords internal +.splitAtLevels = function(obj, currentLevel, levels, by) { + if (isPop(obj)) { + # A Pop at the currentLevel: split only if currentLevel requested (or Inf) + if (any(is.infinite(levels)) || currentLevel %in% levels) { + return(.splitPop(obj, by)) + } else { + return(obj) + } + } + + if (isMultiPop(obj)) { + # For each child: if child is MultiPop, its children are one level deeper; + # if child is Pop, it sits at currentLevel. + obj@pops = lapply(obj@pops, function(child) { + if (isMultiPop(child)) { + .splitAtLevels(child, currentLevel = currentLevel + 1L, levels = levels, by = by) + } else { + .splitAtLevels(child, currentLevel = currentLevel, levels = levels, by = by) + } + }) + validObject(obj) + return(obj) + } +} + +#' Helper function to recursively split a Pop object +#' +#' @param pop A \code{\link{Pop-class}} object. +#' @param by A list of grouping specs (vectors or functions). Each element is +#' applied in order to produce nested splits. +#' +#' @return A \code{\link{MultiPop-class}} produced by recursively applying \code{by}. +#' +#' @keywords internal +.splitPop = function(pop, by) { + if (length(by) == 0) { + return(pop) + } + + f = by[[1]] + groups = if (is.function(f)) f(pop) else f + + if (!is.atomic(groups)) { + stop("Grouping spec must be an atomic vector or a function returning one.") + } + if (any(is.na(groups))) { + stop("Grouping vector contains NA values.") + } + + n = nInd(pop) + if (length(groups) == 1L) { + groups = rep(groups, n) + } + if (length(groups) != n) { + stop( + "Grouping vector length (", length(groups), + ") must equal `nInd(pop)` (", n, "), or be length 1." + ) + } + + popList = split(pop, groups) + mp = newEmptyMultiPop() + mp@pops = lapply(popList, .splitPop, by = by[-1]) + return(mp) +} \ No newline at end of file diff --git a/man/dot-splitAtLevels.Rd b/man/dot-splitAtLevels.Rd new file mode 100644 index 00000000..616238e2 --- /dev/null +++ b/man/dot-splitAtLevels.Rd @@ -0,0 +1,24 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/mergePops.R +\name{.splitAtLevels} +\alias{.splitAtLevels} +\title{Helper to apply splitting only at requested nesting levels} +\usage{ +.splitAtLevels(obj, currentLevel, levels, by) +} +\arguments{ +\item{obj}{A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object.} + +\item{currentLevel}{Integer; current depth during recursion (top-level = 1).} + +\item{levels}{Integer vector or \code{Inf}; levels at which to apply splits.} + +\item{by}{A list of grouping specs (vectors or functions).} +} +\value{ +The input object with splits applied at the requested levels. +} +\description{ +Helper to apply splitting only at requested nesting levels +} +\keyword{internal} diff --git a/man/dot-splitPop.Rd b/man/dot-splitPop.Rd new file mode 100644 index 00000000..17e5520f --- /dev/null +++ b/man/dot-splitPop.Rd @@ -0,0 +1,21 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/mergePops.R +\name{.splitPop} +\alias{.splitPop} +\title{Helper function to recursively split a Pop object} +\usage{ +.splitPop(pop, by) +} +\arguments{ +\item{pop}{A \code{\link{Pop-class}} object.} + +\item{by}{A list of grouping specs (vectors or functions). Each element is +applied in order to produce nested splits.} +} +\value{ +A \code{\link{MultiPop-class}} produced by recursively applying \code{by}. +} +\description{ +Helper function to recursively split a Pop object +} +\keyword{internal} diff --git a/man/splitPop.Rd b/man/splitPop.Rd new file mode 100644 index 00000000..0fa9d21d --- /dev/null +++ b/man/splitPop.Rd @@ -0,0 +1,67 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/mergePops.R +\name{splitPop} +\alias{splitPop} +\title{Split Pop or MultiPop} +\usage{ +splitPop(x, by, level = Inf) +} +\arguments{ +\item{x}{a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object} + +\item{by}{a vector, function, or list of vectors/functions defining +groupings. Functions are called with a \code{Pop} and must return an +atomic vector of length \code{nInd(pop)} (or 1).} + +\item{level}{A positive integer, a vector of positive integers, or +\code{Inf}. Only relevant when \code{x} is a \code{\link{MultiPop-class}} +object. If \code{level = Inf}, split any \code{\link{Pop-class}} within +\code{x} according to \code{by}. If \code{level = c(a, b)}, only +\code{Pop-class} objects at levels \code{a} and \code{b} (top \code{level=1}) +are split. An error is raised if any requested level is deeper than the +object's maximum nesting depth.} +} +\value{ +Returns a \code{\link{MultiPop-class}} object. +} +\description{ +Split a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object +into a \code{MultiPop} at specified nesting \code{level}(s) using one +or more grouping specifications (\code{by}). Grouping specs can be +atomic vectors (length \code{nInd(x)} or 1) or functions that return +such vectors. If \code{by} is a list, each element is applied +recursively to create nested \code{MultiPop} objects. +} +\examples{ +# Create founder haplotypes +founderPop = quickHaplo(nInd = 10, nChr = 1, segSites = 10) + +# Set simulation parameters +SP = SimParam$new(founderPop) +\dontshow{SP$nThreads = 1L} +SP$addTraitA(10) + +# Create population +pop = newPop(founderPop, simParam = SP) + +# Split pop into groups A/B deterministically +mp1 = splitPop(pop, by = rep(c("A", "B"), length.out = nInd(pop))) +mp1 + +# Nested split: First using a random grouping vector (three groups), then by family +mp2 = splitPop( + pop, + by = list( + sample(LETTERS[1:3], nInd(pop), replace = TRUE), + function(x) paste(x@mother, x@father, sep = "_") + ) +) +mp2 + +# When x is a MultiPop, control which nesting levels are split +mp_nested = newMultiPop(pop[1:3], newMultiPop(pop[4:6], pop[7:9])) +splitPop(mp_nested, + by = function(p) rep(c("A","B"), length.out = nInd(p)), + level = 1) + +} diff --git a/tests/testthat/test-mergePops.R b/tests/testthat/test-mergePops.R index 21d420d0..e92eab31 100644 --- a/tests/testthat/test-mergePops.R +++ b/tests/testthat/test-mergePops.R @@ -191,3 +191,146 @@ test_that("cMultiPop_mergeMultiPops_and_flattenMultiPop", { newMultiPop(pop[8:9], pop[10:11], pop[12])) expect_equal(flattenMultiPop(mp3, level = 3), mp3) }) + +test_that("splitPop", { + # Create founder haplotypes and set simulation parameters + founderPop = quickHaplo(nInd = 12, nChr = 1, segSites = 10) + SP = SimParam$new(founderPop) + SP$nThreads = 1L + SP$addTraitA(10) + # Create population + pop = newPop(founderPop, simParam = SP) + + # splitPop basic behavior + by1 = sample(LETTERS[1:3], nInd(pop), replace = TRUE) + mp1 = splitPop(pop, by = by1) + + expect_true(isMultiPop(mp1)) + expect_length(mp1@pops, length(unique(by1))) + expect_identical(.depthMultiPop(mp1), 1L) + + expect_identical( + splitPop(pop, by = "A"), + splitPop(pop, by = rep("A", length(pop))) + ) + + # splitPop recursive behavior + mp2 = splitPop( + pop, + by = list( + sample(LETTERS[1:2], nInd(pop), replace = TRUE), + function(x) getFam(x, famType = "B") + ) + ) + + expect_true(isMultiPop(mp2)) + expect_true(all(vapply(mp2@pops, isMultiPop, logical(1)))) + expect_identical(.depthMultiPop(mp2), 2L) + + # Error handling + expect_error( + splitPop(newEmptyPop(ploidy = 2L, simParam = SP), by = list()), + "`by` must have at least one grouping spec.", + fixed = TRUE + ) + + expect_error( + splitPop(double(), by = 1:5), + "`x` must be a Pop or MultiPop object", + fixed = TRUE + ) + + expect_error( + splitPop(pop, by = numeric(0)), + paste0("Grouping vector length (", length(numeric(0)), + ") must equal `nInd(pop)` (", length(pop), "), or be length 1."), + fixed = TRUE + ) + + expect_error( + splitPop(pop, by = list(NA_real_)), + "Grouping vector contains NA values.", + fixed = TRUE + ) + + expect_error( + splitPop(pop, by = list(1:4)), + paste0("Grouping vector length (", length(list(1:4)[[1]]), + ") must equal `nInd(pop)` (", length(pop), "), or be length 1."), + fixed = TRUE + ) + + expect_error( + splitPop(pop, by = list(list(1, 2, 3))), + "Grouping spec must be an atomic vector or a function returning one.", + fixed = TRUE + ) + + # splitPop works for MultiPop objects + mp3 = splitPop(mp2, by = function(x) sample(LETTERS[1:2], length(x), replace = TRUE)) + + expect_true(is(mp1, "MultiPop")) + expect_true(all(vapply(mp3@pops, isMultiPop, logical(1)))) + expect_identical(.depthMultiPop(mp3), 3L) + + # Level control tests (deterministic splitter) + mp_nested = newMultiPop(pop[1:3], newMultiPop(pop[4:6], pop[7:9])) + splitter = function(p) rep(c("A", "B"), length.out = nInd(p)) + + # level = 1: only top-level Pop nodes are split + res1 = splitPop(mp_nested, by = splitter, level = 1) + expect_identical(res1[[1]], splitPop(mp_nested[[1]], by = splitter)) + expect_identical(res1[[2]], mp_nested[[2]]) + + # level = 2: only Pop nodes at depth 2 are split + res2 = splitPop(mp_nested, by = splitter, level = 2) + expect_identical(res2[[1]], mp_nested[[1]]) + expect_identical(res2[[2]], splitPop(mp_nested[[2]], by = splitter)) + + # level = 1:2: Pop nodes at depth 1 and 2 are split + res12 = splitPop(mp_nested, by = splitter, level = 1:2) + expect_identical(res12[[1]], splitPop(mp_nested[[1]], by = splitter)) + expect_identical(res12[[2]], splitPop(mp_nested[[2]], by = splitter)) + + # level = Inf: split every Pop node encountered + resInf = splitPop(mp_nested, by = splitter, level = Inf) + expect_identical(resInf[[1]], splitPop(mp_nested[[1]], by = splitter)) + expect_identical(resInf[[2]], splitPop(mp_nested[[2]], by = splitter)) + + # Error handling for level argument + expect_error( + splitPop(mp_nested, by = splitter, level = -1L), + "`level` must be a positive integer or Inf", + fixed = TRUE + ) + + expect_error( + splitPop(mp_nested, by = splitter, level = 4L), + "requested `level` exceeds max depth of `x` (2)", + fixed = TRUE + ) + + expect_error( + splitPop(mp_nested, by = splitter, level = c(1, NA)), + "`level` must be a numeric vector of positive integers (no NA)", + fixed = TRUE + ) + + expect_error( + splitPop(mp_nested, by = splitter, level = c(1, Inf)), + "cannot mix Inf with integer levels", + fixed = TRUE + ) + + expect_error( + splitPop(mp_nested, by = splitter, level = c(1, 1.5)), + "`level` must be a numeric vector of positive integers (no NA)", + fixed = TRUE + ) + + expect_error( + splitPop(mp_nested, by = splitter, level = c(1, 5)), + "requested level(s) exceed max depth of `x` (2)", + fixed = TRUE + ) +}) From 10583dfe4697f12a4cb1eb4b21c414a0bf8066bf Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Thu, 4 Jun 2026 14:56:00 +0200 Subject: [PATCH 2/9] Update `flattenMultiPop` to preserve names --- R/Class-Pop.R | 4 + R/mergePops.R | 148 +++++++++++++++++++++++++++----- man/dot-flattenMultiPop.Rd | 2 +- man/flattenMultiPop.Rd | 48 +++++++++-- man/unnameMultiPop.Rd | 1 + tests/testthat/test-mergePops.R | 39 ++++++++- 6 files changed, 209 insertions(+), 33 deletions(-) diff --git a/R/Class-Pop.R b/R/Class-Pop.R index ebf398c3..a27c60c4 100644 --- a/R/Class-Pop.R +++ b/R/Class-Pop.R @@ -1044,6 +1044,7 @@ setValidity("MultiPop",function(object){ errors = character() # Check that all populations are valid for(i in seq_len(length(object@pops))){ + # TODO: Validate names before returning multiPops if(!validObject(object@pops[[i]]) & (is(object@pops[[i]], "Pop") | is(object@pops[[i]],"MultiPop"))){ @@ -1433,6 +1434,7 @@ newEmptyMultiPop = function(){ #' #' # Set simulation parameters #' SP = SimParam$new(founderPop) +#' \dontshow{SP$nThreads = 1L} #' #' # Create population #' pop = newPop(founderPop, simParam = SP) @@ -1458,6 +1460,8 @@ newEmptyMultiPop = function(){ #' #' @export unnameMultiPop = function(x, level = Inf) { + stopifnot(isMultiPop(x)) + # Get max depth of nesting in MultiPop md = .depthMultiPop(x) diff --git a/R/mergePops.R b/R/mergePops.R index ecccad14..815ec241 100644 --- a/R/mergePops.R +++ b/R/mergePops.R @@ -215,6 +215,15 @@ mergePops = function(popList){ #' @param x A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object. #' @param level Integer scalar >= 1. Number of \code{MultiPop} levels to #' preserve. +#' @param preserveNames Character scalar controlling how to handle names when +#' flattening. One of: \cr +#' - \code{"auto"} (default): Keep names only if all are present and +#' unique after flattening; otherwise drop them. \cr +#' - \code{"concatenate"}: Prefix child names with parent names using +#' underscore separator. Keep names only if resulting names are unique. \cr +#' - \code{"force"}: Always keep concatenated names, making them unique +#' via \code{\link{make.unique}} if duplicates are found (with a warning). \cr +#' - \code{"none"}: Drop all names. #' #' @details #' The \code{level} argument controls how many levels of nesting are preserved. \cr @@ -226,11 +235,20 @@ mergePops = function(popList){ #' #' If \code{level} is greater than or equal to the nesting depth, the original #' object is returned unchanged. +#' +#' The \code{preserveNames} argument allows control over name preservation during +#' flattening. This is useful when you want to track the source of flattened +#' populations through their hierarchical names. The default \code{"auto"} mode +#' is conservative: names are kept only when they are meaningful (all present +#' and unique). Use \code{"concatenate"} to always build hierarchical names +#' (when possible), \code{"force"} to guarantee names with uniqueness enforcement, +#' or \code{"none"} to explicitly discard all names. #' #' @return If \code{x} is a \code{\link{Pop-class}}, the same \code{x} #' object is returned. Otherwise a \code{\link{MultiPop-class}} is returned #' whose \code{x@pops} slot contains \code{\link{Pop-class}} (and possibly -#' \code{\link{MultiPop-class}}) objects flattened according to \code{level}. +#' \code{\link{MultiPop-class}}) objects flattened according to \code{level} +#' and \code{preserveNames}. #' #' @seealso \code{\link{mergeMultiPops}} and \code{\link{mergePops}} #' @@ -246,32 +264,51 @@ mergePops = function(popList){ #' pop = newPop(founderPop, simParam=SP) #' #' # Create a multi-population with down to level 3 nesting -#' mp_nested = newMultiPop(pop[1:2], -#' newMultiPop(pop[3:4], -#' newMultiPop(pop[5:7], pop[8:12]))) +#' mp_nested = newMultiPop(pop1 = pop[1:2], +#' mp1 = newMultiPop(pop2 = pop[3:4], +#' mp2 = newMultiPop(pop3 = pop[5:7], pop4 = pop[8:12]))) #' mp_nested #' #' # Completely flatten to a single top-level MultiPop +#' # With default "auto" mode: names are kept if unique #' flattenMultiPop(mp_nested) #' -#' # Preserve two levels of nesting -#' flattenMultiPop(mp_nested, level=2) +#' # With "concatenate": build hierarchical names +#' flattenMultiPop(mp_nested, preserveNames = "concatenate") +#' +#' # With "force": guarantee unique names with suffixes +#' # (not needed here but useful if duplicates were present) +#' flattenMultiPop(mp_nested, preserveNames = "force") +#' +#' # With "none": drop all names +#' flattenMultiPop(mp_nested, preserveNames = "none") +#' +#' # Preserve two levels of nesting (level 2 names are kept) +#' flattenMultiPop(mp_nested, level = 2) #' #' @export -flattenMultiPop = function(x, level=1) { - if (isPop(x)) return(x) +flattenMultiPop = function(x, level = 1, + preserveNames = c("auto", "concatenate", "force", "none")) { + preserveNames = match.arg(preserveNames) + if (isPop(x)) { + return(x) + } stopifnot(isMultiPop(x)) - multi = which(sapply(x@pops, isMultiPop)) - while (level > 1) { - level = level - 1 - for (i in multi) { - x@pops[[i]] = flattenMultiPop(x@pops[[i]], level = level) + + if (level > 1L) { + for (i in seq_along(x@pops)) { + if (isMultiPop(x@pops[[i]])) { + x@pops[[i]] = flattenMultiPop( + x@pops[[i]], + level = level - 1L, + preserveNames = preserveNames + ) + } } - multiPop = do.call(newMultiPop, x@pops) - validObject(multiPop) - return(multiPop) + validObject(x) + return(x) } - flatPopList = .flattenMultiPop(x) + flatPopList = .flattenMultiPop(x, preserveNames = preserveNames) multiPop = do.call(newMultiPop, flatPopList) validObject(multiPop) return(multiPop) @@ -282,18 +319,85 @@ flattenMultiPop = function(x, level=1) { #' @param mp \code{\link{MultiPop-class}} object #' #' @keywords internal -.flattenMultiPop = function(mp) { +.flattenMultiPop = function(mp, preserveNames = c("auto", "concatenate", "force", "none")) { + preserveNames = match.arg(preserveNames) + + if (.depthMultiPop(mp) == 1L) { + res = mp@pops + if (preserveNames == "none") { + names(res) = NULL + } + return(res) + } + + nm = names(mp@pops) popList = list() - for (item in mp@pops) { + nameVec = character() + + for (i in seq_along(mp@pops)) { + item = mp@pops[[i]] + label = if (!is.null(nm) && nzchar(nm[i])) nm[i] else as.character(i) + if (isPop(item)) { popList = c(popList, list(item)) + nameVec = c(nameVec, label) } else if (isMultiPop(item)) { - popList = c(popList, .flattenMultiPop(item)) + child = .flattenMultiPop(item, preserveNames = preserveNames) + childNames = names(child) + + if (preserveNames %in% c("concatenate", "force")) { + if (is.null(childNames)) { + childNames = as.character(seq_along(child)) + } + childNames = paste(label, childNames, sep = "_") + names(child) = childNames + } + + popList = c(popList, child) + nameVec = c( + nameVec, + if (is.null(names(child))) { + rep(NA_character_, length(child)) + } else { + names(child) + } + ) } } - return(popList) -} + if (preserveNames == "none") { + names(popList) = NULL + return(popList) + } + + # decide whether to keep/transform names + if (preserveNames == "force") { + # nm_final = nameVec + nameVec[is.na(nameVec) | !nzchar(nameVec)] = as.character(which( + is.na(nameVec) | !nzchar(nameVec) + )) + if (any(duplicated(nameVec))) { + warning( + "Duplicate names found in 'force' mode. Making names unique by appending suffixes.", + call. = FALSE + ) + } + names(popList) = make.unique(nameVec) + return(popList) + } + + # auto or concatenate: keep only if no NA and no duplicates + if (length(nameVec) > 0L && + !any(is.na(nameVec)) && + !any(duplicated(nameVec)) + ) { + names(popList) = nameVec + } else { + names(popList) = NULL + } + + popList +} #' @title Merge Pop and MultiPop objects #' #' @description diff --git a/man/dot-flattenMultiPop.Rd b/man/dot-flattenMultiPop.Rd index 1a545491..13bfaf55 100644 --- a/man/dot-flattenMultiPop.Rd +++ b/man/dot-flattenMultiPop.Rd @@ -4,7 +4,7 @@ \alias{.flattenMultiPop} \title{Helper function to recursively extract Pop objects from a MultiPop} \usage{ -.flattenMultiPop(mp) +.flattenMultiPop(mp, preserveNames = c("auto", "concatenate", "force", "none")) } \arguments{ \item{mp}{\code{\link{MultiPop-class}} object} diff --git a/man/flattenMultiPop.Rd b/man/flattenMultiPop.Rd index 01d55d87..a0a213b2 100644 --- a/man/flattenMultiPop.Rd +++ b/man/flattenMultiPop.Rd @@ -4,19 +4,34 @@ \alias{flattenMultiPop} \title{Flatten a MultiPop object to a specified depth} \usage{ -flattenMultiPop(x, level = 1) +flattenMultiPop( + x, + level = 1, + preserveNames = c("auto", "concatenate", "force", "none") +) } \arguments{ \item{x}{A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object.} \item{level}{Integer scalar >= 1. Number of \code{MultiPop} levels to preserve.} + +\item{preserveNames}{Character scalar controlling how to handle names when +flattening. One of: \cr +- \code{"auto"} (default): Keep names only if all are present and + unique after flattening; otherwise drop them. \cr +- \code{"concatenate"}: Prefix child names with parent names using + underscore separator. Keep names only if resulting names are unique. \cr +- \code{"force"}: Always keep concatenated names, making them unique + via \code{\link{make.unique}} if duplicates are found (with a warning). \cr +- \code{"none"}: Drop all names.} } \value{ If \code{x} is a \code{\link{Pop-class}}, the same \code{x} object is returned. Otherwise a \code{\link{MultiPop-class}} is returned whose \code{x@pops} slot contains \code{\link{Pop-class}} (and possibly -\code{\link{MultiPop-class}}) objects flattened according to \code{level}. +\code{\link{MultiPop-class}}) objects flattened according to \code{level} +and \code{preserveNames}. } \description{ Recursively flatten a \code{\link{MultiPop-class}} object into a shallower @@ -33,6 +48,14 @@ The \code{level} argument controls how many levels of nesting are preserved. \c If \code{level} is greater than or equal to the nesting depth, the original object is returned unchanged. + +The \code{preserveNames} argument allows control over name preservation during +flattening. This is useful when you want to track the source of flattened +populations through their hierarchical names. The default \code{"auto"} mode +is conservative: names are kept only when they are meaningful (all present +and unique). Use \code{"concatenate"} to always build hierarchical names +(when possible), \code{"force"} to guarantee names with uniqueness enforcement, +or \code{"none"} to explicitly discard all names. } \examples{ # Create founder haplotypes @@ -46,16 +69,27 @@ SP = SimParam$new(founderPop) pop = newPop(founderPop, simParam=SP) # Create a multi-population with down to level 3 nesting -mp_nested = newMultiPop(pop[1:2], - newMultiPop(pop[3:4], - newMultiPop(pop[5:7], pop[8:12]))) +mp_nested = newMultiPop(pop1 = pop[1:2], + mp1 = newMultiPop(pop2 = pop[3:4], + mp2 = newMultiPop(pop3 = pop[5:7], pop4 = pop[8:12]))) mp_nested # Completely flatten to a single top-level MultiPop +# With default "auto" mode: names are kept if unique flattenMultiPop(mp_nested) -# Preserve two levels of nesting -flattenMultiPop(mp_nested, level=2) +# With "concatenate": build hierarchical names +flattenMultiPop(mp_nested, preserveNames = "concatenate") + +# With "force": guarantee unique names with suffixes +# (not needed here but useful if duplicates were present) +flattenMultiPop(mp_nested, preserveNames = "force") + +# With "none": drop all names +flattenMultiPop(mp_nested, preserveNames = "none") + +# Preserve two levels of nesting (level 2 names are kept) +flattenMultiPop(mp_nested, level = 2) } \seealso{ diff --git a/man/unnameMultiPop.Rd b/man/unnameMultiPop.Rd index cd6d2d27..d7487d81 100644 --- a/man/unnameMultiPop.Rd +++ b/man/unnameMultiPop.Rd @@ -39,6 +39,7 @@ founderPop = quickHaplo(nInd = 10, nChr = 1, segSites = 10) # Set simulation parameters SP = SimParam$new(founderPop) +\dontshow{SP$nThreads = 1L} # Create population pop = newPop(founderPop, simParam = SP) diff --git a/tests/testthat/test-mergePops.R b/tests/testthat/test-mergePops.R index e92eab31..f89f70b8 100644 --- a/tests/testthat/test-mergePops.R +++ b/tests/testthat/test-mergePops.R @@ -159,12 +159,12 @@ test_that("cMultiPop_mergeMultiPops_and_flattenMultiPop", { # Nothing to flatten, so the same input should be returned expect_identical(flattenMultiPop(mp1, level = 0), mp1) - expect_identical(flattenMultiPop(mp1, level = 1), mp1) - expect_identical(flattenMultiPop(mp1, level = 2), mp1) - tmp = flattenMultiPop(mp1) + tmp = flattenMultiPop(mp1, level = 1) + expect_identical(tmp, mp1) expect_identical(tmp[[1]], pop[1:2]) expect_identical(tmp[[2]], pop[3:5]) tmp = flattenMultiPop(mp1, level = 2) + expect_identical(tmp, mp1) expect_identical(tmp[[1]], pop[1:2]) expect_identical(tmp[[2]], pop[3:5]) @@ -190,6 +190,39 @@ test_that("cMultiPop_mergeMultiPops_and_flattenMultiPop", { expect_identical(tmp[[3]], newMultiPop(pop[8:9], pop[10:11], pop[12])) expect_equal(flattenMultiPop(mp3, level = 3), mp3) + + # Flattening with preserveNames + mp_named = newMultiPop(A = pop[1:2], B = pop[3:5], + C = newMultiPop(D = pop[6:7], E = pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "auto"), + newMultiPop(A = pop[1:2], B = pop[3:5], D = pop[6:7], E = pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "concatenate"), + newMultiPop(A = pop[1:2], B = pop[3:5], C_D = pop[6:7], C_E = pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "force"), + newMultiPop(A = pop[1:2], B = pop[3:5], C_D = pop[6:7], C_E = pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "none"), + newMultiPop(pop[1:2], pop[3:5], pop[6:7], pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 2, preserveNames = "none"), + newMultiPop(A = pop[1:2], B = pop[3:5], + C = newMultiPop(pop[6:7], pop[8:10]))) + + mp_named = newMultiPop(A = pop[1:2], B = pop[3:5], + C = newMultiPop(D = pop[6:7], D = pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "auto"), + newMultiPop(pop[1:2], pop[3:5], pop[6:7], pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "concatenate"), + newMultiPop(pop[1:2], pop[3:5], pop[6:7], pop[8:10])) + + expect_identical( + expect_warning(flattenMultiPop(mp_named, level = 1, preserveNames = "force"), + "Duplicate names found in 'force' mode. Making names unique by appending suffixes.", + fixed = TRUE), + newMultiPop(A = pop[1:2], B = pop[3:5], C_D = pop[6:7], C_D.1 = pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 1, preserveNames = "none"), + newMultiPop(pop[1:2], pop[3:5], pop[6:7], pop[8:10])) + expect_identical(flattenMultiPop(mp_named, level = 2, preserveNames = "none"), + newMultiPop(A = pop[1:2], B = pop[3:5], + C = newMultiPop(pop[6:7], pop[8:10]))) }) test_that("splitPop", { From af2a49ed97a8e843efd39add6d5fe86d5bd1f871 Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Thu, 4 Jun 2026 18:02:54 +0200 Subject: [PATCH 3/9] Update `calcPopValue` to support simplification of output matrices by `level` --- NAMESPACE | 1 + R/selection.R | 259 ++++++++++++++++++++++++-------- man/calcPopValue.Rd | 96 ++++++++---- man/dot-collectLeafPaths.Rd | 17 +++ man/dot-formatCalcPopSource.Rd | 19 +++ tests/testthat/test-selection.R | 137 +++++++++++++++++ 6 files changed, 438 insertions(+), 91 deletions(-) create mode 100644 man/dot-collectLeafPaths.Rd create mode 100644 man/dot-formatCalcPopSource.Rd diff --git a/NAMESPACE b/NAMESPACE index 5fd04848..0a287a00 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -19,6 +19,7 @@ export(attrition) export(bv) export(cChr) export(calcGCA) +export(calcPopValue) export(dd) export(doubleGenome) export(ebv) diff --git a/R/selection.R b/R/selection.R index 9f466320..0006bbd2 100644 --- a/R/selection.R +++ b/R/selection.R @@ -85,36 +85,78 @@ getResponse = function(pop,trait,use,simParam=NULL,nThreads=NULL,...){ return(response) } -#' Calculate summary response from a population +#' @title Calculate values from Pop or MultiPop objects #' -#' Calculates a summary response from a \code{\link{Pop-class}} or -#' \code{\link{MultiPop-class}} object. For \code{MultiPop} objects, -#' the function is applied recursively to all populations. +#' @description +#' Apply a summary function to each \code{\link{Pop-class}} in a +#' \code{\link{MultiPop-class}} and optionally simplify the output to a +#' requested level of nesting. #' -#' @param x a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object. -#' @param trait the trait for selection. Either a number indicating -#' a single trait or a character for a trait name, or a function -#' returning a vector of values. -#' @param use the selection criterion. Either a character ("rand", "gv", -#' "ebv", or "pheno") or a custom function. -#' @param FUN a summary function to be applied to the vector or matrix of -#' values. Default is \code{mean}. -#' @param FUN.ARGS a list of additional arguments passed to \code{FUN}. -#' @param returnList logical to return a list (when \code{TRUE}) -#' or a vector (when \code{FALSE} and the return object can be -#' simplified to a vector). -#' @param simParam an object of \code{\link{SimParam}}. -#' @param ... additional arguments passed to \code{trait} or \code{use} -#' when they are custom functions. +#' @param x A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object. +#' @param FUN Function to apply to each \code{Pop}. Must accept a \code{Pop} +#' as its first argument. +#' @param simplify Logical. If \code{TRUE}, simplify the output by flattening +#' the \code{MultiPop} in \code{x} to the requested \code{level} using +#' \code{\link{flattenMultiPop}}. The output matrices from each \code{Pop} +#' are combined with \code{\link{rbind}}, and a \code{"source"} attribute +#' is added to indicate the origin of each row. +#' @param level Integer scalar >= 1. Number of \code{MultiPop} levels to +#' preserve when \code{simplify=TRUE}. Passed to \code{\link{flattenMultiPop}}. +#' Ignored if \code{simplify=FALSE}. +#' @param simParam an object of class \code{\link{SimParam}}. If \code{NULL}, +#' the function uses the object named \code{SP} from the global environment. +#' @param ... Additional arguments passed to \code{FUN}. #' -#' @keywords internal +#' @details +#' The \code{level} argument controls the depth of nesting retained when +#' \code{simplify=TRUE}. If \code{level} exceeds the depth of \code{x}, the +#' output structure is returned unchanged. The \code{"source"} attribute is +#' a data frame with columns \code{level1}, \code{level2}, etc., indicating +#' the origin of each row in the simplified output. +#' +#' @return If \code{x} is a \code{Pop}, returns the value of \code{FUN(x, ...)}. +#' Otherwise returns a list (or a simplified matrix) of results, optionally with +#' a \code{"source"} attribute. +#' +#' @seealso \code{\link{MultiPop-class}}, \code{\link{flattenMultiPop}} +#' +#' @examples +#' # Create founder haplotypes +#' founderPop = quickHaplo(nInd=6, nChr=1, segSites=10) +#' +#' # Set simulation parameters +#' SP = SimParam$new(founderPop) +#' \dontshow{SP$nThreads = 1L} +#' SP$addTraitA(10, mean = c(0, 0), var = c(1, 1)) +#' SP$setVarE(h2 = c(0.6, 0.4)) +#' +#' # Create population +#' pop = newPop(founderPop, simParam=SP) +#' pop2 = randCross(pop, nCrosses = 3, nProgeny = 4) +#' +#' # Create a nested MultiPop +#' mp1 = splitPop( +#' pop2, +#' by = list( +#' function(x) rep(LETTERS[1:2], length.out = length(x)), +#' function(x) paste(x@mother, x@father, sep = "_") +#' ) +#' ) +#' +#' # Calculate phenotypes at different levels +#' calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 2) +#' calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 1) +#' +#' # Custom function returning a summary matrix +#' calcPopValue(mp1, FUN = function(x) colMeans(pheno(x)), +#' simplify = TRUE, level = 1) +#' +#' @export calcPopValue = function( x, - trait = 1, - use = 'pheno', - FUN = mean, - FUN.ARGS = list(), - returnList = FALSE, + FUN, + simplify = FALSE, + level = 1, simParam = NULL, ... ) { @@ -122,46 +164,72 @@ calcPopValue = function( simParam = get("SP", envir = .GlobalEnv) } - if (isMultiPop(x)) { - popValueList = lapply(x@pops, function(x) { - calcPopValue( - x, - trait = trait, - use = use, - FUN = FUN, - FUN.ARGS = FUN.ARGS, - returnList = returnList, - simParam = simParam, - ... - ) - }) + dots = list(...) + .level_offset = dots[[".level_offset"]] + if (is.null(.level_offset)) { + .level_offset = 0L + } + dots[[".level_offset"]] = NULL - if (returnList || any(sapply(popValueList, length) != 1)) { - return(popValueList) - } else { - popValue = do.call('c', popValueList) - return(unname(popValue)) - } - } else if (isPop(x)) { - # TODO: Update support for use='bv' with genParamPop() - if (is.character(use)) { - if (use == 'bv') { - stop("use='bv' is not currently supported for populations") - } + if (isPop(x)) { + FUN.ARGS = c(list(x), dots) + return(do.call(FUN, FUN.ARGS)) + } + + if (!isMultiPop(x)) { + stop("`x` must be a Pop or MultiPop object.") + } + + source = NULL + if (simplify && level == 1L) { + source = .collectLeafPaths(x) + } + + if (simplify) { + x = flattenMultiPop(x, level = level) + } + + popValueList = lapply( + x@pops, + function(pop) { + FUN.ARGS = c( + list( + pop, + FUN = FUN, + simplify = simplify, + level = level - 1L, + simParam = simParam, + .level_offset = .level_offset + 1L + ), + dots + ) + do.call(calcPopValue, FUN.ARGS) } - response = getResponse( - x, - trait = trait, - use = use, - simParam = simParam, - ... + ) + + if (simplify && level == 1L) { + popValue = do.call(rbind, popValueList) + src = .formatCalcPopSource( + source, + popValueList, + level_offset = .level_offset ) - if (is.matrix(response)) { - stopifnot(ncol(response) == 1) + + if (NROW(src) != NROW(popValue)) { + warning( + "The number of rows in the source data frame does not match the number of rows in the output matrix.", + call. = FALSE + ) } - FUN.ARGS = append(FUN.ARGS, list(response), after = 0) - return(do.call(FUN, FUN.ARGS)) + attr(popValue, "source") = src + return(popValue) } + + if (!is.null(names(x@pops))) { + names(popValueList) = names(x@pops) + } + + return(popValueList) } #' Identify candidate individuals @@ -905,3 +973,74 @@ selectPop = function( return(x[take[0:nPop]]) } + +#' Helper function to collect leaf paths in a \code{MultiPop} +#' +#' @param x \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object. +#' @param path Character vector of labels representing the current path. +#' +#' @keywords internal +.collectLeafPaths = function(x, path = list()) { + if (isPop(x)) { + return(list(path)) + } + stopifnot(isMultiPop(x)) + + out = list() + nm = names(x@pops) + + for (i in seq_along(x@pops)) { + child = x@pops[[i]] + if (!is.null(nm) && length(nm) >= i && !is.na(nm[i]) && nzchar(nm[i])) { + label = nm[i] + } else { + label = i + } + childPath = c(path, list(label)) + out = c(out, .collectLeafPaths(child, path = childPath)) + } + + return(out) +} + +#' Helper function to format the source attribute for \code{calcPopValue} +#' +#' @param paths List of character vectors representing leaf paths. +#' @param values List of values returned from \code{FUN}. +#' @param level_offset Integer scalar used to label source levels. +#' +#' @keywords internal +.formatCalcPopSource = function(paths, values, level_offset = 0L) { + if (length(paths) == 0L) { + return(NULL) + } + + maxDepth = max(lengths(paths)) + pathList = lapply(paths, function(p) { + pad = rep(list(NA), maxDepth - length(p)) + c(p, pad) + }) + + pathDf = data.frame(sapply(pathList, `[[`, 1)) + for (i in seq_len(maxDepth)[-1]) { + pathDf = cbind(pathDf, sapply(pathList, `[[`, i)) + } + colnames(pathDf) = paste0("level", level_offset + seq_len(maxDepth)) + + for (j in seq_len(ncol(pathDf))) { + col = pathDf[[j]] + first = col[!is.na(col)][1] + if (!is.null(first) && is.numeric(first)) { + pathDf[[j]] = as.integer(col) + } else { + pathDf[[j]] = as.character(col) + } + } + + nRows = vapply(values, NROW, integer(1L)) + idx = rep(seq_along(nRows), nRows) + pathDf = pathDf[idx, , drop = FALSE] + rownames(pathDf) = NULL + + return(pathDf) +} \ No newline at end of file diff --git a/man/calcPopValue.Rd b/man/calcPopValue.Rd index 4d4ac004..1560b75f 100644 --- a/man/calcPopValue.Rd +++ b/man/calcPopValue.Rd @@ -2,46 +2,80 @@ % Please edit documentation in R/selection.R \name{calcPopValue} \alias{calcPopValue} -\title{Calculate summary response from a population} +\title{Calculate values from Pop or MultiPop objects} \usage{ -calcPopValue( - x, - trait = 1, - use = "pheno", - FUN = mean, - FUN.ARGS = list(), - returnList = FALSE, - simParam = NULL, - ... -) +calcPopValue(x, FUN, simplify = FALSE, level = 1, simParam = NULL, ...) } \arguments{ -\item{x}{a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object.} +\item{x}{A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object.} + +\item{FUN}{Function to apply to each \code{Pop}. Must accept a \code{Pop} +as its first argument.} + +\item{simplify}{Logical. If \code{TRUE}, simplify the output by flattening +the \code{MultiPop} in \code{x} to the requested \code{level} using +\code{\link{flattenMultiPop}}. The output matrices from each \code{Pop} + are combined with \code{\link{rbind}}, and a \code{"source"} attribute +is added to indicate the origin of each row.} -\item{trait}{the trait for selection. Either a number indicating -a single trait or a character for a trait name, or a function -returning a vector of values.} +\item{level}{Integer scalar >= 1. Number of \code{MultiPop} levels to +preserve when \code{simplify=TRUE}. Passed to \code{\link{flattenMultiPop}}. +Ignored if \code{simplify=FALSE}.} + +\item{simParam}{an object of class \code{\link{SimParam}}. If \code{NULL}, +the function uses the object named \code{SP} from the global environment.} + +\item{...}{Additional arguments passed to \code{FUN}.} +} +\value{ +If \code{x} is a \code{Pop}, returns the value of \code{FUN(x, ...)}. +Otherwise returns a list (or a simplified matrix) of results, optionally with +a \code{"source"} attribute. +} +\description{ +Apply a summary function to each \code{\link{Pop-class}} in a +\code{\link{MultiPop-class}} and optionally simplify the output to a +requested level of nesting. +} +\details{ +The \code{level} argument controls the depth of nesting retained when +\code{simplify=TRUE}. If \code{level} exceeds the depth of \code{x}, the +output structure is returned unchanged. The \code{"source"} attribute is +a data frame with columns \code{level1}, \code{level2}, etc., indicating +the origin of each row in the simplified output. +} +\examples{ +# Create founder haplotypes +founderPop = quickHaplo(nInd=6, nChr=1, segSites=10) -\item{use}{the selection criterion. Either a character ("rand", "gv", -"ebv", or "pheno") or a custom function.} +# Set simulation parameters +SP = SimParam$new(founderPop) +\dontshow{SP$nThreads = 1L} +SP$addTraitA(10, mean = c(0, 0), var = c(1, 1)) +SP$setVarE(h2 = c(0.6, 0.4)) -\item{FUN}{a summary function to be applied to the vector or matrix of -values. Default is \code{mean}.} +# Create population +pop = newPop(founderPop, simParam=SP) +pop2 = randCross(pop, nCrosses = 3, nProgeny = 4) -\item{FUN.ARGS}{a list of additional arguments passed to \code{FUN}.} +# Create a nested MultiPop +mp1 = splitPop( + pop2, + by = list( + function(x) rep(LETTERS[1:2], length.out = length(x)), + function(x) paste(x@mother, x@father, sep = "_") + ) +) -\item{returnList}{logical to return a list (when \code{TRUE}) -or a vector (when \code{FALSE} and the return object can be -simplified to a vector).} +# Calculate phenotypes at different levels +calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 2) +calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 1) -\item{simParam}{an object of \code{\link{SimParam}}.} +# Custom function returning a summary matrix +calcPopValue(mp1, FUN = function(x) colMeans(pheno(x)), + simplify = TRUE, level = 1) -\item{...}{additional arguments passed to \code{trait} or \code{use} -when they are custom functions.} } -\description{ -Calculates a summary response from a \code{\link{Pop-class}} or -\code{\link{MultiPop-class}} object. For \code{MultiPop} objects, -the function is applied recursively to all populations. +\seealso{ +\code{\link{MultiPop-class}}, \code{\link{flattenMultiPop}} } -\keyword{internal} diff --git a/man/dot-collectLeafPaths.Rd b/man/dot-collectLeafPaths.Rd new file mode 100644 index 00000000..313e3eea --- /dev/null +++ b/man/dot-collectLeafPaths.Rd @@ -0,0 +1,17 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/selection.R +\name{.collectLeafPaths} +\alias{.collectLeafPaths} +\title{Helper function to collect leaf paths in a \code{MultiPop}} +\usage{ +.collectLeafPaths(x, path = list()) +} +\arguments{ +\item{x}{\code{\link{Pop-class}} or \code{\link{MultiPop-class}} object.} + +\item{path}{Character vector of labels representing the current path.} +} +\description{ +Helper function to collect leaf paths in a \code{MultiPop} +} +\keyword{internal} diff --git a/man/dot-formatCalcPopSource.Rd b/man/dot-formatCalcPopSource.Rd new file mode 100644 index 00000000..283bb7a7 --- /dev/null +++ b/man/dot-formatCalcPopSource.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/selection.R +\name{.formatCalcPopSource} +\alias{.formatCalcPopSource} +\title{Helper function to format the source attribute for \code{calcPopValue}} +\usage{ +.formatCalcPopSource(paths, values, level_offset = 0L) +} +\arguments{ +\item{paths}{List of character vectors representing leaf paths.} + +\item{values}{List of values returned from \code{FUN}.} + +\item{level_offset}{Integer scalar used to label source levels.} +} +\description{ +Helper function to format the source attribute for \code{calcPopValue} +} +\keyword{internal} diff --git a/tests/testthat/test-selection.R b/tests/testthat/test-selection.R index 767e4d41..18262d55 100644 --- a/tests/testthat/test-selection.R +++ b/tests/testthat/test-selection.R @@ -182,3 +182,140 @@ test_that("selectPop_and_calcPopValue",{ selectPop(mp3[[3]][[2]], nPop = 1, level = 1, simParam = SP) ) }) + +# tests/testthat/test-calcPopValue.R +test_that("calcPopValue", { + founderPop = quickHaplo(nInd = 14, nChr = 1, segSites = 10) + SP = SimParam$new(founderPop) + SP$nThreads = 1L + SP$addTraitA(10, mean = c(0, 0), var = c(1, 1)) + SP$setVarE(h2 = c(0.6, 0.4)) + pop = newPop(founderPop, simParam = SP) + + # Pop case + expect_identical( + calcPopValue(pop, FUN = pheno, simplify = FALSE, simParam = SP), + pheno(pop) + ) + expect_identical( + calcPopValue( + pop, + FUN = function(x, ...) cor(x@pheno, ...), + simplify = TRUE, + simParam = SP, + method = "kendall" + ), + cor(pheno(pop), method = "kendall") + ) + + # Nested MultiPop case + mp1 = splitPop( + randCross(pop, nCrosses = 3, nProgeny = 4, simParam = SP), + by = list( + function(x) rep(LETTERS[1:2], length.out = length(x)), + function(x) getFam(x, famType = "B") + ) + ) + + # simplify=FALSE preserves nesting and names + outF = calcPopValue(mp1, FUN = pheno, simplify = FALSE, simParam = SP) + expect_true(is.list(outF)) + expect_identical(names(outF), names(mp1)) + expect_true(is.list(outF$A)) + expect_identical(outF$A, lapply(mp1$A@pops, pheno)) + expect_identical(outF$B, lapply(mp1$B@pops, pheno)) + + # level > depth returns unchanged structure + expect_identical( + outF, + calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 3, simParam = SP) + ) + + # simplify=TRUE level=2 attaches source with level2 + out2 = calcPopValue( + mp1, + FUN = pheno, + simplify = TRUE, + level = 2, + simParam = SP + ) + expect_true(is.list(out2)) + expect_identical(names(out2), names(mp1)) + + outA = calcPopValue(mp1$A, FUN = pheno, simplify = TRUE, simParam = SP) + expect_identical(attributes(out2$A)$source[[1]], attributes(outA)$source[[1]]) + expect_equal(nrow(attributes(outA)$source), nrow(outA)) + + outB = calcPopValue(mp1$B, FUN = pheno, simplify = TRUE, simParam = SP) + expect_identical(attributes(out2$B)$source[[1]], attributes(outB)$source[[1]]) + expect_equal(nrow(attributes(outB)$source), nrow(outB)) + + expect_equal(rbind(outA), do.call('rbind', lapply(mp1$A@pops, pheno))) + expect_equal(rbind(outB), do.call('rbind', lapply(mp1$B@pops, pheno))) + + # simplify=TRUE level=1 attaches level1/level2 source + out1 = calcPopValue( + mp1, + FUN = pheno, + simplify = TRUE, + level = 1, + simParam = SP + ) + expect_true(is.matrix(out1)) + expect_equal(nrow(out1), nInd(mergeMultiPops(mp1))) + expect_identical( + rbind(attributes(out2$A)$source, attributes(out2$B)$source), + attributes(out1)$source[, 2, drop = FALSE] + ) + expect_equal(rbind(out1), rbind(outA, outB)) + + # Ragged and nested MultiPop + mp3 = newMultiPop( + p1 = pop[1:2], + mp1 = newMultiPop( + p2 = pop[3:4], + mp2 = newMultiPop(p3 = pop[5:6], p4 = pop[7:8]) + ), + mp3 = newMultiPop( + p5 = pop[9:10], + mp4 = newMultiPop(p6 = pop[11:12], p7 = pop[13:14]) + ) + ) + + # simplify=FALSE preserves nesting and names + outF = calcPopValue(mp3, FUN = pheno, simplify = FALSE, simParam = SP) + expect_true(is.list(outF)) + expect_identical(names(outF), names(mp3)) + expect_identical(names(outF$mp1), names(mp3$mp1)) + expect_identical(names(outF$mp3), names(mp3$mp3)) + expect_identical(outF$mp1$mp2, lapply(mp3$mp1$mp2@pops, pheno)) + expect_identical(outF$mp3$mp4, lapply(mp3$mp3$mp4@pops, pheno)) + + # simplify=TRUE level=1 attaches level1/level2/level3 source + out1 = calcPopValue( + mp3, + FUN = pheno, + simplify = TRUE, + level = 1, + simParam = SP + ) + expect_true(is.matrix(out1)) + expect_equal(nrow(out1), nInd(mergeMultiPops(mp3))) + + # Error handling + expect_error( + calcPopValue("not_a_pop", FUN = pheno, simplify = FALSE, simParam = SP), + "`x` must be a Pop or MultiPop object.", + fixed = TRUE + ) + expect_error( + calcPopValue(1:5, FUN = pheno, simplify = FALSE, simParam = SP), + "`x` must be a Pop or MultiPop object.", + fixed = TRUE + ) + expect_error( + calcPopValue(list(pop), FUN = pheno, simplify = FALSE, simParam = SP), + "`x` must be a Pop or MultiPop object.", + fixed = TRUE + ) +}) From 08274262c9f70661c57c6b3ce53d9b214860f940 Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Thu, 4 Jun 2026 18:12:43 +0200 Subject: [PATCH 4/9] Update `selectPop` to use the updated `calcPopValue` --- NEWS.md | 6 +++- R/selection.R | 52 ++++++++++++++++++++------------- tests/testthat/test-selection.R | 31 ++++++++++++-------- 3 files changed, 56 insertions(+), 33 deletions(-) diff --git a/NEWS.md b/NEWS.md index 350f14d7..9c038384 100644 --- a/NEWS.md +++ b/NEWS.md @@ -12,12 +12,16 @@ * Added `flattenMultiPop()` function to reduce the depth of a nested MultiPop object up to a defined `level`. -* Added `selectPop()` function to select a subset of populations from a MultiPop. This function uses the new internal `calcPopValue()` function to calculate a summary value for each population, which is then used for selection. +* Added `calcPopValue()` function to calculate summary values for each Pop in a MultiPop. + +* Added `selectPop()` function to select a subset of populations from a MultiPop. This function uses the new `calcPopValue()` function. * Added `unnameMultiPop()` function to remove names from a MultiPop at specific `level`(s). * Added `newEmptyMultiPop()` function to create an empty MultiPop. It has the same outcome as `newMultiPop()` with no arguments. To make this behavior consistent with `newEmptyPop()`, `newPop()` now also supports creating an empty population with defined ploidy. +* Added `splitPop()` function to divide a Pop or MultiPop object into multiple subpopulations at specified level(s). + * Update `mergePops()` to handle nested MultiPop objects. * Added asLogNormal() function. diff --git a/R/selection.R b/R/selection.R index 0006bbd2..d4797156 100644 --- a/R/selection.R +++ b/R/selection.R @@ -897,6 +897,7 @@ selectPop = function( FUN = mean, FUN.ARGS = list() ) { + stopifnot(nPop >= 0) if (is.null(simParam)) { simParam = get("SP", envir = .GlobalEnv) @@ -907,17 +908,22 @@ selectPop = function( } stopifnot(isMultiPop(x)) - multi = which(sapply(unname(x@pops), isMultiPop)) + if (length(level) != 1L || !is.numeric(level) || is.na(level) || + level < 1 || level != as.integer(level)) { + stop("`level` must be a positive integer") + } - if (level > 1 & identical(multi, integer(0))) { - stop(paste( - "The MultiPop object does not contain other MultiPop objects", - "at this level. You may want to decrease the value of 'level'" - )) + md = .depthMultiPop(x) + if (md < level) { + stop(sprintf("requested `level` exceeds max depth of `x` (%d)", md)) } - while (level > 1) { - level = level - 1 + is_multi = vapply(unname(x@pops), isMultiPop, logical(1L)) + is_pop = vapply(unname(x@pops), isPop, logical(1L)) + multi = which(is_multi) + + while (level > 1L) { + level = level - 1L for (i in multi) { x@pops[[i]] = selectPop( x = x[[i]], @@ -932,8 +938,7 @@ selectPop = function( ... ) } - multiPop = do.call(newMultiPop, x@pops) - return(multiPop) + return(x) } if (!identical(multi, integer(0))) { @@ -946,7 +951,7 @@ selectPop = function( )) } - eligible = which(sapply(x@pops, isPop)) + eligible = which(is_pop) if (length(eligible) < nPop) { nPop = length(eligible) @@ -957,21 +962,28 @@ selectPop = function( ) } - popValues = calcPopValue( + if (is.character(use) && use == 'bv') { + stop("use='bv' is not currently supported for populations") + } + + response = calcPopValue( x, - trait = trait, - use = use, - FUN = FUN, - returnList = FALSE, - FUN.ARGS = FUN.ARGS, - simParam = simParam, - ... + FUN = function(pop) { + getResponse(pop = pop, trait = trait, use = use, simParam = simParam, ...) + }, + simplify = FALSE, + level = 1L, + simParam = simParam ) + popValues = vapply(response, function(res) { + do.call(FUN, c(list(res), FUN.ARGS)) + }, numeric(1L)) + take = order(popValues, decreasing = selectTop) take = take[take %in% eligible] - return(x[take[0:nPop]]) + return(x[take[seq_len(nPop)]]) } #' Helper function to collect leaf paths in a \code{MultiPop} diff --git a/tests/testthat/test-selection.R b/tests/testthat/test-selection.R index 18262d55..a7b6059b 100644 --- a/tests/testthat/test-selection.R +++ b/tests/testthat/test-selection.R @@ -37,8 +37,8 @@ test_that("selectInd_and_getResponse",{ expect_equal(pop5@id, pop6@id) }) -test_that("selectPop_and_calcPopValue",{ - founderPop = quickHaplo(nInd=100, nChr=1, segSites=10) +test_that("selectPop",{ + founderPop = quickHaplo(nInd=14, nChr=1, segSites=10) SP = SimParam$new(founderPop) SP$nThreads = 1L SP$addTraitA(10, mean = c(0, 0), var = c(1, 1)) @@ -70,6 +70,14 @@ test_that("selectPop_and_calcPopValue",{ paste("Suitable candidate populations smaller than nPop, returning", length(mp1), "populations")) + # Invalid level + expect_error(selectPop(mp1, nPop = 2, level = 0.5, simParam = SP), + "`level` must be a positive integer", + fixed=TRUE) + expect_error(selectPop(mp1, nPop = 2, level = 1:2, simParam = SP), + "`level` must be a positive integer", + fixed=TRUE) + # TODO: Update support for use='bv' with genParamPop() # bv is not currently supported expect_error(selectPop(mp1, nPop = 2, use = 'bv', simParam = SP), @@ -78,8 +86,7 @@ test_that("selectPop_and_calcPopValue",{ # Selecting nested populations from a non-nested object expect_error(selectPop(mp1, nPop = 1, level = 3, use = 'pheno', simParam = SP), - paste("The MultiPop object does not contain other MultiPop objects", - "at this level. You may want to decrease the value of 'level'"), + "requested `level` exceeds max depth of `x` (1)", fixed=TRUE) # Use a custom trait obtained by summing up the two available traits @@ -119,9 +126,9 @@ test_that("selectPop_and_calcPopValue",{ # MultiPop with 1 nested object - mp2 = newMultiPop(pop[1:20], pop[21:40], - newMultiPop(pop[41:60], - newMultiPop(pop[61:80], pop[81:100]))) + mp2 = newMultiPop(pop[1:2], pop[3:4], + newMultiPop(pop[5:6], + newMultiPop(pop[7:8], pop[9:10]))) # setPheno for all traits mp2 = setPheno(mp2, varE = c(1,1), simParam = SP) @@ -146,11 +153,11 @@ test_that("selectPop_and_calcPopValue",{ # MultiPop with >1 nested object - mp3 = newMultiPop(pop[1:20], - newMultiPop(pop[21:30], - newMultiPop(pop[31:35], pop[36:40])), - newMultiPop(pop[41:60], - newMultiPop(pop[61:80], pop[81:100]))) + mp3 = newMultiPop(pop[1:2], + newMultiPop(pop[3:4], + newMultiPop(pop[5:6], pop[7:8])), + newMultiPop(pop[9:10], + newMultiPop(pop[11:12], pop[13:14]))) mp3 = setPheno(mp3, varE = c(1,1), simParam = SP) From 9c298a3d00c73f2c978770f7ff5d00f18a46c440 Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Fri, 5 Jun 2026 15:16:34 +0200 Subject: [PATCH 5/9] Fix bug in `calcPopValue` where the nrows of source and simplified output were different --- R/selection.R | 41 ++++++++++++++++++++++----------- man/calcPopValue.Rd | 4 ++-- man/dot-formatCalcPopSource.Rd | 4 ++-- tests/testthat/test-selection.R | 7 +++++- 4 files changed, 38 insertions(+), 18 deletions(-) diff --git a/R/selection.R b/R/selection.R index d4797156..d163ed1f 100644 --- a/R/selection.R +++ b/R/selection.R @@ -103,7 +103,7 @@ getResponse = function(pop,trait,use,simParam=NULL,nThreads=NULL,...){ #' @param level Integer scalar >= 1. Number of \code{MultiPop} levels to #' preserve when \code{simplify=TRUE}. Passed to \code{\link{flattenMultiPop}}. #' Ignored if \code{simplify=FALSE}. -#' @param simParam an object of class \code{\link{SimParam}}. If \code{NULL}, +#' @param simParam an object of class \code{\link{SimParam}}. If \code{NULL}, #' the function uses the object named \code{SP} from the global environment. #' @param ... Additional arguments passed to \code{FUN}. #' @@ -148,7 +148,7 @@ getResponse = function(pop,trait,use,simParam=NULL,nThreads=NULL,...){ #' calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 1) #' #' # Custom function returning a summary matrix -#' calcPopValue(mp1, FUN = function(x) colMeans(pheno(x)), +#' calcPopValue(mp1, FUN = function(x) colMeans(pheno(x)), #' simplify = TRUE, level = 1) #' #' @export @@ -208,20 +208,36 @@ calcPopValue = function( ) if (simplify && level == 1L) { - popValue = do.call(rbind, popValueList) - src = .formatCalcPopSource( - source, + nRows = vapply( popValueList, - level_offset = .level_offset + function(v) { + if (is.null(v)) { + return(0L) + } + if (is.data.frame(v) || is.matrix(v) || is.array(v)) { + return(dim(v)[1]) # number of rows + } + if (is.atomic(v) && is.null(dim(v))) { + return(1L) # vectors are one row + } + return(NA_integer_) + }, + integer(1L) ) - - if (NROW(src) != NROW(popValue)) { + if (any(is.na(nRows))) { warning( - "The number of rows in the source data frame does not match the number of rows in the output matrix.", + "Some values returned by FUN have unsupported types for simplification. Returning list output.", call. = FALSE ) + return(popValueList) } - attr(popValue, "source") = src + + popValue = do.call(rbind, popValueList) + attr(popValue, "source") = .formatCalcPopSource( + paths = source, + nRows = nRows, + level_offset = .level_offset + ) return(popValue) } @@ -1018,11 +1034,11 @@ selectPop = function( #' Helper function to format the source attribute for \code{calcPopValue} #' #' @param paths List of character vectors representing leaf paths. -#' @param values List of values returned from \code{FUN}. +#' @param nRows Vector of row counts for each path. #' @param level_offset Integer scalar used to label source levels. #' #' @keywords internal -.formatCalcPopSource = function(paths, values, level_offset = 0L) { +.formatCalcPopSource = function(paths, nRows, level_offset = 0L) { if (length(paths) == 0L) { return(NULL) } @@ -1049,7 +1065,6 @@ selectPop = function( } } - nRows = vapply(values, NROW, integer(1L)) idx = rep(seq_along(nRows), nRows) pathDf = pathDf[idx, , drop = FALSE] rownames(pathDf) = NULL diff --git a/man/calcPopValue.Rd b/man/calcPopValue.Rd index 1560b75f..81ef5eb0 100644 --- a/man/calcPopValue.Rd +++ b/man/calcPopValue.Rd @@ -22,7 +22,7 @@ is added to indicate the origin of each row.} preserve when \code{simplify=TRUE}. Passed to \code{\link{flattenMultiPop}}. Ignored if \code{simplify=FALSE}.} -\item{simParam}{an object of class \code{\link{SimParam}}. If \code{NULL}, +\item{simParam}{an object of class \code{\link{SimParam}}. If \code{NULL}, the function uses the object named \code{SP} from the global environment.} \item{...}{Additional arguments passed to \code{FUN}.} @@ -72,7 +72,7 @@ calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 2) calcPopValue(mp1, FUN = pheno, simplify = TRUE, level = 1) # Custom function returning a summary matrix -calcPopValue(mp1, FUN = function(x) colMeans(pheno(x)), +calcPopValue(mp1, FUN = function(x) colMeans(pheno(x)), simplify = TRUE, level = 1) } diff --git a/man/dot-formatCalcPopSource.Rd b/man/dot-formatCalcPopSource.Rd index 283bb7a7..7665788f 100644 --- a/man/dot-formatCalcPopSource.Rd +++ b/man/dot-formatCalcPopSource.Rd @@ -4,12 +4,12 @@ \alias{.formatCalcPopSource} \title{Helper function to format the source attribute for \code{calcPopValue}} \usage{ -.formatCalcPopSource(paths, values, level_offset = 0L) +.formatCalcPopSource(paths, nRows, level_offset = 0L) } \arguments{ \item{paths}{List of character vectors representing leaf paths.} -\item{values}{List of values returned from \code{FUN}.} +\item{nRows}{Vector of row counts for each path.} \item{level_offset}{Integer scalar used to label source levels.} } diff --git a/tests/testthat/test-selection.R b/tests/testthat/test-selection.R index a7b6059b..dfa67845 100644 --- a/tests/testthat/test-selection.R +++ b/tests/testthat/test-selection.R @@ -309,7 +309,12 @@ test_that("calcPopValue", { expect_true(is.matrix(out1)) expect_equal(nrow(out1), nInd(mergeMultiPops(mp3))) - # Error handling + # Error and warning handling + expect_warning( + calcPopValue(mp1, FUN = \(x) list(x@pheno), simplify = TRUE, simParam = SP), + "Some values returned by FUN have unsupported types for simplification. Returning list output.", + fixed = TRUE + ) expect_error( calcPopValue("not_a_pop", FUN = pheno, simplify = FALSE, simParam = SP), "`x` must be a Pop or MultiPop object.", From 53b14dcdb1092d31960cd5e34293e5312eea0389 Mon Sep 17 00:00:00 2001 From: darizasu Date: Fri, 5 Jun 2026 13:20:01 +0000 Subject: [PATCH 6/9] Update documentation --- DESCRIPTION | 2 +- man/AlphaSimR-package.Rd | 5 + man/SimParam.Rd | 1929 ++++++++++++++++++-------------------- 3 files changed, 893 insertions(+), 1043 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 1daa5b74..d95e0ee0 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -43,7 +43,7 @@ Depends: R (>= 4.0.0) Imports: Rcpp (>= 0.12.7), Rdpack, methods, R6 RdMacros: Rdpack LinkingTo: Rcpp, RcppArmadillo (>= 0.7.500.0.0), BH, dqrng (>= 0.4.1) -RoxygenNote: 7.3.3 Suggests: knitr, rmarkdown, testthat VignetteBuilder: knitr NeedsCompilation: true +Config/roxygen2/version: 8.0.0 diff --git a/man/AlphaSimR-package.Rd b/man/AlphaSimR-package.Rd index a4b2a8fe..fb8e0658 100644 --- a/man/AlphaSimR-package.Rd +++ b/man/AlphaSimR-package.Rd @@ -35,6 +35,11 @@ Useful links: \author{ \strong{Maintainer}: Chris Gaynor \email{gaynor.robert@hotmail.com} (\href{https://orcid.org/0000-0003-0558-6656}{ORCID}) +Authors: +\itemize{ + \item Chris Gaynor \email{gaynor.robert@hotmail.com} (\href{https://orcid.org/0000-0003-0558-6656}{ORCID}) +} + Other contributors: \itemize{ \item Gregor Gorjanc (\href{https://orcid.org/0000-0001-8008-2787}{ORCID}) [contributor] diff --git a/man/SimParam.Rd b/man/SimParam.Rd index fa63e257..c06b5b74 100644 --- a/man/SimParam.Rd +++ b/man/SimParam.Rd @@ -17,7 +17,7 @@ new trait values. \examples{ ## ------------------------------------------------ -## Method `SimParam$new` +## Method `SimParam$new()` ## ------------------------------------------------ #Create founder haplotypes @@ -27,7 +27,7 @@ founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) SP = SimParam$new(founderPop) ## ------------------------------------------------ -## Method `SimParam$setTrackPed` +## Method `SimParam$setTrackPed()` ## ------------------------------------------------ #Create founder haplotypes @@ -39,7 +39,7 @@ SP = SimParam$new(founderPop) SP$setTrackPed(TRUE) ## ------------------------------------------------ -## Method `SimParam$setTrackRec` +## Method `SimParam$setTrackRec()` ## ------------------------------------------------ #Create founder haplotypes @@ -51,7 +51,7 @@ SP = SimParam$new(founderPop) SP$setTrackRec(TRUE) ## ------------------------------------------------ -## Method `SimParam$resetPed` +## Method `SimParam$resetPed()` ## ------------------------------------------------ #Create founder haplotypes @@ -71,7 +71,7 @@ pop2 = newPop(founderPop, simParam=SP) pop2@id # 1:10 ## ------------------------------------------------ -## Method `SimParam$restrSegSites` +## Method `SimParam$restrSegSites()` ## ------------------------------------------------ #Create founder haplotypes @@ -83,7 +83,7 @@ SP = SimParam$new(founderPop) SP$restrSegSites(minQtlPerChr=5, minSnpPerChr=5) ## ------------------------------------------------ -## Method `SimParam$setSexes` +## Method `SimParam$setSexes()` ## ------------------------------------------------ #Create founder haplotypes @@ -95,7 +95,7 @@ SP = SimParam$new(founderPop) SP$setSexes("yes_sys") ## ------------------------------------------------ -## Method `SimParam$addSnpChip` +## Method `SimParam$addSnpChip()` ## ------------------------------------------------ #Create founder haplotypes @@ -107,7 +107,7 @@ SP = SimParam$new(founderPop) SP$addSnpChip(10) ## ------------------------------------------------ -## Method `SimParam$addSnpChipByName` +## Method `SimParam$addSnpChipByName()` ## ------------------------------------------------ #Create founder haplotypes @@ -119,7 +119,7 @@ SP = SimParam$new(founderPop) SP$addSnpChipByName(c("1_1","1_3")) ## ------------------------------------------------ -## Method `SimParam$addTraitA` +## Method `SimParam$addTraitA()` ## ------------------------------------------------ #Create founder haplotypes @@ -131,7 +131,7 @@ SP = SimParam$new(founderPop) SP$addTraitA(10) ## ------------------------------------------------ -## Method `SimParam$addTraitAD` +## Method `SimParam$addTraitAD()` ## ------------------------------------------------ #Create founder haplotypes @@ -143,7 +143,7 @@ SP = SimParam$new(founderPop) SP$addTraitAD(10, meanDD=0.5) ## ------------------------------------------------ -## Method `SimParam$altAddTraitAD` +## Method `SimParam$altAddTraitAD()` ## ------------------------------------------------ #Create founder haplotypes @@ -155,7 +155,7 @@ SP = SimParam$new(founderPop) SP$altAddTraitAD(nQtlPerChr=10, mean=0, varA=1, varD=0.05, inbrDepr=0.2) ## ------------------------------------------------ -## Method `SimParam$addTraitAG` +## Method `SimParam$addTraitAG()` ## ------------------------------------------------ #Create founder haplotypes @@ -167,7 +167,7 @@ SP = SimParam$new(founderPop) SP$addTraitAG(10, varGxE=2) ## ------------------------------------------------ -## Method `SimParam$addTraitADG` +## Method `SimParam$addTraitADG()` ## ------------------------------------------------ #Create founder haplotypes @@ -179,7 +179,7 @@ SP = SimParam$new(founderPop) SP$addTraitADG(10, meanDD=0.5, varGxE=2) ## ------------------------------------------------ -## Method `SimParam$addTraitAE` +## Method `SimParam$addTraitAE()` ## ------------------------------------------------ #Create founder haplotypes @@ -191,7 +191,7 @@ SP = SimParam$new(founderPop) SP$addTraitAE(10, relAA=0.1) ## ------------------------------------------------ -## Method `SimParam$addTraitADE` +## Method `SimParam$addTraitADE()` ## ------------------------------------------------ #Create founder haplotypes @@ -203,7 +203,7 @@ SP = SimParam$new(founderPop) SP$addTraitADE(10) ## ------------------------------------------------ -## Method `SimParam$addTraitAEG` +## Method `SimParam$addTraitAEG()` ## ------------------------------------------------ #Create founder haplotypes @@ -215,7 +215,7 @@ SP = SimParam$new(founderPop) SP$addTraitAEG(10, varGxE=2) ## ------------------------------------------------ -## Method `SimParam$addTraitADEG` +## Method `SimParam$addTraitADEG()` ## ------------------------------------------------ #Create founder haplotypes @@ -227,7 +227,7 @@ SP = SimParam$new(founderPop) SP$addTraitADEG(10, meanDD=0.5, varGxE=2) ## ------------------------------------------------ -## Method `SimParam$setVarE` +## Method `SimParam$setVarE()` ## ------------------------------------------------ #Create founder haplotypes @@ -240,7 +240,7 @@ SP$addTraitA(10) SP$setVarE(h2=0.5) ## ------------------------------------------------ -## Method `SimParam$setCorE` +## Method `SimParam$setCorE()` ## ------------------------------------------------ #Create founder haplotypes @@ -255,7 +255,7 @@ E = 0.5*diag(2)+0.5 #Positively correlated error SP$setCorE(E) ## ------------------------------------------------ -## Method `SimParam$rescaleTraits` +## Method `SimParam$rescaleTraits()` ## ------------------------------------------------ #Create founder haplotypes @@ -277,7 +277,7 @@ pop = resetPop(pop, simParam=SP) meanG(pop) ## ------------------------------------------------ -## Method `SimParam$setRecombRatio` +## Method `SimParam$setRecombRatio()` ## ------------------------------------------------ #Create founder haplotypes @@ -289,17 +289,17 @@ SP = SimParam$new(founderPop) SP$setRecombRatio(2) #Twice as much recombination in females } \section{Public fields}{ -\if{html}{\out{
}} -\describe{ -\item{\code{snpChips}}{list of SNP chips} + \if{html}{\out{
}} + \describe{ + \item{\code{snpChips}}{list of SNP chips} -\item{\code{invalidQtl}}{list of segregating sites that aren't valid QTL} + \item{\code{invalidQtl}}{list of segregating sites that aren't valid QTL} -\item{\code{invalidSnp}}{list of segregating sites that aren't valid SNP} + \item{\code{invalidSnp}}{list of segregating sites that aren't valid SNP} -\item{\code{founderPop}}{founder population used for variance scaling} + \item{\code{founderPop}}{founder population used for variance scaling} -\item{\code{finalizePop}}{function applied to newly created populations. + \item{\code{finalizePop}}{function applied to newly created populations. It is run as the last step of creating a new population and provides a way to automatically modify populations. The function must satisfy these four requirements: @@ -311,10 +311,10 @@ The function must satisfy these four requirements: 4) The return is a \code{\link{Pop-class}} object. See \code{\link{asCategorical}} for an example.} -\item{\code{allowEmptyPop}}{if true, population arguments with nInd=0 will + \item{\code{allowEmptyPop}}{if true, population arguments with nInd=0 will return an empty population with a warning instead of an error.} -\item{\code{finalizePheno}}{function applied to newly generated phenotype values. + \item{\code{finalizePheno}}{function applied to newly generated phenotype values. It is run as the last step of generating new phenotype values and provides a way to automatically modify phenotype values. While \code{finalizePop} function can be used to a similar effect, @@ -336,135 +336,135 @@ The function must satisfy these five requirements: row order. See \code{\link{asCategorical}} for an example.} -\item{\code{v}}{the crossover interference parameter for a gamma model of + \item{\code{v}}{the crossover interference parameter for a gamma model of recombination. A value of 1 indicates no crossover interference (e.g. Haldane mapping function). A value of 2.6 approximates the degree of crossover interference implied by the Kosambi mapping function. (default is 2.6)} -\item{\code{p}}{the proportion of crossovers coming from a non-interfering + \item{\code{p}}{the proportion of crossovers coming from a non-interfering pathway. (default is 0)} -\item{\code{quadProb}}{the probability of quadrivalent pairing in an + \item{\code{quadProb}}{the probability of quadrivalent pairing in an autopolyploid. (default is 0)} -} -\if{html}{\out{
}} + } + \if{html}{\out{
}} } \section{Active bindings}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nThreads}}{number of threads used with OpenMP (when available). + \if{html}{\out{
}} + \describe{ + \item{\code{nThreads}}{number of threads used with OpenMP (when available). Assign \code{NULL} to reset to \code{getNumThreads()}. See \code{vignette("parallelization", package="AlphaSimR")} for setup details.} -\item{\code{traitNames}}{vector of trait names} + \item{\code{traitNames}}{vector of trait names} -\item{\code{snpChipNames}}{vector of chip names} + \item{\code{snpChipNames}}{vector of chip names} -\item{\code{traits}}{list of traits} + \item{\code{traits}}{list of traits} -\item{\code{nChr}}{number of chromosomes} + \item{\code{nChr}}{number of chromosomes} -\item{\code{nTraits}}{number of traits} + \item{\code{nTraits}}{number of traits} -\item{\code{nSnpChips}}{number of SNP chips} + \item{\code{nSnpChips}}{number of SNP chips} -\item{\code{segSites}}{segregating sites per chromosome} + \item{\code{segSites}}{segregating sites per chromosome} -\item{\code{sexes}}{sexes used for mating} + \item{\code{sexes}}{sexes used for mating} -\item{\code{sepMap}}{are there seperate genetic maps for + \item{\code{sepMap}}{are there seperate genetic maps for males and females} -\item{\code{genMap}}{list of chromosome genetic maps} + \item{\code{genMap}}{list of chromosome genetic maps} -\item{\code{femaleMap}}{list of chromosome genetic maps for + \item{\code{femaleMap}}{list of chromosome genetic maps for females} -\item{\code{maleMap}}{list of chromosome genetic maps for + \item{\code{maleMap}}{list of chromosome genetic maps for males} -\item{\code{centromere}}{position of centromeres genetic map} + \item{\code{centromere}}{position of centromeres genetic map} -\item{\code{femaleCentromere}}{position of centromeres on female + \item{\code{femaleCentromere}}{position of centromeres on female genetic map} -\item{\code{maleCentromere}}{position of centromeres on male + \item{\code{maleCentromere}}{position of centromeres on male genetic map} -\item{\code{lastId}}{last ID number assigned} + \item{\code{lastId}}{last ID number assigned} -\item{\code{isTrackPed}}{is pedigree being tracked} + \item{\code{isTrackPed}}{is pedigree being tracked} -\item{\code{pedigree}}{pedigree matrix for all individuals} + \item{\code{pedigree}}{pedigree matrix for all individuals} -\item{\code{isTrackRec}}{is recombination being tracked} + \item{\code{isTrackRec}}{is recombination being tracked} -\item{\code{recHist}}{list of historic recombination events} + \item{\code{recHist}}{list of historic recombination events} -\item{\code{haplotypes}}{list of computed IBD haplotypes} + \item{\code{haplotypes}}{list of computed IBD haplotypes} -\item{\code{varA}}{additive genetic variance in founderPop} + \item{\code{varA}}{additive genetic variance in founderPop} -\item{\code{varG}}{total genetic variance in founderPop} + \item{\code{varG}}{total genetic variance in founderPop} -\item{\code{varE}}{default error variance} + \item{\code{varE}}{default error variance} -\item{\code{version}}{the version of AlphaSimR used to generate this object} + \item{\code{version}}{the version of AlphaSimR used to generate this object} -\item{\code{activeQtl}}{a LociMap representing all active QTL in simulation} + \item{\code{activeQtl}}{a LociMap representing all active QTL in simulation} -\item{\code{qtlIndex}}{a list of vectors giving trait specific QTL indices + \item{\code{qtlIndex}}{a list of vectors giving trait specific QTL indices relative to all active QTL} -} -\if{html}{\out{
}} + } + \if{html}{\out{
}} } \section{Methods}{ \subsection{Public methods}{ -\itemize{ -\item \href{#method-SimParam-new}{\code{SimParam$new()}} -\item \href{#method-SimParam-setTrackPed}{\code{SimParam$setTrackPed()}} -\item \href{#method-SimParam-setTrackRec}{\code{SimParam$setTrackRec()}} -\item \href{#method-SimParam-resetPed}{\code{SimParam$resetPed()}} -\item \href{#method-SimParam-restrSegSites}{\code{SimParam$restrSegSites()}} -\item \href{#method-SimParam-setSexes}{\code{SimParam$setSexes()}} -\item \href{#method-SimParam-setFounderHap}{\code{SimParam$setFounderHap()}} -\item \href{#method-SimParam-addSnpChip}{\code{SimParam$addSnpChip()}} -\item \href{#method-SimParam-addSnpChipByName}{\code{SimParam$addSnpChipByName()}} -\item \href{#method-SimParam-addStructuredSnpChip}{\code{SimParam$addStructuredSnpChip()}} -\item \href{#method-SimParam-addTraitA}{\code{SimParam$addTraitA()}} -\item \href{#method-SimParam-addTraitAD}{\code{SimParam$addTraitAD()}} -\item \href{#method-SimParam-altAddTraitAD}{\code{SimParam$altAddTraitAD()}} -\item \href{#method-SimParam-addTraitAG}{\code{SimParam$addTraitAG()}} -\item \href{#method-SimParam-addTraitADG}{\code{SimParam$addTraitADG()}} -\item \href{#method-SimParam-addTraitAE}{\code{SimParam$addTraitAE()}} -\item \href{#method-SimParam-addTraitADE}{\code{SimParam$addTraitADE()}} -\item \href{#method-SimParam-addTraitAEG}{\code{SimParam$addTraitAEG()}} -\item \href{#method-SimParam-addTraitADEG}{\code{SimParam$addTraitADEG()}} -\item \href{#method-SimParam-manAddTrait}{\code{SimParam$manAddTrait()}} -\item \href{#method-SimParam-importTrait}{\code{SimParam$importTrait()}} -\item \href{#method-SimParam-switchTrait}{\code{SimParam$switchTrait()}} -\item \href{#method-SimParam-removeTrait}{\code{SimParam$removeTrait()}} -\item \href{#method-SimParam-setVarE}{\code{SimParam$setVarE()}} -\item \href{#method-SimParam-setCorE}{\code{SimParam$setCorE()}} -\item \href{#method-SimParam-rescaleTraits}{\code{SimParam$rescaleTraits()}} -\item \href{#method-SimParam-setRecombRatio}{\code{SimParam$setRecombRatio()}} -\item \href{#method-SimParam-switchGenMap}{\code{SimParam$switchGenMap()}} -\item \href{#method-SimParam-switchFemaleMap}{\code{SimParam$switchFemaleMap()}} -\item \href{#method-SimParam-switchMaleMap}{\code{SimParam$switchMaleMap()}} -\item \href{#method-SimParam-addToRec}{\code{SimParam$addToRec()}} -\item \href{#method-SimParam-ibdHaplo}{\code{SimParam$ibdHaplo()}} -\item \href{#method-SimParam-updateLastId}{\code{SimParam$updateLastId()}} -\item \href{#method-SimParam-addToPed}{\code{SimParam$addToPed()}} -\item \href{#method-SimParam-clone}{\code{SimParam$clone()}} -} + \itemize{ + \item \href{#method-SimParam-initialize}{\code{SimParam$new()}} + \item \href{#method-SimParam-setTrackPed}{\code{SimParam$setTrackPed()}} + \item \href{#method-SimParam-setTrackRec}{\code{SimParam$setTrackRec()}} + \item \href{#method-SimParam-resetPed}{\code{SimParam$resetPed()}} + \item \href{#method-SimParam-restrSegSites}{\code{SimParam$restrSegSites()}} + \item \href{#method-SimParam-setSexes}{\code{SimParam$setSexes()}} + \item \href{#method-SimParam-setFounderHap}{\code{SimParam$setFounderHap()}} + \item \href{#method-SimParam-addSnpChip}{\code{SimParam$addSnpChip()}} + \item \href{#method-SimParam-addSnpChipByName}{\code{SimParam$addSnpChipByName()}} + \item \href{#method-SimParam-addStructuredSnpChip}{\code{SimParam$addStructuredSnpChip()}} + \item \href{#method-SimParam-addTraitA}{\code{SimParam$addTraitA()}} + \item \href{#method-SimParam-addTraitAD}{\code{SimParam$addTraitAD()}} + \item \href{#method-SimParam-altAddTraitAD}{\code{SimParam$altAddTraitAD()}} + \item \href{#method-SimParam-addTraitAG}{\code{SimParam$addTraitAG()}} + \item \href{#method-SimParam-addTraitADG}{\code{SimParam$addTraitADG()}} + \item \href{#method-SimParam-addTraitAE}{\code{SimParam$addTraitAE()}} + \item \href{#method-SimParam-addTraitADE}{\code{SimParam$addTraitADE()}} + \item \href{#method-SimParam-addTraitAEG}{\code{SimParam$addTraitAEG()}} + \item \href{#method-SimParam-addTraitADEG}{\code{SimParam$addTraitADEG()}} + \item \href{#method-SimParam-manAddTrait}{\code{SimParam$manAddTrait()}} + \item \href{#method-SimParam-importTrait}{\code{SimParam$importTrait()}} + \item \href{#method-SimParam-switchTrait}{\code{SimParam$switchTrait()}} + \item \href{#method-SimParam-removeTrait}{\code{SimParam$removeTrait()}} + \item \href{#method-SimParam-setVarE}{\code{SimParam$setVarE()}} + \item \href{#method-SimParam-setCorE}{\code{SimParam$setCorE()}} + \item \href{#method-SimParam-rescaleTraits}{\code{SimParam$rescaleTraits()}} + \item \href{#method-SimParam-setRecombRatio}{\code{SimParam$setRecombRatio()}} + \item \href{#method-SimParam-switchGenMap}{\code{SimParam$switchGenMap()}} + \item \href{#method-SimParam-switchFemaleMap}{\code{SimParam$switchFemaleMap()}} + \item \href{#method-SimParam-switchMaleMap}{\code{SimParam$switchMaleMap()}} + \item \href{#method-SimParam-addToRec}{\code{SimParam$addToRec()}} + \item \href{#method-SimParam-ibdHaplo}{\code{SimParam$ibdHaplo()}} + \item \href{#method-SimParam-updateLastId}{\code{SimParam$updateLastId()}} + \item \href{#method-SimParam-addToPed}{\code{SimParam$addToPed()}} + \item \href{#method-SimParam-clone}{\code{SimParam$clone()}} + } } \if{html}{\out{
}} -\if{html}{\out{}} -\if{latex}{\out{\hypertarget{method-SimParam-new}{}}} -\subsection{Method \code{new()}}{ -Starts the process of building a new simulation +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-SimParam-initialize}{}}} +\subsection{\code{SimParam$new()}}{ + Starts the process of building a new simulation by creating a new \code{SimParam} object and assigning a founder population to the class. It is recommended that you save the object with the name \code{SP}, because subsequent functions will @@ -472,136 +472,132 @@ check your global environment for an object of this name if their \code{simParam} arguments are \code{NULL}. This allows you to call these functions without explicitly supplying a \code{simParam} argument with every call. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$new(founderPop)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{founderPop}}{an object of \code{\link{MapPop-class}}} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$new(founderPop)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{founderPop}}{an object of \code{\link{MapPop-class}}} + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setTrackPed}{}}} -\subsection{Method \code{setTrackPed()}}{ -Sets pedigree tracking for the simulation. +\subsection{\code{SimParam$setTrackPed()}}{ + Sets pedigree tracking for the simulation. By default pedigree tracking is turned off. When turned on, the pedigree of all individuals created will be tracked, except those created by \code{\link{hybridCross}}. Turning off pedigree tracking will turn off recombination tracking if it is turned on. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setTrackPed(isTrackPed, force = FALSE)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{isTrackPed}}{should pedigree tracking be on.} - -\item{\code{force}}{should the check for a running simulation be + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setTrackPed(isTrackPed, force = FALSE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{isTrackPed}}{should pedigree tracking be on.} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$setTrackPed(TRUE) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setTrackRec}{}}} -\subsection{Method \code{setTrackRec()}}{ -Sets recombination tracking for the simulation. +\subsection{\code{SimParam$setTrackRec()}}{ + Sets recombination tracking for the simulation. By default recombination tracking is turned off. When turned on recombination tracking will also turn on pedigree tracking. Recombination tracking keeps records of all individuals created, except those created by \code{\link{hybridCross}}, because their pedigree is not tracked. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setTrackRec(isTrackRec, force = FALSE)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{isTrackRec}}{should recombination tracking be on.} - -\item{\code{force}}{should the check for a running simulation be + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setTrackRec(isTrackRec, force = FALSE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{isTrackRec}}{should recombination tracking be on.} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$setTrackRec(TRUE) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-resetPed}{}}} -\subsection{Method \code{resetPed()}}{ -Resets the internal lastId, the pedigree +\subsection{\code{SimParam$resetPed()}}{ + Resets the internal lastId, the pedigree and recombination tracking (if in use) to the supplied lastId. Be careful using this function because it may introduce a bug if you use individuals from the deleted portion of the pedigree. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$resetPed(lastId = 0L)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{lastId}}{last ID to include in pedigree} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$resetPed(lastId = 0L)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{lastId}}{last ID to include in pedigree} + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} #Create population pop = newPop(founderPop, simParam=SP) @@ -612,16 +608,15 @@ SP$resetPed() pop2 = newPop(founderPop, simParam=SP) pop2@id # 1:10 } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-restrSegSites}{}}} -\subsection{Method \code{restrSegSites()}}{ -Sets restrictions on which segregating sites can +\subsection{\code{SimParam$restrSegSites()}}{ + Sets restrictions on which segregating sites can serve as a SNP and/or QTL. The default behavior of AlphaSimR is to randomly sample QTL or SNP from all eligible sites and then mark the sampled sites ineligible to be sampled as the other @@ -641,59 +636,53 @@ overlap=FALSE to preallocate sites as QTL and SNP respectively. This option is useful when simulating multiple traits and/or SNP chips, because it can be used to guarantee that enough eligible sites are available when running addTrait and or addSnpChip functions. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$restrSegSites( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$restrSegSites( minQtlPerChr = NULL, minSnpPerChr = NULL, excludeQtl = NULL, excludeSnp = NULL, overlap = FALSE, minSnpFreq = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{minQtlPerChr}}{the minimum number of segregating sites for +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{minQtlPerChr}}{the minimum number of segregating sites for QTLs. Can be a single value or a vector values for each chromosome.} - -\item{\code{minSnpPerChr}}{the minimum number of segregating sites for SNPs. + \item{\code{minSnpPerChr}}{the minimum number of segregating sites for SNPs. Can be a single value or a vector values for each chromosome.} - -\item{\code{excludeQtl}}{an optional vector of segregating site names to + \item{\code{excludeQtl}}{an optional vector of segregating site names to exclude from consideration as a viable QTL.} - -\item{\code{excludeSnp}}{an optional vector of segregating site names to + \item{\code{excludeSnp}}{an optional vector of segregating site names to exclude from consideration as a viable SNP.} - -\item{\code{overlap}}{should SNP and QTL sites be allowed to overlap.} - -\item{\code{minSnpFreq}}{minimum allowable frequency for SNP loci. + \item{\code{overlap}}{should SNP and QTL sites be allowed to overlap.} + \item{\code{minSnpFreq}}{minimum allowable frequency for SNP loci. No minimum SNP frequency is used if value is NULL.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$restrSegSites(minQtlPerChr=5, minSnpPerChr=5) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setSexes}{}}} -\subsection{Method \code{setSexes()}}{ -Changes how sexes are determined in the simulation. +\subsection{\code{SimParam$setSexes()}}{ + Changes how sexes are determined in the simulation. The default sexes is "no", indicating all individuals are hermaphrodites. To add sexes to the simulation, run this function with "yes_sys" or "yes_rand". The value "yes_sys" will systematically assign @@ -701,165 +690,160 @@ sexes to newly created individuals as first male and then female. Populations with an odd number of individuals will have one more male than female. The value "yes_rand" will randomly assign a sex to each individual. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setSexes(sexes, force = FALSE)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{sexes}}{acceptable value are "no", "yes_sys", or + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setSexes(sexes, force = FALSE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{sexes}}{acceptable value are "no", "yes_sys", or "yes_rand"} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$setSexes("yes_sys") } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setFounderHap}{}}} -\subsection{Method \code{setFounderHap()}}{ -Allows for the manual setting of founder haplotypes. This functionality +\subsection{\code{SimParam$setFounderHap()}}{ + Allows for the manual setting of founder haplotypes. This functionality is not fully documented, because it is still experimental. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setFounderHap(hapMap)}\if{html}{\out{
}} + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setFounderHap(hapMap)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{hapMap}}{a list of founder haplotypes} + } + \if{html}{\out{
}} + } } -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{hapMap}}{a list of founder haplotypes} -} -\if{html}{\out{
}} -} -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addSnpChip}{}}} -\subsection{Method \code{addSnpChip()}}{ -Randomly assigns eligible SNPs to a SNP chip -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addSnpChip(nSnpPerChr, minSnpFreq = NULL, refPop = NULL, name = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nSnpPerChr}}{number of SNPs per chromosome. +\subsection{\code{SimParam$addSnpChip()}}{ + Randomly assigns eligible SNPs to a SNP chip + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addSnpChip(nSnpPerChr, minSnpFreq = NULL, refPop = NULL, name = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nSnpPerChr}}{number of SNPs per chromosome. Can be a single value or nChr values.} - -\item{\code{minSnpFreq}}{minimum allowable frequency for SNP loci. + \item{\code{minSnpFreq}}{minimum allowable frequency for SNP loci. If NULL, no minimum frequency is used.} - -\item{\code{refPop}}{reference population for calculating SNP + \item{\code{refPop}}{reference population for calculating SNP frequency. If NULL, the founder population is used.} - -\item{\code{name}}{optional name for chip} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + \item{\code{name}}{optional name for chip} + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addSnpChip(10) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addSnpChipByName}{}}} -\subsection{Method \code{addSnpChipByName()}}{ -Assigns SNPs to a SNP chip by supplying marker names. This function does +\subsection{\code{SimParam$addSnpChipByName()}}{ + Assigns SNPs to a SNP chip by supplying marker names. This function does check against excluded SNPs and will not add the SNPs to the list of excluded QTL for the purpose of avoiding overlap between SNPs and QTL. Excluding these SNPs from being used as QTL can be accomplished using the excludeQtl argument in SimParam's restrSegSites function. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addSnpChipByName(markers, name = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{markers}}{a vector of names for the markers} - -\item{\code{name}}{optional name for chip} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addSnpChipByName(markers, name = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{markers}}{a vector of names for the markers} + \item{\code{name}}{optional name for chip} + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addSnpChipByName(c("1_1","1_3")) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addStructuredSnpChip}{}}} -\subsection{Method \code{addStructuredSnpChip()}}{ -Randomly selects the number of snps in structure and then +\subsection{\code{SimParam$addStructuredSnpChip()}}{ + Randomly selects the number of snps in structure and then assigns them to chips based on structure -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addStructuredSnpChip(nSnpPerChr, structure, force = FALSE)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nSnpPerChr}}{number of SNPs per chromosome. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addStructuredSnpChip(nSnpPerChr, structure, force = FALSE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nSnpPerChr}}{number of SNPs per chromosome. Can be a single value or nChr values.} - -\item{\code{structure}}{a matrix. Rows are snp chips, columns are chips. + \item{\code{structure}}{a matrix. Rows are snp chips, columns are chips. If value is true then that snp is on that chip.} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitA}{}}} -\subsection{Method \code{addTraitA()}}{ -Randomly assigns eligible QTLs for one or more additive traits. +\subsection{\code{SimParam$addTraitA()}}{ + Randomly assigns eligible QTLs for one or more additive traits. If simulating more than one trait, all traits will be pleiotropic with correlated additive effects. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitA( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitA( nQtlPerChr, mean = 0, var = 1, @@ -869,60 +853,51 @@ with correlated additive effects. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitA(10) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitAD}{}}} -\subsection{Method \code{addTraitAD()}}{ -Randomly assigns eligible QTLs for one or more traits with dominance. +\subsection{\code{SimParam$addTraitAD()}}{ + Randomly assigns eligible QTLs for one or more traits with dominance. If simulating more than one trait, all traits will be pleiotropic with correlated effects. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitAD( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitAD( nQtlPerChr, mean = 0, var = 1, @@ -936,70 +911,57 @@ with correlated effects. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{meanDD}}{mean dominance degree} - -\item{\code{varDD}}{variance of dominance degree} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corDD}}{a matrix of correlations between dominance degrees} - -\item{\code{useVarA}}{tune according to additive genetic variance if true. If +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{meanDD}}{mean dominance degree} + \item{\code{varDD}}{variance of dominance degree} + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corDD}}{a matrix of correlations between dominance degrees} + \item{\code{useVarA}}{tune according to additive genetic variance if true. If FALSE, tuning is performed according to total genetic variance.} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitAD(10, meanDD=0.5) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-altAddTraitAD}{}}} -\subsection{Method \code{altAddTraitAD()}}{ -An alternative method for adding a trait with additive and dominance effects +\subsection{\code{SimParam$altAddTraitAD()}}{ + An alternative method for adding a trait with additive and dominance effects to an AlphaSimR simulation. The function attempts to create a trait matching user defined values for number of QTL, inbreeding depression, additive genetic variance and dominance genetic variance. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$altAddTraitAD( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$altAddTraitAD( nQtlPerChr, mean = 0, varA = 1, @@ -1011,41 +973,31 @@ variance and dominance genetic variance. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{desired mean of the trait} - -\item{\code{varA}}{desired additive variance} - -\item{\code{varD}}{desired dominance variance} - -\item{\code{inbrDepr}}{desired inbreeding depression, see details} - -\item{\code{limMeanDD}}{limits for meanDD, see details} - -\item{\code{limVarDD}}{limits for varDD, see details} - -\item{\code{silent}}{should summary details be printed to the console} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{mean}}{desired mean of the trait} + \item{\code{varA}}{desired additive variance} + \item{\code{varD}}{desired dominance variance} + \item{\code{inbrDepr}}{desired inbreeding depression, see details} + \item{\code{limMeanDD}}{limits for meanDD, see details} + \item{\code{limVarDD}}{limits for varDD, see details} + \item{\code{silent}}{should summary details be printed to the console} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Details}{ -This function will always add a trait to 'SimParam', unless an error occurs + } + \if{html}{\out{
}} + } + \subsection{Details}{ + This function will always add a trait to 'SimParam', unless an error occurs with picking QTLs. The resulting trait will always have the desired mean and additive genetic variance. However, it may not have the desired values for inbreeding depression and dominance variance. Thus, it is strongly recommended @@ -1078,32 +1030,30 @@ Summary information on this trait is printed to the console when silent=FALSE. The summary information reports the inbreeding depression and dominance variance for the population as well as the dominance degree mean and variance applied to the trait. -} - -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$altAddTraitAD(nQtlPerChr=10, mean=0, varA=1, varD=0.05, inbrDepr=0.2) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitAG}{}}} -\subsection{Method \code{addTraitAG()}}{ -Randomly assigns eligible QTLs for one or more additive GxE traits. +\subsection{\code{SimParam$addTraitAG()}}{ + Randomly assigns eligible QTLs for one or more additive GxE traits. If simulating more than one trait, all traits will be pleiotropic with correlated effects. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitAG( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitAG( nQtlPerChr, mean = 0, var = 1, @@ -1116,64 +1066,52 @@ with correlated effects. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} - -\item{\code{varEnv}}{a vector of environmental variances for one or more traits} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corGxE}}{a matrix of correlations between GxE effects} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} + \item{\code{varEnv}}{a vector of environmental variances for one or more traits} + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corGxE}}{a matrix of correlations between GxE effects} + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitAG(10, varGxE=2) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitADG}{}}} -\subsection{Method \code{addTraitADG()}}{ -Randomly assigns eligible QTLs for a trait with dominance and GxE. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitADG( +\subsection{\code{SimParam$addTraitADG()}}{ + Randomly assigns eligible QTLs for a trait with dominance and GxE. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitADG( nQtlPerChr, mean = 0, var = 1, @@ -1190,75 +1128,59 @@ Randomly assigns eligible QTLs for a trait with dominance and GxE. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{varEnv}}{a vector of environmental variances for one or more traits} - -\item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} - -\item{\code{meanDD}}{mean dominance degree} - -\item{\code{varDD}}{variance of dominance degree} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corDD}}{a matrix of correlations between dominance degrees} - -\item{\code{corGxE}}{a matrix of correlations between GxE effects} - -\item{\code{useVarA}}{tune according to additive genetic variance if true} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{varEnv}}{a vector of environmental variances for one or more traits} + \item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} + \item{\code{meanDD}}{mean dominance degree} + \item{\code{varDD}}{variance of dominance degree} + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corDD}}{a matrix of correlations between dominance degrees} + \item{\code{corGxE}}{a matrix of correlations between GxE effects} + \item{\code{useVarA}}{tune according to additive genetic variance if true} + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitADG(10, meanDD=0.5, varGxE=2) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitAE}{}}} -\subsection{Method \code{addTraitAE()}}{ -Randomly assigns eligible QTLs for one or more additive and epistasis +\subsection{\code{SimParam$addTraitAE()}}{ + Randomly assigns eligible QTLs for one or more additive and epistasis traits. If simulating more than one trait, all traits will be pleiotropic with correlated additive effects. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitAE( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitAE( nQtlPerChr, mean = 0, var = 1, @@ -1271,68 +1193,56 @@ with correlated additive effects. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{relAA}}{the relative value of additive-by-additive variance compared +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{relAA}}{the relative value of additive-by-additive variance compared to additive variance in a diploid organism with allele frequency 0.5} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} - -\item{\code{useVarA}}{tune according to additive genetic variance if true. If + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} + \item{\code{useVarA}}{tune according to additive genetic variance if true. If FALSE, tuning is performed according to total genetic variance.} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitAE(10, relAA=0.1) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitADE}{}}} -\subsection{Method \code{addTraitADE()}}{ -Randomly assigns eligible QTLs for one or more traits with dominance and +\subsection{\code{SimParam$addTraitADE()}}{ + Randomly assigns eligible QTLs for one or more traits with dominance and epistasis. If simulating more than one trait, all traits will be pleiotropic with correlated effects. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitADE( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitADE( nQtlPerChr, mean = 0, var = 1, @@ -1348,74 +1258,59 @@ with correlated effects. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{meanDD}}{mean dominance degree} - -\item{\code{varDD}}{variance of dominance degree} - -\item{\code{relAA}}{the relative value of additive-by-additive variance compared +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{meanDD}}{mean dominance degree} + \item{\code{varDD}}{variance of dominance degree} + \item{\code{relAA}}{the relative value of additive-by-additive variance compared to additive variance in a diploid organism with allele frequency 0.5} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corDD}}{a matrix of correlations between dominance degrees} - -\item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} - -\item{\code{useVarA}}{tune according to additive genetic variance if true. If + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corDD}}{a matrix of correlations between dominance degrees} + \item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} + \item{\code{useVarA}}{tune according to additive genetic variance if true. If FALSE, tuning is performed according to total genetic variance.} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitADE(10) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitAEG}{}}} -\subsection{Method \code{addTraitAEG()}}{ -Randomly assigns eligible QTLs for one or more additive and epistasis +\subsection{\code{SimParam$addTraitAEG()}}{ + Randomly assigns eligible QTLs for one or more additive and epistasis GxE traits. If simulating more than one trait, all traits will be pleiotropic with correlated effects. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitAEG( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitAEG( nQtlPerChr, mean = 0, var = 1, @@ -1431,73 +1326,58 @@ with correlated effects. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{relAA}}{the relative value of additive-by-additive variance compared +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{relAA}}{the relative value of additive-by-additive variance compared to additive variance in a diploid organism with allele frequency 0.5} - -\item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} - -\item{\code{varEnv}}{a vector of environmental variances for one or more traits} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} - -\item{\code{corGxE}}{a matrix of correlations between GxE effects} - -\item{\code{useVarA}}{tune according to additive genetic variance if true. If + \item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} + \item{\code{varEnv}}{a vector of environmental variances for one or more traits} + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} + \item{\code{corGxE}}{a matrix of correlations between GxE effects} + \item{\code{useVarA}}{tune according to additive genetic variance if true. If FALSE, tuning is performed according to total genetic variance.} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitAEG(10, varGxE=2) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addTraitADEG}{}}} -\subsection{Method \code{addTraitADEG()}}{ -Randomly assigns eligible QTLs for a trait with dominance, +\subsection{\code{SimParam$addTraitADEG()}}{ + Randomly assigns eligible QTLs for a trait with dominance, epistasis and GxE. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addTraitADEG( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addTraitADEG( nQtlPerChr, mean = 0, var = 1, @@ -1516,111 +1396,92 @@ epistasis and GxE. force = FALSE, name = NULL, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{nQtlPerChr}}{number of QTLs per chromosome. Can be a single value or nChr values.} - -\item{\code{mean}}{a vector of desired mean genetic values for one or more traits} - -\item{\code{var}}{a vector of desired genetic variances for one or more traits} - -\item{\code{varEnv}}{a vector of environmental variances for one or more traits} - -\item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} - -\item{\code{meanDD}}{mean dominance degree} - -\item{\code{varDD}}{variance of dominance degree} - -\item{\code{relAA}}{the relative value of additive-by-additive variance compared + \item{\code{mean}}{a vector of desired mean genetic values for one or more traits} + \item{\code{var}}{a vector of desired genetic variances for one or more traits} + \item{\code{varEnv}}{a vector of environmental variances for one or more traits} + \item{\code{varGxE}}{a vector of total genotype-by-environment variances for the traits} + \item{\code{meanDD}}{mean dominance degree} + \item{\code{varDD}}{variance of dominance degree} + \item{\code{relAA}}{the relative value of additive-by-additive variance compared to additive variance in a diploid organism with allele frequency 0.5} - -\item{\code{corA}}{a matrix of correlations between additive effects} - -\item{\code{corDD}}{a matrix of correlations between dominance degrees} - -\item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} - -\item{\code{corGxE}}{a matrix of correlations between GxE effects} - -\item{\code{useVarA}}{tune according to additive genetic variance if true} - -\item{\code{gamma}}{should a gamma distribution be used instead of normal} - -\item{\code{shape}}{the shape parameter for the gamma distribution + \item{\code{corA}}{a matrix of correlations between additive effects} + \item{\code{corDD}}{a matrix of correlations between dominance degrees} + \item{\code{corAA}}{a matrix of correlations between additive-by-additive effects} + \item{\code{corGxE}}{a matrix of correlations between GxE effects} + \item{\code{useVarA}}{tune according to additive genetic variance if true} + \item{\code{gamma}}{should a gamma distribution be used instead of normal} + \item{\code{shape}}{the shape parameter for the gamma distribution (the rate/scale parameter of the gamma distribution is accounted for via the desired level of genetic variance, the var argument)} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing.} - -\item{\code{name}}{optional name for trait(s)} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{name}}{optional name for trait(s)} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitADEG(10, meanDD=0.5, varGxE=2) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-manAddTrait}{}}} -\subsection{Method \code{manAddTrait()}}{ -Manually add a new trait to the simulation. Trait must +\subsection{\code{SimParam$manAddTrait()}}{ + Manually add a new trait to the simulation. Trait must be formatted as a \code{\link{LociMap-class}}. If the trait is not already formatted, consider using importTrait. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$manAddTrait(lociMap, varE = NA_real_, force = FALSE, nThreads = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{lociMap}}{a new object descended from + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$manAddTrait(lociMap, varE = NA_real_, force = FALSE, nThreads = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{lociMap}}{a new object descended from \code{\link{LociMap-class}}} - -\item{\code{varE}}{default error variance for phenotype, optional} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{varE}}{default error variance for phenotype, optional} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-importTrait}{}}} -\subsection{Method \code{importTrait()}}{ -Manually add a new trait(s) to the simulation. Unlike the +\subsection{\code{SimParam$importTrait()}}{ + Manually add a new trait(s) to the simulation. Unlike the manAddTrait function, this function does not require formatting the trait as a \code{\link{LociMap-class}}. The formatting is performed automatically for the user, with more user friendly data.frames or matrices taken as inputs. This function only works for A and AD trait types. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$importTrait( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$importTrait( markerNames, addEff, domEff = NULL, @@ -1629,211 +1490,194 @@ inputs. This function only works for A and AD trait types. varE = NULL, force = FALSE, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{markerNames}}{a vector of names for the QTL} - -\item{\code{addEff}}{a matrix of additive effects (nLoci x nTraits). +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{markerNames}}{a vector of names for the QTL} + \item{\code{addEff}}{a matrix of additive effects (nLoci x nTraits). Alternatively, a vector of length nLoci can be supplied for a single trait.} - -\item{\code{domEff}}{optional dominance effects for each locus} - -\item{\code{intercept}}{optional intercepts for each trait} - -\item{\code{name}}{optional name(s) for the trait(s)} - -\item{\code{varE}}{default error variance for phenotype, optional} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{domEff}}{optional dominance effects for each locus} + \item{\code{intercept}}{optional intercepts for each trait} + \item{\code{name}}{optional name(s) for the trait(s)} + \item{\code{varE}}{default error variance for phenotype, optional} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-switchTrait}{}}} -\subsection{Method \code{switchTrait()}}{ -Switch a trait in the simulation. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$switchTrait( +\subsection{\code{SimParam$switchTrait()}}{ + Switch a trait in the simulation. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$switchTrait( traitPos, lociMap, varE = NA_real_, force = FALSE, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{traitPos}}{an integer indicate which trait to switch} - -\item{\code{lociMap}}{a new object descended from +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{traitPos}}{an integer indicate which trait to switch} + \item{\code{lociMap}}{a new object descended from \code{\link{LociMap-class}}} - -\item{\code{varE}}{default error variance for phenotype, optional} - -\item{\code{force}}{should the check for a running simulation be + \item{\code{varE}}{default error variance for phenotype, optional} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-removeTrait}{}}} -\subsection{Method \code{removeTrait()}}{ -Remove a trait from the simulation -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$removeTrait(traits, force = FALSE)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{traits}}{an integer vector indicating which traits to remove} - -\item{\code{force}}{should the check for a running simulation be +\subsection{\code{SimParam$removeTrait()}}{ + Remove a trait from the simulation + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$removeTrait(traits, force = FALSE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{traits}}{an integer vector indicating which traits to remove} + \item{\code{force}}{should the check for a running simulation be ignored. Only set to TRUE if you know what you are doing} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setVarE}{}}} -\subsection{Method \code{setVarE()}}{ -Defines a default values for error +\subsection{\code{SimParam$setVarE()}}{ + Defines a default values for error variances used in \code{\link{setPheno}}. These defaults will be used to automatically generate phenotypes when new populations are created. See the details section of \code{\link{setPheno}} for more information about each arguments and how they should be used. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setVarE(h2 = NULL, H2 = NULL, varE = NULL, corE = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{h2}}{a vector of desired narrow-sense heritabilities} - -\item{\code{H2}}{a vector of desired broad-sense heritabilities} - -\item{\code{varE}}{a vector or matrix of error variances} - -\item{\code{corE}}{an optional matrix of error correlations} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setVarE(h2 = NULL, H2 = NULL, varE = NULL, corE = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{h2}}{a vector of desired narrow-sense heritabilities} + \item{\code{H2}}{a vector of desired broad-sense heritabilities} + \item{\code{varE}}{a vector or matrix of error variances} + \item{\code{corE}}{an optional matrix of error correlations} + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitA(10) SP$setVarE(h2=0.5) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setCorE}{}}} -\subsection{Method \code{setCorE()}}{ -Defines a correlation structure for default +\subsection{\code{SimParam$setCorE()}}{ + Defines a correlation structure for default error variances. You must call \code{setVarE} first to define the default error variances. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setCorE(corE)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{corE}}{a correlation matrix for the error variances} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setCorE(corE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{corE}}{a correlation matrix for the error variances} + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitA(10, mean=c(0,0), var=c(1,1), corA=diag(2)) SP$setVarE(varE=c(1,1)) E = 0.5*diag(2)+0.5 #Positively correlated error SP$setCorE(E) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-rescaleTraits}{}}} -\subsection{Method \code{rescaleTraits()}}{ -Linearly scales all traits to achieve desired +\subsection{\code{SimParam$rescaleTraits()}}{ + Linearly scales all traits to achieve desired values of means and variances in the founder population. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$rescaleTraits( + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$rescaleTraits( mean = 0, var = 1, varEnv = 0, varGxE = 1e-06, useVarA = TRUE, nThreads = NULL -)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{mean}}{a vector of new trait means} - -\item{\code{var}}{a vector of new trait variances} - -\item{\code{varEnv}}{a vector of new environmental variances} - -\item{\code{varGxE}}{a vector of new GxE variances} - -\item{\code{useVarA}}{tune according to additive genetic variance if true} - -\item{\code{nThreads}}{number of threads to use if OpenMP is available. +)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{mean}}{a vector of new trait means} + \item{\code{var}}{a vector of new trait variances} + \item{\code{varEnv}}{a vector of new environmental variances} + \item{\code{varGxE}}{a vector of new GxE variances} + \item{\code{useVarA}}{tune according to additive genetic variance if true} + \item{\code{nThreads}}{number of threads to use if OpenMP is available. If \code{NULL}, the number is obtained from \code{self$nThreads}.} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$addTraitA(10) #Create population @@ -1846,218 +1690,219 @@ SP$rescaleTraits(mean=1) pop = resetPop(pop, simParam=SP) meanG(pop) } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-setRecombRatio}{}}} -\subsection{Method \code{setRecombRatio()}}{ -Set the relative recombination rates between males +\subsection{\code{SimParam$setRecombRatio()}}{ + Set the relative recombination rates between males and females. This allows for sex-specific recombination rates, under the assumption of equivalent recombination landscapes. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$setRecombRatio(femaleRatio)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{femaleRatio}}{relative ratio of recombination in females compared to + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$setRecombRatio(femaleRatio)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{femaleRatio}}{relative ratio of recombination in females compared to males. A value of 2 indicate twice as much recombination in females. The value must be greater than 0. (default is 1)} -} -\if{html}{\out{
}} -} -\subsection{Examples}{ -\if{html}{\out{
}} -\preformatted{#Create founder haplotypes + } + \if{html}{\out{
}} + } + \subsection{Examples}{ + \if{html}{\out{
}} + \preformatted{#Create founder haplotypes founderPop = quickHaplo(nInd=10, nChr=1, segSites=10) #Set simulation parameters SP = SimParam$new(founderPop) -\dontshow{SP$nThreads = 1L} SP$setRecombRatio(2) #Twice as much recombination in females } -\if{html}{\out{
}} - + \if{html}{\out{
}} + } } -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-switchGenMap}{}}} -\subsection{Method \code{switchGenMap()}}{ -Replaces existing genetic map. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$switchGenMap(genMap, centromere = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{genMap}}{a list of length nChr containing +\subsection{\code{SimParam$switchGenMap()}}{ + Replaces existing genetic map. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$switchGenMap(genMap, centromere = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{genMap}}{a list of length nChr containing numeric vectors for the position of each segregating site on a chromosome.} - -\item{\code{centromere}}{a numeric vector of centromere + \item{\code{centromere}}{a numeric vector of centromere positions. If NULL, the centromere are assumed to be metacentric.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-switchFemaleMap}{}}} -\subsection{Method \code{switchFemaleMap()}}{ -Replaces existing female genetic map. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$switchFemaleMap(genMap, centromere = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{genMap}}{a list of length nChr containing +\subsection{\code{SimParam$switchFemaleMap()}}{ + Replaces existing female genetic map. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$switchFemaleMap(genMap, centromere = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{genMap}}{a list of length nChr containing numeric vectors for the position of each segregating site on a chromosome.} - -\item{\code{centromere}}{a numeric vector of centromere + \item{\code{centromere}}{a numeric vector of centromere positions. If NULL, the centromere are assumed to be metacentric.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-switchMaleMap}{}}} -\subsection{Method \code{switchMaleMap()}}{ -Replaces existing male genetic map. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$switchMaleMap(genMap, centromere = NULL)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{genMap}}{a list of length nChr containing +\subsection{\code{SimParam$switchMaleMap()}}{ + Replaces existing male genetic map. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$switchMaleMap(genMap, centromere = NULL)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{genMap}}{a list of length nChr containing numeric vectors for the position of each segregating site on a chromosome.} - -\item{\code{centromere}}{a numeric vector of centromere + \item{\code{centromere}}{a numeric vector of centromere positions. If NULL, the centromere are assumed to be metacentric.} + } + \if{html}{\out{
}} + } } -\if{html}{\out{
}} -} -} + \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addToRec}{}}} -\subsection{Method \code{addToRec()}}{ -For internal use only. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addToRec(lastId, id, mother, father, isDH, hist, ploidy)}\if{html}{\out{
}} +\subsection{\code{SimParam$addToRec()}}{ + For internal use only. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addToRec(lastId, id, mother, father, isDH, hist, ploidy)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{lastId}}{ID of last individual} + \item{\code{id}}{the name of each individual} + \item{\code{mother}}{vector of mother iids} + \item{\code{father}}{vector of father iids} + \item{\code{isDH}}{indicator for DH lines} + \item{\code{hist}}{new recombination history} + \item{\code{ploidy}}{ploidy level} + } + \if{html}{\out{
}} + } } -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{lastId}}{ID of last individual} - -\item{\code{id}}{the name of each individual} - -\item{\code{mother}}{vector of mother iids} - -\item{\code{father}}{vector of father iids} - -\item{\code{isDH}}{indicator for DH lines} - -\item{\code{hist}}{new recombination history} - -\item{\code{ploidy}}{ploidy level} -} -\if{html}{\out{
}} -} -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-ibdHaplo}{}}} -\subsection{Method \code{ibdHaplo()}}{ -For internal use only. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$ibdHaplo(iid)}\if{html}{\out{
}} +\subsection{\code{SimParam$ibdHaplo()}}{ + For internal use only. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$ibdHaplo(iid)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{iid}}{internal ID} + } + \if{html}{\out{
}} + } } -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{iid}}{internal ID} -} -\if{html}{\out{
}} -} -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-updateLastId}{}}} -\subsection{Method \code{updateLastId()}}{ -For internal use only. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$updateLastId(lastId)}\if{html}{\out{
}} +\subsection{\code{SimParam$updateLastId()}}{ + For internal use only. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$updateLastId(lastId)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{lastId}}{last ID assigned} + } + \if{html}{\out{
}} + } } -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{lastId}}{last ID assigned} -} -\if{html}{\out{
}} -} -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-addToPed}{}}} -\subsection{Method \code{addToPed()}}{ -For internal use only. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$addToPed(lastId, id, mother, father, isDH)}\if{html}{\out{
}} +\subsection{\code{SimParam$addToPed()}}{ + For internal use only. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$addToPed(lastId, id, mother, father, isDH)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{lastId}}{ID of last individual} + \item{\code{id}}{the name of each individual} + \item{\code{mother}}{vector of mother iids} + \item{\code{father}}{vector of father iids} + \item{\code{isDH}}{indicator for DH lines} + } + \if{html}{\out{
}} + } } -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{lastId}}{ID of last individual} - -\item{\code{id}}{the name of each individual} - -\item{\code{mother}}{vector of mother iids} - -\item{\code{father}}{vector of father iids} - -\item{\code{isDH}}{indicator for DH lines} -} -\if{html}{\out{
}} -} -} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-SimParam-clone}{}}} -\subsection{Method \code{clone()}}{ -The objects of this class are cloneable with this method. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{SimParam$clone(deep = FALSE)}\if{html}{\out{
}} +\subsection{\code{SimParam$clone()}}{ + The objects of this class are cloneable with this method. + \subsection{Usage}{ + \if{html}{\out{
}} + \preformatted{SimParam$clone(deep = FALSE)} + \if{html}{\out{
}} + } + \subsection{Arguments}{ + \if{html}{\out{
}} + \describe{ + \item{\code{deep}}{Whether to make a deep clone.} + } + \if{html}{\out{
}} + } } -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{deep}}{Whether to make a deep clone.} -} -\if{html}{\out{
}} -} -} } From ffd088a13eab45cdbc599e5d30fd9618b90d6edd Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Wed, 8 Jul 2026 16:58:11 +0200 Subject: [PATCH 7/9] Update `splitPop` to use data frames as grouping specs --- R/mergePops.R | 129 +++++++++++++++++++++++++++++++------------- man/dot-splitPop.Rd | 12 ++++- man/splitPop.Rd | 74 ++++++++++++++++++------- 3 files changed, 159 insertions(+), 56 deletions(-) diff --git a/R/mergePops.R b/R/mergePops.R index 815ec241..b8f24924 100644 --- a/R/mergePops.R +++ b/R/mergePops.R @@ -527,35 +527,55 @@ mergeMultiPops = function(..., level=0){ #' #' @description #' Split a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object -#' into a \code{MultiPop} at specified nesting \code{level}(s) using one -#' or more grouping specifications (\code{by}). Grouping specs can be -#' atomic vectors (length \code{nInd(x)} or 1) or functions that return -#' such vectors. If \code{by} is a list, each element is applied -#' recursively to create nested \code{MultiPop} objects. +#' into a \code{MultiPop} using one or more grouping specifications passed +#' through \code{by}. #' #' @param x a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object -#' @param by a vector, function, or list of vectors/functions defining -#' groupings. Functions are called with a \code{Pop} and must return an -#' atomic vector of length \code{nInd(pop)} (or 1). -#' @param level A positive integer, a vector of positive integers, or +#' @param by a vector, function, list, or data frame defining groupings +#' @param level a positive integer, a vector of positive integers, or #' \code{Inf}. Only relevant when \code{x} is a \code{\link{MultiPop-class}} -#' object. If \code{level = Inf}, split any \code{\link{Pop-class}} within -#' \code{x} according to \code{by}. If \code{level = c(a, b)}, only -#' \code{Pop-class} objects at levels \code{a} and \code{b} (top \code{level=1}) -#' are split. An error is raised if any requested level is deeper than the -#' object's maximum nesting depth. -#' +#' object. See Details. +#' +#' @details +#' The \code{by} argument can take several forms: \cr +#' - Atomic vectors: used directly as grouping labels. The +#' vector is passed to \code{\link[base]{split}} and may be recycled as +#' needed. If the vector length is not a multiple of the population size, +#' \code{\link[base]{split}} issues a warning. \cr +#' - Functions: called on each \code{Pop} object and must return +#' an atomic vector of grouping labels. \cr +#' - Lists: each element must be a vector or function. Elements +#' are applied recursively to create nested \code{\link{MultiPop-class}} +#' objects. \cr +#' - Data frames: each column defines one nesting level. In this +#' case, \code{x} must be a \code{\link{Pop-class}} object. Rows of the data +#' frame must correspond to individuals in \code{x}. If row names are present +#' and match \code{x@id}, they are used to align rows to individuals; +#' otherwise rows are assumed to be in population order and a warning is +#' issued. \code{NA} values are allowed in grouping columns and can be used +#' to represent uneven nesting structures. If all grouping values for a +#' subpopulation are \code{NA} at a given level, splitting stops for that +#' branch. +#' +#' The \code{level} argument is only relevant when \code{x} is a +#' \code{\link{MultiPop-class}} object. If \code{level = Inf}, all +#' \code{\link{Pop-class}} objects contained in \code{x} are split according +#' to \code{by}. If \code{level = c(a, b)}, only \code{Pop-class} objects at +#' nesting levels \code{a} and \code{b} (with the top level being \code{1}) are +#' split. An error is raised if any requested level exceeds the maximum nesting +#' depth of \code{x}. +#' #' @return Returns a \code{\link{MultiPop-class}} object. #' #' @examples #' # Create founder haplotypes #' founderPop = quickHaplo(nInd = 10, nChr = 1, segSites = 10) -#' +#' #' # Set simulation parameters #' SP = SimParam$new(founderPop) #' \dontshow{SP$nThreads = 1L} #' SP$addTraitA(10) -#' +#' #' # Create population #' pop = newPop(founderPop, simParam = SP) #' @@ -563,7 +583,23 @@ mergeMultiPops = function(..., level=0){ #' mp1 = splitPop(pop, by = rep(c("A", "B"), length.out = nInd(pop))) #' mp1 #' -#' # Nested split: First using a random grouping vector (three groups), then by family +#' # Split using a data frame, with rows aligned by individual ID +#' by_df = data.frame( +#' level1 = sample(LETTERS[1:3], nInd(pop), replace = TRUE), +#' level2 = sample(1:2, nInd(pop), replace = TRUE), +#' row.names = pop@id +#' ) +#' splitPop(pop, by = by_df) +#' +#' # Uneven nested structure using NA values in lower levels +#' by_df2 = data.frame( +#' level1 = c(rep("A", 4), rep("B", 6)), +#' level2 = c(rep(NA_character_, 4), rep(c("C", "D"), each = 3)), +#' row.names = pop@id +#' ) +#' splitPop(pop, by = by_df2) +#' +#' # Nested split: first by a random grouping vector, then by family #' mp2 = splitPop( #' pop, #' by = list( @@ -572,12 +608,14 @@ mergeMultiPops = function(..., level=0){ #' ) #' ) #' mp2 -#' +#' #' # When x is a MultiPop, control which nesting levels are split #' mp_nested = newMultiPop(pop[1:3], newMultiPop(pop[4:6], pop[7:9])) -#' splitPop(mp_nested, -#' by = function(p) rep(c("A","B"), length.out = nInd(p)), -#' level = 1) +#' splitPop( +#' mp_nested, +#' by = function(p) rep(c("A", "B"), length.out = nInd(p)), +#' level = 1 +#' ) #' #' @export splitPop = function(x, by, level = Inf) { @@ -592,6 +630,20 @@ splitPop = function(x, by, level = Inf) { if (length(by) == 0) { stop("`by` must have at least one grouping spec.") } + if (is.data.frame(by) && isPop(x)){ + if (.row_names_info(by) > 0) { + if (setequal(rownames(by), x@id)){ + by = by[x@id, , drop = FALSE] + } else { + warning("Row names of data frame `by` don't match `x@id`. Mapping rows by order instead of names.") + } + } else { + warning("Mapping rows of data frame (`by`) to individuals' identifiers (`x@id`) by order.\n", + "Consider setting row names of `by` to match `x@id` for clarity.") + } + res = .splitPop(x, by) + return(res) + } # Get max depth of nesting in MultiPop md = ifelse(isMultiPop(x), .depthMultiPop(x), 1L) @@ -669,8 +721,16 @@ splitPop = function(x, by, level = Inf) { #' Helper function to recursively split a Pop object #' #' @param pop A \code{\link{Pop-class}} object. -#' @param by A list of grouping specs (vectors or functions). Each element is -#' applied in order to produce nested splits. +#' @param by A grouping specification. This may be: \cr +#' - an atomic vector. \cr +#' - a function returning an atomic vector. \cr +#' - a list of vectors/functions for recursive nested splitting. \cr +#' - a data frame whose columns define nested grouping levels. +#' +#' For data frames, rows are aligned to \code{pop@id} when possible, otherwise +#' they are assumed to already be in population order. During recursion, the +#' data frame is subset to each child population so that row alignment is +#' preserved across levels. #' #' @return A \code{\link{MultiPop-class}} produced by recursively applying \code{by}. #' @@ -686,23 +746,20 @@ splitPop = function(x, by, level = Inf) { if (!is.atomic(groups)) { stop("Grouping spec must be an atomic vector or a function returning one.") } + if (all(is.na(groups))) { + return(pop) + } if (any(is.na(groups))) { stop("Grouping vector contains NA values.") } - n = nInd(pop) - if (length(groups) == 1L) { - groups = rep(groups, n) - } - if (length(groups) != n) { - stop( - "Grouping vector length (", length(groups), - ") must equal `nInd(pop)` (", n, "), or be length 1." - ) - } - popList = split(pop, groups) mp = newEmptyMultiPop() - mp@pops = lapply(popList, .splitPop, by = by[-1]) + if (is.data.frame(by)){ + groups = lapply(split(by, groups), `[`, -1) + mp@pops = mapply(FUN = .splitPop, pop = popList, by = groups) + } else { + mp@pops = lapply(popList, .splitPop, by = by[-1]) + } return(mp) } \ No newline at end of file diff --git a/man/dot-splitPop.Rd b/man/dot-splitPop.Rd index 17e5520f..f69a27f0 100644 --- a/man/dot-splitPop.Rd +++ b/man/dot-splitPop.Rd @@ -9,8 +9,16 @@ \arguments{ \item{pop}{A \code{\link{Pop-class}} object.} -\item{by}{A list of grouping specs (vectors or functions). Each element is -applied in order to produce nested splits.} +\item{by}{A grouping specification. This may be: \cr +- an atomic vector. \cr +- a function returning an atomic vector. \cr +- a list of vectors/functions for recursive nested splitting. \cr +- a data frame whose columns define nested grouping levels. + +For data frames, rows are aligned to \code{pop@id} when possible, otherwise +they are assumed to already be in population order. During recursion, the +data frame is subset to each child population so that row alignment is +preserved across levels.} } \value{ A \code{\link{MultiPop-class}} produced by recursively applying \code{by}. diff --git a/man/splitPop.Rd b/man/splitPop.Rd index 0fa9d21d..64643579 100644 --- a/man/splitPop.Rd +++ b/man/splitPop.Rd @@ -9,28 +9,48 @@ splitPop(x, by, level = Inf) \arguments{ \item{x}{a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object} -\item{by}{a vector, function, or list of vectors/functions defining -groupings. Functions are called with a \code{Pop} and must return an -atomic vector of length \code{nInd(pop)} (or 1).} +\item{by}{a vector, function, list, or data frame defining groupings} -\item{level}{A positive integer, a vector of positive integers, or +\item{level}{a positive integer, a vector of positive integers, or \code{Inf}. Only relevant when \code{x} is a \code{\link{MultiPop-class}} -object. If \code{level = Inf}, split any \code{\link{Pop-class}} within -\code{x} according to \code{by}. If \code{level = c(a, b)}, only -\code{Pop-class} objects at levels \code{a} and \code{b} (top \code{level=1}) -are split. An error is raised if any requested level is deeper than the -object's maximum nesting depth.} +object. See Details.} } \value{ Returns a \code{\link{MultiPop-class}} object. } \description{ Split a \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object -into a \code{MultiPop} at specified nesting \code{level}(s) using one -or more grouping specifications (\code{by}). Grouping specs can be -atomic vectors (length \code{nInd(x)} or 1) or functions that return -such vectors. If \code{by} is a list, each element is applied -recursively to create nested \code{MultiPop} objects. +into a \code{MultiPop} using one or more grouping specifications passed +through \code{by}. +} +\details{ +The \code{by} argument can take several forms: \cr +- Atomic vectors: used directly as grouping labels. The + vector is passed to \code{\link[base]{split}} and may be recycled as + needed. If the vector length is not a multiple of the population size, + \code{\link[base]{split}} issues a warning. \cr +- Functions: called on each \code{Pop} object and must return + an atomic vector of grouping labels. \cr +- Lists: each element must be a vector or function. Elements + are applied recursively to create nested \code{\link{MultiPop-class}} + objects. \cr +- Data frames: each column defines one nesting level. In this + case, \code{x} must be a \code{\link{Pop-class}} object. Rows of the data + frame must correspond to individuals in \code{x}. If row names are present + and match \code{x@id}, they are used to align rows to individuals; + otherwise rows are assumed to be in population order and a warning is + issued. \code{NA} values are allowed in grouping columns and can be used + to represent uneven nesting structures. If all grouping values for a + subpopulation are \code{NA} at a given level, splitting stops for that + branch. + +The \code{level} argument is only relevant when \code{x} is a +\code{\link{MultiPop-class}} object. If \code{level = Inf}, all +\code{\link{Pop-class}} objects contained in \code{x} are split according +to \code{by}. If \code{level = c(a, b)}, only \code{Pop-class} objects at +nesting levels \code{a} and \code{b} (with the top level being \code{1}) are +split. An error is raised if any requested level exceeds the maximum nesting +depth of \code{x}. } \examples{ # Create founder haplotypes @@ -48,7 +68,23 @@ pop = newPop(founderPop, simParam = SP) mp1 = splitPop(pop, by = rep(c("A", "B"), length.out = nInd(pop))) mp1 -# Nested split: First using a random grouping vector (three groups), then by family +# Split using a data frame, with rows aligned by individual ID +by_df = data.frame( + level1 = sample(LETTERS[1:3], nInd(pop), replace = TRUE), + level2 = sample(1:2, nInd(pop), replace = TRUE), + row.names = pop@id +) +splitPop(pop, by = by_df) + +# Uneven nested structure using NA values in lower levels +by_df2 = data.frame( + level1 = c(rep("A", 4), rep("B", 6)), + level2 = c(rep(NA_character_, 4), rep(c("C", "D"), each = 3)), + row.names = pop@id +) +splitPop(pop, by = by_df2) + +# Nested split: first by a random grouping vector, then by family mp2 = splitPop( pop, by = list( @@ -60,8 +96,10 @@ mp2 # When x is a MultiPop, control which nesting levels are split mp_nested = newMultiPop(pop[1:3], newMultiPop(pop[4:6], pop[7:9])) -splitPop(mp_nested, - by = function(p) rep(c("A","B"), length.out = nInd(p)), - level = 1) +splitPop( + mp_nested, + by = function(p) rep(c("A", "B"), length.out = nInd(p)), + level = 1 +) } From 36ac4e3742f1bfc56b2c5a73ad3904e9970bc834 Mon Sep 17 00:00:00 2001 From: Gregor Gorjanc Date: Thu, 9 Jul 2026 10:15:59 +0200 Subject: [PATCH 8/9] Incorporate suggestions from code review --- R/Class-Pop.R | 2 +- R/mergePops.R | 18 ++++++---- R/selection.R | 2 +- man/calcPopValue.Rd | 2 +- man/dot-flattenMultiPop.Rd | 2 +- man/flattenMultiPop.Rd | 9 ++--- tests/testthat/test-mergePops.R | 59 +++++++++++++++++++++++++++------ 7 files changed, 66 insertions(+), 28 deletions(-) diff --git a/R/Class-Pop.R b/R/Class-Pop.R index a27c60c4..77b5a18d 100644 --- a/R/Class-Pop.R +++ b/R/Class-Pop.R @@ -1077,7 +1077,7 @@ setMethod("show", indexLabel = paste0("[[", idx, "]] ") # Dont print name label if it's empty - if (nameLabel %in% c(" \"NA\" - ", " \"\" - ")) { + if (nameLabel %in% c(" \"NA\" - ", " \"\" - ", NA_character_)) { nameLabel = "" } diff --git a/R/mergePops.R b/R/mergePops.R index b8f24924..f95e282e 100644 --- a/R/mergePops.R +++ b/R/mergePops.R @@ -247,8 +247,7 @@ mergePops = function(popList){ #' @return If \code{x} is a \code{\link{Pop-class}}, the same \code{x} #' object is returned. Otherwise a \code{\link{MultiPop-class}} is returned #' whose \code{x@pops} slot contains \code{\link{Pop-class}} (and possibly -#' \code{\link{MultiPop-class}}) objects flattened according to \code{level} -#' and \code{preserveNames}. +#' \code{\link{MultiPop-class}}) objects flattened according to \code{level}. #' #' @seealso \code{\link{mergeMultiPops}} and \code{\link{mergePops}} #' @@ -287,9 +286,11 @@ mergePops = function(popList){ #' flattenMultiPop(mp_nested, level = 2) #' #' @export -flattenMultiPop = function(x, level = 1, - preserveNames = c("auto", "concatenate", "force", "none")) { - preserveNames = match.arg(preserveNames) +flattenMultiPop = function(x, level = 1, preserveNames = "auto") { + preserveNames = match.arg( + arg = preserveNames, + choices = c("auto", "concatenate", "force", "none") + ) if (isPop(x)) { return(x) } @@ -319,8 +320,11 @@ flattenMultiPop = function(x, level = 1, #' @param mp \code{\link{MultiPop-class}} object #' #' @keywords internal -.flattenMultiPop = function(mp, preserveNames = c("auto", "concatenate", "force", "none")) { - preserveNames = match.arg(preserveNames) +.flattenMultiPop = function(mp, preserveNames = "auto") { + preserveNames = match.arg( + arg = preserveNames, + choices = c("auto", "concatenate", "force", "none") + ) if (.depthMultiPop(mp) == 1L) { res = mp@pops diff --git a/R/selection.R b/R/selection.R index d163ed1f..9951ce6f 100644 --- a/R/selection.R +++ b/R/selection.R @@ -98,7 +98,7 @@ getResponse = function(pop,trait,use,simParam=NULL,nThreads=NULL,...){ #' @param simplify Logical. If \code{TRUE}, simplify the output by flattening #' the \code{MultiPop} in \code{x} to the requested \code{level} using #' \code{\link{flattenMultiPop}}. The output matrices from each \code{Pop} -#' are combined with \code{\link{rbind}}, and a \code{"source"} attribute +#' are combined with \code{\link{rbind}}, and a \code{"source"} attribute #' is added to indicate the origin of each row. #' @param level Integer scalar >= 1. Number of \code{MultiPop} levels to #' preserve when \code{simplify=TRUE}. Passed to \code{\link{flattenMultiPop}}. diff --git a/man/calcPopValue.Rd b/man/calcPopValue.Rd index 81ef5eb0..ba355188 100644 --- a/man/calcPopValue.Rd +++ b/man/calcPopValue.Rd @@ -15,7 +15,7 @@ as its first argument.} \item{simplify}{Logical. If \code{TRUE}, simplify the output by flattening the \code{MultiPop} in \code{x} to the requested \code{level} using \code{\link{flattenMultiPop}}. The output matrices from each \code{Pop} - are combined with \code{\link{rbind}}, and a \code{"source"} attribute +are combined with \code{\link{rbind}}, and a \code{"source"} attribute is added to indicate the origin of each row.} \item{level}{Integer scalar >= 1. Number of \code{MultiPop} levels to diff --git a/man/dot-flattenMultiPop.Rd b/man/dot-flattenMultiPop.Rd index 13bfaf55..5b61f463 100644 --- a/man/dot-flattenMultiPop.Rd +++ b/man/dot-flattenMultiPop.Rd @@ -4,7 +4,7 @@ \alias{.flattenMultiPop} \title{Helper function to recursively extract Pop objects from a MultiPop} \usage{ -.flattenMultiPop(mp, preserveNames = c("auto", "concatenate", "force", "none")) +.flattenMultiPop(mp, preserveNames = "auto") } \arguments{ \item{mp}{\code{\link{MultiPop-class}} object} diff --git a/man/flattenMultiPop.Rd b/man/flattenMultiPop.Rd index a0a213b2..a365eb2c 100644 --- a/man/flattenMultiPop.Rd +++ b/man/flattenMultiPop.Rd @@ -4,11 +4,7 @@ \alias{flattenMultiPop} \title{Flatten a MultiPop object to a specified depth} \usage{ -flattenMultiPop( - x, - level = 1, - preserveNames = c("auto", "concatenate", "force", "none") -) +flattenMultiPop(x, level = 1, preserveNames = "auto") } \arguments{ \item{x}{A \code{\link{Pop-class}} or \code{\link{MultiPop-class}} object.} @@ -30,8 +26,7 @@ flattening. One of: \cr If \code{x} is a \code{\link{Pop-class}}, the same \code{x} object is returned. Otherwise a \code{\link{MultiPop-class}} is returned whose \code{x@pops} slot contains \code{\link{Pop-class}} (and possibly -\code{\link{MultiPop-class}}) objects flattened according to \code{level} -and \code{preserveNames}. +\code{\link{MultiPop-class}}) objects flattened according to \code{level}. } \description{ Recursively flatten a \code{\link{MultiPop-class}} object into a shallower diff --git a/tests/testthat/test-mergePops.R b/tests/testthat/test-mergePops.R index f89f70b8..2a3df306 100644 --- a/tests/testthat/test-mergePops.R +++ b/tests/testthat/test-mergePops.R @@ -273,24 +273,22 @@ test_that("splitPop", { fixed = TRUE ) - expect_error( - splitPop(pop, by = numeric(0)), - paste0("Grouping vector length (", length(numeric(0)), - ") must equal `nInd(pop)` (", length(pop), "), or be length 1."), + expect_warning( + splitPop(pop, by = numeric(7)), + "data length is not a multiple of split variable", fixed = TRUE ) expect_error( - splitPop(pop, by = list(NA_real_)), + splitPop(pop, by = list(c(NA_real_, 0))), "Grouping vector contains NA values.", fixed = TRUE ) - expect_error( - splitPop(pop, by = list(1:4)), - paste0("Grouping vector length (", length(list(1:4)[[1]]), - ") must equal `nInd(pop)` (", length(pop), "), or be length 1."), - fixed = TRUE + # split.default() recycling behavior + expect_identical( + splitPop(pop, by = 1:4), + splitPop(pop, by = rep(1:4, length.out = nInd(pop))) ) expect_error( @@ -366,4 +364,45 @@ test_that("splitPop", { "requested level(s) exceed max depth of `x` (2)", fixed = TRUE ) + + # splitPop using a data frame for grouping + by1 = data.frame(level1 = by1) + expect_warning( + expect_identical(mp1, splitPop(pop, by = by1)), + "Mapping rows of data frame (`by`) to individuals' identifiers (`x@id`) by order.\nConsider setting row names of `by` to match `x@id` for clarity.", + fixed = TRUE + ) + + rownames(by1) = letters[pop@iid] + expect_warning( + splitPop(pop, by = data.frame(level1 = by1)), + "Row names of data frame `by` don't match `x@id`. Mapping rows by order instead of names.", + fixed = TRUE + ) + + # Specify row names of the data frame to match `x@id` + by_df = data.frame(level1 = by1, row.names = pop@id) + expect_identical(mp1, splitPop(pop, by = by_df)) + + # splitPop can create nested MultiPop objects using a data frame with multiple columns + nms = sapply(mp2@pops, lengths) + by_df = data.frame(level1 = rep(names(nms), nms), + level2 = "0_0", + row.names = mergePops(mp2)@id) + expect_identical(mp2, splitPop(pop, by = by_df)) + + # A data frame with NA values in the grouping columns creates an uneven nested structure + by_df = data.frame(level1 = c(rep("A", 4), rep("B", 8)), + level2 = c(rep(NA_character_, 4), rep(LETTERS[3:4], each = 4)), + level3 = c(rep(NA_character_, 8), rep("E", 4)), + row.names = pop@id) + + mp4 = newMultiPop(A = pop[1:4], + B = newMultiPop( + C = pop[5:8], + D = newMultiPop(E = pop[9:12]))) + + expect_identical(mp4, splitPop(pop, by = by_df)) + + }) From 11d6d9d8152a6e210c76dd1b9e0b927b63042611 Mon Sep 17 00:00:00 2001 From: Daniel Ariza Suarez Date: Thu, 9 Jul 2026 11:04:50 +0200 Subject: [PATCH 9/9] Add prototypes for `genParamPop` and `calcGenParamPop` --- NAMESPACE | 1 + R/popSummary.R | 380 +++++++++++++++++++++++++++++++++ man/calcGenParamPop_R.Rd | 41 ++++ man/genParamPop.Rd | 44 ++++ tests/testthat/test-genParam.R | 146 +++++++++++++ 5 files changed, 612 insertions(+) create mode 100644 man/calcGenParamPop_R.Rd create mode 100644 man/genParamPop.Rd diff --git a/NAMESPACE b/NAMESPACE index 0a287a00..16b9817a 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -28,6 +28,7 @@ export(editGenomeTopQtl) export(fastRRBLUP) export(flattenMultiPop) export(genParam) +export(genParamPop) export(genicVarA) export(genicVarAA) export(genicVarD) diff --git a/R/popSummary.R b/R/popSummary.R index bee70277..a2eeebb5 100644 --- a/R/popSummary.R +++ b/R/popSummary.R @@ -902,3 +902,383 @@ mendelianSampling = function(pop, parents = NULL, mothers = NULL, fathers = NULL nInd = function(pop){ pop@nInd } + +#' @title Calculate trait parameters for populations +#' +#' @description +#' Internal helper used by \code{\link{genParamPop}} to compute +#' population summaries for a single trait. Written in R for prototyping +#' and discussion. +#' +#' @param multiPop an object of \code{\link{MultiPop-class}}. +#' @param trait a trait object from \code{simParam$traits}. +#' @param FUN function used to summarize values within each population +#' (e.g. \code{mean}, \code{sum}). +#' @param simParam an object of class \code{\link{SimParam}}. If +#' \code{NULL}, the function uses the object named \code{SP} from the +#' global environment. +#' @param ... additional arguments passed to \code{FUN}. +#' +#' @return +#' A list with per-population summaries for one trait, including: +#' \itemize{ +#' \item \code{gv}, \code{bv}, and (when available) \code{dd} +#' \item \code{gv_a}, and (when available) \code{gv_d} +#' \item \code{mu}, \code{gv_mu}, \code{alpha} +#' \item currently implemented variance components (e.g. \code{genicVarA}) +#' } +#' +#' @details +#' Prototype implementation under development. Intended for exploratory +#' population-level summaries and API discussion. +#' +#' @keywords internal +calcGenParamPop_R = function(multiPop, trait, FUN = mean, simParam=NULL, ...) { + + if(is.null(simParam)){ + simParam = get("SP",envir=.GlobalEnv) + } + + # Extract ploidy from population information + ploidy = multiPop@pops[[1]]@ploidy + + # Extract trait information + lociLoc = trait@lociLoc + lociPerChr = trait@lociPerChr + a = trait@addEff # Additive effects + nLoci = length(a) + intercept = trait@intercept + nPop = length(multiPop) + nIndPop = lengths(multiPop) + interceptPop = rep(0, nPop) + + + # Genotype dosage values (0, 1, 2 for diploids) + x = 0:ploidy + + # Additive coding: scales genotype to [-1, 1] for diploids + # For diploid: xa = (x - 1) * 1 = c(-1, 0, 1) + xa = (x - ploidy/2) * (2/ploidy) + + # Allele substitution effects + alpha = numeric(nLoci) + + # Initialize output vectors + gv_a = bvMat = vector(mode = 'double', length = nPop) + + # Additive genic variance (sum over loci, assuming no LD) + genicVarA = 0 # Additive genic variance (observed frequencies) + + # Check if trait has dominance effects + hasD = methods::.hasSlot(trait, "domEff") + if (hasD) { + d = trait@domEff # Dominance effects + + # Initialize output vectors + gv_d = ddMat = gv_a + + # Dominance coding: measures deviation from additivity + # For diploid: xd = x*(2-x)*(1/2)^2 = c(0, 0.5, 0) -> heterozygote = 0.5 + xd = x * (ploidy - x) * (2/ploidy)^2 + + # Dominance genic variance (sum over loci, assuming no LD) + genicVarD = 0 # Dominance genic variance (observed frequencies) + } + + + # Loop over each locus + for (loc in seq_len(nLoci)) { + + # Additive genetic value: aEff = xa * a[i] + aEff = xa * a[loc] + + # Initialize pop-specific vectors + genoMu = gv_a_pop = gv_d_pop = double() + + if (hasD) { + # Dominance genetic value: dEff = xd * d[i] + dEff = xd * d[loc] + # Total genetic value at this locus + gv = aEff + dEff + # Initialize pop-specific vectors + gv_d_pop = gv_a_pop + } else { + gv = aEff + } + + # Loop over each pop + for (p in seq_along(multiPop)) { + pop = multiPop@pops[[p]] + + # Get genotype matrix (individuals x QTL) + # Each cell contains the dosage (0 to ploidy) + lociMap = selectLoci(chr=NULL, lociPerChr, lociLoc) + genoMat = getGeno(pop@geno, lociMap$lociPerChr, lociMap$lociLoc, + simParam$nThreads) + genoMat = convToImat(genoMat) + rownames(genoMat) = pop@id + colnames(genoMat) = getLociNames(lociMap$lociPerChr, lociMap$lociLoc, + simParam$genMap) + + # Mean allele dosage for each population + genoMu = c(genoMu, mean(genoMat[,loc])) + + # Vectorize operations over individuals + geno = genoMat[,loc] + 1 # +1 for R indexing (genotype 0 -> index 1) + + # Append average additive genetic values per population + # gv_a_pop = c(gv_a_pop, mean(aEff[geno])) + gv_a_pop = c(gv_a_pop, FUN(aEff[geno], ...)) + + if (hasD) { + # Append average dominance genetic values per population + # gv_d_pop = c(gv_d_pop, mean(dEff[geno])) + gv_d_pop = c(gv_d_pop, FUN(dEff[geno], ...)) + } + } + + # Cumulative additive genetic values across loci + gv_a = gv_a + gv_a_pop + + if (hasD) { + # Cumulative dominance genetic values across loci + gv_d = gv_d + gv_d_pop + # Total genetic values per population + gv_t_pop = gv_a_pop + gv_d_pop + + } else { + # Total genetic values per population + gv_t_pop = gv_a_pop + } + + + # Allele substitution effect for this locus + alpha[loc] = cov(genoMu, gv_t_pop) / var(genoMu) + # Breeding values per population at this locus + bv = (genoMu - mean(genoMu)) * alpha[loc] + # Cumulative breeding values across loci + bvMat = bvMat + bv + + # Var(BV) assuming no LD between loci + # Using observed frequencies + # TODO: popVar(bv) -> genicVarA + genicVarA = genicVarA + popVar(matrix(bv))[1,1] + + if (hasD) { + # Dominance values per population at this locus + dd = gv_t_pop - mean(gv_t_pop) - bv + # Cumulative dominance deviations across loci + ddMat = ddMat + dd + # Var(BV) assuming no LD between loci + # Using observed frequencies + # TODO: popVar(dd) -> genicVarD + genicVarD = genicVarD + popVar(matrix(dd))[1,1] + } + } + + # Calculate population-specific intercepts + for (p in seq_along(multiPop)) { + interceptPop[p] = FUN(rep(intercept, nIndPop[p]), ...) + } + + if (hasD) { + # Cumulative total genetic values across loci + gv_t = gv_a + gv_d + result = list( + gv = gv_t + interceptPop, # Total genetic value + bv = bvMat, # Breeding values + dd = ddMat, # Dominance deviations + genicVarA = genicVarA, # Additive genic variance + varA = popVar(matrix(bvMat))[1,1], # Additive genetic variance + genicVarD = genicVarD, # Dominance genic variance + varD = popVar(matrix(ddMat))[1,1], # Dominance genetic variance + varG = popVar(matrix(gv_t))[1,1], # Total genetic variance + # genicVarA2 = genicVarA2, # Additive genic variance (HWE) + # genicVarD2 = genicVarD2, # Dominance genic variance (HWE) + mu = mean(gv_t + interceptPop) , # Population mean + # mu_HWE = mu_HWE + intercept, # Population mean under HWE + gv_a = gv_a, # Additive genetic values + gv_d = gv_d, # Dominance genetic values + gv_mu = interceptPop, # Intercept + alpha = alpha # Allele substitution effects + # alpha_HW = alpha_HW # Allele substitution effects (HWE) + ) + } else { + result = list( + gv = gv_a + interceptPop, + bv = bvMat, + genicVarA = genicVarA, + # genicVarA2 = genicVarA2, + mu = mean(gv_a) + interceptPop, + # mu_HWE = mu_HWE + intercept, + gv_a = gv_a, + gv_mu = interceptPop, + alpha = alpha + # alpha_HW = alpha_HW + ) + } + + return(result) +} + +#' @title Summarize genetic parameters across populations +#' +#' @description +#' Computes trait-wise population summaries of genetic values and related +#' variance components for a \code{\link{MultiPop-class}} object. +#' This is a prototype implementation intended for development and feedback. +#' +#' @param multiPop an object of \code{\link{MultiPop-class}}. +#' @param FUN function used to summarize values within each population +#' (e.g. \code{mean}, \code{sum}). +#' @param simParam an object of class \code{\link{SimParam}}. If +#' \code{NULL}, the function uses the object named \code{SP} from the +#' global environment. +#' @param nThreads number of threads to use if OpenMP is available. +#' If \code{NULL}, the number is obtained from \code{simParam$nThreads}. +#' @param ... additional arguments passed to \code{FUN}. +#' +#' @return +#' A list containing trait-level and population-level summaries, including: +#' \itemize{ +#' \item \code{varA}, \code{varD}, \code{varAA}, \code{varG} +#' \item \code{mu} +#' \item \code{gv}, \code{bv}, \code{dd}, \code{aa} +#' \item \code{gv_mu}, \code{gv_a}, \code{gv_d}, \code{gv_aa} +#' \item \code{alpha} +#' \item selected covariance and genic-variance components currently available +#' } +#' +#' @details +#' This function currently wraps \code{calcGenParamPop_R} trait-by-trait. +#' Behavior for non-linear summary functions depends on \code{FUN} and should +#' be interpreted with care. +#' +#' @export +genParamPop = function(multiPop, FUN = mean, simParam = NULL, nThreads = NULL, ...){ + if(is.null(simParam)){ + simParam = get("SP",envir=.GlobalEnv) + } + if(is.null(nThreads)){ + nThreads = simParam$nThreads + }else{ + nThreads = as.integer(nThreads) + } + + nPop = length(multiPop) + nTraits = simParam$nTraits + traitNames = simParam$traitNames + + # Blank nPop multiPop nTrait matrices + gv = matrix(NA_real_, nrow=nPop, ncol=nTraits) + colnames(gv) = traitNames + bv = dd = aa = gv_a = gv_d = gv_aa = gv_mu = gv + + # Blank nTrait vectors + genicVarA = rep(NA_real_, nTraits) + names(genicVarA) = traitNames + genicVarD = genicVarAA = covA_HW = covD_HW = covAA_HW = + covG_HW = mu = mu_HW = covAAA_L = covDAA_L = + covAD_L = genicVarA + + # Average effect of an allele substitution + alpha = vector("list", length=nTraits) + names(alpha) = traitNames + # alpha_HW = alpha + + #Loop through trait calculations + for(i in seq_len(nTraits)){ + trait = simParam$traits[[i]] + tmp = calcGenParamPop_R(multiPop = multiPop, trait = trait, FUN = FUN, simParam = simParam, ...) + # genicVarA[i] = tmp$genicVarA2 + # covA_HW[i] = tmp$genicVarA-tmp$genicVarA2 + gv[,i] = tmp$gv + bv[,i] = tmp$bv + mu[i] = tmp$mu + # mu_HW[i] = tmp$mu_HWE + gv_a[,i] = tmp$gv_a + gv_mu[,i] = tmp$gv_mu + if(.hasSlot(trait,"domEff")){ + # genicVarD[i] = tmp$genicVarD2 + # covD_HW[i] = tmp$genicVarD-tmp$genicVarD2 + dd[,i] = tmp$dd + gv_d[,i] = tmp$gv_d + }else{ + # genicVarD[i] = 0 + # covD_HW[i] = 0 + dd[,i] = rep(0,nPop) + gv_d[,i] = rep(0,nPop) + } + if(.hasSlot(trait,"epiEff")){ + genicVarAA[i] = tmp$genicVarAA2 + covAA_HW[i] = tmp$genicVarAA-tmp$genicVarAA2 + aa[,i] = tmp$aa + gv_aa[,i] = tmp$gv_aa + }else{ + genicVarAA[i] = 0 + covAA_HW[i] = 0 + aa[,i] = rep(0,nPop) + gv_aa[,i] = rep(0,nPop) + } + if(nPop==1){ + covAD_L[i] = 0 + covAAA_L[i] = 0 + covDAA_L[i] = 0 + } else { + covAD_L[i] = popVar(cbind(bv[,i],dd[,i]))[1,2] + covAAA_L[i] = popVar(cbind(bv[,i],aa[,i]))[1,2] + covDAA_L[i] = popVar(cbind(dd[,i],aa[,i]))[1,2] + } + alpha[[i]] = tmp$alpha + # alpha_HW[[i]] = tmp$alpha_HW + } + + varA = popVar(bv) + rownames(varA) = colnames(varA) = traitNames + + varD = popVar(dd) + rownames(varD) = colnames(varD) = traitNames + + varAA = popVar(aa) + rownames(varAA) = colnames(varAA) = traitNames + + varG = popVar(gv) + rownames(varG) = colnames(varG) = traitNames + + # genicVarG = genicVarA + genicVarD + genicVarAA + # covG_HW = covA_HW + covD_HW + covAA_HW + + output = list(varA=varA, + varD=varD, + varAA=varAA, + varG=varG, + # genicVarA=genicVarA, + # genicVarD=genicVarD, + genicVarAA=genicVarAA, + # genicVarG=genicVarG, + # # covA_HW=covA_HW, + # covD_HW=covD_HW, + covAA_HW=covAA_HW, + covG_HW=covG_HW, + # covA_L=diag(varA)-genicVarA-covA_HW, + # covD_L=diag(varD)-genicVarD-covD_HW, + covAA_L=diag(varAA)-genicVarAA-covAA_HW, + covAD_L=covAD_L, + covAAA_L=covAAA_L, + covDAA_L=covDAA_L, + # covG_L=diag(varG)-genicVarG-covG_HW, + mu=mu, + # mu_HW=mu_HW, + gv=gv, + bv=bv, + dd=dd, + aa=aa, + gv_mu=gv_mu, + gv_a=gv_a, + gv_d=gv_d, + gv_aa=gv_aa, + alpha=alpha#, + # alpha_HW=alpha_HW + ) + return(output) +} diff --git a/man/calcGenParamPop_R.Rd b/man/calcGenParamPop_R.Rd new file mode 100644 index 00000000..d467a051 --- /dev/null +++ b/man/calcGenParamPop_R.Rd @@ -0,0 +1,41 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/popSummary.R +\name{calcGenParamPop_R} +\alias{calcGenParamPop_R} +\title{Calculate trait parameters for populations} +\usage{ +calcGenParamPop_R(multiPop, trait, FUN = mean, simParam = NULL, ...) +} +\arguments{ +\item{multiPop}{an object of \code{\link{MultiPop-class}}.} + +\item{trait}{a trait object from \code{simParam$traits}.} + +\item{FUN}{function used to summarize values within each population +(e.g. \code{mean}, \code{sum}).} + +\item{simParam}{an object of class \code{\link{SimParam}}. If +\code{NULL}, the function uses the object named \code{SP} from the +global environment.} + +\item{...}{additional arguments passed to \code{FUN}.} +} +\value{ +A list with per-population summaries for one trait, including: +\itemize{ + \item \code{gv}, \code{bv}, and (when available) \code{dd} + \item \code{gv_a}, and (when available) \code{gv_d} + \item \code{mu}, \code{gv_mu}, \code{alpha} + \item currently implemented variance components (e.g. \code{genicVarA}) +} +} +\description{ +Internal helper used by \code{\link{genParamPop}} to compute +population summaries for a single trait. Written in R for prototyping +and discussion. +} +\details{ +Prototype implementation under development. Intended for exploratory +population-level summaries and API discussion. +} +\keyword{internal} diff --git a/man/genParamPop.Rd b/man/genParamPop.Rd new file mode 100644 index 00000000..cdf3e825 --- /dev/null +++ b/man/genParamPop.Rd @@ -0,0 +1,44 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/popSummary.R +\name{genParamPop} +\alias{genParamPop} +\title{Summarize genetic parameters across populations} +\usage{ +genParamPop(multiPop, FUN = mean, simParam = NULL, nThreads = NULL, ...) +} +\arguments{ +\item{multiPop}{an object of \code{\link{MultiPop-class}}.} + +\item{FUN}{function used to summarize values within each population +(e.g. \code{mean}, \code{sum}).} + +\item{simParam}{an object of class \code{\link{SimParam}}. If +\code{NULL}, the function uses the object named \code{SP} from the +global environment.} + +\item{nThreads}{number of threads to use if OpenMP is available. +If \code{NULL}, the number is obtained from \code{simParam$nThreads}.} + +\item{...}{additional arguments passed to \code{FUN}.} +} +\value{ +A list containing trait-level and population-level summaries, including: +\itemize{ + \item \code{varA}, \code{varD}, \code{varAA}, \code{varG} + \item \code{mu} + \item \code{gv}, \code{bv}, \code{dd}, \code{aa} + \item \code{gv_mu}, \code{gv_a}, \code{gv_d}, \code{gv_aa} + \item \code{alpha} + \item selected covariance and genic-variance components currently available +} +} +\description{ +Computes trait-wise population summaries of genetic values and related +variance components for a \code{\link{MultiPop-class}} object. +This is a prototype implementation intended for development and feedback. +} +\details{ +This function currently wraps \code{calcGenParamPop_R} trait-by-trait. +Behavior for non-linear summary functions depends on \code{FUN} and should +be interpreted with care. +} diff --git a/tests/testthat/test-genParam.R b/tests/testthat/test-genParam.R index 7cc367ef..a7ec32a4 100644 --- a/tests/testthat/test-genParam.R +++ b/tests/testthat/test-genParam.R @@ -1044,3 +1044,149 @@ test_that("genParam", { # Multiple locus traits # TODO }) + +test_that("genParamPop", { + + # Recycling example from `genParam` test above, but using `genParamPop` instead + # Haplotypes (a bit more than usual so get desired allele and genotype freqs) + haplo = matrix(data = 0, nrow = 200, ncol = 3) + + # Locus 1 - in Hardy-Weinberg equilibrium (=random mating, no inbreeding) + # Aiming for q = 0.1 as in Falconer (1996) example 7.1 + haplo[1:2, 1] = 0 # 1 individual with 0/0 + haplo[3:38, 1] = rep(c(0, 1), times = 18) # 18 individuals with 0/1 + haplo[39:200, 1] = 1 # 81 individuals with 1/1 + + # Locus 2 - in Hardy-Weinberg equilibrium (=random mating, no inbreeding) + # Aiming for q = 0.4 as in Falconer (1996) example 7.1 + haplo[1:32, 2] = 0 # 16 individuals with 0/0 + haplo[33:128, 2] = rep(c(0, 1), times = 48) # 48 individuals with 0/1 + haplo[129:200, 2] = 1 # 36 individuals with 1/1 + + # Locus 3 - not in Hardy-Weinberg equilibrium (=non-random mating, inbreeding) + # Aiming for q = 0.4 as in Falconer (1996) example 7.1 + haplo[1:68, 3] = 0 # 34 individuals with 0/0 + haplo[69:92, 3] = rep(c(0, 1), times = 12) # 12 individuals with 0/1 + haplo[93:200, 3] = 1 # 54 individuals with 1/1 + + nInd = nrow(haplo) / 2 + nLoc = ncol(haplo) + colnames(haplo) = letters[1:nLoc] + + genMap = data.frame( + markerName = letters[1:nLoc], + chromosome = c(1, 1, 1), + position = c(0, 0.2, 0.4) + ) + + ped = data.frame( + id = as.character(1:nInd), + mother = rep(0, nInd), + father = rep(0, nInd) + ) + + founderPop = importHaplo( + haplo = haplo, + genMap = genMap, + ploidy = 2L, + ped = ped + ) + + SP = SimParam$new(founderPop = founderPop) + SP$nThreads = 1L + + a = c(4, 4, 4) + d = c(2, 2, 2) + SP$importTrait( + markerNames = "a", + addEff = a[1], + domEff = d[1], + # markerNames = letters[1:nLoc], + # addEff = a, + # domEff = d, + intercept = 10, + name = "Falconer7.1_q=0.1_HWE" + ) + SP$importTrait( + markerNames = "b", + addEff = a[2], + domEff = d[2], + intercept = 10, + name = "Falconer7.1_q=0.4_HWE" + ) + SP$importTrait( + markerNames = "c", + addEff = a[3], + domEff = d[3], + intercept = 10, + name = "Falconer7.1_q=0.4_not_HWE" + ) + SP$importTrait( + markerNames = letters[1:3], + addEff = a, + domEff = d, + intercept = 10, + name = "Multi-locus" + ) + + # Create a new shuffled population + pop = newPop(founderPop, simParam = SP)[sample(1:100, 100, replace = F)] + multiPop = split(pop, 1:100) + multiPop = do.call(newMultiPop, unname(multiPop)) + + ans1_C = calcGenParam(trait = SP$traits[[1]], pop, SP$nThreads) + ans1_R = calcGenParamPop_R(trait = SP$traits[[1]], multiPop = multiPop, simParam = SP) + + expect_equal(ans1_R$gv, ans1_C$gv[,1]) + expect_equal(ans1_R$bv, ans1_C$bv[,1]) + expect_equal(ans1_R$dd, ans1_C$dd[,1]) + expect_equal(ans1_R$genicVarA, ans1_C$genicVarA) + expect_equal(ans1_R$genicVarD, ans1_C$genicVarD) + expect_equal(ans1_R$mu, ans1_C$mu) + expect_equal(ans1_R$gv_a, ans1_C$gv_a[,1]) + expect_equal(ans1_R$gv_d, ans1_C$gv_d[,1]) + expect_equal(mean(ans1_R$gv_mu), ans1_C$gv_mu) + expect_equal(ans1_R$alpha, ans1_C$alpha[,1]) + + ans2_C = calcGenParam(trait = SP$traits[[2]], pop, SP$nThreads) + ans2_R = calcGenParamPop_R(trait = SP$traits[[2]], multiPop = multiPop, simParam = SP) + + expect_equal(ans2_R$gv, ans2_C$gv[,1]) + expect_equal(ans2_R$bv, ans2_C$bv[,1]) + expect_equal(ans2_R$dd, ans2_C$dd[,1]) + expect_equal(ans2_R$genicVarA, ans2_C$genicVarA) + expect_equal(ans2_R$genicVarD, ans2_C$genicVarD) + expect_equal(ans2_R$mu, ans2_C$mu) + expect_equal(ans2_R$gv_a, ans2_C$gv_a[,1]) + expect_equal(ans2_R$gv_d, ans2_C$gv_d[,1]) + expect_equal(mean(ans2_R$gv_mu), ans2_C$gv_mu) + expect_equal(ans2_R$alpha, ans2_C$alpha[,1]) + + ans3_C = calcGenParam(trait = SP$traits[[3]], pop, SP$nThreads) + ans3_R = calcGenParamPop_R(trait = SP$traits[[3]], multiPop = multiPop, simParam = SP) + + expect_equal(ans3_R$gv, ans3_C$gv[,1]) + expect_equal(ans3_R$bv, ans3_C$bv[,1]) + expect_equal(ans3_R$dd, ans3_C$dd[,1]) + expect_equal(ans3_R$genicVarA, ans3_C$genicVarA) + expect_equal(ans3_R$genicVarD, ans3_C$genicVarD) + expect_equal(ans3_R$mu, ans3_C$mu) + expect_equal(ans3_R$gv_a, ans3_C$gv_a[,1]) + expect_equal(ans3_R$gv_d, ans3_C$gv_d[,1]) + expect_equal(mean(ans3_R$gv_mu), ans3_C$gv_mu) + expect_equal(ans3_R$alpha, ans3_C$alpha[,1]) + + ans4_C = calcGenParam(trait = SP$traits[[4]], pop, SP$nThreads) + ans4_R = calcGenParamPop_R(trait = SP$traits[[4]], multiPop = multiPop, simParam = SP) + + expect_equal(ans4_R$gv, ans4_C$gv[,1]) + expect_equal(ans4_R$bv, ans4_C$bv[,1]) + expect_equal(ans4_R$dd, ans4_C$dd[,1]) + expect_equal(ans4_R$genicVarA, ans4_C$genicVarA) + expect_equal(ans4_R$genicVarD, ans4_C$genicVarD) + expect_equal(ans4_R$mu, ans4_C$mu) + expect_equal(ans4_R$gv_a, ans4_C$gv_a[,1]) + expect_equal(ans4_R$gv_d, ans4_C$gv_d[,1]) + expect_equal(mean(ans4_R$gv_mu), ans4_C$gv_mu) + expect_equal(ans4_R$alpha, ans4_C$alpha[,1]) +})