Final installation of methylseq analysis.

diff --git a/htswanalysis/reference_data/qPCR_Tests/RNA b/htswanalysis/reference_data/qPCR_Tests/RNA
new file mode 100755 (executable)
index 0000000..de12e2c
--- /dev/null
+++ b/htswanalysis/reference_data/qPCR_Tests/RNA
@@ -0,0 +1,6 @@
+1      chr19   18251563        18253432        1
+2      chr12   6517076 6517400 1
+3      chr7    5534457 5534627 1
+4      chrX    153412800       153413450       1
+5      chr22   39902800        39904900        1
+6      chr6    44327720        44327950        1

diff --git a/htswanalysis/scripts/ConfigureTasks.pm b/htswanalysis/scripts/ConfigureTasks.pm
index 8fe851ea12611df2f123b8d9cf8d753aea6e927e..d4348c8b25b10273ac9a0141740909de24d384ea 100755 (executable)
--- a/htswanalysis/scripts/ConfigureTasks.pm
+++ b/htswanalysis/scripts/ConfigureTasks.pm
@@ -52,8 +52,7 @@ my $SUMMARIZE_PROJECT="$root_dir/scripts/SummarizeProject2.pm";
 ###
 ### Important: get the xml file form the server.
 ###
-#my $projects_file = getProjectsXML($parm);
-my $projects_file = "/Users/Data/Projects.xml";
+my $projects_file = getProjectsXML($parm);
 
 my $xmldoc = XML::Simple->new();
 
@@ -101,7 +100,7 @@ close(MAKE);
 
 print STDERR "Project makefile compete. Starting to build tasks\n";
 
