src/jalview/io/HMMFile.java

   1 package jalview.io;
   2
   3 import jalview.api.AlignExportSettingI;
   4 import jalview.api.AlignmentViewPanel;
   5 import jalview.datamodel.HMMNode;
   6 import jalview.datamodel.HiddenMarkovModel;
   7 import jalview.datamodel.SequenceI;
   8
   9 import java.io.BufferedReader;
  10 import java.io.IOException;
  11 import java.util.ArrayList;
  12 import java.util.List;
  13 import java.util.Scanner;
  14
  15
  16 /**
  17  * Adds capability to read in and write out HMMER3 files. .
  18  *
  19  *
  20  * @author TZVanaalten
  21  *
  22  */
  23 public class HMMFile extends AlignFile
  24         implements AlignmentFileReaderI, AlignmentFileWriterI
  25 {
  26   private static final int NUMBER_OF_TRANSITIONS = 7;
  27
  28   private static final String SPACE = " ";
  29
  30   private static final String COMPO = "COMPO";
  31
  32   private static final String EMPTY = "";
  33
  34   /*
  35    * guide line added to an output HMMER file, purely for readability
  36    */
  37   private static final String TRANSITIONTYPELINE = "            m->m     m->i     m->d     i->m     i->i     d->m     d->d";
  38
  39   private static String NL = "\n";
  40
  41   private HiddenMarkovModel hmm;
  42
  43   // number of symbols in the alphabet used in the hidden Markov model
  44   int numberOfSymbols;
  45
  46   /**
  47    * Parses immediately.
  48    *
  49    * @param inFile
  50    * @param type
  51    * @throws IOException
  52    */
  53   public HMMFile(String inFile, DataSourceType type) throws IOException
  54   {
  55     super(inFile, type);
  56   }
  57
  58   /**
  59    * Parses immediately.
  60    *
  61    * @param source
  62    * @throws IOException
  63    */
  64   public HMMFile(FileParse source) throws IOException
  65   {
  66     super(source);
  67   }
  68
  69   /**
  70    * Default constructor, do not use!
  71    */
  72   public HMMFile()
  73   {
  74
  75   }
  76
  77   /**
  78    * Constructor for HMMFile used for exporting.
  79    *
  80    * @param hmm
  81    * @param exportImmediately
  82    */
  83   public HMMFile(HiddenMarkovModel markov)
  84   {
  85     hmm = markov;
  86   }
  87
  88   /**
  89    * For testing, do not use.
  90    *
  91    * @param br
  92    */
  93   HMMFile(BufferedReader br)
  94   {
  95     dataIn = br;
  96   }
  97
  98   /**
  99    * Returns the HMM produced by reading in a HMMER3 file.
 100    *
 101    * @return
 102    */
 103   public HiddenMarkovModel getHMM()
 104   {
 105     return hmm;
 106   }
 107
 108   /**
 109    * Sets the HMM used in this file.
 110    *
 111    * @param model
 112    */
 113   public void setHMM(HiddenMarkovModel model)
 114   {
 115     this.hmm = model;
 116   }
 117
 118   /**
 119    * Gets the name of the hidden Markov model.
 120    *
 121    * @return
 122    */
 123   public String getName()
 124   {
 125     return hmm.getName();
 126   }
 127
 128   /**
 129    * Reads the data from HMM file into the HMM field on this object.
 130    *
 131    * @throws IOException
 132    */
 133   @Override
 134   public void parse() throws IOException
 135   {
 136     try
 137     {
 138       hmm = new HiddenMarkovModel();
 139       parseFileProperties(dataIn);
 140       parseModel(dataIn);
 141     } catch (Exception e)
 142     {
 143       e.printStackTrace();
 144     }
 145   }
 146
 147   /**
 148    * Reads the data from HMM file into the HMM field on this object.
 149    *
 150    * @throws IOException
 151    */
 152
 153   public void parse(BufferedReader br) throws IOException
 154   {
 155     hmm = new HiddenMarkovModel();
 156     parseFileProperties(br);
 157     parseModel(br);
 158   }
 159
 160
 161
 162   /**
 163    * Imports the file properties from a HMMER3 file.
 164    *
 165    * @param input
 166    *          The buffered reader used to read in the file.
 167    * @throws IOException
 168    */
 169   void parseFileProperties(BufferedReader input) throws IOException
 170   {
 171     boolean readingFile = true;
 172     hmm.setFileHeader(input.readLine());
 173     String line = input.readLine();
 174     while (readingFile)
 175     {
 176       if (line != null)
 177       {
 178         Scanner parser = new Scanner(line);
 179         String next = parser.next();
 180         if ("HMM".equals(next)) // indicates start of HMM data (end of file
 181                               // properties)
 182         {
 183           readingFile = false;
 184           fillSymbols(parser);
 185           numberOfSymbols = hmm.getNumberOfSymbols();
 186         }
 187         else if ("STATS".equals(next))
 188         {
 189           parser.next();
 190           String key;
 191           String value;
 192           key = parser.next();
 193           value = parser.next() + SPACE + SPACE + parser.next();
 194           hmm.addFileProperty(key, value);
 195         }
 196         else
 197         {
 198           String key = next;
 199           String value = parser.next();
 200           while (parser.hasNext())
 201           {
 202             value = value + SPACE + parser.next();
 203           }
 204           hmm.addFileProperty(key, value);
 205         }
 206         parser.close();
 207       }
 208       line = input.readLine();
 209       if (line == null)
 210       {
 211         readingFile = false;
 212       }
 213     }
 214
 215   }
 216
 217   /**
 218    * Parses the model data from the HMMER3 file
 219    *
 220    * @param input
 221    *          The buffered reader used to read the file.
 222    * @throws IOException
 223    */
 224   void parseModel(BufferedReader input) throws IOException
 225   {
 226     boolean first = true;
 227     String line = input.readLine();
 228     while (!"//".equals(line))
 229     {
 230       HMMNode node = new HMMNode();
 231       hmm.getNodes().add(node);
 232       Scanner matchReader = new Scanner(line);
 233       String next = matchReader.next();
 234       if (next.equals(COMPO) || !first)
 235       {
 236         // stores match emission line in list
 237         List<Double> matches = new ArrayList<>();
 238         matches = fillList(matchReader, numberOfSymbols);
 239         node.setMatchEmissions(matches);
 240         if (!first)
 241         {
 242           int column = parseAnnotations(matchReader, node);
 243           hmm.setAlignmentColumn(node, column - 1);
 244         }
 245       }
 246       matchReader.close();
 247       // stores insert emission line in list
 248       line = input.readLine();
 249       Scanner insertReader = new Scanner(line);
 250       List<Double> inserts = new ArrayList<>();
 251       inserts = fillList(insertReader, numberOfSymbols);
 252       node.setInsertEmissions(inserts);
 253       insertReader.close();
 254
 255       // stores state transition line in list
 256       line = input.readLine();
 257       Scanner transitionReader = new Scanner(line);
 258       List<Double> transitions = new ArrayList<>();
 259       transitions = fillList(transitionReader, NUMBER_OF_TRANSITIONS);
 260       node.setStateTransitions(transitions);
 261       transitionReader.close();
 262       line = input.readLine();
 263
 264       first = false;
 265     }
 266   }
 267
 268   /**
 269    * Parses the annotations on the match emission line and add them to the node.
 270    * (See p109 of the HMMER User Guide (V3.1b2) for the specification.) Returns
 271    * the alignment column number (base 1) that the node maps to, if provided,
 272    * else zero.
 273    *
 274    * @param scanner
 275    * @param node
 276    */
 277   int parseAnnotations(Scanner scanner, HMMNode node)
 278   {
 279     /*
 280      * map from hmm node to alignment column index, if provided
 281      * HMM counts columns from 1, convert to base 0 for Jalview
 282      */
 283     int column = 0;
 284     if (hmm.mapIsActive() && scanner.hasNext())
 285     {
 286       column = scanner.nextInt();
 287       node.setAlignmentColumn(column - 1);
 288     }
 289     else
 290     {
 291       scanner.next();
 292     }
 293
 294     /*
 295      * hmm consensus residue if provided, else -
 296      */
 297     if (scanner.hasNext())
 298     {
 299       char consensusR;
 300       consensusR = charValue(scanner.next());
 301       node.setConsensusResidue(consensusR);
 302     }
 303
 304     /*
 305      * RF reference annotation, if provided, else -
 306      */
 307     if (scanner.hasNext())
 308     {
 309       char reference;
 310       reference = charValue(scanner.next());
 311       node.setReferenceAnnotation(reference);
 312     }
 313
 314     /*
 315      * 'm' for masked position, if provided, else -
 316      */
 317     if (scanner.hasNext())
 318     {
 319       char value;
 320       value = charValue(scanner.next());
 321       node.setMaskValue(value);
 322     }
 323
 324     /*
 325      * structure consensus symbol, if provided, else -
 326      */
 327     if (scanner.hasNext())
 328     {
 329       char consensusS;
 330       consensusS = charValue(scanner.next());
 331       node.setConsensusStructure(consensusS);
 332     }
 333
 334     return column;
 335   }
 336
 337   /**
 338    * Fills a list of doubles from an input line
 339    *
 340    * @param input
 341    *          The scanner for the line containing the data to be transferred to
 342    *          the list.
 343    * @param numberOfElements
 344    *          The number of elements in the list to be filled.
 345    * @return filled list Returns the list of doubles.
 346    * @throws IOException
 347    */
 348   static List<Double> fillList(Scanner input,
 349           int numberOfElements) throws IOException
 350   {
 351     List<Double> list = new ArrayList<>();
 352     for (int i = 0; i < numberOfElements; i++)
 353     {
 354
 355       String next = input.next();
 356       if (next.contains("*")) // state transitions to or from delete states
 357                               // occasionally have values of -infinity. These
 358                               // values are represented by an * in the .hmm
 359                               // file.
 360       {
 361         list.add(Double.NEGATIVE_INFINITY);
 362       }
 363       else
 364       {
 365         double prob = Double.valueOf(next);
 366         prob = Math.pow(Math.E, -prob);
 367         list.add(prob);
 368       }
 369     }
 370     if (list.size() < numberOfElements)
 371     {
 372       throw new IOException("Incomplete data");
 373     }
 374     return list;
 375   }
 376
 377   /**
 378    * Returns a string to be added to the StringBuilder containing the entire
 379    * output String.
 380    *
 381    * @param initialColumnSeparation
 382    *          The initial whitespace separation between the left side of the
 383    *          file and first character.
 384    * @param columnSeparation
 385    *          The separation between subsequent data entries.
 386    * @param data
 387    *          The list fo data to be added to the String.
 388    * @return
 389    */
 390   String addData(int initialColumnSeparation,
 391           int columnSeparation, List<String> data)
 392   {
 393     String line = EMPTY;
 394     int index = 0;
 395     for (String value : data)
 396     {
 397       if (index == 0)
 398       {
 399         line += String.format("%" + initialColumnSeparation + "s", value);
 400       }
 401       else
 402       {
 403         line += String.format("%" + columnSeparation + "s", value);
 404       }
 405       index++;
 406     }
 407     return line;
 408   }
 409
 410   /**
 411    * Converts list of characters into a list of Strings.
 412    *
 413    * @param list
 414    * @return Returns the list of Strings.
 415    */
 416   List<String> charListToStringList(List<Character> list)
 417   {
 418     List<String> strList = new ArrayList<>();
 419     for (char value : list)
 420     {
 421       String strValue = Character.toString(value);
 422       strList.add(strValue);
 423     }
 424     return strList;
 425   }
 426
 427   /**
 428    * Converts a list of doubles into a list of Strings, rounded to the nearest
 429    * 5th decimal place.
 430    *
 431    * @param list
 432    * @param noOfDecimals
 433    * @return
 434    */
 435   List<String> doubleListToStringList(List<Double> list)
 436   {
 437     List<String> strList = new ArrayList<>();
 438     for (double value : list)
 439     {
 440       String strValue;
 441       if (value > 0)
 442       {
 443         strValue = String.format("%.5f", value);
 444
 445       }
 446       else if (value == -0.00000d)
 447       {
 448         strValue = "0.00000";
 449       }
 450       else
 451       {
 452         strValue = "*";
 453       }
 454
 455       strList.add(strValue);
 456     }
 457     return strList;
 458   }
 459
 460   /**
 461    * Converts a primitive array of Strings to a list of Strings.
 462    *
 463    * @param array
 464    * @return
 465    */
 466   List<String> stringArrayToStringList(String[] array)
 467   {
 468     List<String> list = new ArrayList<>();
 469     for (String value : array)
 470     {
 471       list.add(value);
 472     }
 473
 474     return list;
 475   }
 476
 477   /**
 478    * Returns a string containing the model data.
 479    */
 480   String getModelAsString()
 481   {
 482     StringBuilder output = new StringBuilder();
 483     String symbolLine = "HMM";
 484     List<Character> charSymbols = hmm.getSymbols();
 485     List<String> strSymbols;
 486     strSymbols = charListToStringList(charSymbols);
 487     symbolLine += addData(11, 9, strSymbols);
 488     output.append(symbolLine);
 489     output.append(NL).append(TRANSITIONTYPELINE);
 490
 491     int length = hmm.getLength();
 492
 493     for (int node = 0; node <= length; node++)
 494     {
 495       String matchLine;
 496       if (node == 0)
 497       {
 498         matchLine = String.format("%7s", "COMPO");
 499       }
 500       else
 501       {
 502         matchLine = String.format("%7s", node);
 503       }
 504
 505       List<String> strMatches;
 506       List<Double> doubleMatches;
 507       doubleMatches = convertListToLogSpace(
 508               hmm.getNode(node).getMatchEmissions());
 509       strMatches = doubleListToStringList(doubleMatches);
 510       matchLine += addData(10, 9, strMatches);
 511
 512
 513       if (node != 0)
 514       {
 515         matchLine += SPACE + (hmm.getNodeAlignmentColumn(node) + 1);
 516         matchLine += SPACE + hmm.getConsensusResidue(node);
 517         matchLine += SPACE + hmm.getReferenceAnnotation(node);
 518         if (hmm.getFileHeader().contains("HMMER3/f"))
 519         {
 520           matchLine += SPACE + hmm.getMaskedValue(node);
 521           matchLine += SPACE + hmm.getConsensusStructure(node);
 522         }
 523
 524       }
 525
 526       output.append(NL).append(matchLine);
 527
 528       String insertLine = EMPTY;
 529       List<String> strInserts;
 530       List<Double> doubleInserts;
 531       doubleInserts = convertListToLogSpace(
 532               hmm.getNode(node).getInsertEmissions());
 533       strInserts = doubleListToStringList(doubleInserts);
 534       insertLine += addData(17, 9, strInserts);
 535
 536       output.append(NL).append(insertLine);
 537
 538       String transitionLine = EMPTY;
 539       List<String> strTransitions;
 540       List<Double> doubleTransitions;
 541       doubleTransitions = convertListToLogSpace(
 542               hmm.getNode(node).getStateTransitions());
 543       strTransitions = doubleListToStringList(doubleTransitions);
 544       transitionLine += addData(17, 9, strTransitions);
 545
 546       output.append(NL).append(transitionLine);
 547     }
 548     return output.toString();
 549   }
 550
 551   /**
 552    * Returns a String containing the HMM file properties
 553    */
 554   String getFilePropertiesAsString()
 555   {
 556     StringBuffer output = new StringBuffer();
 557     String line;
 558
 559     output.append(hmm.getFileHeader());
 560
 561     line = String.format("%-5s %1s", "NAME", hmm.getName());
 562     output.append(NL + line);
 563
 564     if (hmm.getAccessionNumber() != null)
 565     {
 566     line = String.format("%-5s %1s", "ACC", hmm.getAccessionNumber());
 567       output.append(NL + line);
 568     }
 569
 570     if (hmm.getDescription() != null)
 571     {
 572     line = String.format("%-5s %1s", "DESC", hmm.getDescription());
 573       output.append(NL + line);
 574     }
 575     line = String.format("%-5s %1s", "LENG", hmm.getLength());
 576     output.append(NL + line);
 577
 578     if (hmm.getMaxInstanceLength() != null)
 579     {
 580     line = String.format("%-5s %1s", "MAXL", hmm.getMaxInstanceLength());
 581       output.append(NL + line);
 582     }
 583     line = String.format("%-5s %1s", "ALPH", hmm.getAlphabetType());
 584     output.append(NL + line);
 585
 586     boolean status;
 587     String statusStr;
 588
 589     status = hmm.referenceAnnotationIsActive();
 590     statusStr = HiddenMarkovModel.findStringFromBoolean(status);
 591     line = String.format("%-5s %1s", "RF",
 592             statusStr);
 593     output.append(NL + line);
 594
 595     status = hmm.maskValueIsActive();
 596     statusStr = HiddenMarkovModel.findStringFromBoolean(status);
 597     line = String.format("%-5s %1s", "MM",
 598             statusStr);
 599     output.append(NL + line);
 600
 601     status = hmm.consensusResidueIsActive();
 602     statusStr = HiddenMarkovModel.findStringFromBoolean(status);
 603     line = String.format("%-5s %1s", "CONS",
 604             statusStr);
 605     output.append(NL + line);
 606
 607     status = hmm.consensusStructureIsActive();
 608     statusStr = HiddenMarkovModel.findStringFromBoolean(status);
 609     line = String.format("%-5s %1s", "CS",
 610             statusStr);
 611     output.append(NL + line);
 612
 613     status = hmm.mapIsActive();
 614     statusStr = HiddenMarkovModel.findStringFromBoolean(status);
 615     line = String.format("%-5s %1s", "MAP",
 616             statusStr);
 617     output.append(NL + line);
 618
 619
 620     if (hmm.getDate() != null)
 621     {
 622     line = String.format("%-5s %1s", "DATE", hmm.getDate());
 623       output.append(NL + line);
 624     }
 625     if (hmm.getNumberOfSequences() != null)
 626     {
 627     line = String.format("%-5s %1s", "NSEQ", hmm.getNumberOfSequences());
 628       output.append(NL + line);
 629     }
 630     if (hmm.getEffectiveNumberOfSequences() != null)
 631     {
 632     line = String.format("%-5s %1s", "EFFN",
 633             hmm.getEffectiveNumberOfSequences());
 634       output.append(NL + line);
 635     }
 636     if (hmm.getCheckSum() != null)
 637     {
 638     line = String.format("%-5s %1s", "CKSUM", hmm.getCheckSum());
 639       output.append(NL + line);
 640     }
 641     if (hmm.getGatheringThreshold() != null)
 642     {
 643     line = String.format("%-5s %1s", "GA", hmm.getGatheringThreshold());
 644       output.append(NL + line);
 645     }
 646
 647     if (hmm.getTrustedCutoff() != null)
 648     {
 649     line = String.format("%-5s %1s", "TC", hmm.getTrustedCutoff());
 650       output.append(NL + line);
 651     }
 652     if (hmm.getNoiseCutoff() != null)
 653     {
 654     line = String.format("%-5s %1s", "NC", hmm.getNoiseCutoff());
 655       output.append(NL + line);
 656     }
 657     if (hmm.getMSV() != null)
 658     {
 659       line = String.format("%-19s %18s", "STATS LOCAL MSV", hmm.getMSV());
 660       output.append(NL + line);
 661
 662       line = String.format("%-19s %18s", "STATS LOCAL VITERBI",
 663               hmm.getViterbi());
 664       output.append(NL + line);
 665
 666       line = String.format("%-19s %18s", "STATS LOCAL FORWARD",
 667               hmm.getForward());
 668       output.append(NL + line);
 669     }
 670     return output.toString();
 671   }
 672
 673
 674   /**
 675    * Returns the char value of a single lettered String.
 676    *
 677    * @param string
 678    * @return
 679    */
 680   char charValue(String string)
 681   {
 682     char character;
 683     character = string.charAt(0);
 684     return character;
 685
 686   }
 687
 688   @Override
 689   public String print(SequenceI[] seqs, boolean jvsuffix)
 690   {
 691     if (seqs[0].getHMM() != null)
 692     {
 693       hmm = seqs[0].getHMM();
 694     }
 695     return print();
 696   }
 697
 698   /**
 699    * Prints the .hmm file to a String.
 700    *
 701    * @return
 702    */
 703   public String print()
 704   {
 705     StringBuffer output = new StringBuffer();
 706     output.append(getFilePropertiesAsString());
 707     output.append(NL);
 708     output.append(getModelAsString());
 709     output.append(NL + "//");
 710     return output.toString();
 711   }
 712
 713   /**
 714    * Converts the probabilities contained in a list into log space.
 715    *
 716    * @param list
 717    */
 718   List<Double> convertListToLogSpace(List<Double> list)
 719   {
 720
 721     List<Double> convertedList = new ArrayList<>();
 722     for (int i = 0; i < list.size(); i++)
 723     {
 724       double prob = list.get(i);
 725       double logProb = -1 * Math.log(prob);
 726
 727       convertedList.add(logProb);
 728     }
 729     return convertedList;
 730
 731
 732   }
 733
 734   /**
 735    * Returns the HMM sequence produced by reading a .hmm file.
 736    */
 737   @Override
 738   public SequenceI[] getSeqsAsArray()
 739   {
 740     SequenceI hmmSeq = hmm.initHMMSequence();
 741     SequenceI[] seq = new SequenceI[1];
 742     seq[0] = hmmSeq;
 743     return seq;
 744
 745   }
 746
 747   /**
 748    * Fills symbol array and adds each symbol to an index lookup
 749    *
 750    * @param parser
 751    *          The scanner scanning the symbol line in the file.
 752    */
 753   public void fillSymbols(Scanner parser)
 754   {
 755     int i = 0;
 756     while (parser.hasNext())
 757     {
 758       String strSymbol = parser.next();
 759       char[] symbol = strSymbol.toCharArray();
 760       hmm.getSymbols().add(symbol[0]);
 761       hmm.setSymbolIndex(symbol[0], i);
 762       i++;
 763     }
 764   }
 765
 766   @Override
 767   public void setNewlineString(String newLine)
 768   {
 769     NL = newLine;
 770   }
 771
 772   @Override
 773   public void setExportSettings(AlignExportSettingI exportSettings)
 774   {
 775
 776   }
 777
 778   @Override
 779   public void configureForView(AlignmentViewPanel viewpanel)
 780   {
 781
 782   }
 783
 784   @Override
 785   public boolean hasWarningMessage()
 786   {
 787     return false;
 788   }
 789
 790   @Override
 791   public String getWarningMessage()
 792   {
 793     return "warning message";
 794   }
 795
 796 }
 797