diff options
| -rw-r--r--[-rwxr-xr-x] | wqflask/utility/Plot.py | 352 | ||||
| -rw-r--r-- | wqflask/wqflask/marker_regression/marker_regression.py | 90 | ||||
| -rw-r--r-- | wqflask/wqflask/marker_regression/marker_regression_gn1.py | 156 | ||||
| -rwxr-xr-x | wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js | 3 | ||||
| -rw-r--r-- | wqflask/wqflask/templates/marker_regression_gn1.html | 17 | ||||
| -rw-r--r-- | wqflask/wqflask/views.py | 2 |
6 files changed, 325 insertions, 295 deletions
diff --git a/wqflask/utility/Plot.py b/wqflask/utility/Plot.py index 51a57a6d..ad11a81e 100755..100644 --- a/wqflask/utility/Plot.py +++ b/wqflask/utility/Plot.py @@ -24,10 +24,9 @@ # # Last updated by GeneNetwork Core Team 2010/10/20 -#import piddle as pid - from __future__ import print_function +import piddle as pid from pprint import pformat as pf print("Lysol") @@ -478,181 +477,180 @@ def plotSecurity(canvas, text="12345"): # parameter: data is either object returned by reaper permutation function (called by MarkerRegressionPage.py) # or the first object returned by direct (pair-scan) permu function (called by DirectPlotPage.py) -#def plotBar(canvas, data, barColor=pid.blue, axesColor=pid.black, labelColor=pid.black, XLabel=None, YLabel=None, title=None, offset= (60, 20, 40, 40), zoom = 1): -# -# xLeftOffset, xRightOffset, yTopOffset, yBottomOffset = offset -# -# plotWidth = canvas.size[0] - xLeftOffset - xRightOffset -# plotHeight = canvas.size[1] - yTopOffset - yBottomOffset -# if plotHeight<=0 or plotWidth<=0: -# return -# -# if len(data) < 2: -# return -# -# max_D = max(data) -# min_D = min(data) -# #add by NL 06-20-2011: fix the error: when max_D is infinite, log function in detScale will go wrong -# if max_D == float('inf') or max_D>webqtlConfig.MAXLRS: -# max_D=webqtlConfig.MAXLRS #maximum LRS value -# -# xLow, xTop, stepX = detScale(min_D, max_D) -# -# #reduce data -# step = ceil((xTop-xLow)/50.0) -# j = xLow -# dataXY = [] -# Count = [] -# while j <= xTop: -# dataXY.append(j) -# Count.append(0) -# j += step -# -# for i, item in enumerate(data): -# if item == float('inf') or item>webqtlConfig.MAXLRS: -# item = webqtlConfig.MAXLRS #maximum LRS value -# j = int((item-xLow)/step) -# Count[j] += 1 -# -# yLow, yTop, stepY=detScale(0,max(Count)) -# -# #draw data -# xScale = plotWidth/(xTop-xLow) -# yScale = plotHeight/(yTop-yLow) -# barWidth = xScale*step -# -# for i, count in enumerate(Count): -# if count: -# xc = (dataXY[i]-xLow)*xScale+xLeftOffset -# yc =-(count-yLow)*yScale+yTopOffset+plotHeight -# canvas.drawRect(xc+2,yc,xc+barWidth-2,yTopOffset+plotHeight,edgeColor=barColor,fillColor=barColor) -# -# #draw drawing region -# canvas.drawRect(xLeftOffset, yTopOffset, xLeftOffset+plotWidth, yTopOffset+plotHeight) -# -# #draw scale -# scaleFont=pid.Font(ttf="cour",size=11,bold=1) -# x=xLow -# for i in range(stepX+1): -# xc=xLeftOffset+(x-xLow)*xScale -# canvas.drawLine(xc,yTopOffset+plotHeight,xc,yTopOffset+plotHeight+5, color=axesColor) -# strX = cformat(d=x, rank=0) -# canvas.drawString(strX,xc-canvas.stringWidth(strX,font=scaleFont)/2,yTopOffset+plotHeight+14,font=scaleFont) -# x+= (xTop - xLow)/stepX -# -# y=yLow -# for i in range(stepY+1): -# yc=yTopOffset+plotHeight-(y-yLow)*yScale -# canvas.drawLine(xLeftOffset,yc,xLeftOffset-5,yc, color=axesColor) -# strY = "%d" %y -# canvas.drawString(strY,xLeftOffset-canvas.stringWidth(strY,font=scaleFont)-6,yc+5,font=scaleFont) -# y+= (yTop - yLow)/stepY -# -# #draw label -# labelFont=pid.Font(ttf="tahoma",size=17,bold=0) -# if XLabel: -# canvas.drawString(XLabel,xLeftOffset+(plotWidth-canvas.stringWidth(XLabel,font=labelFont))/2.0, -# yTopOffset+plotHeight+yBottomOffset-10,font=labelFont,color=labelColor) -# -# if YLabel: -# canvas.drawString(YLabel, 19, yTopOffset+plotHeight-(plotHeight-canvas.stringWidth(YLabel,font=labelFont))/2.0, -# font=labelFont,color=labelColor,angle=90) -# -# labelFont=pid.Font(ttf="verdana",size=16,bold=0) -# if title: -# canvas.drawString(title,xLeftOffset+(plotWidth-canvas.stringWidth(title,font=labelFont))/2.0, -# 20,font=labelFont,color=labelColor) -# -#def plotBarText(canvas, data, label, variance=None, barColor=pid.blue, axesColor=pid.black, labelColor=pid.black, XLabel=None, YLabel=None, title=None, sLabel = None, offset= (80, 20, 40, 100), barSpace = 2, zoom = 1): -# xLeftOffset, xRightOffset, yTopOffset, yBottomOffset = offset -# plotWidth = canvas.size[0] - xLeftOffset - xRightOffset -# plotHeight = canvas.size[1] - yTopOffset - yBottomOffset -# if plotHeight<=0 or plotWidth<=0: -# return -# -# NNN = len(data) -# if NNN < 2 or NNN != len(label): -# return -# if variance and len(variance)!=NNN: -# variance = [] -# -# Y2 = data[:] -# if variance: -# for i in range(NNN): -# if variance[i]: -# Y2 += [data[i]-variance[i]] -# -# #Y axis -# YLow, YTop, stepY = detScale(min(Y2), max(Y2)) -# YScale = plotHeight/(YTop - YLow) -# -# if YLow < 0 and YTop > 0: -# drawZero = 1 -# else: -# drawZero = 0 -# -# #X axis -# X = range(NNN) -# Xll= 0 -# Xur= NNN-1 -# -# -# if drawZero: -# YZero = yTopOffset+plotHeight-YScale*(0-YLow) -# canvas.drawLine(xLeftOffset, YZero, xLeftOffset+plotWidth, YZero) -# else: -# YZero = yTopOffset+plotHeight -# #draw data -# spaceWidth = barSpace -# if spaceWidth < 1: -# spaceWidth = 1 -# barWidth = int((plotWidth - (NNN-1.0)*spaceWidth)/NNN) -# -# xc= xLeftOffset -# scaleFont=pid.Font(ttf="verdana",size=11,bold=0) -# for i in range(NNN): -# yc = yTopOffset+plotHeight-(data[i]-YLow)*YScale -# canvas.drawRect(xc,YZero,xc+barWidth-1, yc, edgeColor=barColor,fillColor=barColor) -# if variance and variance[i]: -# varlen = variance[i]*YScale -# if yc-varlen < yTopOffset: -# topYd = yTopOffset -# else: -# topYd = yc-varlen -# canvas.drawLine(xc+barWidth/2-2,yc-varlen,xc+barWidth/2+2,yc-varlen,color=pid.red) -# canvas.drawLine(xc+barWidth/2,yc+varlen,xc+barWidth/2,topYd,color=pid.red) -# canvas.drawLine(xc+barWidth/2-2,yc+varlen,xc+barWidth/2+2,yc+varlen,color=pid.red) -# strX = label[i] -# canvas.drawString(strX,xc+barWidth/2.0+2,yTopOffset+plotHeight+2+canvas.stringWidth(strX,font=scaleFont),font=scaleFont,angle=90) -# xc += barWidth + spaceWidth -# -# #draw drawing region -# canvas.drawRect(xLeftOffset, yTopOffset, xLeftOffset+plotWidth, yTopOffset+plotHeight) -# -# #draw Y scale -# scaleFont=pid.Font(ttf="cour",size=16,bold=1) -# y=YLow -# for i in range(stepY+1): -# yc=yTopOffset+plotHeight-(y-YLow)*YScale -# canvas.drawLine(xLeftOffset,yc,xLeftOffset-5,yc, color=axesColor) -# strY = cformat(d=y, rank=0) -# canvas.drawString(strY,xLeftOffset-canvas.stringWidth(strY,font=scaleFont)-6,yc+5,font=scaleFont) -# y+= (YTop - YLow)/stepY -# -# #draw label -# labelFont=pid.Font(ttf="verdana",size=17,bold=0) -# if XLabel: -# canvas.drawString(XLabel,xLeftOffset+(plotWidth-canvas.stringWidth(XLabel,font=labelFont))/2.0,yTopOffset+plotHeight+65,font=labelFont,color=labelColor) -# -# if YLabel: -# canvas.drawString(YLabel,xLeftOffset-50, yTopOffset+plotHeight-(plotHeight-canvas.stringWidth(YLabel,font=labelFont))/2.0,font=labelFont,color=labelColor,angle=90) -# -# labelFont=pid.Font(ttf="verdana",size=18,bold=0) -# if title: -# canvas.drawString(title,xLeftOffset,yTopOffset-15,font=labelFont,color=labelColor) -# -# return -# +def plotBar(canvas, data, barColor=pid.blue, axesColor=pid.black, labelColor=pid.black, XLabel=None, YLabel=None, title=None, offset= (60, 20, 40, 40), zoom = 1): + xLeftOffset, xRightOffset, yTopOffset, yBottomOffset = offset + + plotWidth = canvas.size[0] - xLeftOffset - xRightOffset + plotHeight = canvas.size[1] - yTopOffset - yBottomOffset + if plotHeight<=0 or plotWidth<=0: + return + + if len(data) < 2: + return + + max_D = max(data) + min_D = min(data) + #add by NL 06-20-2011: fix the error: when max_D is infinite, log function in detScale will go wrong + if max_D == float('inf') or max_D>webqtlConfig.MAXLRS: + max_D=webqtlConfig.MAXLRS #maximum LRS value + + xLow, xTop, stepX = detScale(min_D, max_D) + + #reduce data + step = ceil((xTop-xLow)/50.0) + j = xLow + dataXY = [] + Count = [] + while j <= xTop: + dataXY.append(j) + Count.append(0) + j += step + + for i, item in enumerate(data): + if item == float('inf') or item>webqtlConfig.MAXLRS: + item = webqtlConfig.MAXLRS #maximum LRS value + j = int((item-xLow)/step) + Count[j] += 1 + + yLow, yTop, stepY=detScale(0,max(Count)) + + #draw data + xScale = plotWidth/(xTop-xLow) + yScale = plotHeight/(yTop-yLow) + barWidth = xScale*step + + for i, count in enumerate(Count): + if count: + xc = (dataXY[i]-xLow)*xScale+xLeftOffset + yc =-(count-yLow)*yScale+yTopOffset+plotHeight + canvas.drawRect(xc+2,yc,xc+barWidth-2,yTopOffset+plotHeight,edgeColor=barColor,fillColor=barColor) + + #draw drawing region + canvas.drawRect(xLeftOffset, yTopOffset, xLeftOffset+plotWidth, yTopOffset+plotHeight) + + #draw scale + scaleFont=pid.Font(ttf="cour",size=11,bold=1) + x=xLow + for i in range(int(stepX)+1): + xc=xLeftOffset+(x-xLow)*xScale + canvas.drawLine(xc,yTopOffset+plotHeight,xc,yTopOffset+plotHeight+5, color=axesColor) + strX = cformat(d=x, rank=0) + canvas.drawString(strX,xc-canvas.stringWidth(strX,font=scaleFont)/2,yTopOffset+plotHeight+14,font=scaleFont) + x+= (xTop - xLow)/stepX + + y=yLow + for i in range(int(stepY)+1): + yc=yTopOffset+plotHeight-(y-yLow)*yScale + canvas.drawLine(xLeftOffset,yc,xLeftOffset-5,yc, color=axesColor) + strY = "%d" %y + canvas.drawString(strY,xLeftOffset-canvas.stringWidth(strY,font=scaleFont)-6,yc+5,font=scaleFont) + y+= (yTop - yLow)/stepY + + #draw label + labelFont=pid.Font(ttf="tahoma",size=17,bold=0) + if XLabel: + canvas.drawString(XLabel,xLeftOffset+(plotWidth-canvas.stringWidth(XLabel,font=labelFont))/2.0, + yTopOffset+plotHeight+yBottomOffset-10,font=labelFont,color=labelColor) + + if YLabel: + canvas.drawString(YLabel, 19, yTopOffset+plotHeight-(plotHeight-canvas.stringWidth(YLabel,font=labelFont))/2.0, + font=labelFont,color=labelColor,angle=90) + + labelFont=pid.Font(ttf="verdana",size=16,bold=0) + if title: + canvas.drawString(title,xLeftOffset+(plotWidth-canvas.stringWidth(title,font=labelFont))/2.0, + 20,font=labelFont,color=labelColor) + +def plotBarText(canvas, data, label, variance=None, barColor=pid.blue, axesColor=pid.black, labelColor=pid.black, XLabel=None, YLabel=None, title=None, sLabel = None, offset= (80, 20, 40, 100), barSpace = 2, zoom = 1): + xLeftOffset, xRightOffset, yTopOffset, yBottomOffset = offset + plotWidth = canvas.size[0] - xLeftOffset - xRightOffset + plotHeight = canvas.size[1] - yTopOffset - yBottomOffset + if plotHeight<=0 or plotWidth<=0: + return + + NNN = len(data) + if NNN < 2 or NNN != len(label): + return + if variance and len(variance)!=NNN: + variance = [] + + Y2 = data[:] + if variance: + for i in range(NNN): + if variance[i]: + Y2 += [data[i]-variance[i]] + + #Y axis + YLow, YTop, stepY = detScale(min(Y2), max(Y2)) + YScale = plotHeight/(YTop - YLow) + + if YLow < 0 and YTop > 0: + drawZero = 1 + else: + drawZero = 0 + + #X axis + X = range(NNN) + Xll= 0 + Xur= NNN-1 + + + if drawZero: + YZero = yTopOffset+plotHeight-YScale*(0-YLow) + canvas.drawLine(xLeftOffset, YZero, xLeftOffset+plotWidth, YZero) + else: + YZero = yTopOffset+plotHeight + #draw data + spaceWidth = barSpace + if spaceWidth < 1: + spaceWidth = 1 + barWidth = int((plotWidth - (NNN-1.0)*spaceWidth)/NNN) + + xc= xLeftOffset + scaleFont=pid.Font(ttf="verdana",size=11,bold=0) + for i in range(NNN): + yc = yTopOffset+plotHeight-(data[i]-YLow)*YScale + canvas.drawRect(xc,YZero,xc+barWidth-1, yc, edgeColor=barColor,fillColor=barColor) + if variance and variance[i]: + varlen = variance[i]*YScale + if yc-varlen < yTopOffset: + topYd = yTopOffset + else: + topYd = yc-varlen + canvas.drawLine(xc+barWidth/2-2,yc-varlen,xc+barWidth/2+2,yc-varlen,color=pid.red) + canvas.drawLine(xc+barWidth/2,yc+varlen,xc+barWidth/2,topYd,color=pid.red) + canvas.drawLine(xc+barWidth/2-2,yc+varlen,xc+barWidth/2+2,yc+varlen,color=pid.red) + strX = label[i] + canvas.drawString(strX,xc+barWidth/2.0+2,yTopOffset+plotHeight+2+canvas.stringWidth(strX,font=scaleFont),font=scaleFont,angle=90) + xc += barWidth + spaceWidth + + #draw drawing region + canvas.drawRect(xLeftOffset, yTopOffset, xLeftOffset+plotWidth, yTopOffset+plotHeight) + + #draw Y scale + scaleFont=pid.Font(ttf="cour",size=16,bold=1) + y=YLow + for i in range(stepY+1): + yc=yTopOffset+plotHeight-(y-YLow)*YScale + canvas.drawLine(xLeftOffset,yc,xLeftOffset-5,yc, color=axesColor) + strY = cformat(d=y, rank=0) + canvas.drawString(strY,xLeftOffset-canvas.stringWidth(strY,font=scaleFont)-6,yc+5,font=scaleFont) + y+= (YTop - YLow)/stepY + + #draw label + labelFont=pid.Font(ttf="verdana",size=17,bold=0) + if XLabel: + canvas.drawString(XLabel,xLeftOffset+(plotWidth-canvas.stringWidth(XLabel,font=labelFont))/2.0,yTopOffset+plotHeight+65,font=labelFont,color=labelColor) + + if YLabel: + canvas.drawString(YLabel,xLeftOffset-50, yTopOffset+plotHeight-(plotHeight-canvas.stringWidth(YLabel,font=labelFont))/2.0,font=labelFont,color=labelColor,angle=90) + + labelFont=pid.Font(ttf="verdana",size=18,bold=0) + if title: + canvas.drawString(title,xLeftOffset,yTopOffset-15,font=labelFont,color=labelColor) + + return + #def plotXY(canvas, dataX, dataY, rank=0, dataLabel=[], plotColor = pid.black, axesColor=pid.black, labelColor=pid.black, lineSize="thin", lineColor=pid.grey, idFont="arial", idColor=pid.blue, idSize="14", symbolColor=pid.black, symbolType="circle", filled="yes", symbolSize="tiny", XLabel=None, YLabel=None, title=None, fitcurve=None, connectdot=1, displayR=None, loadingPlot = 0, offset= (80, 20, 40, 60), zoom = 1, specialCases=[], showLabel = 1, bufferSpace = 15): # 'displayR : correlation scatter plot, loadings : loading plot' # diff --git a/wqflask/wqflask/marker_regression/marker_regression.py b/wqflask/wqflask/marker_regression/marker_regression.py index f08e47ca..112dcd5c 100644 --- a/wqflask/wqflask/marker_regression/marker_regression.py +++ b/wqflask/wqflask/marker_regression/marker_regression.py @@ -79,10 +79,8 @@ class MarkerRegression(object): self.pair_scan = False # Initializing this since it is checked in views to determine which template to use self.score_type = "LRS" #ZS: LRS or LOD self.mapping_scale = "physic" - self.num_perm = 0 self.bootstrap_results = [] - #ZS: This is passed to GN1 code for single chr mapping self.selected_chr = -1 if "selected_chr" in start_vars: @@ -99,6 +97,12 @@ class MarkerRegression(object): if "haplotypeAnalystCheck" in start_vars: self.haplotypeAnalystCheck = start_vars['haplotypeAnalystCheck'] if "startMb" in start_vars: #ZS: This is to ensure showGenes, Legend, etc are checked the first time you open the mapping page, since startMb will only not be set during the first load + if "permCheck" in start_vars: + self.permCheck = "ON" + else: + self.permCheck = False + self.num_perm = int(start_vars['num_perm']) + if "showGenes" in start_vars: self.showGenes = start_vars['showGenes'] else: @@ -108,6 +112,15 @@ class MarkerRegression(object): else: self.viewLegend = False else: + try: + if int(start_vars['num_perm']) > 0: + self.num_perm = int(start_vars['num_perm']) + else: + self.num_perm = 0 + except: + self.num_perm = 0 + + self.permCheck = "ON" self.showGenes = "ON" self.viewLegend = "ON" @@ -120,45 +133,32 @@ class MarkerRegression(object): with Bench("Getting markers from csv"): marker_obs = get_markers_from_csv(included_markers, p_values, self.dataset.group.name) results = marker_obs - #self.dataset.group.get_specified_markers(markers = included_markers) - #self.dataset.group.markers.add_pvalues(p_values) - #results = self.dataset.group.markers.markers elif self.mapping_method == "rqtl_plink": results = self.run_rqtl_plink() elif self.mapping_method == "rqtl_geno": self.score_type = "LOD" self.mapping_scale = "morgan" - if start_vars['num_perm'] == "": - self.num_perm = 0 - else: - self.num_perm = start_vars['num_perm'] self.control_marker = start_vars['control_marker'] self.do_control = start_vars['do_control'] self.method = start_vars['mapmethod_rqtl_geno'] self.model = start_vars['mapmodel_rqtl_geno'] - if start_vars['pair_scan'] == "true": self.pair_scan = True - results = self.run_rqtl_geno() - elif self.mapping_method == "reaper": - if start_vars['num_perm'] == "": - self.num_perm = 0 - else: - self.num_perm = int(start_vars['num_perm']) - + elif self.mapping_method == "reaper": if "startMb" in start_vars: #ZS: Check if first time page loaded, so it can default to ON - if "bootCheck" in start_vars: - self.bootCheck = "ON" - self.num_bootstrap = int(start_vars['num_bootstrap']) - else: - self.bootCheck = False - self.num_bootstrap = int(start_vars['num_bootstrap']) if "additiveCheck" in start_vars: self.additiveCheck = start_vars['additiveCheck'] else: self.additiveCheck = False + + if "bootCheck" in start_vars: + self.bootCheck = "ON" + else: + self.bootCheck = False + self.num_bootstrap = int(start_vars['num_bootstrap']) else: + self.additiveCheck = "ON" try: if int(start_vars['num_bootstrap']) > 0: self.bootCheck = "ON" @@ -166,12 +166,9 @@ class MarkerRegression(object): else: self.bootCheck = False self.num_bootstrap = 0 - #ZS: If some string that can't be converted to int is input for num_bootstrap except: - self.num_bootstrap = 0 self.bootCheck = False - - self.additiveCheck = "ON" + self.num_bootstrap = 0 self.control_marker = start_vars['control_marker'] self.do_control = start_vars['do_control'] @@ -180,10 +177,8 @@ class MarkerRegression(object): results = self.run_plink() elif self.mapping_method == "pylmm": print("RUNNING PYLMM") - self.num_perm = start_vars['num_perm'] - if self.num_perm != "": - if int(self.num_perm) > 0: - self.run_permutations(str(temp_uuid)) + if self.num_perm > 0: + self.run_permutations(str(temp_uuid)) results = self.gen_data(str(temp_uuid)) else: print("RUNNING NOTHING") @@ -449,11 +444,11 @@ class MarkerRegression(object): else: print("No covariates"); result_data_frame = scanone(cross_object, pheno = "the_pheno", model=self.model, method=self.method) - if int(self.num_perm) > 0: # Do permutation (if requested by user) + if self.num_perm > 0 and self.permCheck == "ON": # Do permutation (if requested by user) if self.do_control == "true": - perm_data_frame = scanone(cross_object, pheno_col = "the_pheno", addcovar = covar, n_perm = int(self.num_perm), model=self.model, method=self.method) + perm_data_frame = scanone(cross_object, pheno_col = "the_pheno", addcovar = covar, n_perm = self.num_perm, model=self.model, method=self.method) else: - perm_data_frame = scanone(cross_object, pheno_col = "the_pheno", n_perm = int(self.num_perm), model=self.model, method=self.method) + perm_data_frame = scanone(cross_object, pheno_col = "the_pheno", n_perm = self.num_perm, model=self.model, method=self.method) self.process_rqtl_perm_results(perm_data_frame) # Functions that sets the thresholds for the webinterface @@ -531,10 +526,11 @@ class MarkerRegression(object): def process_rqtl_perm_results(self, results): perm_vals = [] - for line in str(results).split("\n")[1:(int(self.num_perm)+1)]: + for line in str(results).split("\n")[1:(self.num_perm+1)]: #print("R/qtl permutation line:", line.split()) perm_vals.append(float(line.split()[1])) + self.perm_output = perm_vals self.suggestive = np.percentile(np.array(perm_vals), 67) self.significant = np.percentile(np.array(perm_vals), 95) @@ -666,8 +662,8 @@ class MarkerRegression(object): trimmed_values.append(values[i]) self.lrs_array = genotype.permutation(strains = trimmed_samples, - trait = trimmed_values, - nperm= self.num_perm) + trait = trimmed_values, + nperm= self.num_perm) self.suggestive = self.lrs_array[int(self.num_perm*0.37-1)] self.significant = self.lrs_array[int(self.num_perm*0.95-1)] @@ -699,17 +695,25 @@ class MarkerRegression(object): control_geno.append(control_geno2[_idx]) self.bootstrap_results = genotype.bootstrap(strains = trimmed_samples, - trait = trimmed_values, - control = control_geno, - nboot = self.num_bootstrap) + trait = trimmed_values, + control = control_geno, + nboot = self.num_bootstrap) else: reaper_results = genotype.regression(strains = trimmed_samples, trait = trimmed_values) if self.bootCheck: self.bootstrap_results = genotype.bootstrap(strains = trimmed_samples, - trait = trimmed_values, - nboot = self.num_bootstrap) + trait = trimmed_values, + nboot = self.num_bootstrap) + + if self.num_perm < 100: + self.suggestive = 0 + self.significant = 0 + else: + self.perm_output = genotype.permutation(strains = trimmed_samples, trait = trimmed_values, nperm=self.num_perm) + self.suggestive = self.perm_output[int(self.num_perm*0.37-1)] + self.significant = self.perm_output[int(self.num_perm*0.95-1)] self.json_data['chr'] = [] self.json_data['pos'] = [] @@ -823,7 +827,7 @@ class MarkerRegression(object): #print("self.num_perm:", self.num_perm) - for permutation in range(int(self.num_perm)): + for permutation in range(self.num_perm): pheno_vector = np.array([val == "x" and np.nan or float(val) for val in self.vals]) np.random.shuffle(pheno_vector) diff --git a/wqflask/wqflask/marker_regression/marker_regression_gn1.py b/wqflask/wqflask/marker_regression/marker_regression_gn1.py index 99f0ac99..b8307aaa 100644 --- a/wqflask/wqflask/marker_regression/marker_regression_gn1.py +++ b/wqflask/wqflask/marker_regression/marker_regression_gn1.py @@ -220,33 +220,32 @@ class MarkerRegression(object): else: self.plotScale = "physic" - #self.plotScale = fd.formdata.getvalue('scale', 'physic') - #if self.plotScale == 'physic' and not fd.genotype.Mbmap: #ZS: Not sure where "Mbmap" is stored, if at all; should be fine without this though - # self.plotScale = 'morgan' - if start_vars['num_perm'] != "": - self.nperm = int(start_vars['num_perm']) - else: - self.nperm = 0 - if 'num_bootstrap' in start_vars.keys(): - self.nboot = int(start_vars['num_bootstrap']) + if 'permCheck' in start_vars.keys(): + self.permChecked = start_vars['permCheck'] else: - self.nboot = 0 - if (start_vars['num_perm'] == "") or (start_vars['num_perm'] < 1): self.permChecked = False + if start_vars['num_perm'] > 0: + self.nperm = int(start_vars['num_perm']) + if self.permChecked: + self.perm_output = start_vars['perm_output'] + self.suggestive = start_vars['suggestive'] + self.significant = start_vars['significant'] else: - self.permChecked = True - #self.permChecked = fd.formdata.getvalue('permCheck', True) - + self.nperm = 0 + if 'bootCheck' in start_vars.keys(): self.bootChecked = start_vars['bootCheck'] else: self.bootChecked = False + if 'num_bootstrap' in start_vars.keys(): + self.nboot = int(start_vars['num_bootstrap']) + else: + self.nboot = 0 if 'bootstrap_results' in start_vars.keys(): self.bootResult = start_vars['bootstrap_results'] else: self.bootResult = [] - if 'do_control' in start_vars.keys(): self.doControl = start_vars['do_control'] else: @@ -593,9 +592,9 @@ class MarkerRegression(object): else: showLocusForm = intImg - if self.permChecked and not self.multipleInterval and 0<self.nperm: - perm_histogram = self.drawPermutationHistogram() - perm_text_file = self.permutationTextFile() + if self.permChecked and self.nperm > 0 and not self.multipleInterval and 0 < self.nperm: + self.perm_filename = self.drawPermutationHistogram() + #perm_text_file = self.permutationTextFile() ################################################################ # footnote goes here @@ -608,12 +607,12 @@ class MarkerRegression(object): if self.traitList and self.traitList[0].dataset and self.traitList[0].dataset.type == 'Geno': btminfo.append(HT.BR(), 'Mapping using genotype data as a trait will result in infinity LRS at one locus. In order to display the result properly, all LRSs higher than 100 are capped at 100.') - if self.permChecked and not self.multipleInterval and 0 < self.nperm: - TD_LR = HT.TD(HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo, HT.P(), perm_histogram, HT.P(), perm_text_file), bgColor='#eeeeee', height = 200) - #TD_LR = HT.TD(HT.Blockquote(topTable), HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo, HT.P(), perm_histogram, HT.P(), perm_text_file), bgColor='#eeeeee', height = 200) - else: - TD_LR = HT.TD(HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo), bgColor='#eeeeee', height = 200) - #TD_LR = HT.TD(HT.Blockquote(topTable), HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo, HT.P(), perm_histogram, HT.P(), perm_text_file), bgColor='#eeeeee', height = 200) + #if self.permChecked and not self.multipleInterval and 0 < self.nperm: + # TD_LR = HT.TD(HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo, HT.P(), perm_histogram, HT.P(), perm_text_file), bgColor='#eeeeee', height = 200) + # #TD_LR = HT.TD(HT.Blockquote(topTable), HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo, HT.P(), perm_histogram, HT.P(), perm_text_file), bgColor='#eeeeee', height = 200) + #else: + TD_LR = HT.TD(HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo), bgColor='#eeeeee', height = 200) + #TD_LR = HT.TD(HT.Blockquote(topTable), HT.Blockquote(gifmap, showLocusForm, HT.P(), btminfo, HT.P(), perm_histogram, HT.P(), perm_text_file), bgColor='#eeeeee', height = 200) if geneTable: @@ -683,7 +682,7 @@ class MarkerRegression(object): fpText.write(time.strftime("Date and Time (US Center): %b %d, %Y at %I.%M %p\n", time.localtime())) fpText.write("Trait ID: %s\n" % self.this_trait.name) fpText.write("Suggestive LRS = %0.2f\n" % self.suggestive) - fpText.write("Significant LRS = %0.2f\n" % self.significance) + fpText.write("Significant LRS = %0.2f\n" % self.significant) """ if self.this_trait.symbol and self.this_trait.chr and self.this_trait.mb: writeSymbol, writeChromosome, writeMb = self.this_trait.symbol, self.this_trait.chr, self.this_trait.mb @@ -720,7 +719,7 @@ class MarkerRegression(object): else: lrs_lod = marker['lod_score'] - P_value = self.calculatePValue(lrs_lod, self.LRSArray) + P_value = self.calculatePValue(lrs_lod, self.perm_output) #if _dominance: # fpText.write("%s\t%s\t%2.3f\t%s\t%2.3f\t%2.3f\t%2.3f\t%2.3f\n" %(qtlresult.locus.chr, \ @@ -1139,7 +1138,7 @@ class MarkerRegression(object): canvas.drawLine(startPosX+54,startPosY,startPosX+67,startPosY,color=self.HAPLOTYPE_RECOMBINATION, width=4) canvas.drawString('Haplotypes (Pat, Mat, Het, Unk)',startPosX+76,startPosY+5,font=labelFont,color=pid.black) - if self.permChecked: + if self.permChecked and self.nperm > 0: startPosY += stepPosY startPosX = xLeftOffset canvas.drawLine(startPosX, startPosY, startPosX + 32, startPosY, color=self.SIGNIFICANT_COLOR, width=self.SIGNIFICANT_WIDTH) @@ -1147,7 +1146,7 @@ class MarkerRegression(object): lod = 1 if self.LRS_LOD == 'LOD': lod = self.LODFACTOR - canvas.drawString('Significant %s = %2.2f' % (self.LRS_LOD, self.significance/lod),xLeftOffset+42,startPosY +5,font=labelFont,color=pid.black) + canvas.drawString('Significant %s = %2.2f' % (self.LRS_LOD, self.significant/lod),xLeftOffset+42,startPosY +5,font=labelFont,color=pid.black) canvas.drawString('Suggestive %s = %2.2f' % (self.LRS_LOD, self.suggestive/lod),xLeftOffset+42,startPosY + 5 +stepPosY,font=labelFont,color=pid.black) @@ -1884,10 +1883,10 @@ class MarkerRegression(object): #LRSMax = max(map(max, self.qtlresults)).lod_score #genotype trait will give infinite LRS LRSMax = min(LRSMax, webqtlConfig.MAXLRS) - if self.permChecked and not self.multipleInterval: - self.significance = min(self.significance, webqtlConfig.MAXLRS) + if self.permChecked and self.nperm > 0 and not self.multipleInterval: + self.significant = min(self.significant, webqtlConfig.MAXLRS) self.suggestive = min(self.suggestive, webqtlConfig.MAXLRS) - LRSMax = max(self.significance, LRSMax) + LRSMax = max(self.significant, LRSMax) else: LRSMax = self.lrsMax*lodm @@ -1922,34 +1921,46 @@ class MarkerRegression(object): #Draw the LRS/LOD Y axis label canvas.drawString(scaleStr, xLeftOffset-4-canvas.stringWidth(scaleStr, font=LRSScaleFont)-5, yLRS+3, font=LRSScaleFont, color=self.LRS_COLOR) - - #"Significant" and "Suggestive" Drawing Routine - # ======= Draw the thick lines for "Significant" and "Suggestive" ===== (crowell: I tried to make the SNPs draw over these lines, but piddle wouldn't have it...) - if self.permChecked and not self.multipleInterval: - significantY = yZero - self.significance*LRSHeightThresh/LRSMax + if self.permChecked and self.nperm > 0 and not self.multipleInterval: + significantY = yZero - self.significant*LRSHeightThresh/LRSMax suggestiveY = yZero - self.suggestive*LRSHeightThresh/LRSMax startPosX = xLeftOffset - for i, _chr in enumerate(self.genotype): - rightEdge = int(startPosX + self.ChrLengthDistList[i]*plotXScale - self.SUGGESTIVE_WIDTH/1.5) - canvas.drawLine(startPosX+self.SUGGESTIVE_WIDTH/1.5, suggestiveY, rightEdge, suggestiveY, color=self.SUGGESTIVE_COLOR, + + #"Significant" and "Suggestive" Drawing Routine + # ======= Draw the thick lines for "Significant" and "Suggestive" ===== (crowell: I tried to make the SNPs draw over these lines, but piddle wouldn't have it...) + + #ZS: I don't know if what I did here with this inner function is clever or overly complicated, but it's the only way I could think of to avoid duplicating the code inside this function + def add_suggestive_significant_lines_and_legend(start_pos_x, chr_length_dist): + rightEdge = int(start_pos_x + chr_length_dist*plotXScale - self.SUGGESTIVE_WIDTH/1.5) + canvas.drawLine(start_pos_x+self.SUGGESTIVE_WIDTH/1.5, suggestiveY, rightEdge, suggestiveY, color=self.SUGGESTIVE_COLOR, width=self.SUGGESTIVE_WIDTH*zoom, clipX=(xLeftOffset, xLeftOffset + plotWidth-2)) - canvas.drawLine(startPosX+self.SUGGESTIVE_WIDTH/1.5, significantY, rightEdge, significantY, color=self.SIGNIFICANT_COLOR, + canvas.drawLine(start_pos_x+self.SUGGESTIVE_WIDTH/1.5, significantY, rightEdge, significantY, color=self.SIGNIFICANT_COLOR, width=self.SIGNIFICANT_WIDTH*zoom, clipX=(xLeftOffset, xLeftOffset + plotWidth-2)) - sugg_coords = "%d, %d, %d, %d" % (startPosX, suggestiveY-2, rightEdge + 2*zoom, suggestiveY+2) - sig_coords = "%d, %d, %d, %d" % (startPosX, significantY-2, rightEdge + 2*zoom, significantY+2) + sugg_coords = "%d, %d, %d, %d" % (start_pos_x, suggestiveY-2, rightEdge + 2*zoom, suggestiveY+2) + sig_coords = "%d, %d, %d, %d" % (start_pos_x, significantY-2, rightEdge + 2*zoom, significantY+2) if self.LRS_LOD == 'LRS': sugg_title = "Suggestive LRS = %0.2f" % self.suggestive - sig_title = "Significant LRS = %0.2f" % self.significance + sig_title = "Significant LRS = %0.2f" % self.significant else: sugg_title = "Suggestive LOD = %0.2f" % (self.suggestive/4.61) - sig_title = "Significant LOD = %0.2f" % (self.significance/4.61) + sig_title = "Significant LOD = %0.2f" % (self.significant/4.61) Areas1 = HT.Area(shape='rect',coords=sugg_coords,title=sugg_title) Areas2 = HT.Area(shape='rect',coords=sig_coords,title=sig_title) gifmap.areas.append(Areas1) gifmap.areas.append(Areas2) - startPosX += (self.ChrLengthDistList[i]+self.GraphInterval)*plotXScale - + start_pos_x += (chr_length_dist+self.GraphInterval)*plotXScale + return start_pos_x + + for i, _chr in enumerate(self.genotype): + if self.selectedChr != -1: + if _chr.name == self.ChrList[self.selectedChr][0]: + startPosX = add_suggestive_significant_lines_and_legend(startPosX, self.ChrLengthDistList[0]) + break + else: + continue + else: + startPosX = add_suggestive_significant_lines_and_legend(startPosX, self.ChrLengthDistList[i]) if self.multipleInterval: lrsEdgeWidth = 1 @@ -2217,17 +2228,17 @@ class MarkerRegression(object): if self.multipleInterval: self.suggestive = 0 - self.significance = 0 + self.significant = 0 if self.selectedChr > -1: self.genotype.chromosome = [self.genotype[self.selectedChr]] else: #single interval mapping try: self.suggestive = float(fd.formdata.getvalue('permSuggestive')) - self.significance = float(fd.formdata.getvalue('permSignificance')) + self.significant = float(fd.formdata.getvalue('permSignificance')) except: self.suggestive = None - self.significance = None + self.significant = None _strains, _vals, _vars = self.traitList[0].exportInformative(weightedRegression) @@ -2258,21 +2269,21 @@ class MarkerRegression(object): return "The control marker you selected is not in the genofile." if weightedRegression: - self.LRSArray = self.genotype.permutation(strains = _strains, trait = _vals, + self.perm_output = self.genotype.permutation(strains = _strains, trait = _vals, variance = _vars, nperm=self.nperm) else: - self.LRSArray = self.genotype.permutation(strains = _strains, trait = _vals, + self.perm_output = self.genotype.permutation(strains = _strains, trait = _vals, nperm=self.nperm) - if self.significance and self.suggestive: + if self.significant and self.suggestive: pass else: if self.nperm < 100: self.suggestive = 0 - self.significance = 0 + self.significant = 0 else: - self.suggestive = self.LRSArray[int(self.nperm*0.37-1)] - self.significance = self.LRSArray[int(self.nperm*0.95-1)] + self.suggestive = self.perm_output[int(self.nperm*0.37-1)] + self.significant = self.perm_output[int(self.nperm*0.95-1)] #calculating bootstrap #from now on, genotype could only contain a single chromosome @@ -2443,7 +2454,7 @@ class MarkerRegression(object): controlsForm.append(controlsTable) controlsForm.append(HT.Input(name="permSuggestive", value=self.suggestive, type="hidden")) - controlsForm.append(HT.Input(name="permSignificance", value=self.significance, type="hidden")) + controlsForm.append(HT.Input(name="permSignificance", value=self.significant, type="hidden")) ## BEGIN HaplotypeAnalyst #### haplotypeAnalystCheck added below ## END HaplotypeAnalyst @@ -2504,26 +2515,27 @@ class MarkerRegression(object): # Permutation Graph ######################################### myCanvas = pid.PILCanvas(size=(400,300)) - #plotBar(myCanvas,10,10,390,290,LRSArray,XLabel='LRS',YLabel='Frequency',title=' Histogram of Permutation Test',identification=fd.identification) - Plot.plotBar(myCanvas, self.LRSArray,XLabel='LRS',YLabel='Frequency',title=' Histogram of Permutation Test') + Plot.plotBar(myCanvas, self.perm_output, XLabel='LRS', YLabel='Frequency', title=' Histogram of Permutation Test') filename= webqtlUtil.genRandStr("Reg_") myCanvas.save(webqtlConfig.IMGDIR+filename, format='gif') - img=HT.Image('/image/'+filename+'.gif',border=0,alt='Histogram of Permutation Test') - + + return filename + + # img=HT.Image('/image/'+filename+'.gif',border=0,alt='Histogram of Permutation Test') - self.suggestive = self.LRSArray[int(self.nperm*0.37-1)] - self.significant = self.LRSArray[int(self.nperm*0.95-1)] - self.highlysignificant = self.LRSArray[int(self.nperm*0.99-1)] + # self.suggestive = self.perm_output[int(self.nperm*0.37-1)] + # self.significant = self.perm_output[int(self.nperm*0.95-1)] + # self.highlysignificant = self.perm_output[int(self.nperm*0.99-1)] - permutationHeading = HT.Paragraph('Histogram of Permutation Test') - permutationHeading.__setattr__("class","title") + # permutationHeading = HT.Paragraph('Histogram of Permutation Test') + # permutationHeading.__setattr__("class","title") - permutation = HT.TableLite() - permutation.append(HT.TR(HT.TD(img)), - HT.TR(HT.TD('')), - HT.TR(HT.TD('Total of %d permutations'%self.nperm))) + # permutation = HT.TableLite() + # permutation.append(HT.TR(HT.TD(img)), + # HT.TR(HT.TD('')), + # HT.TR(HT.TD('Total of %d permutations'%self.nperm))) - return permutation + # return permutation def permutationTextFile(self): filename= webqtlUtil.genRandStr("Reg_") @@ -2531,14 +2543,14 @@ class MarkerRegression(object): fpText.write('Suggestive LRS (p = 0.63) = %3.2f\n'%self.suggestive) fpText.write('Significant LRS (p = 0.05) = %3.2f\n'%self.significant) fpText.write('Highly Significant LRS (p = 0.01) = %3.2f\n\n'%self.highlysignificant) - fpText.write('%s Permutations\n\n' % str(len(self.LRSArray))) + fpText.write('%s Permutations\n\n' % str(len(self.perm_output))) LRSInfo =HT.Paragraph(' Suggestive LRS = %3.2f\n'%self.suggestive, HT.BR(), ' Significant LRS =%3.2f\n'%self.significant, HT.BR(), ' Highly Significant LRS =%3.2f\n' % self.highlysignificant) - for lrs_value in self.LRSArray: + for lrs_value in self.perm_output: fpText.write(str(lrs_value) + "\n") textUrl = HT.Href(text = 'Download Permutation Results', url= '/tmp/'+filename+'.txt', target = "_blank", Class='fs12 fwn') diff --git a/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js b/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js index 68744aa9..519d1304 100755 --- a/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js +++ b/wqflask/wqflask/static/new/javascript/show_trait_mapping_tools.js @@ -228,9 +228,10 @@ //$("#progress_bar_container").modal(); url = "/marker_regression"; $('input[name=method]').val("reaper"); - $('input[name=manhattan_plot]').val($('input[name=manhattan_plot_reaper]:checked').val()); + $('input[name=num_perm]').val($('input[name=num_perm_reaper]').val()); $('input[name=control_marker]').val($('input[name=control_reaper]').val()); $('input[name=do_control]').val($('input[name=do_control_reaper]:checked').val()); + $('input[name=manhattan_plot]').val($('input[name=manhattan_plot_reaper]:checked').val()); $('input[name=mapping_display_all]').val($('input[name=display_all_reaper]')); $('input[name=suggestive]').val($('input[name=suggestive_reaper]')); form_data = $('#trait_data_form').serialize(); diff --git a/wqflask/wqflask/templates/marker_regression_gn1.html b/wqflask/wqflask/templates/marker_regression_gn1.html index 88fb0e89..8ec8a104 100644 --- a/wqflask/wqflask/templates/marker_regression_gn1.html +++ b/wqflask/wqflask/templates/marker_regression_gn1.html @@ -21,7 +21,7 @@ <input type="hidden" name="maf"> <input type="hidden" name="selected_chr" value="{{ selectedChr }}"> <input type="hidden" name="manhattan_plot"> - <input type="hidden" name="num_perm"> + <input type="hidden" name="num_perm" value="{{ nperm }}"> <input type="hidden" name="num_bootstrap" value="{{ nboot }}"> <input type="hidden" name="do_control" value="{{ doControl }}"> <input type="hidden" name="control_marker" value="{{ controlLocus }}"> @@ -69,8 +69,13 @@ </table> </div> <div class="col-xs-4" style="padding: 0px;"> - {% if mapping_method == "reaper" %} + {% if (mapping_method == "reaper" or mapping_method == "rqtl_geno") and nperm > 0 %} + <input type="checkbox" name="permCheck" class="checkbox" style="display: inline; margin-top: 0px;" {% if permChecked|upper == "ON" %}value="ON" checked{% endif %}> <span style="font-size: 12px;">Permutation Test <br> + {% endif %} + {% if mapping_method == "reaper" and nboot > 0 %} <input type="checkbox" name="bootCheck" class="checkbox" style="display: inline; margin-top: 0px;" {% if bootChecked|upper == "ON" %}value="ON" checked{% endif %}> <span style="font-size: 12px;">Bootstrap Test <br> + {% endif %} + {% if mapping_method == "reaper" %} <input type="checkbox" name="additiveCheck" class="checkbox" style="display: inline; margin-top: 0px;" {% if additiveChecked|upper == "ON" %}value="ON" checked{% endif %}> <span style="font-size: 12px;">Allele Effects<br> {% endif %} <input type="checkbox" name="showGenes" class="checkbox" style="display: inline; margin-top: 0px;" {% if geneChecked|upper == "ON" %}value="ON" checked{% endif %}> <span style="font-size: 12px;">Gene Track </span> <span style="color:red;">*</span><br> @@ -98,6 +103,14 @@ <div class="qtlcharts"> {{ gifmap|safe }} <img src="/static/output/{{ filename }}.jpeg" usemap="#WebQTLImageMap"> + {% if additiveChecked|upper == "ON" %} + <br> + <span style="white-space: nowrap;">A positive additive coefficient (green line) indicates that {{ dataset.group.parlist[1] }} alleles increase trait values. In contrast, a negative additive coefficient (orange line) indicates that {{ dataset.group.parlist[0] }} alleles increase trait values.</span> + {% endif %} + {% if permChecked|upper == "ON" %} + <br><br> + <img src="/static/output/{{ perm_filename }}.gif"> + {% endif %} </div> </div> {% if mapping_method != "gemma" %} diff --git a/wqflask/wqflask/views.py b/wqflask/wqflask/views.py index 9d0fd131..acfc5fc6 100644 --- a/wqflask/wqflask/views.py +++ b/wqflask/wqflask/views.py @@ -347,8 +347,10 @@ def marker_regression_page(): 'mapping_scale', 'score_type', 'suggestive', + 'significant', 'num_perm', 'permCheck', + 'perm_output', 'num_bootstrap', 'bootCheck', 'bootstrap_results', |