From 9a4433a7fe6bab8298f58c42f54731d7ade18d4e Mon Sep 17 00:00:00 2001 From: zsloan Date: Thu, 12 Mar 2020 12:53:20 -0500 Subject: Made some changes to how the mapping Y axis is determined/drawn; left in the commented-out old code until we're sure this works well --- .../marker_regression/display_mapping_results.py | 52 +++++++++++++++------- 1 file changed, 36 insertions(+), 16 deletions(-) diff --git a/wqflask/wqflask/marker_regression/display_mapping_results.py b/wqflask/wqflask/marker_regression/display_mapping_results.py index e74f35f5..66b5e25e 100644 --- a/wqflask/wqflask/marker_regression/display_mapping_results.py +++ b/wqflask/wqflask/marker_regression/display_mapping_results.py @@ -1654,10 +1654,6 @@ class DisplayMappingResults(object): INTERCROSS = (self.genotype.type=="intercross") - LRSHeightThresh = drawAreaHeight - AdditiveHeightThresh = drawAreaHeight/2 - DominanceHeightThresh = drawAreaHeight/2 - #draw the LRS scale #We first determine whether or not we are using a sliding scale. #If so, we need to compute the maximum LRS value to determine where the max y-value should be, and call this LRS_LOD_Max. @@ -1694,7 +1690,9 @@ class DisplayMappingResults(object): pass if self.permChecked and self.nperm > 0 and not self.multipleInterval: - LRS_LOD_Max = max(self.significant, LRS_LOD_Max) + if self.significant > LRS_LOD_Max: + LRS_LOD_Max = self.significant * 1.1 + #LRS_LOD_Max = max(self.significant, LRS_LOD_Max) else: LRS_LOD_Max = 1.15*LRS_LOD_Max @@ -1715,6 +1713,12 @@ class DisplayMappingResults(object): LRSLODFont=pid.Font(ttf="verdana", size=18*zoom*1.5, bold=0) yZero = yTopOffset + plotHeight + #LRSHeightThresh = drawAreaHeight + #AdditiveHeightThresh = drawAreaHeight/2 + #DominanceHeightThresh = drawAreaHeight/2 + LRSHeightThresh = (yZero - yTopOffset + 30*(zoom - 1)) + AdditiveHeightThresh = LRSHeightThresh/2 + DominanceHeightThresh = LRSHeightThresh/2 if LRS_LOD_Max > 100: LRSScale = 20.0 @@ -1728,7 +1732,7 @@ class DisplayMappingResults(object): LRSAxisList = Plot.frange(LRSScale, LRS_LOD_Max, LRSScale) #make sure the user's value appears on the y-axis #update by NL 6-21-2011: round the LOD value to 100 when LRS_LOD_Max is equal to 460 - LRSAxisList.append(round(LRS_LOD_Max)) + LRSAxisList.append(ceil(LRS_LOD_Max)) #ZS: Convert to int if all axis values are whole numbers all_int = True @@ -1751,7 +1755,9 @@ class DisplayMappingResults(object): for item in LRSAxisList: if LRS_LOD_Max == 0.0: LRS_LOD_Max = 0.000001 - yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh + yTopOffset + 30*(zoom - 1) + #yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh + yLRS = yZero - (item/LRSAxisList[-1]) * LRSHeightThresh canvas.drawLine(xLeftOffset, yLRS, xLeftOffset - 4, yLRS, color=self.LRS_COLOR, width=1*zoom) if all_int: scaleStr = "%d" % item @@ -1761,8 +1767,10 @@ class DisplayMappingResults(object): canvas.drawString(scaleStr, xLeftOffset-4-canvas.stringWidth(scaleStr, font=LRSScaleFont)-5, yLRS+3, font=LRSScaleFont, color=self.LRS_COLOR) if self.permChecked and self.nperm > 0 and not self.multipleInterval: - significantY = yZero - self.significant*LRSHeightThresh/LRS_LOD_Max - suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRS_LOD_Max + # significantY = yZero - self.significant*LRSHeightThresh/LRS_LOD_Max + # suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRS_LOD_Max + significantY = yZero - self.significant*LRSHeightThresh/LRSAxisList[-1] + suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRSAxisList[-1] startPosX = xLeftOffset #"Significant" and "Suggestive" Drawing Routine @@ -1878,25 +1886,37 @@ class DisplayMappingResults(object): # updated by NL 06-18-2011: # fix the over limit LRS graph issue since genotype trait may give infinite LRS; # for any lrs is over than 460(LRS max in this system), it will be reset to 460 + + yLRS = yZero - (item/LRSAxisList[-1]) * (yZero - yTopOffset + 30*(zoom - 1)) + #yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh + + if 'lrs_value' in qtlresult: if self.LRS_LOD == "LOD" or self.LRS_LOD == "-log(p)": if qtlresult['lrs_value'] > 460 or qtlresult['lrs_value']=='inf': - Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) + Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRSAxisList[-1]*self.LODFACTOR) + #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) else: - Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) + Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRSAxisList[-1]*self.LODFACTOR) + #Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) else: if qtlresult['lrs_value'] > 460 or qtlresult['lrs_value']=='inf': - Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max + Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRSAxisList[-1] + #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max else: - Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRS_LOD_Max + Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRSAxisList[-1] + #Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRS_LOD_Max else: if qtlresult['lod_score'] > 100 or qtlresult['lod_score']=='inf': - Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max + Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRSAxisList[-1] + #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max else: if self.LRS_LOD == "LRS": - Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRS_LOD_Max + Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRSAxisList[-1] + #Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRS_LOD_Max else: - Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRS_LOD_Max + Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRSAxisList[-1] + #Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRS_LOD_Max if self.manhattan_plot == True: if self.selectedChr == -1 and (previous_chr_as_int % 2 == 1): -- cgit v1.2.3 From bb1f854ca7acdd4b234c9744725fce50bf0ea2d8 Mon Sep 17 00:00:00 2001 From: zsloan Date: Thu, 12 Mar 2020 15:03:15 -0500 Subject: Added back a library that is apparently used for the Copy To Other Collection feature --- wqflask/wqflask/templates/collections/view.html | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/wqflask/wqflask/templates/collections/view.html b/wqflask/wqflask/templates/collections/view.html index e5188ed8..97e46756 100644 --- a/wqflask/wqflask/templates/collections/view.html +++ b/wqflask/wqflask/templates/collections/view.html @@ -26,7 +26,7 @@ @@ -163,6 +163,7 @@ {% block js %} + -- cgit v1.2.3 From f0f97202cf77a4459dd7ca525e94e34a56ddb147 Mon Sep 17 00:00:00 2001 From: zsloan Date: Sun, 15 Mar 2020 15:03:15 -0500 Subject: Fixed digits issue, though still need to figure out how to increase gap between y axis and title --- wqflask/wqflask/static/new/javascript/show_trait.js | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/wqflask/wqflask/static/new/javascript/show_trait.js b/wqflask/wqflask/static/new/javascript/show_trait.js index 819c4d12..0162f858 100644 --- a/wqflask/wqflask/static/new/javascript/show_trait.js +++ b/wqflask/wqflask/static/new/javascript/show_trait.js @@ -1017,10 +1017,12 @@ get_bar_range = function(sample_vals, sample_errors = null){ root.chart_range = get_bar_range(get_sample_vals(sample_lists[0]), get_sample_errors(sample_lists[0])[0]) val_range = root.chart_range[1] - root.chart_range[0] -if (val_range < 5){ - tick_digits = '.1f' +if (val_range < 0.05){ + tick_digits = '.3f' } else if (val_range < 0.5) { tick_digits = '.2f' +} else if (val_range < 5){ + tick_digits = '.1f' } else { tick_digits = 'f' } @@ -1045,7 +1047,7 @@ if (js_data.num_values < 256) { }, }, yaxis: { - title: js_data.unit_type, + title: "" + js_data.unit_type + "", range: root.chart_range, titlefont: { size: 16 -- cgit v1.2.3 From f97da31add711dcc0e75797b206821ca2d05c105 Mon Sep 17 00:00:00 2001 From: zsloan Date: Thu, 19 Mar 2020 11:46:49 -0500 Subject: Made a minor change just to help qtlreaper deal with genotype files with positions in bases instead of megabases --- wqflask/wqflask/marker_regression/qtlreaper_mapping.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/wqflask/wqflask/marker_regression/qtlreaper_mapping.py b/wqflask/wqflask/marker_regression/qtlreaper_mapping.py index 41764c1c..0c560582 100644 --- a/wqflask/wqflask/marker_regression/qtlreaper_mapping.py +++ b/wqflask/wqflask/marker_regression/qtlreaper_mapping.py @@ -97,7 +97,10 @@ def parse_reaper_output(gwa_filename, permu_filename, bootstrap_filename): except: marker['chr'] = line.split("\t")[2] marker['cM'] = float(line.split("\t")[3]) - marker['Mb'] = float(line.split("\t")[4]) + if float(line.split("\t")[4]) > 1000: + marker['Mb'] = float(line.split("\t")[4])/1000000 + else: + marker['Mb'] = float(line.split("\t")[4]) if float(line.split("\t")[7]) != 1: marker['p_value'] = float(line.split("\t")[7]) marker['lrs_value'] = float(line.split("\t")[5]) -- cgit v1.2.3 From 2511934d7259bd9f377ec9585bcf285a46c4de17 Mon Sep 17 00:00:00 2001 From: zsloan Date: Wed, 25 Mar 2020 13:53:20 -0500 Subject: Added fix for issue when creating collections (related to hmac import change) --- wqflask/wqflask/collect.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/wqflask/wqflask/collect.py b/wqflask/wqflask/collect.py index 2a94e63d..74eb869f 100644 --- a/wqflask/wqflask/collect.py +++ b/wqflask/wqflask/collect.py @@ -43,10 +43,10 @@ def process_traits(unprocessed_traits): unprocessed_traits = unprocessed_traits.split(",") traits = set() for trait in unprocessed_traits: - data, _separator, hmac = trait.rpartition(':') + data, _separator, the_hmac = trait.rpartition(':') data = data.strip() if g.user_session.logged_in: - assert hmac == hmac.data_hmac(data), "Data tampering?" + assert the_hmac == hmac.hmac_creation(data), "Data tampering?" traits.add(str(data)) return traits -- cgit v1.2.3 From 0e240cbfc6034b65b42c2d21f59ae2663dd2f6ec Mon Sep 17 00:00:00 2001 From: zsloan Date: Fri, 27 Mar 2020 16:32:08 -0500 Subject: Added integration for using RData files with R/qtl, plus some other changes to how it decides which mapping methods to display on the trait page --- wqflask/base/data_set.py | 11 ++- wqflask/utility/gen_geno_ob.py | 91 +++++++++++++++------- .../marker_regression/display_mapping_results.py | 86 ++++++++++---------- wqflask/wqflask/marker_regression/rqtl_mapping.py | 40 +++++++--- wqflask/wqflask/marker_regression/run_mapping.py | 3 +- wqflask/wqflask/templates/collections/view.html | 4 +- .../templates/show_trait_mapping_tools.html | 39 +++++----- 7 files changed, 168 insertions(+), 106 deletions(-) diff --git a/wqflask/base/data_set.py b/wqflask/base/data_set.py index 7fe9a8ac..abfdd277 100644 --- a/wqflask/base/data_set.py +++ b/wqflask/base/data_set.py @@ -307,9 +307,11 @@ class DatasetGroup(object): mapping_id = g.db.execute("select MappingMethodId from InbredSet where Name= '%s'" % self.name).fetchone()[0] if mapping_id == "1": - mapping_names = ["QTLReaper", "R/qtl"] + mapping_names = ["GEMMA", "QTLReaper", "R/qtl"] elif mapping_id == "2": - mapping_names = ["GEMMA"] + mapping_names = [] + elif mapping_id == "3": + mapping_names = ["R/qtl"] elif mapping_id == "4": mapping_names = ["GEMMA", "PLINK"] else: @@ -392,7 +394,10 @@ class DatasetGroup(object): # reaper barfs on unicode filenames, so here we ensure it's a string if self.genofile: - full_filename = str(locate(self.genofile, 'genotype')) + if "RData" in self.genofile: #ZS: This is a temporary fix; I need to change the way the JSON files that point to multiple genotype files are structured to point to other file types like RData + full_filename = str(locate(self.genofile.split(".")[0] + ".geno", 'genotype')) + else: + full_filename = str(locate(self.genofile, 'genotype')) else: full_filename = str(locate(self.name + '.geno', 'genotype')) diff --git a/wqflask/utility/gen_geno_ob.py b/wqflask/utility/gen_geno_ob.py index 44e2722f..aa5b27c4 100644 --- a/wqflask/utility/gen_geno_ob.py +++ b/wqflask/utility/gen_geno_ob.py @@ -1,5 +1,8 @@ from __future__ import absolute_import, division, print_function +import utility.logger +logger = utility.logger.getLogger(__name__ ) + class genotype(object): """ Replacement for reaper.Dataset so we can remove qtlreaper use while still generating mapping output figure @@ -34,8 +37,34 @@ class genotype(object): def __len__(self): return len(self.chromosomes) - def read_file(self, filename): + def read_rdata_output(self, qtl_results): + #ZS: This is necessary because R/qtl requires centimorgan marker positions, which it normally gets from the .geno file, but that doesn't exist for HET3-ITP (which only has RData), so it needs to read in the marker cM positions from the results + self.chromosomes = [] #ZS: Overwriting since the .geno file's contents are just placeholders + + this_chr = "" #ZS: This is so it can track when the chromosome changes as it iterates through markers + chr_ob = None + for marker in qtl_results: + locus = Locus(self) + if str(marker['chr']) != this_chr: + if this_chr != "": + self.chromosomes.append(chr_ob) + this_chr = str(marker['chr']) + chr_ob = Chr(this_chr, self) + if 'chr' in marker: + locus.chr = str(marker['chr']) + if 'name' in marker: + locus.name = marker['name'] + if 'Mb' in marker: + locus.Mb = marker['Mb'] + if 'cM' in marker: + locus.cM = marker['cM'] + chr_ob.loci.append(locus) + + self.chromosomes.append(chr_ob) + + return self + def read_file(self, filename): with open(filename, 'r') as geno_file: lines = geno_file.readlines() @@ -109,33 +138,39 @@ class Chr(object): return len(self.loci) def add_marker(self, marker_row): - self.loci.append(Locus(marker_row, self.geno_ob)) + self.loci.append(Locus(self.geno_ob, marker_row)) class Locus(object): - def __init__(self, marker_row, geno_ob): - self.chr = marker_row[0] - self.name = marker_row[1] - try: - self.cM = float(marker_row[geno_ob.cm_column]) - except: - self.cM = float(marker_row[geno_ob.mb_column]) if geno_ob.mb_exists else 0 - self.Mb = float(marker_row[geno_ob.mb_column]) if geno_ob.mb_exists else None - - geno_table = { - geno_ob.mat: -1, - geno_ob.pat: 1, - geno_ob.het: 0, - geno_ob.unk: "U" - } - + def __init__(self, geno_ob, marker_row = None): + self.chr = None + self.name = None + self.cM = None + self.Mb = None self.genotype = [] - if geno_ob.mb_exists: - start_pos = 4 - else: - start_pos = 3 - - for allele in marker_row[start_pos:]: - if allele in geno_table.keys(): - self.genotype.append(geno_table[allele]) - else: #ZS: Some genotype appears that isn't specified in the metadata, make it unknown - self.genotype.append("U") \ No newline at end of file + if marker_row: + self.chr = marker_row[0] + self.name = marker_row[1] + try: + self.cM = float(marker_row[geno_ob.cm_column]) + except: + self.cM = float(marker_row[geno_ob.mb_column]) if geno_ob.mb_exists else 0 + self.Mb = float(marker_row[geno_ob.mb_column]) if geno_ob.mb_exists else None + + geno_table = { + geno_ob.mat: -1, + geno_ob.pat: 1, + geno_ob.het: 0, + geno_ob.unk: "U" + } + + self.genotype = [] + if geno_ob.mb_exists: + start_pos = 4 + else: + start_pos = 3 + + for allele in marker_row[start_pos:]: + if allele in geno_table.keys(): + self.genotype.append(geno_table[allele]) + else: #ZS: Some genotype appears that isn't specified in the metadata, make it unknown + self.genotype.append("U") \ No newline at end of file diff --git a/wqflask/wqflask/marker_regression/display_mapping_results.py b/wqflask/wqflask/marker_regression/display_mapping_results.py index 66b5e25e..cf4508dd 100644 --- a/wqflask/wqflask/marker_regression/display_mapping_results.py +++ b/wqflask/wqflask/marker_regression/display_mapping_results.py @@ -247,11 +247,29 @@ class DisplayMappingResults(object): self.strainlist = start_vars['samples'] + self.traitList = [] + thisTrait = start_vars['this_trait'] + self.traitList.append(thisTrait) + + ################################################################ + # Calculations QTL goes here + ################################################################ + self.multipleInterval = len(self.traitList) > 1 + self.qtlresults = start_vars['qtl_results'] + + if self.multipleInterval: + self.colorCollection = Plot.colorSpectrum(len(self.qtlresults)) + else: + self.colorCollection = [self.LRS_COLOR] + if self.mapping_method == "reaper" and self.manhattan_plot != True: self.genotype = self.dataset.group.read_genotype_file(use_reaper=True) else: self.genotype = self.dataset.group.read_genotype_file() + if self.mapping_method == "rqtl_geno" and self.genotype.filler == True: + self.genotype = self.genotype.read_rdata_output(self.qtlresults) + #Darwing Options try: if self.selectedChr > -1: @@ -346,10 +364,6 @@ class DisplayMappingResults(object): else: self.GraphInterval = self.cMGraphInterval #cM - self.traitList = [] - thisTrait = start_vars['this_trait'] - self.traitList.append(thisTrait) - ## BEGIN HaplotypeAnalyst ## count the amount of individuals to be plotted, and increase self.graphHeight if self.haplotypeAnalystChecked and self.selectedChr > -1: @@ -371,16 +385,7 @@ class DisplayMappingResults(object): self.graphHeight = self.graphHeight + 2 * (self.NR_INDIVIDUALS+10) * self.EACH_GENE_HEIGHT ## END HaplotypeAnalyst - ################################################################ - # Calculations QTL goes here - ################################################################ - self.multipleInterval = len(self.traitList) > 1 - self.qtlresults = start_vars['qtl_results'] - if self.multipleInterval: - self.colorCollection = Plot.colorSpectrum(len(self.qtlresults)) - else: - self.colorCollection = [self.LRS_COLOR] ######################### @@ -1713,12 +1718,12 @@ class DisplayMappingResults(object): LRSLODFont=pid.Font(ttf="verdana", size=18*zoom*1.5, bold=0) yZero = yTopOffset + plotHeight - #LRSHeightThresh = drawAreaHeight - #AdditiveHeightThresh = drawAreaHeight/2 - #DominanceHeightThresh = drawAreaHeight/2 - LRSHeightThresh = (yZero - yTopOffset + 30*(zoom - 1)) - AdditiveHeightThresh = LRSHeightThresh/2 - DominanceHeightThresh = LRSHeightThresh/2 + LRSHeightThresh = drawAreaHeight + AdditiveHeightThresh = drawAreaHeight/2 + DominanceHeightThresh = drawAreaHeight/2 + # LRSHeightThresh = (yZero - yTopOffset + 30*(zoom - 1)) + # AdditiveHeightThresh = LRSHeightThresh/2 + # DominanceHeightThresh = LRSHeightThresh/2 if LRS_LOD_Max > 100: LRSScale = 20.0 @@ -1756,8 +1761,8 @@ class DisplayMappingResults(object): if LRS_LOD_Max == 0.0: LRS_LOD_Max = 0.000001 yTopOffset + 30*(zoom - 1) - #yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh - yLRS = yZero - (item/LRSAxisList[-1]) * LRSHeightThresh + yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh + #yLRS = yZero - (item/LRSAxisList[-1]) * LRSHeightThresh canvas.drawLine(xLeftOffset, yLRS, xLeftOffset - 4, yLRS, color=self.LRS_COLOR, width=1*zoom) if all_int: scaleStr = "%d" % item @@ -1767,10 +1772,10 @@ class DisplayMappingResults(object): canvas.drawString(scaleStr, xLeftOffset-4-canvas.stringWidth(scaleStr, font=LRSScaleFont)-5, yLRS+3, font=LRSScaleFont, color=self.LRS_COLOR) if self.permChecked and self.nperm > 0 and not self.multipleInterval: - # significantY = yZero - self.significant*LRSHeightThresh/LRS_LOD_Max - # suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRS_LOD_Max - significantY = yZero - self.significant*LRSHeightThresh/LRSAxisList[-1] - suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRSAxisList[-1] + significantY = yZero - self.significant*LRSHeightThresh/LRS_LOD_Max + suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRS_LOD_Max + # significantY = yZero - self.significant*LRSHeightThresh/LRSAxisList[-1] + # suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRSAxisList[-1] startPosX = xLeftOffset #"Significant" and "Suggestive" Drawing Routine @@ -1887,36 +1892,35 @@ class DisplayMappingResults(object): # fix the over limit LRS graph issue since genotype trait may give infinite LRS; # for any lrs is over than 460(LRS max in this system), it will be reset to 460 - yLRS = yZero - (item/LRSAxisList[-1]) * (yZero - yTopOffset + 30*(zoom - 1)) - #yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh + yLRS = yZero - (item/LRS_LOD_Max) * LRSHeightThresh if 'lrs_value' in qtlresult: if self.LRS_LOD == "LOD" or self.LRS_LOD == "-log(p)": if qtlresult['lrs_value'] > 460 or qtlresult['lrs_value']=='inf': - Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRSAxisList[-1]*self.LODFACTOR) - #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) + #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRSAxisList[-1]*self.LODFACTOR) + Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) else: - Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRSAxisList[-1]*self.LODFACTOR) - #Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) + #Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRSAxisList[-1]*self.LODFACTOR) + Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/(LRS_LOD_Max*self.LODFACTOR) else: if qtlresult['lrs_value'] > 460 or qtlresult['lrs_value']=='inf': - Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRSAxisList[-1] - #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max + #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRSAxisList[-1] + Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max else: - Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRSAxisList[-1] - #Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRS_LOD_Max + #Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRSAxisList[-1] + Yc = yZero - qtlresult['lrs_value']*LRSHeightThresh/LRS_LOD_Max else: if qtlresult['lod_score'] > 100 or qtlresult['lod_score']=='inf': - Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRSAxisList[-1] - #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max + #Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRSAxisList[-1] + Yc = yZero - webqtlConfig.MAXLRS*LRSHeightThresh/LRS_LOD_Max else: if self.LRS_LOD == "LRS": - Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRSAxisList[-1] - #Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRS_LOD_Max + #Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRSAxisList[-1] + Yc = yZero - qtlresult['lod_score']*self.LODFACTOR*LRSHeightThresh/LRS_LOD_Max else: - Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRSAxisList[-1] - #Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRS_LOD_Max + #Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRSAxisList[-1] + Yc = yZero - qtlresult['lod_score']*LRSHeightThresh/LRS_LOD_Max if self.manhattan_plot == True: if self.selectedChr == -1 and (previous_chr_as_int % 2 == 1): diff --git a/wqflask/wqflask/marker_regression/rqtl_mapping.py b/wqflask/wqflask/marker_regression/rqtl_mapping.py index 41d67012..2e3ea406 100644 --- a/wqflask/wqflask/marker_regression/rqtl_mapping.py +++ b/wqflask/wqflask/marker_regression/rqtl_mapping.py @@ -9,8 +9,6 @@ import utility.logger logger = utility.logger.getLogger(__name__ ) def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, control_marker, manhattan_plot, pair_scan): - geno_to_rqtl_function(dataset) - ## Get pointers to some common R functions r_library = ro.r["library"] # Map the library function r_c = ro.r["c"] # Map the c function @@ -21,16 +19,24 @@ def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, print(r_library("qtl")) # Load R/qtl ## Get pointers to some R/qtl functions - scanone = ro.r["scanone"] # Map the scanone function - scantwo = ro.r["scantwo"] # Map the scantwo function - calc_genoprob = ro.r["calc.genoprob"] # Map the calc.genoprob function - GENOtoCSVR = ro.r["GENOtoCSVR"] # Map the local GENOtoCSVR function + scanone = ro.r["scanone"] # Map the scanone function + scantwo = ro.r["scantwo"] # Map the scantwo function + calc_genoprob = ro.r["calc.genoprob"] # Map the calc.genoprob function crossname = dataset.group.name - genofilelocation = locate(crossname + ".geno", "genotype") - crossfilelocation = TMPDIR + crossname + ".cross" - - cross_object = GENOtoCSVR(genofilelocation, crossfilelocation) # TODO: Add the SEX if that is available + try: + generate_cross_from_rdata(dataset) + read_cross_from_rdata = ro.r["generate_cross_from_rdata"] # Map the local read_cross_from_rdata function + genofilelocation = locate(crossname + ".RData", "genotype/rdata") + cross_object = read_cross_from_rdata(genofilelocation) # Map the local GENOtoCSVR function + except: + generate_cross_from_geno(dataset) + GENOtoCSVR = ro.r["GENOtoCSVR"] # Map the local GENOtoCSVR function + crossfilelocation = TMPDIR + crossname + ".cross" + genofilelocation = locate(crossname + ".geno", "genotype") + + GENOtoCSVR = ro.r["GENOtoCSVR"] # Map the local GENOtoCSVR function + cross_object = GENOtoCSVR(genofilelocation, crossfilelocation) # TODO: Add the SEX if that is available if manhattan_plot: cross_object = calc_genoprob(cross_object) @@ -71,7 +77,18 @@ def run_rqtl_geno(vals, dataset, method, model, permCheck, num_perm, do_control, else: return process_rqtl_results(result_data_frame) -def geno_to_rqtl_function(dataset): # TODO: Need to figure out why some genofiles have the wrong format and don't convert properly +def generate_cross_from_rdata(dataset): + rdata_location = locate(dataset.group.name + ".RData", "genotype/rdata") + ro.r(""" + generate_cross_from_rdata <- function(filename = '%s') { + load(file=filename) + cross = cunique + return(cross) + } + """ % (rdata_location)) + +def generate_cross_from_geno(dataset): # TODO: Need to figure out why some genofiles have the wrong format and don't convert properly + ro.r(""" trim <- function( x ) { gsub("(^[[:space:]]+|[[:space:]]+$)", "", x) } @@ -170,6 +187,7 @@ def process_rqtl_results(result): # TODO: how to make this a one liner an marker = {} marker['name'] = result.rownames[i] marker['chr'] = output[i][0] + marker['cM'] = output[i][1] marker['Mb'] = output[i][1] marker['lod_score'] = output[i][2] qtl_results.append(marker) diff --git a/wqflask/wqflask/marker_regression/run_mapping.py b/wqflask/wqflask/marker_regression/run_mapping.py index ba41ffc3..f03b046e 100644 --- a/wqflask/wqflask/marker_regression/run_mapping.py +++ b/wqflask/wqflask/marker_regression/run_mapping.py @@ -81,7 +81,6 @@ class RunMapping(object): self.vals.append(value) else: self.samples = [] - for sample in self.dataset.group.samplelist: # sample is actually the name of an individual if (len(genofile_samplelist) == 0) or (sample in genofile_samplelist): in_trait_data = False @@ -186,7 +185,7 @@ class RunMapping(object): self.showGenes = "ON" self.viewLegend = "ON" - self.dataset.group.get_markers() + #self.dataset.group.get_markers() if self.mapping_method == "gemma": self.first_run = True self.output_files = None diff --git a/wqflask/wqflask/templates/collections/view.html b/wqflask/wqflask/templates/collections/view.html index 97e46756..1be6539d 100644 --- a/wqflask/wqflask/templates/collections/view.html +++ b/wqflask/wqflask/templates/collections/view.html @@ -161,10 +161,10 @@ {% endblock %} {% block js %} + + - - diff --git a/wqflask/wqflask/templates/show_trait_mapping_tools.html b/wqflask/wqflask/templates/show_trait_mapping_tools.html index a806a8b3..ad7412b2 100644 --- a/wqflask/wqflask/templates/show_trait_mapping_tools.html +++ b/wqflask/wqflask/templates/show_trait_mapping_tools.html @@ -4,25 +4,18 @@
-
+ {% for mapping_method in dataset.group.mapping_names %} + {% if mapping_method == "GEMMA" %} +
@@ -102,7 +97,7 @@
- {% if dataset.group.mapping_id == "1" %} + {% elif mapping_method == "QTLReaper" %}
@@ -217,7 +212,8 @@
-
+ {% elif mapping_method == "R/qtl" %} +
@@ -341,19 +337,24 @@
{% endif %} + {% endfor %}
+ {% for mapping_method in dataset.group.mapping_names %} + {% if mapping_method == "GEMMA" %}
GEMMA
Maps traits with correction for kinship among samples using a linear mixed model method, and also allows users to fit multiple covariates such as sex, age, treatment, and genetic markers (PMID: 2453419, and GitHub code). GEMMA incorporates the Leave One Chromosome Out (LOCO) method to ensure that the correction for kinship does not remove useful genetic variance near each marker. Markers can be filtered to include only those with minor allele frequencies (MAF) above a threshold. The default MAF is 0.05.
- {% if dataset.group.mapping_id == "1" %} + {% elif mapping_method == "R/qtl" %}
R/qtl
Major upgrade of R/qtl that supports most experimental populations including those with complex admixture and two or more parental lines as well as large omic data sets (PMID: 30591514). Both R/qtl and R/qtl2 are available as stand-alone R packages (R suite).
+ {% elif mapping_method == "QTLReaper" %}
Haley-Knott Regression
Fast linear mapping method (PMID 16718932) works well with F2 intercrosses and backcrosses, but that is not recommended for complex or admixed populations (e.g., GWAS or heterogeneous stock studies) or for advanced intercrosses, recombinant inbred families, or diallel crosses. Interactive plots in GeneNetwork have relied on the fast HK mapping for two decades and we still use this method for mapping omics data sets and computing genome-wide permutation threshold (QTL Reaper code).
{% endif %} + {% endfor %}