Singlevariant testing¶
Linear models¶

limix.qtl.
qtl_test_lm
(snps, pheno, covs=None, test=’lrt’, verbose=None)[source]¶ Wrapper function for univariate singlevariant association testing using a linear model.
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 covs (ndarray, optional) – (N, D) ndarray of D covariates for N individuals.
By default,
covs
is a (N, 1) array of ones.  test ({'lrt', 'f'}, optional) – test statistic. ‘lrt’ for likelihood ratio test (default) or ‘f’ for Ftest.
 verbose (bool, optional) – if True, details such as runtime as displayed.
Returns: LIMIX LMM object
Return type: Examples
>>> from numpy.random import RandomState >>> from limix.qtl import qtl_test_lm >>> from numpy import set_printoptions >>> set_printoptions(4) >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, 1) >>> >>> lm = qtl_test_lm(snps, pheno) >>> print(lm.getPv()[:,:4]) [[ 0.8796 0.5065 0.5666 0.6016]]
Linear mixed models¶

limix.qtl.
qtl_test_lmm
(snps, pheno, K=None, covs=None, test=’lrt’, NumIntervalsDelta0=100, NumIntervalsDeltaAlt=100, searchDelta=False, verbose=None)[source]¶ Wrapper function for univariate singlevariant association testing using a linear mixed model.
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 K (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients. If not provided, then standard linear regression is considered.
 covs (ndarray, optional) – (N, D) ndarray of D covariates for N individuals.
By default,
covs
is a (N, 1) array of ones.  test ({'lrt', 'f'}, optional) – test statistic. ‘lrt’ for likelihood ratio test (default) or ‘f’ for Ftest.
 NumIntervalsDelta0 (int, optional) – number of steps for delta optimization on the null model.
By default
NumIntervalsDelta0
is 100.  NumIntervalsDeltaAlt (int, optional) – number of steps for delta optimization on the alternative model.
Requires
searchDelta=True
to have an effect.  searchDelta (bool, optional) – if True, delta optimization on the alternative model is carried out.
By default
searchDelta
is False.  verbose (bool, optional) – if True, details such as runtime as displayed.
Returns: LIMIX LMM object
Return type: Examples
>>> from numpy.random import RandomState >>> from numpy import dot >>> from limix.qtl import qtl_test_lmm >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, 1) >>> W = random.randn(N, 10) >>> kinship = dot(W, W.T) / float(10) >>> >>> lmm = qtl_test_lmm(snps, pheno, kinship) >>> print(lmm.getPv()[:,:4]) [[ 0.8571 0.4668 0.5872 0.5589]]

limix.qtl.
qtl_test_interaction_lmm
(snps, pheno, Inter, Inter0=None, covs=None, K=None, test=’lrt’)[source]¶ Wrapper function for singlevariant interaction test using LMMs.
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 Inter (ndarray, optional) – (N, I) ndarray of I interaction variables to be tested for N individuals.
 Inter0 (ndarray, optional) – (N, I0) ndarray of I0 interaction variables to be included for N individuals to be included both in the alternative and in the null model. By default Inter0 is a (N, 1) array of ones.
 covs (ndarray, optional) – (N, D) ndarray of D covariates for N individuals.
By default,
covs
is a (N, 1) array of ones.  K (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients (optional) If not provided, then standard linear regression is considered.
 test ({'lrt', 'f'}, optional) – test statistic. ‘lrt’ for likelihood ratio test (default) or ‘f’ for Ftest.
Returns: limix CInteractLMM object
Return type: limix_legacy.deprecated.CInteractLMM()
Examples
>>> from numpy.random import RandomState >>> from numpy import dot, eye, ones, concatenate >>> from limix.qtl import qtl_test_interaction_lmm >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, 1) >>> W = random.randn(N, 10) >>> kinship = dot(W, W.T) / float(10) >>> kinship+= 1e4 * eye(N) >>> >>> #environments to include in null and alt models >>> Inter0 = ones((N,1)) >>> >>> #environments to include in alt model >>> Inter = random.randn(N,2) >>> >>> #environments are also included as fixed effect covariates >>> mean = ones([N, 1]) >>> covs = concatenate([mean, Inter], 1) >>> >>> #interaction test >>> lmi = qtl_test_interaction_lmm(snps, pheno, Inter=Inter, ... Inter0=Inter0) >>> pvi = lmi.getPv() >>> >>> print(pvi.shape) (1, 1000) >>> print(pvi[:,:4]) [[ 0.8179 0.2185 0.3338 0.077 ]]

