Extending The Theory of Junctions

Description

The theory of junctions is extended by this package by including the effect of a finite number of recombination sites along the chromosome. The package provides functions to calculate the estimated number of junctions, depending on the time since the onset of hybridization, population size, number of recombination sites, initial heterozygosity and the number of crossovers per meiosis.

Details

This package provides individual based simulations in order to simulate the accumulation of junctions over time, both for chromosomes with a finite and an infinite number of recombination sites. Furthermore, the package provides mathematical tools to verify the outcomes of the individual based simulations.

Update version 2.1.1 : fixed typo in one of the examples

Update version 2.1.0 : updated tbb::task_scheduler_init to tbb::global_control

Update version 2.0.2 : simplified some tests

Update version 2.0 : merged many functions with similar functionality, added vignette that provides overview of all functionality.

Update version 1.9 : added c++ versions of the unphased and phased likelihoods.

Update version 1.8 : added multithreading using the TBB library.

Update version 1.7 : further improved the recombination function following Hanno Hildenbrandt's suggestions

Update version 1.6 : improved the recombination function to run twice as fast

Update version 1.5.1: added option to track the true number of junctions

Update version 1.5: added support for inferring the time since admixture based on phased and unphased data. Also included are simulation functions to simulate appropriate data (e.g. phased and unphased).

Update version 1.4: added support for estimating the number of junctions, and simulating the number of junctions, under a backcrossing scheme, using the code supplied in Lavretsky et al. 2019.

Update version 1.3: added support for estimating the time since admixture using unphased data.
Update version 1.3: added individual based simulations returning phased and unphased data.
Update version 1.3: Updated entire package to Roxygen.

Update version 1.2: added support for estimating the expected number of junctions for arbitrarily distributed markers.

Update version 1.1: updated underlying random number generator for picking recombination sites. The previous generator had limited precision, which could generate duplicate recombination sites. This update fixes that

Author(s)

Maintainer: Thijs Janzen <thijsjanzen@gmail.com>

References

Janzen, T. , Nolte, A. W. and Traulsen, A. (2018), The breakdown of genomic ancestry blocks in hybrid lineages given a finite number of recombination sites. Evolution, 72: 735-750. doi:10.1111/evo.13436

Lavretsky, P, Janzen, T. and McCracken, KG. (2019) Identifying hybrids & the genomics of hybridization: Mallards & American black ducks of Eastern North America. Ecology and Evolution 9: 3470-3490. doi:10.1002/ece3.4981

Calculate the limit of the number of junctions

Description

Calculate the average number of junctions after an infinite number of generations, provided information on the initial heterozygosity, population size and the number of generations.

Usage

calc_k(N = Inf, R = Inf, H_0 = 0.5, C = 1)

Arguments

Value

The number of junctions for at time = infinity

Examples

k <-  calc_k(N = 100, R = 1000, H_0 = 0.5, C = 1)

Function to calculate the maximum accurate time

Description

Function that calculates the maximum time after hybridization after which the number of junctions can still be reliably used to estimate the onset of hybridization. This is following equation 15 in Janzen et al. 2018.

Usage

calculate_mat(N = Inf, R = Inf, H_0 = 0.5, C = 1)

Arguments

Value

The maximum accurate time

Examples

calculate_mat(N = Inf, R = 1000, H_0 = 0.5, C = 1)

Estimate the time since the onset of hybridization, using the number of junctions

Description

Estimate the time since the onset of hybridization, following equation 14 in Janzen et al. 2018

Usage

estimate_time(J = NA, N = Inf, R = Inf, H_0 = 0.5, C = 1)

Arguments

Value

The number of generations passed since the onset of hybridization

Examples

J <- number_of_junctions(N = 100, R = 1000, H_0 = 0.5, C = 1, t = 200)
estimate_time(J = J, N = 100, R = 1000, H_0 = 0.5, C = 1)
# should be 200 again

estimates the time since admixture, given diploid ancestry data.

Description

Calculates the time since admixture, given unphased ancestry data.

Usage

