R/simulate_recombination.R

In sams25/rcombinator_old: Tools to Simulate the Evolution of Genetic Sequences

#' Simulates the evolution of a family of sequences by taking into
#' consideration both point mutation and burst events
#' @param sequences The original sequences to begin with
#' @param num_jumps The number of time steps to simulate for
#' @param transition_matrix The transition matrix for the point mutations
#' @param burst_threshhold The probability of a burst even occurring
#' @param burst_mean The mean (for a Poisson distribution) increase in copy
#' number during a burst event
#' @return A list of sequences that have evolved over time with both point
#' mutation and recombination
#' @export
simulate_recombination <- function(sequences,
                                   num_jumps,
                                   transition_matrix = get_k80_matrix(10^6),
                                   burst_threshold = 0.2,
                                   burst_mean = 1,
                                   flip_probability = 0.3,
                                   max_copy_number = 30,
                                   dist_method = "levenshtein",
                                   num_dist_matrices = 0,
                                   options = NULL,
                                   save_sequences = FALSE,
                                   save_file = NULL)
{
  if(options == "expensive")
  {
    sim_data <- simulate_recombination_expensive(sequences, num_jumps, transition_matrix,
                                                 burst_threshold, burst_mean, flip_probability)
  }
  else if(options == "cheap")
  {
    sim_data <- simulate_recombination_cheap(sequences, num_jumps, transition_matrix,
                                             burst_threshold, burst_mean, flip_probability,
                                             max_copy_number, dist_method)
  }
  else if(options == "summary")
  {
    sim_data <- simulate_recombination_summary(sequences, num_jumps, transition_matrix,
                                             burst_threshold, burst_mean, flip_probability,
                                             max_copy_number, dist_method)
  }
  else if(options == "families")
  {
    sim_data <- simulate_recombination_families(sequences, num_jumps, transition_matrix,
                                             burst_threshold, burst_mean, flip_probability,
                                             max_copy_number, dist_method, num_dist_matrices)
  }
  else
  {
    print("Invalid options passed onto function")
  }
  if(!is.null(save_file))
  {
    save(sim_data, file = paste(save_file, ".RData", sep = ""))
  }
  return(sim_data)
}