-#`cd $data_dir/Projects && make -j 4 -f TaskMakefile > make_tasks.log 2> make_tasks.err`;
+`cd $data_dir/Projects && make -j 4 -f TaskMakefile > make_tasks.log 2> make_tasks.err`;
 
 sub writeProject {
   my $project = shift;
@@ -187,7 +186,7 @@ sub WriteMethylseqTasks {
          $seqcheck .= "\tif [ ! -e $root_dir/reference_data/".$genome."_methylseq_regions ]; then $root_dir/scripts/analys_track_main.py updsts $task \"No methylseq regions available for $genome.\"; fi;\n"; 
 
       my $cmd  = "$outfile: $data_dir/Libraries/$msp1.txt $data_dir/Libraries/$hpa2.txt\n";
-         $cmd .= "\t$METHYLSEQDIR/Methylseq_Analysis $data_dir/Libraries/$msp1.txt $data_dir/Libraries/$hpa2.txt $root_dir/reference_data/".$genome."_methylseq_regions > \$@\n";
+         $cmd .= "\t$METHYLSEQDIR/Methylseq_Analysis $data_dir/Libraries/$msp1.txt $data_dir/Libraries/$hpa2.txt $root_dir/reference_data/".$genome."_methylseq_regions $name > \$@\n";
       writeTask($mseq, $root, $outfile, $cmd, $seqcheck);
 
       $tasks .= $task." ";
@@ -381,8 +380,8 @@ sub WritePeakCallingTasks {
        $cmd .= "\n.PRECIOUS: ".$name.'_peaks.xls '.$name.'_peaks.bed '.$name."_peaks.fasta\n\n";
 
         $cmd .= "\n".$name."_peaks.xls: $data_dir/Libraries/$signal.txt $data_dir/Libraries/$bg.txt\n";
-        $cmd .= "\tcat $data_dir/Libraries/$signal.txt | $root_dir/scripts/align_to_bed.pm | grep -v 'contam' | grep -v 'humRibosomal' > $signal.bed\n";
-        $cmd .= "\tcat $data_dir/Libraries/$bg.txt | $root_dir/scripts/align_to_bed.pm  | grep -v 'contam' | grep -v 'humRibosomal' > $bg.bed\n";
+        $cmd .= "\tcat $data_dir/Libraries/$signal.txt | $root_dir/scripts/align_to_bed.pm | grep -v 'splice' | grep -v 'contam' | grep -v 'humRibosomal' > $signal.bed\n";
+        $cmd .= "\tcat $data_dir/Libraries/$bg.txt | $root_dir/scripts/align_to_bed.pm  | grep -v 'splice' | grep -v 'contam' | grep -v 'humRibosomal' > $bg.bed\n";
         $cmd .= "\t$MACS -t $signal.bed -c ./$bg.bed --name=$name --pvalue=1e-10 --mfold=20 > $name.log 2> $name.err\n";
        $cmd .= "\trm -f $signal.bed $bg.bed\n";
        $cmd .= "\trm -f Background.ELAND.pos $name.ELAND.pos $name.R0.*.bar Background.ELAND.sorted $name.ELAND.sorted\n";

diff --git a/htswanalysis/scripts/SummarizeProject2.pm b/htswanalysis/scripts/SummarizeProject2.pm
index 25ffa0273497b5a05919370da182e92c564c87a3..1ba8334d3d446f3bde32199ec14d5d91c70d65b4 100755 (executable)
--- a/htswanalysis/scripts/SummarizeProject2.pm
+++ b/htswanalysis/scripts/SummarizeProject2.pm
@@ -7,7 +7,7 @@ my $library_info = shift;
 my %libs; 
 my %quest;
 my $expxml = XML::Simple->new();
-my $xml = $expxml->XMLin("Project.xml", ForceArray => [ qw(ComparePeakCalls CompareLibraries PeakCalling ProfileReads qPCR MotifFinding)]);
+my $xml = $expxml->XMLin("Project.xml", ForceArray => [ qw(ComparePeakCalls CompareLibraries PeakCalling ProfileReads qPCR MotifFinding Methylseq)]);
 my @peak_calling = ();
 
 my $project_name = $xml->{Name};
@@ -76,8 +76,10 @@ $qPCR_Summary .= "</TABLE>\n";
 
 my $methylseq_summary = "";
 if(exists($xml->{Methylseq})) {
-  $methylseq_sumary "<TABLE BORDER=1>\n";
+  $methylseq_summary .= "<H2>Methylseq Analysis</H2>\n";
+  $methylseq_summary  .= "<TABLE BORDER=1>\n";
   $methylseq_summary .= "<TR>\n";
+  $methylseq_summary .= "<TD><EM>Task</EM></TD>\n";
   $methylseq_summary .= "<TD><EM>msp1</EM></TD>\n";
   $methylseq_summary .= "<TD><EM>hpa2</EM></TD>\n";
   $methylseq_summary .= "<TD><EM>Assayed</EM></TD>\n";
@@ -93,7 +95,7 @@ if(exists($xml->{Methylseq})) {
     my $msp1   = $xml->{Methylseq}->[$i]->{Msp1}->{Library};
     my $hpa2   = $xml->{Methylseq}->[$i]->{Hpa2}->{Library};
 
-    $methylseq_summary .= `$root_dir/scripts/SummarizeMethylseq.pm $name $task $msp1 $hpa2 $genome`;
+    $methylseq_summary .= `$root_dir/scripts/Summarize_Methylseq.pm $name $task ../../Tasks $msp1 $hpa2 $genome`;
 
   }
   $methylseq_summary .= "</TABLE>\n";
@@ -112,7 +114,7 @@ if(exists($xml->{ProfileReads})) {
     $libs{$lib} = 0;
 
     if($genome eq "scer") { $genome = "sacCer1"; }
-    $profile_summary .= "<TR><TD>$name ($lib)</TD><TD><A HREF=../../Tasks/$task/$lib.wig.gz>$lib.wig.gz</A><BR><A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$genome&position=chr1:1000-2000&hgt.customText=http://171.65.76.194/Tasks/$task/$lib.wig.gz TARGET=\"_blank\">View in UCSC Genome Browser</A></TD>\n";
+    $profile_summary .= "<TR><TD>$name ($lib)</TD><TD><A HREF=../../Tasks/$task/$lib.wig.gz>$lib.wig.gz</A><BR><A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$genome&position=chr1:1000-2000&hgt.customText=http://illumina-mac.stanford.edu/Tasks/$task/$lib.wig.gz TARGET=\"_blank\">View in UCSC Genome Browser</A></TD>\n";
     $profile_summary .= "<TD><A HREF=../../Tasks/$task/$lib.profile.gif><IMG WIDTH=300 SRC=../../Tasks/$task/$lib.profile.gif></A></TD></TR>\n";
   }
   $profile_summary .= "</TABLE>\n";
@@ -223,7 +225,7 @@ for(@peak_calling) {
   my $color = "#EEEEFF";
   print "<TR BGCOLOR=$color><TD><B>$hash{Name}</B></TD>";
   print "<TD>$hash{Caller}</TD><TD>$hash{Summary}</TD>\n";
-  print "<TD><A HREF=$hash{outfile}>BED file</A><BR><A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$hash{Genome}&hgt.customText=http://171.65.76.194/Tasks/$hash{Task}/$hash{outfile}>View in Genome Browser</A></TD><TD><A HREF=$hash{fasta}>FASTA</A></TD>\n";
+  print "<TD><A HREF=$hash{outfile}>BED file</A><BR><A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$hash{Genome}&hgt.customText=http://illumina-mac.stanford.edu/Tasks/$hash{Task}/$hash{outfile}>View in Genome Browser</A></TD><TD><A HREF=$hash{fasta}>FASTA</A></TD>\n";
   if(exists($hash{primer_design})) { print "<TD><A HREF=$hash{primer_design}>Validation Primers</A></TD>\n"; }
   print "</TR>\n";
 }

diff --git a/htswanalysis/scripts/Summarize_Methylseq.pm b/htswanalysis/scripts/Summarize_Methylseq.pm
index 044f6c911cc1797c9d4bda9056bbf3d171d62a0a..d5f84efaae9e7fd96e4d7306efdb30aa487f8967 100755 (executable)
--- a/htswanalysis/scripts/Summarize_Methylseq.pm
+++ b/htswanalysis/scripts/Summarize_Methylseq.pm
@@ -7,21 +7,33 @@ my $name = shift;
 if(!defined($name)) { print STDERR "Error: No name provided to summarize methylseq\n"; exit(1); }
 my $task = shift;
 if(!defined($task)) { print STDERR "Error: No task provided to summarize methylseq\n"; exit(1); }
+my $tasks_root = shift;
+if(!defined($tasks_root)) { print STDERR "Error: No root task directory provided to summarize methylseq\n"; exit(1); }
 my $msp1_lib = shift;
 if(!defined($msp1_lib)) { print STDERR "Error: No msp1 library provided to summarize methylseq\n"; exit(1); }
 my $hpa2_lib = shift;
 if(!defined($hpa2_lib)) { print STDERR "Error: No hpa2 library provided to summarize methylseq\n"; exit(1); }
 my $genome = shift;
 
-my $assayed_bed = "Assayed.bed";
-my $methylated_bed = "Methylated.bed";
+my $methylseq_bed = "MethylseqCalls.bed";
 my $methylseq_calls = "MethylseqCalls.tab";
 my $summary = "methylseq_summary";
-if(!defined($summary)) { return; }
 
-open(FILE,$summary);
+if(! -e $tasks_root."/".$task."/".$summary) {
+  printf "<TR><TD BGCOLOR=#FFCCCC>Processing...</TD></TR>\n";
+  exit;
+}
 
-my $tag_count = <FILE>; chomp $tag_count;
+open(FILE,$tasks_root."/".$task."/".$summary);
+
+my $tag_count = <FILE>;
+
+if(!defined($tag_count)) {
+  printf "<TR><TD BGCOLOR=#FFCCCC>Processing...</TD></TR>\n";
+  exit;
+}
+
+chomp $tag_count;
 $tag_count =~ /^(\d+) MSP1 tags, (\d+) HPA2 tags$/;
 my($msp_count,$hpa_count) = ($1,$2);
 
@@ -37,6 +49,7 @@ my $error_count = <FILE>; chomp $error_count;
    $error_count =~ /^(\d+) of \d+ regions with error/;
    my $n_error = $1;
 
+
 my $pcnt_assayed = $n_assayed / $n_regions;
 my $pcnt_methylated = $n_methylated / $n_assayed;
 my $pcnt_error = $n_error / $n_regions;
@@ -57,10 +70,8 @@ elsif($pcnt_error > 0.01) { $error_color = "4444FF"; }
 else {$error_color= "66FF66"; }
 
 print "<TR>";
-printf "<TD>$msp1_lib</TD><TD>$hpa2_lib</TD><TD>$n_assayed (%.4f)</TD><TD>$n_methylated (%0.4f)</TD><TD>$n_error (%0.4f)</TD>",$pcnt_assayed,$pcnt_methylated,$pcnt_error;
-print "<TD><A HREF=../../Tasks/$task/$methylseq_calls>Calls</A><BR>\n";
-print "<A HREF=../../Tasks/$task/Assayed.bed>Assayed Regions</A>";
-print "<A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$genome&position=chr1:1000-2000&hgt.customText=http://illumina-mac.stanford.edu/Tasks/$task/Assayed.bed TARGET=\"_blank\">(Genome Browser)</A></BR>\n";
-print "<A HREF=../../Tasks/$task/Methylated.bed>Methylated Regions</A> ";
-print "<A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$genome&position=chr1:1000-2000&hgt.customText=http://illumina-mac.stanford.edu/Tasks/$task/Methylated.bed TARGET=\"_blank\">(Genome Browser)</A></TD>\n";
+printf "<TD>$name</TD><TD>$msp1_lib</TD><TD>$hpa2_lib</TD><TD BGCOLOR=#$assayed_color>$n_assayed (%.4f)</TD><TD BGCOLOR=#$methylated_color>$n_methylated (%0.4f)</TD><TD BGCOLOR=$error_color>$n_error (%0.4f)</TD>",$pcnt_assayed,$pcnt_methylated,$pcnt_error;
+print "<TD><A HREF=../../Tasks/$task/$methylseq_calls>Results</A><BR>\n";
+print "<A HREF=http://genome.ucsc.edu/cgi-bin/hgTracks?db=$genome&position=chr1:1000-2000&hgt.customText=http://illumina-mac.stanford.edu/Tasks/$task/$methylseq_bed TARGET=\"_blank\">(Genome Browser)</A></BR>\n";
+print "</TD>";
 print "</TR>";