Description Usage Arguments Details Value Author(s) References See Also Examples
Marker data is simulated for a specified person and several markers using paramlinkConditionalOne
for each marker.
The only difference between FamiliasConditional
and paramlinkConditional
is that the former uses
the C implementation of FamiliasPosterior
for likelihood calculation while the latter is based on the R
implementation in paramlink
.
1 2 3 4 | paramlinkConditional(Nsim = 5, datamatrix, loci, pedigrees,
truePed = 1, ref=NULL, available = NULL, prior=NULL, seed = NULL)
paramlinkConditionalOne(Nsim = 5, mark = 1, ref = 2, datamatrix, loci,
pedigrees, truePed = 1, available = NULL, prior = NULL, seed = NULL)
|
Nsim |
Integer. Number of simulations. |
ref |
Integer Denominator of LR. |
datamatrix |
A data frame or a matrix. The row names must be the names of the persons you have data for. The columns contain the alleles, two columns for each marker, in the same order used in the loci list. |
loci |
A list of |
pedigrees |
An list with elements of type |
truePed |
Integer. Index of pedigree from which marker data are simulated. |
available |
A character giving the name of the person to be simulated or the integer ID. |
prior |
Double vector. The prior on pedigrees. |
seed |
Integer used to fix simulations. |
mark |
Integer. Index of |
Marker data is simulated for a specified person and one specified marker
using paramlinkConditional
which calls markerSim
.
The marker data is then loaded into a datamatrix and likelihoods calculated using 'Familias'.
LR.All.Markers |
One LR for each simulation for each pedigree. |
lik.All.Markers |
One likehood for each simulation for each pedigree. |
LR.Per.Marker |
One LR for each simulation and marker for each pedigree. |
lik.Per.Marker |
Onelikelihood for each simulation and marker for each pedigree. |
first.Sim |
NULL |
Thore Egeland <Thore.Egeland@gmail.com> and Magnus Dehli Vigeland
Kling et al. (2017)
See also FamiliasConditional
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 | # Example
data(grandmother)
pedigrees = grandmother$pedigrees
datamatrix = grandmother$datamatrix
loci = grandmother$loci
persons = rownames(datamatrix)
Nsim = 5 #Increase to 1000
res1 = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 1, available = "GS", ref=2, seed=17)
LR.H1 = table(res1$LR.All.Markers[,1])/Nsim
# LR(H1) distribution, agrees well with theory:
# Pr(LR=0.5|H1) = 0.49005, Pr(LR=25.5|H1) = 0.50490, Pr(LR=50.5|H1) = 0.00505
# Simulate from unrelated alternative:
res2 = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 2, available = "GS", ref=2, seed=17)
LR.H2 = table(res2$LR.All.Markers[,1])/Nsim
# Try mutation
p = as.double(loci[[1]]$alleles)
loci = list(FamiliasLocus(p, 1:3, "L1", MutationModel = "Proportional",
MutationRate=0.005))
x = Familias2linkdat(pedigrees, datamatrix, loci)
res1 = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 1, available = "GS", ref=2, seed=17)
# Simulate father instead
res1 = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 1, available = "FAT", ref=2, seed=17)
LR.H1 = table(res1$LR.All.Markers[,1])/Nsim
res2 = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 2, available = "FAT", ref=2, seed=17)
LR.H2 = table(res2$LR.All.Markers[,1])/Nsim
### Example
data(Demo3Markers)
pedigrees = Demo3Markers$pedigrees
datamatrix = Demo3Markers$datamatrix
loci = Demo3Markers$loci
persons = rownames(datamatrix)
Nsim = 5
res1 = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 1, available = "Mother", ref=2, seed=177)
res1[[1]][,1] #Always the same LR! Why?:
x = Familias2linkdat(pedigrees, datamatrix,loci)
P.mother.H1 = oneMarkerDistribution(x[[1]],3,3, verbose=FALSE)
P.mother.H2 = oneMarkerDistribution(x[[2]][[1]],3,3, verbose=FALSE)
round(P.mother.H1/P.mother.H2,10) == 1
# The probability distribution of the mother is the same for both hypotheses
# and therefore we always get the same LR.
# Example
data(symmetric)
pedigrees = symmetric$pedigrees
datamatrix = symmetric$datamatrix
loci = symmetric$loci
persons = rownames(datamatrix)
res1 = paramlinkConditional(Nsim = 2, datamatrix, loci, pedigrees,
truePed = 1, available = NULL, ref=2, seed=17)
# Without mutation, all LRs 1. With mutation:
res1$LR.All.Markers
# Example Mariana's F21 example
## Not run: #Takes a few minutes; compares 'Familias' and 'paramlink'.
#Results and paramlink execution times
data(F21)
pedigrees = F21$pedigrees
datamatrix = F21$datamatrix
loci = F21$loci
persons = rownames(datamatrix)
Nsim = 1000
start.time <- Sys.time()
res1.paramlink = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 1, available = "Missing Person", ref=2, seed=17)
res2.paramlink = paramlinkConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 2, available = "Missing Person", ref=2, seed=17)
end.time <- Sys.time()
paramlink.time <- end.time - start.time
start.time <- Sys.time()
res1.familias = FamiliasConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 1, available = "Missing Person", ref=2, seed=17)
res2.familias = FamiliasConditional(Nsim = Nsim, datamatrix, loci, pedigrees,
truePed = 2, available = "Missing Person", ref=2, seed=17)
end.time <- Sys.time()
familias.time <- end.time - start.time
familias.time - paramlink.time #around -0.7 mins
LR1.familias = res1.familias[[1]][,1]
LR1.paramlink = res1.paramlink[[1]][,1]
aa=cbind(LR1.familias, LR1.paramlink)
foo=apply(aa,2, quantile)
foo = LR1.familias - LR1.paramlink
max(foo/(0.5*(LR1.familias+LR1.paramlink)))
LR1.familias = res2.familias[[1]][,1]
LR1.paramlink = res2.paramlink[[1]][,1]
aa=cbind(LR1.familias, LR1.paramlink)
foo=apply(aa,2, quantile)
foo = LR1.familias - LR1.paramlink
max(foo)
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.