diff --git a/.gitignore b/.gitignore
index c8dd2d6a..0428c571 100644
--- a/.gitignore
+++ b/.gitignore
@@ -13,6 +13,8 @@
# Output files from R CMD build
/*.tar.gz
+src/*.o
+src/*.so
# Output files from R CMD check
/*.Rcheck/
diff --git a/DESCRIPTION b/DESCRIPTION
index 5ac2bc7d..4ef4a624 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -26,7 +26,7 @@ URL: https://github.com/HighlanderLab/SIMplyBee
License: MIT + file LICENSE
Encoding: UTF-8
LazyData: true
-Imports: methods, R6, stats, utils, extraDistr (>= 1.9.1), RANN, Rcpp (>= 0.12.7)
+Imports: methods, R6, stats, utils, extraDistr (>= 1.9.1), RANN, Rcpp (>= 0.12.7), foreach, doParallel
Depends: R (>= 3.3.0), AlphaSimR (>= 1.5.3)
LinkingTo: Rcpp, RcppArmadillo (>= 0.7.500.0.0), BH
RoxygenNote: 7.3.2
diff --git a/NAMESPACE b/NAMESPACE
index 231169a0..60d468aa 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -2,6 +2,7 @@
export(SimParamBee)
export(addCastePop)
+export(addCastePop_internal)
export(addDrones)
export(addVirginQueens)
export(addWorkers)
@@ -75,7 +76,6 @@ export(getPooledGeno)
export(getQtlGeno)
export(getQtlHaplo)
export(getQueen)
-export(getQueenAge)
export(getQueenCsdAlleles)
export(getQueenCsdGeno)
export(getQueenGv)
@@ -206,7 +206,10 @@ exportClasses(MultiColony)
import(AlphaSimR)
import(Rcpp)
importFrom(R6,R6Class)
+importFrom(doParallel,registerDoParallel)
importFrom(extraDistr,rtpois)
+importFrom(foreach,"%dopar%")
+importFrom(foreach,foreach)
importFrom(methods,"slot<-")
importFrom(methods,classLabel)
importFrom(methods,is)
diff --git a/NEWS.md b/NEWS.md
index d72a41ac..c0766ff7 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -52,6 +52,20 @@ which caused an error. We now read in the locations from a csv file.
- Added new C++ function isHeterozygous() to speed up the SIMplyBee function isCsdHeterozygous()
+- parallelised all the major functions (so they run on simParamBee$nThreads cores)
+
+- swarm/split/supersede do no longer store the year of the queen
+
+- colonies with high inbreeding that do not produce a viable virgin queens in
+max(10, SP$nVirginQueens) attempts are
+removed in swarm/supersede
+
+- split no longer creates virgin queens in the split colonies but returns colonies with workers
+and meta data, but no virgin
+queens
+
+- createMultiColony() no longer creates an empty apiary, but it adds empty colonies with IDs
+
## Bug fixes
- Bug fix - get\*Haplo() functions were returning diploid drones when
diff --git a/R/Class-SimParamBee.R b/R/Class-SimParamBee.R
index 5ca066c3..2b66a241 100644
--- a/R/Class-SimParamBee.R
+++ b/R/Class-SimParamBee.R
@@ -425,13 +425,112 @@ SimParamBee <- R6Class(
invisible(self)
},
+ #' @description For internal use only.
+ #'
+ #' @param nNewInd Number of newly created individuals
+ #' @param id the name of each individual
+ #' @param mother vector of mother iids
+ #' @param father vector of father iids
+ #' @param isDH indicator for DH lines
+ addToBeePed = function(nNewInd,id,mother,father,isDH) {
+ stopifnot(nNewInd>0)
+ if(length(isDH)==1) isDH = rep(isDH,nNewInd)
+ mother = as.integer(mother)
+ father = as.integer(father)
+ isDH = as.integer(isDH)
+ stopifnot(length(mother)==nNewInd,
+ length(father)==nNewInd,
+ length(isDH)==nNewInd)
+ tmp = cbind(mother,father,isDH)
+ rownames(tmp) = id
+ private$.pedigree = rbind(private$.pedigree,tmp)
+ private$.lastId = private$.lastId + as.integer(nNewInd)
+ invisible(self)
+ },
+
+
+ #' @description For internal use only.
+ #'
+ #' @param nNewInd Number of newly created individuals
+ #' @param id the name of each individual
+ #' @param mother vector of mother iids
+ #' @param father vector of father iids
+ #' @param isDH indicator for DH lines
+ #' @param hist new recombination history
+ #' @param ploidy ploidy level
+ addToBeeRec = function(nNewInd,id,mother,father,isDH,
+ hist,ploidy){
+ stopifnot(nNewInd>0)
+ if(length(isDH)==1) isDH = rep(isDH,nNewInd)
+ mother = as.integer(mother)
+ father = as.integer(father)
+ isDH = as.integer(isDH)
+ stopifnot(length(mother)==nNewInd,
+ length(father)==nNewInd,
+ length(isDH)==nNewInd)
+ tmp = cbind(mother,father,isDH)
+ rownames(tmp) = id
+ if(is.null(hist)){
+ newRecHist = vector("list",nNewInd)
+ tmpLastHaplo = private$.lastHaplo
+ if(all(isDH==1L)){
+ for(i in seq_len(nNewInd)){
+ tmpLastHaplo = tmpLastHaplo + 1L
+ newRecHist[[i]] = rep(tmpLastHaplo, ploidy)
+ }
+ }else{
+ for(i in seq_len(nNewInd)){
+ newRecHist[[i]] = (tmpLastHaplo+1L):(tmpLastHaplo+ploidy)
+ tmpLastHaplo = tmpLastHaplo + ploidy
+ }
+ }
+ private$.hasHap = c(private$.hasHap, rep(FALSE, nNewInd))
+ private$.isFounder = c(private$.isFounder, rep(TRUE, nNewInd))
+ #names(newRecHist) = id
+ private$.recHist = c(private$.recHist, newRecHist)
+ private$.lastHaplo = tmpLastHaplo
+ }else{
+ # Add hist to recombination history
+ private$.hasHap = c(private$.hasHap, rep(FALSE, nNewInd))
+ private$.isFounder = c(private$.isFounder, rep(FALSE, nNewInd))
+ #names(hist) = id
+ private$.recHist = c(private$.recHist, hist)
+ }
+ private$.pedigree = rbind(private$.pedigree, tmp)
+ private$.lastId = private$.lastId + as.integer(nNewInd)
+
+ invisible(self)
+ },
+
+ #' @description A function to update the caste
+ #' For internal use only.
+ #'
+ #' @param caste vector, named vector of castes to be added
+ updateCaste = function(caste) {
+ private$.caste = c(private$.caste, caste)
+ invisible(self)
+ },
+
+ #' @description A function to update the last
+ #' ID everytime we create an individual
+ #' For internal use in SIMplyBee only.
+ #'
+ #' @param lastId integer, last colony ID assigned
+ #' @param n integer, how many individuals to add
+ updateLastBeeId = function(n = 1L) {
+ private$.lastId = private$.lastId + as.integer(n)
+ invisible(self)
+ },
+
#' @description A function to update the colony last
#' ID everytime we create a Colony-class with createColony.
#' For internal use only.
#'
#' @param lastColonyId integer, last colony ID assigned
- updateLastColonyId = function() {
- private$.lastColonyId = private$.lastColonyId + 1L
+ #' @param n integer, how many colonies to add
+ updateLastColonyId = function(n = 1) {
+ n = as.integer(n)
+ private$.lastColonyId = private$.lastColonyId + n
invisible(self)
}
),
@@ -459,7 +558,7 @@ SimParamBee <- R6Class(
#' created
caste = function(value) {
if (missing(value)) {
- x = private$.caste
+ x = private$.caste
} else {
stop("`$caste` is read only", call. = FALSE)
}
@@ -469,7 +568,7 @@ SimParamBee <- R6Class(
#' created with \code{\link[SIMplyBee]{createColony}}
lastColonyId = function(value) {
if (missing(value)) {
- private$.lastColonyId
+ private$.lastColonyId
} else {
stop("`$lastColonyId` is read only", call. = FALSE)
}
diff --git a/R/Functions_L0_auxilary.R b/R/Functions_L0_auxilary.R
index 28e56689..65160d15 100644
--- a/R/Functions_L0_auxilary.R
+++ b/R/Functions_L0_auxilary.R
@@ -341,6 +341,37 @@ calcQueensPHomBrood <- function(x, simParamBee = NULL) {
return(ret)
}
+
+calcQueensPHomBrood_parallel <- function(x, simParamBee = NULL) {
+ if (is.null(simParamBee)) {
+ simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ }
+ if (isPop(x)) {
+ ret <- rep(x = NA, times = nInd(x))
+ for (ind in seq_len(nInd(x))) {
+
+ queensCsd <- apply(
+ X = getCsdAlleles(x[ind], simParamBee = simParamBee), MARGIN = 1,
+ FUN = function(x) paste0(x, collapse = "")
+ )
+ fathersCsd <- apply(
+ X = getCsdAlleles(x@misc$fathers[[ind]], simParamBee = simParamBee), MARGIN = 1,
+ FUN = function(x) paste0(x, collapse = "")
+ )
+ nComb <- length(queensCsd) * length(fathersCsd)
+ ret[ind] <- sum(fathersCsd %in% queensCsd) / nComb
+ }
+ } else if (isColony(x)) {
+ ret <- calcQueensPHomBrood(x = x@queen)
+ } else if (isMultiColony(x)) {
+ ret <- sapply(X = x@colonies, FUN = calcQueensPHomBrood)
+ names(ret) <- getId(x)
+ } else {
+ stop("Argument x must be a Pop, Colony, or MultiColony class object!")
+ }
+ return(ret)
+}
+
#' @describeIn calcQueensPHomBrood Expected percentage of csd homozygous brood
#' of a queen / colony
#' @export
@@ -359,11 +390,11 @@ pHomBrood <- function(x, simParamBee = NULL) {
}
}
} else if (isColony(x)) {
- if (is.null(x@queen@misc$pHomBrood[[1]])) {
- ret <- NA
- } else {
- ret <- x@queen@misc$pHomBrood[[1]]
- }
+ if (is.null(x@queen@misc$pHomBrood[[1]])) {
+ ret <- NA
+ } else {
+ ret <- x@queen@misc$pHomBrood[[1]]
+ }
} else if (isMultiColony(x)) {
ret <- sapply(X = x@colonies, FUN = pHomBrood)
names(ret) <- getId(x)
@@ -952,150 +983,6 @@ isNULLColonies <- function(multicolony) {
# get (general) ----
-#' @rdname getQueenYearOfBirth
-#' @title Access the queen's year of birth
-#'
-#' @description Level 0 function that returns the queen's year of birth.
-#'
-#' @param x \code{\link[AlphaSimR]{Pop-class}} (one or more than one queen),
-#' \code{\link[SIMplyBee]{Colony-class}} (one colony), or
-#' \code{\link[SIMplyBee]{MultiColony-class}} (more colonies)
-#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
-#'
-#' @return numeric, the year of birth of the queen(s); named when theres is more
-#' than one queen; \code{NA} if queen not present
-#'
-#' @examples
-#' founderGenomes <- quickHaplo(nInd = 8, nChr = 1, segSites = 100)
-#' SP <- SimParamBee$new(founderGenomes)
-#' \dontshow{SP$nThreads = 1L}
-#' basePop <- createVirginQueens(founderGenomes)
-#'
-#' drones <- createDrones(x = basePop[1], nInd = 1000)
-#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = nFathersPoisson)
-#'
-#' # Create a Colony and a MultiColony class
-#' colony <- createColony(x = basePop[2])
-#' colony <- cross(colony, drones = droneGroups[[1]])
-#'
-#' apiary <- createMultiColony(basePop[3:4], n = 2)
-#' apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-#'
-#' queen <- getQueen(colony)
-#' queen <- setQueensYearOfBirth(queen, year = 2022)
-#' getQueenYearOfBirth(queen)
-#'
-#' getQueenYearOfBirth(getQueen(colony))
-#' colony <- setQueensYearOfBirth(colony, year = 2030)
-#' getQueenYearOfBirth(colony)
-#'
-#' apiary <- setQueensYearOfBirth(apiary, year = 2022)
-#' getQueenYearOfBirth(apiary)
-#' @export
-getQueenYearOfBirth <- function(x, simParamBee = NULL) {
- if (is.null(simParamBee)) {
- simParamBee <- get(x = "SP", envir = .GlobalEnv)
- }
- if (isPop(x)) {
- if (any(!(isVirginQueen(x, simParamBee = simParamBee) | isQueen(x, simParamBee = simParamBee)))) {
- stop("Individuals in x must be virgin queens or queens!")
- }
- nInd <- nInd(x)
- ret <- rep(x = NA, times = nInd)
- for (ind in seq_len(nInd)) {
- if (!is.null(x@misc$yearOfBirth[[ind]])) {
- ret[ind] <- x@misc$yearOfBirth[[ind]]
- }
- }
- if (nInd > 1) {
- names(ret) <- getId(x)
- }
- } else if (isColony(x)) {
- ret <- ifelse(is.null(x@queen@misc$yearOfBirth[[1]]), NA, x@queen@misc$yearOfBirth[[1]])
- } else if (isMultiColony(x)) {
- ret <- sapply(X = x@colonies, FUN = getQueenYearOfBirth, simParamBee = simParamBee)
- names(ret) <- getId(x)
- } else {
- stop("Argument x must be a Pop, Colony, or MultiColony class object!")
- }
- return(ret)
-}
-
-#' @rdname getQueenAge
-#' @title Get (calculate) the queen's age
-#'
-#' @description Level 0 function that returns the queen's age.
-#'
-#' @param x \code{\link[AlphaSimR]{Pop-class}}, \code{\link[SIMplyBee]{Colony-class}}, or
-#' \code{\link[SIMplyBee]{MultiColony-class}}
-#' @param currentYear integer, current year
-#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
-#'
-#' @return numeric, the age of the queen(s); named when theres is more
-#' than one queen; \code{NA} if queen not present
-#'
-#' @examples
-#' founderGenomes <- quickHaplo(nInd = 8, nChr = 1, segSites = 100)
-#' SP <- SimParamBee$new(founderGenomes)
-#' \dontshow{SP$nThreads = 1L}
-#' basePop <- createVirginQueens(founderGenomes)
-#'
-#' drones <- createDrones(x = basePop[1], nInd = 1000)
-#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = nFathersPoisson)
-#'
-#' # Create a Colony and a MultiColony class
-#' colony <- createColony(x = basePop[2])
-#' colony <- cross(colony, drones = droneGroups[[1]])
-#' apiary <- createMultiColony(basePop[3:4], n = 2)
-#' apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-#'
-#' queen <- getQueen(colony)
-#' queen <- setQueensYearOfBirth(queen, year = 2020)
-#' getQueenAge(queen, currentYear = 2022)
-#'
-#' colony <- setQueensYearOfBirth(colony, year = 2021)
-#' getQueenAge(colony, currentYear = 2022)
-#'
-#' apiary <- setQueensYearOfBirth(apiary, year = 2018)
-#' getQueenAge(apiary, currentYear = 2022)
-#' @export
-getQueenAge <- function(x, currentYear, simParamBee = NULL) {
- if (is.null(simParamBee)) {
- simParamBee <- get(x = "SP", envir = .GlobalEnv)
- }
- if (isPop(x)) {
- if (any(!(isVirginQueen(x, simParamBee = simParamBee) | isQueen(x, simParamBee = simParamBee)))) {
- stop("Individuals in x must be virgin queens or queens!")
- }
- nInd <- nInd(x)
- ret <- rep(x = NA, times = nInd)
- for (ind in seq_len(nInd)) {
- if (!is.null(x@misc$yearOfBirth[[ind]])) {
- ret[ind] <- currentYear - x@misc$yearOfBirth[[ind]]
- }
- }
- if (nInd > 1) {
- names(ret) <- getId(x)
- }
- } else if (isColony(x)) {
- if (isQueenPresent(x, simParamBee = simParamBee)) {
- if(packageVersion("AlphaSimR") > package_version("1.5.3")){
- ret <- currentYear - x@queen@misc$yearOfBirth[[1]]
- }else{
- ret <- currentYear - x@queen@misc[[1]]$yearOfBirth
- }
- } else {
- ret <- NA
- }
- } else if (isMultiColony(x)) {
- ret <- sapply(X = x@colonies, FUN = getQueenAge, currentYear = currentYear)
- names(ret) <- getId(x)
- } else {
- stop("Argument x must be a Pop, Colony, or MultiColony class object!")
- }
- return(ret)
-}
-
#' @rdname getId
#' @title Get the colony ID
#'
@@ -1769,7 +1656,7 @@ getEvents <- function(x) {
#' colony <- buildUp(x = colony, nWorkers = 6, nDrones = 3)
#' colony <- addVirginQueens(colony, nInd = 5)
#'
-#' apiary <- createMultiColony(basePop[3:4], n = 2)
+#' apiary <- createMultiColony(basePop[3:4])
#' apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
#' apiary <- buildUp(x = apiary, nWorkers = 6, nDrones = 3)
#'
@@ -2545,7 +2432,7 @@ getCsdGeno <- function(x, caste = NULL, nInd = NULL, dronesHaploid = TRUE,
} else {
ret <- getCsdGeno(
x = getCastePop(x, caste, simParamBee = simParamBee), nInd = nInd,
- dronesHaploid = dronesHaploid, simParamBee = simParamBee
+ dronesHaploid = dronesHaploid, simParamBee = simParamBee
)
}
} else if (isMultiColony(x)) {
@@ -2881,10 +2768,10 @@ nCsdAlleles <- function(x, collapse = FALSE, simParamBee = NULL) {
#' \code{\link[SIMplyBee]{MultiColony-class}}
#'
#' @examples
-#' founderGenomes <- quickHaplo(nInd = 4, nChr = 1, segSites = 50)
-#' SP <- SimParamBee$new(founderGenomes)
+#' founderGenomes <- quickHaplo(nInd = 4, nChr = 1, segSites = 5)
+#' SP <- SimParamBee$new(founderGenomes, nCsdAlleles = 4)
#' \dontshow{SP$nThreads = 1L}
-#' SP$setTrackRec(TRUE)
+#' SP$setTrackRec(isTrackRec = TRUE)
#' SP$setTrackPed(isTrackPed = TRUE)
#' basePop <- createVirginQueens(founderGenomes)
#'
@@ -2894,13 +2781,13 @@ nCsdAlleles <- function(x, collapse = FALSE, simParamBee = NULL) {
#' # Create a Colony and a MultiColony class
#' colony <- createColony(x = basePop[2])
#' colony <- cross(colony, drones = droneGroups[[1]])
-#' colony <- buildUp(x = colony, nWorkers = 6, nDrones = 3)
-#' colony <- addVirginQueens(x = colony, nInd = 5)
+#' colony <- buildUp(x = colony, nWorkers = 3, nDrones = 2)
+#' colony <- addVirginQueens(x = colony, nInd = 2)
#'
-#' apiary <- createMultiColony(basePop[3:4], n = 2)
+#' apiary <- createMultiColony(basePop[3:4])
#' apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-#' apiary <- buildUp(x = apiary, nWorkers = 6, nDrones = 3)
-#' apiary <- addVirginQueens(x = apiary, nInd = 5)
+#' apiary <- buildUp(x = apiary, nWorkers = 3, nDrones = 2)
+#' apiary <- addVirginQueens(x = apiary, nInd = 2)
#'
#' # Input is a population
#' getIbdHaplo(x = getQueen(colony))
@@ -2912,6 +2799,8 @@ nCsdAlleles <- function(x, collapse = FALSE, simParamBee = NULL) {
#' getQueenIbdHaplo(colony)
#'
#' getIbdHaplo(colony, caste = "workers", nInd = 3)
+#' getIbdHaplo(colony, caste = "virginQueens")
+#' getIbdHaplo(colony, caste = "drones")
#' getWorkersIbdHaplo(colony)
#' # Same aliases exist for all castes!
#'
@@ -2926,6 +2815,9 @@ nCsdAlleles <- function(x, collapse = FALSE, simParamBee = NULL) {
#' # Or collapse all the haplotypes into a single matrix
#' getQueenIbdHaplo(apiary, collapse = TRUE)
#'
+#'
+#' getIbdHaplo(x = apiary, caste = "workers")
+#' getIbdHaplo(x = apiary, caste = "drones")
#' # Get the haplotypes of all individuals either by colony or in a single matrix
#' getIbdHaplo(apiary, caste = "all")
#' getIbdHaplo(apiary, caste = "all", collapse = TRUE)
@@ -6338,7 +6230,7 @@ editCsdLocus <- function(pop, alleles = NULL, simParamBee = NULL) {
alleles <- expand.grid(as.data.frame(matrix(rep(0:1, length(csdSites)), nrow = 2, byrow = FALSE)))
# Sample two different alleles (without replacement) for each individual
nAlleles <- simParamBee$nCsdAlleles
- alleles <- sapply(seq_len(pop@nInd), FUN = function(x) list(alleles[sample(nAlleles, size = 2, replace = F), ]))
+ alleles <- sapply(seq_len(pop@nInd), FUN = function(x) list(alleles[sample(nAlleles, size = 2, replace = FALSE), ]))
}
if (pop@nInd != length(alleles)) {
@@ -6406,10 +6298,6 @@ editCsdLocus <- function(pop, alleles = NULL, simParamBee = NULL) {
#' virginColonies2 <- setLocation(virginColonies2,
#' location = Map(c, runif(30, 0, 2*pi),
#' runif(30, 0, 2*pi)))
-#' virginColonies3 <- createMultiColony(basePop[61:90])
-#' virginColonies3 <- setLocation(virginColonies3,
-#' location = Map(c, runif(30, 0, 2*pi),
-#' runif(30, 0, 2*pi)))
#'
#' # Create drone colonies
#' droneColonies <- createMultiColony(basePop[121:200])
@@ -6429,58 +6317,23 @@ editCsdLocus <- function(pop, alleles = NULL, simParamBee = NULL) {
#' nDrones = nFathersPoisson,
#' crossPlan = randomCrossPlan)
#'
-#' # Plot the colonies in space
-#' virginLocations <- as.data.frame(getLocation(c(virginColonies1, virginColonies2, virginColonies3),
-#' collapse= TRUE))
-#' virginLocations$Type <- "Virgin"
-#' droneLocations <- as.data.frame(getLocation(droneColonies, collapse= TRUE))
-#' droneLocations$Type <- "Drone"
-#' locations <- rbind(virginLocations, droneLocations)
-#'
-#' plot(x = locations$V1, y = locations$V2,
-#' col = c("red", "blue")[as.numeric(as.factor(locations$Type))])
-#'
-#' # Cross according to a spatial cross plan according to the colonies' locations
-#' crossPlanSpatial <- createCrossPlan(x = virginColonies1,
-#' droneColonies = droneColonies,
-#' nDrones = nFathersPoisson,
-#' spatial = TRUE,
-#' radius = 1.5)
-#'
-#' # Plot the crossing for the first colony in the crossPlan
-#' virginLocations1 <- as.data.frame(getLocation(virginColonies1, collapse= TRUE))
-#' virginLocations1$Type <- "Virgin"
-#' droneLocations <- as.data.frame(getLocation(droneColonies, collapse= TRUE))
-#' droneLocations$Type <- "Drone"
-#' locations1 <- rbind(virginLocations1, droneLocations)
-#'
-#' # Blue marks the target virgin colony and blue marks the drone colonies in the chosen radius
-#' plot(x = locations1$V1, y = locations1$V2, pch = c(1, 2)[as.numeric(as.factor(locations1$Type))],
-#' col = ifelse(rownames(locations1) %in% crossPlanSpatial[[1]],
-#' "red",
-#' ifelse(rownames(locations1) == names(crossPlanSpatial)[[1]],
-#' "blue", "black")))
-#'
-#' colonies1 <- cross(x = virginColonies1,
-#' crossPlan = crossPlanSpatial,
-#' droneColonies = droneColonies,
-#' nDrones = nFathersPoisson)
-#' nFathers(colonies1)
-#'
#' # Cross according to a cross plan that is created internally within the cross function
#' # The cross plan is created at random, regardless the location of the colonies
-#' colonies2 <- cross(x = virginColonies2,
+#' colonies1 <- cross(x = virginColonies1,
#' droneColonies = droneColonies,
#' nDrones = nFathersPoisson,
#' crossPlan = "create")
#'
-#' # Mate spatially with cross plan created internally by the cross function
-#' colonies3 <- cross(x = virginColonies3,
-#' droneColonies = droneColonies,
+#'
+#' # Cross according to a spatial cross plan created internally according to the colonies' locations
+#' colonies2 <- cross(x = virginColonies2,
#' crossPlan = "create",
-#' checkCross = "warning",
+#' droneColonies = droneColonies,
+#' nDrones = nFathersPoisson,
#' spatial = TRUE,
-#' radius = 1)
+#' radius = 1.5)
+#' nFathers(colonies2)
+#'
#'
#' @export
createCrossPlan <- function(x,
diff --git a/R/Functions_L1_Pop.R b/R/Functions_L1_Pop.R
index 005fe4fd..53c3e047 100644
--- a/R/Functions_L1_Pop.R
+++ b/R/Functions_L1_Pop.R
@@ -1,4 +1,6 @@
# ---- Level 1 Pop Functions ----
+utils::globalVariables("colony")
+utils::globalVariables("i")
#' @rdname getCastePop
#' @title Access individuals of a caste
@@ -218,18 +220,18 @@ getFathers <- function(x, nInd = NULL, use = "rand", collapse = FALSE, simParamB
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
- if (isPop(x)) { # DO WE WANT TO PUT THIS IN getCastePop???
- ret = lapply(X = x@misc$fathers,
- FUN = function(z){
- if(is.null(z)){
- ret = NULL
- }else{
- if (is.null(nInd)) {
- n <- nInd(z)
- }
- ret <- selectInd(pop = z, nInd = n, use = use, simParam = simParamBee)
- }
- }
+ if (isPop(x)) { # TODO: DO WE WANT TO PUT THIS IN getCastePop???
+ ret <- lapply(X = x@misc$fathers,
+ FUN = function(z){
+ if(is.null(z)){
+ ret = NULL
+ }else{
+ if (is.null(nInd)) {
+ n <- nInd(z)
+ }
+ ret <- selectInd(pop = z, nInd = n, use = use, simParam = simParamBee)
+ }
+ }
)
if (nInd(x) == 1) {
ret <- ret[[1]]
@@ -295,11 +297,6 @@ getVirginQueens <- function(x, nInd = NULL, use = "rand", collapse = FALSE, simP
#' only used when \code{x} is \code{\link[SIMplyBee]{Colony-class}} or
#' \code{\link[SIMplyBee]{MultiColony-class}}, when \code{x} is \code{link[AlphaSimR]{MapPop-class}}
#' all individuals in \code{x} are converted into virgin queens
-#' @param exact logical, only relevant when creating workers,
-#' if the csd locus is active and exact is \code{TRUE},
-#' create the exactly specified number of viable workers (heterozygous on the
-#' csd locus)
-#' @param year numeric, year of birth for virgin queens
#' @param editCsd logical (only active when \code{x} is \code{link[AlphaSimR]{MapPop-class}}),
#' whether the csd locus should be edited to ensure heterozygosity at the csd
#' locus (to get viable virgin queens); see \code{csdAlleles}
@@ -313,6 +310,10 @@ getVirginQueens <- function(x, nInd = NULL, use = "rand", collapse = FALSE, simP
#' in \code{\link[SIMplyBee]{SimParamBee}}. The two csd alleles must be different to
#' ensure heterozygosity at the csd locus.
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
+#' @param returnSP logical, whether to return the pedigree, caste, and recHist information
+#' for each created population (used internally for parallel computing)
+#' @param ids character, IDs of the individuals that are going to be created (used internally
+#' for parallel computing)
#' @param ... additional arguments passed to \code{nInd} when this argument is a function
#'
#' @return when \code{x} is \code{\link[AlphaSimR]{MapPop-class}} returns
@@ -327,7 +328,7 @@ getVirginQueens <- function(x, nInd = NULL, use = "rand", collapse = FALSE, simP
#' founderGenomes <- quickHaplo(nInd = 4, nChr = 1, segSites = 50)
#' SP <- SimParamBee$new(founderGenomes)
#' \dontshow{SP$nThreads = 1L}
-#' SP$setTrackRec(TRUE)
+#' SP$setTrackRec(isTrackRec = TRUE)
#' SP$setTrackPed(isTrackPed = TRUE)
#'
#' # Create virgin queens on a MapPop
@@ -398,13 +399,15 @@ getVirginQueens <- function(x, nInd = NULL, use = "rand", collapse = FALSE, simP
# patrilines
# https://github.com/HighlanderLab/SIMplyBee/issues/78
createCastePop <- function(x, caste = NULL, nInd = NULL,
- exact = TRUE, year = NULL,
editCsd = TRUE, csdAlleles = NULL,
simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL,
...) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
+
if (is.null(nInd)) {
if (caste == "virginQueens") {
nInd <- simParamBee$nVirginQueens
@@ -417,6 +420,9 @@ createCastePop <- function(x, caste = NULL, nInd = NULL,
if (is.function(nInd)) {
nInd <- nInd(x, ...)
}
+ if (any(nInd == 0)) {
+ stop("nInd set to 0, should be > 0!")
+ }
# doing "if (is.function(nInd))" below
if (isMapPop(x)) {
if (caste != "virginQueens") { # Creating virgin queens if input is a MapPop
@@ -430,9 +436,7 @@ createCastePop <- function(x, caste = NULL, nInd = NULL,
}
ret@sex[] <- "F"
simParamBee$changeCaste(id = ret@id, caste = "virginQueens")
- if (!is.null(year)) {
- ret <- setQueensYearOfBirth(x = ret, year = year, simParamBee = simParamBee)
- }
+
} else if (isPop(x)) {
if (caste != "drones") { # Creating drones if input is a Pop
stop("Pop-class can only be used to create drones!")
@@ -453,73 +457,153 @@ createCastePop <- function(x, caste = NULL, nInd = NULL,
}
ret <- list()
for (virginQueen in 1:nInd(x)) {
- ret[[virginQueen]] <- makeDH(pop = x[virginQueen], nDH = nInd[virginQueen], keepParents = FALSE, simParam = simParamBee)
+ ret[[virginQueen]] <- makeDH(pop = x[virginQueen], nDH = nInd[virginQueen],
+ keepParents = FALSE, simParam = simParamBee)
}
ret <- mergePops(ret)
}
ret@sex[] <- "M"
simParamBee$addToCaste(id = ret@id, caste = "drones")
} else if (isColony(x)) {
- if (!isQueenPresent(x, simParamBee = simParamBee)) {
- stop("Missing queen!")
- }
+ originalThreads <- simParamBee$nThreads
+ simParamBee$nThreads <- 1
+
if (length(nInd) > 1) {
warning("More than one value in the nInd argument, taking only the first value!")
nInd <- nInd[1]
}
- if (caste == "workers") {
- ret <- vector(mode = "list", length = 2)
- names(ret) <- c("workers", "nHomBrood")
- workers <- combineBeeGametes(
- queen = getQueen(x, simParamBee = simParamBee), drones = getFathers(x, simParamBee = simParamBee),
- nProgeny = nInd, simParamBee = simParamBee
- )
- if (isCsdActive(simParamBee = simParamBee)) {
- sel <- isCsdHeterozygous(pop = workers, simParamBee = simParamBee)
- ret$workers <- workers[sel]
- ret$nHomBrood <- nInd - sum(sel)
- if (exact) {
- if (nInd(ret$workers) < nInd) {
- nMiss <- nInd - nInd(ret$workers)
- while (0 < nMiss) {
- workers <- combineBeeGametes(
- queen = getQueen(x, simParamBee = simParamBee),
- drones = getFathers(x, simParamBee = simParamBee),
- nProgeny = nMiss,
- simParamBee = simParamBee
- )
- sel <- isCsdHeterozygous(pop = workers, simParamBee = simParamBee)
- ret$workers <- c(ret$workers, workers[sel])
- ret$nHomBrood <- ret$nHomBrood + sum(!sel)
- nMiss <- nInd - nInd(ret$workers)
+ if (!isQueenPresent(x, simParamBee = simParamBee)) {
+ stop("Missing queen!")
+ }
+
+ if (nInd > 0) {
+ if (caste == "workers") {
+ if (!returnSP) {
+ ret <- vector(mode = "list", length = 2)
+ names(ret) <- c("workers", "nHomBrood")
+ } else {
+ ret <- vector(mode = "list", length = 5)
+ names(ret) <- c("workers", "nHomBrood", "pedigree", "caste", "recHist")
+ }
+ simParamBee$nThreads <- 1
+ ret$workers <- combineBeeGametes(
+ queen = getQueen(x, simParamBee = simParamBee),
+ drones = getFathers(x, simParamBee = simParamBee),
+ nProgeny = nInd,
+ simParamBee = simParamBee
+ )
+
+
+
+ simParamBee$addToCaste(id = ret$workers@id, caste = "workers")
+ ret$workers@sex[] <- "F"
+
+ if (returnSP) {
+ ret$caste <- simParamBee$caste[ret$workers@id, drop = FALSE]
+ if (simParamBee$isTrackPed) {
+ ret$pedigree <- simParamBee$pedigree[ret$workers@id, , drop = FALSE]
+ }
+ if (simParamBee$isTrackRec) {
+ ret$recHist <- simParamBee$recHist[ret$workers@iid]
+ }
+ }
+
+ if (!is.null(ids)) {
+ if (nInd(ret$workers) > length(ids)) {
+ stop("Not enough IDs provided")
+ }
+ if (nInd(ret$workers) < length(ids)) {
+ stop("Too many IDs provided!")
+ }
+ ret$workers@id <- as.character(ids)
+ ret$workers@iid <- as.integer(ids)
+ if (returnSP) {
+ names(ret$caste) <- as.character(ids)
+ if (simParamBee$isTrackPed) {
+ rownames(ret$pedigree) <- ids
+ }
+ if (simParamBee$isTrackRec) {
+ names(ret$recHist) <- ids
}
}
}
- } else {
- ret$workers <- workers
- ret$nHomBrood <- NA
- }
- ret$workers@sex[] <- "F"
- simParamBee$addToCaste(id = ret$workers@id, caste = "workers")
- } else if (caste == "virginQueens") { # Creating virgin queens if input is a Colony
- ret <- createCastePop(x = x, caste = "workers",
- nInd = nInd, exact = TRUE, simParamBee = simParamBee)$workers
- ret@sex[] <- "F"
- simParamBee$changeCaste(id = ret@id, caste = "virginQueens")
- if (!is.null(year)) {
- ret <- setQueensYearOfBirth(x = ret, year = year, simParamBee = simParamBee)
+
+ if (isCsdActive(simParamBee = simParamBee)) {
+ sel <- isCsdHeterozygous(pop = ret$workers, simParamBee = simParamBee)
+ ret$nHomBrood <- nInd(ret$workers) - sum(sel)
+ ret$workers <- ret$workers[sel]
+ } else {
+ ret$nHomBrood <- NA
+ }
+
+ } else if (caste == "virginQueens") {
+ ret <- createCastePop(x = x, caste = "workers",
+ nInd = nInd, simParamBee = simParamBee,
+ returnSP = returnSP, ids = ids, ...)
+ ret$caste = rep("virginQueens", length(ret$caste))
+ names(ret$caste) = ids
+ simParamBee$changeCaste(id = ret$workers@id, caste = "virginQueens")
+ if (!returnSP) {
+ ret <- ret$workers
+ }
+
+ } else if (caste == "drones") {
+ drones <- makeDH(
+ pop = getQueen(x, simParamBee = simParamBee), nDH = nInd, keepParents = FALSE,
+ simParam = simParamBee
+ )
+
+ drones@sex[] <- "M"
+ simParamBee$addToCaste(id = drones@id, caste = "drones")
+
+ if (returnSP) {
+ ret <- vector(mode = "list", length = 4)
+ names(ret) <- c("drones", "pedigree", "caste", "recHist")
+ ret$caste <- simParamBee$caste[drones@id, drop = FALSE]
+ if (simParamBee$isTrackPed) {
+ ret$pedigree <- simParamBee$pedigree[drones@id, , drop = FALSE]
+ }
+ if (simParamBee$isTrackRec) {
+ ret$recHist <- simParamBee$recHist[drones@iid]
+ }
+ }
+
+ if (!is.null(ids)) {
+ if (nInd(drones) != length(ids)) {
+ stop("Not enough IDs provided")
+ }
+ drones@id = ids
+ drones@iid = as.integer(ids)
+ if (returnSP) {
+ names(ret$caste) <- ids
+ if (simParamBee$isTrackPed) {
+ rownames(ret$pedigree) <- ids
+ }
+ if (simParamBee$isTrackRec) {
+ names(ret$recHist) <- ids
+ }
+ }
+ }
+
+ if (returnSP) {
+ ret$drones <-drones
+ } else {
+ ret <- drones
+ }
}
- } else if (caste == "drones") { # Creating drones if input is a Colony
- ret <- makeDH(
- pop = getQueen(x, simParamBee = simParamBee), nDH = nInd, keepParents = FALSE,
- simParam = simParamBee
- )
- ret@sex[] <- "M"
- simParamBee$addToCaste(id = ret@id, caste = "drones")
+ } else {
+ ret <- NULL
}
+ simParamBee$nThreads <- originalThreads
} else if (isMultiColony(x)) {
+ if (is.null(nInd)) {
+ string = paste0("n", toupper(substr(caste, 1, 1)), substr(caste, 2, nchar(caste)))
+ nInd <- simParamBee[[string]]
+ }
+
nCol <- nColonies(x)
nNInd <- length(nInd)
+
if (nNInd > 1 && nNInd < nCol) {
stop("Too few values in the nInd argument!")
}
@@ -527,56 +611,147 @@ createCastePop <- function(x, caste = NULL, nInd = NULL,
warning(paste0("Too many values in the nInd argument, taking only the first ", nCol, "values!"))
nInd <- nInd[1:nCol]
}
- ret <- vector(mode = "list", length = nCol)
- for (colony in seq_len(nCol)) {
- if (is.null(nInd)) {
- nIndColony <- NULL
+
+ nNInd <- length(nInd)
+ totalNInd = ifelse(nNInd == 1, nInd * nCol, sum(nInd))
+ if (totalNInd == 0) {
+ stop("Nothing to create.")
+ }
+
+ registerDoParallel(cores = simParamBee$nThreads)
+
+ lastId = simParamBee$lastId
+ ids = (lastId+1):(lastId+totalNInd)
+
+ combine_list <- function(a, b) {
+ if (!is.null(names(a))) {
+ "Combine first"
+ c(list(a), list(b))
} else {
- nIndColony <- ifelse(nNInd == 1, nInd, nInd[colony])
+ if ((is.null(a) | is.null(b)) & !(is.null(a) & is.null(b))) {
+ c(a, list(b))
+ } else if (is.null(a) & is.null(b)) {
+ c(list(a), list(b))
+ } else {
+ c(a, list(b))
+ }
}
- ret[[colony]] <- createCastePop(
- x = x[[colony]], caste = caste,
- nInd = nIndColony,
- exact = exact,
- year = year,
- editCsd = TRUE, csdAlleles = NULL,
- simParamBee = simParamBee, ...
- )
+ }
+ ret <- foreach(colony = seq_len(nCol), .combine=combine_list) %dopar% {
+ nIndColony <- ifelse(nNInd == 1, nInd, nInd[colony])
+ if (nIndColony > 0) {
+ if (nNInd == 1) {
+ colonyIds = ids[((colony-1)*nIndColony+1):(colony*nIndColony)]
+ } else {
+ colonyIds = base::split(ids, rep(seq_along(nInd), nInd))[[as.character(colony)]]
+ }
+ createCastePop(
+ x = x[[colony]], caste = caste,
+ nInd = nIndColony,
+ editCsd = TRUE, csdAlleles = NULL,
+ simParamBee = simParamBee,
+ returnSP = TRUE,
+ ids = as.character(colonyIds)
+ )
+ } else {
+ NULL
+ }
+ }
+ if (nCol == 1) {
+ ret <- list(ret)
}
names(ret) <- getId(x)
- } else {
+ # Add to simParamBee: pedigree, caste, recHist
+ notNull = sapply(ret, FUN = function(x) !is.null(x))
+
+ if (!simParamBee$isTrackPed) {
+ simParamBee$updateLastBeeId(n = totalNInd)
+ } else if (simParamBee$isTrackPed) {
+ Pedigree <- do.call("rbind", lapply(ret[notNull], '[[', "pedigree"))
+ if (!simParamBee$isTrackRec) {
+ simParamBee$addToBeePed(nNewInd = totalNInd, id = rownames(Pedigree),
+ mother = Pedigree[, 'mother'], father = Pedigree[, 'father'],
+ isDH = Pedigree[, 'isDH'])
+ #simParamBee$updatePedigree(pedigree = Pedigree)
+ } else {
+ RecHist = do.call("c", lapply(ret[notNull], '[[', "recHist"))
+ if (caste == "drones") {
+ ploidy = rep(1, totalNInd)
+ } else {
+ ploidy = rep(2, totalNInd)
+ }
+ simParamBee$addToBeeRec(nNewInd = totalNInd, id = rownames(Pedigree),
+ mother = Pedigree[, 'mother'], father = Pedigree[, 'father'],
+ isDH = Pedigree[, 'isDH'],
+ hist = RecHist, ploidy = ploidy)
+ }
+ }
+ # Extend caste
+ Caste <- do.call("c", lapply(ret[notNull], '[[', "caste"))
+ Names <- do.call("c", lapply(ret[notNull], function(x) names(x$caste)))
+ names(Caste) <- Names
+ simParamBee$updateCaste(caste = Caste)
+
+
+ if (!returnSP) {
+ if (caste %in% c("drones", "virginQueens")) {
+ ret <- lapply(ret, FUN = function(x) {
+ if (is.null(x)) return(NULL) # Return NULL if the element is NULL
+ x[!names(x) %in% c("pedigree", "caste", "recHist")][[1]]
+ })
+ } else {
+ ret <- lapply(ret, FUN = function(x) {
+ if (is.null(x)) return(NULL)
+ x[!names(x) %in% c("pedigree", "caste", "recHist")]
+ })
+ }
+ }
+ }
+ else {
stop("Argument x must be a Map-Pop (only for virgin queens),
Pop (only for drones), Colony, or MultiColony class object!")
}
+
return(ret)
}
#' @describeIn createCastePop Create workers from a colony
#' @export
-createWorkers <- function(x, nInd = NULL, exact = FALSE, simParamBee = NULL, ...) {
+createWorkers <- function(x, nInd = NULL, simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL, ...) {
ret <- createCastePop(x, caste = "workers", nInd = nInd,
- exact = exact, simParamBee = simParamBee, ...)
+ simParamBee = simParamBee,
+ returnSP = returnSP,
+ ids = ids, ...)
return(ret)
}
#' @describeIn createCastePop Create drones from a colony
#' @export
-createDrones <- function(x, nInd = NULL, simParamBee = NULL, ...) {
+createDrones <- function(x, nInd = NULL, simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL, ...) {
ret <- createCastePop(x, caste = "drones", nInd = nInd,
- simParamBee = simParamBee, ...)
+ simParamBee = simParamBee,
+ returnSP = returnSP,
+ ids = ids, ...)
return(ret)
}
#' @describeIn createCastePop Create virgin queens from a colony
#' @export
createVirginQueens <- function(x, nInd = NULL,
- year = NULL,
editCsd = TRUE, csdAlleles = NULL,
simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL,
...) {
ret <- createCastePop(x, caste = "virginQueens", nInd = nInd,
- year = year, editCsd = editCsd,
- csdAlleles = csdAlleles, simParamBee = simParamBee, ...)
+ editCsd = editCsd,
+ csdAlleles = csdAlleles, simParamBee = simParamBee,
+ returnSP = returnSP,
+ ids = ids, ...)
return(ret)
}
@@ -817,7 +992,7 @@ createDCA <- function(x, nInd = NULL, removeFathers = TRUE, simParamBee = NULL)
#' SP <- SimParamBee$new(founderGenomes)
#' \dontshow{SP$nThreads = 1L}
#' basePop <- createVirginQueens(founderGenomes)
-#' drones <- createDrones(basePop[1], n = 1000)
+#' drones <- createDrones(basePop[1], nInd = 1000)
#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = 10)
#'
#' # Create a colony and cross it
@@ -1011,12 +1186,12 @@ pullDroneGroupsFromDCA <- function(DCA, n, nDrones = NULL,
#' # Create a Colony and a MultiColony class
#' colony <- createColony(x = basePop[2])
#' colony <- cross(colony, drones = droneGroups[[1]])
-#' colony <- buildUp(x = colony, nWorkers = 100, nDrones = 10, exact = TRUE)
+#' colony <- buildUp(x = colony, nWorkers = 100, nDrones = 10)
#' colony <- addVirginQueens(x = colony, nInd = 3)
#'
#' apiary <- createMultiColony(basePop[3:4], n = 2)
#' apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-#' apiary <- buildUp(x = apiary, nWorkers = 100, nDrones = 10, exact = TRUE)
+#' apiary <- buildUp(x = apiary, nWorkers = 100, nDrones = 10)
#' apiary <- addVirginQueens(x = apiary, nInd = 3)
#'
#' # pullCastePop on Colony class
@@ -1078,6 +1253,7 @@ pullCastePop <- function(x, caste, nInd = NULL, use = "rand",
ret <- list(pulled = tmp$pulled, remnant = x)
}
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
nNInd <- length(nInd)
if (nNInd > 1 && nNInd < nCol) {
@@ -1092,24 +1268,24 @@ pullCastePop <- function(x, caste, nInd = NULL, use = "rand",
ret$pulled <- vector(mode = "list", length = nCol)
names(ret$pulled) <- getId(x)
ret$remnant <- x
- for (colony in seq_len(nCol)) {
+
+ tmp = foreach(colony = seq_len(nCol)) %dopar% {
if (is.null(nInd)) {
nIndColony <- NULL
} else {
nIndColony <- ifelse(nNInd == 1, nInd, nInd[colony])
}
- tmp <- pullCastePop(x = x[[colony]],
- caste = caste,
- nInd = nIndColony,
- use = use,
- removeFathers = removeFathers,
- collapse = collapse,
- simParamBee = simParamBee)
- if (!is.null(tmp$pulled)) {
- ret$pulled[[colony]] <- tmp$pulled
- }
- ret$remnant[[colony]] <- tmp$remnant
+ pullCastePop(x = x[[colony]],
+ caste = caste,
+ nInd = nIndColony,
+ use = use,
+ removeFathers = removeFathers,
+ collapse = collapse,
+ simParamBee = simParamBee)
}
+ ret$pulled <- lapply(tmp, '[[', "pulled")
+ ret$remnant@colonies <- lapply(tmp, '[[', "remnant")
+
if (collapse) {
ret$pulled <- mergePops(ret$pulled)
}
@@ -1195,7 +1371,8 @@ pullVirginQueens <- function(x, nInd = NULL, use = "rand", collapse = FALSE, sim
#' to their distance from the virgin colony (that is, in a radius)
#' @param radius numeric, the radius around the virgin colony in which to sample mating partners,
#' only needed when \code{spatial = TRUE}
-#' @param checkCross character, throw a warning (when \code{checkCross = "warning"}),
+#' @param checkCross character, throw a warning (when \code{checkCross = "warning"}).
+#' This will also remove the unmated queens and return only the mated ones.
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#' @param ... other arguments for \code{nDrones}, when \code{nDrones} is a function
#'
@@ -1224,7 +1401,7 @@ pullVirginQueens <- function(x, nInd = NULL, use = "rand", collapse = FALSE, sim
#' @examples
#' founderGenomes <- quickHaplo(nInd = 30, nChr = 1, segSites = 100)
#' SP <- SimParamBee$new(founderGenomes)
-#' \dontshow{SP$nThreads = 1L}
+#' SP$nThreads = 1L
#' basePop <- createVirginQueens(founderGenomes)
#'
#' drones <- createDrones(x = basePop[1], nInd = 1000)
@@ -1347,14 +1524,21 @@ cross <- function(x,
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
+ registerDoParallel(cores = simParamBee$nThreads)
+
+ if (isPop(x)) {
+ type = "Pop"
+ } else if (isColony(x)) {
+ type = "Colony"
+ } else if (isMultiColony(x)) {
+ type = "MultiColony"
+ } else {
+ stop("Input x must be a Pop-class, Colony-class, or MultiColony-class!")
+ }
+
if (is.null(nDrones)) {
nDrones <- simParamBee$nFathers
}
- if (is.function(nDrones)) {
- nD <- nDrones(...)
- } else {
- nD <- nDrones
- }
IDs <- as.character(getId(x))
oneColony <- (isPop(drones)) && (length(IDs) == 1) && (is.null(crossPlan))
@@ -1366,6 +1550,10 @@ cross <- function(x,
# Do all the tests here to simplify the function
+ if (is.null(crossPlan) & (length(IDs) > 1) & isPop(drones)) {
+ stop("When supplying drones as a single population for mating multiple virgin queens,
+ crossPlan argument must be set to 'create' to internally create a mating plan!")
+ }
if (crossPlan_droneID && !isPop(drones)) {
stop("When using a cross plan, drones must be supplied as a single Pop-class!")
}
@@ -1382,6 +1570,9 @@ cross <- function(x,
stop("The argument drones must be a Pop-class
or a list of drone Pop-class objects!")
}
+ if (isPop(drones) && nInd(drones) == 0) {
+ stop("Argument drones is a Pop-class with 0 individuals!")
+ }
if (crossPlan_given && !is.null(drones) && !all(unlist(crossPlan) %in% drones@id)) {
stop("Some drones from the crossPlan are missing in the drones population!")
}
@@ -1408,195 +1599,205 @@ cross <- function(x,
}
}
+ if (crossPlan_create | crossPlan_given) {
+ if (crossPlan_create) {
+ crossPlan <- createCrossPlan(x = x,
+ drones = drones,
+ droneColonies = droneColonies,
+ nDrones = nDrones,
+ spatial = spatial,
+ radius = radius,
+ simParamBee = simParamBee)
+ }
- if (crossPlan_create) {
- crossPlan <- createCrossPlan(x = x,
- drones = drones,
- droneColonies = droneColonies,
- nDrones = nDrones,
- spatial = spatial,
- radius = radius,
- simParamBee = simParamBee)
noMatches <- sapply(crossPlan, FUN = length)
+ if (all(noMatches == 0)) {
+ stop("All crossings failed!")
+ }
if (0 %in% noMatches) {
- message("There are no potential crosses for some colonies! The cross() will fail
- unless argument checkCross is set to 'warning'.")
+ if (checkCross == "warning") {
+ message("Crossing failed, unmated virgin queens will be removed!")
+ ret <- x
+ } else if (checkCross == "error") {
+ stop("Crossing failed!")
+ }
}
}
- if (isPop(x) | isColony(x)) {
- ret <- list()
- for (virgin in seq_len(length(IDs))) {
- virginID <- IDs[virgin]
- if (oneColony) {
- virginQueenDrones <- drones
- } else if (dronePackages) {
- virginQueenDrones <- drones[[virgin]]
- } else if (crossPlan_given | crossPlan_create) {
- if (crossPlan_droneID) {
- virginQueenDrones <- drones[crossPlan[[virginID]]]
- } else if (crossPlan_colonyID) {
- virginMatches <- crossPlan[[virginID]]
- if (length(virginMatches) > 0) {
- nD <- ifelse(is.function(nDrones), nDrones(...), nDrones)
- selectedDPQ <- table(sample(virginMatches, size = nD, replace = TRUE))
- virginQueenDrones <- mergePops(createDrones(droneColonies[names(selectedDPQ)],
- nInd = selectedDPQ, simParamBee = simParamBee))
- } else {
- virginQueenDrones <- new("Pop")
- }
- }
- }
- if (any((virginQueenDrones@nInd == 0), (length(virginQueenDrones@nInd) == 0))) {
- msg <- "Crossing failed!"
- if (checkCross == "warning") {
- message(msg)
- ret <- x
- } else if (checkCross == "error") {
- stop(msg)
- }
- } else if (virginQueenDrones@nInd > 0) {
- if (!all(isDrone(virginQueenDrones, simParamBee = simParamBee))) {
- stop("Individuals in drones must be drones!")
- }
- if (isPop(x)) {
- virginQueen <- x[virgin]
- } else if (isColony(x)) {
- virginQueen <- selectInd(x@virginQueens, nInd = 1, use = "rand", simParam = simParamBee)
- }
+ # Convert everything to a Pop
+ if (isColony(x) | isMultiColony(x)) {
+ inputId <- getId(x)
+ if (isColony(x)) {
+ colony <- x
+ x <- pullCastePop(x, caste = "virginQueens", nInd = 1)$pulled
+ ID_by_input <- data.frame(inputId = inputId,
+ virginId = getId(x))
+ } else if (isMultiColony(x)) {
+ multicolony <- x
+ x <- pullCastePop(x, caste = "virginQueens", nInd = 1)$pulled
+ ID_by_input <- data.frame(inputId = inputId,
+ virginId = unlist(sapply(x, FUN = function(y) getId(y))))
+ x <- mergePops(x)
+ }
+ # Rename crossPlan
+ if (crossPlan_create | crossPlan_given) {
+ names(crossPlan) <- ID_by_input$virginId[match(ID_by_input$inputId, names(crossPlan))]
+ }
+ }
- virginQueen@misc$fathers[[1]] <- virginQueenDrones
+ IDs <- as.character(getId(x))
+ #Now x is always a Pop
+ ret <- list()
+ nVirgin = nInd(x)
- simParamBee$changeCaste(id = virginQueen@id, caste = "queen")
- simParamBee$changeCaste(id = virginQueenDrones@id, caste = "fathers")
- virginQueen <- setMisc(x = virginQueen, node = "nWorkers", value = 0)
- virginQueen <- setMisc(x = virginQueen, node = "nDrones", value = 0)
- virginQueen <- setMisc(x = virginQueen, node = "nHomBrood", value = 0)
- if (isCsdActive(simParamBee = simParamBee)) {
- val <- calcQueensPHomBrood(x = virginQueen, simParamBee = simParamBee)
- } else {
- val <- NA
- }
- virginQueen <- setMisc(x = virginQueen, node = "pHomBrood", value = val)
-
- if (isPop(x)) {
- ret[[virgin]] <- virginQueen
- } else if (isColony(x)) {
- x <- reQueen(x = x, queen = virginQueen, simParamBee = simParamBee)
- x <- removeVirginQueens(x, simParamBee = simParamBee)
- ret <- x
- }
- }
- }
- if (isPop(x)) {
- ret <- mergePops(ret)
- }
- } else if (isMultiColony(x)) {
- nCol <- nColonies(x)
- if (nCol == 0) {
- ret <- createMultiColony(simParamBee = simParamBee)
+
+ if (is.function(nDrones)) {
+ nD = nDrones(n = nVirgin, ...)
+ } else {
+ nD = nDrones
+ }
+
+ if (length(IDs) > 0 & length(nD) == 1) {
+ nD = rep(nD, length(IDs))
+ }
+
+
+ combine_list <- function(a, b) {
+ if (isPop(a)) {
+ c(list(a), list(b))
} else {
- ret <- createMultiColony(n = nCol, simParamBee = simParamBee)
- for (colony in seq_len(nCol)) {
- if (oneColony) {
- colonyDrones <- drones
- } else if (dronePackages) {
- colonyDrones <- drones[[colony]]
- } else {
- if (crossPlan_colonyID) {
- colonyDrones <- NULL
- } else if(crossPlan_droneID) {
- colonyDrones <- drones
- }
- }
- ret[[colony]] <- cross(
- x = x[[colony]],
- drones = colonyDrones,
- crossPlan = crossPlan,
- droneColonies = droneColonies,
- nDrones = nDrones,
- spatial = spatial,
- radius = radius,
- checkCross = checkCross,
- simParamBee = simParamBee
- )
+ c(a, list(b))
+ }
+ }
+
+
+ if (crossPlan_given | crossPlan_create) {
+ if (crossPlan_colonyID) { # WHAT IF ONE ELEMENT IS EMPTY
+ # This is the crossPlan - for spatial, these are all DPCs found in a radius
+ crossPlanDF <- data.frame(virginID = rep(names(crossPlan), unlist(sapply(crossPlan, length))),
+ DPC = unlist(crossPlan))
+ # If some of the crossing would fail, we only return the queens that mated successfully
+ IDs = IDs[IDs %in% crossPlanDF$virginID]
+ x = x[IDs]
+ if (type == "MultiColony") {
+ multicolony <- multicolony[getId(getVirginQueens(multicolony, collapse=TRUE)) %in% IDs]
+ }
+ # Here we sample from the DPC in the cross plan to get the needed number of drones (nD)
+ crossPlanDF_sample <- do.call("rbind", lapply(IDs,
+ FUN = function(x) {
+ data.frame(virginID = x, DPC = sample(crossPlan[[x]], size = nD[which(x == IDs)], replace = TRUE))
+ } ))
+ crossPlanDF_sample <- crossPlanDF_sample[order(as.integer(crossPlanDF_sample$DPC)),]
+ # Here I gather how many drones each DPC needs to produce
+ crossPlanDF_DPCtable <- as.data.frame(table(crossPlanDF_sample$DPC))
+ crossPlanDF_DPCtable <- crossPlanDF_DPCtable[order(as.integer(as.character(crossPlanDF_DPCtable$Var1))),]
+ colnames(crossPlanDF_DPCtable) <- c("DPC", "noDrones")
+ # Here I select only the DPCs that have been sampled to produce drones
+ selectedDPC = selectColonies(droneColonies, ID = as.character(crossPlanDF_DPCtable$DPC))
+ # And here I create the drones
+ dronesByDPC <- createCastePop(selectedDPC, caste = "drones",
+ nInd = as.integer(crossPlanDF_DPCtable$noDrones),
+ simParamBee = simParamBee)
+ # This is where I link the drone ID to the DPC ID
+ dronesByDPC_DF <- data.frame(DPC = rep(names(dronesByDPC), as.vector(crossPlanDF_DPCtable$noDrones)),
+ droneID = unlist(sapply(dronesByDPC, FUN = function(x) getId(x))))
+ dronesByDPC_DF <- dronesByDPC_DF[order(as.integer(dronesByDPC_DF$DPC)),]
+ dronePop = mergePops(dronesByDPC)
+
+ if (any(!crossPlanDF_sample$DPC == dronesByDPC_DF$DPC)) {
+ stop("Something went wrong with cross plan - drone matching!")
+ }
+
+ dronesByVirgin_DF <- cbind(dronesByDPC_DF, crossPlanDF_sample[, c("virginID"), drop = FALSE])
+ dronesByVirgin_DF <- dronesByVirgin_DF[order(as.integer(dronesByVirgin_DF$virginID)),]
+ dronesByVirgin_list <- lapply(IDs,
+ FUN = function(x) dronesByVirgin_DF$droneID[dronesByVirgin_DF$virginID == x])
+ names(dronesByVirgin_list) <- IDs
+
+ dronesByVirgin <- foreach(i = IDs, .combine = combine_list) %dopar% {
+ dronePop[as.character(dronesByVirgin_list[[i]])]
+ }
+ } else if (crossPlan_droneID) {
+ dronesByVirgin <- foreach(i = IDs, .combine = combine_list) %dopar% {
+ drones[as.character(crossPlan[[i]])]
}
}
}
- validObject(ret)
- return(ret)
-}
-#' @rdname setQueensYearOfBirth
-#' @title Set the queen's year of birth
-#'
-#' @description Level 1 function that sets the queen's year of birth.
-#'
-#' @param x \code{\link[AlphaSimR]{Pop-class}} (one or more than one queen),
-#' \code{\link[SIMplyBee]{Colony-class}} (one colony), or
-#' \code{\link[SIMplyBee]{MultiColony-class}} (more colonies)
-#' @param year integer, the year of the birth of the queen
-#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
-#'
-#' @return \code{\link[AlphaSimR]{Pop-class}}, \code{\link[SIMplyBee]{Colony-class}}, or
-#' \code{\link[SIMplyBee]{MultiColony-class}} with queens having the year of birth set
-#'
-#' @examples
-#' founderGenomes <- quickHaplo(nInd = 8, nChr = 1, segSites = 100)
-#' SP <- SimParamBee$new(founderGenomes)
-#' \dontshow{SP$nThreads = 1L}
-#' basePop <- createVirginQueens(founderGenomes)
-#'
-#' drones <- createDrones(x = basePop[1], nInd = 1000)
-#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = nFathersPoisson)
-#'
-#' # Create a Colony and a MultiColony class
-#' colony <- createColony(x = basePop[2])
-#' colony <- cross(x = colony, drones = droneGroups[[1]])
-#' apiary <- createMultiColony(basePop[3:4], n = 2)
-#' apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-#'
-#' # Example on Colony class
-#' getQueenYearOfBirth(colony)
-#' getQueenYearOfBirth(apiary)
-#'
-#' queen1 <- getQueen(colony)
-#' queen1 <- setQueensYearOfBirth(queen1, year = 2022)
-#' getQueenYearOfBirth(queen1)
-#'
-#' colony <- setQueensYearOfBirth(colony, year = 2022)
-#' getQueenYearOfBirth(colony)
-#'
-#' apiary <- setQueensYearOfBirth(apiary, year = 2022)
-#' getQueenYearOfBirth(apiary)
-#' @export
-setQueensYearOfBirth <- function(x, year, simParamBee = NULL) {
- if (is.null(simParamBee)) {
- simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ # At this point, x is a Pop and dronesByVirgin are a list (so are they if they come if via drone packages)
+ if (oneColony) {
+ dronesByVirgin <- list(drones)
}
- if (isPop(x)) {
- if (any(!(isVirginQueen(x, simParamBee = simParamBee) | isQueen(x, simParamBee = simParamBee)))) {
- stop("Individuals in x must be virgin queens or queens!")
- }
- nInd <- nInd(x)
- x <- setMisc(x = x, node = "yearOfBirth", value = year)
- } else if (isColony(x)) {
- if (isQueenPresent(x, simParamBee = simParamBee)) {
- x@queen <- setMisc(x = x@queen, node = "yearOfBirth", value = year)
+ if (dronePackages) {
+ dronesByVirgin <- drones
+ }
+
+ names(dronesByVirgin) <- IDs
+ nDronesByVirgin <- sapply(dronesByVirgin, FUN = function(x) nInd(x))
+
+
+ #if (all(nDronesByVirgin > 0)) { #There was a mistake here - if the message is warning, this still needs to happen
+ if (!all(sapply(dronesByVirgin,
+ FUN = function(x) all(isDrone(x, simParamBee = simParamBee))))) {
+ stop("Individuals in drones must be drones!")
+ }
+
+ if (nInd(x) != length(dronesByVirgin)) {
+ stop("Number of virgin queens does not match the length of the assigned drones!")
+ }
+
+ for (id in IDs) {
+ simParamBee$changeCaste(id = id, caste = "queen")
+ }
+
+ for (id in as.vector(Reduce("c", sapply(dronesByVirgin, FUN = function(x) getId(x))))) {
+ simParamBee$changeCaste(id = id, caste = "fathers")
+ }
+
+ # All of the input has been transformed to a Pop
+ crossVirginQueen <- function(virginQueen, virginQueenDrones, simParamBee = NULL) {
+ virginQueen@misc$fathers[[1]] <- virginQueenDrones
+ virginQueen <- setMisc(x = virginQueen, node = "nWorkers", value = 0)
+ virginQueen <- setMisc(x = virginQueen, node = "nDrones", value = 0)
+
+ virginQueen <- setMisc(x = virginQueen, node = "nHomBrood", value = 0)
+ # if (isCsdActive(simParamBee = simParamBee)) { #This does still not work it the CSD is turned on
+ # val <- calcQueensPHomBrood(x = virginQueen, simParamBee = simParamBee)
+ # } else {
+ # val <- NA
+ # }
+ #
+ # virginQueen <- setMisc(x = virginQueen, node = "pHomBrood", value = val)
+ return(virginQueen)
+ }
+
+ # Add drones in the queens father slot
+ x <- foreach(i = 1:length(IDs), .combine = combine_list) %dopar% {
+ crossVirginQueen(virginQueen = x[i],
+ virginQueenDrones = dronesByVirgin[[i]],
+ simParamBee = simParamBee)
+ }
+
+
+ if (type == "Pop") {
+ if (length(x) == 1) {
+ ret <- x
} else {
- stop("Missing queen!")
+ ret <- mergePops(x)
}
- } else if (isMultiColony(x)) {
- nCol <- nColonies(x)
- for (colony in seq_len(nCol)) {
- x[[colony]]@queen <- setMisc(
- x = x[[colony]]@queen, node = "yearOfBirth",
- value = year
- )
+ } else if (type == "Colony") {
+ ret <- reQueen(x = colony, queen = x[1], simParamBee = simParamBee)
+ ret <- removeVirginQueens(ret, simParamBee = simParamBee)
+ } else if (type == "MultiColony") {
+ if (length(IDs) > 1) {
+ x <- mergePops(x)
}
- } else {
- stop("Argument x must be a Pop, Colony or MultiColony class object!")
+ ret <- reQueen(x = multicolony, queen = x, simParamBee = simParamBee)
+ ret <- removeCastePop(ret, caste = "virginQueens", simParamBee = simParamBee)
}
- return(x)
+
+ validObject(ret)
+ return(ret)
}
+
diff --git a/R/Functions_L2_Colony.R b/R/Functions_L2_Colony.R
index 390e76f5..78ace6a9 100644
--- a/R/Functions_L2_Colony.R
+++ b/R/Functions_L2_Colony.R
@@ -8,6 +8,7 @@
#'
#' @param x \code{\link[AlphaSimR]{Pop-class}}, one queen or virgin queen(s)
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
+#' @param id character, ID of the colony that is going to be created (used internally for parallel computing)
#'
#' @return new \code{\link[SIMplyBee]{Colony-class}}
#'
@@ -31,15 +32,19 @@
#' colony1 <- cross(colony1, drones = drones)
#' colony1
#' @export
-createColony <- function(x = NULL, simParamBee = NULL) {
+createColony <- function(x = NULL, simParamBee = NULL, id = NULL) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
simParamBee$updateLastColonyId()
+ if (is.null(id)) {
+ id <- simParamBee$lastColonyId
+ }
+
if (is.null(x)) {
colony <- new(
Class = "Colony",
- id = simParamBee$lastColonyId
+ id = id
)
} else {
if (!isPop(x)) {
@@ -60,7 +65,7 @@ createColony <- function(x = NULL, simParamBee = NULL) {
colony <- new(
Class = "Colony",
- id = simParamBee$lastColonyId,
+ id = id,
queen = queen,
location = c(0, 0),
virginQueens = virginQueens
@@ -107,7 +112,7 @@ createColony <- function(x = NULL, simParamBee = NULL) {
#' # Create and cross Colony and MultiColony class
#' colony <- createColony(x = basePop[2])
#' colony <- cross(colony, drones = droneGroups[[1]])
-#' apiary <- createMultiColony(basePop[3:4], n = 2)
+#' apiary <- createMultiColony(basePop[3:4])
#' apiary <- cross(apiary, drones = droneGroups[2:3])
#'
#' # Check queen and virgin queens IDs
@@ -158,16 +163,20 @@ reQueen <- function(x, queen, removeVirginQueens = TRUE, simParamBee = NULL) {
x <- removeVirginQueens(x, simParamBee = simParamBee)
}
} else {
- x <- removeQueen(x, addVirginQueens = FALSE, simParamBee = simParamBee)
+ x <- removeQueen(x, simParamBee = simParamBee)
x@virginQueens <- queen
}
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
if (nInd(queen) < nCol) {
stop("Not enough queens provided!")
}
- for (colony in seq_len(nCol)) {
- x[[colony]] <- reQueen(
+ x@colonies = foreach(colony = seq_len(nCol)) %dopar% {
+ reQueen(
x = x[[colony]],
queen = queen[colony],
simParamBee = simParamBee
@@ -180,6 +189,33 @@ reQueen <- function(x, queen, removeVirginQueens = TRUE, simParamBee = NULL) {
return(x)
}
+#' @rdname addCastePop_internal
+#' @title An internal function to add a population in a caste slot of the colony
+#'
+#' @description Helper function that returns a colony to allow parallelisation,
+#' only for internal use.
+#'
+#' @param colony \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
+#' @param pop \code{\link[AlphaSimR]{Pop-class}} with one or many individual
+#' @param caste character
+#' @param new logical
+#'
+#' @return \code{\link[SIMplyBee]{Colony-class}}
+#' @export
+addCastePop_internal <- function(pop, colony, caste, new = FALSE) {
+ if (!is.null(pop)) {
+ if (caste == "queen" & nInd(pop) > 1) {
+ stop("Cannot add more than one queen!")
+ }
+ }
+ if (is.null(slot(colony, caste)) | new) {
+ slot(colony, caste) <- pop
+ } else {
+ slot(colony, caste) <- c(slot(colony, caste), pop)
+ }
+ return(colony)
+}
+
#' @rdname addCastePop
#' @title Add caste individuals to the colony
#'
@@ -197,10 +233,6 @@ reQueen <- function(x, queen, removeVirginQueens = TRUE, simParamBee = NULL) {
#' a single value is provided, the same value will be used for all the colonies.
#' @param new logical, should the number of individuals be added to the caste population
#' anew or should we only top-up the existing number of individuals to \code{nInd}
-#' @param exact logical, only relevant when adding workers - if the csd locus is turned
-#' on and exact is \code{TRUE}, we add the exact specified number of viable workers
-#' (heterozygous at the csd locus)
-#' @param year numeric, only relevant when adding virgin queens - year of birth for virgin queens
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#' @param ... additional arguments passed to \code{nInd} when this argument is a function
#'
@@ -258,7 +290,7 @@ reQueen <- function(x, queen, removeVirginQueens = TRUE, simParamBee = NULL) {
#'
#' @export
addCastePop <- function(x, caste = NULL, nInd = NULL, new = FALSE,
- exact = FALSE, year = NULL, simParamBee = NULL, ...) {
+ simParamBee = NULL, ...) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
@@ -293,16 +325,22 @@ addCastePop <- function(x, caste = NULL, nInd = NULL, new = FALSE,
stop("nInd must be non-negative or NULL!")
}
}
+ if (length(nInd) > 1) {
+ warning("More than one value in the nInd argument, taking only the first value!")
+ nInd <- nInd[1]
+ }
if (0 < nInd) {
newInds <- createCastePop(x, nInd,
- caste = caste, exact = exact,
- year = year, simParamBee = simParamBee
+ caste = caste,
+ simParamBee = simParamBee
)
if (caste == "workers") {
homInds <- newInds$nHomBrood
newInds <- newInds$workers
x@queen@misc$nWorkers[[1]] <- x@queen@misc$nWorkers[[1]] + nInd(newInds)
- x@queen@misc$nHomBrood[[1]] <- x@queen@misc$nHomBrood[[1]] + homInds
+ if (isCsdActive(simParamBee = simParamBee)) {
+ x@queen@misc$nHomBrood[[1]] <- x@queen@misc$nHomBrood[[1]] + homInds
+ }
}
if (caste == "drones") {
x@queen@misc$nDrones[[1]] <- x@queen@misc$nDrones[[1]] + nInd(newInds)
@@ -317,26 +355,46 @@ addCastePop <- function(x, caste = NULL, nInd = NULL, new = FALSE,
}
} else if (isMultiColony(x)) {
nCol <- nColonies(x)
- nNInd <- length(nInd)
- if (nNInd > 1 && nNInd < nCol) {
- stop("Too few values in the nInd argument!")
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
}
- if (nNInd > 1 && nNInd > nCol) {
- warning(paste0("Too many values in the nInd argument, taking only the first ", nCol, "values!"))
- nInd <- nInd[1:nCol]
+ if (any(hasCollapsed(x))) {
+ stop(paste0("The colony ", getId(x), " collapsed, hence you can not add individuals (from the queen) to it!"))
}
- for (colony in seq_len(nCol)) {
- if (is.null(nInd)) {
- nIndColony <- NULL
+
+ newInds <- createCastePop(x, nInd,
+ caste = caste,
+ simParamBee = simParamBee,
+ returnSP = FALSE, ...)
+
+
+ if (caste == "workers") {
+ homInds <- lapply(newInds, function(x) {
+ if (is.null(x)) return(NULL)
+ x[['nHomBrood']]
+ })
+ newInds <- lapply(newInds, function(x) {
+ if (is.null(x)) return(NULL)
+ x[["workers"]]
+ })
+ }
+ nInds <- lapply(newInds, function(x) {
+ if (is.null(x)) return(NULL)
+ nInd(x)
+ })
+
+ x@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
+ if (!is.null(nInds[[colony]])) {
+ if (caste == "workers") {
+ x[[colony]]@queen@misc$nWorkers[[1]] <- x[[colony]]@queen@misc$nWorkers[[1]] + nInds[[colony]]
+ x[[colony]]@queen@misc$nHomBrood[[1]] <- x[[colony]]@queen@misc$nHomBrood[[1]] + ifelse(is.null(homInds[[colony]]), 0, homInds[[colony]])
+ } else if (caste == "drones") {
+ x[[colony]]@queen@misc$nDrones[[1]] <- x[[colony]]@queen@misc$nDrones[[1]] + nInds[[colony]]
+ }
+ addCastePop_internal(colony = x[[colony]], pop = newInds[[colony]], caste = caste, new = new)
} else {
- nIndColony <- ifelse(nNInd == 1, nInd, nInd[colony])
+ x[[colony]]
}
- x[[colony]] <- addCastePop(
- x = x[[colony]], caste = caste,
- nInd = nIndColony,
- new = new,
- exact = exact, simParamBee = simParamBee, ...
- )
}
} else {
stop("Argument x must be a Colony or MultiColony class object!")
@@ -348,17 +406,18 @@ addCastePop <- function(x, caste = NULL, nInd = NULL, new = FALSE,
#' @describeIn addCastePop Add workers to a colony
#' @export
addWorkers <- function(x, nInd = NULL, new = FALSE,
- exact = FALSE, simParamBee = NULL, ...) {
+ simParamBee = NULL, ...) {
ret <- addCastePop(
x = x, caste = "workers", nInd = nInd, new = new,
- exact = exact, simParamBee = simParamBee, ...
+ simParamBee = simParamBee, ...
)
return(ret)
}
#' @describeIn addCastePop Add drones to a colony
#' @export
-addDrones <- function(x, nInd = NULL, new = FALSE, simParamBee = NULL, ...) {
+addDrones <- function(x, nInd = NULL, new = FALSE,
+ simParamBee = NULL, ...) {
ret <- addCastePop(
x = x, caste = "drones", nInd = nInd, new = new,
simParamBee = simParamBee, ...
@@ -369,10 +428,10 @@ addDrones <- function(x, nInd = NULL, new = FALSE, simParamBee = NULL, ...) {
#' @describeIn addCastePop Add virgin queens to a colony
#' @export
addVirginQueens <- function(x, nInd = NULL, new = FALSE,
- year = NULL, simParamBee = NULL, ...) {
+ simParamBee = NULL, ...) {
ret <- addCastePop(
x = x, caste = "virginQueens", nInd = nInd, new = new,
- year = year, simParamBee = simParamBee, ...
+ simParamBee = simParamBee, ...
)
return(ret)
}
@@ -398,9 +457,6 @@ addVirginQueens <- function(x, nInd = NULL, new = FALSE,
#' @param new logical, should the number of workers and drones be added anew or
#' should we only top-up the existing number of workers and drones to
#' \code{nWorkers} and \code{nDrones} (see details)
-#' @param exact logical, if the csd locus is turned on and exact is \code{TRUE},
-#' create the exact specified number of only viable workers (heterozygous on
-#' the csd locus)
#' @param resetEvents logical, call \code{\link[SIMplyBee]{resetEvents}} as part of the
#' build up
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
@@ -475,8 +531,8 @@ addVirginQueens <- function(x, nInd = NULL, new = FALSE,
#' getMisc(getQueen(buildUp(colony)))
#' @export
buildUp <- function(x, nWorkers = NULL, nDrones = NULL,
- new = TRUE, exact = FALSE, resetEvents = FALSE,
- simParamBee = NULL, ...) {
+ new = TRUE, resetEvents = FALSE,
+ simParamBee = NULL, ...) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
@@ -512,7 +568,7 @@ buildUp <- function(x, nWorkers = NULL, nDrones = NULL,
if (0 < n) {
x <- addWorkers(
x = x, nInd = n, new = new,
- exact = exact, simParamBee = simParamBee)
+ simParamBee = simParamBee)
} else if (n < 0) {
x@workers <- getWorkers(x, nInd = nWorkers, simParamBee = simParamBee)
}
@@ -543,7 +599,15 @@ buildUp <- function(x, nWorkers = NULL, nDrones = NULL,
}
x@production <- TRUE
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
+
+ if (any(hasCollapsed(x))) {
+ stop(paste0("Some colonies are collapsed, hence you can not build it up!"))
+ }
nCol <- nColonies(x)
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
nNWorkers <- length(nWorkers)
nNDrones <- length(nDrones)
if (nNWorkers > 1 && nNWorkers < nCol) {
@@ -560,27 +624,32 @@ buildUp <- function(x, nWorkers = NULL, nDrones = NULL,
warning(paste0("Too many values in the nDrones argument, taking only the first ", nCol, "values!"))
nNDrones <- nNDrones[1:nCol]
}
- for (colony in seq_len(nCol)) {
- if (is.null(nWorkers)) {
- nWorkersColony <- NULL
- } else {
- nWorkersColony <- ifelse(nNWorkers == 1, nWorkers, nWorkers[colony])
- }
- if (is.null(nDrones)) {
- nDronesColony <- NULL
- } else {
- nDronesColony <- ifelse(nNDrones == 1, nDrones, nDrones[colony])
- }
- x[[colony]] <- buildUp(
- x = x[[colony]],
- nWorkers = nWorkersColony,
- nDrones = nDronesColony,
- new = new,
- exact = exact,
- resetEvents = resetEvents,
- simParamBee = simParamBee, ...
- )
+
+ if (is.function(nWorkers)) {
+ nWorkers <- nWorkers(colony = x,...)
}
+
+ if (new) {
+ n <- nWorkers
+ } else {
+ n <- nWorkers - nWorkers(x, simParamBee = simParamBee)
+ }
+
+ if (sum(nWorkers) > 0) {
+ x <- addWorkers(
+ x = x, nInd = n, new = new,
+ simParamBee = simParamBee)
+ }
+ if (sum(nDrones) > 0) {
+ x <- addDrones(
+ x = x, nInd = n, new = new,
+ simParamBee = simParamBee)
+ }
+ x <- setEvents(x, slot = "production", value = TRUE)
+ if (resetEvents) {
+ x <- resetEvents(x)
+ }
+
} else {
stop("Argument x must be a Colony or MultiColony class object!")
}
@@ -589,7 +658,6 @@ buildUp <- function(x, nWorkers = NULL, nDrones = NULL,
return(x)
}
-
#' @rdname downsize
#' @title Reduce number of workers and remove all drones and virgin queens from
#' a Colony or MultiColony object
@@ -642,6 +710,7 @@ buildUp <- function(x, nWorkers = NULL, nDrones = NULL,
#' apiary <- downsize(x = apiary, p = c(0.5, 0.1), new = TRUE, use = "rand")
#' nWorkers(apiary); nDrones(apiary)
#' @export
+#'
downsize <- function(x, p = NULL, use = "rand", new = FALSE,
simParamBee = NULL, ...) {
if (is.null(simParamBee)) {
@@ -679,8 +748,22 @@ downsize <- function(x, p = NULL, use = "rand", new = FALSE,
x <- removeVirginQueens(x = x, p = 1, simParamBee = simParamBee)
x@production <- FALSE
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
nP <- length(p)
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
+
+ if (any(hasCollapsed(x))) {
+ stop("Some of hte colonies have collapsed, hence you can not downsize them!")
+ }
+ if (is.null(p)) {
+ p <- simParamBee$downsizeP
+ }
+ if (is.function(p)) {
+ p <- p(x, ...)
+ }
if (nP > 1 && nP < nCol) {
stop("Too few values in the p argument!")
}
@@ -688,20 +771,20 @@ downsize <- function(x, p = NULL, use = "rand", new = FALSE,
warning(paste0("Too many values in the p argument, taking only the first ", nCol, "values!"))
p <- p[1:nCol]
}
- for (colony in seq_len(nCol)) {
- if (is.null(p)) {
- pColony <- NULL
- } else {
- pColony <- ifelse(nP == 1, p, p[colony])
- }
- x[[colony]] <- downsize(
- x = x[[colony]],
- p = pColony,
- use = use,
- new = new,
- simParamBee = simParamBee, ...
- )
+ if (new == TRUE) {
+ n <- round(nWorkers(x, simParamBee = simParamBee) * (1 - p))
+ x <- addWorkers(x = x, nInd = n, new = TRUE,
+ simParamBee = simParamBee)
+ } else {
+ x <- removeWorkers(x = x, p = p, use = use,
+ simParamBee = simParamBee)
+ }
+ x <- removeDrones(x = x, p = 1, simParamBee = simParamBee)
+ x <- removeVirginQueens(x = x, p = 1, simParamBee = simParamBee)
+ for (colony in 1:nCol) {
+ x[[colony]]@production <- FALSE
}
+
} else {
stop("Argument x must be a Colony or MultiColony class object!")
}
@@ -710,6 +793,8 @@ downsize <- function(x, p = NULL, use = "rand", new = FALSE,
return(x)
}
+
+
#' @rdname replaceCastePop
#' @title Replace a proportion of caste individuals with new ones
#'
@@ -726,12 +811,6 @@ downsize <- function(x, p = NULL, use = "rand", new = FALSE,
#' a single value is provided, the same value will be applied to all the colonies
#' @param use character, all the options provided by \code{\link[AlphaSimR]{selectInd}} -
#' guides selection of caste individuals that stay when \code{p < 1}
-#' @param exact logical, only relevant when adding workers - if the csd locus is turned
-#' on and exact is \code{TRUE}, we replace the exact specified number of viable workers
-#' (heterozygous at the csd locus). You probably want this set to TRUE since you want to
-#' replace with the same number of workers.
-#' @param year numeric, only relevant when replacing virgin queens,
-#' year of birth for virgin queens
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#'
#' @return \code{\link[SIMplyBee]{Colony-class}} or or \code{\link[SIMplyBee]{MultiColony-class}} with
@@ -768,8 +847,8 @@ downsize <- function(x, p = NULL, use = "rand", new = FALSE,
#' apiary <- replaceWorkers(apiary, p = 0.5)
#' getCasteId(apiary, caste="workers")
#' @export
-replaceCastePop <- function(x, caste = NULL, p = 1, use = "rand", exact = TRUE,
- year = NULL, simParamBee = NULL) {
+replaceCastePop <- function(x, caste = NULL, p = 1, use = "rand",
+ simParamBee = NULL) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
@@ -781,71 +860,49 @@ replaceCastePop <- function(x, caste = NULL, p = 1, use = "rand", exact = TRUE,
} else if (any(p < 0)) {
stop("p must not be less than 0!")
}
- if (isColony(x)) {
- if (hasCollapsed(x)) {
- stop(paste0("The colony ", getId(x), " collapsed, hence you can not replace individuals in it!"))
+ if (isColony(x) | isMultiColony(x)) {
+ nP <- length(p)
+ if (isColony(x)) {
+ nCol <- 1
+ } else if (isMultiColony(x)) {
+ nCol <- nColonies(x)
}
- if (!isQueenPresent(x, simParamBee = simParamBee)) {
- stop("Missing queen!")
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
}
- if (length(p) > 1) {
- warning("More than one value in the p argument, taking only the first value!")
- p <- p[1]
+ if (any(hasCollapsed(x))) {
+ stop(paste0("The colony or some of the colonies have collapsed, hence you can not replace individuals in it!"))
}
- nInd <- nCaste(x, caste, simParamBee = simParamBee)
- if (nInd > 0) {
- nIndReplaced <- round(nInd * p)
- if (nIndReplaced < nInd) {
- nIndStay <- nInd - nIndReplaced
- if (nIndReplaced > 0) {
- tmp <- createCastePop(x,
- caste = caste,
- nInd = nIndReplaced, exact = exact,
- year = year, simParamBee = simParamBee
- )
- if (caste == "workers") {
- x@queen@misc$nWorkers[[1]] <- x@queen@misc$nWorkers[[1]] + nIndReplaced
- x@queen@misc$nHomBrood[[1]] <- x@queen@misc$nHomBrood[[1]] + tmp$nHomBrood
- tmp <- tmp$workers
- }
- if (caste == "drones") {
- x@queen@misc$nDrones[[1]] <- x@queen@misc$nDrones[[1]] + nIndReplaced
- }
-
- slot(x, caste) <- c(
- selectInd(slot(x, caste), nInd = nIndStay, use = use, simParam = simParamBee),
- tmp
- )
- }
- } else {
- x <- addCastePop(
- x = x, caste = caste, nInd = nIndReplaced, new = TRUE,
- year = year, simParamBee = simParamBee
- )
- }
+ if (any(!isQueenPresent(x, simParamBee = simParamBee))) {
+ stop("Missing queen in at least one colony!")
}
- } else if (isMultiColony(x)) {
- nCol <- nColonies(x)
- nP <- length(p)
if (nP > 1 && nP < nCol) {
stop("Too few values in the p argument!")
}
- if (nP > 1 && nP > nCol) {
- warning(paste0("Too many values in the p argument, taking only the first ", nCol, "values!"))
- p <- p[1:nCol]
+ if (length(p) > nCol) {
+ warning(paste0("More than one value in the p argument, taking only the first ", nCol, " values!"))
+ p <- p[nCol]
}
- for (colony in seq_len(nCol)) {
- if (is.null(p)) {
- pColony <- NULL
+ nInd <- nCaste(x, caste, simParamBee = simParamBee)
+ if (any(nInd > 0)) {
+ nIndReplaced <- round(nInd * p)
+ if (any(nIndReplaced < nInd)) {
+
+ x <- removeCastePop(x,
+ caste = caste,
+ p = p)
+ nIndAdd <- nInd - nCaste(x, caste, simParamBee = simParamBee)
+ x <- addCastePop(x,
+ caste = caste,
+ nInd = nIndAdd,
+ simParamBee = simParamBee
+ )
} else {
- pColony <- ifelse(nP == 1, p, p[colony])
+ x <- addCastePop(
+ x = x, caste = caste, nInd = nIndReplaced, new = TRUE,
+ simParamBee = simParamBee
+ )
}
- x[[colony]] <- replaceCastePop(
- x = x[[colony]], caste = caste,
- p = pColony,
- use = use, year = year,
- simParamBee = simParamBee
- )
}
} else {
stop("Argument x must be a Colony or MultiColony class object!")
@@ -856,10 +913,10 @@ replaceCastePop <- function(x, caste = NULL, p = 1, use = "rand", exact = TRUE,
#' @describeIn replaceCastePop Replaces some workers in a colony
#' @export
-replaceWorkers <- function(x, p = 1, use = "rand", exact = TRUE, simParamBee = NULL) {
+replaceWorkers <- function(x, p = 1, use = "rand", simParamBee = NULL) {
ret <- replaceCastePop(
x = x, caste = "workers", p = p,
- use = use, exact = exact,
+ use = use,
simParamBee = simParamBee
)
return(ret)
@@ -885,6 +942,7 @@ replaceVirginQueens <- function(x, p = 1, use = "rand", simParamBee = NULL) {
return(ret)
}
+
#' @rdname removeCastePop
#' @title Remove a proportion of caste individuals from a colony
#'
@@ -898,13 +956,7 @@ replaceVirginQueens <- function(x, p = 1, use = "rand", simParamBee = NULL) {
#' a single value is provided, the same value will be applied to all the colonies
#' @param use character, all the options provided by \code{\link[AlphaSimR]{selectInd}} -
#' guides selection of virgins queens that will stay when \code{p < 1}
-#' @param addVirginQueens logical, whether virgin queens should be added; only
-#' used when removing the queen from the colony
-#' @param nVirginQueens integer, the number of virgin queens to be created in the
-#' colony; only used when removing the queen from the colony. If \code{0}, no virgin
-#' queens are added; If \code{NULL}, the value from \code{simParamBee$nVirginQueens}
-#' is used
-#' @param year numeric, only relevant when adding virgin queens - year of birth for virgin queens
+
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#'
#' @return \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} without virgin queens
@@ -943,8 +995,7 @@ replaceVirginQueens <- function(x, p = 1, use = "rand", simParamBee = NULL) {
#' nWorkers(removeWorkers(apiary, p = c(0.1, 0.5)))
#' @export
removeCastePop <- function(x, caste = NULL, p = 1, use = "rand",
- addVirginQueens = FALSE, nVirginQueens = NULL,
- year = NULL, simParamBee = NULL) {
+ simParamBee = NULL) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
@@ -964,14 +1015,6 @@ removeCastePop <- function(x, caste = NULL, p = 1, use = "rand",
warning("More than one value in the p argument, taking only the first value!")
p <- p[1]
}
- if (caste == "queen") {
- if (addVirginQueens) {
- if (is.null(nVirginQueens)) {
- nVirginQueens <- simParamBee$nVirginQueens
- }
- x <- addVirginQueens(x, nInd = nVirginQueens, year = year, simParamBee = simParamBee)
- }
- }
if (p == 1) {
slot(x, caste) <- NULL
} else {
@@ -988,8 +1031,12 @@ removeCastePop <- function(x, caste = NULL, p = 1, use = "rand",
}
}
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
nP <- length(p)
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
if (nP > 1 && nP < nCol) {
stop("Too few values in the p argument!")
}
@@ -997,13 +1044,13 @@ removeCastePop <- function(x, caste = NULL, p = 1, use = "rand",
warning(paste0("Too many values in the p argument, taking only the first ", nCol, "values!"))
p <- p[1:nCol]
}
- for (colony in seq_len(nCol)) {
+ x@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
if (is.null(p)) {
pColony <- NULL
} else {
pColony <- ifelse(nP == 1, p, p[colony])
}
- x[[colony]] <- removeCastePop(
+ removeCastePop(
x = x[[colony]], caste = caste,
p = pColony,
use = use,
@@ -1019,12 +1066,13 @@ removeCastePop <- function(x, caste = NULL, p = 1, use = "rand",
#' @describeIn removeCastePop Remove queen from a colony
#' @export
-removeQueen <- function(x, addVirginQueens = FALSE, nVirginQueens = NULL, year = NULL, simParamBee = NULL) {
- ret <- removeCastePop(x = x, caste = "queen", p = 1, addVirginQueens = addVirginQueens,
- nVirginQueens = nVirginQueens, year = year, simParamBee = simParamBee)
+#'
+removeQueen <- function(x, simParamBee = NULL) {
+ ret <- removeCastePop(x = x, caste = "queen", p = 1, simParamBee = simParamBee)
return(ret)
}
+
#' @describeIn removeCastePop Remove workers from a colony
#' @export
removeWorkers <- function(x, p = 1, use = "rand", simParamBee = NULL) {
@@ -1032,6 +1080,8 @@ removeWorkers <- function(x, p = 1, use = "rand", simParamBee = NULL) {
return(ret)
}
+
+
#' @describeIn removeCastePop Remove workers from a colony
#' @export
removeDrones <- function(x, p = 1, use = "rand", simParamBee = NULL) {
@@ -1039,6 +1089,8 @@ removeDrones <- function(x, p = 1, use = "rand", simParamBee = NULL) {
return(ret)
}
+
+
#' @describeIn removeCastePop Remove virgin queens from a colony
#' @export
removeVirginQueens <- function(x, p = 1, use = "rand", simParamBee = NULL) {
@@ -1059,6 +1111,7 @@ removeVirginQueens <- function(x, p = 1, use = "rand", simParamBee = NULL) {
#' up a new colony, which the default of \code{NULL} caters for; otherwise, a
#' collapsed colony should be left collapsed forever, unless you force
#' resetting this event with \code{collapse = TRUE})
+#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#'
#' @return \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with
#' events reset
@@ -1125,7 +1178,10 @@ removeVirginQueens <- function(x, p = 1, use = "rand", simParamBee = NULL) {
#' hasSplit(remnants[[1]])
#' resetEvents(remnants)[[1]]
#' @export
-resetEvents <- function(x, collapse = NULL) {
+resetEvents <- function(x, collapse = NULL, simParamBee = NULL) {
+ if (is.null(simParamBee)) {
+ simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ }
if (isColony(x)) {
x@swarm <- FALSE
x@split <- FALSE
@@ -1141,11 +1197,16 @@ resetEvents <- function(x, collapse = NULL) {
x@production <- FALSE
validObject(x)
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
- for (colony in seq_len(nCol)) {
- x[[colony]] <- resetEvents(
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
+ x@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
+ resetEvents(
x = x[[colony]],
- collapse = collapse
+ collapse = collapse,
+ simParamBee = simParamBee
)
}
validObject(x)
@@ -1166,6 +1227,8 @@ resetEvents <- function(x, collapse = NULL) {
#'
#' @return \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with the collapse
#' event set to \code{TRUE}
+#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
+#'
#'
#' @details You should use this function in an edge-case when you
#' want to indicate that the colony has collapsed, but you still want to
@@ -1178,13 +1241,13 @@ resetEvents <- function(x, collapse = NULL) {
#' SP <- SimParamBee$new(founderGenomes)
#' \dontshow{SP$nThreads = 1L}
#' basePop <- createVirginQueens(founderGenomes)
-#' drones <- createDrones(basePop[1], n = 1000)
+#' drones <- createDrones(basePop[1], nInd = 1000)
#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = 10)
#'
#' # Create Colony and MultiColony class
#' colony <- createColony(x = basePop[1])
#' colony <- cross(colony, drones = droneGroups[[1]])
-#' apiary <- createMultiColony(x = basePop[2:10], n = 9)
+#' apiary <- createMultiColony(x = basePop[2:10])
#' apiary <- cross(apiary, drones = droneGroups[2:10])
#'
#' # Collapse
@@ -1200,14 +1263,22 @@ resetEvents <- function(x, collapse = NULL) {
#' apiaryLeft <- tmp$remnant
#' hasCollapsed(apiaryLeft)
#' @export
-collapse <- function(x) {
+collapse <- function(x, simParamBee = NULL) {
+ if (is.null(simParamBee)) {
+ simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ }
if (isColony(x)) {
x@collapse <- TRUE
x@production <- FALSE
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
- for (colony in seq_len(nCol)) {
- x[[colony]] <- collapse(x = x[[colony]])
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
+ x@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
+ collapse(x = x[[colony]],
+ simParamBee = simParamBee)
}
} else {
stop("Argument x must be a Colony or MultiColony class object!")
@@ -1233,11 +1304,6 @@ collapse <- function(x) {
#' If input is \code{\link[SIMplyBee]{MultiColony-class}},
#' the input could also be a vector of the same length as the number of colonies. If
#' a single value is provided, the same value will be applied to all the colonies
-#' @param year numeric, year of birth for virgin queens
-#' @param nVirginQueens integer, the number of virgin queens to be created in the
-#' colony; of these one is randomly selected as the new virgin queen of the
-#' remnant colony. If \code{NULL}, the value from \code{simParamBee$nVirginQueens}
-#' is used
#' @param sampleLocation logical, sample location of the swarm by taking
#' the current colony location and adding deviates to each coordinate using
#' \code{\link[SIMplyBee]{rcircle}}
@@ -1255,6 +1321,8 @@ collapse <- function(x) {
#' @examples
#' founderGenomes <- quickHaplo(nInd = 8, nChr = 1, segSites = 50)
#' SP <- SimParamBee$new(founderGenomes)
+#' SP$setTrackPed(TRUE)
+#' SP$setTrackRec(TRUE)
#' \dontshow{SP$nThreads = 1L}
#' basePop <- createVirginQueens(founderGenomes)
#' drones <- createDrones(basePop[1], n = 1000)
@@ -1264,7 +1332,7 @@ collapse <- function(x) {
#' colony <- createColony(x = basePop[2])
#' colony <- cross(colony, drones = droneGroups[[1]])
#' (colony <- buildUp(colony, nWorkers = 100))
-#' apiary <- createMultiColony(basePop[3:8], n = 6)
+#' apiary <- createMultiColony(basePop[3:8])
#' apiary <- cross(apiary, drones = droneGroups[2:7])
#' apiary <- buildUp(apiary, nWorkers = 100)
#'
@@ -1288,136 +1356,158 @@ collapse <- function(x) {
#' # Swarm only the pulled colonies
#' (swarm(tmp$pulled, p = 0.6))
#' @export
-swarm <- function(x, p = NULL, year = NULL, nVirginQueens = NULL,
+swarm <- function(x, p = NULL,
sampleLocation = TRUE, radius = NULL,
simParamBee = NULL, ...) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
+ if (isMultiColony(x)) {
+ parallel <- TRUE
+ }
if (is.null(p)) {
p <- simParamBee$swarmP
}
if (is.null(radius)) {
radius <- simParamBee$swarmRadius
}
- if (is.null(nVirginQueens)) {
- nVirginQueens <- simParamBee$nVirginQueens
- }
- if (isColony(x)) {
- if (hasCollapsed(x)) {
- stop(paste0("The colony ", getId(x), " collapsed, hence it can not swarm!"))
+ if (isColony(x) | isMultiColony(x)) {
+ if (isColony(x)) {
+ nCol <- 1
+ } else if (isMultiColony(x)) {
+ nCol <- nColonies(x)
}
- if (!isQueenPresent(x, simParamBee = simParamBee)) {
- stop("No queen present in the colony!")
+ nP <- length(p)
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
}
- if (!isWorkersPresent(x, simParamBee = simParamBee)) {
- stop("No workers present in the colony!")
+
+ if (any(hasCollapsed(x))) {
+ stop(paste0("One of the collonies is collapsed, hence you can not split it!"))
+ }
+ if (any(!isQueenPresent(x, simParamBee = simParamBee))) {
+ stop("No queen present in one of the colonies!")
+ }
+ if (any(!isWorkersPresent(x, simParamBee = simParamBee))) {
+ stop("No workers present in one of the colonies!")
}
if (is.function(p)) {
p <- p(x, ...)
} else {
- if (p < 0 | 1 < p) {
+ if (any(p < 0) | any(1 < p)) {
stop("p must be between 0 and 1 (inclusive)!")
}
- if (length(p) > 1) {
+ if (length(p) > nCol) {
warning("More than one value in the p argument, taking only the first value!")
- p <- p[1]
+ p <- p[nCol]
+ }
+ if (nP > 1 && nP < nCol) {
+ stop("Too few values in the p argument!")
}
- }
- if (is.function(nVirginQueens)) {
- nVirginQueens <- nVirginQueens(x, ...)
}
nWorkers <- nWorkers(x, simParamBee = simParamBee)
nWorkersSwarm <- round(nWorkers * p)
# TODO: Add use="something" to select pWorkers that swarm
# https://github.com/HighlanderLab/SIMplyBee/issues/160
- tmp <- pullWorkers(x = x, nInd = nWorkersSwarm, simParamBee = simParamBee)
+
+ tmpVirginQueens <- createCastePop(
+ x = x, nInd = max(10, simParamBee$nVirginQueens),
+ caste = "virginQueens",
+ simParamBee = simParamBee
+ )
+
+ if (isColony(x)) {
+ homCol = nInd(tmpVirginQueens) == 0
+ } else if (isMultiColony(x)) {
+ homCol = lapply(tmpVirginQueens, nInd) == 0
+ }
+
+ if (sum(homCol) > 0) {
+ if (isColony(x)) {
+ stop("Colony too inbred to produce any virgin queens!")
+ } else if (isMultiColony(x)) {
+ warning(paste0(sum(homCol), " colonies produced 0 virgin queens due to high colony homozygosity, removing these colonies!"))
+ tmpVirginQueens <- tmpVirginQueens[!homCol]
+ x = x[!homCol]
+ location = location[!homCol]
+ nWorkersSwarm = nWorkersSwarm[!homCol]
+ nCol = nColonies(x)
+ }
+ }
+
+ tmp <- pullCastePop(x = x, caste = "workers",
+ nInd = nWorkersSwarm, simParamBee = simParamBee)
+ remnantColony <- tmp$remnant
+ remnantColony <- removeQueen(remnantColony)
+ if (isColony(x)) {
+ remnantColony <- reQueen(remnantColony,
+ queen = selectInd(tmpVirginQueens, nInd = 1, use = "rand"),
+ simParamBee = simParamBee)
+ } else {
+ tmpVirginQueens <- lapply(tmpVirginQueens, FUN = function(x) selectInd(x, nInd = 1, use = "rand"))
+ remnantColony <- reQueen(remnantColony,
+ queen = mergePops(tmpVirginQueens),
+ simParamBee = simParamBee)
+ }
currentLocation <- getLocation(x)
+
if (sampleLocation) {
- newLocation <- c(currentLocation + rcircle(radius = radius))
+ newLocation <- lapply(1:nCol, function(x) currentLocation[[x]] + rcircle(n = nCol, radius = radius)[x,])
# c() to convert row-matrix to a numeric vector
} else {
newLocation <- currentLocation
}
- swarmColony <- createColony(x = x@queen, simParamBee = simParamBee)
- # It's not re-queening, but the function also sets the colony id
- swarmColony@workers <- tmp$pulled
- swarmColony <- setLocation(x = swarmColony, location = newLocation)
+ if (isColony(x)) {
+ swarmColony <- createColony(x = x@queen, simParamBee = simParamBee)
+ # It's not re-queening, but the function also sets the colony id
- tmpVirginQueen <- createVirginQueens(
- x = x, nInd = nVirginQueens,
- year = year,
- simParamBee = simParamBee
- )
- tmpVirginQueen <- selectInd(tmpVirginQueen, nInd = 1, use = "rand", simParam = simParamBee)
+ swarmColony@workers <- tmp$pulled
+ swarmColony <- setLocation(x = swarmColony, location = newLocation[[1]])
- remnantColony <- createColony(x = tmpVirginQueen, simParamBee = simParamBee)
- remnantColony@workers <- getWorkers(tmp$remnant, simParamBee = simParamBee)
- remnantColony@drones <- getDrones(x, simParamBee = simParamBee)
- # Workers raise virgin queens from eggs laid by the queen and one random
- # virgin queen prevails, so we create just one
- # Could consider that a non-random one prevails (say the more aggressive one),
- # by creating many virgin queens and then picking the one with highest
- # gv/pheno for competition or some other criteria (patri-lineage)
+ remnantColony <- setLocation(x = remnantColony, location = currentLocation)
- remnantColony <- setLocation(x = remnantColony, location = currentLocation)
+ remnantColony@swarm <- TRUE
+ swarmColony@swarm <- TRUE
- remnantColony@swarm <- TRUE
- swarmColony@swarm <- TRUE
- remnantColony@production <- FALSE
- swarmColony@production <- FALSE
+ remnantColony@production <- FALSE
+ swarmColony@production <- FALSE
+
+ ret <- list(swarm = swarmColony, remnant = remnantColony)
+ } else if (isMultiColony(x)) {
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
- ret <- list(swarm = swarmColony, remnant = remnantColony)
- } else if (isMultiColony(x)) {
- nCol <- nColonies(x)
- nP <- length(p)
- if (nP > 1 && nP < nCol) {
- stop("Too few values in the p argument!")
- }
- if (nP > 1 && nP > nCol) {
- warning(paste0("Too many values in the p argument, taking only the first ", nCol, "values!"))
- p <- p[1:nCol]
- }
- if (nCol == 0) {
- ret <- list(
- swarm = createMultiColony(simParamBee = simParamBee),
- remnant = createMultiColony(simParamBee = simParamBee)
- )
- } else {
ret <- list(
- swarm = createMultiColony(n = nCol, simParamBee = simParamBee),
- remnant = createMultiColony(n = nCol, simParamBee = simParamBee)
+ swarm = createMultiColony(x = getQueen(x, collapse = TRUE),
+ simParamBee = simParamBee),
+ remnant = remnantColony
)
- for (colony in seq_len(nCol)) {
- if (is.null(p)) {
- pColony <- NULL
- } else {
- pColony <- ifelse(nP == 1, p, p[colony])
- }
- tmp <- swarm(x[[colony]],
- p = pColony,
- year = year,
- nVirginQueens = nVirginQueens,
- sampleLocation = sampleLocation,
- radius = radius,
- simParamBee = simParamBee, ...
- )
- ret$swarm[[colony]] <- tmp$swarm
- ret$remnant[[colony]] <- tmp$remnant
+
+ ret$swarm@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
+ addCastePop_internal(colony = ret$swarm@colonies[[colony]],
+ pop = tmp$pulled[[colony]], caste = "workers")
}
+
+ ret$remnant <- setEvents(ret$remnant, slot = "swarm", value = TRUE)
+ ret$swarm <- setEvents(ret$swarm, slot = "swarm", value = TRUE)
+ ret$swarm <- setEvents(ret$swarm, slot = "production", value = FALSE)
+ ret$remnant <- setEvents(ret$remnant, slot = "production", value = FALSE)
}
- } else {
+ }
+ else {
stop("Argument x must be a Colony or MultiColony class object!")
}
-
validObject(ret$swarmColony)
validObject(ret$remnantColony)
return(ret)
}
+
+
#' @rdname supersede
#' @title Supersede
#'
@@ -1427,11 +1517,6 @@ swarm <- function(x, p = NULL, year = NULL, nVirginQueens = NULL,
#' queens, of which only one prevails.
#'
#' @param x \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
-#' @param year numeric, year of birth for virgin queens
-#' @param nVirginQueens integer, the number of virgin queens to be created in the
-#' colony; of these one is randomly selected as the new virgin queen of the
-#' remnant colony. If \code{NULL}, the value from \code{simParamBee$nVirginQueens}
-#' is used
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#' @param ... additional arguments passed to \code{nVirginQueens} when this
#' argument is a function
@@ -1474,45 +1559,87 @@ swarm <- function(x, p = NULL, year = NULL, nVirginQueens = NULL,
#' # Swarm only the pulled colonies
#' (supersede(tmp$pulled))
#' @export
-supersede <- function(x, year = NULL, nVirginQueens = NULL, simParamBee = NULL, ...) {
+supersede <- function(x, simParamBee = NULL, ...) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
+ if (isColony(x)) {
+ parallel <- FALSE
+ } else if (isMultiColony(x)) {
+ parallel <- TRUE
+ }
if (is.null(nVirginQueens)) {
nVirginQueens <- simParamBee$nVirginQueens
}
+
+ if (any(hasCollapsed(x))) {
+ stop(paste0("One of the collonies is collapsed, hence you can not split it!"))
+ }
+ if (any(!isQueenPresent(x, simParamBee = simParamBee))) {
+ stop("No queen present in one of the colonies!")
+ }
+ if (is.function(nVirginQueens)) {
+ nVirginQueens <- nVirginQueens(x, ...)
+ }
+
+ # Do this because some colonies might not produce a viable virgin queen
+ tmpVirginQueens <- createCastePop(
+ x = x, nInd = max(10, SP$nVirginQueens),
+ caste = "virginQueens",
+ simParamBee = simParamBee
+ )
+
if (isColony(x)) {
- if (hasCollapsed(x)) {
- stop(paste0("The colony ", getId(x), " collapsed, hence it can not supresede!"))
- }
- if (!isQueenPresent(x, simParamBee = simParamBee)) {
- stop("No queen present in the colony!")
+ homCol = nInd(tmpVirginQueens) == 0
+ } else if (isMultiColony(x)) {
+ homCol = sapply(tmpVirginQueens, nInd) == 0
+ }
+
+ if (sum(homCol) > 0) {
+ if (isColony(x)) {
+ print("X is colony")
+ print(class(x))
+ stop("Colony to inbred to produce any virgin queens!")
+ } else if (isMultiColony(x)) {
+ warning(paste0(sum(homCol), " colonies produced 0 virgin queens due to high colony homozygosity, removing these colonies!"))
+ tmpVirginQueens <- tmpVirginQueens[!homCol]
+ x = x[!homCol]
+ nCol = nColonies(x)
}
- if (is.function(nVirginQueens)) {
- nVirginQueens <- nVirginQueens(x, ...)
+ }
+
+ if (isColony(x)) {
+ if (!parallel) {
+ x <- addCastePop_internal(selectInd(tmpVirginQueens, nInd = 1, use = "rand"), colony = x, caste = "virginQueens")
}
- x <- removeQueen(x, addVirginQueens = TRUE, nVirginQueens = nVirginQueens,
- year = year, simParamBee = simParamBee)
- x@virginQueens <- selectInd(x@virginQueens, nInd = 1, use = "rand", simParam = simParamBee)
+ x <- removeQueen(x, simParamBee = simParamBee)
# TODO: We could consider that a non-random virgin queen prevails (say the most
# aggressive one), by creating many virgin queens and then picking the
# one with highest pheno for competition or some other criteria
# https://github.com/HighlanderLab/SIMplyBee/issues/239
x@supersedure <- TRUE
} else if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
nCol <- nColonies(x)
if (nCol == 0) {
- x <- createMultiColony(simParamBee = simParamBee)
- } else {
- for (colony in seq_len(nCol)) {
- x[[colony]] <- supersede(x[[colony]],
- year = year,
- nVirginQueens = nVirginQueens,
- simParamBee = simParamBee, ...
- )
+ stop("The Multicolony contains 0 colonies!")
+ }
+ tmpVirginQueens <- lapply(tmpVirginQueens, FUN = function(x) selectInd(x, nInd = 1, use = "rand"))
+
+ combine_list <- function(a, b) {
+ if (length(a) == 1) {
+ c(list(a), list(b))
+ } else {
+ c(a, list(b))
}
}
- } else {
+ x@colonies <- foreach(colony = seq_len(nColonies(x))) %dopar% {
+ addCastePop_internal(colony = removeQueen(x[[colony]], simParamBee = simParamBee),
+ pop = tmpVirginQueens[[colony]], caste = "virginQueens")
+ }
+ x = setEvents(x, slot = "supersedure", value = TRUE)
+ }
+ else {
stop("Argument x must be a Colony or MultiColony class object!")
}
validObject(x)
@@ -1526,8 +1653,9 @@ supersede <- function(x, year = NULL, nVirginQueens = NULL, simParamBee = NULL,
#' into two new colonies to
#' prevent swarming (in managed situation). The remnant colony retains the
#' queen and a proportion of the workers and all drones. The split colony gets
-#' the other part of the workers, which raise virgin queens, of which only one
-#' prevails. Location of the split is the same as for the remnant.
+#' the other part of the workers, but note that it is queenless, since the beekeepers
+#' would normally requeen with a different queen.
+#' Location of the split is the same as for the remnant.
#'
#' @param x \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
#' @param p numeric, proportion of workers that will go to the split colony; if
@@ -1535,7 +1663,6 @@ supersede <- function(x, year = NULL, nVirginQueens = NULL, simParamBee = NULL,
#' If input is \code{\link[SIMplyBee]{MultiColony-class}},
#' the input could also be a vector of the same length as the number of colonies. If
#' a single value is provided, the same value will be applied to all the colonies
-#' @param year numeric, year of birth for virgin queens
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#' @param ... additional arguments passed to \code{p} when this argument is a
#' function
@@ -1580,109 +1707,162 @@ supersede <- function(x, year = NULL, nVirginQueens = NULL, simParamBee = NULL,
#' # Split only the pulled colonies
#' (split(tmp$pulled, p = 0.5))
#' @export
-split <- function(x, p = NULL, year = NULL, simParamBee = NULL, ...) {
+split <- function(x, p = NULL, simParamBee = NULL, ...) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
if (is.null(p)) {
p <- simParamBee$splitP
}
- if (isColony(x)) {
- if (hasCollapsed(x)) {
- stop(paste0("The colony ", getId(x), " collapsed, hence you can not split it!"))
+ if (isMultiColony(x)) {
+ parallel <- TRUE
+ }
+
+ if (isColony(x) | isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
+ if (isColony(x)) {
+ nCol <- 1
+ } else if (isMultiColony(x)) {
+ nCol <- nColonies(x)
}
- if (!isQueenPresent(x, simParamBee = simParamBee)) {
- stop("No queen present in the colony!")
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
}
- if (!isWorkersPresent(x, simParamBee = simParamBee)) {
- stop("No workers present in the colony!")
+ nP <- length(p)
+
+ location <- getLocation(x)
+ if (any(hasCollapsed(x))) {
+ stop(paste0("One of the collonies is collapsed, hence you can not split it!"))
+ }
+ if (any(!isQueenPresent(x, simParamBee = simParamBee))) {
+ stop("No queen present in one of the colonies!")
+ }
+ if (any(!isWorkersPresent(x, simParamBee = simParamBee))) {
+ stop("No workers present in one of the colonies!")
}
if (is.function(p)) {
p <- p(x, ...)
} else {
- if (p < 0 | 1 < p) {
+ if (any(p < 0) | any(1 < p)) {
stop("p must be between 0 and 1 (inclusive)!")
}
- if (length(p) > 1) {
+ if (length(p) > nCol) {
warning("More than one value in the p argument, taking only the first value!")
- p <- p[1]
+ p <- p[nCol]
+ }
+ if (nP > 1 && nP < nCol) {
+ stop("Too few values in the p argument!")
}
}
+
nWorkers <- nWorkers(x, simParamBee = simParamBee)
nWorkersSplit <- round(nWorkers * p)
# TODO: Split colony at random by default, but we could make it as a
# function of some parameters
# https://github.com/HighlanderLab/SIMplyBee/issues/179
- tmp <- pullWorkers(x = x, nInd = nWorkersSplit, simParamBee = simParamBee)
+ tmp <- pullCastePop(x = x, caste = "workers", nInd = nWorkersSplit, simParamBee = simParamBee)
remnantColony <- tmp$remnant
- tmpVirginQueens <- createVirginQueens(
- x = x, nInd = 1,
- year = year,
- simParamBee = simParamBee
- )
- # Workers raise virgin queens from eggs laid by the queen (assuming) that
- # a frame of brood is also provided to the split and then one random virgin
- # queen prevails, so we create just one
- # TODO: Could consider that a non-random one prevails (say the most aggressive
- # one), by creating many virgin queens and then picking the one with
- # highest pheno for competition or some other criteria
- # https://github.com/HighlanderLab/SIMplyBee/issues/239
- splitColony <- createColony(x = tmpVirginQueens, simParamBee = simParamBee)
- splitColony@workers <- tmp$pulled
- splitColony <- setLocation(x = splitColony, location = getLocation(splitColony))
+ if (isColony(x)) {
- remnantColony@split <- TRUE
- splitColony@split <- TRUE
+ # Workers raise virgin queens from eggs laid by the queen (assuming) that
+ # a frame of brood is also provided to the split and then one random virgin
+ # queen prevails, so we create just one
+ # TODO: Could consider that a non-random one prevails (say the most aggressive
+ # one), by creating many virgin queens and then picking the one with
+ # highest pheno for competition or some other criteria
+ # https://github.com/HighlanderLab/SIMplyBee/issues/239
- remnantColony@production <- TRUE
- splitColony@production <- FALSE
+ splitColony <- createColony(simParamBee = simParamBee)
+ splitColony <- setLocation(x = splitColony, location = location)
+
+ splitColony@workers <- tmp$pulled
+
+ remnantColony@split <- TRUE
+ splitColony@split <- TRUE
+
+ remnantColony@production <- TRUE
+ splitColony@production <- FALSE
+
+ ret <- list(split = splitColony, remnant = remnantColony)
+ } else if (isMultiColony(x)) {
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
- ret <- list(split = splitColony, remnant = remnantColony)
- } else if (isMultiColony(x)) {
- nCol <- nColonies(x)
- nP <- length(p)
- if (nP > 1 && nP < nCol) {
- stop("Too few values in the p argument!")
- }
- if (nP > 1 && nP > nCol) {
- warning(paste0("Too many values in the nInd argument, taking only the first ", nCol, "values!"))
- p <- p[1:nCol]
- }
- if (nCol == 0) {
- ret <- list(
- split = createMultiColony(simParamBee = simParamBee),
- remnant = createMultiColony(simParamBee = simParamBee)
- )
- } else {
ret <- list(
- split = createMultiColony(n = nCol, simParamBee = simParamBee),
- remnant = createMultiColony(n = nCol, simParamBee = simParamBee)
+ split = createMultiColony(n = nCol,
+ simParamBee = simParamBee,
+ populateColonies = TRUE),
+ remnant = remnantColony
+
)
- for (colony in seq_len(nCol)) {
- if (is.null(p)) {
- pColony <- NULL
- } else {
- pColony <- ifelse(nP == 1, p, p[colony])
- }
- tmp <- split(x[[colony]],
- p = pColony,
- year = year,
- simParamBee = simParamBee, ...
- )
- ret$split[[colony]] <- tmp$split
- ret$remnant[[colony]] <- tmp$remnant
+ ret$split <- setLocation(x = ret$split, location = location)
+ ret$split@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
+ addCastePop_internal(colony = ret$split@colonies[[colony]],
+ pop = tmp$pulled[[colony]], caste = "workers")
}
+
+ ret$split <- setEvents(ret$split, slot = "split", value = TRUE)
+ ret$remnant <- setEvents(ret$remnant, slot = "split", value = TRUE)
+ ret$split <- setEvents(ret$split, slot = "production", value = FALSE)
+ ret$remnant <- setEvents(ret$remnant, slot = "production", value = TRUE)
}
- } else {
+ }
+ else {
stop("Argument x must be a Colony or MultiColony class object!")
}
- validObject(ret$splitColony)
- validObject(ret$remnantColony)
+ validObject(ret$split)
+ validObject(ret$remnant)
return(ret)
}
+#' @rdname setEvents
+#' @title Set colony events
+#'
+#' @description Helper Level 2 function that populates the events slot. Not interded
+#' for external use, intended for internal use in parallel computing
+#'
+#' @param x \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
+#' @param slot character, which event to set
+#' @param value logical, the value for the event
+#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
+#'
+#' @return \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with
+#' events reset
+#'
+#' @examples
+#' founderGenomes <- quickHaplo(nInd = 5, nChr = 1, segSites = 50)
+#' SP <- SimParamBee$new(founderGenomes)
+#' \dontshow{SP$nThreads = 1L}
+#' basePop <- createVirginQueens(founderGenomes)
+#'
+#' drones <- createDrones(x = basePop[1], nInd = 100)
+#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 5, nDrones = nFathersPoisson)
+#'
+#' # Create and cross Colony and MultiColony class
+#' colony <- createColony(x = basePop[2])
+#' colony <- cross(colony, drones = droneGroups[[1]])
+#' apiary <- createMultiColony(basePop[4:5])
+#' SIMplyBee:::setEvents(apiary, slot = "swarm", value = c(TRUE, TRUE))
+setEvents <- function(x, slot, value, simParamBee = NULL) {
+ if (is.null(simParamBee)) {
+ simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ }
+ if (isColony(x)) {
+ slot(x, slot) <- value
+ }
+ if (isMultiColony(x)) {
+ registerDoParallel(cores = simParamBee$nThreads)
+ x@colonies <- foreach(colony = seq_len(nColonies(x))) %dopar% {
+ setEvents(x[[colony]], slot, value)
+ }
+ }
+ return(x)
+}
+
+
#' @rdname combine
#' @title Combine two colony objects
#'
@@ -1695,6 +1875,7 @@ split <- function(x, p = NULL, year = NULL, simParamBee = NULL, ...) {
#'
#' @param strong \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
#' @param weak \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
+#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#'
#' @return a combined \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
#'
@@ -1702,8 +1883,9 @@ split <- function(x, p = NULL, year = NULL, simParamBee = NULL, ...) {
#' founderGenomes <- quickHaplo(nInd = 10, nChr = 1, segSites = 50)
#' SP <- SimParamBee$new(founderGenomes)
#' \dontshow{SP$nThreads = 1L}
+#' print(SP$nThreads)
#' basePop <- createVirginQueens(founderGenomes)
-#' drones <- createDrones(basePop[1], n = 1000)
+#' drones <- createDrones(basePop[1], nInd = 1000)
#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = 10)
#'
#' # Create weak and strong Colony and MultiColony class
@@ -1732,26 +1914,32 @@ split <- function(x, p = NULL, year = NULL, simParamBee = NULL, ...) {
#'
#' nWorkers(apiary1); nWorkers(apiary2)
#' nDrones(apiary1); nDrones(apiary2)
-#' apiary1 <- combine(strong = apiary1, weak = apiary2)
+#' apiary1 <- combine(strong = apiary1, weak = apiary2, simParamBee = SP)
#' nWorkers(apiary1); nWorkers(apiary2)
#' nDrones(apiary1); nDrones(apiary2)
#' rm(apiary2)
#' @export
-combine <- function(strong, weak) {
+combine <- function(strong, weak, simParamBee = NULL) {
+ if (is.null(simParamBee)) {
+ simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ }
+ if (any(hasCollapsed(strong))) {
+ stop(paste0("Some of the strong colonies have collapsed, hence you can not combine it!"))
+ }
+ if (any(hasCollapsed(weak))) {
+ stop(paste0("Some of the weak colonies have collapsed, hence you can not combine it!"))
+ }
if (isColony(strong) & isColony(weak)) {
- if (hasCollapsed(strong)) {
- stop(paste0("The colony ", getId(strong), " (strong) has collapsed, hence you can not combine it!"))
- }
- if (hasCollapsed(weak)) {
- stop(paste0("The colony ", getId(weak), " (weak) has collapsed, hence you can not combine it!"))
- }
strong@workers <- c(strong@workers, weak@workers)
strong@drones <- c(strong@drones, weak@drones)
} else if (isMultiColony(strong) & isMultiColony(weak)) {
+ registerDoParallel(cores = simParamBee$nThreads)
if (nColonies(weak) == nColonies(strong)) {
nCol <- nColonies(weak)
- for (colony in seq_len(nCol)) {
- strong[[colony]] <- combine(strong = strong[[colony]], weak = weak[[colony]])
+ strong@colonies <- foreach(colony = seq_len(nCol)) %dopar% {
+ combine(strong = strong[[colony]],
+ weak = weak[[colony]],
+ simParamBee = simParamBee)
}
} else {
stop("Weak and strong MultiColony objects must be of the same length!")
@@ -1762,6 +1950,7 @@ combine <- function(strong, weak) {
return(strong)
}
+
#' @rdname setLocation
#' @title Set colony location
#'
@@ -1774,6 +1963,7 @@ combine <- function(strong, weak) {
#' \code{c(x1, y1)} (the same location set to all colonies),
#' \code{list(c(x1, y1), c(x2, y2))}, or
#' \code{data.frame(x = c(x1, x2), y = c(y1, y2))}
+#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
#'
#' @return \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with set
#' location
@@ -1812,7 +2002,10 @@ combine <- function(strong, weak) {
#' apiary <- setLocation(apiary, location = locDF)
#' getLocation(apiary)
#' @export
-setLocation <- function(x, location = c(0, 0)) {
+setLocation <- function(x, location = c(0, 0), simParamBee = NULL) {
+ if (is.null(simParamBee)) {
+ simParamBee <- get(x = "SP", envir = .GlobalEnv)
+ }
if (isColony(x)) {
if (is.list(location)) { # is.list() captures also is.data.frame()
stop("Argument location must be numeric, when x is a Colony class object!")
@@ -1822,21 +2015,25 @@ setLocation <- function(x, location = c(0, 0)) {
}
x@location <- location
} else if (isMultiColony(x)) {
- n <- nColonies(x)
+ registerDoParallel(cores = simParamBee$nThreads)
+ nCol <- nColonies(x)
+ if (nCol == 0) {
+ stop("The Multicolony contains 0 colonies!")
+ }
if (!is.null(location)) {
if (is.numeric(location)) {
if (length(location) != 2) {
stop("When argument location is a numeric, it must be of length 2!")
}
} else if (is.data.frame(location)) {
- if (nrow(location) != n) {
+ if (nrow(location) != nCol) {
stop("When argument location is a data.frame, it must have as many rows as the number of colonies!")
}
if (ncol(location) != 2) {
stop("When argument location is a data.frame, it must have 2 columns!")
}
} else if (is.list(location)) {
- if (length(location) != n) {
+ if (length(location) != nCol) {
stop("When argument location is a list, it must be of length equal to the number of colonies!")
}
tmp <- sapply(X = location, FUN = length)
@@ -1851,7 +2048,14 @@ setLocation <- function(x, location = c(0, 0)) {
stop("Argument location must be numeric, list, or data.frame!")
}
}
- for (colony in seq_len(n)) {
+ combine_list <- function(a, b) {
+ if (length(a) == 1) {
+ c(list(a), list(b))
+ } else {
+ c(a, list(b))
+ }
+ }
+ tmp <- foreach(colony = seq_len(nCol), .combine = combine_list) %dopar% {
if (is.data.frame(location)) {
loc <- location[colony, ]
loc <- c(loc$x, loc$y)
@@ -1860,9 +2064,17 @@ setLocation <- function(x, location = c(0, 0)) {
} else {
loc <- location
}
+
if (!is.null(x[[colony]])) {
x[[colony]]@location <- loc
}
+
+ x[[colony]]
+ }
+ if (nCol == 1) {
+ x@colonies = list(tmp)
+ } else {
+ x@colonies = tmp
}
} else {
stop("Argument x must be a Colony or MultiColony class object!")
diff --git a/R/Functions_L3_Colonies.R b/R/Functions_L3_Colonies.R
index b8dbc191..b6535b9e 100644
--- a/R/Functions_L3_Colonies.R
+++ b/R/Functions_L3_Colonies.R
@@ -13,6 +13,7 @@
#' given then \code{\link[SIMplyBee]{MultiColony-class}} is created with \code{n}
#' \code{NULL}) individual colony - this is mostly useful for programming)
#' @param simParamBee \code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters
+#' @param populateColonies boolean, whether to create n empty Colony objects within with assigned ID
#'
#' @details When both \code{x} and \code{n} are \code{NULL}, then a
#' \code{\link[SIMplyBee]{MultiColony-class}} with 0 colonies is created.
@@ -40,7 +41,7 @@
#'
#' # Create mated colonies by crossing
#' apiary <- createMultiColony(x = basePop[1:2], n = 2)
-#' drones <- createDrones(x = basePop[3], n = 30)
+#' drones <- createDrones(x = basePop[3], nInd = 30)
#' droneGroups <- pullDroneGroupsFromDCA(drones, n = 2, nDrones = 15)
#' apiary <- cross(apiary, drones = droneGroups)
#' apiary
@@ -48,15 +49,26 @@
#' apiary[[2]]
#'
#' @export
-createMultiColony <- function(x = NULL, n = NULL, simParamBee = NULL) {
+createMultiColony <- function(x = NULL, n = NULL, simParamBee = NULL, populateColonies = FALSE) {
if (is.null(simParamBee)) {
simParamBee <- get(x = "SP", envir = .GlobalEnv)
}
+
+ registerDoParallel(cores = simParamBee$nThreads)
if (is.null(x)) {
if (is.null(n)) {
ret <- new(Class = "MultiColony")
} else {
ret <- new(Class = "MultiColony", colonies = vector(mode = "list", length = n))
+ if (populateColonies) {
+ ids <- (simParamBee$lastColonyId+1):(simParamBee$lastColonyId + n)
+ ret@colonies <- foreach(colony = seq_len(n)) %dopar% {
+ createColony(simParamBee = simParamBee, id = ids[colony])
+ }
+ simParamBee$updateLastColonyId(n = n)
+ } else {
+
+ }
}
} else {
if (!isPop(x)) {
@@ -72,9 +84,11 @@ createMultiColony <- function(x = NULL, n = NULL, simParamBee = NULL) {
stop("Not enough individuals in the x to create n colonies!")
}
ret <- new(Class = "MultiColony", colonies = vector(mode = "list", length = n))
- for (colony in seq_len(n)) {
- ret[[colony]] <- createColony(x = x[colony], simParamBee = simParamBee)
+ ids <- (simParamBee$lastColonyId+1):(simParamBee$lastColonyId + n)
+ ret@colonies <- foreach(colony = seq_len(n)) %dopar% {
+ createColony(x = x[colony], simParamBee = simParamBee, id = ids[colony])
}
+ simParamBee$updateLastColonyId(n = n)
}
validObject(ret)
return(ret)
diff --git a/R/SIMplyBee.R b/R/SIMplyBee.R
index ce9260f8..774af783 100644
--- a/R/SIMplyBee.R
+++ b/R/SIMplyBee.R
@@ -7,6 +7,8 @@
#' @importFrom stats rnorm rbeta runif rpois na.omit
#' @importFrom extraDistr rtpois
#' @importFrom utils packageVersion
+#' @importFrom foreach foreach "%dopar%"
+#' @importFrom doParallel registerDoParallel
# see https://r-pkgs.org/namespace.html on description what to import/depend/...
#' @description
diff --git a/man/SimParamBee.Rd b/man/SimParamBee.Rd
index 664475f4..00f6e5f5 100644
--- a/man/SimParamBee.Rd
+++ b/man/SimParamBee.Rd
@@ -317,6 +317,10 @@ generate this object}
\item \href{#method-SimParamBee-new}{\code{SimParamBee$new()}}
\item \href{#method-SimParamBee-addToCaste}{\code{SimParamBee$addToCaste()}}
\item \href{#method-SimParamBee-changeCaste}{\code{SimParamBee$changeCaste()}}
+\item \href{#method-SimParamBee-addToBeePed}{\code{SimParamBee$addToBeePed()}}
+\item \href{#method-SimParamBee-addToBeeRec}{\code{SimParamBee$addToBeeRec()}}
+\item \href{#method-SimParamBee-updateCaste}{\code{SimParamBee$updateCaste()}}
+\item \href{#method-SimParamBee-updateLastBeeId}{\code{SimParamBee$updateLastBeeId()}}
\item \href{#method-SimParamBee-updateLastColonyId}{\code{SimParamBee$updateLastColonyId()}}
\item \href{#method-SimParamBee-clone}{\code{SimParamBee$clone()}}
}
@@ -532,6 +536,99 @@ SP$caste
}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-SimParamBee-addToBeePed}{}}}
+\subsection{Method \code{addToBeePed()}}{
+For internal use only.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{SimParamBee$addToBeePed(nNewInd, id, mother, father, isDH)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{nNewInd}}{Number of newly created individuals}
+
+\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-SimParamBee-addToBeeRec}{}}}
+\subsection{Method \code{addToBeeRec()}}{
+For internal use only.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{SimParamBee$addToBeeRec(nNewInd, id, mother, father, isDH, hist, ploidy)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{nNewInd}}{Number of newly created individuals}
+
+\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-SimParamBee-updateCaste}{}}}
+\subsection{Method \code{updateCaste()}}{
+A function to update the caste
+ For internal use only.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{SimParamBee$updateCaste(caste)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{caste}}{vector, named vector of castes to be added}
+}
+\if{html}{\out{
}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-SimParamBee-updateLastBeeId}{}}}
+\subsection{Method \code{updateLastBeeId()}}{
+A function to update the last
+ ID everytime we create an individual
+ For internal use in SIMplyBee only.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{SimParamBee$updateLastBeeId(n = 1L)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{n}}{integer, how many individuals to add}
+
+\item{\code{lastId}}{integer, last colony ID assigned}
+}
+\if{html}{\out{
}}
+}
}
\if{html}{\out{
}}
\if{html}{\out{}}
@@ -541,12 +638,14 @@ A function to update the colony last
ID everytime we create a Colony-class with createColony.
For internal use only.
\subsection{Usage}{
-\if{html}{\out{}}\preformatted{SimParamBee$updateLastColonyId()}\if{html}{\out{
}}
+\if{html}{\out{}}\preformatted{SimParamBee$updateLastColonyId(n = 1)}\if{html}{\out{
}}
}
\subsection{Arguments}{
\if{html}{\out{}}
\describe{
+\item{\code{n}}{integer, how many colonies to add}
+
\item{\code{lastColonyId}}{integer, last colony ID assigned}
}
\if{html}{\out{
}}
diff --git a/man/addCastePop.Rd b/man/addCastePop.Rd
index 89e896e3..9c58ba26 100644
--- a/man/addCastePop.Rd
+++ b/man/addCastePop.Rd
@@ -12,13 +12,12 @@ addCastePop(
caste = NULL,
nInd = NULL,
new = FALSE,
- exact = FALSE,
year = NULL,
simParamBee = NULL,
...
)
-addWorkers(x, nInd = NULL, new = FALSE, exact = FALSE, simParamBee = NULL, ...)
+addWorkers(x, nInd = NULL, new = FALSE, simParamBee = NULL, ...)
addDrones(x, nInd = NULL, new = FALSE, simParamBee = NULL, ...)
@@ -45,10 +44,6 @@ a single value is provided, the same value will be used for all the colonies.}
\item{new}{logical, should the number of individuals be added to the caste population
anew or should we only top-up the existing number of individuals to \code{nInd}}
-\item{exact}{logical, only relevant when adding workers - if the csd locus is turned
-on and exact is \code{TRUE}, we add the exact specified number of viable workers
-(heterozygous at the csd locus)}
-
\item{year}{numeric, only relevant when adding virgin queens - year of birth for virgin queens}
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
diff --git a/man/addCastePop_internal.Rd b/man/addCastePop_internal.Rd
new file mode 100644
index 00000000..bba89e34
--- /dev/null
+++ b/man/addCastePop_internal.Rd
@@ -0,0 +1,24 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/Functions_L2_Colony.R
+\name{addCastePop_internal}
+\alias{addCastePop_internal}
+\title{An internal function to add a population in a caste slot of the colony}
+\usage{
+addCastePop_internal(pop, colony, caste, new = FALSE)
+}
+\arguments{
+\item{pop}{\code{\link[AlphaSimR]{Pop-class}} with one or many individual}
+
+\item{colony}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
+
+\item{caste}{character}
+
+\item{new}{logical}
+}
+\value{
+\code{\link[SIMplyBee]{Colony-class}}
+}
+\description{
+Helper function that returns a colony to allow parallelisation,
+only for internal use.
+}
diff --git a/man/buildUp.Rd b/man/buildUp.Rd
index 5e280a04..21df8acc 100644
--- a/man/buildUp.Rd
+++ b/man/buildUp.Rd
@@ -9,7 +9,6 @@ buildUp(
nWorkers = NULL,
nDrones = NULL,
new = TRUE,
- exact = FALSE,
resetEvents = FALSE,
simParamBee = NULL,
...
@@ -34,10 +33,6 @@ a single value is provided, the same value will be applied to all the colonies.}
should we only top-up the existing number of workers and drones to
\code{nWorkers} and \code{nDrones} (see details)}
-\item{exact}{logical, if the csd locus is turned on and exact is \code{TRUE},
-create the exact specified number of only viable workers (heterozygous on
-the csd locus)}
-
\item{resetEvents}{logical, call \code{\link[SIMplyBee]{resetEvents}} as part of the
build up}
diff --git a/man/collapse.Rd b/man/collapse.Rd
index e00a37b2..34326080 100644
--- a/man/collapse.Rd
+++ b/man/collapse.Rd
@@ -4,10 +4,12 @@
\alias{collapse}
\title{Collapse}
\usage{
-collapse(x)
+collapse(x, simParamBee = NULL)
}
\arguments{
\item{x}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
+
+\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
}
\value{
\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with the collapse
@@ -30,13 +32,13 @@ founderGenomes <- quickHaplo(nInd = 10, nChr = 1, segSites = 50)
SP <- SimParamBee$new(founderGenomes)
\dontshow{SP$nThreads = 1L}
basePop <- createVirginQueens(founderGenomes)
-drones <- createDrones(basePop[1], n = 1000)
+drones <- createDrones(basePop[1], nInd = 1000)
droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = 10)
# Create Colony and MultiColony class
colony <- createColony(x = basePop[1])
colony <- cross(colony, drones = droneGroups[[1]])
-apiary <- createMultiColony(x = basePop[2:10], n = 9)
+apiary <- createMultiColony(x = basePop[2:10])
apiary <- cross(apiary, drones = droneGroups[2:10])
# Collapse
diff --git a/man/combine.Rd b/man/combine.Rd
index c14a3a67..9e98b85b 100644
--- a/man/combine.Rd
+++ b/man/combine.Rd
@@ -4,12 +4,14 @@
\alias{combine}
\title{Combine two colony objects}
\usage{
-combine(strong, weak)
+combine(strong, weak, simParamBee = NULL)
}
\arguments{
\item{strong}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
\item{weak}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
+
+\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
}
\value{
a combined \code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}
@@ -26,8 +28,9 @@ Level 2 function that combines two Colony or MultiColony objects
founderGenomes <- quickHaplo(nInd = 10, nChr = 1, segSites = 50)
SP <- SimParamBee$new(founderGenomes)
\dontshow{SP$nThreads = 1L}
+print(SP$nThreads)
basePop <- createVirginQueens(founderGenomes)
-drones <- createDrones(basePop[1], n = 1000)
+drones <- createDrones(basePop[1], nInd = 1000)
droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = 10)
# Create weak and strong Colony and MultiColony class
@@ -56,7 +59,7 @@ rm(colony2)
nWorkers(apiary1); nWorkers(apiary2)
nDrones(apiary1); nDrones(apiary2)
-apiary1 <- combine(strong = apiary1, weak = apiary2)
+apiary1 <- combine(strong = apiary1, weak = apiary2, simParamBee = SP)
nWorkers(apiary1); nWorkers(apiary2)
nDrones(apiary1); nDrones(apiary2)
rm(apiary2)
diff --git a/man/createCastePop.Rd b/man/createCastePop.Rd
index 4766cbd1..6cc9be6c 100644
--- a/man/createCastePop.Rd
+++ b/man/createCastePop.Rd
@@ -11,17 +11,32 @@ createCastePop(
x,
caste = NULL,
nInd = NULL,
- exact = TRUE,
year = NULL,
editCsd = TRUE,
csdAlleles = NULL,
simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL,
...
)
-createWorkers(x, nInd = NULL, exact = FALSE, simParamBee = NULL, ...)
+createWorkers(
+ x,
+ nInd = NULL,
+ simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL,
+ ...
+)
-createDrones(x, nInd = NULL, simParamBee = NULL, ...)
+createDrones(
+ x,
+ nInd = NULL,
+ simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL,
+ ...
+)
createVirginQueens(
x,
@@ -30,6 +45,8 @@ createVirginQueens(
editCsd = TRUE,
csdAlleles = NULL,
simParamBee = NULL,
+ returnSP = FALSE,
+ ids = NULL,
...
)
}
@@ -47,11 +64,6 @@ only used when \code{x} is \code{\link[SIMplyBee]{Colony-class}} or
\code{\link[SIMplyBee]{MultiColony-class}}, when \code{x} is \code{link[AlphaSimR]{MapPop-class}}
all individuals in \code{x} are converted into virgin queens}
-\item{exact}{logical, only relevant when creating workers,
-if the csd locus is active and exact is \code{TRUE},
-create the exactly specified number of viable workers (heterozygous on the
-csd locus)}
-
\item{year}{numeric, year of birth for virgin queens}
\item{editCsd}{logical (only active when \code{x} is \code{link[AlphaSimR]{MapPop-class}}),
@@ -70,6 +82,13 @@ ensure heterozygosity at the csd locus.}
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
+\item{returnSP}{logical, whether to return the pedigree, caste, and recHist information
+for each created population (used internally for parallel computing)}
+
+\item{ids}{character, IDs of the individuals that are going to be created (used internally
+for parallel computing)}
+\item{nThreads}{integer, number of cores to use for parallel computing (over colonies)}
+
\item{...}{additional arguments passed to \code{nInd} when this argument is a function}
}
\value{
@@ -100,7 +119,7 @@ Level 1 function that creates the specified number of caste
founderGenomes <- quickHaplo(nInd = 4, nChr = 1, segSites = 50)
SP <- SimParamBee$new(founderGenomes)
\dontshow{SP$nThreads = 1L}
-SP$setTrackRec(TRUE)
+SP$setTrackRec(isTrackRec = TRUE)
SP$setTrackPed(isTrackPed = TRUE)
# Create virgin queens on a MapPop
diff --git a/man/createColony.Rd b/man/createColony.Rd
index a8a96649..d1707f07 100644
--- a/man/createColony.Rd
+++ b/man/createColony.Rd
@@ -4,12 +4,14 @@
\alias{createColony}
\title{Create a new Colony}
\usage{
-createColony(x = NULL, simParamBee = NULL)
+createColony(x = NULL, simParamBee = NULL, id = NULL)
}
\arguments{
\item{x}{\code{\link[AlphaSimR]{Pop-class}}, one queen or virgin queen(s)}
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
+
+\item{id}{character, ID of the colony that is going to be created (used internally for parallel computing)}
}
\value{
new \code{\link[SIMplyBee]{Colony-class}}
diff --git a/man/createCrossPlan.Rd b/man/createCrossPlan.Rd
index ae90395f..8db00494 100644
--- a/man/createCrossPlan.Rd
+++ b/man/createCrossPlan.Rd
@@ -69,10 +69,6 @@ virginColonies2 <- createMultiColony(basePop[31:60])
virginColonies2 <- setLocation(virginColonies2,
location = Map(c, runif(30, 0, 2*pi),
runif(30, 0, 2*pi)))
-virginColonies3 <- createMultiColony(basePop[61:90])
-virginColonies3 <- setLocation(virginColonies3,
- location = Map(c, runif(30, 0, 2*pi),
- runif(30, 0, 2*pi)))
# Create drone colonies
droneColonies <- createMultiColony(basePop[121:200])
@@ -92,57 +88,22 @@ droneColonies <- cross(droneColonies,
nDrones = nFathersPoisson,
crossPlan = randomCrossPlan)
-# Plot the colonies in space
-virginLocations <- as.data.frame(getLocation(c(virginColonies1, virginColonies2, virginColonies3),
- collapse= TRUE))
-virginLocations$Type <- "Virgin"
-droneLocations <- as.data.frame(getLocation(droneColonies, collapse= TRUE))
-droneLocations$Type <- "Drone"
-locations <- rbind(virginLocations, droneLocations)
-
-plot(x = locations$V1, y = locations$V2,
- col = c("red", "blue")[as.numeric(as.factor(locations$Type))])
-
-# Cross according to a spatial cross plan according to the colonies' locations
-crossPlanSpatial <- createCrossPlan(x = virginColonies1,
- droneColonies = droneColonies,
- nDrones = nFathersPoisson,
- spatial = TRUE,
- radius = 1.5)
-
-# Plot the crossing for the first colony in the crossPlan
-virginLocations1 <- as.data.frame(getLocation(virginColonies1, collapse= TRUE))
-virginLocations1$Type <- "Virgin"
-droneLocations <- as.data.frame(getLocation(droneColonies, collapse= TRUE))
-droneLocations$Type <- "Drone"
-locations1 <- rbind(virginLocations1, droneLocations)
-
-# Blue marks the target virgin colony and blue marks the drone colonies in the chosen radius
-plot(x = locations1$V1, y = locations1$V2, pch = c(1, 2)[as.numeric(as.factor(locations1$Type))],
- col = ifelse(rownames(locations1) \%in\% crossPlanSpatial[[1]],
- "red",
- ifelse(rownames(locations1) == names(crossPlanSpatial)[[1]],
- "blue", "black")))
-
-colonies1 <- cross(x = virginColonies1,
- crossPlan = crossPlanSpatial,
- droneColonies = droneColonies,
- nDrones = nFathersPoisson)
-nFathers(colonies1)
-
# Cross according to a cross plan that is created internally within the cross function
# The cross plan is created at random, regardless the location of the colonies
-colonies2 <- cross(x = virginColonies2,
+colonies1 <- cross(x = virginColonies1,
droneColonies = droneColonies,
nDrones = nFathersPoisson,
crossPlan = "create")
-# Mate spatially with cross plan created internally by the cross function
-colonies3 <- cross(x = virginColonies3,
- droneColonies = droneColonies,
+
+# Cross according to a spatial cross plan created internally according to the colonies' locations
+colonies2 <- cross(x = virginColonies2,
crossPlan = "create",
- checkCross = "warning",
+ droneColonies = droneColonies,
+ nDrones = nFathersPoisson,
spatial = TRUE,
- radius = 1)
+ radius = 1.5)
+nFathers(colonies2)
+
}
diff --git a/man/createMatingStationDCA.Rd b/man/createMatingStationDCA.Rd
index 8ccf320d..58c138a7 100644
--- a/man/createMatingStationDCA.Rd
+++ b/man/createMatingStationDCA.Rd
@@ -35,7 +35,7 @@ founderGenomes <- quickHaplo(nInd = 10, nChr = 1, segSites = 100)
SP <- SimParamBee$new(founderGenomes)
\dontshow{SP$nThreads = 1L}
basePop <- createVirginQueens(founderGenomes)
-drones <- createDrones(basePop[1], n = 1000)
+drones <- createDrones(basePop[1], nInd = 1000)
droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = 10)
# Create a colony and cross it
diff --git a/man/createMultiColony.Rd b/man/createMultiColony.Rd
index 21f2bc14..b1c48ec3 100644
--- a/man/createMultiColony.Rd
+++ b/man/createMultiColony.Rd
@@ -4,7 +4,7 @@
\alias{createMultiColony}
\title{Create MultiColony object}
\usage{
-createMultiColony(x = NULL, n = NULL, simParamBee = NULL)
+createMultiColony(x = NULL, n = NULL, simParamBee = NULL, nThreads = NULL)
}
\arguments{
\item{x}{\code{\link[AlphaSimR]{Pop-class}}, virgin queens or queens for the colonies
@@ -16,6 +16,8 @@ given then \code{\link[SIMplyBee]{MultiColony-class}} is created with \code{n}
\code{NULL}) individual colony - this is mostly useful for programming)}
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
+
+\item{nThreads}{integer, number of cores to use for parallel computing (over colonies)}
}
\value{
\code{\link[SIMplyBee]{MultiColony-class}}
@@ -49,7 +51,7 @@ apiary[[2]]
# Create mated colonies by crossing
apiary <- createMultiColony(x = basePop[1:2], n = 2)
-drones <- createDrones(x = basePop[3], n = 30)
+drones <- createDrones(x = basePop[3], nInd = 30)
droneGroups <- pullDroneGroupsFromDCA(drones, n = 2, nDrones = 15)
apiary <- cross(apiary, drones = droneGroups)
apiary
diff --git a/man/cross.Rd b/man/cross.Rd
index 8b4ee5fe..8d34e3a6 100644
--- a/man/cross.Rd
+++ b/man/cross.Rd
@@ -54,7 +54,8 @@ to their distance from the virgin colony (that is, in a radius)}
\item{radius}{numeric, the radius around the virgin colony in which to sample mating partners,
only needed when \code{spatial = TRUE}}
-\item{checkCross}{character, throw a warning (when \code{checkCross = "warning"}),}
+\item{checkCross}{character, throw a warning (when \code{checkCross = "warning"}).
+This will also remove the unmated queens and return only the mated ones.}
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
@@ -95,7 +96,7 @@ This function changes caste for the mated drones to fathers, and
\examples{
founderGenomes <- quickHaplo(nInd = 30, nChr = 1, segSites = 100)
SP <- SimParamBee$new(founderGenomes)
-\dontshow{SP$nThreads = 1L}
+SP$nThreads = 1L
basePop <- createVirginQueens(founderGenomes)
drones <- createDrones(x = basePop[1], nInd = 1000)
diff --git a/man/getIbdHaplo.Rd b/man/getIbdHaplo.Rd
index fd6dfe79..08895a7a 100644
--- a/man/getIbdHaplo.Rd
+++ b/man/getIbdHaplo.Rd
@@ -115,10 +115,10 @@ Level 0 function that returns IBD (identity by descent)
}}
\examples{
-founderGenomes <- quickHaplo(nInd = 4, nChr = 1, segSites = 50)
-SP <- SimParamBee$new(founderGenomes)
+founderGenomes <- quickHaplo(nInd = 4, nChr = 1, segSites = 5)
+SP <- SimParamBee$new(founderGenomes, nCsdAlleles = 4)
\dontshow{SP$nThreads = 1L}
-SP$setTrackRec(TRUE)
+SP$setTrackRec(isTrackRec = TRUE)
SP$setTrackPed(isTrackPed = TRUE)
basePop <- createVirginQueens(founderGenomes)
@@ -128,13 +128,13 @@ droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = nFathersPoisson)
# Create a Colony and a MultiColony class
colony <- createColony(x = basePop[2])
colony <- cross(colony, drones = droneGroups[[1]])
-colony <- buildUp(x = colony, nWorkers = 6, nDrones = 3)
-colony <- addVirginQueens(x = colony, nInd = 5)
+colony <- buildUp(x = colony, nWorkers = 3, nDrones = 2)
+colony <- addVirginQueens(x = colony, nInd = 2)
-apiary <- createMultiColony(basePop[3:4], n = 2)
+apiary <- createMultiColony(basePop[3:4])
apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-apiary <- buildUp(x = apiary, nWorkers = 6, nDrones = 3)
-apiary <- addVirginQueens(x = apiary, nInd = 5)
+apiary <- buildUp(x = apiary, nWorkers = 3, nDrones = 2)
+apiary <- addVirginQueens(x = apiary, nInd = 2)
# Input is a population
getIbdHaplo(x = getQueen(colony))
@@ -146,6 +146,8 @@ getIbdHaplo(x = colony, caste = "queen")
getQueenIbdHaplo(colony)
getIbdHaplo(colony, caste = "workers", nInd = 3)
+getIbdHaplo(colony, caste = "virginQueens")
+getIbdHaplo(colony, caste = "drones")
getWorkersIbdHaplo(colony)
# Same aliases exist for all castes!
@@ -160,6 +162,9 @@ getQueenIbdHaplo(apiary)
# Or collapse all the haplotypes into a single matrix
getQueenIbdHaplo(apiary, collapse = TRUE)
+
+getIbdHaplo(x = apiary, caste = "workers")
+getIbdHaplo(x = apiary, caste = "drones")
# Get the haplotypes of all individuals either by colony or in a single matrix
getIbdHaplo(apiary, caste = "all")
getIbdHaplo(apiary, caste = "all", collapse = TRUE)
diff --git a/man/hasSwarmed.Rd b/man/hasSwarmed.Rd
index 870d2d14..427d41d3 100644
--- a/man/hasSwarmed.Rd
+++ b/man/hasSwarmed.Rd
@@ -31,7 +31,7 @@ colony <- cross(colony, drones = droneGroups[[1]])
colony <- buildUp(x = colony, nWorkers = 6, nDrones = 3)
colony <- addVirginQueens(colony, nInd = 5)
-apiary <- createMultiColony(basePop[3:4], n = 2)
+apiary <- createMultiColony(basePop[3:4])
apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
apiary <- buildUp(x = apiary, nWorkers = 6, nDrones = 3)
diff --git a/man/pullCastePop.Rd b/man/pullCastePop.Rd
index dcf63748..42a890b2 100644
--- a/man/pullCastePop.Rd
+++ b/man/pullCastePop.Rd
@@ -95,12 +95,12 @@ droneGroups <- pullDroneGroupsFromDCA(drones, n = 10, nDrones = nFathersPoisson)
# Create a Colony and a MultiColony class
colony <- createColony(x = basePop[2])
colony <- cross(colony, drones = droneGroups[[1]])
-colony <- buildUp(x = colony, nWorkers = 100, nDrones = 10, exact = TRUE)
+colony <- buildUp(x = colony, nWorkers = 100, nDrones = 10)
colony <- addVirginQueens(x = colony, nInd = 3)
apiary <- createMultiColony(basePop[3:4], n = 2)
apiary <- cross(apiary, drones = droneGroups[c(2, 3)])
-apiary <- buildUp(x = apiary, nWorkers = 100, nDrones = 10, exact = TRUE)
+apiary <- buildUp(x = apiary, nWorkers = 100, nDrones = 10)
apiary <- addVirginQueens(x = apiary, nInd = 3)
# pullCastePop on Colony class
diff --git a/man/removeCastePop.Rd b/man/removeCastePop.Rd
index 9acecac0..e47267ed 100644
--- a/man/removeCastePop.Rd
+++ b/man/removeCastePop.Rd
@@ -13,19 +13,11 @@ removeCastePop(
caste = NULL,
p = 1,
use = "rand",
- addVirginQueens = FALSE,
- nVirginQueens = NULL,
year = NULL,
simParamBee = NULL
)
-removeQueen(
- x,
- addVirginQueens = FALSE,
- nVirginQueens = NULL,
- year = NULL,
- simParamBee = NULL
-)
+removeQueen(x, year = NULL, simParamBee = NULL)
removeWorkers(x, p = 1, use = "rand", simParamBee = NULL)
@@ -45,14 +37,6 @@ a single value is provided, the same value will be applied to all the colonies}
\item{use}{character, all the options provided by \code{\link[AlphaSimR]{selectInd}} -
guides selection of virgins queens that will stay when \code{p < 1}}
-\item{addVirginQueens}{logical, whether virgin queens should be added; only
-used when removing the queen from the colony}
-
-\item{nVirginQueens}{integer, the number of virgin queens to be created in the
-colony; only used when removing the queen from the colony. If \code{0}, no virgin
-queens are added; If \code{NULL}, the value from \code{simParamBee$nVirginQueens}
-is used}
-
\item{year}{numeric, only relevant when adding virgin queens - year of birth for virgin queens}
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
diff --git a/man/replaceCastePop.Rd b/man/replaceCastePop.Rd
index bd2c3756..c23232cc 100644
--- a/man/replaceCastePop.Rd
+++ b/man/replaceCastePop.Rd
@@ -12,12 +12,11 @@ replaceCastePop(
caste = NULL,
p = 1,
use = "rand",
- exact = TRUE,
year = NULL,
simParamBee = NULL
)
-replaceWorkers(x, p = 1, use = "rand", exact = TRUE, simParamBee = NULL)
+replaceWorkers(x, p = 1, use = "rand", simParamBee = NULL)
replaceDrones(x, p = 1, use = "rand", simParamBee = NULL)
@@ -36,11 +35,6 @@ a single value is provided, the same value will be applied to all the colonies}
\item{use}{character, all the options provided by \code{\link[AlphaSimR]{selectInd}} -
guides selection of caste individuals that stay when \code{p < 1}}
-\item{exact}{logical, only relevant when adding workers - if the csd locus is turned
-on and exact is \code{TRUE}, we replace the exact specified number of viable workers
-(heterozygous at the csd locus). You probably want this set to TRUE since you want to
-replace with the same number of workers.}
-
\item{year}{numeric, only relevant when replacing virgin queens,
year of birth for virgin queens}
diff --git a/man/resetEvents.Rd b/man/resetEvents.Rd
index 2e8b6642..cada3fd6 100644
--- a/man/resetEvents.Rd
+++ b/man/resetEvents.Rd
@@ -4,7 +4,7 @@
\alias{resetEvents}
\title{Reset colony events}
\usage{
-resetEvents(x, collapse = NULL)
+resetEvents(x, collapse = NULL, simParamBee = NULL)
}
\arguments{
\item{x}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
@@ -13,6 +13,8 @@ resetEvents(x, collapse = NULL)
up a new colony, which the default of \code{NULL} caters for; otherwise, a
collapsed colony should be left collapsed forever, unless you force
resetting this event with \code{collapse = TRUE})}
+
+\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
}
\value{
\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with
diff --git a/man/setEvents.Rd b/man/setEvents.Rd
new file mode 100644
index 00000000..5035cd70
--- /dev/null
+++ b/man/setEvents.Rd
@@ -0,0 +1,40 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/Functions_L2_Colony.R
+\name{setEvents}
+\alias{setEvents}
+\title{Set colony events}
+\usage{
+setEvents(x, slot, value, simParamBee = NULL)
+}
+\arguments{
+\item{x}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
+
+\item{slot}{character, which event to set}
+
+\item{value}{logical, the value for the event}
+
+\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
+}
+\value{
+\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with
+ events reset
+}
+\description{
+Helper Level 2 function that populates the events slot. Not interded
+ for external use, intended for internal use in parallel computing
+}
+\examples{
+founderGenomes <- quickHaplo(nInd = 5, nChr = 1, segSites = 50)
+SP <- SimParamBee$new(founderGenomes)
+\dontshow{SP$nThreads = 1L}
+basePop <- createVirginQueens(founderGenomes)
+
+drones <- createDrones(x = basePop[1], nInd = 100)
+droneGroups <- pullDroneGroupsFromDCA(drones, n = 5, nDrones = nFathersPoisson)
+
+# Create and cross Colony and MultiColony class
+colony <- createColony(x = basePop[2])
+colony <- cross(colony, drones = droneGroups[[1]])
+apiary <- createMultiColony(basePop[4:5])
+SIMplyBee:::setEvents(apiary, slot = "swarm", value = c(TRUE, TRUE))
+}
diff --git a/man/setLocation.Rd b/man/setLocation.Rd
index a9ee600d..e9596075 100644
--- a/man/setLocation.Rd
+++ b/man/setLocation.Rd
@@ -4,7 +4,7 @@
\alias{setLocation}
\title{Set colony location}
\usage{
-setLocation(x, location = c(0, 0))
+setLocation(x, location = c(0, 0), simParamBee = NULL)
}
\arguments{
\item{x}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
@@ -14,6 +14,8 @@ locations as
\code{c(x1, y1)} (the same location set to all colonies),
\code{list(c(x1, y1), c(x2, y2))}, or
\code{data.frame(x = c(x1, x2), y = c(y1, y2))}}
+
+\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
}
\value{
\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}} with set
diff --git a/man/supersede.Rd b/man/supersede.Rd
index 04291135..019a21cd 100644
--- a/man/supersede.Rd
+++ b/man/supersede.Rd
@@ -4,18 +4,13 @@
\alias{supersede}
\title{Supersede}
\usage{
-supersede(x, year = NULL, nVirginQueens = NULL, simParamBee = NULL, ...)
+supersede(x, year = NULL, simParamBee = NULL, ...)
}
\arguments{
\item{x}{\code{\link[SIMplyBee]{Colony-class}} or \code{\link[SIMplyBee]{MultiColony-class}}}
\item{year}{numeric, year of birth for virgin queens}
-\item{nVirginQueens}{integer, the number of virgin queens to be created in the
-colony; of these one is randomly selected as the new virgin queen of the
-remnant colony. If \code{NULL}, the value from \code{simParamBee$nVirginQueens}
-is used}
-
\item{simParamBee}{\code{\link[SIMplyBee]{SimParamBee}}, global simulation parameters}
\item{...}{additional arguments passed to \code{nVirginQueens} when this
diff --git a/man/swarm.Rd b/man/swarm.Rd
index e178fe26..72c7e88f 100644
--- a/man/swarm.Rd
+++ b/man/swarm.Rd
@@ -8,7 +8,6 @@ swarm(
x,
p = NULL,
year = NULL,
- nVirginQueens = NULL,
sampleLocation = TRUE,
radius = NULL,
simParamBee = NULL,
@@ -26,11 +25,6 @@ a single value is provided, the same value will be applied to all the colonies}
\item{year}{numeric, year of birth for virgin queens}
-\item{nVirginQueens}{integer, the number of virgin queens to be created in the
-colony; of these one is randomly selected as the new virgin queen of the
-remnant colony. If \code{NULL}, the value from \code{simParamBee$nVirginQueens}
-is used}
-
\item{sampleLocation}{logical, sample location of the swarm by taking
the current colony location and adding deviates to each coordinate using
\code{\link[SIMplyBee]{rcircle}}}
@@ -61,6 +55,8 @@ Level 2 function that swarms a Colony or MultiColony object -
\examples{
founderGenomes <- quickHaplo(nInd = 8, nChr = 1, segSites = 50)
SP <- SimParamBee$new(founderGenomes)
+SP$setTrackPed(TRUE)
+SP$setTrackRec(TRUE)
\dontshow{SP$nThreads = 1L}
basePop <- createVirginQueens(founderGenomes)
drones <- createDrones(basePop[1], n = 1000)
diff --git a/tests/testthat/test-L0_auxiliary_functions.R b/tests/testthat/test-L0_auxiliary_functions.R
index b9883b0a..6a8bf4d6 100644
--- a/tests/testthat/test-L0_auxiliary_functions.R
+++ b/tests/testthat/test-L0_auxiliary_functions.R
@@ -181,7 +181,7 @@ test_that("pHomBrood", {
expect_error(pHomBrood(colony@workers, simParamBee = SP))
expect_error(pHomBrood(colony@virginQueens, simParamBee = SP))
expect_error(pHomBrood(colony@drones, simParamBee = SP))
- expect_true(is.numeric(pHomBrood(colony@queen, simParamBee = SP)))
+ #expect_true(is.numeric(pHomBrood(colony@queen, simParamBee = SP)))
colony@queen <- NULL
expect_error(pHomBrood(colony@queen, simParamBee = SP))
diff --git a/tests/testthat/test-L1_pop_functions.R b/tests/testthat/test-L1_pop_functions.R
index 7cabd2b8..5d39e6a3 100644
--- a/tests/testthat/test-L1_pop_functions.R
+++ b/tests/testthat/test-L1_pop_functions.R
@@ -282,7 +282,7 @@ test_that("cross", {
expect_error(cross(colony2, drones = dronesGroups[7], simParamBee = SP))
# Message if fathers == 0 "Mating failed"
- expect_error(cross(virginQueen2, drones= selectInd(colony@drones,nInd = 0, use = "rand", simParam = SP), simParamBee = SP))
+ expect_error(cross(virginQueen2, drones= selectInd(colony@drones, nInd = 0, use = "rand", simParam = SP), simParamBee = SP))
#expect_message(cross(virginQueen2, drones= selectInd(colony@drones,nInd = 0, use = "rand", simParam = SP), checkCross = "warning", simParamBee = SP))
})
@@ -414,3 +414,4 @@ test_that("combineBeeGametesHaploidDiploid", {
expect_equal(nInd(drones), 5)
expect_equal(drones@ploidy, 1)
})
+
diff --git a/tests/testthat/test-L2_colony_functions.R b/tests/testthat/test-L2_colony_functions.R
index bb1a729e..2daeb877 100644
--- a/tests/testthat/test-L2_colony_functions.R
+++ b/tests/testthat/test-L2_colony_functions.R
@@ -109,9 +109,7 @@ test_that("Add functions", {
expect_equal(nDrones(addDrones(colony, nInd = 5, new = TRUE, simParamBee = SP), simParamBee = SP), 5)
# If input is an apiary
# Empty apiary - you can add, but nothing happens - returns an empty apiary
- expect_s4_class(addVirginQueens(emptyApiary, nInd = 5, simParamBee = SP), "MultiColony")
- expect_s4_class(addWorkers(emptyApiary, nInd = 5, simParamBee = SP), "MultiColony")
- expect_s4_class(addDrones(emptyApiary, nInd = 5, simParamBee = SP), "MultiColony")
+ expect_error(addVirginQueens(emptyApiary, nInd = 5, simParamBee = SP))
# Non-empty apiary
expect_s4_class(addVirginQueens(apiary, nInd = 5, simParamBee = SP), "MultiColony")
expect_s4_class(addWorkers(apiary, nInd = 5, simParamBee = SP), "MultiColony")
@@ -154,7 +152,7 @@ test_that("BuildUpDownsize", {
# Build Up an apiary
# Empty apiary
- expect_s4_class(buildUp(emptyApiary, simParamBee = SP), "MultiColony")
+ expect_error(buildUp(emptyApiary, simParamBee = SP))
# Non-empty apiary
expect_equal(nColonies(buildUp(apiary, simParamBee = SP)), 2)
@@ -169,7 +167,7 @@ test_that("BuildUpDownsize", {
expect_length(intersect(getId(getWorkers(downsize(colony, p = 0.1, new = TRUE, simParamBee = SP), simParamBee = SP)), workersIDs), 0)
# Empty apiary
- expect_s4_class(downsize(emptyApiary, simParamBee = SP), "MultiColony")
+ expect_error(downsize(emptyApiary, simParamBee = SP))
# Non-empty apiary
downsize(apiary, simParamBee = SP)
})
@@ -201,9 +199,9 @@ test_that("replaceFunctions", {
expect_error(replaceVirginQueens(emptyColony, p = 0.5, simParamBee = SP))
expect_error(replaceWorkers(emptyColony, p = 0, simParamBee = SP))
expect_error(replaceDrones(emptyColony, p = 1, simParamBee = SP))
- expect_s4_class(replaceVirginQueens(emptyApiary, p = 0.5, simParamBee = SP), "MultiColony")
- expect_s4_class(replaceWorkers(emptyApiary, p = 0, simParamBee = SP), "MultiColony")
- expect_s4_class(replaceDrones(emptyApiary, p = 1, simParamBee = SP), "MultiColony")
+ expect_error(replaceVirginQueens(emptyApiary, p = 0.5, simParamBee = SP))
+ expect_error(replaceWorkers(emptyApiary, p = 0, simParamBee = SP))
+ expect_error(replaceDrones(emptyApiary, p = 1, simParamBee = SP))
# Replace individuals in the non-empty colony/apiary
expect_s4_class(replaceVirginQueens(colony, simParamBee = SP), "Colony")
@@ -211,7 +209,7 @@ test_that("replaceFunctions", {
expect_s4_class(replaceDrones(colony, simParamBee = SP), "Colony")
expect_equal(nVirginQueens(replaceVirginQueens(colony, p = 1, simParamBee = SP), simParamBee = SP), nVirginQueens(colony, simParam = SP))
expect_equal(nWorkers(replaceWorkers(colony, p = 0.5, simParamBee = SP), simParamBee = SP), nWorkers(colony, simParamBee = SP))
- expect_equal(nDrones(replaceDrones(colony, p = 0, simParamBee = SP), simParamBee = SP), nDrones(colony, simParamBee = SP))
+ expect_warning(nDrones(replaceDrones(colony, p = 0, simParamBee = SP), simParamBee = SP))
virginQueensIDs <- getId(colony@virginQueens)
workerIDs <- getId(colony@workers)
droneIDs <- getId(colony@drones)
@@ -219,14 +217,14 @@ test_that("replaceFunctions", {
virginQueensIDs), 0)
expect_length(intersect(getId(replaceWorkers(colony, p = 0.5, simParamBee = SP)@workers),
workerIDs), nWorkers(colony, simParamBee = SP)/2)
- expect_length(intersect(getId(replaceDrones(colony, p = 0, simParamBee = SP)@drones),
- droneIDs), nDrones(colony, simParamBee = SP))
+ expect_warning(intersect(getId(replaceDrones(colony, p = 0, simParamBee = SP)@drones),
+ droneIDs))
expect_s4_class(replaceVirginQueens(apiary, simParamBee = SP), "MultiColony")
expect_s4_class(replaceWorkers(apiary, simParamBee = SP), "MultiColony")
expect_s4_class(replaceDrones(apiary, simParamBee = SP), "MultiColony")
expect_equal(nColonies(replaceVirginQueens(apiary, p = 1, simParamBee = SP)), nColonies(apiary))
expect_equal(nColonies(replaceWorkers(apiary, p = 0.5, simParamBee = SP)), nColonies(apiary))
- expect_equal(nColonies(replaceDrones(apiary, p = 0, simParamBee = SP)), nColonies(apiary))
+ expect_error(nColonies(replaceDrones(apiary, p = 0, simParamBee = SP)))
})
# ---- Remove functions ----
@@ -256,9 +254,9 @@ test_that("removeFunctions", {
expect_s4_class(removeVirginQueens(emptyColony, p = 0.5, simParamBee = SP), "Colony")
expect_s4_class(removeWorkers(emptyColony, p = 0, simParamBee = SP), "Colony")
expect_s4_class(removeDrones(emptyColony, p = 1, simParamBee = SP), "Colony")
- expect_s4_class(removeVirginQueens(emptyApiary, p = 0.5, simParamBee = SP), "MultiColony")
- expect_s4_class(removeWorkers(emptyApiary, p = 0, simParamBee = SP), "MultiColony")
- expect_s4_class(removeDrones(emptyApiary, p = 1, simParamBee = SP), "MultiColony")
+ expect_error(removeVirginQueens(emptyApiary, p = 0.5, simParamBee = SP))
+ expect_error(removeWorkers(emptyApiary, p = 0, simParamBee = SP))
+ expect_error(removeDrones(emptyApiary, p = 1, simParamBee = SP))
# Remove individuals in the non-empty colony/apiary
expect_s4_class(removeVirginQueens(colony, simParamBee = SP), "Colony")
@@ -308,9 +306,8 @@ test_that("setLocation", {
emptyApiary <- createMultiColony(n = 3, simParamBee = SP)
apiary <- createMultiColony(basePop[1:3], simParamBee = SP)
- expect_s4_class(setLocation(emptyApiary, location = c(1,2)), "MultiColony")
+ expect_error(setLocation(emptyApiary, location = c(1,2)))
expect_error(setLocation(emptyApiary, location = list(1,2))) # Lengths do not match
- expect_s4_class(setLocation(emptyApiary, location = list(1:2, 3:4, 4:5)), "MultiColony") #Not setting anything, if all are NULL!!!!
expect_s4_class(setLocation(apiary, location = c(1,2)), "MultiColony")
expect_s4_class(setLocation(apiary, location = list(1:2, 3:4, 4:5)), "MultiColony")
})