src/jalview/io/StockholmFile.java

   1 /*\r
   2  * Jalview - A Sequence Alignment Editor and Viewer (Version 2.7)\r
   3  * Copyright (C) 2011 J Procter, AM Waterhouse, G Barton, M Clamp, S Searle\r
   4  * \r
   5  * This file is part of Jalview.\r
   6  * \r
   7  * Jalview is free software: you can redistribute it and/or\r
   8  * modify it under the terms of the GNU General Public License \r
   9  * as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.\r
  10  * \r
  11  * Jalview is distributed in the hope that it will be useful, but \r
  12  * WITHOUT ANY WARRANTY; without even the implied warranty \r
  13  * of MERCHANTABILITY or FITNESS FOR A PARTICULAR \r
  14  * PURPOSE.  See the GNU General Public License for more details.\r
  15  * \r
  16  * You should have received a copy of the GNU General Public License along with Jalview.  If not, see <http://www.gnu.org/licenses/>.\r
  17  */\r
  18 /*\r
  19  * This extension was written by Benjamin Schuster-Boeckler at sanger.ac.uk\r
  20  */\r
  21 package jalview.io;\r
  22 \r
  23 import java.io.*;\r
  24 import java.util.*;\r
  25 \r
  26 import com.stevesoft.pat.*;\r
  27 import jalview.datamodel.*;\r
  28 \r
  29 // import org.apache.log4j.*;\r
  30 \r
  31 /**\r
  32  * This class is supposed to parse a Stockholm format file into Jalview There\r
  33  * are TODOs in this class: we do not know what the database source and version\r
  34  * is for the file when parsing the #GS= AC tag which associates accessions with\r
  35  * sequences. Database references are also not parsed correctly: a separate\r
  36  * reference string parser must be added to parse the database reference form\r
  37  * into Jalview's local representation.\r
  38  * \r
  39  * @author bsb at sanger.ac.uk\r
  40  * @version 0.3 + jalview mods\r
  41  * \r
  42  */\r
  43 public class StockholmFile extends AlignFile\r
  44 {\r
  45   // static Logger logger = Logger.getLogger("jalview.io.StockholmFile");\r
  46 \r
  47   public StockholmFile()\r
  48   {\r
  49   }\r
  50 \r
  51   public StockholmFile(String inFile, String type) throws IOException\r
  52   {\r
  53     super(inFile, type);\r
  54   }\r
  55 \r
  56   public StockholmFile(FileParse source) throws IOException\r
  57   {\r
  58     super(source);\r
  59   }\r
  60 \r
  61   public void initData()\r
  62   {\r
  63     super.initData();\r
  64   }\r
  65 \r
  66   /**\r
  67    * Parse a file in Stockholm format into Jalview's data model. The file has to\r
  68    * be passed at construction time\r
  69    * \r
  70    * @throws IOException\r
  71    *           If there is an error with the input file\r
  72    */\r
  73   public void parse() throws IOException\r
  74   {\r
  75     StringBuffer treeString = new StringBuffer();\r
  76     String treeName = null;\r
  77     // --------------- Variable Definitions -------------------\r
  78     String line;\r
  79     String version;\r
  80     // String id;\r
  81     Hashtable seqAnn = new Hashtable(); // Sequence related annotations\r
  82     Hashtable seqs = new Hashtable();\r
  83     Regex p, r, rend, s, x;\r
  84 \r
  85     // ------------------ Parsing File ----------------------\r
  86     // First, we have to check that this file has STOCKHOLM format, i.e. the\r
  87     // first line must match\r
  88     r = new Regex("# STOCKHOLM ([\\d\\.]+)");\r
  89     if (!r.search(nextLine()))\r
  90     {\r
  91       throw new IOException(\r
  92               "This file is not in valid STOCKHOLM format: First line does not contain '# STOCKHOLM'");\r
  93     }\r
  94     else\r
  95     {\r
  96       version = r.stringMatched(1);\r
  97       // logger.debug("Stockholm version: " + version);\r
  98     }\r
  99 \r
 100     // We define some Regexes here that will be used regularily later\r
 101     rend = new Regex("^\\s*\\/\\/"); // Find the end of an alignment\r
 102     p = new Regex("(\\S+)\\/(\\d+)\\-(\\d+)"); // split sequence id in\r
 103     // id/from/to\r
 104     s = new Regex("(\\S+)\\s+(\\S*)\\s+(.*)"); // Parses annotation subtype\r
 105     r = new Regex("#=(G[FSRC]?)\\s+(.*)"); // Finds any annotation line\r
 106     x = new Regex("(\\S+)\\s+(\\S+)"); // split id from sequence\r
 107 \r
 108     rend.optimize();\r
 109     p.optimize();\r
 110     s.optimize();\r
 111     r.optimize();\r
 112     x.optimize();\r
 113 \r
 114     while ((line = nextLine()) != null)\r
 115     {\r
 116       if (line.length() == 0)\r
 117       {\r
 118         continue;\r
 119       }\r
 120       if (rend.search(line))\r
 121       {\r
 122         // End of the alignment, pass stuff back\r
 123 \r
 124         this.noSeqs = seqs.size();\r
 125         // logger.debug("Number of sequences: " + this.noSeqs);\r
 126         Enumeration accs = seqs.keys();\r
 127         while (accs.hasMoreElements())\r
 128         {\r
 129           String acc = (String) accs.nextElement();\r
 130           // logger.debug("Processing sequence " + acc);\r
 131           String seq = (String) seqs.remove(acc);\r
 132           if (maxLength < seq.length())\r
 133           {\r
 134             maxLength = seq.length();\r
 135           }\r
 136           int start = 1;\r
 137           int end = -1;\r
 138           String sid = acc;\r
 139           // Retrieve hash of annotations for this accession\r
 140           Hashtable accAnnotations = null;\r
 141 \r
 142           if (seqAnn != null && seqAnn.containsKey(acc))\r
 143           {\r
 144             accAnnotations = (Hashtable) seqAnn.remove(acc);\r
 145           }\r
 146 \r
 147           // Split accession in id and from/to\r
 148           if (p.search(acc))\r
 149           {\r
 150             sid = p.stringMatched(1);\r
 151             start = Integer.parseInt(p.stringMatched(2));\r
 152             end = Integer.parseInt(p.stringMatched(3));\r
 153           }\r
 154           // logger.debug(sid + ", " + start + ", " + end);\r
 155 \r
 156           Sequence seqO = new Sequence(sid, seq, start, end);\r
 157           // Add Description (if any)\r
 158           if (accAnnotations != null && accAnnotations.containsKey("DE"))\r
 159           {\r
 160             String desc = (String) accAnnotations.get("DE");\r
 161             seqO.setDescription((desc == null) ? "" : desc);\r
 162           }\r
 163           // Add DB References (if any)\r
 164           if (accAnnotations != null && accAnnotations.containsKey("DR"))\r
 165           {\r
 166             String dbr = (String) accAnnotations.get("DR");\r
 167             if (dbr != null && dbr.indexOf(";") > -1)\r
 168             {\r
 169               String src = dbr.substring(0, dbr.indexOf(";"));\r
 170               String acn = dbr.substring(dbr.indexOf(";") + 1);\r
 171               jalview.util.DBRefUtils.parseToDbRef(seqO, src, "0", acn);\r
 172               // seqO.addDBRef(dbref);\r
 173             }\r
 174           }\r
 175           Hashtable features = null;\r
 176           // We need to adjust the positions of all features to account for gaps\r
 177           try\r
 178           {\r
 179             features = (Hashtable) accAnnotations.remove("features");\r
 180           } catch (java.lang.NullPointerException e)\r
 181           {\r
 182             // loggerwarn("Getting Features for " + acc + ": " +\r
 183             // e.getMessage());\r
 184             // continue;\r
 185           }\r
 186           // if we have features\r
 187           if (features != null)\r
 188           {\r
 189             int posmap[] = seqO.findPositionMap();\r
 190             Enumeration i = features.keys();\r
 191             while (i.hasMoreElements())\r
 192             {\r
 193               // TODO: parse out secondary structure annotation as annotation\r
 194               // row\r
 195               // TODO: parse out scores as annotation row\r
 196               // TODO: map coding region to core jalview feature types\r
 197               String type = i.nextElement().toString();\r
 198               Hashtable content = (Hashtable) features.remove(type);\r
 199               Enumeration j = content.keys();\r
 200               while (j.hasMoreElements())\r
 201               {\r
 202                 String desc = j.nextElement().toString();\r
 203                 String ns = content.get(desc).toString();\r
 204                 char[] byChar = ns.toCharArray();\r
 205                 for (int k = 0; k < byChar.length; k++)\r
 206                 {\r
 207                   char c = byChar[k];\r
 208                   if (!(c == ' ' || c == '_' || c == '-' || c == '.')) // PFAM\r
 209                   // uses\r
 210                   // '.'\r
 211                   // for\r
 212                   // feature\r
 213                   // background\r
 214                   {\r
 215                     int new_pos = posmap[k]; // look up nearest seqeunce\r
 216                     // position to this column\r
 217                     SequenceFeature feat = new SequenceFeature(type, desc,\r
 218                             new_pos, new_pos, 0f, null);\r
 219 \r
 220                     seqO.addSequenceFeature(feat);\r
 221                   }\r
 222                 }\r
 223               }\r
 224 \r
 225             }\r
 226 \r
 227           }\r
 228           // garbage collect\r
 229 \r
 230           // logger.debug("Adding seq " + acc + " from " + start + " to " + end\r
 231           // + ": " + seq);\r
 232           this.seqs.addElement(seqO);\r
 233         }\r
 234         return; // finished parsing this segment of source\r
 235       }\r
 236       else if (!r.search(line))\r
 237       {\r
 238         // System.err.println("Found sequence line: " + line);\r
 239 \r
 240         // Split sequence in sequence and accession parts\r
 241         if (!x.search(line))\r
 242         {\r
 243           // logger.error("Could not parse sequence line: " + line);\r
 244           throw new IOException("Could not parse sequence line: " + line);\r
 245         }\r
 246         String ns = (String) seqs.get(x.stringMatched(1));\r
 247         if (ns == null)\r
 248         {\r
 249           ns = "";\r
 250         }\r
 251         ns += x.stringMatched(2);\r
 252 \r
 253         seqs.put(x.stringMatched(1), ns);\r
 254       }\r
 255       else\r
 256       {\r
 257         String annType = r.stringMatched(1);\r
 258         String annContent = r.stringMatched(2);\r
 259 \r
 260         // System.err.println("type:" + annType + " content: " + annContent);\r
 261 \r
 262         if (annType.equals("GF"))\r
 263         {\r
 264           /*\r
 265            * Generic per-File annotation, free text Magic features: #=GF NH\r
 266            * <tree in New Hampshire eXtended format> #=GF TN <Unique identifier\r
 267            * for the next tree> Pfam descriptions: 7. DESCRIPTION OF FIELDS\r
 268            * \r
 269            * Compulsory fields: ------------------\r
 270            * \r
 271            * AC Accession number: Accession number in form PFxxxxx.version or\r
 272            * PBxxxxxx. ID Identification: One word name for family. DE\r
 273            * Definition: Short description of family. AU Author: Authors of the\r
 274            * entry. SE Source of seed: The source suggesting the seed members\r
 275            * belong to one family. GA Gathering method: Search threshold to\r
 276            * build the full alignment. TC Trusted Cutoff: Lowest sequence score\r
 277            * and domain score of match in the full alignment. NC Noise Cutoff:\r
 278            * Highest sequence score and domain score of match not in full\r
 279            * alignment. TP Type: Type of family -- presently Family, Domain,\r
 280            * Motif or Repeat. SQ Sequence: Number of sequences in alignment. AM\r
 281            * Alignment Method The order ls and fs hits are aligned to the model\r
 282            * to build the full align. // End of alignment.\r
 283            * \r
 284            * Optional fields: ----------------\r
 285            * \r
 286            * DC Database Comment: Comment about database reference. DR Database\r
 287            * Reference: Reference to external database. RC Reference Comment:\r
 288            * Comment about literature reference. RN Reference Number: Reference\r
 289            * Number. RM Reference Medline: Eight digit medline UI number. RT\r
 290            * Reference Title: Reference Title. RA Reference Author: Reference\r
 291            * Author RL Reference Location: Journal location. PI Previous\r
 292            * identifier: Record of all previous ID lines. KW Keywords: Keywords.\r
 293            * CC Comment: Comments. NE Pfam accession: Indicates a nested domain.\r
 294            * NL Location: Location of nested domains - sequence ID, start and\r
 295            * end of insert.\r
 296            * \r
 297            * Obsolete fields: ----------- AL Alignment method of seed: The\r
 298            * method used to align the seed members.\r
 299            */\r
 300           // Let's save the annotations, maybe we'll be able to do something\r
 301           // with them later...\r
 302           Regex an = new Regex("(\\w+)\\s*(.*)");\r
 303           if (an.search(annContent))\r
 304           {\r
 305             if (an.stringMatched(1).equals("NH"))\r
 306             {\r
 307               treeString.append(an.stringMatched(2));\r
 308             }\r
 309             else if (an.stringMatched(1).equals("TN"))\r
 310             {\r
 311               if (treeString.length() > 0)\r
 312               {\r
 313                 if (treeName == null)\r
 314                 {\r
 315                   treeName = "Tree " + (getTreeCount() + 1);\r
 316                 }\r
 317                 addNewickTree(treeName, treeString.toString());\r
 318               }\r
 319               treeName = an.stringMatched(2);\r
 320               treeString = new StringBuffer();\r
 321             }\r
 322             setAlignmentProperty(an.stringMatched(1), an.stringMatched(2));\r
 323           }\r
 324         }\r
 325         else if (annType.equals("GS"))\r
 326         {\r
 327           // Generic per-Sequence annotation, free text\r
 328           /*\r
 329            * Pfam uses these features: Feature Description ---------------------\r
 330            * ----------- AC <accession> ACcession number DE <freetext>\r
 331            * DEscription DR <db>; <accession>; Database Reference OS <organism>\r
 332            * OrganiSm (species) OC <clade> Organism Classification (clade, etc.)\r
 333            * LO <look> Look (Color, etc.)\r
 334            */\r
 335           if (s.search(annContent))\r
 336           {\r
 337             String acc = s.stringMatched(1);\r
 338             String type = s.stringMatched(2);\r
 339             String content = s.stringMatched(3);\r
 340             // TODO: store DR in a vector.\r
 341             // TODO: store AC according to generic file db annotation.\r
 342             Hashtable ann;\r
 343             if (seqAnn.containsKey(acc))\r
 344             {\r
 345               ann = (Hashtable) seqAnn.get(acc);\r
 346             }\r
 347             else\r
 348             {\r
 349               ann = new Hashtable();\r
 350             }\r
 351             ann.put(type, content);\r
 352             seqAnn.put(acc, ann);\r
 353           }\r
 354           else\r
 355           {\r
 356             throw new IOException("Error parsing " + line);\r
 357           }\r
 358         }\r
 359         else if (annType.equals("GC"))\r
 360         {\r
 361           // Generic per-Column annotation, exactly 1 char per column\r
 362           // always need a label.\r
 363           if (x.search(annContent))\r
 364           {\r
 365             // parse out and create alignment annotation directly.\r
 366             parseAnnotationRow(annotations, x.stringMatched(1),\r
 367                     x.stringMatched(2));\r
 368           }\r
 369         }\r
 370         else if (annType.equals("GR"))\r
 371         {\r
 372           // Generic per-Sequence AND per-Column markup, exactly 1 char per\r
 373           // column\r
 374           /*\r
 375            * Feature Description Markup letters ------- -----------\r
 376            * -------------- SS Secondary Structure [HGIEBTSCX] SA Surface\r
 377            * Accessibility [0-9X] (0=0%-10%; ...; 9=90%-100%) TM TransMembrane\r
 378            * [Mio] PP Posterior Probability [0-9*] (0=0.00-0.05; 1=0.05-0.15;\r
 379            * *=0.95-1.00) LI LIgand binding [*] AS Active Site [*] IN INtron (in\r
 380            * or after) [0-2]\r
 381            */\r
 382           if (s.search(annContent))\r
 383           {\r
 384             String acc = s.stringMatched(1);\r
 385             String type = s.stringMatched(2);\r
 386             String seq = new String(s.stringMatched(3));\r
 387             String description = null;\r
 388             // Check for additional information about the current annotation\r
 389             // We use a simple string tokenizer here for speed\r
 390             StringTokenizer sep = new StringTokenizer(seq, " \t");\r
 391             description = sep.nextToken();\r
 392             if (sep.hasMoreTokens())\r
 393             {\r
 394               seq = sep.nextToken();\r
 395             }\r
 396             else\r
 397             {\r
 398               seq = description;\r
 399               description = new String();\r
 400             }\r
 401             // sequence id with from-to fields\r
 402 \r
 403             Hashtable ann;\r
 404             // Get an object with all the annotations for this sequence\r
 405             if (seqAnn.containsKey(acc))\r
 406             {\r
 407               // logger.debug("Found annotations for " + acc);\r
 408               ann = (Hashtable) seqAnn.get(acc);\r
 409             }\r
 410             else\r
 411             {\r
 412               // logger.debug("Creating new annotations holder for " + acc);\r
 413               ann = new Hashtable();\r
 414               seqAnn.put(acc, ann);\r
 415             }\r
 416 \r
 417             Hashtable features;\r
 418             // Get an object with all the content for an annotation\r
 419             if (ann.containsKey("features"))\r
 420             {\r
 421               // logger.debug("Found features for " + acc);\r
 422               features = (Hashtable) ann.get("features");\r
 423             }\r
 424             else\r
 425             {\r
 426               // logger.debug("Creating new features holder for " + acc);\r
 427               features = new Hashtable();\r
 428               ann.put("features", features);\r
 429             }\r
 430 \r
 431             Hashtable content;\r
 432             if (features.containsKey(this.id2type(type)))\r
 433             {\r
 434               // logger.debug("Found content for " + this.id2type(type));\r
 435               content = (Hashtable) features.get(this.id2type(type));\r
 436             }\r
 437             else\r
 438             {\r
 439               // logger.debug("Creating new content holder for " +\r
 440               // this.id2type(type));\r
 441               content = new Hashtable();\r
 442               features.put(this.id2type(type), content);\r
 443             }\r
 444             String ns = (String) content.get(description);\r
 445             if (ns == null)\r
 446             {\r
 447               ns = "";\r
 448             }\r
 449             ns += seq;\r
 450             content.put(description, ns);\r
 451           }\r
 452           else\r
 453           {\r
 454             System.err\r
 455                     .println("Warning - couldn't parse sequence annotation row line:\n"\r
 456                             + line);\r
 457             // throw new IOException("Error parsing " + line);\r
 458           }\r
 459         }\r
 460         else\r
 461         {\r
 462           throw new IOException("Unknown annotation detected: " + annType\r
 463                   + " " + annContent);\r
 464         }\r
 465       }\r
 466     }\r
 467     if (treeString.length() > 0)\r
 468     {\r
 469       if (treeName == null)\r
 470       {\r
 471         treeName = "Tree " + (1 + getTreeCount());\r
 472       }\r
 473       addNewickTree(treeName, treeString.toString());\r
 474     }\r
 475   }\r
 476 \r
 477   private AlignmentAnnotation parseAnnotationRow(Vector annotation,\r
 478           String label, String annots)\r
 479   {\r
 480     String type = (label.indexOf("_cons") == label.length() - 5) ? label\r
 481             .substring(0, label.length() - 5) : label;\r
 482     boolean ss = false;\r
 483     type = id2type(type);\r
 484     if (type.equals("secondary structure"))\r
 485     {\r
 486       ss = true;\r
 487     }\r
 488     // decide on secondary structure or not.\r
 489     Annotation[] els = new Annotation[annots.length()];\r
 490     for (int i = 0; i < annots.length(); i++)\r
 491     {\r
 492       String pos = annots.substring(i, i + 1);\r
 493       Annotation ann;\r
 494       ann = new Annotation(pos, "", ' ', 0f); // 0f is 'valid' null - will not\r
 495       // be written out\r
 496       if (ss)\r
 497       {\r
 498         ann.secondaryStructure = jalview.schemes.ResidueProperties\r
 499                 .getDssp3state(pos).charAt(0);\r
 500         if (ann.secondaryStructure == pos.charAt(0) || pos.charAt(0) == 'C')\r
 501         {\r
 502           ann.displayCharacter = ""; // null; // " ";\r
 503         }\r
 504         else\r
 505         {\r
 506           ann.displayCharacter = " " + ann.displayCharacter;\r
 507         }\r
 508       }\r
 509 \r
 510       els[i] = ann;\r
 511     }\r
 512     AlignmentAnnotation annot = null;\r
 513     Enumeration e = annotation.elements();\r
 514     while (e.hasMoreElements())\r
 515     {\r
 516       annot = (AlignmentAnnotation) e.nextElement();\r
 517       if (annot.label.equals(type))\r
 518         break;\r
 519       annot = null;\r
 520     }\r
 521     if (annot == null)\r
 522     {\r
 523       annot = new AlignmentAnnotation(type, type, els);\r
 524       annotation.addElement(annot);\r
 525     }\r
 526     else\r
 527     {\r
 528       Annotation[] anns = new Annotation[annot.annotations.length\r
 529               + els.length];\r
 530       System.arraycopy(annot.annotations, 0, anns, 0,\r
 531               annot.annotations.length);\r
 532       System.arraycopy(els, 0, anns, annot.annotations.length, els.length);\r
 533       annot.annotations = anns;\r
 534     }\r
 535     return annot;\r
 536   }\r
 537 \r
 538   public static String print(SequenceI[] s)\r
 539   {\r
 540     return "not yet implemented";\r
 541   }\r
 542 \r
 543   public String print()\r
 544   {\r
 545     return print(getSeqsAsArray());\r
 546   }\r
 547 \r
 548   private static Hashtable typeIds = null;\r
 549   static\r
 550   {\r
 551     if (typeIds == null)\r
 552     {\r
 553       typeIds = new Hashtable();\r
 554       typeIds.put("SS", "secondary structure");\r
 555       typeIds.put("SA", "surface accessibility");\r
 556       typeIds.put("TM", "transmembrane");\r
 557       typeIds.put("PP", "posterior probability");\r
 558       typeIds.put("LI", "ligand binding");\r
 559       typeIds.put("AS", "active site");\r
 560       typeIds.put("IN", "intron");\r
 561       typeIds.put("IR", "interacting residue");\r
 562       typeIds.put("AC", "accession");\r
 563       typeIds.put("OS", "organism");\r
 564       typeIds.put("CL", "class");\r
 565       typeIds.put("DE", "description");\r
 566       typeIds.put("DR", "reference");\r
 567       typeIds.put("LO", "look");\r
 568       typeIds.put("RF", "reference positions");\r
 569 \r
 570     }\r
 571   }\r
 572 \r
 573   private String id2type(String id)\r
 574   {\r
 575     if (typeIds.containsKey(id))\r
 576     {\r
 577       return (String) typeIds.get(id);\r
 578     }\r
 579     System.err.println("Warning : Unknown Stockholm annotation type code "\r
 580             + id);\r
 581     return id;\r
 582   }\r
 583 }\r