inprogress
[jalview.git] / forester / java / src / org / forester / application / msa_compactor.java
index 996ccf3..19d5caf 100644 (file)
@@ -63,17 +63,18 @@ public class msa_compactor {
     final static private String       STEP_FOR_DIAGNOSTICS_OPTION            = "sd";
     final static private String       MIN_LENGTH_OPTION                      = "ml";
     final static private String       GAP_RATIO_LENGTH_OPTION                = "gr";
-    final static private String       REPORT_ALN_MEAN_IDENTITY               = "q";
+    final static private String       REPORT_ENTROPY                         = "e";
     final static private String       OUTPUT_FORMAT_PHYLIP_OPTION            = "p";
     final static private String       OUTPUT_REMOVED_SEQS_OPTION             = "ro";
     final static private String       MAFFT_OPTIONS                          = "mo";
+    final static private String       PERFORM_PHYLOGENETIC_INFERENCE         = "t";
     //        
     final static private String       PATH_TO_MAFFT_OPTION                   = "mafft";
     final static private String       DO_NOT_NORMALIZE_FOR_EFF_LENGTH_OPTION = "nn";
     final static private String       PRG_NAME                               = "msa_compactor";
     final static private String       PRG_DESC                               = "multiple sequence aligment compactor";
-    final static private String       PRG_VERSION                            = "0.2";
-    final static private String       PRG_DATE                               = "140428";
+    final static private String       PRG_VERSION                            = "0.3";
+    final static private String       PRG_DATE                               = "140508";
     final static private String       E_MAIL                                 = "czmasek@sanfordburham.org";
     final static private String       WWW                                    = "https://sites.google.com/site/cmzmasek/home/software/forester";
 
@@ -100,10 +101,11 @@ public class msa_compactor {
             int step_for_diagnostics = 1;
             int min_length = -1;
             double gap_ratio = -1;
-            boolean report_aln_mean_identity = false;
+            boolean report_entropy = false;
             MSA_FORMAT output_format = MSA_FORMAT.FASTA;
             File removed_seqs_out_base = null;
             String mafft_options = "--auto";
+            boolean perform_phylogenetic_inference = false;
             final List<String> allowed_options = new ArrayList<String>();
             allowed_options.add( REMOVE_WORST_OFFENDERS_OPTION );
             allowed_options.add( AV_GAPINESS_OPTION );
@@ -115,10 +117,11 @@ public class msa_compactor {
             allowed_options.add( STEP_FOR_DIAGNOSTICS_OPTION );
             allowed_options.add( MIN_LENGTH_OPTION );
             allowed_options.add( GAP_RATIO_LENGTH_OPTION );
-            allowed_options.add( REPORT_ALN_MEAN_IDENTITY );
+            allowed_options.add( REPORT_ENTROPY );
             allowed_options.add( OUTPUT_FORMAT_PHYLIP_OPTION );
             allowed_options.add( OUTPUT_REMOVED_SEQS_OPTION );
             allowed_options.add( MAFFT_OPTIONS );
+            allowed_options.add( PERFORM_PHYLOGENETIC_INFERENCE );
             final String dissallowed_options = cla.validateAllowedOptionsAsString( allowed_options );
             if ( dissallowed_options.length() > 0 ) {
                 ForesterUtil.fatalError( PRG_NAME, "unknown option(s): " + dissallowed_options );
@@ -132,9 +135,10 @@ public class msa_compactor {
                 msa = DeleteableMsa.createInstance( GeneralMsaParser.parse( is ) );
             }
             final DescriptiveStatistics initial_msa_stats = MsaMethods.calculateEffectiveLengthStatistics( msa );
-            //System.out.println( initial_msa_stats.toString() );
-            if ( ( cla.isOptionSet( LENGTH_OPTION ) || cla.isOptionSet( REMOVE_WORST_OFFENDERS_OPTION ) || cla
-                    .isOptionSet( AV_GAPINESS_OPTION ) ) && ( out == null ) ) {
+            final boolean chart_only = ( !cla.isOptionSet( LENGTH_OPTION ) )
+                    && ( !cla.isOptionSet( REMOVE_WORST_OFFENDERS_OPTION ) )
+                    && ( !cla.isOptionSet( AV_GAPINESS_OPTION ) );
+            if ( !chart_only && ( out == null ) ) {
                 ForesterUtil.fatalError( PRG_NAME, "outfile file missing" );
             }
             if ( cla.isOptionSet( REMOVE_WORST_OFFENDERS_OPTION ) ) {
@@ -211,8 +215,8 @@ public class msa_compactor {
                     ForesterUtil.fatalError( PRG_NAME, "gap ratio is out of range: " + gap_ratio );
                 }
             }
-            if ( cla.isOptionSet( REPORT_ALN_MEAN_IDENTITY ) ) {
-                report_aln_mean_identity = true;
+            if ( cla.isOptionSet( REPORT_ENTROPY ) ) {
+                report_entropy = true;
             }
             if ( cla.isOptionSet( OUTPUT_FORMAT_PHYLIP_OPTION ) ) {
                 output_format = MSA_FORMAT.PHYLIP;
@@ -229,14 +233,33 @@ public class msa_compactor {
                 if ( cla.isOptionSet( MAFFT_OPTIONS ) ) {
                     mafft_options = cla.getOptionValueAsCleanString( MAFFT_OPTIONS );
                     if ( ForesterUtil.isEmpty( mafft_options ) || ( mafft_options.length() < 3 ) ) {
-                        ForesterUtil.fatalError( PRG_NAME, "gap ratio is out of range: " + gap_ratio );
+                        ForesterUtil.fatalError( PRG_NAME, "illegal or empty MAFFT options: " + mafft_options );
                     }
                 }
             }
-            if ( ( !cla.isOptionSet( LENGTH_OPTION ) && !cla.isOptionSet( REMOVE_WORST_OFFENDERS_OPTION ) && !cla
-                    .isOptionSet( AV_GAPINESS_OPTION ) ) && ( ( out != null ) || ( removed_seqs_out_base != null ) ) ) {
-                ForesterUtil.fatalError( PRG_NAME,
-                                         "chart only, no outfile(s) produced (no need to indicate output file(s))" );
+            else if ( cla.isOptionSet( MAFFT_OPTIONS ) ) {
+                ForesterUtil.fatalError( PRG_NAME, "no need to indicate MAFFT options without realigning" );
+            }
+            if ( cla.isOptionSet( PERFORM_PHYLOGENETIC_INFERENCE ) ) {
+                perform_phylogenetic_inference = true;
+            }
+            if ( chart_only ) {
+                if ( ( out != null ) || ( removed_seqs_out_base != null ) ) {
+                    ForesterUtil
+                            .fatalError( PRG_NAME,
+                                         "chart only, no outfile(s) produced, thus no need to indicate output file(s)" );
+                }
+                if ( !realign && cla.isOptionSet( STEP_OPTION ) ) {
+                    ForesterUtil.fatalError( PRG_NAME,
+                                             "chart only, no re-aligning, thus no need to use step for output and re-aligning; use -"
+                                                     + STEP_FOR_DIAGNOSTICS_OPTION + " instead" );
+                }
+            }
+            if ( perform_phylogenetic_inference ) {
+                if ( step_for_diagnostics != 1 ) {
+                    ForesterUtil.fatalError( PRG_NAME,
+                                             "step for diagnostics reports needs to be set to 1 for tree calculation" );
+                }
             }
             ForesterUtil.printProgramInformation( PRG_NAME,
                                                   PRG_DESC,
@@ -253,7 +276,7 @@ public class msa_compactor {
                     + NF_1.format( initial_msa_stats.arithmeticMean() ) );
             System.out.println( "  Max sequence length                : " + ( ( int ) initial_msa_stats.getMax() ) );
             System.out.println( "  Min sequence length                : " + ( ( int ) initial_msa_stats.getMin() ) );
-            if ( out != null ) {
+            if ( !chart_only ) {
                 System.out.println( "Output                               : " + out );
             }
             else {
@@ -265,22 +288,35 @@ public class msa_compactor {
             if ( worst_remove > 0 ) {
                 System.out.println( "Number of worst offenders to remove  : " + worst_remove );
             }
-            else if ( av_gap > 0 ) {
+            if ( av_gap > 0 ) {
                 System.out.println( "Target gap-ratio                     : " + av_gap );
             }
-            else if ( length > 0 ) {
+            if ( length > 0 ) {
                 System.out.println( "Target MSA length                    : " + length );
             }
-            else {
-                System.out.println( "Chart and diagnostics only           : true" );
+            if ( min_length > 1 ) {
+                System.out.println( "Minimal effective sequence length    : " + min_length );
+            }
+            if ( gap_ratio > -1 ) {
+                System.out.println( "Maximum allowed gap ratio per column : " + gap_ratio );
             }
             if ( ( out != null ) || ( removed_seqs_out_base != null ) ) {
                 System.out.println( "Output format                        : "
                         + ( output_format == MSA_FORMAT.FASTA ? "fasta" : "phylip" ) );
             }
-            System.out.println( "Step for output and re-aligning      : " + step );
+            if ( chart_only && !realign ) {
+                System.out.println( "Step for output and re-aligning      : n/a" );
+            }
+            else {
+                if ( chart_only ) {
+                    System.out.println( "Step for re-aligning                 : " + step );
+                }
+                else {
+                    System.out.println( "Step for output and re-aligning      : " + step );
+                }
+            }
             System.out.println( "Step for diagnostics reports         : " + step_for_diagnostics );
-            System.out.println( "Calculate mean identity              : " + report_aln_mean_identity );
+            System.out.println( "Calculate normalized Shannon Entropy : " + report_entropy );
             if ( !norm ) {
                 System.out.println( "Normalize                            : " + norm );
             }
@@ -288,87 +324,46 @@ public class msa_compactor {
             if ( realign ) {
                 System.out.println( "MAFFT options                        : " + mafft_options );
             }
-            if ( min_length > -1 ) {
-                System.out.println( "Minimal effective sequence length    : " + min_length );
-            }
-            if ( gap_ratio > -1 ) {
-                System.out.println( "Maximum allowed gap ratio per column : " + gap_ratio );
-            }
+            System.out.println( "Simple tree (Kimura distances, NJ)   : " + perform_phylogenetic_inference );
             System.out.println();
-            //
-            //
             final int initial_number_of_seqs = msa.getNumberOfSequences();
             List<MsaProperties> msa_props = null;
-            if ( worst_remove > 0 ) {
-                final MsaCompactor mc = new MsaCompactor( msa );
-                mc.setRealign( realign );
+            final MsaCompactor mc = new MsaCompactor( msa );
+            mc.setInfileName( in.getName() );
+            mc.setNorm( norm );
+            mc.setRealign( realign );
+            if ( realign ) {
+                mc.setPathToMafft( path_to_mafft );
+                mc.setMafftOptions( mafft_options );
+            }
+            mc.setStep( step );
+            mc.setStepForDiagnostics( step_for_diagnostics );
+            mc.setCalculateNormalizedShannonEntropy( report_entropy );
+            mc.setPeformPhylogenticInference( perform_phylogenetic_inference );
+            if ( ( worst_remove > 0 ) || ( av_gap > 0 ) || ( length > 0 ) ) {
                 mc.setOutputFormat( output_format );
-                if ( realign ) {
-                    mc.setPathToMafft( path_to_mafft );
-                    mc.setMafftOptions( mafft_options );
-                }
-                mc.setNorm( norm );
-                mc.setReportAlnMeanIdentity( report_aln_mean_identity );
                 mc.setOutFileBase( out );
                 if ( removed_seqs_out_base != null ) {
                     mc.setRemovedSeqsOutBase( removed_seqs_out_base );
                 }
-                mc.setStep( step );
-                mc.setStepForDiagnostics( step_for_diagnostics );
+            }
+            if ( min_length > 1 ) {
+                mc.removeSequencesByMinimalLength( min_length );
+                mc.writeMsa( new File( "removed" ) );
+            }
+            if ( worst_remove > 0 ) {
                 msa_props = mc.removeWorstOffenders( worst_remove );
             }
             else if ( av_gap > 0 ) {
-                final MsaCompactor mc = new MsaCompactor( msa );
-                mc.setRealign( realign );
-                mc.setOutputFormat( output_format );
-                if ( realign ) {
-                    mc.setPathToMafft( path_to_mafft );
-                    mc.setMafftOptions( mafft_options );
-                }
-                mc.setNorm( norm );
-                mc.setReportAlnMeanIdentity( report_aln_mean_identity );
-                mc.setOutFileBase( out );
-                if ( removed_seqs_out_base != null ) {
-                    mc.setRemovedSeqsOutBase( removed_seqs_out_base );
-                }
-                mc.setStep( step );
-                mc.setStepForDiagnostics( step_for_diagnostics );
                 msa_props = mc.removeViaGapAverage( av_gap );
             }
             else if ( length > 0 ) {
-                final MsaCompactor mc = new MsaCompactor( msa );
-                mc.setRealign( realign );
-                mc.setOutputFormat( output_format );
-                if ( realign ) {
-                    mc.setPathToMafft( path_to_mafft );
-                    mc.setMafftOptions( mafft_options );
-                }
-                mc.setNorm( norm );
-                mc.setReportAlnMeanIdentity( report_aln_mean_identity );
-                mc.setOutFileBase( out );
-                if ( removed_seqs_out_base != null ) {
-                    mc.setRemovedSeqsOutBase( removed_seqs_out_base );
-                }
-                mc.setStep( step );
-                mc.setStepForDiagnostics( step_for_diagnostics );
                 msa_props = mc.removeViaLength( length );
             }
             else {
-                //MsaCompactor.chart( msa, step, realign, norm, path_to_mafft );
-                final MsaCompactor mc = new MsaCompactor( msa );
-                mc.setRealign( realign );
-                if ( realign ) {
-                    mc.setPathToMafft( path_to_mafft );
-                    mc.setMafftOptions( mafft_options );
-                }
-                mc.setNorm( norm );
-                mc.setReportAlnMeanIdentity( report_aln_mean_identity );
-                mc.setOutFileBase( out );
-                mc.setStep( step );
-                mc.setStepForDiagnostics( step_for_diagnostics );
                 msa_props = mc.chart( step, realign, norm );
             }
-            Chart.display( msa_props, initial_number_of_seqs, report_aln_mean_identity, in.toString() );
+            Chart.display( msa_props, initial_number_of_seqs, report_entropy, in.getName() );
         }
         catch ( final IllegalArgumentException iae ) {
             //  iae.printStackTrace(); //TODO remove me
@@ -428,17 +423,17 @@ public class msa_compactor {
         System.out.println( "   -" + STEP_OPTION + "=<integer>   step for output and re-aligning (default: 1)" );
         System.out.println( "   -" + STEP_FOR_DIAGNOSTICS_OPTION
                 + "=<integer>  step for diagnostics reports (default: 1)" );
-        //   System.out.println( "   -" + MIN_LENGTH_OPTION
-        //           + "=<integer>  minimal effecive sequence length (for deleting of shorter sequences)" );
-        //   System.out.println( "   -" + GAP_RATIO_LENGTH_OPTION
-        //           + "=<decimal>  maximal allowed gap ratio per column (for deleting of columms) (0.0-1.0)" );
-        System.out
-                .println( "   -"
-                        + REPORT_ALN_MEAN_IDENTITY
-                        + "             to calculate mean MSA column identity (\"MSA quality\") (not recommended for very large alignments)" );
+        System.out.println( "   -" + REPORT_ENTROPY
+                + "             to calculate normalized Shannon Entropy (not recommended for very large alignments)" );
         System.out.println( "   -" + OUTPUT_FORMAT_PHYLIP_OPTION
                 + "             to write output alignments in phylip format instead of fasta" );
         System.out.println( "   -" + OUTPUT_REMOVED_SEQS_OPTION + "=<file>     to output the removed sequences" );
+        System.out.println( "   -" + MIN_LENGTH_OPTION
+                + "=<integer>  minimal effecive sequence length (for deleting of shorter sequences)" );
+        System.out.println( "   -" + GAP_RATIO_LENGTH_OPTION
+                + "=<decimal>  maximal allowed gap ratio per column (for deleting of columms) (0.0-1.0)" );
+        System.out.println( "   -" + PERFORM_PHYLOGENETIC_INFERENCE
+                + "             to calculate a simple phylogenetic tree (Kimura distances, NJ)" );
         System.out.println();
         System.out.println();
         System.out.println();