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