src/jalview/io/HMMFile.java

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