forester/java/src/org/forester/msa_compactor/MsaCompactor.java

   1 // $Id:
   2 // FORESTER -- software libraries and applications
   3 // for evolutionary biology research and applications.
   4 //
   5 // Copyright (C) 2014 Christian M. Zmasek
   6 // Copyright (C) 2014 Sanford-Burnham Medical Research Institute
   7 // All rights reserved
   8 //
   9 // This library is free software; you can redistribute it and/or
  10 // modify it under the terms of the GNU Lesser General Public
  11 // License as published by the Free Software Foundation; either
  12 // version 2.1 of the License, or (at your option) any later version.
  13 //
  14 // This library is distributed in the hope that it will be useful,
  15 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  16 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
  17 // Lesser General Public License for more details.
  18 //
  19 // You should have received a copy of the GNU Lesser General Public
  20 // License along with this library; if not, write to the Free Software
  21 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA
  22 //
  23 // WWW: https://sites.google.com/site/cmzmasek/home/software/forester
  24
  25 package org.forester.msa_compactor;
  26
  27 import java.awt.Color;
  28 import java.io.File;
  29 import java.io.IOException;
  30 import java.io.Writer;
  31 import java.math.RoundingMode;
  32 import java.text.DecimalFormat;
  33 import java.text.NumberFormat;
  34 import java.util.ArrayList;
  35 import java.util.Arrays;
  36 import java.util.List;
  37 import java.util.SortedSet;
  38 import java.util.TreeSet;
  39
  40 import org.forester.archaeopteryx.Archaeopteryx;
  41 import org.forester.archaeopteryx.Configuration;
  42 import org.forester.evoinference.distance.NeighborJoiningF;
  43 import org.forester.evoinference.distance.PairwiseDistanceCalculator;
  44 import org.forester.evoinference.distance.PairwiseDistanceCalculator.PWD_DISTANCE_METHOD;
  45 import org.forester.evoinference.matrix.distance.BasicSymmetricalDistanceMatrix;
  46 import org.forester.evoinference.tools.BootstrapResampler;
  47 import org.forester.io.parsers.nhx.NHXParser.TAXONOMY_EXTRACTION;
  48 import org.forester.io.parsers.phyloxml.PhyloXmlDataFormatException;
  49 import org.forester.io.parsers.util.ParserUtils;
  50 import org.forester.io.writers.SequenceWriter;
  51 import org.forester.io.writers.SequenceWriter.SEQ_FORMAT;
  52 import org.forester.msa.DeleteableMsa;
  53 import org.forester.msa.Mafft;
  54 import org.forester.msa.Msa;
  55 import org.forester.msa.Msa.MSA_FORMAT;
  56 import org.forester.msa.MsaInferrer;
  57 import org.forester.msa.MsaMethods;
  58 import org.forester.msa.ResampleableMsa;
  59 import org.forester.phylogeny.Phylogeny;
  60 import org.forester.phylogeny.PhylogenyMethods;
  61 import org.forester.phylogeny.PhylogenyMethods.DESCENDANT_SORT_PRIORITY;
  62 import org.forester.phylogeny.PhylogenyNode;
  63 import org.forester.phylogeny.data.NodeVisualData;
  64 import org.forester.phylogeny.data.NodeVisualData.NodeFill;
  65 import org.forester.phylogeny.data.NodeVisualData.NodeShape;
  66 import org.forester.phylogeny.iterators.PhylogenyNodeIterator;
  67 import org.forester.sequence.Sequence;
  68 import org.forester.tools.ConfidenceAssessor;
  69 import org.forester.util.BasicDescriptiveStatistics;
  70 import org.forester.util.ForesterUtil;
  71
  72 public class MsaCompactor {
  73
  74     final private static NumberFormat NF_3                       = new DecimalFormat( "#.###" );
  75     final private static NumberFormat NF_4                       = new DecimalFormat( "#.####" );
  76     private boolean                   _calculate_shannon_entropy = false;
  77     //
  78     private String                    _infile_name               = null;
  79     private final short               _longest_id_length;
  80     //
  81     private String                    _maffts_opts               = "--auto";
  82     private DeleteableMsa             _msa                       = null;
  83     private boolean                   _norm                      = true;
  84     private File                      _out_file_base             = null;
  85     private MSA_FORMAT                _output_format             = MSA_FORMAT.FASTA;
  86     private String                    _path_to_mafft             = null;
  87     private boolean                   _phylogentic_inference     = false;
  88     //
  89     private boolean                   _realign                   = false;
  90     private final SortedSet<String>   _removed_seq_ids;
  91     private final ArrayList<Sequence> _removed_seqs;
  92     private File                      _removed_seqs_out_base     = null;
  93     private int                       _step                      = -1;
  94     private int                       _step_for_diagnostics      = -1;
  95     static {
  96         NF_4.setRoundingMode( RoundingMode.HALF_UP );
  97         NF_3.setRoundingMode( RoundingMode.HALF_UP );
  98     }
  99
 100     public MsaCompactor( final DeleteableMsa msa ) {
 101         _msa = msa;
 102         _removed_seq_ids = new TreeSet<String>();
 103         _longest_id_length = _msa.determineMaxIdLength();
 104         _removed_seqs = new ArrayList<Sequence>();
 105     }
 106
 107     public final Phylogeny calcTree() {
 108         final Phylogeny phy = inferNJphylogeny( PWD_DISTANCE_METHOD.KIMURA_DISTANCE, _msa, false, "" );
 109         PhylogenyMethods.midpointRoot( phy );
 110         PhylogenyMethods.orderAppearance( phy.getRoot(), true, true, DESCENDANT_SORT_PRIORITY.NODE_NAME );
 111         final boolean x = PhylogenyMethods.extractFastaInformation( phy );
 112         if ( !x ) {
 113             final PhylogenyNodeIterator it = phy.iteratorExternalForward();
 114             while ( it.hasNext() ) {
 115                 final PhylogenyNode n = it.next();
 116                 final String name = n.getName().trim();
 117                 if ( !ForesterUtil.isEmpty( name ) ) {
 118                     try {
 119                         ParserUtils.extractTaxonomyDataFromNodeName( n, TAXONOMY_EXTRACTION.AGGRESSIVE );
 120                     }
 121                     catch ( final PhyloXmlDataFormatException e ) {
 122                         // Ignore.
 123                     }
 124                 }
 125             }
 126         }
 127         return phy;
 128     }
 129
 130     public final List<MsaProperties> chart( final int step, final boolean realign, final boolean norm )
 131             throws IOException, InterruptedException {
 132         final GapContribution stats[] = calcGapContribtionsStats( norm );
 133         final List<String> to_remove_ids = new ArrayList<String>();
 134         final List<MsaProperties> msa_props = new ArrayList<MsaProperties>();
 135         for( final GapContribution gap_gontribution : stats ) {
 136             to_remove_ids.add( gap_gontribution.getId() );
 137         }
 138         Phylogeny phy = null;
 139         if ( _phylogentic_inference ) {
 140             System.out.println( "calculating phylogentic tree..." );
 141             System.out.println();
 142             phy = calcTree();
 143         }
 144         if ( !_realign ) {
 145             _step = -1;
 146         }
 147         int x = ForesterUtil.roundToInt( _msa.getNumberOfSequences() / 10.0 );
 148         if ( x < 2 ) {
 149             x = 2;
 150         }
 151         MsaProperties msa_prop = new MsaProperties( _msa, "", _calculate_shannon_entropy );
 152         msa_props.add( msa_prop );
 153         printTableHeader();
 154         printMsaProperties( msa_prop );
 155         System.out.println();
 156         int i = 0;
 157         while ( _msa.getNumberOfSequences() > x ) {
 158             final String id = to_remove_ids.get( i );
 159             _msa.deleteRow( id, false );
 160             if ( realign && isPrintMsaStatsWriteOutfileAndRealign( i ) ) {
 161                 removeGapColumns();
 162                 realignWithMafft();
 163                 msa_prop = new MsaProperties( _msa, id, _calculate_shannon_entropy );
 164                 msa_props.add( msa_prop );
 165                 printMsaProperties( msa_prop );
 166                 System.out.print( "(realigned)" );
 167                 System.out.println();
 168             }
 169             else if ( isPrintMsaStats( i ) ) {
 170                 removeGapColumns();
 171                 msa_prop = new MsaProperties( _msa, id, _calculate_shannon_entropy );
 172                 msa_props.add( msa_prop );
 173                 printMsaProperties( msa_prop );
 174                 System.out.println();
 175             }
 176             ++i;
 177         }
 178         if ( _phylogentic_inference ) {
 179             decorateTree( phy, msa_props, true );
 180             displayTree( phy );
 181         }
 182         return msa_props;
 183     }
 184
 185     public final void decorateTree( final Phylogeny phy, final List<MsaProperties> msa_props, final boolean chart_only ) {
 186         final BasicDescriptiveStatistics length_stats = new BasicDescriptiveStatistics();
 187         for( int i = 0; i < msa_props.size(); ++i ) {
 188             final MsaProperties msa_prop = msa_props.get( i );
 189             final String id = msa_prop.getRemovedSeq();
 190             if ( !ForesterUtil.isEmpty( id ) ) {
 191                 length_stats.addValue( msa_prop.getLength() );
 192             }
 193         }
 194         final double mean = length_stats.arithmeticMean();
 195         final double min = length_stats.getMin();
 196         final double max = length_stats.getMax();
 197         final Color min_color = new Color( 0, 255, 0 );
 198         final Color max_color = new Color( 255, 0, 0 );
 199         final Color mean_color = new Color( 255, 255, 0 );
 200         final PhylogenyNodeIterator it = phy.iteratorExternalForward();
 201         if ( chart_only ) {
 202             while ( it.hasNext() ) {
 203                 final NodeVisualData vis = new NodeVisualData();
 204                 vis.setFillType( NodeFill.SOLID );
 205                 vis.setShape( NodeShape.RECTANGLE );
 206                 vis.setNodeColor( min_color );
 207                 it.next().getNodeData().setNodeVisualData( vis );
 208             }
 209         }
 210         for( int i = 0; i < msa_props.size(); ++i ) {
 211             final MsaProperties msa_prop = msa_props.get( i );
 212             final String id = msa_prop.getRemovedSeq();
 213             if ( !ForesterUtil.isEmpty( id ) ) {
 214                 final PhylogenyNode n = phy.getNode( id );
 215                 n.setName( n.getName() + " [" + i + "]" );
 216                 if ( !chart_only ) {
 217                     final NodeVisualData vis = new NodeVisualData();
 218                     vis.setFillType( NodeFill.SOLID );
 219                     vis.setShape( NodeShape.RECTANGLE );
 220                     vis.setNodeColor( ForesterUtil.calcColor( msa_prop.getLength(), min, max, mean_color, max_color ) );
 221                     n.getNodeData().setNodeVisualData( vis );
 222                 }
 223                 else {
 224                     n.getNodeData()
 225                             .getNodeVisualData()
 226                             .setNodeColor( ForesterUtil.calcColor( msa_prop.getLength(),
 227                                                                    min,
 228                                                                    max,
 229                                                                    mean,
 230                                                                    min_color,
 231                                                                    max_color,
 232                                                                    mean_color ) );
 233                 }
 234             }
 235         }
 236     }
 237
 238     final public void deleteGapColumns( final double max_allowed_gap_ratio ) {
 239         _msa.deleteGapColumns( max_allowed_gap_ratio );
 240     }
 241
 242     public final void displayTree( final Phylogeny phy ) {
 243         final Configuration config = new Configuration();
 244         config.setDisplayAsPhylogram( true );
 245         config.setUseStyle( true );
 246         config.setDisplayTaxonomyCode( false );
 247         config.setDisplayTaxonomyCommonNames( false );
 248         config.setDisplayTaxonomyScientificNames( false );
 249         config.setDisplaySequenceNames( false );
 250         config.setDisplaySequenceSymbols( false );
 251         config.setDisplayGeneNames( false );
 252         config.setShowScale( true );
 253         config.setAddTaxonomyImagesCB( false );
 254         config.setBaseFontSize( 9 );
 255         config.setBaseFontFamilyName( "Arial" );
 256         Archaeopteryx.createApplication( phy, config, _infile_name );
 257     }
 258
 259     final public Msa getMsa() {
 260         return _msa;
 261     }
 262
 263     public final void removeSequencesByMinimalLength( final int min_effective_length ) {
 264         printMsaProperties( new MsaProperties( _msa, "", _calculate_shannon_entropy ) );
 265         System.out.println();
 266         _msa = DeleteableMsa.createInstance( MsaMethods.removeSequencesByMinimalLength( _msa, min_effective_length ) );
 267         removeGapColumns();
 268         printMsaProperties( new MsaProperties( _msa, "", _calculate_shannon_entropy ) );
 269         System.out.println();
 270     }
 271
 272     public final List<MsaProperties> removeViaGapAverage( final double mean_gapiness ) throws IOException,
 273             InterruptedException {
 274         final GapContribution stats[] = calcGapContribtionsStats( _norm );
 275         final List<String> to_remove_ids = new ArrayList<String>();
 276         final List<MsaProperties> msa_props = new ArrayList<MsaProperties>();
 277         for( final GapContribution gap_gontribution : stats ) {
 278             to_remove_ids.add( gap_gontribution.getId() );
 279         }
 280         Phylogeny phy = null;
 281         if ( _phylogentic_inference ) {
 282             System.out.println( "calculating phylogentic tree..." );
 283             System.out.println();
 284             phy = calcTree();
 285         }
 286         printTableHeader();
 287         MsaProperties msa_prop = new MsaProperties( _msa, "", _calculate_shannon_entropy );
 288         msa_props.add( msa_prop );
 289         printMsaProperties( msa_prop );
 290         System.out.println();
 291         int i = 0;
 292         while ( MsaMethods.calcGapRatio( _msa ) > mean_gapiness ) {
 293             final String id = to_remove_ids.get( i );
 294             _removed_seq_ids.add( id );
 295             final Sequence deleted = _msa.deleteRow( id, true );
 296             _removed_seqs.add( deleted );
 297             removeGapColumns();
 298             if ( isPrintMsaStatsWriteOutfileAndRealign( i ) || ( MsaMethods.calcGapRatio( _msa ) <= mean_gapiness ) ) {
 299                 msa_prop = printMsaStatsWriteOutfileAndRealign( _realign, id );
 300                 msa_props.add( msa_prop );
 301                 System.out.println();
 302             }
 303             else if ( isPrintMsaStats( i ) ) {
 304                 msa_prop = new MsaProperties( _msa, id, _calculate_shannon_entropy );
 305                 msa_props.add( msa_prop );
 306                 printMsaProperties( msa_prop );
 307                 System.out.println();
 308             }
 309             ++i;
 310         }
 311         if ( _removed_seqs_out_base != null ) {
 312             final String msg = writeAndAlignRemovedSeqs();
 313             System.out.println();
 314             System.out.println( msg );
 315         }
 316         if ( _phylogentic_inference ) {
 317             decorateTree( phy, msa_props, false );
 318             displayTree( phy );
 319         }
 320         return msa_props;
 321     }
 322
 323     public List<MsaProperties> removeViaLength( final int length ) throws IOException, InterruptedException {
 324         final GapContribution stats[] = calcGapContribtionsStats( _norm );
 325         final List<String> to_remove_ids = new ArrayList<String>();
 326         final List<MsaProperties> msa_props = new ArrayList<MsaProperties>();
 327         for( final GapContribution gap_gontribution : stats ) {
 328             to_remove_ids.add( gap_gontribution.getId() );
 329         }
 330         Phylogeny phy = null;
 331         if ( _phylogentic_inference ) {
 332             System.out.println( "calculating phylogentic tree..." );
 333             System.out.println();
 334             phy = calcTree();
 335         }
 336         printTableHeader();
 337         MsaProperties msa_prop = new MsaProperties( _msa, "", _calculate_shannon_entropy );
 338         msa_props.add( msa_prop );
 339         printMsaProperties( msa_prop );
 340         System.out.println();
 341         int i = 0;
 342         while ( _msa.getLength() > length ) {
 343             final String id = to_remove_ids.get( i );
 344             _removed_seq_ids.add( id );
 345             final Sequence deleted = _msa.deleteRow( id, true );
 346             _removed_seqs.add( deleted );
 347             removeGapColumns();
 348             if ( isPrintMsaStatsWriteOutfileAndRealign( i ) || ( _msa.getLength() <= length ) ) {
 349                 msa_prop = printMsaStatsWriteOutfileAndRealign( _realign, id );
 350                 msa_props.add( msa_prop );
 351                 System.out.println();
 352             }
 353             else if ( isPrintMsaStats( i ) ) {
 354                 msa_prop = new MsaProperties( _msa, id, _calculate_shannon_entropy );
 355                 printMsaProperties( msa_prop );
 356                 msa_props.add( msa_prop );
 357                 System.out.println();
 358             }
 359             ++i;
 360         }
 361         if ( _removed_seqs_out_base != null ) {
 362             final String msg = writeAndAlignRemovedSeqs();
 363             System.out.println();
 364             System.out.println( msg );
 365         }
 366         if ( _phylogentic_inference ) {
 367             decorateTree( phy, msa_props, false );
 368             displayTree( phy );
 369         }
 370         return msa_props;
 371     }
 372
 373     public final List<MsaProperties> removeWorstOffenders( final int to_remove ) throws IOException,
 374             InterruptedException {
 375         final GapContribution stats[] = calcGapContribtionsStats( _norm );
 376         final List<String> to_remove_ids = new ArrayList<String>();
 377         final List<MsaProperties> msa_props = new ArrayList<MsaProperties>();
 378         for( int j = 0; j < to_remove; ++j ) {
 379             to_remove_ids.add( stats[ j ].getId() );
 380         }
 381         Phylogeny phy = null;
 382         if ( _phylogentic_inference ) {
 383             System.out.println( "calculating phylogentic tree..." );
 384             System.out.println();
 385             phy = calcTree();
 386         }
 387         printTableHeader();
 388         MsaProperties msa_prop = new MsaProperties( _msa, "", _calculate_shannon_entropy );
 389         msa_props.add( msa_prop );
 390         printMsaProperties( msa_prop );
 391         System.out.println();
 392         for( int i = 0; i < to_remove_ids.size(); ++i ) {
 393             final String id = to_remove_ids.get( i );
 394             _removed_seq_ids.add( id );
 395             final Sequence deleted = _msa.deleteRow( id, true );
 396             _removed_seqs.add( deleted );
 397             removeGapColumns();
 398             if ( isPrintMsaStatsWriteOutfileAndRealign( i ) || ( i == ( to_remove_ids.size() - 1 ) ) ) {
 399                 msa_prop = printMsaStatsWriteOutfileAndRealign( _realign, id );
 400                 msa_props.add( msa_prop );
 401                 System.out.println();
 402             }
 403             else if ( isPrintMsaStats( i ) ) {
 404                 msa_prop = new MsaProperties( _msa, id, _calculate_shannon_entropy );
 405                 msa_props.add( msa_prop );
 406                 printMsaProperties( msa_prop );
 407                 System.out.println();
 408             }
 409         }
 410         if ( _removed_seqs_out_base != null ) {
 411             final String msg = writeAndAlignRemovedSeqs();
 412             System.out.println();
 413             System.out.println( msg );
 414         }
 415         if ( _phylogentic_inference ) {
 416             decorateTree( phy, msa_props, false );
 417             displayTree( phy );
 418         }
 419         return msa_props;
 420     }
 421
 422     public final void setCalculateNormalizedShannonEntropy( final boolean calculate_shannon_entropy ) {
 423         _calculate_shannon_entropy = calculate_shannon_entropy;
 424     }
 425
 426     public void setInfileName( final String infile_name ) {
 427         _infile_name = infile_name;
 428     }
 429
 430     public final void setMafftOptions( final String maffts_opts ) {
 431         _maffts_opts = maffts_opts;
 432     }
 433
 434     public final void setNorm( final boolean norm ) {
 435         _norm = norm;
 436     }
 437
 438     final public void setOutFileBase( final File out_file_base ) {
 439         _out_file_base = out_file_base;
 440     }
 441
 442     public final void setOutputFormat( final MSA_FORMAT output_format ) {
 443         _output_format = output_format;
 444     }
 445
 446     public void setPathToMafft( final String path_to_mafft ) {
 447         _path_to_mafft = path_to_mafft;
 448     }
 449
 450     public void setPeformPhylogenticInference( final boolean phylogentic_inference ) {
 451         _phylogentic_inference = phylogentic_inference;
 452     }
 453
 454     public final void setRealign( final boolean realign ) {
 455         _realign = realign;
 456     }
 457
 458     public final void setRemovedSeqsOutBase( final File removed_seqs_out_base ) {
 459         _removed_seqs_out_base = removed_seqs_out_base;
 460     }
 461
 462     public final void setStep( final int step ) {
 463         _step = step;
 464     }
 465
 466     public final void setStepForDiagnostics( final int step_for_diagnostics ) {
 467         _step_for_diagnostics = step_for_diagnostics;
 468     }
 469
 470     final public String writeAndAlignRemovedSeqs() throws IOException, InterruptedException {
 471         final StringBuilder msg = new StringBuilder();
 472         final String n = _removed_seqs_out_base + "_" + _removed_seqs.size() + ".fasta";
 473         SequenceWriter.writeSeqs( _removed_seqs, new File( n ), SEQ_FORMAT.FASTA, 100 );
 474         msg.append( "wrote " + _removed_seqs.size() + " removed sequences to " + "\"" + n + "\"" );
 475         if ( _realign ) {
 476             final MsaInferrer mafft = Mafft.createInstance( _path_to_mafft );
 477             final List<String> opts = new ArrayList<String>();
 478             for( final String o : _maffts_opts.split( "\\s" ) ) {
 479                 opts.add( o );
 480             }
 481             final Msa removed_msa = mafft.infer( _removed_seqs, opts );
 482             final Double gr = MsaMethods.calcGapRatio( removed_msa );
 483             String s = _removed_seqs_out_base + "_" + removed_msa.getNumberOfSequences() + "_"
 484                     + removed_msa.getLength() + "_" + ForesterUtil.roundToInt( gr * 100 );
 485             final String suffix = obtainSuffix();
 486             s += suffix;
 487             writeMsa( removed_msa, s, _output_format );
 488             msg.append( ", and as MSA of length " + removed_msa.getLength() + " to \"" + s + "\"" );
 489         }
 490         return msg.toString();
 491     }
 492
 493     final public String writeMsa( final File outfile ) throws IOException {
 494         final Double gr = MsaMethods.calcGapRatio( _msa );
 495         final String s = outfile + "_" + _msa.getNumberOfSequences() + "_" + _msa.getLength() + "_"
 496                 + ForesterUtil.roundToInt( gr * 100 );
 497         writeMsa( _msa, s + obtainSuffix(), _output_format );
 498         return s;
 499     }
 500
 501     final int calcNonGapResidues( final Sequence seq ) {
 502         int ng = 0;
 503         for( int i = 0; i < seq.getLength(); ++i ) {
 504             if ( !seq.isGapAt( i ) ) {
 505                 ++ng;
 506             }
 507         }
 508         return ng;
 509     }
 510
 511     private final GapContribution[] calcGapContribtions( final boolean normalize_for_effective_seq_length ) {
 512         final double gappiness[] = calcGappiness();
 513         final GapContribution stats[] = new GapContribution[ _msa.getNumberOfSequences() ];
 514         for( int row = 0; row < _msa.getNumberOfSequences(); ++row ) {
 515             stats[ row ] = new GapContribution( _msa.getIdentifier( row ) );
 516             for( int col = 0; col < _msa.getLength(); ++col ) {
 517                 if ( !_msa.isGapAt( row, col ) ) {
 518                     stats[ row ].addToValue( gappiness[ col ] );
 519                 }
 520             }
 521             if ( normalize_for_effective_seq_length ) {
 522                 stats[ row ].divideValue( calcNonGapResidues( _msa.getSequence( row ) ) );
 523             }
 524             else {
 525                 stats[ row ].divideValue( _msa.getLength() );
 526             }
 527         }
 528         return stats;
 529     }
 530
 531     final private GapContribution[] calcGapContribtionsStats( final boolean norm ) {
 532         final GapContribution stats[] = calcGapContribtions( norm );
 533         Arrays.sort( stats );
 534         return stats;
 535     }
 536
 537     private final double[] calcGappiness() {
 538         final int l = _msa.getLength();
 539         final double gappiness[] = new double[ l ];
 540         final int seqs = _msa.getNumberOfSequences();
 541         for( int i = 0; i < l; ++i ) {
 542             gappiness[ i ] = ( double ) MsaMethods.calcGapSumPerColumn( _msa, i ) / seqs;
 543         }
 544         return gappiness;
 545     }
 546
 547     private final Phylogeny inferNJphylogeny( final PWD_DISTANCE_METHOD pwd_distance_method,
 548                                               final Msa msa,
 549                                               final boolean write_matrix,
 550                                               final String matrix_name ) {
 551         BasicSymmetricalDistanceMatrix m = null;
 552         switch ( pwd_distance_method ) {
 553             case KIMURA_DISTANCE:
 554                 m = PairwiseDistanceCalculator.calcKimuraDistances( msa );
 555                 break;
 556             case POISSON_DISTANCE:
 557                 m = PairwiseDistanceCalculator.calcPoissonDistances( msa );
 558                 break;
 559             case FRACTIONAL_DISSIMILARITY:
 560                 m = PairwiseDistanceCalculator.calcFractionalDissimilarities( msa );
 561                 break;
 562             default:
 563                 throw new IllegalArgumentException( "invalid pwd method" );
 564         }
 565         if ( write_matrix ) {
 566             try {
 567                 m.write( ForesterUtil.createBufferedWriter( matrix_name ) );
 568             }
 569             catch ( final IOException e ) {
 570                 e.printStackTrace();
 571             }
 572         }
 573         final NeighborJoiningF nj = NeighborJoiningF.createInstance( false, 5 );
 574         final Phylogeny phy = nj.execute( m );
 575         return phy;
 576     }
 577
 578     private final boolean isPrintMsaStats( final int i ) {
 579         return ( ( ( _step == 1 ) && ( _step_for_diagnostics == 1 ) ) || ( ( _step_for_diagnostics > 0 ) && ( ( ( i + 1 ) % _step_for_diagnostics ) == 0 ) ) );
 580     }
 581
 582     private final boolean isPrintMsaStatsWriteOutfileAndRealign( final int i ) {
 583         return ( ( ( _step == 1 ) && ( _step_for_diagnostics == 1 ) ) || ( ( _step > 0 ) && ( ( ( i + 1 ) % _step ) == 0 ) ) );
 584     }
 585
 586     private final StringBuilder msaPropertiesAsSB( final MsaProperties msa_properties ) {
 587         final StringBuilder sb = new StringBuilder();
 588         sb.append( msa_properties.getNumberOfSequences() );
 589         sb.append( "\t" );
 590         sb.append( msa_properties.getLength() );
 591         sb.append( "\t" );
 592         sb.append( NF_4.format( msa_properties.getGapRatio() ) );
 593         if ( _calculate_shannon_entropy ) {
 594             sb.append( "\t" );
 595             sb.append( NF_4.format( msa_properties.getEntropy7() ) );
 596             sb.append( "\t" );
 597             sb.append( NF_4.format( msa_properties.getEntropy21() ) );
 598         }
 599         return sb;
 600     }
 601
 602     private String obtainSuffix() {
 603         if ( _output_format == MSA_FORMAT.FASTA ) {
 604             return ".fasta";
 605         }
 606         else if ( _output_format == MSA_FORMAT.PHYLIP ) {
 607             return ".aln";
 608         }
 609         return "";
 610     }
 611
 612     private final Phylogeny pi( final String matrix, final int boostrap ) {
 613         final Phylogeny master_phy = inferNJphylogeny( PWD_DISTANCE_METHOD.KIMURA_DISTANCE, _msa, true, matrix );
 614         final int seed = 15;
 615         final int n = 100;
 616         final ResampleableMsa resampleable_msa = new ResampleableMsa( _msa );
 617         final int[][] resampled_column_positions = BootstrapResampler.createResampledColumnPositions( _msa.getLength(),
 618                                                                                                       n,
 619                                                                                                       seed );
 620         final Phylogeny[] eval_phys = new Phylogeny[ n ];
 621         for( int i = 0; i < n; ++i ) {
 622             resampleable_msa.resample( resampled_column_positions[ i ] );
 623             eval_phys[ i ] = inferNJphylogeny( PWD_DISTANCE_METHOD.KIMURA_DISTANCE, resampleable_msa, false, null );
 624         }
 625         ConfidenceAssessor.evaluate( "bootstrap", eval_phys, master_phy, true, 1 );
 626         PhylogenyMethods.extractFastaInformation( master_phy );
 627         return master_phy;
 628     }
 629
 630     private final void printMsaProperties( final MsaProperties msa_properties ) {
 631         if ( ( _step == 1 ) || ( _step_for_diagnostics == 1 ) ) {
 632             System.out.print( ForesterUtil.pad( msa_properties.getRemovedSeq(), _longest_id_length, ' ', false ) );
 633             System.out.print( "\t" );
 634         }
 635         System.out.print( msaPropertiesAsSB( msa_properties ) );
 636         System.out.print( "\t" );
 637     }
 638
 639     final private MsaProperties printMsaStatsWriteOutfileAndRealign( final boolean realign, final String id )
 640             throws IOException, InterruptedException {
 641         if ( realign ) {
 642             realignWithMafft();
 643         }
 644         final MsaProperties msa_prop = new MsaProperties( _msa, id, _calculate_shannon_entropy );
 645         printMsaProperties( msa_prop );
 646         final String s = writeOutfile();
 647         System.out.print( "-> " + s + ( realign ? "\t(realigned)" : "" ) );
 648         return msa_prop;
 649     }
 650
 651     private final void printTableHeader() {
 652         if ( ( _step == 1 ) || ( _step_for_diagnostics == 1 ) ) {
 653             System.out.print( ForesterUtil.pad( "Id", _longest_id_length, ' ', false ) );
 654             System.out.print( "\t" );
 655         }
 656         System.out.print( "Seqs" );
 657         System.out.print( "\t" );
 658         System.out.print( "Length" );
 659         System.out.print( "\t" );
 660         System.out.print( "Gaps" );
 661         System.out.print( "\t" );
 662         if ( _calculate_shannon_entropy ) {
 663             System.out.print( "entn7" );
 664             System.out.print( "\t" );
 665             System.out.print( "entn21" );
 666             System.out.print( "\t" );
 667         }
 668         System.out.println();
 669     }
 670
 671     final private void realignWithMafft() throws IOException, InterruptedException {
 672         final MsaInferrer mafft = Mafft.createInstance( _path_to_mafft );
 673         final List<String> opts = new ArrayList<String>();
 674         for( final String o : _maffts_opts.split( "\\s" ) ) {
 675             opts.add( o );
 676         }
 677         _msa = DeleteableMsa.createInstance( mafft.infer( _msa.asSequenceList(), opts ) );
 678     }
 679
 680     final private void removeGapColumns() {
 681         _msa.deleteGapOnlyColumns();
 682     }
 683
 684     private final String writeOutfile() throws IOException {
 685         final String s = writeMsa( _out_file_base );
 686         return s;
 687     }
 688
 689     // Returns null if not path found.
 690     final public static String guessPathToMafft() {
 691         String path;
 692         if ( ForesterUtil.OS_NAME.toLowerCase().indexOf( "win" ) >= 0 ) {
 693             path = "C:\\Program Files\\mafft-win\\mafft.bat";
 694             if ( MsaInferrer.isInstalled( path ) ) {
 695                 return path;
 696             }
 697         }
 698         path = "/home/czmasek/SOFTWARE/MSA/MAFFT/mafft-7.130-without-extensions/scripts/mafft";
 699         if ( MsaInferrer.isInstalled( path ) ) {
 700             return path;
 701         }
 702         path = "/usr/local/bin/mafft";
 703         if ( MsaInferrer.isInstalled( path ) ) {
 704             return path;
 705         }
 706         path = "/usr/bin/mafft";
 707         if ( MsaInferrer.isInstalled( path ) ) {
 708             return path;
 709         }
 710         path = "/bin/mafft";
 711         if ( MsaInferrer.isInstalled( path ) ) {
 712             return path;
 713         }
 714         path = "mafft";
 715         if ( MsaInferrer.isInstalled( path ) ) {
 716             return path;
 717         }
 718         return null;
 719     }
 720
 721     final private static void writeMsa( final Msa msa, final String outfile, final MSA_FORMAT format )
 722             throws IOException {
 723         final Writer w = ForesterUtil.createBufferedWriter( outfile );
 724         msa.write( w, format );
 725         w.close();
 726     }
 727 }