evolve: Evolution of genetic and genomic landscapes

Description Usage Arguments Details Value References See Also Examples

View source: R/glads_funcs.R

Description

This function simulates the evolution of individual-based populations forward in time.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
evolve(
  x,
  time,
  type = c("constant", "dynamic", "additive", "custom"),
  recombination = c("map", "average"),
  recom.rate,
  loci.pos = NULL,
  chromo_mb = NULL,
  init.sex = NULL,
  migration.rate = NULL,
  mutation.rate = NULL,
  param.z = NULL,
  param.w = NULL,
  fun = c(phenotype = NULL, fitness = NULL)
)

Arguments

x

List of objects of class "struct" with the initial genetic structure of populations.

time

Number of simulated generations.

type

Type of simulated evolution. Current options are 'constant', 'dynamic', 'additive', and 'custom' (see Details).

recombination

Type of recombination between homologous chromosomes. Current options are 'map' and 'average' (see Details).

recom.rate

A numerical value for the recombination type 'average', or a vector of nl - 1 elements, with the recombination rate between neighbouring loci for type 'map' (nl: number of loci) (see Details).

loci.pos

A vector with the position of loci. The default value is NULL, but it is required for the recombination type 'average'.

chromo_mb

A numerical value with the size of the simulated chromosome (in megabase). The default value is NULL, but it is required for the recombination type 'average'.

init.sex

A list of vectors defining the sex of the initial individuals in the populations (1: females and 2: males). The default value is NULL and assigns the sex of the first generation randomly.

migration.rate

A single value setting the migration rate between all populations, or a squared matrix of order equal to the number of populations, with the migration rate between them (migration['from', 'to']). It is ignored for single population simulations. The default value is NULL with no migration between populations.

mutation.rate

A numerical value setting the mutation rate per site. It is currently restricted to biallelic SNPs (genetic structures with values 1 or 2). The default value is NULL.

param.z

A list with the parameter values for the function computing phenotypes. The list of parameters for each population should be included as a list of list (see Example)

param.w

A list with the parameter values for the fitness function. The list of parameters for each population should be included as a list of list (see Example)

fun

A character vector with the names of the custom phenotype and fitness functions. The default value is NULL, but it is required for the 'custom' type of evolution

Details

This function returns a list of populations composed of two-dimensional arrays that represents a pair of homologous chromosomes. Rows represent individuals and columns the different loci. Each element of the array is an integer defining the copy of a given allele at a given locus.

Different types of evolution are available for simulations:

The recombination between homologous chromosomes are either of type 'map' or 'average'. The first case needs a vector with the recombination rate (ρ) between neigbour loci of length equal to nl - 1 (nl: number of loci). The probability of having a crossover (1) or not (0) is uniformly distributed at a rate defined by the value of ρ between loci (i.e. positions with a probability smaller than ρ recombine). The uniform distribution allows each position with the same values of ρ to have an equal chance of crossover across all iterations. There is no recombination between homologous chromosomes when ρ = 0, both loci are completely linked (e.g. within an inversion or situated close to centromeres), while with a value of ρ = 0.5, the recombination rate is completely random (i.e. both loci are very distant on the same chromosome or are located on different chromosomes). A value of ρ < 0.5 means the loci are physically linked. The second case, when recombination is of type 'average', a numerical value with the average recombination rate per base pair should be supplied, the loci position and the size of the chromosome in megabase are also required. The crossover points are exponentially distributed as a Poisson process (see Example).

Value

A list of objects of class "struct" or array.

References

Quilodrán, C. S., Ruegg, K., Sendell-Price, A. T., Anderson, E., Coulson, T. and Clegg, S. (2020). The many population genetic and demographic routes to islands of genomic divergence. Methods in Ecology and Evolution 11(1):6-21.. doi: 10.1111/2041-210X.13324.

See Also

initial.struct

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
## Not run: 
## We start with a population of 20 individuals and 10 biallelic SNPs
initial.population.size <- 20
n.loci <- 10
n.alleles.per.locus <- 2

start1 <- initial.struct(initial.population.size,n.loci,n.alleles.per.locus)

##################################
# Type of evolution = "constant" #
##################################

## We set a recombination map and the number of generations to be simulated
recom.map <- rep(0.5, n.loci-1) #all loci are independent
n.gens <- 10
pop <- evolve(list(start1), n.gens, "constant", "map", recom.map)

## A similar simulation but with recombination type "average".
## We need to specify the position of loci and the size of the chromosome (MB)
loci.pos<- sample(80000: 12000000, 10) #random loci positions
chromo_mb<-12000000 	 #chromosome size
crossover <- 1/100000000.0 #average recombination rate (1 cM/MB)
pop <- evolve(list(start1), n.gens, "constant", "average", crossover, loci.pos, chromo_mb)

