Final of this initial set of 3 quantitative trait simulations. The first reason I set out to do these simulations is that none of the figures in text books were very good (okay these aren’t pretty but at least I understand them), and also so that they could be used for demonstration for my grad. student class. This year I may try and get my big undergrad. class to download and run R (though this may not work).

Parents are generated from the population as before [1,2], with an additive genetic variance of 1, and a user specified environmental variance.

In this sim. there is a single generation of truncation selection, where individuals in the top SEL% of the phenotypic distribution are chosen to breed. The simulation produces 3 histograms of the parental distribution before and after selection (vertical bars show the mean) and the distribution of phenotypes in the offspring generation.

Ideas:
Height loci?
Allow user to specify alternative forms to selection (i.e. linear, quadratic).
Correlated traits
Multiple generations of selection, showing: movement of trait towards optimum; change in freq. of individual loci resulting in loss/fixation of genetic variation over generations; decline in heritability.

This slideshow requires JavaScript.

I usually show this one after the parent-offspring regression. If I were clever I should also could a version with the parent-child distributions as side-panels on the regression graph. But that would mean leaning how to draw a histogram on its side (sure it must a simple flag in par or some such thing).

one.gen.sel<-function(L=1000,environ.var,sel){
##Quantitative genetics sims
allele.freq<-0.5   ###each locus is assumed to have the same allele frequencies. This is just to simplify the coding, in reality these results work when each locus has its own frequency (and the coding wouldn't be too much harder). 
 

Num_inds=10000

##MAKE A MUM
## For each mother, at each locus we draw an allele (either 0 or 1) from the population allele frequency. 
##We do this twice for each mother two represent the two haplotypes in the mother 
mum.hap.1<-replicate(Num_inds, rbinom(L,1,allele.freq) )
mum.hap.2<-replicate(Num_inds, rbinom(L,1,allele.freq) )
##type mum.hap.1[,1] to see the 1st mothers 1st haplotype

##Each mothers genotype at each locus is either 0,1,2
mum.geno<-mum.hap.1+mum.hap.2

additive.genetic<-colSums(mum.geno)
mean.genetic<-mean(additive.genetic)
genetic.var<-sd(additive.genetic)

additive.genetic<-additive.genetic / sd(additive.genetic)
mum.pheno<- additive.genetic + rnorm(Num_inds,sd=sqrt(environ.var))
mum.pheno<-mum.pheno-mean(mum.pheno)



###FAMILIES


##MAKE A DAD (same code as make a mum, only said in a deeper voice)
dad.hap.1<-replicate(Num_inds, rbinom(L,1,allele.freq) )
dad.hap.2<-replicate(Num_inds, rbinom(L,1,allele.freq) )
dad.geno<-dad.hap.1+dad.hap.2


additive.genetic<-colSums(dad.geno)
additive.genetic<-additive.genetic / sd(additive.genetic)
dad.pheno<- additive.genetic + rnorm(Num_inds,sd=sqrt(environ.var))
dad.pheno<-dad.pheno-mean(dad.pheno)

### Make a child
child.geno<-dad.hap.1+mum.hap.1 ##1/2 from mum 1/2 from dad

additive.genetic<-colSums(child.geno)
additive.genetic<-additive.genetic / sd(additive.genetic)
child.pheno<- additive.genetic + rnorm(Num_inds,sd=sqrt(environ.var))
child.pheno<-child.pheno-mean(child.pheno)



 ##Selection of top sel% of individuals
 
top.sel.per.mums<- mum.pheno>quantile(mum.pheno,p=1-sel) 
top.sel.per.dads<- dad.pheno>quantile(dad.pheno,p=1-sel)
 
 
child.geno<-dad.hap.1[,top.sel.per.dads]+mum.hap.1[,top.sel.per.mums] ##1/2 from mum 1/2 from dad

additive.genetic<-(colSums(child.geno)-mean.genetic)
additive.genetic<-additive.genetic/genetic.var
child.pheno<- additive.genetic + rnorm(length(child.geno),sd=sqrt(environ.var))

layout(1:3)
my.lim<-quantile(c(mum.pheno,dad.pheno),p=c(0.01,0.99))
my.lim[2]<-quantile(child.pheno,p=c(0.99))

hist(c(mum.pheno,dad.pheno),breaks=100,xlim=my.lim,xlab="Phenotype",main=paste("Phenotype distribution before selection,Mean=0, Vg=1 Ve=",environ.var,"Taking top",round(100*sel),"%")); 
hist(c(mum.pheno[top.sel.per.mums],dad.pheno[top.sel.per.dads]),breaks=100,xlim=my.lim,xlab="Phenotype",main=paste("Phenotype distribution after selection, parental mean",format(mean(c(mum.pheno[top.sel.per.mums],dad.pheno[top.sel.per.dads])),dig=4))); 
abline(v= mean(c(mum.pheno[top.sel.per.mums],dad.pheno[top.sel.per.dads])),col="red",lwd=3)
hist(child.pheno,xlim=my.lim,breaks=100,xlab="Phenotype",main=paste("Phenotype distribution in the children Mean in children = ",format(mean(child.pheno),dig=4))); 

abline(v= mean(child.pheno),col="red",lwd=3)
##Mean phenotype after selection
cat("Selected parental mean",mean(c(mum.pheno[top.sel.per.mums],dad.pheno[top.sel.per.dads])),"\n")
##Mean child phenotype
cat("Mean in children = ",mean(child.pheno),"\n")
}
for(sel in c(0.1,0.4))
for(environ.var in c(.01,1,50)){
png(file=paste("One_generation_selection_varg_1_vare_",format(environ.var,dig=1),"sel_",format(sel,dig=1),".png",sep=""))
one.gen.sel(L=50,environ.var=environ.var,sel=sel) ##high herit.
dev.off()
}
 

This work is licensed under a Creative Commons Attribution 3.0 Unported License.