limix.qtl.
qtl_test_lmm_kronecker
(snps, phenos, covs=None, Acovs=None, Asnps=None, K1r=None, K1c=None, K2r=None, K2c=None, NumIntervalsDelta0=100, NumIntervalsDeltaAlt=100, searchDelta=False)[source]¶ Wrapper function for singlevariant multitrait association testing class using linear mixed models.
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 covs (ndarray or list of ndarrays, optional) – sample design matrices of the different fixed effect terms.
In the case of a single fixed effect term,
covs
can be specified as a (N, D) ndarray of D covariates for N individuals. In this caseAcovs
should also be a (K, P) ndarray, where K is the rank of the trait design and P is the number of traits. For multiple terms,covs
is a list of ndarray, where each elementcovs[i]
has shape (N, Di). In this case,Acovs
should also be list with the same number of ndarray elements. Each elementAcovs[i]
has shape (Ki, P), where P is the number of analyzed phenotypes. By default covs is (N, 1) ndarray of ones andAcovs[i]
is a (P, P) identity ndarray.  Acovs (ndarray or list of ndarrays, optional) – trait design matrices of the different fixed effect terms. See covs for more info.
 Asnps (ndarray, optional) – (K, P) trait design of rank K of the snp fixed effect on the
P traits. By default
Asnps
is (1, P) ndrray of ones.  K1r (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients.
If not provided,
K1r = numpy.dot(snps, snps.T)
is assumed.  K1c (ndarray, optional) – (P, P) ndarray of LMMcovariance/kinship coefficients. If either K1c or K2c is not provided, K1c is estimated by (restricted) maximum likelihood from the model without the snp effect.
 K2r (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients.
By default
K2r
is a (N, N) identity ndarray.  K2c (ndarray, optional) – (P, P) ndarray of residual covariances. If either K1c or K2c is not provided, K1c is estimated by (restricted) maximum likelihood from the model without the snp effect.
 NumIntervalsDelta0 (int, optional) – number of steps for delta optimization on the null model.
By default
NumIntervalsDelta0
is 100.  NumIntervalsDeltaAlt (int, optional) – number of steps for delta optimization on the alternative model.
Requires
searchDelta=True
to have an effect.  searchDelta (bool, optional) – if True, delta optimization on the alternative model is carried out.
By default
searchDelta
is False.
Returns:  tuple containing:
 lmm (:class:`limix_legacy.deprecated.CKroneckerLMM()`): CKroneckerLMM limix_legacy object
 pv (ndarray): Pvalues for all SNPs from liklelihood ratio test
Return type: (tuple)
Examples
Example of multitrait singlevariant association testing using a linear mixed model and specifying a common effects.
>>> from numpy.random import RandomState >>> from numpy import dot, eye, ones >>> from limix.qtl import qtl_test_lmm_kronecker >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> P = 2 >>> >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, P) >>> W = random.randn(N, 10) >>> kinship = dot(W, W.T) / float(10) >>> kinship+= 1e4 * eye(N) >>> >>> #common effect test >>> Asnps = ones((1, P)) >>> lmm, pv = qtl_test_lmm_kronecker(snps, pheno, K1r=kinship, ... Asnps=Asnps) >>> beta = lmm.getBetaSNP() >>> >>> print(pv.shape) (1, 1000) >>> print(pv[:,:4]) [[ 0.8722 0.6509 0.3593 0.6816]] >>> >>> print(beta.shape) (1, 1000) >>> print(beta[:,:4]) [[ 0.0281 0.0758 0.1458 0.0701]]
Example showing how to set an any effect test. For more information on effect designs see limix_tutorials.
>>> #any effect test >>> Asnps = eye(P) >>> lmm, pv = qtl_test_lmm_kronecker(snps, pheno, K1r=kinship, ... Asnps=Asnps) >>> beta = lmm.getBetaSNP() >>> >>> print(pv.shape) (1, 1000) >>> print(pv[:,:4]) [[ 0.9838 0.8664 0.3555 0.5244]] >>> >>> print(beta.shape) (2, 1000) >>> print(beta[:,:4]) [[ 0.0413 0.0312 0.0177 0.098 ] [ 0.0139 0.1234 0.3204 0.2495]]

limix.qtl.
qtl_test_interaction_lmm_kronecker
(snps, phenos, covs=None, Acovs=None, Asnps1=None, Asnps0=None, K1r=None, K1c=None, K2r=None, K2c=None, NumIntervalsDelta0=100, NumIntervalsDeltaAlt=100, searchDelta=False, return_lmm=False)[source]¶ Wrapper function for singlevariant multitrait interaction test
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 covs (ndarray or list of ndarrays, optional) – sample design matrices of the different fixed effect terms. See limix.qtl.qtl_test_lmm_kronecker for more info.
 Acovs (ndarray or list of ndarrays, optional) – trait design matrices of the different fixed effect terms. See limix.qtl.qtl_test_lmm_kronecker for more info.
 Asnps1 (ndarray or list of ndarray, optional) – (K, P) trait design of rank K of the snp fixed effect
on the P traits for the alternative model.
Asnps1
is tested against the null trait designAsnps0
. Multiple designs can be tested vs the null trait design If a list of ndarray trait designs is specified forAsnps1
, the multiple designs are tested vs the null trait designAsnps0
onebyone. By defaultAsnps1
is a (P, P) identity ndarray.  Asnps0 (ndarray, optional) – (K0, P) trait design of rank K of the snp fixed effect
on the P traits for the null model.
By default
Asnps0
is (1, P) ndrray of ones.  K1r (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients.
If not provided,
K1r = numpy.dot(snps, snps.T)
is assumed.  K1c (ndarray, optional) – (P, P) ndarray of LMMcovariance/kinship coefficients. If either K1c or K2c is not provided, K1c is estimated by (restricted) maximum likelihood from the model without the snp effect.
 K2r (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients.
By default
K2r
is a (N, N) identity ndarray.  K2c (ndarray, optional) – (P, P) ndarray of residual covariances. If either K1c or K2c is not provided, K1c is estimated by (restricted) maximum likelihood from the model without the snp effect.
 NumIntervalsDelta0 (int, optional) – number of steps for delta optimization on the null model.
By default
NumIntervalsDelta0
is 100.  NumIntervalsDeltaAlt (int, optional) – number of steps for delta optimization on the alternative model.
Requires
searchDelta=True
to have an effect.  searchDelta (bool, optional) – if True, delta optimization on the alternative model is carried out.
By default
searchDelta
is False.
Returns:  tuple containing:
 pv (ndarray):
Pvalues of the interaction test (
Asnps1
vsAsnps0
) IfAsnps1
is a list of trait designs, pv contains the P values for the different trait designs (as rows).  pv0 (ndarray):
Pvalues of the association test
Asnps0
vs0
.  pvAlt (ndarray):
Pvalues of the association test
Asnps1
vs0
. IfAsnps1
is a list of trait designs, pvAlt contains the P values for the different trait designs (as rows).
 pv (ndarray):
Pvalues of the interaction test (
Return type: (tuple)
Examples
Example of how to perform a specific, common and any effect test. For more information see limix_tutorials.
>>> from numpy.random import RandomState >>> from numpy import dot, eye, ones >>> from limix.qtl import qtl_test_interaction_lmm_kronecker >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> P = 2 >>> >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, P) >>> W = random.randn(N, 10) >>> kinship = dot(W, W.T) / float(10) >>> kinship+= 1e4 * eye(N) >>> >>> #interaction effect test >>> Asnps0 = ones((1,P)) >>> Asnps1 = eye(P) >>> pv, pv0, pvAlt = qtl_test_interaction_lmm_kronecker(snps, pheno, ... K1r=kinship, ... Asnps0=Asnps0, ... Asnps1=Asnps1) >>> #interaction P value >>> print(pv.shape) (1, 1000) >>> print(pv[:,:4]) [[ 0.9347 0.7744 0.2678 0.2893]] >>> >>> #common effect P value >>> print(pv0.shape) (1, 1000) >>> print(pv0[:,:4]) [[ 0.8722 0.6509 0.3593 0.6816]] >>> >>> #any effect P value >>> print(pvAlt.shape) (1, 1000) >>> print(pvAlt[:,:4]) [[ 0.9838 0.8664 0.3555 0.5244]]

limix.qtl.
forward_lmm
(snps, pheno, K=None, covs=None, qvalues=False, threshold=5e08, maxiter=2, test=’lrt’, verbose=None, **kw_args)[source]¶ Wrapper function for univariate singlevariant test with forward selection
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 K (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients (optional) If not provided, then standard linear regression is considered.
 covs (ndarray, optional) – (N, D) ndarray of D covariates for N individuals.
By default,
covs
is a (N, 1) array of ones.  qvalues (bool, optional) – if True,
threshold
is set to Storey qvalues to control FDR. (default is False)  threshold (float, optional) – Pvalue thrashold for inclusion in forward selection (default 5e8).
If
qvalues=True
, the threshold is set to to Storey qvalues.  maxiter (int, optional) – maximum number of interaction scans. First scan is without inclusion, so maxiter1 inclusions can be performed. (default 2)
 test ({'lrt', 'f'}, optional) – test statistic. ‘lrt’ for likelihood ratio test (default) or ‘f’ for Ftest.
 verbose (bool, optional) – print verbose output. (default is False)
Returns:  tuple containing:
 lmm (:class:`limix.qtl.LMM`): LIMIX LMM object
 res (dict): {iadded, pvadded, pvall}. iadded is an ndarray of indices of SNPs included in order of inclusion, pvadded is an ndarray of Pvalues obtained by the included SNPs in iteration and pvall is a (Nadded, S) ndarray of Pvalues for all iterations
Return type: (tuple)
Examples
>>> from numpy.random import RandomState >>> from numpy import dot, eye, ones >>> from limix.qtl import forward_lmm >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, 1) >>> W = random.randn(N, 10) >>> kinship = dot(W, W.T) / float(10) >>> kinship+= 1e4 * eye(N) >>> >>> #forward lmm >>> lmm, res = forward_lmm(snps, pheno, K=kinship, threshold=1.) >>> >>> print(res['pvall'].shape) (2, 1000) >>> print(res['pvall'][:,:4]) [[ 0.8571 0.4668 0.5872 0.5589] [ 0.77 0.4226 0.6165 0.8727]]

class
limix.qtl.
LMM
(snps, pheno, K=None, covs=None, test=’lrt’, NumIntervalsDelta0=100, NumIntervalsDeltaAlt=100, searchDelta=False, verbose=None)[source]¶ Basic class for univariate singlevariant association testing with LMMs.
Parameters:  snps (ndarray) – (N, S) ndarray of S SNPs for N individuals.
 pheno (ndarray) – (N, P) ndarray of P phenotype sfor N individuals. If phenotypes have missing values, then the subset of individuals used for each phenotype column will be subsetted.
 K (ndarray, optional) – (N, N) ndarray of LMMcovariance/kinship coefficients (optional) If not provided, then standard linear regression is considered.
 covs (ndarray, optional) – (N, D) ndarray of D covariates for N individuals.
By default,
covs
is a (N, 1) array of ones.  test ({'lrt', 'f'}, optional) – test statistic. ‘lrt’ for likelihood ratio test (default) or ‘f’ for Ftest.
 NumIntervalsDelta0 (int, optional) – number of steps for delta optimization on the null model.
By default
NumIntervalsDelta0
is 100.  NumIntervalsDeltaAlt (int, optional) – number of steps for delta optimization on the alternative model.
Requires
searchDelta=True
to have an effect.  searchDelta (bool, optional) – if True, delta optimization on the alternative model is carried out.
By default
searchDelta
is False.  verbose (bool, optional) – if True, details such as runtime as displayed.
Examples
Example of singlevariant association testing using a linear mixed model.
>>> from numpy.random import RandomState >>> from limix.qtl import LMM >>> from numpy import dot >>> from numpy import set_printoptions >>> set_printoptions(4) >>> random = RandomState(1) >>> >>> N = 100 >>> S = 1000 >>> >>> # generate data >>> snps = (random.rand(N, S) < 0.2).astype(float) >>> pheno = random.randn(N, 1) >>> W = random.randn(N, 10) >>> kinship = dot(W, W.T) / float(10) >>> >>> # run singlevariant associaiton testing with LMM >>> lmm = LMM(snps, pheno, kinship) >>> pv = lmm.getPv() >>> beta = lmm.getBetaSNP() >>> beta_ste = lmm.getBetaSNPste() >>> >>> print(pv.shape) (1, 1000) >>> print(pv[:,:4]) [[ 0.8571 0.4668 0.5872 0.5589]] >>> print(beta[:,:4]) [[ 0.0436 0.1695 0.12 0.1388]] >>> print(beta_ste[:,:4]) [[ 0.2422 0.2329 0.221 0.2375]]
As shown in the exmaple below, the method can be applied directly on multiple phenotypes. Each phenotype is indipendently tested against all variants onebyone.
>>> random = RandomState(1) >>> >>> P = 3 >>> phenos = random.randn(N, P) >>> >>> lmm = LMM(snps, phenos, kinship) >>> pv = lmm.getPv() >>> >>> print(pv.shape) (3, 1000) >>> print(pv[:,:4]) [[ 0.4712 0.694 0.3369 0.5221] [ 0.7336 0.8799 0.8109 0.1071] [ 0.0662 0.9203 0.2873 0.8268]]
Generalised linear mixed models¶

limix.qtl.
qtl_test_glmm
(snps, pheno, lik, K, covs=None, test=’lrt’, NumIntervalsDeltaAlt=100, searchDelta=False, verbose=True)[source]¶ Wrapper function for univariate singlevariant association testing using a generalised linear mixed model.
Parameters:  snps (array_like) – N individuals by S SNPs.
 pheno (tuple, array_like) – Either a tuple of two arrays of N individuals each (Binomial phenotypes) or an array of N individuals (Poisson or Bernoulli phenotypes). It does not support missing values yet.
 lik ({'bernoulli', 'binomial', 'poisson'}) – Sample likelihood describing the residual distribution.
 K (array_like) – N by N covariance matrix (e.g., kinship coefficients).
 covs (array_like, optional) – N individuals by D covariates.
By default,
covs
is a (N, 1) array of ones.  test ({'lrt'}, optional) – Likelihood ratio test (default).
 NumIntervalsDeltaAlt (int, optional) – number of steps for delta optimization on the alternative model.
Requires
searchDelta=True
to have an effect.  searchDelta (bool, optional) – if
True
, delta optimization on the alternative model is carried out. By defaultsearchDelta
isFalse
.  verbose (bool, optional) – if
True
, details such as runtime are displayed.
Returns: LIMIX LMM object
Return type: Examples
>>> from numpy import dot, exp, sqrt >>> from numpy.random import RandomState >>> from limix.qtl import qtl_test_glmm >>> >>> random = RandomState(0) >>> >>> G = random.randn(250, 500) / sqrt(500) >>> beta = 0.01 * random.randn(500) >>> >>> z = dot(G, beta) + 0.1 * random.randn(250) >>> z += dot(G[:, 0], 1) # causal SNP >>> >>> y = random.poisson(exp(z)) >>> >>> candidates = G[:, :5] >>> K = dot(G[:, 5:], G[:, 5:].T) >>> lm = qtl_test_glmm(candidates, y, 'poisson', K, verbose=False) >>> >>> print(lm.getPv()) [[ 0.0694 0.3336 0.5899 0.7388 0.7796]]