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