# We include a mutation rate for the simulation of biallelic loci
pop <- evolve(list(start1), n.gens, "constant", "average", crossover,
              loci.pos, chromo_mb, mutation.rate = 0.0001)

# A second population is included to incorporate the effect of migration
start2 <- initial.struct(initial.population.size,n.loci,n.alleles.per.locus)
pop <- evolve(list(start1, start2), n.gens, "constant", "average", crossover,
              loci.pos, chromo_mb, migration.rate = 0.01)

# The sex of individuals may be set for the starting generations.
sex.start1 <- sample(1:2, initial.population.size, replace=T)
sex.start2 <- sample(1:2, initial.population.size, replace=T)
init.sex <- list(sex.start1, sex.start2)
pop <- evolve(list(start1, start2), n.gens, "constant", "average", crossover,
              loci.pos, chromo_mb, init.sex = init.sex)


##################################
# Type of evolution = "dynamic"  #
##################################

## We set the parameters for the computation of phenotypes and for the
## fitness function of the type 'dynamic'
sex.ratio <- 0.5
mean.fitness <- 3 #mean number of offsprings per breeding pair
d.d <- 0.01 #density-dependent demographic effect

set.seed(1) #setting the seed for reproducible random numbers
pop <- evolve(list(start1), n.gens, "dynamic", "map", recom.map, param.z = list(list(sex.ratio)),
              param.w = list(list(mean.fitness, d.d)))

##################################
# Type of evolution = "additive" #
##################################

## We set the parameters for the computation of phenotypes and for
## the fitness function of the type 'additive'
sex.ratio <- 0.5
fitness.pos <- 1:n.loci #all simulated loci are additive
add.loci <- n.loci

# next line is breeding value of additive loci
bvs <- t(array( seq(0,1, length = n.alleles.per.locus) ,c(n.alleles.per.locus, add.loci)))

e.v <- 0.01 					# stochastic environmental variant
param.z <- list(sex.ratio, fitness.pos, bvs, add.loci, e.v)

b0 <- 6 # maximum number of offspring
b1 <- 1 # phenotypic optima
b2 <- 0.75 # variance of the Gaussian curve
b3 <- 0.01 # intensity of the density-dependence
d.v = 1 # stochastic demographic variant
param.w <- list(b0,b1,b2,b3, d.v, add.loci)

set.seed(1) #setting the seed for reproducible random numbers
pop<-evolve(list(start1), n.gens, "additive", "map", recom.map,
            param.z = list(param.z), param.w = list(param.w))

##################################
#  Type of evolution = "custom"  #
##################################

## We set a custom 'phenotype' and 'fitness' function.
## In this example, the custom functions are very similar to the default 'additive' ones.

phenotype2 <- function(struct, ...){
 pop.struct <- struct[ , fitness.pos, ]
 temp <- dim(pop.struct)
 mat <- matrix(1:temp[2],temp[1],temp[2],byrow=TRUE)

 # the next line gives an array that gives the locus index at each position in pop.struct
 loci.n <- array(mat,c(temp[1],temp[2],2))

 new.n <- (pop.struct-1)*add.loci+loci.n
 bvv <- as.vector(bvs) # turn to bvs
 outp <- array(bvv[new.n],c(temp[1],temp[2],2))
 z <- apply(outp,1,sum)-add.loci+rnorm(temp[1],0,e.v)
 sex <- rbinom(nrow(struct), 1, sex.ratio)+1
 return(cbind(z = z, sex = sex))
}

fitness2 <- function(z, sex, n, ...){
 ng <- n.loci
 a <- b0
 b <- b1
 c <- b2
 w <- round( a*exp(-((z-b*ng)^2)/(2*(c*ng)^2) ) - b3*n + rnorm(n,0,d.v) , 0)
 w <- ifelse(w<0,0,w)
 return(w)
}

set.seed(1) #setting the seed for reproducible random numbers
pop<-evolve(list(start1), n.gens, "custom", "map", recom.map,
           param.z =list(sex.ratio = sex.ratio, fitness.pos = fitness.pos,
                         bvs = bvs, add.loci = add.loci, e.v = e.v),
           param.w = list(b0 = b0, b1 = b1, b2 = b2, b3 = b3, d.v = d.v, n.loci = add.loci),
           fun=c("phenotype2", "fitness2"))


## End(Not run)

eriqande/glads documentation built on April 24, 2020, 11:53 a.m.