src/jalview/hmmer/HMMSearch.java

   1 package jalview.hmmer;
   2
   3 import jalview.datamodel.Alignment;
   4 import jalview.datamodel.AlignmentAnnotation;
   5 import jalview.datamodel.AlignmentI;
   6 import jalview.datamodel.Annotation;
   7 import jalview.datamodel.HiddenMarkovModel;
   8 import jalview.datamodel.SequenceI;
   9 import jalview.gui.AlignFrame;
  10 import jalview.gui.JvOptionPane;
  11 import jalview.io.DataSourceType;
  12 import jalview.io.FileParse;
  13 import jalview.io.StockholmFile;
  14 import jalview.util.FileUtils;
  15 import jalview.util.MessageManager;
  16 import jalview.ws.params.ArgumentI;
  17 import jalview.ws.params.simple.BooleanOption;
  18
  19 import java.io.BufferedReader;
  20 import java.io.File;
  21 import java.io.FileReader;
  22 import java.io.IOException;
  23 import java.util.ArrayList;
  24 import java.util.List;
  25 import java.util.Scanner;
  26
  27 import javax.swing.JOptionPane;
  28
  29 public class HMMSearch extends HmmerCommand
  30 {
  31   static final String HMMSEARCH = "hmmsearch";
  32
  33   boolean realign = false;
  34
  35   boolean trim = false;
  36
  37   int seqsToReturn = Integer.MAX_VALUE;
  38
  39   SequenceI[] seqs;
  40
  41   /**
  42    * Constructor for the HMMSearchThread
  43    *
  44    * @param af
  45    */
  46   public HMMSearch(AlignFrame af, List<ArgumentI> args)
  47   {
  48     super(af, args);
  49   }
  50
  51   /**
  52    * Runs the HMMSearchThread: the data on the alignment or group is exported,
  53    * then the command is executed in the command line and then the data is
  54    * imported and displayed in a new frame. Call this method directly to execute
  55    * synchronously, or via start() in a new Thread for asynchronously.
  56    */
  57   @Override
  58   public void run()
  59   {
  60     HiddenMarkovModel hmm = af.getSelectedHMM();
  61     if (hmm == null)
  62     {
  63       JOptionPane.showMessageDialog(af,
  64               MessageManager.getString("warn.no_selected_hmm"));
  65       return;
  66     }
  67
  68     SequenceI hmmSeq = af.getSelectedHMMSequence();
  69     long msgId = System.currentTimeMillis();
  70     af.setProgressBar(MessageManager.getString("status.running_hmmsearch"),
  71             msgId);
  72
  73     try
  74     {
  75       File hmmFile = FileUtils.createTempFile("hmm", ".hmm");
  76       File hitsAlignmentFile = FileUtils.createTempFile("hitAlignment",
  77               ".sto");
  78       File searchOutputFile = FileUtils.createTempFile("searchOutput",
  79               ".sto");
  80
  81       exportHmm(hmm, hmmFile.getAbsoluteFile());
  82
  83       boolean ran = runCommand(searchOutputFile, hitsAlignmentFile, hmmFile);
  84       if (!ran)
  85       {
  86         JvOptionPane.showInternalMessageDialog(af,
  87                 MessageManager.getString("warn.hmmsearch_failed"));
  88         return;
  89       }
  90
  91       importData(hmmSeq, hitsAlignmentFile, hmmFile, searchOutputFile);
  92       // TODO make realignment of search results a step at this level
  93       // and make it conditional on this.realign
  94     } catch (IOException | InterruptedException e)
  95     {
  96       e.printStackTrace();
  97     }
  98     finally
  99     {
 100       af.setProgressBar("", msgId);
 101     }
 102   }
 103
 104   /**
 105    * Executes an hmmsearch with the given hmm as input. The database to be
 106    * searched is a local file as specified by the 'Database' parameter, or the
 107    * current alignment (written to file) if none is specified.
 108    *
 109    * @param searchOutputFile
 110    * @param hitsAlignmentFile
 111    * @param hmmFile
 112    *
 113    * @return
 114    * @throws IOException
 115    */
 116   private boolean runCommand(File searchOutputFile, File hitsAlignmentFile,
 117           File hmmFile) throws IOException
 118   {
 119     String command = getCommandPath(HMMSEARCH);
 120     if (command == null)
 121     {
 122       return false;
 123     }
 124
 125     List<String> args = new ArrayList<>();
 126     args.add(command);
 127     args.add("-o");
 128     args.add(getFilePath(searchOutputFile));
 129     args.add("-A");
 130     args.add(getFilePath(hitsAlignmentFile));
 131
 132     boolean dbFound = false;
 133     String dbPath = "";
 134     File databaseFile = null;
 135
 136     if (params != null)
 137     {
 138       for (ArgumentI arg : params)
 139       {
 140         String name = arg.getName();
 141         if (MessageManager.getString("label.number_of_results")
 142                 .equals(name))
 143         {
 144           seqsToReturn = Integer.parseInt(arg.getValue());
 145         }
 146         else if (MessageManager.getString("label.auto_align_seqs")
 147                 .equals(name))
 148         {
 149           realign = true; // TODO: not used
 150         }
 151         else if (MessageManager.getString("label.use_accessions")
 152                 .equals(name))
 153         {
 154           args.add("--acc");
 155         }
 156         else if (MessageManager.getString("label.seq_e_value").equals(name))
 157         {
 158           args.add("--incE");
 159           args.add(arg.getValue());
 160         }
 161         else if (MessageManager.getString("label.seq_score").equals(name))
 162         {
 163           args.add("-incT");
 164           args.add(arg.getValue());
 165         }
 166         else if (MessageManager.getString("label.dom_e_value_desc")
 167                 .equals(name))
 168         {
 169           args.add("--incdomE");
 170           args.add(arg.getValue());
 171         }
 172         else if (MessageManager.getString("label.dom_score").equals(name))
 173         {
 174           args.add("--incdomT");
 175           args.add(arg.getValue());
 176         }
 177         else if (MessageManager.getString("label.trim_termini")
 178                 .equals(name))
 179         {
 180           trim = true;
 181         }
 182         else if (MessageManager.getString("label.database").equals(name))
 183         {
 184           dbFound = true;
 185           dbPath = arg.getValue();
 186           if (!MessageManager.getString("label.this_alignment")
 187                   .equals(dbPath))
 188           {
 189             databaseFile = new File(dbPath);
 190           }
 191         }
 192       }
 193     }
 194
 195     if (!dbFound || MessageManager.getString("label.this_alignment")
 196             .equals(dbPath))
 197     {
 198       /*
 199        * no external database specified for search, so
 200        * export current alignment as 'database' to search
 201        */
 202       databaseFile = FileUtils.createTempFile("database", ".sto");
 203       AlignmentI al = af.getViewport().getAlignment();
 204       AlignmentI copy = new Alignment(al);
 205       SequenceI hmms = copy.getHmmConsensus();
 206       if (hmms != null)
 207       {
 208         copy.deleteSequence(hmms);
 209       }
 210       exportStockholm(copy.getSequencesArray(), databaseFile, null);
 211       // StockholmFile stoFile = new StockholmFile(copy);
 212       // stoFile.setSeqs(copy.getSequencesArray());
 213       // String alignmentString = stoFile.print();
 214       // PrintWriter writer = new PrintWriter(databaseFile);
 215       // writer.print(alignmentString);
 216       // writer.close();
 217     }
 218
 219     args.add(getFilePath(hmmFile));
 220     args.add(getFilePath(databaseFile));
 221
 222     return runCommand(args);
 223   }
 224
 225   /**
 226    * Imports the data from the temporary file to which the output of hmmsearch
 227    * is directed.
 228    *
 229    * @param hmmSeq
 230    */
 231   private void importData(SequenceI hmmSeq, File inputAlignmentTemp,
 232           File hmmTemp, File searchOutputFile)
 233           throws IOException, InterruptedException
 234   {
 235     BufferedReader br = new BufferedReader(
 236             new FileReader(inputAlignmentTemp));
 237     try
 238     {
 239       if (br.readLine() == null)
 240       {
 241         JOptionPane.showMessageDialog(af,
 242                 MessageManager.getString("label.no_sequences_found"));
 243         return;
 244       }
 245       StockholmFile file = new StockholmFile(new FileParse(
 246               inputAlignmentTemp.getAbsolutePath(), DataSourceType.FILE));
 247       seqs = file.getSeqsAsArray();
 248
 249       readTable(searchOutputFile);
 250
 251       int seqCount = Math.min(seqs.length, seqsToReturn);
 252       SequenceI[] hmmAndSeqs = new SequenceI[seqCount + 1];
 253       hmmAndSeqs[0] = hmmSeq;
 254       System.arraycopy(seqs, 0, hmmAndSeqs, 1, seqCount);
 255
 256       AlignmentI alignment = new Alignment(hmmAndSeqs);
 257       AlignFrame frame = new AlignFrame(alignment, 1, 1);
 258       frame.setSelectedHMMSequence(hmmSeq);
 259       List<ArgumentI> alignArgs = new ArrayList<>();
 260       if (trim)
 261       {
 262         alignArgs.add(new BooleanOption(
 263                 MessageManager.getString("label.trim_termini"),
 264                 MessageManager.getString("label.trim_termini_desc"), true,
 265                 true, true, null));
 266       }
 267       HmmerCommand hmmalign = new HMMAlign(frame, alignArgs);
 268       hmmalign.run();
 269       frame = null;
 270       hmmTemp.delete();
 271       inputAlignmentTemp.delete();
 272       searchOutputFile.delete();
 273     } finally
 274     {
 275       if (br != null)
 276       {
 277         br.close();
 278       }
 279     }
 280   }
 281
 282   void readTable(File inputTableTemp) throws IOException
 283   {
 284     BufferedReader br = new BufferedReader(new FileReader(inputTableTemp));
 285     String line = "";
 286     while (!line.startsWith("Query:"))
 287     {
 288       line = br.readLine();
 289     }
 290     for (int i = 0; i < 5; i++)
 291     {
 292       line = br.readLine();
 293     }
 294
 295     int index = 0;
 296     while (!"  ------ inclusion threshold ------".equals(line)
 297             && !"".equals(line))
 298     {
 299       Scanner scanner = new Scanner(line);
 300
 301       String str = scanner.next(); // full sequence eValue score
 302       float eValue = Float.parseFloat(str);
 303       int seqLength = seqs[index].getLength();
 304       Annotation[] annots = new Annotation[seqLength];
 305       for (int j = 0; j < seqLength; j++)
 306       {
 307         annots[j] = new Annotation(eValue);
 308       }
 309       AlignmentAnnotation annot = new AlignmentAnnotation("E-value",
 310               "Score", annots);
 311       annot.setScore(Double.parseDouble(str));
 312       annot.setSequenceRef(seqs[index]);
 313       seqs[index].addAlignmentAnnotation(annot);
 314
 315       scanner.close();
 316       line = br.readLine();
 317       index++;
 318     }
 319
 320     br.close();
 321   }
 322
 323 }