View source: R/exclusionPower.R
exclusionPower | R Documentation |
Computes the power (of a single marker, or for a collection of markers) of excluding a claimed relationship, given the true relationship.
exclusionPower(
claimPed,
truePed,
ids,
markers = NULL,
source = "claim",
disableMutations = NA,
exactMaxL = Inf,
nsim = 1000,
seed = NULL,
alleles = NULL,
afreq = NULL,
knownGenotypes = NULL,
Xchrom = FALSE,
plot = FALSE,
plotMarkers = NULL,
verbose = TRUE
)
claimPed |
A |
truePed |
A |
ids |
Individuals available for genotyping. |
markers |
A vector indicating the names or indices of markers attached
to the source pedigree. If NULL (default), then all markers attached to the
source pedigree are used. If |
source |
Either "claim" (default) or "true", deciding which pedigree is used as source for marker data. |
disableMutations |
This parameter determines how mutation models are treated. Possible values are as follows:
|
exactMaxL |
A positive integer, or |
nsim |
A positive integer; the number of simulations used for markers
whose number of alleles exceeds |
seed |
An integer seed for the random number generator (optional). |
alleles , afreq , Xchrom |
If these are given, they are used (together with
|
knownGenotypes |
A list of triplets |
plot |
Either a logical or the character "plotOnly". If the latter, a plot is drawn, but no further computations are done. |
plotMarkers |
A vector of marker names or indices whose genotypes are to be included in the plot. |
verbose |
A logical. |
This function implements the formula for exclusion power as defined and discussed in (Egeland et al., 2014).
It should be noted that claimPed
and truePed
may be any (lists of)
pedigrees, as long as they both contain the individuals specified by ids
.
In particular, either alternative may have inbred founders (with the same or
different coefficients), but this must be set individually for each.
If plot = "plotOnly"
, the function returns NULL after producing the
plot.
Otherwise, the function returns an EPresult
object, which is essentially
a list with the following entries:
EPperMarker
: A numeric vector containing the exclusion power of each
marker. If the known genotypes of a marker are incompatible with the true
pedigree, the corresponding entry is NA
.
EPtotal
: The total exclusion power, computed as 1 - prod(1 - EPperMarker, na.rm = TRUE)
.
expectedMismatch
: The expected number of markers giving exclusion,
computed as sum(EPperMarker, na.rm = TRUE)
.
distribMismatch
: The probability distribution of the number of markers
giving exclusion. This is given as a numeric vector of length n+1
, where
n
is the number of nonzero elements of EPperMarker
. The vector has
names 0:n
.
time
: The total computation time.
params
: A list containing the (processed) parameters ids
, markers
and disableMutations
.
Magnus Dehli Vigeland
T. Egeland, N. Pinto and M.D. Vigeland, A general approach to power calculation for relationship testing. Forensic Science International: Genetics 9 (2014): 186-190. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1016/j.fsigen.2013.05.001")}
############################################
### A standard case paternity case:
### Compute the power of exclusion when the claimed father is in fact
### unrelated to the child.
############################################
# Claim: 'AF' is the father of 'CH'
claim = nuclearPed(father = "AF", children = "CH")
# Attach two (empty) markers
claim = claim |>
addMarker(alleles = 1:2) |>
addMarker(alleles = 1:3)
# Truth: 'AF' and 'CH' are unrelated
true = singletons(c("AF", "CH"))
# EP when both are available for genotyping
exclusionPower(claim, true, ids = c("AF", "CH"))
# EP when the child is typed; homozygous 1/1 at both markers
claim2 = claim |>
setGenotype(marker = 1:2, id = "CH", geno = "1/1")
exclusionPower(claim2, true, ids = "AF")
############################################
### Two females claim to be mother and daughter, but are in reality sisters.
### We compute the power of various markers to reject the claim.
############################################
ids = c("A", "B")
claim = nuclearPed(father = "NN", mother = "A", children = "B", sex = 2)
true = nuclearPed(children = ids, sex = 2)
# SNP with MAF = 0.1:
PE1 = exclusionPower(claimPed = claim, truePed = true, ids = ids,
alleles = 1:2, afreq = c(0.9, 0.1))
stopifnot(round(PE1$EPtotal, 5) == 0.00405)
# Tetra-allelic marker with one major allele:
PE2 = exclusionPower(claimPed = claim, truePed = true, ids = ids,
alleles = 1:4, afreq = c(0.7, 0.1, 0.1, 0.1))
stopifnot(round(PE2$EPtotal, 5) == 0.03090)
### How does the power change if the true pedigree is inbred?
trueLOOP = halfSibPed(sex2 = 2) |> addChildren(4, 5, ids = ids)
# SNP with MAF = 0.1:
PE3 = exclusionPower(claimPed = claim, truePed = trueLOOP, ids = ids,
alleles = 1:2, afreq = c(0.9, 0.1))
# Power almost doubled compared with PE1
stopifnot(round(PE3$EPtotal, 5) == 0.00765)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.