src/jalview/analysis/CrossRef.java

   1 package jalview.analysis;\r
   2 \r
   3 import java.util.Enumeration;\r
   4 import java.util.Vector;\r
   5 import java.util.Hashtable;\r
   6 \r
   7 import jalview.datamodel.AlignedCodonFrame;\r
   8 import jalview.datamodel.Alignment;\r
   9 import jalview.datamodel.AlignmentI;\r
  10 import jalview.datamodel.DBRefSource;\r
  11 import jalview.datamodel.DBRefEntry;\r
  12 import jalview.datamodel.Sequence;\r
  13 import jalview.datamodel.SequenceI;\r
  14 import jalview.ws.SequenceFetcher;\r
  15 import jalview.ws.seqfetcher.ASequenceFetcher;\r
  16 \r
  17 /**\r
  18  * Functions for cross-referencing sequence databases. user must first specify\r
  19  * if cross-referencing from protein or dna (set dna==true)\r
  20  * \r
  21  * @author JimP\r
  22  * \r
  23  */\r
  24 public class CrossRef\r
  25 {\r
  26   /**\r
  27    * get the DNA or protein references for a protein or dna sequence\r
  28    * \r
  29    * @param dna\r
  30    * @param rfs\r
  31    * @return\r
  32    */\r
  33   public static DBRefEntry[] findXDbRefs(boolean dna, DBRefEntry[] rfs)\r
  34   {\r
  35     if (dna)\r
  36     {\r
  37       rfs = jalview.util.DBRefUtils.selectRefs(rfs, DBRefSource.PROTEINDBS);\r
  38     }\r
  39     else\r
  40     {\r
  41       rfs = jalview.util.DBRefUtils.selectRefs(rfs,\r
  42               DBRefSource.DNACODINGDBS); // could attempt to find other cross\r
  43       // refs and return here - ie PDB xrefs\r
  44       // (not dna, not protein seq)\r
  45     }\r
  46     return rfs;\r
  47   }\r
  48 \r
  49   public static Hashtable classifyDbRefs(DBRefEntry[] rfs)\r
  50   {\r
  51     Hashtable classes = new Hashtable();\r
  52     classes.put(DBRefSource.PROTEINDBS, jalview.util.DBRefUtils.selectRefs(\r
  53             rfs, DBRefSource.PROTEINDBS));\r
  54     classes.put(DBRefSource.DNACODINGDBS, jalview.util.DBRefUtils\r
  55             .selectRefs(rfs, DBRefSource.DNACODINGDBS));\r
  56     classes.put(DBRefSource.DOMAINDBS, jalview.util.DBRefUtils.selectRefs(\r
  57             rfs, DBRefSource.DOMAINDBS));\r
  58     // classes.put(OTHER, )\r
  59     return classes;\r
  60   }\r
  61 \r
  62   /**\r
  63    * @param dna\r
  64    *                true if seqs are DNA seqs\r
  65    * @param seqs\r
  66    * @return a list of sequence database cross reference source types\r
  67    */\r
  68   public static String[] findSequenceXrefTypes(boolean dna, SequenceI[] seqs)\r
  69   {\r
  70     return findSequenceXrefTypes(dna, seqs, null);\r
  71   }\r
  72 \r
  73   /**\r
  74    * Indirect references are references from other sequences from the dataset to\r
  75    * any of the direct DBRefEntrys on the given sequences.\r
  76    * \r
  77    * @param dna\r
  78    *                true if seqs are DNA seqs\r
  79    * @param seqs\r
  80    * @return a list of sequence database cross reference source types\r
  81    */\r
  82   public static String[] findSequenceXrefTypes(boolean dna,\r
  83           SequenceI[] seqs, AlignmentI dataset)\r
  84   {\r
  85     String[] dbrefs = null;\r
  86     Vector refs = new Vector();\r
  87     for (int s = 0; s < seqs.length; s++)\r
  88     {\r
  89       SequenceI dss = seqs[s];\r
  90       while (dss.getDatasetSequence() != null)\r
  91       {\r
  92         dss = dss.getDatasetSequence();\r
  93       }\r
  94       DBRefEntry[] rfs = findXDbRefs(dna, dss.getDBRef());\r
  95       for (int r = 0; rfs != null && r < rfs.length; r++)\r
  96       {\r
  97         if (!refs.contains(rfs[r].getSource()))\r
  98         {\r
  99           refs.addElement(rfs[r].getSource());\r
 100         }\r
 101       }\r
 102       if (dataset != null)\r
 103       {\r
 104         // search for references to this sequence's direct references.\r
 105         DBRefEntry[] lrfs = CrossRef.findXDbRefs(!dna, seqs[s].getDBRef());\r
 106         Vector rseqs = new Vector();\r
 107         CrossRef.searchDatasetXrefs(seqs[s], !dna, lrfs, dataset, rseqs,\r
 108                 null); // don't need to specify codon frame for mapping here\r
 109         Enumeration lr = rseqs.elements();\r
 110         while (lr.hasMoreElements())\r
 111         {\r
 112           SequenceI rs = (SequenceI) lr.nextElement();\r
 113           DBRefEntry[] xrs = findXDbRefs(dna, rs.getDBRef());\r
 114           for (int r = 0; rfs != null && r < rfs.length; r++)\r
 115           {\r
 116             if (!refs.contains(rfs[r].getSource()))\r
 117             {\r
 118               refs.addElement(rfs[r].getSource());\r
 119             }\r
 120           }\r
 121         }\r
 122       }\r
 123     }\r
 124     if (refs.size() > 0)\r
 125     {\r
 126       dbrefs = new String[refs.size()];\r
 127       refs.copyInto(dbrefs);\r
 128     }\r
 129     return dbrefs;\r
 130   }\r
 131 \r
 132   /*\r
 133    * if (dna) { if (rfs[r].hasMap()) { // most likely this is a protein cross\r
 134    * reference if (!refs.contains(rfs[r].getSource())) {\r
 135    * refs.addElement(rfs[r].getSource()); } } }\r
 136    */\r
 137   public static boolean hasCdnaMap(SequenceI[] seqs)\r
 138   {\r
 139     String[] reftypes = findSequenceXrefTypes(false, seqs);\r
 140     for (int s = 0; s < reftypes.length; s++)\r
 141     {\r
 142       if (reftypes.equals(DBRefSource.EMBLCDS))\r
 143       {\r
 144         return true;\r
 145         // no map\r
 146       }\r
 147     }\r
 148     return false;\r
 149   }\r
 150 \r
 151   public static SequenceI[] getCdnaMap(SequenceI[] seqs)\r
 152   {\r
 153     Vector cseqs = new Vector();\r
 154     for (int s = 0; s < seqs.length; s++)\r
 155     {\r
 156       DBRefEntry[] cdna = findXDbRefs(true, seqs[s].getDBRef());\r
 157       for (int c = 0; c < cdna.length; c++)\r
 158       {\r
 159         if (cdna[c].getSource().equals(DBRefSource.EMBLCDS))\r
 160         {\r
 161           // retrieve CDS dataset sequences\r
 162           // need global dataset sequence retriever/resolver to reuse refs\r
 163           // and construct Mapping entry.\r
 164           // insert gaps in CDS according to peptide gaps.\r
 165           // add gapped sequence to cseqs\r
 166         }\r
 167       }\r
 168     }\r
 169     if (cseqs.size() > 0)\r
 170     {\r
 171       SequenceI[] rsqs = new SequenceI[cseqs.size()];\r
 172       cseqs.copyInto(rsqs);\r
 173       return rsqs;\r
 174     }\r
 175     return null;\r
 176 \r
 177   }\r
 178 \r
 179   /**\r
 180    * \r
 181    * @param dna\r
 182    * @param seqs\r
 183    * @return\r
 184    */\r
 185   public static Alignment findXrefSequences(SequenceI[] seqs, boolean dna,\r
 186           String source)\r
 187   {\r
 188     return findXrefSequences(seqs, dna, source, null);\r
 189   }\r
 190 \r
 191   /**\r
 192    * \r
 193    * @param seqs\r
 194    * @param dna\r
 195    * @param source\r
 196    * @param dataset\r
 197    *                alignment to search for product sequences.\r
 198    * @return products (as dataset sequences)\r
 199    */\r
 200   public static Alignment findXrefSequences(SequenceI[] seqs, boolean dna,\r
 201           String source, AlignmentI dataset)\r
 202   {\r
 203     Vector rseqs = new Vector();\r
 204     Alignment ral = null;\r
 205     AlignedCodonFrame cf = new AlignedCodonFrame(0); // nominal width\r
 206     for (int s = 0; s < seqs.length; s++)\r
 207     {\r
 208       SequenceI dss = seqs[s];\r
 209       while (dss.getDatasetSequence() != null)\r
 210       {\r
 211         dss = dss.getDatasetSequence();\r
 212       }\r
 213       boolean found = false;\r
 214       DBRefEntry[] xrfs = CrossRef.findXDbRefs(dna, dss.getDBRef());\r
 215       if ((xrfs == null || xrfs.length == 0) && dataset != null)\r
 216       {\r
 217         System.out.println("Attempting to find ds Xrefs refs.");\r
 218         DBRefEntry[] lrfs = CrossRef.findXDbRefs(!dna, seqs[s].getDBRef()); // less\r
 219         // ambiguous\r
 220         // would\r
 221         // be a\r
 222         // 'find\r
 223         // primary\r
 224         // dbRefEntry'\r
 225         // method.\r
 226         // filter for desired source xref here\r
 227         found = CrossRef.searchDatasetXrefs(dss, !dna, lrfs, dataset,\r
 228                 rseqs, cf);\r
 229       }\r
 230       for (int r = 0; xrfs != null && r < xrfs.length; r++)\r
 231       {\r
 232         if (source != null && !source.equals(xrfs[r].getSource()))\r
 233           continue;\r
 234         if (xrfs[r].hasMap())\r
 235         {\r
 236           if (xrfs[r].getMap().getTo() != null)\r
 237           {\r
 238             Sequence rsq = new Sequence(xrfs[r].getMap().getTo());\r
 239             rseqs.addElement(rsq);\r
 240             if (xrfs[r].getMap().getMap().getFromRatio() != xrfs[r]\r
 241                     .getMap().getMap().getToRatio())\r
 242             {\r
 243               // get sense of map correct for adding to product alignment.\r
 244               if (dna)\r
 245               {\r
 246                 // map is from dna seq to a protein product\r
 247                 cf.addMap(dss, rsq, xrfs[r].getMap().getMap());\r
 248               }\r
 249               else\r
 250               {\r
 251                 // map should be from protein seq to its coding dna\r
 252                 cf.addMap(rsq, dss, xrfs[r].getMap().getMap().getInverse());\r
 253               }\r
 254             }\r
 255             found = true;\r
 256           }\r
 257         }\r
 258         if (!found)\r
 259         {\r
 260           // do a bit more work - search for sequences with references matching\r
 261           // xrefs on this sequence.\r
 262           if (dataset != null)\r
 263           {\r
 264             found |= searchDataset(dss, xrfs[r], dataset, rseqs, cf);\r
 265             if (found)\r
 266               xrfs[r] = null; // we've recovered seqs for this one.\r
 267           }\r
 268         }\r
 269       }\r
 270       if (!found)\r
 271       {\r
 272         if (xrfs != null && xrfs.length > 0)\r
 273         {\r
 274           // Try and get the sequence reference...\r
 275           /*\r
 276            * Ideal world - we ask for a sequence fetcher implementation here if\r
 277            * (jalview.io.RunTimeEnvironment.getSequenceFetcher()) (\r
 278            */\r
 279           ASequenceFetcher sftch = new SequenceFetcher();\r
 280           SequenceI[] retrieved = null;\r
 281           int l = xrfs.length;\r
 282           for (int r = 0; r < xrfs.length; r++)\r
 283           {\r
 284             // filter out any irrelevant or irretrievable references\r
 285             if (xrfs[r] == null\r
 286                     || ((source != null && !source.equals(xrfs[r]\r
 287                             .getSource())) || !sftch.isFetchable(xrfs[r]\r
 288                             .getSource())))\r
 289             {\r
 290               l--;\r
 291               xrfs[r] = null;\r
 292             }\r
 293           }\r
 294           if (l > 0)\r
 295           {\r
 296             System.out\r
 297                     .println("Attempting to retrieve cross referenced sequences.");\r
 298             DBRefEntry[] t = new DBRefEntry[l];\r
 299             l = 0;\r
 300             for (int r = 0; r < xrfs.length; r++)\r
 301             {\r
 302               if (xrfs[r] != null)\r
 303                 t[l++] = xrfs[r];\r
 304             }\r
 305             xrfs = t;\r
 306             try\r
 307             {\r
 308               retrieved = sftch.getSequences(xrfs); // problem here is we don't know which of xrfs resulted in which retrieved element\r
 309             } catch (Exception e)\r
 310             {\r
 311               System.err\r
 312                       .println("Problem whilst retrieving cross references for Sequence : "\r
 313                               + seqs[s].getName());\r
 314               e.printStackTrace();\r
 315             }\r
 316             if (retrieved != null)\r
 317             {\r
 318               for (int rs = 0; rs < retrieved.length; rs++)\r
 319               {\r
 320                 // TODO: examine each sequence for 'redundancy'\r
 321                 jalview.datamodel.DBRefEntry[] dbr = retrieved[rs].getDBRef();\r
 322                 if (dbr != null && dbr.length > 0)\r
 323                 {\r
 324                   for (int di = 0; di < dbr.length; di++)\r
 325                   {\r
 326                     // find any entry where we should put in the sequence being cross-referenced into the map\r
 327                     jalview.datamodel.Mapping map = dbr[di].getMap();\r
 328                     if (map != null)\r
 329                     {\r
 330                       if (map.getTo() != null && map.getMap() != null)\r
 331                       {\r
 332                         // should search the local dataset to find any existing candidates for To !\r
 333                         try\r
 334                         {\r
 335                           // compare ms with dss and replace with dss in mapping if map is congruent\r
 336                           SequenceI ms = map.getTo();\r
 337                           int sf = map.getMap().getToLowest();\r
 338                           int st = map.getMap().getToHighest();\r
 339                           SequenceI mappedrg = ms.getSubSequence(sf, st);\r
 340                           SequenceI loc = dss.getSubSequence(sf, st);\r
 341                           if (mappedrg.getLength()>0 && mappedrg.getSequenceAsString().equals(\r
 342                                   loc.getSequenceAsString()))\r
 343                           {\r
 344                             System.err\r
 345                                     .println("Mapping updated for retrieved crossreference");\r
 346                             // method to update all refs of existing To on retrieved sequence with dss and merge any props on To onto dss.\r
 347                             map.setTo(dss);\r
 348                           }\r
 349                         } catch (Exception e)\r
 350                         {\r
 351                           System.err\r
 352                                   .println("Exception when consolidating Mapped sequence set...");\r
 353                           e.printStackTrace(System.err);\r
 354                         }\r
 355                       }\r
 356                     }\r
 357                   }\r
 358                 }\r
 359                 retrieved[rs].updatePDBIds();\r
 360                 rseqs.addElement(retrieved[rs]);\r
 361               }\r
 362             }\r
 363           }\r
 364         }\r
 365       }\r
 366     }\r
 367     if (rseqs.size() > 0)\r
 368     {\r
 369       SequenceI[] rsqs = new SequenceI[rseqs.size()];\r
 370       rseqs.copyInto(rsqs);\r
 371       ral = new Alignment(rsqs);\r
 372       if (cf != null && cf.getProtMappings() != null)\r
 373       {\r
 374         ral.addCodonFrame(cf);\r
 375       }\r
 376     }\r
 377     return ral;\r
 378   }\r
 379 \r
 380   /**\r
 381    * find references to lrfs in the cross-reference set of each sequence in\r
 382    * dataset (that is not equal to sequenceI) Identifies matching DBRefEntry\r
 383    * based on source and accession string only - Map and Version are nulled.\r
 384    * \r
 385    * @param sequenceI\r
 386    * @param lrfs\r
 387    * @param dataset\r
 388    * @param rseqs\r
 389    * @return true if matches were found.\r
 390    */\r
 391   private static boolean searchDatasetXrefs(SequenceI sequenceI,\r
 392           boolean dna, DBRefEntry[] lrfs, AlignmentI dataset, Vector rseqs,\r
 393           AlignedCodonFrame cf)\r
 394   {\r
 395     boolean found = false;\r
 396     if (lrfs == null)\r
 397       return false;\r
 398     for (int i = 0; i < lrfs.length; i++)\r
 399     {\r
 400       DBRefEntry xref = new DBRefEntry(lrfs[i]);\r
 401       // add in wildcards\r
 402       xref.setVersion(null);\r
 403       xref.setMap(null);\r
 404       found = searchDataset(sequenceI, xref, dataset, rseqs, cf, false, dna);\r
 405     }\r
 406     return found;\r
 407   }\r
 408 \r
 409   /**\r
 410    * search a given sequence dataset for references matching cross-references to\r
 411    * the given sequence\r
 412    * \r
 413    * @param sequenceI\r
 414    * @param xrf\r
 415    * @param dataset\r
 416    * @param rseqs\r
 417    *                set of unique sequences\r
 418    * @param cf\r
 419    * @return true if one or more unique sequences were found and added\r
 420    */\r
 421   public static boolean searchDataset(SequenceI sequenceI, DBRefEntry xrf,\r
 422           AlignmentI dataset, Vector rseqs, AlignedCodonFrame cf)\r
 423   {\r
 424     return searchDataset(sequenceI, xrf, dataset, rseqs, cf, true, false);\r
 425   }\r
 426 \r
 427   /**\r
 428    * TODO: generalise to different protein classifications Search dataset for\r
 429    * DBRefEntrys matching the given one (xrf) and add the associated sequence to\r
 430    * rseq.\r
 431    * \r
 432    * @param sequenceI\r
 433    * @param xrf\r
 434    * @param dataset\r
 435    * @param rseqs\r
 436    * @param direct -\r
 437    *                search all references or only subset\r
 438    * @param dna\r
 439    *                search dna or protein xrefs (if direct=false)\r
 440    * @return true if relationship found and sequence added.\r
 441    */\r
 442   public static boolean searchDataset(SequenceI sequenceI, DBRefEntry xrf,\r
 443           AlignmentI dataset, Vector rseqs, AlignedCodonFrame cf,\r
 444           boolean direct, boolean dna)\r
 445   {\r
 446     boolean found = false;\r
 447     if (dataset == null)\r
 448       return false;\r
 449     if (dataset.getSequences() == null)\r
 450     {\r
 451       System.err.println("Empty dataset sequence set - NO VECTOR");\r
 452       return false;\r
 453     }\r
 454     Enumeration e = dataset.getSequences().elements();\r
 455     while (e.hasMoreElements())\r
 456     {\r
 457       SequenceI nxt = (SequenceI) e.nextElement();\r
 458       if (nxt != null)\r
 459       {\r
 460         if (nxt.getDatasetSequence() != null)\r
 461         {\r
 462           System.err\r
 463                   .println("Implementation warning: getProducts passed a dataset alignment without dataset sequences in it!");\r
 464         }\r
 465         if (nxt != sequenceI && nxt != sequenceI.getDatasetSequence())\r
 466         {\r
 467           // look for direct or indirect references in common\r
 468           DBRefEntry[] poss = null, cands = null;\r
 469           if (direct)\r
 470           {\r
 471             cands = jalview.util.DBRefUtils.searchRefs(poss = nxt\r
 472                     .getDBRef(), xrf);\r
 473           }\r
 474           else\r
 475           {\r
 476             cands = jalview.util.DBRefUtils.searchRefs(poss = CrossRef\r
 477                     .findXDbRefs(dna, nxt.getDBRef()), xrf);\r
 478           }\r
 479           if (cands != null)\r
 480           {\r
 481             if (!rseqs.contains(nxt))\r
 482             {\r
 483               rseqs.addElement(nxt);\r
 484               boolean foundmap = cf != null; // don't search if we aren't given\r
 485               // a codon map object\r
 486               for (int r = 0; foundmap && r < cands.length; r++)\r
 487               {\r
 488                 if (cands[r].hasMap())\r
 489                 {\r
 490                   if (cands[r].getMap().getTo() != null\r
 491                           && cands[r].getMap().getMap().getFromRatio() != cands[r]\r
 492                                   .getMap().getMap().getToRatio())\r
 493                   {\r
 494                     foundmap = true;\r
 495                     // get sense of map correct for adding to product alignment.\r
 496                     if (dna)\r
 497                     {\r
 498                       // map is from dna seq to a protein product\r
 499                       cf.addMap(sequenceI, nxt, cands[r].getMap().getMap());\r
 500                     }\r
 501                     else\r
 502                     {\r
 503                       // map should be from protein seq to its coding dna\r
 504                       cf.addMap(nxt, sequenceI, cands[r].getMap().getMap()\r
 505                               .getInverse());\r
 506                     }\r
 507                   }\r
 508                 }\r
 509               }\r
 510               // TODO: add mapping between sequences if necessary\r
 511               found = true;\r
 512             }\r
 513           }\r
 514 \r
 515         }\r
 516       }\r
 517     }\r
 518     return found;\r
 519   }\r
 520 \r
 521   /**\r
 522    * precalculate different products that can be found for seqs in dataset and\r
 523    * return them.\r
 524    * \r
 525    * @param dna\r
 526    * @param seqs\r
 527    * @param dataset\r
 528    * @param fake -\r
 529    *                don't actually build lists - just get types\r
 530    * @return public static Object[] buildXProductsList(boolean dna, SequenceI[]\r
 531    *         seqs, AlignmentI dataset, boolean fake) { String types[] =\r
 532    *         jalview.analysis.CrossRef.findSequenceXrefTypes( dna, seqs,\r
 533    *         dataset); if (types != null) { System.out.println("Xref Types for:\r
 534    *         "+(dna ? "dna" : "prot")); for (int t = 0; t < types.length; t++) {\r
 535    *         System.out.println("Type: " + types[t]); SequenceI[] prod =\r
 536    *         jalview.analysis.CrossRef.findXrefSequences(seqs, dna, types[t]);\r
 537    *         System.out.println("Found " + ((prod == null) ? "no" : "" +\r
 538    *         prod.length) + " products"); if (prod!=null) { for (int p=0; p<prod.length;\r
 539    *         p++) { System.out.println("Prod "+p+":\r
 540    *         "+prod[p].getDisplayId(true)); } } }\r
 541    *  } else { System.out.println("Trying getProducts for\r
 542    * "+al.getSequenceAt(0).getDisplayId(true)); System.out.println("Search DS\r
 543    * Xref for: "+(dna ? "dna" : "prot")); // have a bash at finding the products\r
 544    * amongst all the retrieved sequences. SequenceI[] prod =\r
 545    * jalview.analysis.CrossRef.findXrefSequences(al .getSequencesArray(), dna,\r
 546    * null, ds); System.out.println("Found " + ((prod == null) ? "no" : "" +\r
 547    * prod.length) + " products"); if (prod!=null) { // select non-equivalent\r
 548    * sequences from dataset list for (int p=0; p<prod.length; p++) {\r
 549    * System.out.println("Prod "+p+": "+prod[p].getDisplayId(true)); } }\r
 550    *  } }\r
 551    */\r
 552 }