src/jalview/ext/so/SequenceOntology.java

   1 /*
   2  * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$)
   3  * Copyright (C) $$Year-Rel$$ The Jalview Authors
   4  *
   5  * This file is part of Jalview.
   6  *
   7  * Jalview is free software: you can redistribute it and/or
   8  * modify it under the terms of the GNU General Public License
   9  * as published by the Free Software Foundation, either version 3
  10  * of the License, or (at your option) any later version.
  11  *
  12  * Jalview is distributed in the hope that it will be useful, but
  13  * WITHOUT ANY WARRANTY; without even the implied warranty
  14  * of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  15  * PURPOSE.  See the GNU General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU General Public License
  18  * along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
  19  * The Jalview Authors are detailed in the 'AUTHORS' file.
  20  */
  21 package jalview.ext.so;
  22
  23 import jalview.datamodel.ontology.OntologyBase;
  24 import jalview.io.gff.SequenceOntologyI;
  25
  26 import java.io.BufferedInputStream;
  27 import java.io.BufferedReader;
  28 import java.io.IOException;
  29 import java.io.InputStream;
  30 import java.io.InputStreamReader;
  31 import java.text.ParseException;
  32 import java.util.ArrayList;
  33 import java.util.Collections;
  34 import java.util.HashMap;
  35 import java.util.HashSet;
  36 import java.util.List;
  37 import java.util.Map;
  38 import java.util.NoSuchElementException;
  39 import java.util.Set;
  40 import java.util.zip.ZipEntry;
  41 import java.util.zip.ZipInputStream;
  42
  43 import org.biojava.nbio.ontology.Ontology;
  44 import org.biojava.nbio.ontology.Synonym;
  45 import org.biojava.nbio.ontology.Term;
  46 import org.biojava.nbio.ontology.Term.Impl;
  47 import org.biojava.nbio.ontology.Triple;
  48 import org.biojava.nbio.ontology.io.OboParser;
  49 import org.biojava.nbio.ontology.utils.Annotation;
  50
  51 /**
  52  * A wrapper class that parses the Sequence Ontology and exposes useful access
  53  * methods. This version uses the BioJava parser.
  54  */
  55 public class SequenceOntology extends OntologyBase
  56         implements SequenceOntologyI
  57 {
  58   /*
  59    * the parsed Ontology data as modelled by BioJava
  60    */
  61   private Ontology ontology;
  62
  63   /*
  64    * the ontology term for the isA relationship
  65    */
  66   private Term isA;
  67
  68   /*
  69    * lookup of terms by user readable name (NB not guaranteed unique)
  70    */
  71   private Map<String, Term> aliases;
  72
  73   /*
  74    * Map where key is a Term and value is a (possibly empty) list of
  75    * all Terms to which the key has an 'isA' relationship, either
  76    * directly or indirectly (A isA B isA C)
  77    */
  78   private Map<Term, List<Term>> termIsA;
  79
  80   private List<String> termsFound;
  81
  82   private List<String> termsNotFound;
  83
  84   /**
  85    * Package private constructor to enforce use of singleton. Parses and caches
  86    * the SO OBO data file.
  87    */
  88   public SequenceOntology()
  89   {
  90     termsFound = new ArrayList<>();
  91     termsNotFound = new ArrayList<>();
  92     aliases = new HashMap<>();
  93     termIsA = new HashMap<>();
  94
  95     loadOntologyZipFile("so-xp-simple.obo");
  96   }
  97
  98   /**
  99    * Loads the given ontology file from a zip file with ".zip" appended
 100    *
 101    * @param ontologyFile
 102    */
 103   protected void loadOntologyZipFile(String ontologyFile)
 104   {
 105     long now = System.currentTimeMillis();
 106     ZipInputStream zipStream = null;
 107     try
 108     {
 109       String zipFile = ontologyFile + ".zip";
 110       InputStream inStream = this.getClass()
 111               .getResourceAsStream("/" + zipFile);
 112       zipStream = new ZipInputStream(new BufferedInputStream(inStream));
 113       ZipEntry entry;
 114       while ((entry = zipStream.getNextEntry()) != null)
 115       {
 116         if (entry.getName().equals(ontologyFile))
 117         {
 118           loadOboFile(zipStream);
 119         }
 120       }
 121       long elapsed = System.currentTimeMillis() - now;
 122       System.out.println("Loaded Sequence Ontology from " + zipFile + " ("
 123               + elapsed + "ms)");
 124     } catch (Exception e)
 125     {
 126       e.printStackTrace();
 127     } finally
 128     {
 129       closeStream(zipStream);
 130     }
 131   }
 132
 133   /**
 134    * Closes the input stream, swallowing all exceptions
 135    *
 136    * @param is
 137    */
 138   protected void closeStream(InputStream is)
 139   {
 140     if (is != null)
 141     {
 142       try
 143       {
 144         is.close();
 145       } catch (IOException e)
 146       {
 147         // ignore
 148       }
 149     }
 150   }
 151
 152   /**
 153    * Reads, parses and stores the OBO file data
 154    *
 155    * @param is
 156    * @throws ParseException
 157    * @throws IOException
 158    */
 159   protected void loadOboFile(InputStream is)
 160           throws ParseException, IOException
 161   {
 162     BufferedReader oboFile = new BufferedReader(new InputStreamReader(is));
 163     OboParser parser = new OboParser();
 164     ontology = parser.parseOBO(oboFile, "SO", "the SO ontology");
 165     isA = ontology.getTerm("is_a");
 166     storeTermAliases();
 167   }
 168
 169   /**
 170    * Stores a lookup table of terms by description or synonym. Note that
 171    * description is not guaranteed unique. Where duplicate descriptions are
 172    * found, try to discard the term that is flagged as obsolete. However we do
 173    * store obsolete terms where there is no duplication of description.
 174    */
 175   protected void storeTermAliases()
 176   {
 177     Set<String> ambiguous = new HashSet<>();
 178
 179     for (Term term : ontology.getTerms())
 180     {
 181       if (term instanceof Impl)
 182       {
 183         boolean newTermIsObsolete = isObsolete(term);
 184         String description = term.getDescription();
 185         if (description != null)
 186         {
 187           description = canonicalise(description);
 188           Term replaced = aliases.get(description);
 189           if (replaced != null)
 190           {
 191             boolean oldTermIsObsolete = isObsolete(replaced);
 192             if (newTermIsObsolete && !oldTermIsObsolete)
 193             {
 194               System.err.println("Ignoring " + term.getName()
 195                       + " as obsolete and duplicated by "
 196                       + replaced.getName());
 197               term = replaced;
 198             }
 199             else if (!newTermIsObsolete && oldTermIsObsolete)
 200             {
 201               System.err.println("Ignoring " + replaced.getName()
 202                       + " as obsolete and duplicated by " + term.getName());
 203             }
 204             else
 205             {
 206               System.err.println("Warning: " + term.getName()
 207                       + " has replaced " + replaced.getName()
 208                       + " for lookup of '" + description + "'");
 209             }
 210           }
 211           aliases.put(description, term);
 212
 213           /*
 214            * also store synonyms if not ambiguous
 215            */
 216           if (!newTermIsObsolete)
 217           {
 218             for (Object syn : term.getSynonyms())
 219             {
 220               String name = ((Synonym) syn).getName();
 221               String synonym = canonicalise(name);
 222               if (aliases.containsKey(synonym))
 223               {
 224                 final Term found = aliases.get(synonym);
 225                 if (found != term)
 226                 {
 227                   /*
 228                    * this alias is ambiguous - matches description,
 229                    * or an alias, of another term
 230                    */
 231                   String msg = String.format(
 232                           "Ambiguous synonym %s for '%s:%s' and '%s:%s'",
 233                           synonym, term.getName(), term.getDescription(),
 234                           found.getName(), found.getDescription());
 235                   System.err.println(msg);
 236
 237                   /*
 238                    * preserve any entry whose canonical description happens to match
 239                    * a synonym (NMD_transcript is a valid description, and also
 240                    * a synonym for NMD_transcript_variant)
 241                    * also preserve a parent (more general) term
 242                    */
 243                   if (synonym.equals(canonicalise(found.getDescription()))
 244                           || termIsA(term, found))
 245                   {
 246                     // leave it alone
 247                   }
 248                   /*
 249                    * replace a specialised term with a more general one
 250                    * with the same alias
 251                    */
 252                   // else if
 253                   // (synonym.equals(canonicalise(term.getDescription())))
 254                   else if (termIsA(found, term))
 255                   {
 256                     aliases.put(synonym, term);
 257                   }
 258                   else
 259                   {
 260                     ambiguous.add(synonym);
 261                   }
 262                 }
 263               }
 264               else
 265               {
 266                 aliases.put(synonym, term);
 267               }
 268             }
 269           }
 270         }
 271       }
 272     }
 273
 274     /*
 275      * remove ambiguous synonyms for safety;
 276      * problem: what if a synonym matches a description?
 277      * only one case found:
 278      * nmd_transcript is synonym for SO:0001621:NMD_transcript_variant
 279      * and also the description for SO:0002114:NMD_transcript
 280      */
 281     for (String syn : ambiguous)
 282     {
 283       aliases.remove(syn);
 284     }
 285   }
 286
 287   /**
 288    * Converts a string to lower case and changes hyphens and spaces to
 289    * underscores
 290    *
 291    * @param s
 292    * @return
 293    */
 294   static String canonicalise(String s)
 295   {
 296     return s == null ? null
 297             : s.toLowerCase().replace('-', '_').replace(' ', '_');
 298   }
 299
 300   /**
 301    * Answers true if the term has property "is_obsolete" with value true, else
 302    * false
 303    *
 304    * @param term
 305    * @return
 306    */
 307   public static boolean isObsolete(Term term)
 308   {
 309     Annotation ann = term.getAnnotation();
 310     if (ann != null)
 311     {
 312       try
 313       {
 314         if (Boolean.TRUE.equals(ann.getProperty("is_obsolete")))
 315         {
 316           return true;
 317         }
 318       } catch (NoSuchElementException e)
 319       {
 320         // fall through to false
 321       }
 322     }
 323     return false;
 324   }
 325
 326   /**
 327    * Test whether the given Sequence Ontology term is nucleotide_match (either
 328    * directly or via is_a relationship)
 329    *
 330    * @param soTerm
 331    *          SO name or description
 332    * @return
 333    */
 334   public boolean isNucleotideMatch(String soTerm)
 335   {
 336     return isA(soTerm, NUCLEOTIDE_MATCH);
 337   }
 338
 339   /**
 340    * Test whether the given Sequence Ontology term is protein_match (either
 341    * directly or via is_a relationship)
 342    *
 343    * @param soTerm
 344    *          SO name or description
 345    * @return
 346    */
 347   public boolean isProteinMatch(String soTerm)
 348   {
 349     return isA(soTerm, PROTEIN_MATCH);
 350   }
 351
 352   /**
 353    * Test whether the given Sequence Ontology term is polypeptide (either
 354    * directly or via is_a relationship)
 355    *
 356    * @param soTerm
 357    *          SO name or description
 358    * @return
 359    */
 360   public boolean isPolypeptide(String soTerm)
 361   {
 362     return isA(soTerm, POLYPEPTIDE);
 363   }
 364
 365   /**
 366    * Returns true if the given term has a (direct or indirect) 'isA'
 367    * relationship with the parent
 368    *
 369    * @param child
 370    * @param parent
 371    * @return
 372    */
 373   @Override
 374   public boolean isA(String child, String parent)
 375   {
 376     if (child == null || parent == null)
 377     {
 378       return false;
 379     }
 380     /*
 381      * optimise trivial checks like isA("CDS", "CDS")
 382      */
 383     if (child.equals(parent))
 384     {
 385       termFound(child);
 386       return true;
 387     }
 388
 389     Term childTerm = getTerm(child);
 390     if (childTerm != null)
 391     {
 392       termFound(child);
 393     }
 394     else
 395     {
 396       termNotFound(child);
 397     }
 398     Term parentTerm = getTerm(parent);
 399
 400     return termIsA(childTerm, parentTerm);
 401   }
 402
 403   /**
 404    * Records a valid term queried for, for reporting purposes
 405    *
 406    * @param term
 407    */
 408   private void termFound(String term)
 409   {
 410     synchronized (termsFound)
 411     {
 412       if (!termsFound.contains(term))
 413       {
 414         termsFound.add(term);
 415       }
 416     }
 417   }
 418
 419   /**
 420    * Records an invalid term queried for, for reporting purposes
 421    *
 422    * @param term
 423    */
 424   private void termNotFound(String term)
 425   {
 426     synchronized (termsNotFound)
 427     {
 428       if (!termsNotFound.contains(term))
 429       {
 430         System.err.println("SO term " + term + " invalid");
 431         termsNotFound.add(term);
 432       }
 433     }
 434   }
 435
 436   /**
 437    * Returns true if the childTerm 'isA' parentTerm (directly or indirectly).
 438    *
 439    * @param childTerm
 440    * @param parentTerm
 441    * @return
 442    */
 443   protected synchronized boolean termIsA(Term childTerm, Term parentTerm)
 444   {
 445     /*
 446      * null term could arise from a misspelled SO description
 447      */
 448     if (childTerm == null || parentTerm == null)
 449     {
 450       return false;
 451     }
 452
 453     /*
 454      * recursive search endpoint:
 455      */
 456     if (childTerm == parentTerm)
 457     {
 458       return true;
 459     }
 460
 461     /*
 462      * lazy initialisation - find all of a term's parents (recursively)
 463      * the first time this is called, and save them in a map.
 464      */
 465     if (!termIsA.containsKey(childTerm))
 466     {
 467       findParents(childTerm);
 468     }
 469
 470     List<Term> parents = termIsA.get(childTerm);
 471     for (Term parent : parents)
 472     {
 473       if (termIsA(parent, parentTerm))
 474       {
 475         /*
 476          * add (great-)grandparents to parents list as they are discovered,
 477          * for faster lookup next time
 478          */
 479         if (!parents.contains(parentTerm))
 480         {
 481           parents.add(parentTerm);
 482         }
 483         return true;
 484       }
 485     }
 486
 487     return false;
 488   }
 489
 490   /**
 491    * Finds all the 'isA' parents of the childTerm and stores them as a (possibly
 492    * empty) list.
 493    *
 494    * @param childTerm
 495    */
 496   protected synchronized void findParents(Term childTerm)
 497   {
 498     List<Term> result = new ArrayList<>();
 499     for (Triple triple : ontology.getTriples(childTerm, null, isA))
 500     {
 501       Term parent = triple.getObject();
 502       result.add(parent);
 503
 504       /*
 505        * and search for the parent's parents recursively
 506        */
 507       findParents(parent);
 508     }
 509     termIsA.put(childTerm, result);
 510   }
 511
 512   /**
 513    * Returns the Term for a given name (e.g. "SO:0000735") or description (e.g.
 514    * "sequence_location"), or alias, or null if not found
 515    *
 516    * @param child
 517    * @return
 518    */
 519   protected Term getTerm(final String nameOrDescription)
 520   {
 521     if (nameOrDescription == null)
 522     {
 523       return null;
 524     }
 525     Term t = aliases.get(canonicalise(nameOrDescription));
 526     if (t == null)
 527     {
 528       try
 529       {
 530         t = ontology.getTerm(nameOrDescription);
 531       } catch (NoSuchElementException e)
 532       {
 533         // not found
 534       }
 535     }
 536     return t;
 537   }
 538
 539   public boolean isSequenceVariant(String soTerm)
 540   {
 541     return isA(soTerm, SEQUENCE_VARIANT);
 542   }
 543
 544   /**
 545    * Sorts (case-insensitive) and returns the list of valid terms queried for
 546    */
 547   @Override
 548   public List<String> termsFound()
 549   {
 550     synchronized (termsFound)
 551     {
 552       Collections.sort(termsFound, String.CASE_INSENSITIVE_ORDER);
 553       return termsFound;
 554     }
 555   }
 556
 557   /**
 558    * Sorts (case-insensitive) and returns the list of invalid terms queried for
 559    */
 560   @Override
 561   public List<String> termsNotFound()
 562   {
 563     synchronized (termsNotFound)
 564     {
 565       Collections.sort(termsNotFound, String.CASE_INSENSITIVE_ORDER);
 566       return termsNotFound;
 567     }
 568   }
 569
 570   /**
 571    * {@inheritDoc}
 572    *
 573    * @throws IllegalStateException
 574    *           if a loop is detected in the ontology
 575    */
 576   @Override
 577   public List<String> getRootParents(final String term)
 578   {
 579     /*
 580      * check in cache first
 581      */
 582     if (rootParents.containsKey(term))
 583     {
 584       return rootParents.get(term);
 585     }
 586     Term t = getTerm(term);
 587     if (t == null)
 588     {
 589       return null;
 590     }
 591
 592     /*
 593      * todo: check for loops using 'seen', allowing for alternate paths e.g.
 594      * stop_gained isA feature_truncation isA feature_variant
 595      * " isA nonsynonymous_variant ... isA geneVariant isA feature_variant
 596      */
 597     List<Term> seen = new ArrayList<>();
 598     List<Term> top = new ArrayList<>();
 599     List<Term> query = new ArrayList<>();
 600     query.add(t);
 601
 602     while (!query.isEmpty())
 603     {
 604       List<Term> nextQuery = new ArrayList<>();
 605       for (Term q : query)
 606       {
 607         Set<Triple> parents = ontology.getTriples(q, null, isA);
 608         if (parents.isEmpty())
 609         {
 610           /*
 611            * q has no parents so is a top level term
 612            */
 613           top.add(q);
 614         }
 615         else
 616         {
 617           /*
 618            * search all parent terms
 619            */
 620           for (Triple triple : parents)
 621           {
 622             Term parent = triple.getObject();
 623             nextQuery.add(parent);
 624           }
 625         }
 626       }
 627       query = nextQuery;
 628     }
 629
 630     List<String> result = new ArrayList<>();
 631     for (Term found : top)
 632     {
 633       String desc = found.getDescription();
 634       if (!result.contains(desc))
 635       {
 636         result.add(desc);
 637       }
 638     }
 639
 640     /*
 641      * save result in cache
 642      */
 643     rootParents.put(term, result);
 644
 645     return result;
 646   }
 647
 648   @Override
 649   public List<String> getParents(String term)
 650   {
 651     List<String> parents = new ArrayList<>();
 652     Term t = getTerm(term);
 653     if (t != null)
 654     {
 655       for (Triple triple : ontology.getTriples(t, null, isA))
 656       {
 657         Term parent = triple.getObject();
 658         parents.add(parent.getDescription());
 659       }
 660     }
 661     return parents;
 662   }
 663
 664   @Override
 665   public boolean isValidTerm(String term)
 666   {
 667     return getTerm(term) != null;
 668   }
 669 }