test allele

doc/api.rst

@@ -48,6 +48,7 @@ Statistics
     limix.stats.linear_kinship
     limix.stats.lrt_pvalues
     limix.stats.multipletests
+    limix.stats.multivariate_normal
     limix.stats.pca
 
 Heritability estimation

limix/__init__.py

@@ -13,7 +13,7 @@
 from ._example import file_example
 from ._testit import test
 
-__version__ = "3.0.dev3"
+__version__ = "3.0.dev4"
 
 
 __all__ = [

limix/_cli/qtl.py

@@ -1,7 +1,6 @@
 import click
 
 from .scan import scan
-from .st_sscan import st_sscan
 
 
 # @click.command(cls=OrderedCommand)
@@ -27,4 +26,3 @@ def qtl(ctx):
 
 
 qtl.add_command(scan)
-qtl.add_command(st_sscan)

limix/stats/_allele.py

@@ -36,9 +36,8 @@ def compute_dosage(X, alt=None):
     X : array_like
         Allele expectations encoded as a variants-by-samples-by-alleles matrix.
     ref : array_like
-        Allele reference of each locus. The allele having the minor allele frequency
-        for the provided ``X`` is used as the reference if `None`. Defaults to
-        `None`.
+        Allele reference of each locus. The allele having the minor allele frequency for
+        the provided ``X`` is used as the reference if `None`. Defaults to ``None``.
 
     Returns
     -------

limix/stats/_random.py

@@ -1,4 +1,22 @@
 def multivariate_normal(random, mean, cov):
+    """
+    Draw random samples from a multivariate normal distribution.
+
+    Parameters
+    ----------
+    random : np.random.RandomState instance
+        Random state.
+    mean : array_like
+        Mean of the n-dimensional distribution.
+    cov : array_like
+        Covariance matrix of the distribution. It must be symmetric and
+        positive-definite for proper sampling.
+
+    Returns
+    -------
+    out : ndarray
+        The drawn sample.
+    """
     from numpy.linalg import cholesky
 
     L = cholesky(cov)

limix/stats/test_allele.py

@@ -0,0 +1,34 @@
+from limix.stats import allele_frequency, compute_dosage, allele_expectation
+from numpy.testing import assert_allclose
+from numpy.random import RandomState
+
+
+def test_allele_frequency():
+    X = RandomState(0).randn(5, 3)
+    assert_allclose(
+        allele_frequency(X),
+        [2.0422233965933585, 0.7940255259577332, 0.5987926640480571],
+    )
+
+
+def test_compute_dosage():
+    X = RandomState(0).randn(2, 3, 4)
+    assert_allclose(
+        compute_dosage(X).tolist(),
+        [
+            [2.240893199201458, -0.1513572082976979, 1.454273506962975],
+            [0.33367432737426683, -0.8540957393017248, -0.7421650204064419],
+        ],
+    )
+
+
+def test_allele_expectation():
+    p = RandomState(0).rand(2, 3)
+    p = p / p.sum(1)[:, None]
+    assert_allclose(
+        allele_expectation(p, 2, 2),
+        [
+            [0.9710998244964055, 1.028900175503594],
+            [0.9374325310073951, 1.0625674689926048],
+        ],
+    )