src/jalview/io/StockholmFile.java

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