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