forester/java/src/org/forester/rio/RIO.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: https://sites.google.com/site/cmzmasek/home/software/forester
  27
  28 package org.forester.rio;
  29
  30 import java.io.File;
  31 import java.io.FileNotFoundException;
  32 import java.io.IOException;
  33 import java.util.ArrayList;
  34 import java.util.Collections;
  35 import java.util.HashMap;
  36 import java.util.HashSet;
  37 import java.util.List;
  38 import java.util.Set;
  39 import java.util.SortedSet;
  40 import java.util.TreeSet;
  41
  42 import org.forester.datastructures.IntMatrix;
  43 import org.forester.io.parsers.PhylogenyParser;
  44 import org.forester.io.parsers.nexus.NexusPhylogeniesParser;
  45 import org.forester.io.parsers.nhx.NHXParser;
  46 import org.forester.io.parsers.nhx.NHXParser.TAXONOMY_EXTRACTION;
  47 import org.forester.io.parsers.util.ParserUtils;
  48 import org.forester.phylogeny.Phylogeny;
  49 import org.forester.phylogeny.PhylogenyMethods;
  50 import org.forester.phylogeny.PhylogenyNode;
  51 import org.forester.phylogeny.data.Taxonomy;
  52 import org.forester.phylogeny.factories.ParserBasedPhylogenyFactory;
  53 import org.forester.phylogeny.factories.PhylogenyFactory;
  54 import org.forester.sdi.GSDI;
  55 import org.forester.sdi.GSDIR;
  56 import org.forester.sdi.SDIException;
  57 import org.forester.sdi.SDIR;
  58 import org.forester.sdi.SDIutil;
  59 import org.forester.sdi.SDIutil.ALGORITHM;
  60 import org.forester.sdi.SDIutil.TaxonomyComparisonBase;
  61 import org.forester.util.BasicDescriptiveStatistics;
  62 import org.forester.util.ForesterUtil;
  63
  64 public final class RIO {
  65
  66     public static final int                  DEFAULT_RANGE = -1;
  67     private Phylogeny[]                      _analyzed_gene_trees;
  68     private List<PhylogenyNode>              _removed_gene_tree_nodes;
  69     private int                              _ext_nodes;
  70     private int                              _int_nodes;
  71     private TaxonomyComparisonBase           _gsdir_tax_comp_base;
  72     private final StringBuilder              _log;
  73     private final BasicDescriptiveStatistics _duplications_stats;
  74     private final boolean                    _produce_log;
  75     private final boolean                    _verbose;
  76     private final REROOTING                  _rerooting;
  77
  78     private RIO( final Phylogeny[] gene_trees,
  79                  final Phylogeny species_tree,
  80                  final ALGORITHM algorithm,
  81                  final REROOTING rerooting,
  82                  final String outgroup,
  83                  int first,
  84                  int last,
  85                  final boolean produce_log,
  86                  final boolean verbose ) throws IOException, SDIException, RIOException {
  87         if ( ( last == DEFAULT_RANGE ) && ( first >= 0 ) ) {
  88             last = gene_trees.length - 1;
  89         }
  90         else if ( ( first == DEFAULT_RANGE ) && ( last >= 0 ) ) {
  91             first = 0;
  92         }
  93         removeSingleDescendentsNodes( species_tree, verbose );
  94         checkPreconditions( gene_trees, species_tree, rerooting, outgroup, first, last );
  95         _produce_log = produce_log;
  96         _verbose = verbose;
  97         _rerooting = rerooting;
  98         _ext_nodes = -1;
  99         _int_nodes = -1;
 100         _log = new StringBuilder();
 101         _gsdir_tax_comp_base = null;
 102         _analyzed_gene_trees = null;
 103         _removed_gene_tree_nodes = null;
 104         _duplications_stats = new BasicDescriptiveStatistics();
 105         inferOrthologs( gene_trees, species_tree, algorithm, outgroup, first, last );
 106     }
 107
 108     public final Phylogeny[] getAnalyzedGeneTrees() {
 109         return _analyzed_gene_trees;
 110     }
 111
 112     public final BasicDescriptiveStatistics getDuplicationsStatistics() {
 113         return _duplications_stats;
 114     }
 115
 116     /**
 117      * Returns the numbers of number of ext nodes in gene trees analyzed (after
 118      * stripping).
 119      *
 120      * @return number of ext nodes in gene trees analyzed (after stripping)
 121      */
 122     public final int getExtNodesOfAnalyzedGeneTrees() {
 123         return _ext_nodes;
 124     }
 125
 126     /**
 127      * Returns the numbers of number of int nodes in gene trees analyzed (after
 128      * stripping).
 129      *
 130      * @return number of int nodes in gene trees analyzed (after stripping)
 131      */
 132     public final int getIntNodesOfAnalyzedGeneTrees() {
 133         return _int_nodes;
 134     }
 135
 136     public final TaxonomyComparisonBase getGSDIRtaxCompBase() {
 137         return _gsdir_tax_comp_base;
 138     }
 139
 140     public final StringBuilder getLog() {
 141         return _log;
 142     }
 143
 144     public final List<PhylogenyNode> getRemovedGeneTreeNodes() {
 145         return _removed_gene_tree_nodes;
 146     }
 147
 148     private final void inferOrthologs( final Phylogeny[] gene_trees,
 149                                        final Phylogeny species_tree,
 150                                        final ALGORITHM algorithm,
 151                                        final String outgroup,
 152                                        final int first,
 153                                        final int last ) throws SDIException, RIOException, FileNotFoundException,
 154             IOException {
 155         if ( algorithm == ALGORITHM.SDIR ) {
 156             // Removes from species_tree all species not found in gene_tree.
 157             PhylogenyMethods.taxonomyBasedDeletionOfExternalNodes( gene_trees[ 0 ], species_tree );
 158             if ( species_tree.isEmpty() ) {
 159                 throw new RIOException( "failed to establish species based mapping between gene and species trees" );
 160             }
 161         }
 162         final Phylogeny[] my_gene_trees;
 163         if ( ( first >= 0 ) && ( last >= first ) && ( last < gene_trees.length ) ) {
 164             my_gene_trees = new Phylogeny[ ( 1 + last ) - first ];
 165             int c = 0;
 166             for( int i = first; i <= last; ++i ) {
 167                 my_gene_trees[ c++ ] = gene_trees[ i ];
 168             }
 169         }
 170         else {
 171             my_gene_trees = gene_trees;
 172         }
 173         if ( log() ) {
 174             preLog( gene_trees, species_tree, algorithm, outgroup, first, last );
 175         }
 176         if ( _verbose && ( my_gene_trees.length >= 4 ) ) {
 177             System.out.println();
 178         }
 179         _analyzed_gene_trees = new Phylogeny[ my_gene_trees.length ];
 180         int gene_tree_ext_nodes = 0;
 181         for( int i = 0; i < my_gene_trees.length; ++i ) {
 182             final Phylogeny gt = my_gene_trees[ i ];
 183             if ( _verbose && ( my_gene_trees.length > 4 ) ) {
 184                 ForesterUtil.updateProgress( ( ( double ) i ) / my_gene_trees.length );
 185             }
 186             if ( i == 0 ) {
 187                 gene_tree_ext_nodes = gt.getNumberOfExternalNodes();
 188             }
 189             else if ( gene_tree_ext_nodes != gt.getNumberOfExternalNodes() ) {
 190                 throw new RIOException( "gene tree #" + ( i + 1 ) + " has a different number of external nodes ("
 191                         + gt.getNumberOfExternalNodes() + ") than the preceding gene trees (" + gene_tree_ext_nodes
 192                         + ")" );
 193             }
 194             if ( algorithm == ALGORITHM.SDIR ) {
 195                 // Removes from gene_tree all species not found in species_tree.
 196                 PhylogenyMethods.taxonomyBasedDeletionOfExternalNodes( species_tree, gt );
 197                 if ( gt.isEmpty() ) {
 198                     throw new RIOException( "failed to establish species based mapping between gene and species trees" );
 199                 }
 200             }
 201             _analyzed_gene_trees[ i ] = performOrthologInference( gt, species_tree, algorithm, outgroup, i );
 202         }
 203         if ( log() ) {
 204             postLog( species_tree );
 205         }
 206         if ( _verbose && ( my_gene_trees.length > 4 ) ) {
 207             System.out.println();
 208             System.out.println();
 209         }
 210     }
 211
 212     private final boolean log() {
 213         return _produce_log;
 214     }
 215
 216     private final void log( final String s ) {
 217         _log.append( s );
 218         _log.append( ForesterUtil.LINE_SEPARATOR );
 219     }
 220
 221     private final void logRemovedGeneTreeNodes() {
 222         log( "Species stripped from gene trees:" );
 223         final SortedSet<String> rn = new TreeSet<String>();
 224         for( final PhylogenyNode n : getRemovedGeneTreeNodes() ) {
 225             final Taxonomy t = n.getNodeData().getTaxonomy();
 226             switch ( getGSDIRtaxCompBase() ) {
 227                 case CODE: {
 228                     rn.add( t.getTaxonomyCode() );
 229                     break;
 230                 }
 231                 case ID: {
 232                     rn.add( t.getIdentifier().toString() );
 233                     break;
 234                 }
 235                 case SCIENTIFIC_NAME: {
 236                     rn.add( t.getScientificName() );
 237                     break;
 238                 }
 239             }
 240         }
 241         for( final String s : rn ) {
 242             log( s );
 243         }
 244         log( "" );
 245     }
 246
 247     private final Phylogeny performOrthologInference( final Phylogeny gene_tree,
 248                                                       final Phylogeny species_tree,
 249                                                       final ALGORITHM algorithm,
 250                                                       final String outgroup,
 251                                                       final int i ) throws SDIException, RIOException {
 252         final Phylogeny assigned_tree;
 253         switch ( algorithm ) {
 254             case SDIR: {
 255                 assigned_tree = performOrthologInferenceBySDI( gene_tree, species_tree );
 256                 break;
 257             }
 258             case GSDIR: {
 259                 assigned_tree = performOrthologInferenceByGSDI( gene_tree, species_tree, outgroup, i );
 260                 break;
 261             }
 262             default: {
 263                 throw new IllegalArgumentException( "illegal algorithm: " + algorithm );
 264             }
 265         }
 266         if ( i == 0 ) {
 267             _ext_nodes = assigned_tree.getNumberOfExternalNodes();
 268             _int_nodes = assigned_tree.getNumberOfInternalNodes();
 269         }
 270         else if ( _ext_nodes != assigned_tree.getNumberOfExternalNodes() ) {
 271             throw new RIOException( "after stripping gene tree #" + ( i + 1 )
 272                     + " has a different number of external nodes (" + assigned_tree.getNumberOfExternalNodes()
 273                     + ") than the preceding gene trees (" + _ext_nodes + ")" );
 274         }
 275         return assigned_tree;
 276     }
 277
 278     private final Phylogeny performOrthologInferenceByGSDI( final Phylogeny gene_tree,
 279                                                             final Phylogeny species_tree,
 280                                                             final String outgroup,
 281                                                             final int i ) throws SDIException, RIOException {
 282         final Phylogeny assigned_tree;
 283         if ( _rerooting == REROOTING.BY_ALGORITHM ) {
 284             final GSDIR gsdir = new GSDIR( gene_tree, species_tree, true, i == 0 );
 285             final List<Phylogeny> assigned_trees = gsdir.getMinDuplicationsSumGeneTrees();
 286             if ( i == 0 ) {
 287                 _removed_gene_tree_nodes = gsdir.getStrippedExternalGeneTreeNodes();
 288                 for( final PhylogenyNode r : _removed_gene_tree_nodes ) {
 289                     if ( !r.getNodeData().isHasTaxonomy() ) {
 290                         throw new RIOException( "node with no (appropriate) taxonomic information found in gene tree #1: "
 291                                 + r.toString() );
 292                     }
 293                 }
 294             }
 295             final List<Integer> shortests = GSDIR.getIndexesOfShortestTree( assigned_trees );
 296             assigned_tree = assigned_trees.get( shortests.get( 0 ) );
 297             if ( log() ) {
 298                 writeStatsToLog( i, gsdir, shortests );
 299             }
 300             if ( i == 0 ) {
 301                 _gsdir_tax_comp_base = gsdir.getTaxCompBase();
 302             }
 303             _duplications_stats.addValue( gsdir.getMinDuplicationsSum() );
 304         }
 305         else {
 306             if ( _rerooting == REROOTING.MIDPOINT ) {
 307                 PhylogenyMethods.midpointRoot( gene_tree );
 308             }
 309             else if ( _rerooting == REROOTING.OUTGROUP ) {
 310                 final PhylogenyNode n = gene_tree.getNode( outgroup );
 311                 gene_tree.reRoot( n );
 312             }
 313             final GSDI gsdi = new GSDI( gene_tree, species_tree, true, true, true );
 314             _removed_gene_tree_nodes = gsdi.getStrippedExternalGeneTreeNodes();
 315             for( final PhylogenyNode r : _removed_gene_tree_nodes ) {
 316                 if ( !r.getNodeData().isHasTaxonomy() ) {
 317                     throw new RIOException( "node with no (appropriate) taxonomic information found in gene tree #1: "
 318                             + r.toString() );
 319                 }
 320             }
 321             assigned_tree = gene_tree;
 322             if ( i == 0 ) {
 323                 _gsdir_tax_comp_base = gsdi.getTaxCompBase();
 324             }
 325             _duplications_stats.addValue( gsdi.getDuplicationsSum() );
 326         }
 327         return assigned_tree;
 328     }
 329
 330     private final Phylogeny performOrthologInferenceBySDI( final Phylogeny gene_tree, final Phylogeny species_tree )
 331             throws SDIException {
 332         final SDIR sdir = new SDIR();
 333         return sdir.infer( gene_tree, species_tree, false, true, true, true, 1 )[ 0 ];
 334     }
 335
 336     private final void postLog( final Phylogeny species_tree ) {
 337         log( "" );
 338         if ( getRemovedGeneTreeNodes().size() > 0 ) {
 339             logRemovedGeneTreeNodes();
 340         }
 341         log( "Species tree external nodes (after stripping)   : " + species_tree.getNumberOfExternalNodes() );
 342         log( "Species tree polytomies (after stripping)       : "
 343                 + PhylogenyMethods.countNumberOfPolytomies( species_tree ) );
 344         log( "Taxonomy linking based on                       : " + getGSDIRtaxCompBase() );
 345         final java.text.DecimalFormat df = new java.text.DecimalFormat( "0.#" );
 346         log( "Gene trees analyzed                             : " + _duplications_stats.getN() );
 347         log( "Mean number of duplications                     : " + df.format( _duplications_stats.arithmeticMean() )
 348                 + " (sd: " + df.format( _duplications_stats.sampleStandardDeviation() ) + ")" + " ("
 349                 + df.format( ( 100.0 * _duplications_stats.arithmeticMean() ) / getIntNodesOfAnalyzedGeneTrees() )
 350                 + "%)" );
 351         if ( _duplications_stats.getN() > 3 ) {
 352             log( "Median number of duplications                   : " + df.format( _duplications_stats.median() )
 353                     + " (" + df.format( ( 100.0 * _duplications_stats.median() ) / getIntNodesOfAnalyzedGeneTrees() )
 354                     + "%)" );
 355         }
 356         log( "Minimum duplications                            : " + ( int ) _duplications_stats.getMin() + " ("
 357                 + df.format( ( 100.0 * _duplications_stats.getMin() ) / getIntNodesOfAnalyzedGeneTrees() ) + "%)" );
 358         log( "Maximum duplications                            : " + ( int ) _duplications_stats.getMax() + " ("
 359                 + df.format( ( 100.0 * _duplications_stats.getMax() ) / getIntNodesOfAnalyzedGeneTrees() ) + "%)" );
 360         log( "Gene tree internal nodes                        : " + getIntNodesOfAnalyzedGeneTrees() );
 361         log( "Gene tree external nodes                        : " + getExtNodesOfAnalyzedGeneTrees() );
 362     }
 363
 364     private final void preLog( final Phylogeny[] gene_trees,
 365                                final Phylogeny species_tree,
 366                                final ALGORITHM algorithm,
 367                                final String outgroup,
 368                                final int first,
 369                                final int last ) {
 370         log( "Number of gene trees (total)                    : " + gene_trees.length );
 371         log( "Algorithm                                       : " + algorithm );
 372         log( "Species tree external nodes (prior to stripping): " + species_tree.getNumberOfExternalNodes() );
 373         log( "Species tree polytomies (prior to stripping)    : "
 374                 + PhylogenyMethods.countNumberOfPolytomies( species_tree ) );
 375         String rs = "";
 376         switch ( _rerooting ) {
 377             case BY_ALGORITHM: {
 378                 rs = "minimizing duplications";
 379                 break;
 380             }
 381             case MIDPOINT: {
 382                 rs = "midpoint";
 383                 break;
 384             }
 385             case OUTGROUP: {
 386                 rs = "outgroup: " + outgroup;
 387                 break;
 388             }
 389             case NONE: {
 390                 rs = "none";
 391                 break;
 392             }
 393         }
 394         log( "Re-rooting                                      : " + rs );
 395         if ( ( first >= 0 ) || ( last >= 0 ) ) {
 396             log( "Gene trees analyzed range                       : " + first + "-" + last );
 397         }
 398         if ( _rerooting == REROOTING.BY_ALGORITHM ) {
 399             writeLogSubHeader();
 400         }
 401     }
 402
 403     private final void writeLogSubHeader() {
 404         _log.append( ForesterUtil.LINE_SEPARATOR );
 405         _log.append( "Some information about duplication numbers in gene trees:" );
 406         _log.append( ForesterUtil.LINE_SEPARATOR );
 407         _log.append( "#" );
 408         _log.append( "\t" );
 409         _log.append( "re-rootings with minimal number of duplications" );
 410         _log.append( "/" );
 411         _log.append( "total root placements" );
 412         _log.append( "\t" );
 413         _log.append( "duplications range" );
 414         _log.append( "\t" );
 415         _log.append( "mininal duplication re-rootings with shortest tree heigth" );
 416         _log.append( ForesterUtil.LINE_SEPARATOR );
 417     }
 418
 419     private final void writeStatsToLog( final int i, final GSDIR gsdir, final List<Integer> shortests ) {
 420         final BasicDescriptiveStatistics stats = gsdir.getDuplicationsSumStats();
 421         _log.append( i );
 422         _log.append( "\t" );
 423         _log.append( gsdir.getMinDuplicationsSumGeneTrees().size() );
 424         _log.append( "/" );
 425         _log.append( stats.getN() );
 426         _log.append( "\t" );
 427         _log.append( ( int ) stats.getMin() );
 428         _log.append( "-" );
 429         _log.append( ( int ) stats.getMax() );
 430         _log.append( "\t" );
 431         _log.append( shortests.size() );
 432         _log.append( ForesterUtil.LINE_SEPARATOR );
 433     }
 434
 435     public final static IntMatrix calculateOrthologTable( final Phylogeny[] analyzed_gene_trees, final boolean sort )
 436             throws RIOException {
 437         final List<String> labels = new ArrayList<String>();
 438         final Set<String> labels_set = new HashSet<String>();
 439         String label;
 440         for( final PhylogenyNode n : analyzed_gene_trees[ 0 ].getExternalNodes() ) {
 441             if ( n.getNodeData().isHasSequence() && !ForesterUtil.isEmpty( n.getNodeData().getSequence().getName() ) ) {
 442                 label = n.getNodeData().getSequence().getName();
 443             }
 444             else if ( n.getNodeData().isHasSequence()
 445                     && !ForesterUtil.isEmpty( n.getNodeData().getSequence().getSymbol() ) ) {
 446                 label = n.getNodeData().getSequence().getSymbol();
 447             }
 448             else if ( !ForesterUtil.isEmpty( n.getName() ) ) {
 449                 label = n.getName();
 450             }
 451             else {
 452                 throw new RIOException( "node " + n + " has no appropriate label" );
 453             }
 454             if ( labels_set.contains( label ) ) {
 455                 throw new RIOException( "label " + label + " is not unique" );
 456             }
 457             labels_set.add( label );
 458             labels.add( label );
 459         }
 460         if ( sort ) {
 461             Collections.sort( labels );
 462         }
 463         final IntMatrix m = new IntMatrix( labels );
 464         int counter = 0;
 465         for( final Phylogeny gt : analyzed_gene_trees ) {
 466             counter++;
 467             PhylogenyMethods.preOrderReId( gt );
 468             final HashMap<String, PhylogenyNode> map = PhylogenyMethods.createNameToExtNodeMap( gt );
 469             for( int x = 0; x < m.size(); ++x ) {
 470                 final String mx = m.getLabel( x );
 471                 final PhylogenyNode nx = map.get( mx );
 472                 if ( nx == null ) {
 473                     throw new RIOException( "node \"" + mx + "\" not present in gene tree #" + counter );
 474                 }
 475                 String my;
 476                 PhylogenyNode ny;
 477                 for( int y = 0; y < m.size(); ++y ) {
 478                     my = m.getLabel( y );
 479                     ny = map.get( my );
 480                     if ( ny == null ) {
 481                         throw new RIOException( "node \"" + my + "\" not present in gene tree #" + counter );
 482                     }
 483                     if ( !PhylogenyMethods.calculateLCAonTreeWithIdsInPreOrder( nx, ny ).isDuplication() ) {
 484                         m.inreaseByOne( x, y );
 485                     }
 486                 }
 487             }
 488         }
 489         return m;
 490     }
 491
 492     public final static RIO executeAnalysis( final File gene_trees_file,
 493                                              final File species_tree_file,
 494                                              final ALGORITHM algorithm,
 495                                              final REROOTING rerooting,
 496                                              final String outgroup,
 497                                              final int first,
 498                                              final int last,
 499                                              final boolean produce_log,
 500                                              final boolean verbose ) throws IOException, SDIException, RIOException {
 501         final Phylogeny[] gene_trees = parseGeneTrees( gene_trees_file );
 502         if ( gene_trees.length < 1 ) {
 503             throw new RIOException( "\"" + gene_trees_file + "\" is devoid of appropriate gene trees" );
 504         }
 505         final Phylogeny species_tree = SDIutil.parseSpeciesTree( gene_trees[ 0 ],
 506                                                                  species_tree_file,
 507                                                                  false,
 508                                                                  true,
 509                                                                  TAXONOMY_EXTRACTION.NO );
 510         return new RIO( gene_trees, species_tree, algorithm, rerooting, outgroup, first, last, produce_log, verbose );
 511     }
 512
 513     public final static RIO executeAnalysis( final File gene_trees_file,
 514                                              final Phylogeny species_tree,
 515                                              final ALGORITHM algorithm,
 516                                              final REROOTING rerooting,
 517                                              final String outgroup,
 518                                              final boolean produce_log,
 519                                              final boolean verbose ) throws IOException, SDIException, RIOException {
 520         return new RIO( parseGeneTrees( gene_trees_file ),
 521                         species_tree,
 522                         algorithm,
 523                         rerooting,
 524                         outgroup,
 525                         DEFAULT_RANGE,
 526                         DEFAULT_RANGE,
 527                         produce_log,
 528                         verbose );
 529     }
 530
 531     public final static RIO executeAnalysis( final File gene_trees_file,
 532                                              final Phylogeny species_tree,
 533                                              final ALGORITHM algorithm,
 534                                              final REROOTING rerooting,
 535                                              final String outgroup,
 536                                              final int first,
 537                                              final int last,
 538                                              final boolean produce_log,
 539                                              final boolean verbose ) throws IOException, SDIException, RIOException {
 540         return new RIO( parseGeneTrees( gene_trees_file ),
 541                         species_tree,
 542                         algorithm,
 543                         rerooting,
 544                         outgroup,
 545                         first,
 546                         last,
 547                         produce_log,
 548                         verbose );
 549     }
 550
 551     public final static RIO executeAnalysis( final Phylogeny[] gene_trees, final Phylogeny species_tree )
 552             throws IOException, SDIException, RIOException {
 553         return new RIO( gene_trees,
 554                         species_tree,
 555                         ALGORITHM.GSDIR,
 556                         REROOTING.BY_ALGORITHM,
 557                         null,
 558                         DEFAULT_RANGE,
 559                         DEFAULT_RANGE,
 560                         false,
 561                         false );
 562     }
 563
 564     public final static RIO executeAnalysis( final Phylogeny[] gene_trees,
 565                                              final Phylogeny species_tree,
 566                                              final ALGORITHM algorithm,
 567                                              final REROOTING rerooting,
 568                                              final String outgroup,
 569                                              final boolean produce_log,
 570                                              final boolean verbose ) throws IOException, SDIException, RIOException {
 571         return new RIO( gene_trees,
 572                         species_tree,
 573                         algorithm,
 574                         rerooting,
 575                         outgroup,
 576                         DEFAULT_RANGE,
 577                         DEFAULT_RANGE,
 578                         produce_log,
 579                         verbose );
 580     }
 581
 582     public final static RIO executeAnalysis( final Phylogeny[] gene_trees,
 583                                              final Phylogeny species_tree,
 584                                              final ALGORITHM algorithm,
 585                                              final REROOTING rerooting,
 586                                              final String outgroup,
 587                                              final int first,
 588                                              final int last,
 589                                              final boolean produce_log,
 590                                              final boolean verbose ) throws IOException, SDIException, RIOException {
 591         return new RIO( gene_trees, species_tree, algorithm, rerooting, outgroup, first, last, produce_log, verbose );
 592     }
 593
 594     private final static void checkPreconditions( final Phylogeny[] gene_trees,
 595                                                   final Phylogeny species_tree,
 596                                                   final REROOTING rerooting,
 597                                                   final String outgroup,
 598                                                   final int first,
 599                                                   final int last ) throws RIOException {
 600         if ( !species_tree.isRooted() ) {
 601             throw new RIOException( "species tree is not rooted" );
 602         }
 603         if ( !( ( last == DEFAULT_RANGE ) && ( first == DEFAULT_RANGE ) )
 604                 && ( ( last < first ) || ( last >= gene_trees.length ) || ( last < 0 ) || ( first < 0 ) ) ) {
 605             throw new RIOException( "attempt to set range (0-based) of gene to analyze to: from " + first + " to "
 606                     + last + " (out of " + gene_trees.length + ")" );
 607         }
 608         if ( ( rerooting == REROOTING.OUTGROUP ) && ForesterUtil.isEmpty( outgroup ) ) {
 609             throw new RIOException( "outgroup not set for midpoint rooting" );
 610         }
 611         if ( ( rerooting != REROOTING.OUTGROUP ) && !ForesterUtil.isEmpty( outgroup ) ) {
 612             throw new RIOException( "outgroup only used for midpoint rooting" );
 613         }
 614         if ( ( rerooting == REROOTING.MIDPOINT )
 615                 && ( PhylogenyMethods.calculateMaxDistanceToRoot( gene_trees[ 0 ] ) <= 0 ) ) {
 616             throw new RIOException( "attempt to use midpoint rooting on gene trees which seem to have no (positive) branch lengths (cladograms)" );
 617         }
 618         if ( rerooting == REROOTING.OUTGROUP ) {
 619             try {
 620                 gene_trees[ 0 ].getNode( outgroup );
 621             }
 622             catch ( final IllegalArgumentException e ) {
 623                 throw new RIOException( "cannot perform re-rooting by outgroup: " + e.getLocalizedMessage() );
 624             }
 625         }
 626     }
 627
 628     private final static Phylogeny[] parseGeneTrees( final File gene_trees_file ) throws FileNotFoundException,
 629             IOException {
 630         final PhylogenyFactory factory = ParserBasedPhylogenyFactory.getInstance();
 631         final PhylogenyParser p = ParserUtils.createParserDependingOnFileType( gene_trees_file, true );
 632         if ( p instanceof NHXParser ) {
 633             final NHXParser nhx = ( NHXParser ) p;
 634             nhx.setReplaceUnderscores( false );
 635             nhx.setIgnoreQuotes( true );
 636             nhx.setTaxonomyExtraction( NHXParser.TAXONOMY_EXTRACTION.YES );
 637         }
 638         else if ( p instanceof NexusPhylogeniesParser ) {
 639             final NexusPhylogeniesParser nex = ( NexusPhylogeniesParser ) p;
 640             nex.setReplaceUnderscores( false );
 641             nex.setIgnoreQuotes( true );
 642             nex.setTaxonomyExtraction( TAXONOMY_EXTRACTION.YES );
 643         }
 644         return factory.create( gene_trees_file, p );
 645     }
 646
 647     private final static void removeSingleDescendentsNodes( final Phylogeny species_tree, final boolean verbose ) {
 648         final int o = PhylogenyMethods.countNumberOfOneDescendantNodes( species_tree );
 649         if ( o > 0 ) {
 650             if ( verbose ) {
 651                 System.out.println( "warning: species tree has " + o
 652                         + " internal nodes with only one descendent which are therefore going to be removed" );
 653             }
 654             PhylogenyMethods.deleteInternalNodesWithOnlyOneDescendent( species_tree );
 655         }
 656     }
 657
 658     public enum REROOTING {
 659         NONE, BY_ALGORITHM, MIDPOINT, OUTGROUP;
 660     }
 661 }