3 import jalview.bin.Cache;
4 import jalview.datamodel.Alignment;
5 import jalview.datamodel.AlignmentAnnotation;
6 import jalview.datamodel.AlignmentI;
7 import jalview.datamodel.Annotation;
8 import jalview.datamodel.HiddenMarkovModel;
9 import jalview.datamodel.SequenceI;
10 import jalview.gui.AlignFrame;
11 import jalview.gui.JvOptionPane;
12 import jalview.io.DataSourceType;
13 import jalview.io.FileParse;
14 import jalview.io.StockholmFile;
15 import jalview.util.FileUtils;
16 import jalview.util.MessageManager;
17 import jalview.ws.params.ArgumentI;
18 import jalview.ws.params.simple.BooleanOption;
19 import jalview.ws.params.simple.Option;
21 import java.io.BufferedReader;
23 import java.io.FileReader;
24 import java.io.IOException;
25 import java.util.ArrayList;
26 import java.util.Collections;
27 import java.util.List;
28 import java.util.Scanner;
30 import javax.swing.JOptionPane;
32 public class HMMSearch extends HmmerCommand
34 static final String HMMSEARCH = "hmmsearch";
37 * constants for i18n lookup of passed parameter names
39 static final String DATABASE_KEY = "label.database";
41 static final String THIS_ALIGNMENT_KEY = "label.this_alignment";
43 static final String USE_ACCESSIONS_KEY = "label.use_accessions";
45 static final String AUTO_ALIGN_SEQS_KEY = "label.auto_align_seqs";
47 static final String NUMBER_OF_RESULTS_KEY = "label.number_of_results";
49 static final String TRIM_TERMINI_KEY = "label.trim_termini";
51 static final String REPORTING_CUTOFF_KEY = "label.reporting_cutoff";
53 static final String CUTOFF_NONE = "None";
55 static final String CUTOFF_SCORE = "Score";
57 static final String CUTOFF_EVALUE = "E-Value";
59 static final String SEQ_EVALUE_KEY = "label.seq_evalue";
61 static final String DOM_EVALUE_KEY = "label.dom_evalue";
63 static final String SEQ_SCORE_KEY = "label.seq_score";
65 static final String DOM_SCORE_KEY = "label.dom_score";
67 boolean realign = false;
71 int seqsToReturn = Integer.MAX_VALUE;
76 * Constructor for the HMMSearchThread
80 public HMMSearch(AlignFrame af, List<ArgumentI> args)
86 * Runs the HMMSearchThread: the data on the alignment or group is exported,
87 * then the command is executed in the command line and then the data is
88 * imported and displayed in a new frame. Call this method directly to execute
89 * synchronously, or via start() in a new Thread for asynchronously.
94 HiddenMarkovModel hmm = getHmmProfile();
97 // shouldn't happen if we got this far
98 Cache.log.error("Error: no hmm for hmmsearch");
102 SequenceI hmmSeq = hmm.getConsensusSequence();// af.getSelectedHMMSequence();
103 long msgId = System.currentTimeMillis();
104 af.setProgressBar(MessageManager.getString("status.running_hmmsearch"),
109 File hmmFile = FileUtils.createTempFile("hmm", ".hmm");
110 File hitsAlignmentFile = FileUtils.createTempFile("hitAlignment",
112 File searchOutputFile = FileUtils.createTempFile("searchOutput",
115 exportHmm(hmm, hmmFile.getAbsoluteFile());
117 boolean ran = runCommand(searchOutputFile, hitsAlignmentFile, hmmFile);
120 JvOptionPane.showInternalMessageDialog(af, MessageManager
121 .formatMessage("warn.command_failed", "hmmsearch"));
125 importData(hmmSeq, hitsAlignmentFile, hmmFile, searchOutputFile);
126 // TODO make realignment of search results a step at this level
127 // and make it conditional on this.realign
128 } catch (IOException | InterruptedException e)
134 af.setProgressBar("", msgId);
139 * Executes an hmmsearch with the given hmm as input. The database to be
140 * searched is a local file as specified by the 'Database' parameter, or the
141 * current alignment (written to file) if none is specified.
143 * @param searchOutputFile
144 * @param hitsAlignmentFile
148 * @throws IOException
150 private boolean runCommand(File searchOutputFile, File hitsAlignmentFile,
151 File hmmFile) throws IOException
153 String command = getCommandPath(HMMSEARCH);
159 List<String> args = new ArrayList<>();
162 args.add(getFilePath(searchOutputFile));
164 args.add(getFilePath(hitsAlignmentFile));
166 boolean dbFound = false;
168 File databaseFile = null;
170 boolean useEvalueCutoff = false;
171 boolean useScoreCutoff = false;
172 String seqEvalueCutoff = null;
173 String domEvalueCutoff = null;
174 String seqScoreCutoff = null;
175 String domScoreCutoff = null;
179 for (ArgumentI arg : params)
181 String name = arg.getName();
182 if (MessageManager.getString(NUMBER_OF_RESULTS_KEY)
185 seqsToReturn = Integer.parseInt(arg.getValue());
187 else if (MessageManager.getString(AUTO_ALIGN_SEQS_KEY)
190 realign = true; // TODO: not used
192 else if (MessageManager.getString(USE_ACCESSIONS_KEY)
197 else if (MessageManager.getString(REPORTING_CUTOFF_KEY)
200 if (CUTOFF_EVALUE.equals(arg.getValue()))
202 useEvalueCutoff = true;
204 else if (CUTOFF_SCORE.equals(arg.getValue()))
206 useScoreCutoff = true;
209 else if (MessageManager.getString(SEQ_EVALUE_KEY).equals(name))
211 seqEvalueCutoff = arg.getValue();
213 else if (MessageManager.getString(SEQ_SCORE_KEY).equals(name))
215 seqScoreCutoff = arg.getValue();
217 else if (MessageManager.getString(DOM_EVALUE_KEY)
220 domEvalueCutoff = arg.getValue();
222 else if (MessageManager.getString(DOM_SCORE_KEY).equals(name))
224 domScoreCutoff = arg.getValue();
226 else if (MessageManager.getString(TRIM_TERMINI_KEY)
231 else if (MessageManager.getString(DATABASE_KEY).equals(name))
234 dbPath = arg.getValue();
235 if (!MessageManager.getString(THIS_ALIGNMENT_KEY)
238 databaseFile = new File(dbPath);
247 args.add(seqEvalueCutoff);
249 args.add(domEvalueCutoff);
251 else if (useScoreCutoff)
254 args.add(seqScoreCutoff);
256 args.add(domScoreCutoff);
259 if (!dbFound || MessageManager.getString(THIS_ALIGNMENT_KEY)
263 * no external database specified for search, so
264 * export current alignment as 'database' to search,
265 * excluding any HMM consensus sequences it contains
267 databaseFile = FileUtils.createTempFile("database", ".sto");
268 AlignmentI al = af.getViewport().getAlignment();
269 AlignmentI copy = new Alignment(al);
270 List<SequenceI> hmms = copy.getHmmSequences();
271 for (SequenceI hmmSeq : hmms)
273 copy.deleteSequence(hmmSeq);
275 exportStockholm(copy.getSequencesArray(), databaseFile, null);
278 args.add(getFilePath(hmmFile));
279 args.add(getFilePath(databaseFile));
281 return runCommand(args);
285 * Imports the data from the temporary file to which the output of hmmsearch
290 private void importData(SequenceI hmmSeq, File inputAlignmentTemp,
291 File hmmTemp, File searchOutputFile)
292 throws IOException, InterruptedException
294 BufferedReader br = new BufferedReader(
295 new FileReader(inputAlignmentTemp));
298 if (br.readLine() == null)
300 JOptionPane.showMessageDialog(af,
301 MessageManager.getString("label.no_sequences_found"));
304 StockholmFile file = new StockholmFile(new FileParse(
305 inputAlignmentTemp.getAbsolutePath(), DataSourceType.FILE));
306 seqs = file.getSeqsAsArray();
308 readTable(searchOutputFile);
310 int seqCount = Math.min(seqs.length, seqsToReturn);
311 SequenceI[] hmmAndSeqs = new SequenceI[seqCount + 1];
312 hmmAndSeqs[0] = hmmSeq;
313 System.arraycopy(seqs, 0, hmmAndSeqs, 1, seqCount);
316 * and align the search results to the HMM profile
318 AlignmentI al = new Alignment(hmmAndSeqs);
319 AlignFrame frame = new AlignFrame(al, 1, 1);
320 List<ArgumentI> alignArgs = new ArrayList<>();
321 String defSeq = hmmSeq.getName();
322 List<String> options = Collections.singletonList(defSeq);
323 Option option = new Option(MessageManager.getString("label.use_hmm"),
324 "", true, defSeq, defSeq, options, null);
325 alignArgs.add(option);
328 alignArgs.add(new BooleanOption(
329 MessageManager.getString(TRIM_TERMINI_KEY),
330 MessageManager.getString("label.trim_termini_desc"), true,
333 HmmerCommand hmmalign = new HMMAlign(frame, alignArgs);
337 inputAlignmentTemp.delete();
338 searchOutputFile.delete();
348 void readTable(File inputTableTemp) throws IOException
350 BufferedReader br = new BufferedReader(new FileReader(inputTableTemp));
352 while (!line.startsWith("Query:"))
354 line = br.readLine();
356 for (int i = 0; i < 5; i++)
358 line = br.readLine();
362 while (!" ------ inclusion threshold ------".equals(line)
365 Scanner scanner = new Scanner(line);
367 String str = scanner.next(); // full sequence eValue score
368 float eValue = Float.parseFloat(str);
369 int seqLength = seqs[index].getLength();
370 Annotation[] annots = new Annotation[seqLength];
371 for (int j = 0; j < seqLength; j++)
373 annots[j] = new Annotation(eValue);
375 AlignmentAnnotation annot = new AlignmentAnnotation("E-value",
377 annot.setScore(Double.parseDouble(str));
378 annot.setSequenceRef(seqs[index]);
379 seqs[index].addAlignmentAnnotation(annot);
382 line = br.readLine();