plotSegmentDistribution {ibdsim2}

Scatter plots of IBD segment distributions

Description

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.

Usage

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

Arguments

...

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.

unit

Length unit; either "cm" (centiMorgan) or "mb" (megabases).

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.

Details

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 \kappa_1 (if type = "ibd1").

Examples

# 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(GR = c(2,7), HS = 4:5, HU = c(4,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")