Please note that full documentation for LDAK-KVIK is provided at www.ldak-kvik.com (this page only provides a summary).
This page explains how to run LDAK-KVIK, our new tool for performing super-fast linear and logistic mixed-model association analysis. Note that when using LDAK-KVIK for single-predictor analysis, there are two steps (Step 1 constructs LOCO PRS, while Step 2 performs the single-SNP association analysis); when using LDAK-KVIK for gene-based analysis, there is also a third step (Step 3 performs the gene-based association analysis).
If running LDAK-KVIK on a cloud server (e.g., those used by UK Biobank, AllofUs or Our Future Health), we recommend first copying the genetic data to the local node; see this GitHub Issue for more details.
If performing a gene-based association analysis, you will require gene annotations. On Resources, you can download gene annotations for the human genome. The page also explains how I created these files, which might be useful if you wish to construct annotations yourself (e.g., for non-human species).
Always read the screen output, which suggests arguments and estimates memory usage.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
LDAK-KVIK Step 1:
The main argument is --kvik-step1 <outfile>.
This requires the options
--bfile/--speed <datastem> or --bgen <datafile> - to specify the genetic data files (see File Formats). Note that it is not possible to use SP or SPED formats; if you have data in these format, you should first convert to SPEED format (see Make Data).
--pheno <phenofile> - to specify phenotypes (in PLINK format). If <phenofile> contains more than one phenotype, specify which should be used with --mpheno <integer> or --pheno-name <string> (the latter requires that <phenofile> has a header row), or use --mpheno ALL to analyse all phenotypes. Note that samples without a phenotype will be excluded (unless using --mpheno ALL, in which case missing phenotypes will be mean-imputed).
You can use --keep <keepfile> and/or --remove <removefile> to restrict to a subset of samples, and --extract <extractfile> and/or --exclude <excludefile> to restrict to a subset of predictors (for more details, see Data Filtering).
You can use --covar <covarfile> or --factors <factorfile> to provide quantitative or categorical covariates (in PLINK format) as fixed effects in the regression (if <covarfile> contains multiple covariates, you can specify a subset using either --covar-names <string> or --covar-numbers <string>).
If the phenotype is binary, you must tell LDAK whether to perform linear regression (--binary NO) or logistic regression (--binary YES).
If you have very many predictors (e.g., over one million), which is common when analysing imputed SNP genotypes, please read the Dense Data Advice below.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
LDAK-KVIK Step 2:
The main argument is --kvik-step2 <outfile>, where the output file name matches that used in Step 1.
This requires the options
--bfile/--gen/--sp/--speed <datastem> or --bgen <datafile> - to specify the genetic data files (see File Formats).
--pheno <phenofile> - to specify phenotypes (in PLINK format). If <phenofile> contains more than one phenotype, specify which should be used with --mpheno <integer>, or use --mpheno ALL to analyse all phenotypes. Note that samples without a phenotype will be excluded (unless using --mpheno ALL, in which case missing phenotypes will be mean-imputed).
You can use --keep <keepfile> and/or --remove <removefile> to restrict to a subset of samples, and --extract <extractfile> and/or --exclude <excludefile> to restrict to a subset of predictors (for more details, see Data Filtering).
You can use --covar <covarfile> or --factors <factorfile> to provide quantitative or categorical covariates (in PLINK format) as fixed effects in the regression (if <covarfile> contains multiple covariates, you can specify a subset using either --covar-names <string> or --covar-numbers <string>).
Please note that the options used in Step 2 must be consistent with those used in Step 1 (e.g., if you used covariates in Step 1, you must also use them in Step 2).
The main output file is <outfile>.step2.assoc, which contains estimated effect sizes, test statistics and p-values for each predictor (note that effect sizes are relative to the A1 alleles). Meanwhile <outfile>.step2.summaries contains summary statistics (in the format required for use with SumHer and MegaPRS).
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
LDAK-KVIK Step 3 (only required for gene-based association testing):
The main argument is --kvik-step3 <outfile>, where the output file name matches that used in Steps 1 & 2.
This requires the options
--bfile/--gen/--sp/--speed <datastem> or --bgen <datafile> - to specify the genetic data files (see File Formats).
--genefile <genefile> - to specify the gene annotations file. The first four columns of <genefile> should specify the name, chromosome, start and end basepairs of each gene. If there is a fifth column, its value should be either + or -, indicating whether the gene is on the forward or backwards strand. Note that the basepairs should provide 0-start, half-open positions (this is the format used by Browser Extensible Data). For example, if Gene ABC is on Chromosome 7 and contains basepairs 1-10, inclusive, the corresponding row of <genefile> would be "ABC 7 0 10" (if providing strand orientation, the row would be either "ABC 7 0 10 +" or "ABC 7 0 10 -").
By default, LDAK will only consider predictors inside genes. To change this, you can use --gene-buffer <integer>, --up-buffer <integer> or --down-buffer <integer>. For example, if you use --gene-buffer 500, each set will contain predictors inside or within 500bp of a gene. Note that you can only use --up-buffer <integer> or --down-buffer <integer> if <genefile> provides the orientation of the gene (i.e., has five columns). If you add --overlap NO, LDAK will ensure that each predictor is only allocated to the nearest gene.
The final results will be saved in <output>.step3.remls.all. For each gene, this file reports the estimated heritability and two p-values (we recommend that you use the second p-value, obtained using permutations). The estimated genetic contributions of all gene will be saved in SP format in the files <output>.step3.prs.all.sp, <output>.step3.prs.all.bim and <output>.step3.prs.all.fam (see Genetic File Formats for details), while the corresponding p-values are saved in <output>.step3.prs.all.pvalues.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
Dense data advice:
Step 1 of LDAK-KVIK can be slow when there are very many predictors. Therefore, if your dataset contains over one million predictors, we recommend you only use a subset of these in Step 1 (e.g., 500,000 predictors). Doing so will have a limited impact on power, but will greatly reduce runtime. Note that you should continue to use all predictors in Step 2.
You may already have a suitable subset of predictors (e.g., you may have a list of directly-genotyped SNPs, or those passing stringent quality control). Otherwise, we you can obtain the SNP subset by performing a moderate thinning of the common predictors (e.g., if analysing SNP data, we suggest identifying SNPs with MAF > 0.01, then filtering so there are no predictors within 100kb with squared correlation above 0.5).
This can be achieved using the main argument --thin-common <output>.
The only required option is
--bfile/--gen/--sp/--speed <datastem> or --bgen <datafile> - to specify the genetic data files (see File Formats).
For details of optional arguments, see Thin Predictors.
When analysing SNP data, this command will produce the file <output>.in, which contains SNPs with MAF > 0.01, thinned so no pair within 100kb has squared correlation above 0.5 (if not using SNP data, the command will instead restrict to predictors with variance > 0.01).
You can then run Step 1 of LDAK-KVIK adding --extract <output>.in.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
Example:
Here we use the binary PLINK files human.bed, human.bim and human.fam, the phenotypes quant.pheno and binary.pheno, and the covariates human.covar from the Test Datasets. For the gene-based analysis, we also use the file anns.txt, which contains annotations for human genes on Chromosomes 21 & 22
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
1 - Linear mixed-model regression (with non-dense data).
We perform single-SNP linear mixed-model regression for the quantitative phenotype quant.pheno by running
./ldak.out --kvik-step1 kvik --bfile human --pheno quant.pheno
followed by
./ldak.out --kvik-step2 kvik --bfile human --pheno quant.pheno
The main results are saved in kvik.step2.assoc. To repeat this analysis including the covariates, we run
./ldak.out --kvik-step1 kvik2 --bfile human --pheno quant.pheno --covar human.covar
followed by
./ldak.out --kvik-step2 kvik2 --bfile human --pheno quant.pheno --covar human.covar
The main results are saved in kvik2.step2.assoc.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
2 - Logistic mixed-model regression (with non-dense data).
We perform single-SNP logistic mixed-model regression for the binary phenotype binary.pheno by running
./ldak.out --kvik-step1 kvik3 --bfile human --pheno binary.pheno --binary YES
followed by
./ldak.out --kvik-step2 kvik3 --bfile human --pheno binary.pheno
The main results are saved in kvik3.step2.assoc.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
3 - Gene-based association analysis (with non-dense data).
We can perform a gene-based analysis for the phenotype quant.pheno by running the following command (note that this analysis uses results from the first part of this example, above).
./ldak.out --kvik-step3 kvik --bfile human --genefile anns.txt
The main results are saved in kvik.step3.remls.all, while estimates of the genetic contributions of each gene are saved in the files kvik.step3.prs.all.sp, kvik.step3.prs.all.bim and kvik.step3.prs.all.fam (the latter files can be used for Clumping). We can find the most significant gene by sorting based on the likelihood ratio test statistic (Column 10).
head -n 1 kvik.step3.remls.all
Gene_Name Gene_Chr Gene_Start Gene_End Length Heritability SD Null_Likelihood Alt_Likelihood LRT_Stat LRT_P_Raw LRT_P_Perm
sort -n -k 10 kvik.step3.remls.all | tail -n 1
CECR2 22 17956628 18037856 20 0.056225 0.024710 -1447.5822 -1432.2956 30.5731 1.6076e-08 5.3785e-08
We see that the most significant gene is CECR2. Its estimated heritability is 0.06 (SD 0.02), while its test statistic is 31. There are two p-value, 1.6e-8 and 5.4e-8. We prefer the second p-value (Column 12), which was obtained via permutations (note that this analysis is stochastic, so results will vary between runs).
Alternatively, we can perform a gene-based analysis for the phenotype binary.pheno by running the following command (this analysis uses results from the second part of this example, above).
./ldak.out --kvik-step3 kvik3 --bfile human --genefile anns.txt
The main results are saved in kvik3.step3.remls.all.
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
4 - Linear mixed-model regression with dense data.
Suppose we had very many SNPs (which is not the case here). We perform single-SNP linear mixed-model regression for the quantitative phenotype quant.pheno by running the following three commands
./ldak.out --thin-common subset --bfile human
./ldak.out --kvik-step1 kvik4 --bfile human --pheno quant.pheno --extract subset.in
./ldak.out --kvik-step2 kvik4 --bfile human --pheno quant.pheno
The first command generates subset.in, a list of thinned common predictors, which are then provided to LDAK in the second command. Note that although only a subset of SNPs are used in Step 1, all SNPs are used in Step 2 (and so the main results file kvik4.step2.assoc contains results for all SNPs).
