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