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