forester/java/src/org/forester/sdi/SDIse.java

   1 // $Id:
   2 // FORESTER -- software libraries and applications
   3 // for evolutionary biology research and applications.
   4 //
   5 // Copyright (C) 2008-2009 Christian M. Zmasek
   6 // Copyright (C) 2008-2009 Burnham Institute for Medical Research
   7 // Copyright (C) 2000-2001 Washington University School of Medicine
   8 // and Howard Hughes Medical Institute
   9 // All rights reserved
  10 //
  11 // This library is free software; you can redistribute it and/or
  12 // modify it under the terms of the GNU Lesser General Public
  13 // License as published by the Free Software Foundation; either
  14 // version 2.1 of the License, or (at your option) any later version.
  15 //
  16 // This library is distributed in the hope that it will be useful,
  17 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  18 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
  19 // Lesser General Public License for more details.
  20 //
  21 // You should have received a copy of the GNU Lesser General Public
  22 // License along with this library; if not, write to the Free Software
  23 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA
  24 //
  25 // Contact: phylosoft @ gmail . com
  26 // WWW: www.phylosoft.org/forester
  27
  28 package org.forester.sdi;
  29
  30 import org.forester.phylogeny.Phylogeny;
  31 import org.forester.phylogeny.PhylogenyMethods;
  32 import org.forester.phylogeny.PhylogenyNode;
  33 import org.forester.phylogeny.data.Event;
  34
  35 /*
  36  * Implements our algorithm for speciation - duplication inference (SDI). <p>
  37  * Reference: </p> <ul> <li>Zmasek, C.M. and Eddy, S.R. (2001) "A simple
  38  * algorithm to infer gene duplication and speciation events on a gene tree".
  39  * Bioinformatics, in press. </ul> <p> The initialization is accomplished by:
  40  * </p> <ul> <li>method "linkExtNodesOfG()" of class SDI: setting the links for
  41  * the external nodes of the gene tree <li>"preorderReID(int)" from class
  42  * Phylogeny: numbering of nodes of the species tree in preorder <li>the
  43  * optional stripping of the species tree is accomplished by method
  44  * "stripTree(Phylogeny,Phylogeny)" of class Phylogeny </ul> <p> The recursion
  45  * part is accomplished by this class' method
  46  * "geneTreePostOrderTraversal(PhylogenyNode)". <p> Requires JDK 1.2 or greater.
  47  *
  48  * @see SDI#linkNodesOfG()
  49  *
  50  * @see Phylogeny#preorderReID(int)
  51  *
  52  * @see
  53  * PhylogenyMethods#taxonomyBasedDeletionOfExternalNodes(Phylogeny,Phylogeny)
  54  *
  55  * @see #geneTreePostOrderTraversal(PhylogenyNode)
  56  *
  57  * @author Christian M. Zmasek
  58  *
  59  * @version 1.102 -- last modified: 10/02/01
  60  */
  61 public class SDIse extends SDI {
  62
  63     /**
  64      * Constructor which sets the gene tree and the species tree to be compared.
  65      * species_tree is the species tree to which the gene tree gene_tree will be
  66      * compared to - with method "infer(boolean)". Both Trees must be completely
  67      * binary and rooted. The actual inference is accomplished with method
  68      * "infer(boolean)". The mapping cost L can then be calculated with method
  69      * "computeMappingCost()".
  70      * <p>
  71      * (Last modified: 01/11/01)
  72      *
  73      * @see #infer(boolean)
  74      * @see SDI#computeMappingCostL()
  75      * @param gene_tree
  76      *            reference to a rooted binary gene Phylogeny to which assign
  77      *            duplication vs speciation, must have species names in the
  78      *            species name fields for all external nodes
  79      * @param species_tree
  80      *            reference to a rooted binary species Phylogeny which might get
  81      *            stripped in the process, must have species names in the
  82      *            species name fields for all external nodes
  83      * @throws SDIException
  84      */
  85     public SDIse( final Phylogeny gene_tree, final Phylogeny species_tree ) throws SDIException {
  86         super( gene_tree, species_tree );
  87         _duplications_sum = 0;
  88         PhylogenyMethods.preOrderReId( getSpeciesTree() );
  89         linkNodesOfG();
  90         geneTreePostOrderTraversal( getGeneTree().getRoot() );
  91     }
  92
  93     // Helper method for updateM( boolean, PhylogenyNode, PhylogenyNode )
  94     // Calculates M for PhylogenyNode n, given that M for the two children
  95     // of n has been calculated.
  96     // (Last modified: 10/02/01)
  97     private void calculateMforNode( final PhylogenyNode n ) {
  98         if ( !n.isExternal() ) {
  99             final boolean was_duplication = n.isDuplication();
 100             PhylogenyNode a = n.getChildNode1().getLink();
 101             PhylogenyNode b = n.getChildNode2().getLink();
 102             while ( a != b ) {
 103                 if ( a.getId() > b.getId() ) {
 104                     a = a.getParent();
 105                 }
 106                 else {
 107                     b = b.getParent();
 108                 }
 109             }
 110             n.setLink( a );
 111             Event event = null;
 112             if ( ( a == n.getChildNode1().getLink() ) || ( a == n.getChildNode2().getLink() ) ) {
 113                 event = Event.createSingleDuplicationEvent();
 114                 if ( !was_duplication ) {
 115                     ++_duplications_sum;
 116                 }
 117             }
 118             else {
 119                 event = Event.createSingleSpeciationEvent();
 120                 if ( was_duplication ) {
 121                     --_duplications_sum;
 122                 }
 123             }
 124             n.getNodeData().setEvent( event );
 125         }
 126     } // calculateMforNode( PhylogenyNode )
 127
 128     /**
 129      * Traverses the subtree of PhylogenyNode g in postorder, calculating the
 130      * mapping function M, and determines which nodes represent speciation
 131      * events and which ones duplication events.
 132      * <p>
 133      * Preconditions: Mapping M for external nodes must have been calculated and
 134      * the species tree must be labelled in preorder.
 135      * <p>
 136      * (Last modified: 01/11/01)
 137      *
 138      * @param g
 139      *            starting node of a gene tree - normally the root
 140      */
 141     void geneTreePostOrderTraversal( final PhylogenyNode g ) {
 142         PhylogenyNode a, b;
 143         if ( !g.isExternal() ) {
 144             geneTreePostOrderTraversal( g.getChildNode( 0 ) );
 145             geneTreePostOrderTraversal( g.getChildNode( 1 ) );
 146             a = g.getChildNode( 0 ).getLink();
 147             b = g.getChildNode( 1 ).getLink();
 148             while ( a != b ) {
 149                 if ( a.getId() > b.getId() ) {
 150                     a = a.getParent();
 151                 }
 152                 else {
 153                     b = b.getParent();
 154                 }
 155             }
 156             g.setLink( a );
 157             // Determines whether dup. or spec.
 158             Event event = null;
 159             if ( ( a == g.getChildNode( 0 ).getLink() ) || ( a == g.getChildNode( 1 ).getLink() ) ) {
 160                 event = Event.createSingleDuplicationEvent();
 161                 ++_duplications_sum;
 162             }
 163             else {
 164                 event = Event.createSingleSpeciationEvent();
 165             }
 166             g.getNodeData().setEvent( event );
 167         }
 168     } // geneTreePostOrderTraversal( PhylogenyNode )
 169
 170     /**
 171      * Updates the mapping function M after the root of the gene tree has been
 172      * moved by one branch. It calculates M for the root of the gene tree and
 173      * one of its two children.
 174      * <p>
 175      * To be used ONLY by method "SDIunrooted.fastInfer(Phylogeny,Phylogeny)".
 176      * <p>
 177      * (Last modfied: 10/02/01)
 178      *
 179      * @param prev_root_was_dup
 180      *            true if the previous root was a duplication, false otherwise
 181      * @param prev_root_c1
 182      *            child 1 of the previous root
 183      * @param prev_root_c2
 184      *            child 2 of the previous root
 185      * @return number of duplications which have been assigned in gene tree
 186      */
 187     int updateM( final boolean prev_root_was_dup, final PhylogenyNode prev_root_c1, final PhylogenyNode prev_root_c2 ) {
 188         final PhylogenyNode root = getGeneTree().getRoot();
 189         if ( ( root.getChildNode1() == prev_root_c1 ) || ( root.getChildNode2() == prev_root_c1 ) ) {
 190             calculateMforNode( prev_root_c1 );
 191         }
 192         else {
 193             calculateMforNode( prev_root_c2 );
 194         }
 195         Event event = null;
 196         if ( prev_root_was_dup ) {
 197             event = Event.createSingleDuplicationEvent();
 198         }
 199         else {
 200             event = Event.createSingleSpeciationEvent();
 201         }
 202         root.getNodeData().setEvent( event );
 203         calculateMforNode( root );
 204         return getDuplicationsSum();
 205     } // updateM( boolean, PhylogenyNode, PhylogenyNode )
 206 } // End of class SDIse.