diff options
-rw-r--r-- | wqflask/wqflask/marker_regression/gemma_mapping.py | 81 | ||||
-rw-r--r-- | wqflask/wqflask/marker_regression/marker_regression.py | 3 | ||||
-rw-r--r-- | wqflask/wqflask/templates/show_trait_mapping_tools.html | 16 | ||||
-rw-r--r-- | wqflask/wqflask/views.py | 2 |
4 files changed, 81 insertions, 21 deletions
diff --git a/wqflask/wqflask/marker_regression/gemma_mapping.py b/wqflask/wqflask/marker_regression/gemma_mapping.py index 31e94266..2c99d1fc 100644 --- a/wqflask/wqflask/marker_regression/gemma_mapping.py +++ b/wqflask/wqflask/marker_regression/gemma_mapping.py @@ -3,12 +3,12 @@ import os, math from base import webqtlConfig from base.trait import GeneralTrait from base.data_set import create_dataset -from utility.tools import flat_files, GEMMA_COMMAND +from utility.tools import flat_files, GEMMA_COMMAND, GEMMA_WRAPPER_COMMAND import utility.logger logger = utility.logger.getLogger(__name__ ) -def run_gemma(this_dataset, samples, vals, covariates, method): +def run_gemma(this_dataset, samples, vals, covariates, method, use_loco): """Generates p-values for each marker using GEMMA""" print("INSIDE GEMMA_MAPPING") @@ -18,7 +18,13 @@ def run_gemma(this_dataset, samples, vals, covariates, method): if not os.path.isfile("{}{}_output.assoc.txt".format(webqtlConfig.GENERATED_IMAGE_DIR, this_dataset.group.name)): open("{}{}_output.assoc.txt".format(webqtlConfig.GENERATED_IMAGE_DIR, this_dataset.group.name), "w+") - logger.debug("COVARIATES_GEMMA:", covariates) + this_chromosomes = this_dataset.species.chromosomes.chromosomes + chr_list_string = "" + for i in range(this_chromosomes): + if i < (len(this_chromosomes) - 1): + chr_list_string += this_chromosomes[i+1].name + "," + else: + chr_list_string += this_chromosomes[i+1].name if covariates != "": gen_covariates_file(this_dataset, covariates) @@ -35,19 +41,38 @@ def run_gemma(this_dataset, samples, vals, covariates, method): # use GEMMA_RUN in the next one, create a unique temp file else: logger.debug("FLAT FILES:", flat_files('mapping')) - #gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.sXX.txt -lmm 1 -maf 0.1 -c %s/%s_covariates.txt -o %s_output' % (flat_files('genotype/bimbam'), - gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.cXX.txt -lmm 1 -maf 0.1 -c %s/%s_covariates.txt -outdir %s -debug -o %s_output' % (flat_files('genotype/bimbam'), - this_dataset.group.name, - flat_files('genotype/bimbam'), - this_dataset.group.name, - flat_files('genotype/bimbam'), - this_dataset.group.name, - flat_files('genotype/bimbam'), - this_dataset.group.name, - flat_files('mapping'), - this_dataset.group.name, - webqtlConfig.GENERATED_IMAGE_DIR, - this_dataset.group.name) + if use_loco: + generate_k_command = GEMMA_WRAPPER_COMMAND + ' --json --loco ' + chr_list_string + ' -- -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -gk -debug > /home/zas1024/tmp/gn2/K.json' % (flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name) + logger.debug("k_command:" + generate_k_command) + os.system(generate_k_command) + + run_loco_command = GEMMA_WRAPPER_COMMAND + ' --json --loco --input /home/zas1024/tmp/gn2/K.json -- -g %s/%s_geno.txt -p %s/%s_pheno.txt -c %s/%s_covariates.txt -a %s/%s_snps.txt -lmm 1 -maf 0.1 -debug > /home/zas1024/tmp/gn2/GWA.json' % (flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('mapping'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name) + os.system(run_loco_command) + else: + gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.cXX.txt -lmm 1 -maf 0.1 -c %s/%s_covariates.txt -outdir %s -debug -o %s_output' % (flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('mapping'), + this_dataset.group.name, + webqtlConfig.GENERATED_IMAGE_DIR, + this_dataset.group.name) else: if method == "gemma": #gemma_command = GEMMA_COMMAND + ' -bfile %s/%s -k %s/%s.sXX.txt -lmm 1 -maf 0.1 -o %s_output' % (flat_files('mapping'), @@ -58,8 +83,28 @@ def run_gemma(this_dataset, samples, vals, covariates, method): webqtlConfig.GENERATED_IMAGE_DIR, this_dataset.group.name) else: - #gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.sXX.txt -lmm 1 -maf 0.1 -o %s_output' % (flat_files('genotype/bimbam'), - gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.cXX.txt -lmm 1 -maf 0.1 -outdir %s -debug -o %s_output' % (flat_files('genotype/bimbam'), + if use_loco: + generate_k_command = GEMMA_WRAPPER_COMMAND + ' --json --loco ' + chr_list_string + ' -- -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -gk -debug > /home/zas1024/tmp/gn2/K.json' % (flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name) + logger.debug("k_command:" + generate_k_command) + os.system(generate_k_command) + + run_loco_command = GEMMA_WRAPPER_COMMAND + ' --json --loco --input /home/zas1024/tmp/gn2/K.json -- -g %s/%s_geno.txt -p %s/%s_pheno.txt -c %s/%s_covariates.txt -a %s/%s_snps.txt -lmm 1 -maf 0.1 -debug > /home/zas1024/tmp/gn2/GWA.json' % (flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name, + flat_files('mapping'), + this_dataset.group.name, + flat_files('genotype/bimbam'), + this_dataset.group.name) + os.system(run_loco_command) + else: + #gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.sXX.txt -lmm 1 -maf 0.1 -o %s_output' % (flat_files('genotype/bimbam'), + gemma_command = GEMMA_COMMAND + ' -g %s/%s_geno.txt -p %s/%s_pheno.txt -a %s/%s_snps.txt -k %s/%s.cXX.txt -lmm 1 -maf 0.1 -outdir %s -debug -o %s_output' % (flat_files('genotype/bimbam'), this_dataset.group.name, flat_files('genotype/bimbam'), this_dataset.group.name, diff --git a/wqflask/wqflask/marker_regression/marker_regression.py b/wqflask/wqflask/marker_regression/marker_regression.py index 55bbacac..24292c35 100644 --- a/wqflask/wqflask/marker_regression/marker_regression.py +++ b/wqflask/wqflask/marker_regression/marker_regression.py @@ -84,6 +84,7 @@ class MarkerRegression(object): self.manhattan_plot = False self.maf = start_vars['maf'] # Minor allele frequency + self.use_loco = start_vars['use_loco'] self.suggestive = "" self.significant = "" self.pair_scan = False # Initializing this since it is checked in views to determine which template to use @@ -158,7 +159,7 @@ class MarkerRegression(object): self.score_type = "-log(p)" self.manhattan_plot = True with Bench("Running GEMMA"): - marker_obs = gemma_mapping.run_gemma(self.dataset, self.samples, self.vals, self.covariates, self.mapping_method) + marker_obs = gemma_mapping.run_gemma(self.dataset, self.samples, self.vals, self.covariates, self.mapping_method, self.use_loco) results = marker_obs elif self.mapping_method == "rqtl_plink": results = self.run_rqtl_plink() diff --git a/wqflask/wqflask/templates/show_trait_mapping_tools.html b/wqflask/wqflask/templates/show_trait_mapping_tools.html index c24a5853..eadc7481 100644 --- a/wqflask/wqflask/templates/show_trait_mapping_tools.html +++ b/wqflask/wqflask/templates/show_trait_mapping_tools.html @@ -302,8 +302,20 @@ <input name="maf_gemma" value="0.01" type="text" class="form-control"> </div> </div> + <div class="mapping_method_fields form-group"> + <label class="col-xs-4 control-label">Use LOCO</label> + <div style="margin-left: 20px;" class="col-xs-4 controls"> + <label class="radio-inline"> + <input type="radio" name="use_loco" value="True"> + Yes + </label> + <label class="radio-inline"> + <input type="radio" name="use_loco" value="False" checked=""> + No + </label> + </div> + </div> </div> - <div style="padding-top: 5px; padding-bottom: 5px; padding-left: 20px;" class="form-horizontal"> <div class="mapping_method_fields form-group"> <button type="button" id="select_covariates" class="btn btn-default"> @@ -395,7 +407,7 @@ <dt>GEMMA</dt> <dd> GEMMA is software implementing the Genome-wide Efficient Mixed Model Association algorithm for a standard linear mixed model for genome-wide association studies (GWAS).<br> - <font style="color: red;">Note: There currently exists an issue where GEMMA can't be run on traits from the same group simultaneously. If you recieve an error, please wait a few minutes and try again.</font> + Note: There currently exists an issue where GEMMA can't be run on traits from the same group simultaneously. If you receive an error, please wait a few minutes and try again. </dd> <dt>PLINK</dt> <dd>PLINK is a free, open-source whole genome association analysis toolset.</dd> diff --git a/wqflask/wqflask/views.py b/wqflask/wqflask/views.py index 2d4fd0f2..f340c954 100644 --- a/wqflask/wqflask/views.py +++ b/wqflask/wqflask/views.py @@ -516,6 +516,7 @@ def loading_page(): 'LRSCheck', 'covariates', 'maf', + 'use_loco', 'manhattan_plot', 'control_marker', 'control_marker_db', @@ -572,6 +573,7 @@ def marker_regression_page(): 'LRSCheck', 'covariates', 'maf', + 'use_loco', 'manhattan_plot', 'control_marker', 'control_marker_db', |