From b1dce9603deddeaa0db2c4867e5a3788b52d6ef0 Mon Sep 17 00:00:00 2001 From: Javier Herrero Date: Mon, 25 Jan 2016 17:31:24 +0000 Subject: [PATCH] Fix mixup between strict and marginal thresholds of significance. Also fixes pink color on dChart --- eForge/ePlot.pm | 33 ++++++++++++++++++--------------- eforge.pl | 19 +++++++++++-------- webserver/cgi-bin/index.pl | 37 ++++++++++++++++++------------------- 3 files changed, 47 insertions(+), 42 deletions(-) diff --git a/eForge/ePlot.pm b/eForge/ePlot.pm index 07f0b04..602ff65 100644 --- a/eForge/ePlot.pm +++ b/eForge/ePlot.pm @@ -71,7 +71,7 @@ This is the original code using standard R plot to generate a static pdf. sub Chart{ print "Making static chart.\n"; - my ($filename, $lab, $resultsdir, $tissues, $cells, $label, $t1, $t2, $data) = @_; + my ($filename, $lab, $resultsdir, $tissues, $cells, $label, $t_marginal, $t_strict, $data) = @_; my $Rdir = $resultsdir; my $chart = "$lab.chart.pdf"; my $rfile = "$Rdir/$lab.chart.R"; @@ -80,17 +80,17 @@ sub Chart{ open my $rfh, ">", $rfile; -#results\$Class<-cut(results\$Pvalue, breaks =c(min(results\$Pvalue), $t1, $t2, max(results\$Pvalue)), labels=FALSE, include.lowest=TRUE) # 99 and 95% CIs 1, 2, 3 -$t1 = sprintf("%.2f", $t1); -$t2 = sprintf("%.2f", $t2); +#results\$Class<-cut(results\$Pvalue, breaks =c(min(results\$Pvalue), $t_marginal, $t_strict, max(results\$Pvalue)), labels=FALSE, include.lowest=TRUE) # 99 and 95% CIs 1, 2, 3 +$t_marginal = sprintf("%.2f", $t_marginal); +$t_strict = sprintf("%.2f", $t_strict); print $rfh "setwd('$Rdir') results<-read.table('$filename', header=TRUE,sep='\\t') -# Class splits the data into non-significant, marginally significant and significant according to $t1 and $t2 (in -log10 scale) -results\$Class <- cut(results\$Pvalue, breaks =c(0, $t2, $t1, 1)/length(unique(results[,'Tissue'])), labels=FALSE, include.lowest=TRUE) +# Class splits the data into non-significant, marginally significant and significant according to $t_marginal and $t_strict (in -log10 scale) +results\$Class <- cut(results\$Pvalue, breaks =c(0, $t_strict, $t_marginal, 1)/length(unique(results[,'Tissue'])), labels=FALSE, include.lowest=TRUE) # Class splits the data into non-significant, marginally significant and significant according to q-value (B-Y FDR adjusted) -results\$Class2 <- cut(results\$Qvalue, breaks =c(0, $t2, $t1, 1), labels=FALSE, include.lowest=TRUE) +results\$Class2 <- cut(results\$Qvalue, breaks =c(0, $t_strict, $t_marginal, 1), labels=FALSE, include.lowest=TRUE) # Re-order the entries according to tissue first and then cell type/line tissue.cell.order <- unique(results[, c('Tissue', 'Cell')]) @@ -131,11 +131,11 @@ mtext(2, text='-log10 binomial p-value', line=2, cex=1.4) # Add legend (internal color first) palette(c('white', '$msig', '$sig')) -legend('topleft', pch=19, legend=c('q < $t2', 'q < $t1', 'non-sig'), col = 3:1, cex=0.8, inset=c(0.001, 0.005), box.col='white', title='FDR q-value', text.col='white', bg='white') +legend('topleft', pch=19, legend=c('q < $t_strict', 'q < $t_marginal', 'non-sig'), col = 3:1, cex=0.8, inset=c(0.001, 0.005), box.col='white', title='FDR q-value', text.col='white', bg='white') # Add contour to the points in the legend palette(c('$ns', '$msig', 'black')) -legend('topleft', pch=1, legend=c('q < $t2', 'q < $t1', 'non-sig'), col = 3:1, cex=0.8, inset=c(0.001, 0.005), box.col='darkgrey', title='FDR q-value') +legend('topleft', pch=1, legend=c('q < $t_strict', 'q < $t_marginal', 'non-sig'), col = 3:1, cex=0.8, inset=c(0.001, 0.005), box.col='darkgrey', title='FDR q-value') palette('default') dev.off() @@ -153,7 +153,7 @@ Make dimple interactive chart. =cut sub dChart{ - my ($filename, $lab, $resultsdir, $data, $label, $t1, $t2, $web) = @_; + my ($filename, $lab, $resultsdir, $data, $label, $t_marginal, $t_strict, $web) = @_; print "Making dChart.\n"; my $chart = "$lab.dchart.html"; @@ -163,19 +163,22 @@ sub dChart{ print $rcfh "setwd(\"$Rdir\") results<-read.table(\"$filename\", header = TRUE, sep=\"\\t\") -# Class splits the data into non-significant, marginally significant and significant according to $t1 and $t2 (in -log10 scale) -results\$Class <- cut(results\$Pvalue, breaks =c(0, $t2, $t1, 1)/length(unique(results[,'Tissue'])), labels=FALSE, include.lowest=TRUE) +# Class splits the data into non-significant, marginally significant and significant according to $t_marginal and $t_strict (in -log10 scale) +results\$Class <- cut(results\$Pvalue, breaks =c(0, $t_strict, $t_marginal, 1)/length(unique(results[,'Tissue'])), labels=FALSE, include.lowest=TRUE) # Class splits the data into non-significant, marginally significant and significant according to q-value (B-Y FDR adjusted) -results\$Class2 <- cut(results\$Qvalue, breaks =c(0, $t2, $t1, 1), labels=FALSE, include.lowest=TRUE) +results\$Class2 <- cut(results\$Qvalue, breaks =c(0, $t_strict, $t_marginal, 1), labels=FALSE, include.lowest=TRUE) color.axis.palette = c(); if (length(which(results\$Class2 == 1)) > 0 ) { color.axis.palette = c('red'); } if (length(which(results\$Class2 == 2)) > 0 ) { - color.axis.palette = c(color.axis.palette, 'pink'); + color.axis.palette = c(color.axis.palette, '#FF82ab'); +} +color.axis.palette = c(color.axis.palette, 'lightblue'); +if (length(color.axis.palette) < 2) { + color.axis.palette = c(color.axis.palette, 'lightblue'); # Add it twice to force the color if only non-significant values } -color.axis.palette = c(color.axis.palette, 'lightblue', 'lightblue'); # Add it twice to force the color if only non-significant values results\$log10pvalue <- -log10(results\$Pvalue) diff --git a/eforge.pl b/eforge.pl index 1a53a9c..de3c034 100644 --- a/eforge.pl +++ b/eforge.pl @@ -338,15 +338,18 @@ =head1 ACKNOWLEDGEMENTS # Define the thresholds to use. -my ($t1, $t2); +my ($t_marginal, $t_strict); if (defined $thresh) { - ($t1, $t2) = split(",", $thresh); - unless (looks_like_number($t1) && looks_like_number($t2)){ - die "You must specify numerical p value thresholds in a comma separated list"; + ($t_marginal, $t_strict) = split(",", $thresh); + unless (looks_like_number($t_marginal) && looks_like_number($t_strict)){ + die "You must specify numerical p value thresholds in a comma separated list\n"; + } + unless ((1 >= $t_marginal) && ($t_marginal >= $t_strict) && ($t_strict >= 0)) { + die "The p-value thresholds should be 1 >= T.marginal >= T.strict >= 0\n"; } } else { - $t1 = 0.05; # set binomial p values, bonferroni is applied later based on number of samples (cells) - $t2 = 0.01; + $t_marginal = 0.05; # set binomial p values, bonferroni is applied later based on number of samples (cells) + $t_strict = 0.01; } # mvps need to come either from a file or a list @@ -576,8 +579,8 @@ =head1 ACKNOWLEDGEMENTS warn "[".scalar(localtime())."] Generating plots...\n"; unless (defined $noplot){ #Plotting and table routines - Chart($filename, $lab, $out_dir, $tissues, $cells, $label, $t1, $t2, $data); # basic pdf plot - dChart($filename, $lab, $out_dir, $data, $label, $t1, $t2, $web); # rCharts Dimple chart + Chart($filename, $lab, $out_dir, $tissues, $cells, $label, $t_marginal, $t_strict, $data); # basic pdf plot + dChart($filename, $lab, $out_dir, $data, $label, $t_marginal, $t_strict, $web); # rCharts Dimple chart table($filename, $lab, $out_dir, $web); # Datatables chart } diff --git a/webserver/cgi-bin/index.pl b/webserver/cgi-bin/index.pl index e078801..2f32255 100755 --- a/webserver/cgi-bin/index.pl +++ b/webserver/cgi-bin/index.pl @@ -277,10 +277,10 @@ sub print_form { textfield('reps', '1000', 10)]), td(["Significance threshold:", ""]), - td(["      Min:", - textfield('thresh1', '0.01', 10)]), - td(["      Max:", - textfield('thresh2', '0.05', 10)]), + td(["      Strict:", + textfield('thresh_strict', '0.01', 10)]), + td(["      Marginal:", + textfield('thresh_marginal', '0.05', 10)]), th({-colspan=>2}, ["
"]), @@ -393,27 +393,26 @@ sub validate_form { } push(@$validated_args, "--reps", $reps); - my $thresh1 = param("thresh1"); - if (!$thresh1 or $thresh1 !~ /^\d+(\.\d*)?$/) { - push(@error_messages, "Min. significance threshold must be a positive number."); - } elsif ($thresh1 >= 1) { - push(@error_messages, "Min. significance threshold must be less than 1."); + my $thresh_strict = param("thresh_strict"); + if (!$thresh_strict or $thresh_strict !~ /^\d+(\.\d*)?$/) { + push(@error_messages, "Strict significance threshold must be a positive number."); + } elsif ($thresh_strict >= 1) { + push(@error_messages, "Strict significance threshold must be less than 1."); } - my $thresh2 = param("thresh2"); - if (!$thresh2 or $thresh2 !~ /^\d+(\.\d*)?$/) { - push(@error_messages, "Max. significance threshold must be a positive number."); - } elsif ($thresh2 >= 1) { - push(@error_messages, "Max. significance threshold must be less than 1."); + my $thresh_marginal = param("thresh_marginal"); + if (!$thresh_marginal or $thresh_marginal !~ /^\d+(\.\d*)?$/) { + push(@error_messages, "Marginal significance threshold must be a positive number."); + } elsif ($thresh_marginal >= 1) { + push(@error_messages, "Marginal significance threshold must be less than 1."); } - if ($thresh1 > $thresh2) { + if ($thresh_strict > $thresh_marginal) { print $q->header; - push(@error_messages, "Min. significance threshold must be less than max. significance". - " threshold."); + push(@error_messages, "Strict significance threshold must be less than marginal one."); return 0; } - push(@$validated_args, "--thresh", "$thresh1,$thresh2"); + push(@$validated_args, "--thresh", "$thresh_marginal,$thresh_strict"); if (@error_messages) { print $q->header; @@ -915,7 +914,7 @@ sub print_help_page {
", "Significance threshold

- Alter the default binomial p value thresholds. (0 < Min < Max < 1) + Alter the default binomial p value thresholds. (0 < Strict < Marginal < 1)
", );