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