estimate_time_diploid(
  ancestry_information,
  analysis_type = "individuals",
  phased = FALSE,
  pop_size = 1000,
  freq_ancestor_1 = 0.5,
  lower_lim = 2,
  upper_lim = 2000,
  num_threads = 1,
  verbose = FALSE
)

Arguments

estimate time using likelihood for a single chromosome

Description

Estimate the time since the onset of hybridization, for a haploid genome

Usage

estimate_time_haploid(
  ancestry_matrix,
  N = 1000,
  freq_ancestor_1 = 0.5,
  lower_lim = 2,
  upper_lim = 1000,
  verbose = FALSE
)

Arguments

Value

The number of generations passed since the onset of hybridization

Estimate the time since the onset of hybridization, using the observed number of junctions, taking into account the distribution of markers on a single chromosome

Description

Estimate the time since the onset of hybridization, following equation 1 in Janzen et al. unpublished

Usage

estimate_time_one_chrom(
  J = NA,
  N = Inf,
  H_0 = 0.5,
  marker_distribution = NA,
  lower_lim = 2,
  upper_lim = 1000
)

Arguments

Value

The number of generations passed since the onset of hybridization

Examples

markers <- seq(from = 0, to = 1, length.out = 100)
J <- number_of_junctions_markers(N = 100, H_0 = 0.5, t = 200,
marker_distribution = markers)
estimate_time_one_chrom(J = J,
                        N = 100,
                        H_0 = 0.5,
                        marker_distribution = markers) #should be 200 again

calculate the log likelihood of observing diploid ancestry data.

Description

Calculates the log likelihood of observing the phased data, given the population size, initial heterozygosity and time since admixture

Usage

log_likelihood_diploid(
  local_anc_matrix,
  pop_size,
  freq_ancestor_1 = 0.5,
  t,
  phased = FALSE,
  num_threads = 1
)

Arguments

Value

log likelihood

log likelihood of the time since admixture for a haploid genome

Description

log likelihood of the time since admixture for a set of single chromosomes (for ex. in Yeast).

Usage

log_likelihood_haploid(ancestry_matrix, N = 1000, freq_ancestor_1 = 0.5, t = 2)

Arguments

Value

loglikelihood

Calculate the average number of junctions

Description

Calculate the average number of junctions in a single chromosome after t generations, provided information on the initial heterozygosity, population size and the number of generations.

Usage

number_of_junctions(N = Inf, R = Inf, H_0 = 0.5, C = 1, t = 100)

Arguments

Value

Estimated number of junctions at time t

Examples

jt <-  number_of_junctions(N = 100, R = 1000, H_0 = 0.5, C = 1, t = 1000)
jt2 <- number_of_junctions(N = 100, R = 1000, H_0 = 0.5, C = 1, t = 0:1000)

Calculate the average number of junctions during backcrossing

Description

Calculate the expected number of junctions after t generations, in a backcrossing mating scheme.

Usage

number_of_junctions_backcross(H_0 = 0.5, C = 1, t = 100)

Arguments

Value

Estimated number of junctions at time t

Examples

jt <-  number_of_junctions_backcross(H_0 = 0.1, C = 1, t = 5)

Calculate the expected number of junctions between two markers separated by a given amount of recombination

Description

Calculate the expected number of junctions after t generations, provided information on the initial heterozygosity, population size, the number of generations since the onset of admixture and the distance between two markers.

Usage

number_of_junctions_di(N = Inf, H_0 = 0.5, t = 100, di = 1e-06)

Arguments

Value

Estimated number of junctions at time t

Examples

number_of_junctions_di(N = 100, H_0 = 0.5, t = 1000, di = 0.01)

Calculate the expected total number of junctions in a chromosome, given the distribution of markers

Description

Calculate the expected number of junctions after t generations, provided information on the initial heterozygosity, population size, the number of generations since the onset of admixture and the distribution of markers.

Usage

number_of_junctions_markers(
  N = Inf,
  H_0 = 0.5,
  t = 100,
  marker_distribution = NA
)

Arguments

Value

