
Online app for IBDsimulations here:ibdsim2-shiny
As of version 2.1.0, the Shiny app is integrated as part of thepackage. To run it locally, simply run the command
ibdsim2::launchApp()The purpose ofibdsim2 is to simulate and analysethe gene flow in pedigrees. In particular, such simulations can be usedto study distributions of chromosomal segments sharedidentical-by-descent (IBD) by pedigree members. In eachmeiosis, the recombination process is simulated using sex specificrecombination rates in the human genome (Halldorsson et al.,2019), or with recombination maps provided by the user. Additionalfeatures include calculation of realised relatedness coefficients,distribution plots of IBD segments, and estimation of two-locusrelatedness coefficients.
ibdsim2 is part of thepedsuite collection ofpackages for pedigree analysis in R. A detailed presentation of thesepackages, including a separate chapter onibdsim2, isavailable in the bookPedigreeanalysis in R (Vigeland, 2021).
To getibdsim2, install from CRAN as follows:
install.packages("ibdsim2")Alternatively, install the latest development version fromGitHub:
# install.packages("remotes")remotes::install_github("magnusdv/ibdsim2")The most important function inibdsim2 isibdsim(), which simulates the recombination process in agiven pedigree. In this example we demonstrate this for in a familyquartet, and show how to visualise the result.
We start by loadingibdsim2.
library(ibdsim2)The main input toibdsim() is a pedigree and arecombination map. In our case we usepedtools::nuclearPed() to create the pedigree, and we loadchromosome 1 of the built-in map of human recombination.
# Pedigree with two siblingsx=nuclearPed(2)# Recombination mapchr1=loadMap("decode19",chrom =1)Now run the simulation! Theseed argument ensuresreproducibility.
sim=ibdsim(x,N =1,map = chr1,seed =1,verbose = F)The output ofibdsim() is a matrix (or a list ofmatrices, ifN > 1). Here are the first few rows of thesimulation we just made:
head(sim)#> chrom startMB endMB startCM endCM 1:p 1:m 2:p 2:m 3:p 3:m 4:p 4:m#> [1,] 1 0.000000 4.647215 0.000000 8.201478 1 2 3 4 1 3 2 3#> [2,] 1 4.647215 9.324570 8.201478 17.302931 1 2 3 4 1 4 2 3#> [3,] 1 9.324570 19.734471 17.302931 39.094957 1 2 3 4 2 4 2 3#> [4,] 1 19.734471 22.758411 39.094957 43.760232 1 2 3 4 2 4 1 3#> [5,] 1 22.758411 41.449745 43.760232 66.991502 1 2 3 4 2 4 1 4#> [6,] 1 41.449745 61.342551 66.991502 86.497704 1 2 3 4 2 4 1 3Each row of the matrix corresponds to a segment of the genome, anddescribes the allelic state of the pedigree in that segment. Eachindividual has two columns, one with the paternal allele (marked by thesuffix “:p”) and one with the maternal (suffix “:m”). The founders (theparents in our case) are assigned alleles 1, 2, 3 and 4.
The functionhaploDraw() interprets the founder allelesas colours and draws the resulting haplotypes onto the pedigree. See?haploDraw for an explanation ofpos and otherarguments.
haploDraw(x, sim,pos =c(2,4,2,4))
In this example we will compare the distributions of counts/lengthsof IBD segments between the following pairwise relationships:
Note that GR and HS have the same relatedness coefficientskappa = (1/2, 1/2, 0), meaning that they are geneticallyindistinguishable in the context of unlinked loci. In contrast, HU haskappa = (3/4, 1/4, 0).
For simplicity we create a pedigree containing all the threerelationships we are interested in.
x=halfSibPed()|>addSon(5)plot(x)
We store the ID labels of the three relationships in a list.
ids=list(GR =c(2,7),HS =4:5,HU =c(4,7))Next, we useibdsim() to produce 500 simulations of theunderlying IBD pattern in the entire pedigree.
s=ibdsim(x,N =500,map ="decode19",seed =1234)#> Simulation parameters:#> Simulations : 500#> Chromosomes : 1-22#> Genome length: 2875 Mb#> 2602.29 cM (male)#> 4180.42 cM (female)#> Recomb model : chi#> Target indivs: 1-7#> Skip recomb : -#> Total time used: 1.71 secsTheplotSegmentDistribution() function, with the optiontype = "ibd1" analyses the IBD segments in each simulation,and makes a nice plot. Note that the names of theids listare used in the legend.
plotSegmentDistribution(s,type ="ibd1",ids = ids,shape =1:3)
We conclude that the three distributions are almost completelydisjoint. Hence the three relationships can typically be distinguishedon the basis of their IBD segments, if these can be determinedaccurately enough.
A Shiny app for visualising IBD distributions is available here:https://magnusdv.shinyapps.io/ibdsim2-shiny/.