Browse Source

added convert_tsv2gff3 function

Robb, Sofia 8 years ago
parent
commit
514453cab1
1 changed files with 106 additions and 1 deletions
  1. 106 1
      api/blast_ui.api.inc

+ 106 - 1
api/blast_ui.api.inc

@@ -141,6 +141,7 @@ function get_BLAST_job($job_id) {
   $job->files->result->xml = $blastjob->result_filestub . '.xml';
   $job->files->result->tsv = $blastjob->result_filestub . '.tsv';
   $job->files->result->html = $blastjob->result_filestub . '.html';
+  $job->files->result->gff = $blastjob->result_filestub . '.gff';
   
   return $job;
 }
@@ -167,6 +168,7 @@ function run_BLAST_tripal_job($program, $query, $database, $output_filestub, $op
   $output_file_xml = $output_filestub . '.xml';
   $output_file_tsv = $output_filestub . '.tsv';
   $output_file_html = $output_filestub . '.html';
+  $output_file_gff = $output_filestub . '.gff';
 
   print "\nExecuting $program\n\n";
   print "Query: $query\n";
@@ -243,6 +245,19 @@ function run_BLAST_tripal_job($program, $query, $database, $output_filestub, $op
     );
   }
 
+  print "\tGFF\n";
+  convert_tsv2gff3($output_file_tsv,$output_file_gff);
+
+  if (!file_exists($output_file_gff)) {
+    tripal_report_error(
+      'blast_ui',
+      TRIPAL_WARNING,
+      "Unable to convert BLAST Tabular Output to GFF Output (%archive => %file).",
+      array('%archive' => $output_file, '%file' => $output_file_gff),
+      array('print' => TRUE)
+    );
+  }
+
   print "\tHTML (includes alignments)\n";
   system("$blast_formatter_command -archive $output_file -outfmt 0 -out $output_file_html -html");
   if (!file_exists($output_file_tsv)) {
@@ -668,4 +683,94 @@ function generate_blast_hit_image($acc = '', $scores, $hits, $tsize, $qsize, $na
   ob_end_clean();
    
   return $b64_img;
-}
+}
+
+/** convert tsv blast output to gff output file
+ *
+ * Created by Sofia Robb
+ * 09/15/2016
+ *
+ * the subject (hit) will be the source feature
+ * the query will be the target
+ *
+ * @param $blast_tsv, $blast_gff
+ *  the name of the blast tsv output file
+ *  the name of the blast gff output file
+ *
+ * @return
+ *   returns nothing but creates the gff file in the same dir as the tsv with gff extension
+ */
+
+function convert_tsv2gff3($blast_tsv,$blast_gff){
+  $gff = fopen($blast_gff,"w");
+  fwrite($gff,"##gff-version 3\n");
+
+  $results = array();
+
+  $tsv = fopen($blast_tsv, "r") or die("Unable to open tsv file!");
+  while(!feof($tsv)) {
+    $line = fgets($tsv);
+    $line = rtrim($line);
+    if (preg_match('/^#/',$line) or preg_match('/^\s*$/',$line)){
+      continue;
+    }
+    //$line has these parts: $queryId, $subjectId, $percIdentity, $alnLength, $mismatchCount,$gapOpenCount, $queryStart, $queryEnd, $subjectStart,$subjectEnd, $eVal, $bitScore
+   $parts = preg_split('/\t/', $line);
+   $hitname = $parts[1];
+   $hspInfo =  "$parts[0],$parts[8],$parts[9],$parts[6],$parts[7]";
+   $eval = $parts[10];
+   $results[$hitname][$hspInfo]=$eval;
+   $IDs = array();
+   $count;
+   $last_q;
+   //need to go thru each line of tsv to find the first and last hsp of a hit.
+   //need to get the smallest and largest coordinate for the parent feature
+   foreach ($results as $s => $hspInfoArray) {
+     $hsp = 0;
+     foreach ($hspInfoArray as $hspInfoStr => $e){
+       list($q,$ss,$se,$qs,$qe) = preg_split('/,/',$hspInfoStr);
+       if (!$last_q or  $q != $last_q){
+          $count++;
+          $hsp=0;
+       }
+       $q_strand = '+';
+       if ($qs > $qe){
+         list($qs,$qe) = array($qe,$qs);
+         $q_strand = '-';
+       }
+       $IDs["$s,$q"]['strand']='+';
+       list($start,$end) = array($ss,$se);
+       if($ss > $se){
+         list($start,$end) = array($se,$ss);
+         $IDs["$s,$q"]['strand']='-';
+       }
+       if (!array_key_exists('SS',$IDs["$s,$q"]) or $ss < $IDs["$s,$q"]['SS']){
+          $IDs["$s,$q"]['SS'] = $ss;
+        }
+        if (!array_key_exists('SE',$IDs["$s,$q"]) or $se > $IDs["$s,$q"]['SE']){
+          $IDs["$s,$q"]['SE'] = $se;
+        }
+        if (!array_key_exists('E',$IDs["$s,$q"]) or $e < $IDs["$s,$q"]['E']){
+          $IDs["$s,$q"]['E'] = $e;
+        }
+        $hsp++;
+        $IDs["$s,$q"]['HSPs'][] = join("\t", array($s, "BLASTRESULT" , "match_part" , $start , $end , $e , $IDs["$s,$q"]['strand'] , '.' , "ID=$s.$count.$hsp;Parent=$s.$count;Target=$q $qs $qe $q_strand"));
+        $last_q = $q;
+
+     }
+   }
+
+  }//end of tsv file while
+  // now can print a parent gff line and all the children
+  // the evalues seem to be sorted properly without actually sorting them, need to make sure this is always true
+  foreach ($IDs as $sq => $value ){
+    list ($s,$q) = preg_split('/,/' , $sq);
+    $evalue =  $IDs["$s,$q"]['E'];
+    $parent =  join ("\t", array($s, "BLASTRESULT" , "match" , $IDs["$s,$q"]['SS'] , $IDs["$s,$q"]['SE'] , $IDs["$s,$q"]['E'] , $IDs["$s,$q"]['strand'] , '.' , "ID=$s.$count;Name=$q($evalue)")) . "\n";
+    $child = join ("\n",$IDs[$sq]['HSPs']) . "\n";
+    fwrite($gff,$parent);
+    fwrite($gff,$child);
+  }
+  fclose($tsv);
+  fclose($gff);
+}