Description Usage Arguments Details Value Author(s) References Examples
For these functions to work, the program MERLIN (see References below) must
be installed and correctly pointed to in the PATH variable. The merlin()
function is a general wrapper which runs MERLIN with the indicated options,
after creating the appropriate input files. For convenience, MERLIN's
"–likelihood" functionality is wrapped in a separate function.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 | merlin(
x,
options,
markers = NULL,
linkageMap = NULL,
verbose = TRUE,
generateFiles = TRUE,
cleanup = TRUE,
dir = tempdir(),
logfile = NULL,
merlinpath = NULL
)
likelihoodMerlin(
x,
markers = NULL,
linkageMap = NULL,
rho = NULL,
logbase = NULL,
options = "--likelihood --bits:100 --megabytes:4000 --quiet",
...
)
checkMerlin()
|
x |
a |
options |
a single string containing all arguments to merlin except for the input file indications. |
markers |
a vector of names or indices of markers attached to |
linkageMap |
a data frame with three columns (chromosome; marker name; centiMorgan position) to be used as the marker map by MERLIN. |
verbose |
a logical. |
generateFiles |
a logical. If TRUE (default), input files to MERLIN
named '_merlin.ped', '_merlin.dat', '_merlin.map', and '_merlin.freq' are
created in the directory indicated by |
cleanup |
a logical. If TRUE (default), the MERLIN input files are deleted after the call to MERLIN. |
dir |
the name of the directory where input files should be written. |
logfile |
a character. If this is given, the MERLIN screen output will be dumped to a file with this name. |
merlinpath |
the path to the folder containing the merlin executables. If the executables are on the system's search path, this can be left as NULL (default). |
rho |
A vector of length one less than the number of markers, specifying the recombination rate between each consecutive pair. |
logbase |
Either NULL (default) or a positive number indicating the
basis for logarithmic output. Typical values are |
... |
Further arguments passed on to |
The merlin()
function creates input files "_merlin.ped", "_merlin.dat",
"_merlin.map" and "_merlin.freq" in the dir
directory, and then runs the
following command through a call to system()
:
1 2 |
likelihoodMerlin()
first runs merlin()
with options = "--likelihood --bits:100 --megabytes:4000 --quiet"
, and then extracts the likelihood
values from the MERLIN output. Note that the output is the total likelihood
including all markers.
For likelihood computations with linked markers, the argument rho
should
indicate the recombination fractions between each consecutive pair of markers
(i.e., rho[i]
is the recombination rate between markers i-1
and i
).
These will be converted to centiMorgan distances using Haldane's map
function, and used to create genetic marker map in a MERLIN-friendly format.
merlin()
returns the screen output of MERLIN invisibly.
likelihoodMerlin()
returns a single number; the total likelihood using
all indicated markers.
checkMerlin()
returns TRUE if MERLIN is installed and available on the
system path, and FALSE otherwise.
Magnus Dehli Vigeland
http://csg.sph.umich.edu/abecasis/Merlin/
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 | if(checkMerlin()) {
### Trivial example for validation
x = nuclearPed(1)
m1 = marker(x, "1" = 1:2) # likelihood = 1/2
m2 = marker(x, "1" = 1, "3" = 1:2) # likelihood = 1/8
x = setMarkers(x, list(m1, m2))
# MERLIN likelihoods
lik1 = likelihoodMerlin(x, markers = 1, verbose = FALSE)
lik2 = likelihoodMerlin(x, markers = 2, verbose = FALSE)
likTot = likelihoodMerlin(x, verbose = FALSE)
stopifnot(all.equal(
round(c(lik1, lik2, likTot), c(3,3,4)), c(1/2, 1/8, 1/16)))
# Example with ped lists
y = list(singleton(1), singleton(2))
y = setMarkers(y, locus = list(alleles = 1:2))
genotype(y[[1]], marker = 1, id = '1') = 1:2
genotype(y[[2]], marker = 1, id = '2') = 1
lik = likelihoodMerlin(y, verbose = FALSE)
stopifnot(all.equal(round(lik, 3), 1/8))
### Linked markers
z = nuclearPed(2)
m = marker(z, geno = c("1/1", "1/2", "1/2", "1/2"))
z = setMarkers(z, list(m, m))
# By MERLIN...
L1 = likelihoodMerlin(z, markers = 1:2, rho = 0.25, verbose = FALSE)
# ...and by pedprobr
L2 = likelihood2(z, marker1 = 1, marker2 = 2, rho = 0.25)
# stopifnot(all.equal(signif(L1, 3), signif(L2, 3)))
}
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.