plotSegmentDistribution: Scatter plots of IBD segment distributions
In ibdsim2: Simulation of Chromosomal Regions Shared by Family Members

Description Usage Arguments Details Examples

View source: R/plotSegmentDistribution.R

Visualise and compare count/length distributions of IBD segments. Two types are currently implemented: Segments of autozygosity (for a single person) and segments with (pairwise) IBD state 1.

plotSegmentDistribution(
  ...,
  type = c("autozygosity", "ibd1"),
  ids = NULL,
  labels = NULL,
  col = NULL,
  shape = 1,
  alpha = 1,
  ellipses = TRUE,
  title = NULL,
  xlab = NULL,
  ylab = NULL,
  legendInside = TRUE
)

`...`	One or several objects of class `genomeSimList`, typically created by `ibdsim()`. They can be entered separately or as a `list`.
`type`	A string indicating which segments should be plotted. Currently, the allowed entries are "autozygosity" and "ibd1".
`ids`	A list of the same length as `...`, where each entry contains one or two ID labels (depending on `type`). By default (NULL), these labels are extracted from the inputs in `...`. Two other short-cuts are possible: If a single vector is given, it is repeated for all pedigrees. Finally, if `ids` is the word "leaves" then `pedtools::leaves()` is used to extract labels in each pedigree.
`labels`	An optional character vector of labels used in the legend. If NULL, the labels are taken from `names(...)`.
`col`	An optional colour vector of the same length as `...`.
`shape`	A vector with point shapes, of the same length as `...`.
`alpha`	A transparency parameter for the scatter points.
`ellipses`	A logical: Should confidence ellipses be added to the plot?
`title, xlab, ylab`	Title and axis labels.
`legendInside`	A logical controlling the legend placement.

This function takes as input one or several complete outputs from the ibdsim(), and produces a scatter plot of the number and average length of IBD segments from each.

Contour curves are added to plot, corresponding to the theoretical/pedigree-based values: either inbreeding coefficients (if type = "autozygosity") or κ_1 (if type = "ibd1").

# Simulation parameters used in the below examples.
map = uniformMap(M = 10)   # recombination map
N = 5                      # number of sims

# For more realistic results, replace with e.g.:
# map = loadMap("decode19")
# N = 1000


#################################################################
# EXAMPLE 1
# Comparison of IBD segment distributions
# between paternal and maternal half siblings.
#################################################################

# Define the pedigrees
xPat = halfSibPed()
xMat = swapSex(xPat, 1)

simPat = ibdsim(xPat, N = N, map = map)
simMat = ibdsim(xMat, N = N, map = map)

# By default, the IBD segments of the "leaves" are computed and plotted
plotSegmentDistribution(simPat, simMat, type = "ibd1", ids = 4:5,
                        labels = c("HSpat", "HSmat"))

#################################################################
# EXAMPLE 2
# Half siblings vs half uncle vs grandparent/grandchild
#################################################################

# Only one pedigree needed here
x = addSon(halfSibPed(), 5)

s = ibdsim(x, N = N, map = map)

# Indicate the pairs explicitly this time.
ids = list(HS = 4:5, HU = c(4,7), GR = c(1,7))

# List names are used as labels in the plot
plotSegmentDistribution(s, type = "ibd1", ids = ids, shape = 1:3)


#################################################################
# EXAMPLE 3
# Comparison of autozygosity distributions in various individuals
# with the same expected inbreeding coefficient (f = 1/8)
#################################################################

G = swapSex(linearPed(2), 5)           # grandfather/granddaughter
G = addChildren(G, 1, 5, 1)
HSpat = swapSex(halfSibPed(), 5)       # paternal half sibs
HSpat = addChildren(HSpat, 4, 5, 1)
HSmat = swapSex(HSpat, 1)              # maternal half sibs
QHFC = quadHalfFirstCousins()          # quad half first cousins
QHFC = addChildren(QHFC, 9, 10, nch = 1)

peds = list(G = G, HSpat = HSpat, HSmat = HSmat, QHFC = QHFC)
plotPedList(peds, newdev = TRUE)
dev.off()

# Simulations
s = lapply(peds, function(p)
  ibdsim(p, N = N, ids = leaves(p), verbose = FALSE, map = map))

# Plot distributions
plotSegmentDistribution(s, type = "autoz", title = "Autozygous segments")