R/simul.pedigree.r

In synbreed: Framework for the Analysis of Genomic Prediction Data using R

#' Simulation of pedigree structure
#'
#' This function can be used to simulate a pedigree for a given number of
#' generations and individuals. Function assumes random mating within
#' generations. Inbred individuals may be generated by chance.
#'
#' If \code{animals=FALSE}, the parents for the current generation will be
#' randomly chosen out of the genotypes in the last generation. If Par1 = Par2,
#' an inbreed is generated. If \code{animal=TRUE}, each ID is either sire or
#' dam. Each ID is progeny of one sire and one dam.
#'
#' @param generations \code{integer}. Number of generations to simulate
#' @param ids \code{integer} or vector of integers. Number of genotypes in each
#' generation. If length equal one, the same number will be replicated and used
#' for each generation.
#' @param animals \code{logical}. Should a pedigree for animals be simulated
#' (no inbreeding)? See 'Details'.
#' @param familySize \code{numeric}. Number of individuals in each full-sib
#' family in the last generation.
#' @return An object of class \code{pedigree} with N=sum(ids) genotypes.
#' @author Valentin Wimmer
#' @seealso \code{\link{simul.phenotype}}, \code{\link{create.pedigree}},
#' \code{\link{plot.pedigree}}
#' @examples
#'
#' # example for plants
#' ped <- simul.pedigree(gener = 4, ids = c(3, 5, 8, 8))
#' plot(ped)
#' # example for animals
#' peda <- simul.pedigree(gener = 4, ids = c(3, 5, 8, 8), animals = TRUE)
#' plot(peda)
#' @export simul.pedigree
simul.pedigree <- function(generations = 2, ids = 4, animals = FALSE, familySize = 1) {

  # if only one value is given, set samenumber for all generations
  if (length(ids) == 1) ids <- rep(ids, times = generations)

  # initialisation
  gener <- rep(1:generations, times = ids)
  ID <- 1:sum(ids)
  Par1 <- Par2 <- rep(0, length(ID))

  # random mating for plants (inbreeds are likely)
  if (!animals) {
    for (i in 2:generations) {
      Par1[gener == i] <- ID[sample(ID[gener == i - 1], size = ids[i], replace = TRUE)]
      Par2[gener == i] <- ID[sample(ID[gener == i - 1], size = ids[i], replace = TRUE)]
    }
    ped <- data.frame(ID = ID, Par1 = Par1, Par2 = Par2, gener = gener - 1, sex = NA)
  }
  # define sire and dams for animals (no inbreeds)
  else {
    # define sex for 1st generation
    # 0 = female
    # 1 = male
    sex <- rep(0, length(ID))
    sex[gener == 1] <- sample(rep(c(0, 1), length = sum(gener == 1)), sum(gener == 1), replace = FALSE)

    for (i in 2:generations) {
      sex[gener == i] <- sample(rep(c(0, 1), length = sum(gener == i)), sum(gener == i), replace = FALSE)

      Par1[gener == i] <- ID[sample(ID[(gener == i - 1) & sex == 1], size = ids[i], replace = TRUE)]
      Par2[gener == i] <- ID[sample(ID[(gener == i - 1) & sex == 0], size = ids[i], replace = TRUE)]
    }
    ped <- data.frame(ID = ID, Par1 = Par1, Par2 = Par2, gener = gener - 1, sex = sex)
  }

  # create full-sib families in the last generation
  pedTemp <- ped[ped$gener == generations - 1, ]
  ped <- ped[ped$gener < generations - 1, ]
  pedFamily <- data.frame(
    ID = rep(1:(nrow(pedTemp) * familySize) + pedTemp$ID[1] - 1),
    Par1 = rep(pedTemp$Par1, each = familySize),
    Par2 = rep(pedTemp$Par2, each = familySize),
    gener = rep(pedTemp$gener, each = familySize),
    sex = rep(pedTemp$sex, each = familySize)
  )
  ped <- rbind(ped, pedFamily)
  ped$ID <- paste("ID", ped$ID, sep = "")
  ped$Par1[ped$gener != 0] <- paste("ID", ped$Par1[ped$gener != 0], sep = "")
  ped$Par2[ped$gener != 0] <- paste("ID", ped$Par2[ped$gener != 0], sep = "")

  class(ped) <- c("pedigree", "data.frame")
  return(ped)
}

simul.pedigree()