relationLR {paramlink} | R Documentation |
Relationship Likelihood Ratio
Description
Computes likelihood for two pedigrees and their ratio, the likelihood ratio (LR).
Usage
relationLR(
ped_numerator,
ped_denominator,
ids,
alleles,
afreq = NULL,
known_genotypes = list(),
loop_breakers = NULL,
Xchrom = FALSE,
plot = TRUE,
title1 = "",
title2 = ""
)
Arguments
ped_numerator |
a |
ped_denominator |
a |
ids |
genotyped individuals. |
alleles |
a numeric or character vector containing marker alleles names |
afreq |
a numerical vector with allele frequencies. An error is given if they don't sum to 1 (rounded to 3 decimals). |
known_genotypes |
list of triplets |
loop_breakers |
Not yet implemented, only default value NULL currently handled |
Xchrom |
a logical: Is the marker on the X chromosome? |
plot |
either a logical or the character 'plot_only', controlling if a plot should be produced. If 'plot_only', a plot is drawn, but no further computations are done. |
title1 |
a character, title of leftmost plot. |
title2 |
a character, title of rightmost plot. |
Details
This function computes the likelihood of two pedigrees (each corresponding to
a hypothesis describing a family relationship). The likelihood ratio is also
reported. Unlike other implementations we are aware of, partial DNA profiles
are allowed here. For instance, if the genotype of a person is reported as
1/0 (0 is 'missing') for a triallelic marker with uniform allele frequencies,
the possible ordered genotypes (1,1), (2,1), (1,2), (1,3) and (3,1) are
treated as equally likely. (For general allele frequencies, genotype
probabilities are obtained by assuming Hardy-Weinberg equilibrium.) A
reasonable future extension would be to allow the user to weigh these
genotypes; typically (1,1) may be more likely than the others. If
plot='plot_only'
, the function returns NULL after producing the plot.
Value
lik.numerator |
likelihood of data given ped_numerator |
lik.denominator |
likelihood of data given ped_denominator |
LR |
likelihood ratio lik.numerator/lik.denominator |
Author(s)
Thore Egeland, Magnus Dehli Vigeland
See Also
Examples
############################################
# A partial DNA profile is obtained from the person
# denoted 4 in the figure produced below
# There are two possibilities:
# H1: 4 is the missing relative of 3 and 6 (as shown to the left) or
# H2: 4 is unrelated to 3 and 6.
############################################
p = c(0.2, 0.8)
alleles = 1:length(p)
g3 = c(1,1); g4 = c(1,0); g6 = c(2,2)
x1 = nuclearPed(2)
x1 = addOffspring(x1, father = 4, sex = 1, noff = 1)
m = marker(x1, 3, g3, 4, g4, 6, g6, alleles = alleles, afreq = p)
x1 = addMarker(x1, m)
x2 = nuclearPed(2)
x2 = addOffspring(x2, father = 4, sex = 1, noff = 1)
m = marker(x2, 3, g3, 6, g6, alleles = alleles, afreq = p)
x2 = addMarker(x2, m)
missing = singleton(4, sex = 1)
m.miss = marker(missing, g4, alleles = alleles, afreq = p)
missing = addMarker(missing, m.miss)
x2 = relabel(x2, c(1:3, 99, 5:6), 1:6)
known = list(c(3, g3), c(4,g4), c(6, g6))
LR = relationLR(x1, list(x2, missing), ids = c(3,4,6),
alleles = alleles, afreq = p, known = known,
title1 = 'H1: Missing person 4 related',
title2 = 'H2:Missing person 4 unrelated')$LR
# Formula:
p = p[1]
LR.a = (1+p)/(2*p*(2-p))
stopifnot(abs(LR - LR.a) < 1e-10)