Estimated number of observed junctions at time t

Examples

markers <- seq(from = 0, to = 1, length.out = 1000)
jt <-  number_of_junctions_markers(N = 100,
                                  H_0 = 0.5,
                                  t = 1000,
                                  marker_distribution = markers)
random_markers <- sort(runif(1000, 0, 1))
jt2 <- number_of_junctions_markers(N = 100,
                                  H_0 = 0.5,
                                  t = 1000,
                                  marker_distribution = random_markers)

Function to simulate data using a back crossing scheme

Description

Individual based simulation of the accumulation of junctions, under a back crossing scheme

Usage

sim_backcrossing(
  population_size = 100,
  freq_ancestor_1 = 0.5,
  total_runtime = 5,
  size_in_morgan = 1,
  number_of_markers = 100,
  seed = 6,
  time_points = -1
)

Arguments

Value

List with five entries: average_junctions: average number of junctions over time, detected_junctions: average number of detected junctions, given the markers. markers: vector with the locations of the molecular markers, junction_distribution: distribution of junctions per time step average_heterozygosity: average heterozygosity.

Examples

sim_backcrossing(population_size = 100,
                       total_runtime = 5,
                       size_in_morgan = 1,
                       number_of_markers = 100,
                       seed = 6,
                       time_points = 1:5)

Individual Based Simulation of the accumulation of junctions

Description

Individual based simulation of the accumulation of junctions for a chromosome with regularly distributed markers.

Usage

sim_fin_chrom(
  pop_size = 100,
  freq_ancestor_1 = 0.5,
  total_runtime = 100,
  morgan = 1,
  seed = 42,
  R = 100
)

Arguments

Value

Examples

sim_fin_chrom(pop_size = 100, freq_ancestor_1 = 0.5,
                   total_runtime = 10, morgan = 1, seed = 42,
                   R = 100)

Individual Based Simulation of the accumulation of junctions

Description

Individual based simulation of the accumulation of junctions for a chromosome with an infinite number of recombination sites.

Usage

sim_inf_chrom(
  pop_size = 100,
  freq_ancestor_1 = 0.5,
  total_runtime = 100,
  morgan = 1,
  markers = -1,
  seed = 42
)

Arguments

Value

Examples

v <- sim_inf_chrom(pop_size = 100, freq_ancestor_1 = 0.5,
                   total_runtime = 10, morgan = 1, markers = 100,
                   seed = 42)
plot(v$avgJunctions, type = "l", xlab = "Generations",
ylab = "Number of Junctions", main = "Example Infinite Chromosome")
lines(v$detectedJunctions, col = "blue")
legend("bottomright", c("Real number","Number detected"),
       lty = 1, col = c("black", "blue"))

Individual Based Simulation of the accumulation of junctions

Description

Individual based simulation of the accumulation of junctions, returning phased and unphased data. Ancestry on both chromosomes of 10 randomly sampled individuals per generations is returned.

Usage

sim_phased_unphased(
  pop_size = 100,
  freq_ancestor_1 = 0.5,
  total_runtime = 100,
  size_in_morgan = 1,
  markers = 100,
  time_points = -1,
  num_threads = 1,
  verbose = FALSE,
  record_true_junctions = FALSE,
  num_indiv_sampled = 10,
  coverage = 1,
  error_rate = 0
)

Arguments

Value

a tibble with five columns: [time, individual, marker location, ancestry chromosome 1, ancestry chromosome 2]

Examples

## Not run: 
sim_phased_unphased(pop_size = 100, freq_ancestor_1 = 0.5,
                    total_runtime = 10, size_in_morgan = 1,
                    markers = 10, time_points = c(0, 5, 10),
                    num_threads = 1)

## End(Not run)

Estimate the error in the time estimate

Description

Calculate the error in the estimate of the onset of hybridization, following Equations 3 & 4 in the Supplementary information of Janzen et al. 2018.

Usage

time_error(t = NA, N = Inf, R = Inf, H_0 = 0.5, C = 1, relative = TRUE)

Arguments

Value

Expected error in the time estimate