X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Fhmmer%2FHMMSearch.java;h=767a56ed4e57995bb98c407c809c08ef6106db43;hb=bd306b69256571ae575562d3079c842fa3be43f5;hp=7413d14ddaf341f0f812f6421a2bd2509ce8752d;hpb=d31b51985d01217340aa5f6470d3fd3c2314e3eb;p=jalview.git diff --git a/src/jalview/hmmer/HMMSearch.java b/src/jalview/hmmer/HMMSearch.java index 7413d14..767a56e 100644 --- a/src/jalview/hmmer/HMMSearch.java +++ b/src/jalview/hmmer/HMMSearch.java @@ -1,6 +1,7 @@ package jalview.hmmer; import jalview.bin.Cache; +import jalview.bin.Console; import jalview.datamodel.Alignment; import jalview.datamodel.AlignmentAnnotation; import jalview.datamodel.AlignmentI; @@ -8,6 +9,7 @@ import jalview.datamodel.Annotation; import jalview.datamodel.HiddenMarkovModel; import jalview.datamodel.SequenceI; import jalview.gui.AlignFrame; +import jalview.gui.Desktop; import jalview.gui.JvOptionPane; import jalview.io.DataSourceType; import jalview.io.FileParse; @@ -25,23 +27,20 @@ import java.io.IOException; import java.util.ArrayList; import java.util.Collections; import java.util.List; -import java.util.Scanner; import javax.swing.JOptionPane; -public class HMMSearch extends HmmerCommand +public class HMMSearch extends Search { - private static final String PARAMNAME_NO_OF_RESULTS = MessageManager.getString("label.number_of_results"); - - static final String HMMSEARCH = "hmmsearch"; boolean realign = false; boolean trim = false; + boolean returnNoOfNewSeqs = false; + int seqsToReturn = Integer.MAX_VALUE; - SequenceI[] seqs; /** * Constructor for the HMMSearchThread @@ -66,13 +65,13 @@ public class HMMSearch extends HmmerCommand if (hmm == null) { // shouldn't happen if we got this far - Cache.log.error("Error: no hmm for hmmsearch"); + Console.error("Error: no hmm for hmmsearch"); return; } - SequenceI hmmSeq = hmm.getConsensusSequence();// af.getSelectedHMMSequence(); + SequenceI hmmSeq = hmm.getConsensusSequence(); long msgId = System.currentTimeMillis(); - af.setProgressBar(MessageManager.getString("status.running_hmmsearch"), + af.setProgressBar(MessageManager.getString("status.running_search"), msgId); try @@ -129,108 +128,15 @@ public class HMMSearch extends HmmerCommand List args = new ArrayList<>(); args.add(command); - args.add("-o"); - args.add(getFilePath(searchOutputFile)); - args.add("-A"); - args.add(getFilePath(hitsAlignmentFile)); - - boolean dbFound = false; - String dbPath = ""; - File databaseFile = null; - - if (params != null) - { - for (ArgumentI arg : params) - { - String name = arg.getName(); - if (MessageManager.getString("label.number_of_results") - .equals(name)) - { - seqsToReturn = Integer.parseInt(arg.getValue()); - } - else if (MessageManager.getString("label.auto_align_seqs") - .equals(name)) - { - realign = true; // TODO: not used - } - else if (MessageManager.getString("label.use_accessions") - .equals(name)) - { - args.add("--acc"); - } - else if (MessageManager.getString("label.seq_e_value").equals(name)) - { - args.add("--incE"); - args.add(arg.getValue()); - } - else if (MessageManager.getString("label.seq_score").equals(name)) - { - args.add("-incT"); - args.add(arg.getValue()); - } - else if (MessageManager.getString("label.dom_e_value_desc") - .equals(name)) - { - args.add("--incdomE"); - args.add(arg.getValue()); - } - else if (MessageManager.getString("label.dom_score").equals(name)) - { - args.add("--incdomT"); - args.add(arg.getValue()); - } - else if (MessageManager.getString("label.trim_termini") - .equals(name)) - { - trim = true; - } - else if (MessageManager.getString("label.database").equals(name)) - { - dbFound = true; - dbPath = arg.getValue(); - if (!MessageManager.getString("label.this_alignment") - .equals(dbPath)) - { - databaseFile = new File(dbPath); - } - } - } - } - - if (!dbFound || MessageManager.getString("label.this_alignment") - .equals(dbPath)) - { - /* - * no external database specified for search, so - * export current alignment as 'database' to search, - * excluding any HMM consensus sequences it contains - */ - databaseFile = FileUtils.createTempFile("database", ".sto"); - AlignmentI al = af.getViewport().getAlignment(); - AlignmentI copy = new Alignment(al); - List hmms = copy.getHmmSequences(); - for (SequenceI hmmSeq : hmms) - { - copy.deleteSequence(hmmSeq); - } - exportStockholm(copy.getSequencesArray(), databaseFile, null); - // StockholmFile stoFile = new StockholmFile(copy); - // stoFile.setSeqs(copy.getSequencesArray()); - // String alignmentString = stoFile.print(); - // PrintWriter writer = new PrintWriter(databaseFile); - // writer.print(alignmentString); - // writer.close(); - } - - args.add(getFilePath(hmmFile)); - args.add(getFilePath(databaseFile)); + buildArguments(args, searchOutputFile, hitsAlignmentFile, hmmFile); return runCommand(args); } + /** * Imports the data from the temporary file to which the output of hmmsearch - * is directed. + * was directed. The results are optionally realigned using hmmalign. * * @param hmmSeq */ @@ -252,34 +158,95 @@ public class HMMSearch extends HmmerCommand inputAlignmentTemp.getAbsolutePath(), DataSourceType.FILE)); seqs = file.getSeqsAsArray(); - readTable(searchOutputFile); + readDomainTable(searchOutputFile, false); + + if (searchAlignment) + { + recoverSequences(sequencesHash, seqs); + } + + // look for PP cons and ref seq in alignment only annotation + AlignmentAnnotation modelpos = null, ppcons = null; + for (AlignmentAnnotation aa : file.getAnnotations()) + { + if (aa.sequenceRef == null) + { + if (aa.label.equals("Reference Positions")) // RF feature type in + // stockholm parser + { + modelpos = aa; + } + if (aa.label.equals("Posterior Probability")) + { + ppcons = aa; + } + } + } + int seqCount = Math.min(seqs.length, seqsToReturn); SequenceI[] hmmAndSeqs = new SequenceI[seqCount + 1]; + hmmSeq = hmmSeq.deriveSequence(); // otherwise all bad things happen hmmAndSeqs[0] = hmmSeq; System.arraycopy(seqs, 0, hmmAndSeqs, 1, seqCount); + if (modelpos != null) + { + // TODO need - get ungapped sequence method + hmmSeq.setSequence( + hmmSeq.getDatasetSequence().getSequenceAsString()); + Annotation[] refpos = modelpos.annotations; + // insert gaps to match with refseq positions + int gc = 0, lcol = 0; + for (int c = 0; c < refpos.length; c++) + { + if (refpos[c] != null && ("x".equals(refpos[c].displayCharacter))) + { + if (gc > 0) + { + hmmSeq.insertCharAt(lcol + 1, gc, '-'); + } + gc = 0; + lcol = c; + } + else + { + gc++; + } + } + } - /* - * and align the search results to the HMM profile - */ - AlignmentI al = new Alignment(hmmAndSeqs); - AlignFrame frame = new AlignFrame(al, 1, 1); - List alignArgs = new ArrayList<>(); - String defSeq = hmmSeq.getName(); - List options = Collections.singletonList(defSeq); - Option option = new Option(MessageManager.getString("label.use_hmm"), - "", true, defSeq, defSeq, options, null); - alignArgs.add(option); - if (trim) + if (realign) + { + realignResults(hmmAndSeqs); + } + else { - alignArgs.add(new BooleanOption( - MessageManager.getString("label.trim_termini"), - MessageManager.getString("label.trim_termini_desc"), true, - true, true, null)); + AlignmentI al = new Alignment(hmmAndSeqs); + if (ppcons != null) + { + al.addAnnotation(ppcons); + } + if (modelpos != null) + { + al.addAnnotation(modelpos); + } + AlignFrame alignFrame = new AlignFrame(al, AlignFrame.DEFAULT_WIDTH, + AlignFrame.DEFAULT_HEIGHT); + String ttl = "hmmSearch of " + databaseName + " using " + + hmmSeq.getName(); + Desktop.addInternalFrame(alignFrame, ttl, AlignFrame.DEFAULT_WIDTH, + AlignFrame.DEFAULT_HEIGHT); + + if (returnNoOfNewSeqs) + { + int nNew = checkForNewSequences(); + JvOptionPane.showMessageDialog(af.alignPanel, nNew + " " + + MessageManager.getString("label.new_returned")); + } + } - HmmerCommand hmmalign = new HMMAlign(frame, alignArgs); - hmmalign.run(); - frame = null; + + hmmTemp.delete(); inputAlignmentTemp.delete(); searchOutputFile.delete(); @@ -292,45 +259,61 @@ public class HMMSearch extends HmmerCommand } } - void readTable(File inputTableTemp) throws IOException + private int checkForNewSequences() { - BufferedReader br = new BufferedReader(new FileReader(inputTableTemp)); - String line = ""; - while (!line.startsWith("Query:")) + int nNew = seqs.length; + + for (SequenceI resultSeq : seqs) { - line = br.readLine(); + for (SequenceI aliSeq : alignment.getSequencesArray()) + { + if (resultSeq.getName().equals(aliSeq.getName())) + { + nNew--; + break; + } + } } - for (int i = 0; i < 5; i++) + + return nNew; + + } + + /** + * Realigns the given sequences using hmmalign, to the HMM profile sequence + * which is the first in the array, and opens the results in a new frame + * + * @param hmmAndSeqs + */ + protected void realignResults(SequenceI[] hmmAndSeqs) + { + /* + * and align the search results to the HMM profile + */ + AlignmentI al = new Alignment(hmmAndSeqs); + AlignFrame frame = new AlignFrame(al, 1, 1); + List alignArgs = new ArrayList<>(); + String alignTo = hmmAndSeqs[0].getName(); + List options = Collections.singletonList(alignTo); + Option option = new Option(MessageManager.getString("label.use_hmm"), + "", true, alignTo, alignTo, options, null); + alignArgs.add(option); + if (trim) { - line = br.readLine(); + alignArgs.add(new BooleanOption( + MessageManager.getString(TRIM_TERMINI_KEY), + MessageManager.getString("label.trim_termini_desc"), true, + true, true, null)); } + HmmerCommand hmmalign = new HMMAlign(frame, alignArgs); + hmmalign.run(); - int index = 0; - while (!" ------ inclusion threshold ------".equals(line) - && !"".equals(line)) + if (returnNoOfNewSeqs) { - Scanner scanner = new Scanner(line); - - String str = scanner.next(); // full sequence eValue score - float eValue = Float.parseFloat(str); - int seqLength = seqs[index].getLength(); - Annotation[] annots = new Annotation[seqLength]; - for (int j = 0; j < seqLength; j++) - { - annots[j] = new Annotation(eValue); - } - AlignmentAnnotation annot = new AlignmentAnnotation("E-value", - "Score", annots); - annot.setScore(Double.parseDouble(str)); - annot.setSequenceRef(seqs[index]); - seqs[index].addAlignmentAnnotation(annot); - - scanner.close(); - line = br.readLine(); - index++; + int nNew = checkForNewSequences(); + JvOptionPane.showMessageDialog(frame.alignPanel, + nNew + " " + MessageManager.getString("label.new_returned")); } - - br.close(); } }