paramlink
linkdat
object x
corresponding to a hypothesis H
. Relevant individuals unrelated to all others, are defined using
singleton
.The likelihood Pr(mixture,Typed contributors,Typed non-contributors|H)=P(R,T,V|H)
is calculated; the notation on the right hand side corresponds to that of Curran, Gill and Bill (2005).
A plot is also produced summarising the essential information.
Compared to previous literature and methods, including a series of papers by Fung and Hu, we generalise calculations to allow for
general, possibly inbred, pedigrees.
Typically calculations are performed for competing hypotheses and the ratio of likelihoods,
the likelihood ratio LR
is calculated and reported. Previous methods have assumed the relationships
between typed contributors to be same for the competing hypotheses. This restriction does not apply for
our approach. The calculation may also be used for identification cases where a mixture and reference samples are available.
Likelihood calculations are performed using the likelihood
of paramlink
.
The function checkInput
checks input to paraMix
.
paraMix(x, R, id.U, id.V = NULL, alleles, afreq = NULL,
Xchrom= FALSE, known_genotypes = list(), loop_breakers =NULL,
eliminate = 0, check = TRUE, plot = TRUE, title= NULL)
checkInput(x, R, id.U, id.V, alleles, all_typed, K, R_not_masked)
linkdat
object, or a list of such (if disconnected), describing the claimed relationship.
FALSE
for autosomal marker.
breakLoops
.
Pr(R,T,V|H)=Pr(R|T,V,H)Pr(T,V|H)= Pr(T,V|H)sum_u Pr(U=u,T,V|H)
where the sum extends over u among persons specified by id.U so that the union of u,T, V is R.
The likelihoohd for each u and the sum is returned.
Assumes alleles to be numbered 1,2,...
famMix
#Example 1: Motivating example Egeland et al. (2013)
require(paramlink)
y1=swapSex(nuclearPed(3),c(3,4))
p=c(0.1,0.2,0.3,0.4)
alleles=1:length(p)
T1=c(1,1)
T2=c(2,2)
R=1:2
known=list(c(3,T1),c(4,T2))
l1=paraMix(y1,R,id.U=5,alleles=alleles,afreq=p,known_genotypes=known)
y2=swapSex(nuclearPed(1),3)
y2=addOffspring(y2,mother=2,noff=1,sex=2)
y2=relabel(y2,c(1:3,6,4),1:5)
l2=paraMix(y2,R,id.U=6,alleles=alleles,afreq=p,known_genotypes=known)
LR1=l1$lik/l2$lik
exact=1/(2*(p[1]+p[2]))
stopifnot(abs(LR1-exact)<10^(-6))
#Example 2. Example 1 in Egeland et al. (2013) based on Fung and Hu (2008)
#Data:
#Mixture 1/2/3
#Suspect=4, genotype 3/3
#Victim=10, genotype 1/2
#H1: Contributors were the suspect and victim (unrelated)
#H2: Contributors were the father of suspect and victim (unrelated)
#H3: Contributors were the brother of suspect and victim (unrelated)
afreq=c(0.044,0.166,0.110,0.680)
alleles=1:length(afreq)
R=1:3 #Mixture
man_ped=nuclearPed(2)
victim = singleton(id=10, sex=2)
known = list(c(4,3,3),c(10,1,2)) #individual 4 is 3/3, and 10 (the victim) is 1/2.
#The likelihoods corresponding to H1,H2 and H3
l1=paraMix(list(man_ped, victim), R, id.U=NULL, id.V=NULL,
alleles=alleles, afreq=afreq, known_genotypes=known)$lik
l2=paraMix(list(man_ped, victim), R, id.U=1, id.V=4,
alleles=alleles, afreq=afreq, known_genotypes=known)$lik
l3=paraMix(list(man_ped, victim), R, id.U=3, id.V=4,
alleles=alleles, afreq=afreq, known_genotypes=known)$lik
LR12=l1/l2
stopifnot(abs(LR12-3.125)<10^(-6))
LR13=l1/l3
stopifnot(abs(LR13- 2.355296)<10^(-6))
Run the code above in your browser using DataLab