forester/ruby/evoruby/lib/evo/apps/msa_processor.rb

   1 #
   2 # = lib/evo/apps/msa_processor.rb - MsaProcessor class
   3 #
   4 # Copyright::  Copyright (C) 2006-2007 Christian M. Zmasek
   5 # License::    GNU Lesser General Public License (LGPL)
   6 #
   7 # $Id: msa_processor.rb,v 1.33 2010/12/13 19:00:10 cmzmasek Exp $
   8 #
   9
  10 require 'date'
  11 require 'set'
  12
  13 require 'lib/evo/util/constants'
  14 require 'lib/evo/util/util'
  15 require 'lib/evo/util/command_line_arguments'
  16 require 'lib/evo/msa/msa_factory'
  17 require 'lib/evo/io/msa_io'
  18 require 'lib/evo/io/writer/phylip_sequential_writer'
  19 require 'lib/evo/io/writer/nexus_writer'
  20 require 'lib/evo/io/writer/fasta_writer'
  21 require 'lib/evo/io/parser/fasta_parser'
  22 require 'lib/evo/io/parser/general_msa_parser'
  23 require 'lib/evo/io/writer/msa_writer'
  24
  25 module Evoruby
  26
  27   class MsaProcessor
  28
  29     PRG_NAME       = "msa_pro"
  30     PRG_DATE       = "2010.03.19"
  31     PRG_DESC       = "processing of multiple sequence alignments"
  32     PRG_VERSION    = "1.05"
  33     COPYRIGHT      = "2008-2010 Christian M Zmasek"
  34     CONTACT        = "phylosoft@gmail.com"
  35     WWW            = "www.phylosoft.org"
  36
  37
  38     NAME_LENGTH_DEFAULT                = 10
  39     WIDTH_DEFAULT_FASTA                = 60
  40     INPUT_TYPE_OPTION                  = "i"
  41     OUTPUT_TYPE_OPTION                 = "o"
  42     MAXIMAL_NAME_LENGTH_OPTION         = "n"
  43     WIDTH_OPTION                       = "w"
  44     CLEAN_UP_SEQ_OPTION                = "c"
  45     REM_RED_OPTION                     = "rem_red"
  46     REMOVE_GAP_COLUMNS_OPTION          = "rgc"
  47     REMOVE_GAP_ONLY_COLUMNS            = "rgoc"
  48     REMOVE_COLUMNS_GAP_RATIO_OPTION    = "rr"
  49     REMOVE_ALL_GAP_CHARACTERS_OPTION   = "rg"
  50     REMOVE_ALL_SEQUENCES_LISTED_OPTION = "r"
  51     KEEP_ONLY_SEQUENCES_LISTED_OPTION  = "k"
  52
  53     KEEP_MATCHING_SEQUENCES_OPTION     = "mk"
  54     REMOVE_MATCHING_SEQUENCES_OPTION   = "mr"
  55
  56     TRIM_OPTION                        = "t"
  57     REMOVE_SEQS_GAP_RATIO_OPTION       = "rsgr"
  58     REMOVE_SEQS_NON_GAP_LENGTH_OPTION  = "rsl"
  59     SPLIT                              = "split"
  60     LOG_SUFFIX                         = "_msa_pro.log"
  61     HELP_OPTION_1                      = "help"
  62     HELP_OPTION_2                      = "h"
  63
  64
  65     def initialize()
  66       @input_format_set = false
  67       @output_format_set = false
  68       @fasta_input      = false
  69       @phylip_input     = true
  70       @name_length      = NAME_LENGTH_DEFAULT
  71       @name_length_set  = false
  72       @width            = WIDTH_DEFAULT_FASTA     # fasta only
  73       @pi_output        = true
  74       @fasta_output     = false
  75       @nexus_output     = false
  76       @clean            = false  # phylip only
  77       @rgc              = false
  78       @rgoc             = false
  79       @rg               = false  # fasta only
  80       @rem_red          = false
  81       @rgr              = -1
  82       @rsgr             = -1
  83       @rsl              = -1
  84       @remove_matching  = nil
  85       @keep_matching    = nil
  86
  87       @seqs_name_file   = nil
  88       @remove_seqs      = false
  89       @keep_seqs        = false
  90       @trim             = false
  91       @split            = -1
  92       @first            = -1
  93       @last             = -1
  94     end
  95
  96
  97     def run()
  98
  99       Util.print_program_information( PRG_NAME,
 100         PRG_VERSION,
 101         PRG_DESC,
 102         PRG_DATE,
 103         COPYRIGHT,
 104         CONTACT,
 105         WWW,
 106         STDOUT )
 107
 108       if ( ARGV == nil || ARGV.length < 1 )
 109         Util.print_message( PRG_NAME, "Illegal number of arguments" )
 110         print_help
 111         exit( -1 )
 112       end
 113
 114       begin
 115         cla = CommandLineArguments.new( ARGV )
 116       rescue ArgumentError => e
 117         Util.fatal_error( PRG_NAME, "Error: " + e.to_s, STDOUT )
 118       end
 119
 120       if ( cla.is_option_set?( HELP_OPTION_1 ) ||
 121            cla.is_option_set?( HELP_OPTION_2 ) )
 122         print_help
 123         exit( 0 )
 124       end
 125
 126       if ( cla.get_number_of_files != 2 || ARGV.length < 2 )
 127         Util.print_message( PRG_NAME, "Illegal number of arguments" )
 128         print_help
 129         exit( -1 )
 130       end
 131
 132       allowed_opts = Array.new
 133       allowed_opts.push( INPUT_TYPE_OPTION )
 134       allowed_opts.push( OUTPUT_TYPE_OPTION )
 135       allowed_opts.push( MAXIMAL_NAME_LENGTH_OPTION )
 136       allowed_opts.push( WIDTH_OPTION )
 137       allowed_opts.push( CLEAN_UP_SEQ_OPTION )
 138       allowed_opts.push( REMOVE_GAP_COLUMNS_OPTION )
 139       allowed_opts.push( REMOVE_GAP_ONLY_COLUMNS )
 140       allowed_opts.push( REMOVE_COLUMNS_GAP_RATIO_OPTION )
 141       allowed_opts.push( REMOVE_ALL_GAP_CHARACTERS_OPTION )
 142       allowed_opts.push( REMOVE_ALL_SEQUENCES_LISTED_OPTION )
 143       allowed_opts.push( KEEP_ONLY_SEQUENCES_LISTED_OPTION )
 144       allowed_opts.push( TRIM_OPTION )
 145       allowed_opts.push( REMOVE_SEQS_GAP_RATIO_OPTION )
 146       allowed_opts.push( REMOVE_SEQS_NON_GAP_LENGTH_OPTION )
 147       allowed_opts.push( SPLIT )
 148       allowed_opts.push( REM_RED_OPTION )
 149       allowed_opts.push( KEEP_MATCHING_SEQUENCES_OPTION )
 150       allowed_opts.push( REMOVE_MATCHING_SEQUENCES_OPTION )
 151
 152       disallowed = cla.validate_allowed_options_as_str( allowed_opts )
 153       if ( disallowed.length > 0 )
 154         Util.fatal_error( PRG_NAME,
 155           "unknown option(s): " + disallowed )
 156       end
 157
 158       input = cla.get_file_name( 0 )
 159       output = cla.get_file_name( 1 )
 160
 161       analyze_command_line( cla )
 162
 163       begin
 164         Util.check_file_for_readability( input )
 165       rescue ArgumentError => e
 166         Util.fatal_error( PRG_NAME, "error: " + e.to_s )
 167       end
 168
 169       begin
 170         Util.check_file_for_writability( output )
 171       rescue ArgumentError => e
 172         Util.fatal_error( PRG_NAME, "error: " + e.to_s )
 173       end
 174
 175       if ( @rg )
 176         set_pi_output( false )
 177         set_fasta_output( true )
 178         set_nexus_output( false )
 179       end
 180
 181       if ( !@input_format_set )
 182         fasta_like = false
 183         begin
 184           fasta_like = Util.looks_like_fasta?( input )
 185         rescue ArgumentError => e
 186           Util.fatal_error( PRG_NAME, "error: " + e.to_s )
 187         end
 188         @fasta_input = fasta_like
 189         @phylip_input = !fasta_like
 190         if ( !@output_format_set )
 191           @fasta_output = fasta_like
 192           @pi_output = !fasta_like
 193           @nexus_output = false
 194         end
 195       end
 196
 197       ld = Constants::LINE_DELIMITER
 198       log = PRG_NAME + " " + PRG_VERSION + " [" + PRG_DATE + "]" + " LOG" + ld
 199       now = DateTime.now
 200       log << "Date/time: " + now.to_s + ld
 201
 202       puts()
 203       puts( "Input alignment  : " + input )
 204       log << "Input alignment  : " + input + ld
 205       puts( "Output alignment : " + output )
 206       log << "Output alignment : " + output + ld
 207       if ( @phylip_input )
 208         puts( "Input is         : Phylip, or something like it" )
 209         log << "Input is         : Phylip, or something like it" + ld
 210       elsif ( @fasta_input )
 211         puts( "Input is         : Fasta" )
 212         log << "Input is         : Fasta" + ld
 213       end
 214       if( @rgr >= 0 )
 215         puts( "Max col gap ratio: " + @rgr.to_s )
 216         log << "Max col gap ratio: " + @rgr.to_s + ld
 217       elsif ( @rgc )
 218         puts( "Remove gap colums" )
 219         log << "Remove gap colums" + ld
 220       elsif( @rgoc )
 221         puts( "Remove gap only colums" )
 222         log << "Remove gap only colums" + ld
 223       end
 224       if ( @clean )
 225         puts( "Clean up         : true" )
 226         log << "Clean up         : true" + ld
 227       end
 228
 229       if ( @pi_output )
 230         puts( "Output is        : Phylip interleaved" )
 231         log << "Output is        : Phylip interleaved" + ld
 232       elsif ( @fasta_output )
 233         puts( "Output is        : Fasta" )
 234         log << "Output is        : Fasta" + ld
 235         if ( @width )
 236           puts( "Width            : " + @width.to_s )
 237           log << "Width            : " + @width.to_s + ld
 238         end
 239         if ( @rg )
 240           puts( "Remove all gap characters (alignment is destroyed)" )
 241           log << "Remove all gap characters (alignment is destroyed)" + ld
 242         end
 243       elsif ( @nexus_output )
 244         puts( "Output is        : Nexus" )
 245         log << "Output is        : Nexus" + ld
 246       end
 247       if ( @name_length_set || !@fasta_output )
 248         puts( "Max name length  : " + @name_length.to_s )
 249         log << "Max name length  : " + @name_length.to_s + ld
 250       end
 251       if( @rsgr >= 0 )
 252         puts( "Remove sequences for which the gap ratio > " + @rsgr.to_s )
 253         log << "Remove sequences for which the gap ratio > " + @rsgr.to_s + ld
 254       end
 255       if( @rsl >= 0 )
 256         puts( "Remove sequences with less than "  + @rsl.to_s + " non-gap characters" )
 257         log << "Remove sequences with less than "  + @rsl.to_s + " non-gap characters" + ld
 258       end
 259       if ( @remove_seqs )
 260         puts( "Remove sequences listed in: " + @seqs_name_file )
 261         log << "Remove sequences listed in: " + @seqs_name_file + ld
 262       elsif ( @keep_seqs )
 263         puts( "Keep only sequences listed in: " + @seqs_name_file )
 264         log << "Keep only sequences listed in: " + @seqs_name_file + ld
 265       end
 266       if ( @trim )
 267         puts( "Keep only columns from: "+ @first.to_s + " to " + @last.to_s )
 268         log << "Keep only columns from: "+ @first.to_s + " to " + @last.to_s + ld
 269       end
 270       if ( @rem_red )
 271         puts( "Remove redundant sequences: true" )
 272         log << "Remove redundant sequences: true" + ld
 273       end
 274       if ( @split > 0 )
 275         puts( "Split            : " + @split.to_s )
 276         log << "Split            : " + @split.to_s + ld
 277       end
 278       puts()
 279
 280       f = MsaFactory.new()
 281
 282       msa = nil
 283
 284       begin
 285         if ( @phylip_input )
 286           msa = f.create_msa_from_file( input, GeneralMsaParser.new() )
 287         elsif ( @fasta_input )
 288           msa = f.create_msa_from_file( input, FastaParser.new() )
 289         end
 290       rescue Exception => e
 291         Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 292       end
 293
 294       if ( msa.is_aligned() )
 295         Util.print_message( PRG_NAME, "Length of original alignment : " + msa.get_length.to_s )
 296         log << "Length of original alignment : " + msa.get_length.to_s + ld
 297       else
 298         Util.print_message( PRG_NAME, "the input is not aligned" )
 299         log << "The input is not aligned" + ld
 300       end
 301
 302       all_names = Set.new()
 303       for i in 0 ... msa.get_number_of_seqs()
 304         current_name = msa.get_sequence( i ).get_name
 305         if all_names.include?( current_name )
 306           Util.print_warning_message( PRG_NAME, "sequence name [" + current_name + "] is not unique" )
 307         else
 308           all_names.add( current_name )
 309         end
 310       end
 311
 312       begin
 313
 314         if ( @remove_seqs || @keep_seqs )
 315           names = Util.file2array( @seqs_name_file, true )
 316           if ( names == nil ||  names.length() < 1 )
 317             error_msg = "file \"" + @seqs_name_file.to_s + "\" appears empty"
 318             Util.fatal_error( PRG_NAME, error_msg )
 319           end
 320
 321           if ( @remove_seqs )
 322             c = 0
 323             for i in 0 ... names.length()
 324               to_delete = msa.find_by_name( names[ i ], true, false )
 325               if ( to_delete.length() < 1 )
 326                 error_msg = "sequence name \"" + names[ i ] + "\" not found"
 327                 Util.fatal_error( PRG_NAME, error_msg )
 328               elsif ( to_delete.length() > 1 )
 329                 error_msg = "sequence name \"" + names[ i ] + "\" is not unique"
 330                 Util.fatal_error( PRG_NAME, error_msg )
 331               else
 332                 msa.remove_sequence!( to_delete[ 0 ] )
 333                 c += 1
 334               end
 335             end
 336             Util.print_message( PRG_NAME, "Removed " + c.to_s + " sequences" )
 337             log <<  "Removed " + c.to_s + " sequences" + ld
 338           elsif ( @keep_seqs )
 339             msa_new = Msa.new()
 340             r = 0
 341             k = 0
 342             for j in 0 ... msa.get_number_of_seqs()
 343               if ( names.include?( msa.get_sequence( j ).get_name() ) )
 344                 msa_new.add_sequence( msa.get_sequence( j ) )
 345                 k += 1
 346               else
 347                 r += 1
 348               end
 349             end
 350             msa = msa_new
 351             Util.print_message( PRG_NAME, "Kept    " + k.to_s + " sequences" )
 352             log << "Kept    " + k.to_s + " sequences" + ld
 353             Util.print_message( PRG_NAME, "Removed " + r.to_s + " sequences" )
 354             log << "removed " + r.to_s + " sequences" + ld
 355           end
 356         end
 357
 358         if ( @trim )
 359           msa.trim!( @first, @last )
 360         end
 361         if( @rgr >= 0 )
 362           msa.remove_gap_columns_w_gap_ratio!( @rgr )
 363         elsif ( @rgc )
 364           msa.remove_gap_columns!()
 365         elsif( @rgoc )
 366           msa.remove_gap_only_columns!()
 367         end
 368         if( @rsgr >= 0 )
 369           n = msa.get_number_of_seqs()
 370           removed = msa.remove_sequences_by_gap_ratio!( @rsgr )
 371           k = msa.get_number_of_seqs()
 372           r = n - k
 373           Util.print_message( PRG_NAME, "Kept    " + k.to_s + " sequences" )
 374           log << "Kept    " + k.to_s + " sequences" + ld
 375           Util.print_message( PRG_NAME, "Removed " + r.to_s + " sequences"  )
 376           log << "Removed " + r.to_s + " sequences:" + ld
 377           removed.each { | seq_name |
 378             log << "         " + seq_name  + ld
 379           }
 380         end
 381         if( @rsl >= 0 )
 382           n = msa.get_number_of_seqs()
 383           removed = msa.remove_sequences_by_non_gap_length!( @rsl )
 384           k = msa.get_number_of_seqs()
 385           r = n - k
 386           Util.print_message( PRG_NAME, "Kept    " + k.to_s + " sequences" )
 387           log << "Kept    " + k.to_s + " sequences" + ld
 388           Util.print_message( PRG_NAME, "Removed " + r.to_s + " sequences" )
 389           log << "Removed " + r.to_s + " sequences:" + ld
 390           removed.each { | seq_name |
 391             log << "         " + seq_name  + ld
 392           }
 393         end
 394         if ( @keep_matching )
 395           n = msa.get_number_of_seqs
 396           to_be_removed = Set.new
 397           for ii in 0 ...  n
 398             seq = msa.get_sequence( ii )
 399             if !seq.get_name.downcase.index( @keep_matching.downcase )
 400               to_be_removed.add( ii )
 401             end
 402           end
 403           to_be_removed_ary = to_be_removed.to_a.sort.reverse
 404           to_be_removed_ary.each { | index |
 405             msa.remove_sequence!( index )
 406           }
 407           # msa = sort( msa )
 408         end
 409         if ( @remove_matching )
 410           n = msa.get_number_of_seqs
 411           to_be_removed = Set.new
 412           for iii in 0 ... n
 413
 414             seq = msa.get_sequence( iii )
 415
 416             if seq.get_name.downcase.index( @remove_matching.downcase )
 417               to_be_removed.add( iii )
 418             end
 419           end
 420           to_be_removed_ary = to_be_removed.to_a.sort.reverse
 421           to_be_removed_ary.each { | index |
 422             msa.remove_sequence!( index )
 423           }
 424           msa = sort( msa )
 425         end
 426
 427
 428
 429         if ( @split > 0 )
 430           begin
 431             msas = msa.split( @split, true )
 432             io = MsaIO.new()
 433             w = MsaWriter
 434             if ( @pi_output )
 435               w = PhylipSequentialWriter.new()
 436               w.clean( @clean )
 437               w.set_max_name_length( @name_length )
 438             elsif( @fasta_output )
 439               w = FastaWriter.new()
 440               w.set_line_width( @width )
 441               if ( @rg )
 442                 w.remove_gap_chars( true )
 443                 Util.print_warning_message( PRG_NAME, "removing gap character, the output is likely to become unaligned" )
 444                 log << "removing gap character, the output is likely to become unaligned" + ld
 445               end
 446               w.clean( @clean )
 447               if ( @name_length_set )
 448                 w.set_max_name_length( @name_length )
 449               end
 450             elsif( @nexus_output )
 451               w = NexusWriter.new()
 452               w.clean( @clean )
 453               w.set_max_name_length( @name_length )
 454             end
 455             i = 0
 456             for m in msas
 457               i = i + 1
 458               io.write_to_file( m, output + "_" + i.to_s, w )
 459             end
 460             Util.print_message( PRG_NAME, "wrote " + msas.length.to_s + " files"  )
 461             log << "wrote " + msas.length.to_s + " files" + ld
 462           rescue Exception => e
 463             Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 464           end
 465
 466         end
 467       rescue Exception => e
 468         Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 469       end
 470
 471       if ( @split <= 0 )
 472
 473         unless ( @rg )
 474           if ( msa.is_aligned() )
 475             Util.print_message( PRG_NAME, "length of processed alignment: " + msa.get_length.to_s )
 476             log <<  "length of processed alignment: " + msa.get_length.to_s + ld
 477           else
 478             Util.print_warning_message( PRG_NAME, "output is not aligned" )
 479             log << "output is not aligned" + ld
 480           end
 481         end
 482
 483         if @rem_red
 484           removed = msa.remove_redundant_sequences!( true, true )
 485           if removed.size > 0
 486             Util.print_message( PRG_NAME, "going to ignore the following " + removed.size.to_s + " redundant sequences:" )
 487             log << "going to ignore the following " + removed.size.to_s + " redundant sequences:" + ld
 488             removed.each { | seq_name |
 489               puts seq_name
 490               log << seq_name + ld
 491             }
 492             Util.print_message( PRG_NAME, "will store " + msa.get_number_of_seqs.to_s + " non-redundant sequences" )
 493             log << "will store " + msa.get_number_of_seqs.to_s + " non-redundant sequences" + ld
 494           end
 495         end
 496
 497         io = MsaIO.new()
 498
 499         w = MsaWriter
 500
 501         if ( @pi_output )
 502           w = PhylipSequentialWriter.new()
 503           w.clean( @clean )
 504           w.set_max_name_length( @name_length )
 505         elsif( @fasta_output )
 506           w = FastaWriter.new()
 507           w.set_line_width( @width )
 508           if ( @rg )
 509             w.remove_gap_chars( true )
 510             Util.print_warning_message( PRG_NAME, "removing gap characters, the output is likely to become unaligned"  )
 511             log << "removing gap character, the output is likely to become unaligned" + ld
 512           end
 513           w.clean( @clean )
 514           if ( @name_length_set )
 515             w.set_max_name_length( @name_length )
 516           end
 517         elsif( @nexus_output )
 518           w = NexusWriter.new()
 519           w.clean( @clean )
 520           w.set_max_name_length( @name_length )
 521         end
 522
 523
 524         begin
 525           io.write_to_file( msa, output, w )
 526         rescue Exception => e
 527           Util.fatal_error( PRG_NAME, "error: " + e.to_s )
 528         end
 529
 530         begin
 531           f = File.open( output + LOG_SUFFIX, 'a' )
 532           f.print( log )
 533           f.close
 534         rescue Exception => e
 535           Util.fatal_error( PRG_NAME, "error: " + e.to_s )
 536         end
 537
 538
 539       end
 540       Util.print_message( PRG_NAME, "OK" )
 541       puts
 542     end
 543
 544
 545     private
 546
 547     def sort( msa )
 548       names = Set.new
 549       for i in 0 ... msa.get_number_of_seqs
 550         name = msa.get_sequence( i ).get_name
 551         names.add( name )
 552       end
 553       sorted_ary = names.to_a.sort
 554       new_msa = Msa.new
 555       sorted_ary.each { | seq_name |
 556         seq = msa.get_sequence( msa.find_by_name( seq_name, true, false )[ 0 ] )
 557         new_msa.add_sequence( seq )
 558       }
 559       new_msa
 560     end
 561
 562     def set_fasta_input( fi = true )
 563       @fasta_input = fi
 564       @input_format_set = true
 565     end
 566     def set_phylip_input( pi = true )
 567       @phylip_input = pi
 568       @input_format_set = true
 569     end
 570     def set_name_length( i )
 571       @name_length = i
 572       @name_length_set = true
 573     end
 574     def set_width( i )
 575       @width = i
 576     end
 577     def set_fasta_output( fo = true )
 578       @fasta_output = fo
 579       @output_format_set = true
 580     end
 581     def set_pi_output( pso = true )
 582       @pi_output = pso
 583       @output_format_set = true
 584     end
 585     def set_nexus_output( nexus = true )
 586       @nexus_output = nexus
 587       @output_format_set = true
 588     end
 589     def set_clean( c = true )
 590       @clean = c
 591     end
 592     def set_remove_gap_columns( rgc = true )
 593       @rgc = rgc
 594     end
 595     def set_remove_gap_only_columns( rgoc = true )
 596       @rgoc = rgoc
 597     end
 598     def set_remove_gaps( rg = true )
 599       @rg = rg
 600     end
 601     def set_remove_gap_ratio( rgr )
 602       @rgr = rgr
 603     end
 604     def set_remove_seqs_gap_ratio( rsgr )
 605       @rsgr = rsgr
 606     end
 607     def set_remove_seqs_min_non_gap_length( rsl )
 608       @rsl = rsl
 609     end
 610     def set_remove_seqs( file )
 611       @seqs_name_file = file
 612       @remove_seqs    = true
 613       @keep_seqs      = false
 614     end
 615     def set_keep_seqs( file )
 616       @seqs_name_file = file
 617       @keep_seqs      = true
 618       @remove_seqs    = false
 619     end
 620     def set_trim( first, last )
 621       @trim            = true
 622       @first           = first
 623       @last            = last
 624     end
 625     def set_remove_matching( remove )
 626       @remove_matching  = remove
 627     end
 628     def set_keep_matching( keep )
 629       @keep_matching = keep
 630     end
 631     def set_rem_red( rr )
 632       @rem_red = rr
 633     end
 634
 635
 636
 637     def set_split( s )
 638       if ( s > 0 )
 639         @split            = s
 640         @clean            = false  # phylip only
 641         @rgc              = false
 642         @rgoc             = false
 643         @rg               = false  # fasta only
 644         @rgr              = -1
 645         @rsgr             = -1
 646         @rsl              = -1
 647         @seqs_name_file   = nil
 648         @remove_seqs      = false
 649         @keep_seqs        = false
 650         @trim             = false
 651         @first            = -1
 652         @last             = -1
 653       end
 654     end
 655     def analyze_command_line( cla )
 656       if ( cla.is_option_set?( INPUT_TYPE_OPTION ) )
 657         begin
 658           type = cla.get_option_value( INPUT_TYPE_OPTION )
 659           if ( type == "p" )
 660             set_phylip_input( true )
 661             set_fasta_input( false )
 662           elsif ( type == "f" )
 663             set_fasta_input( true )
 664             set_phylip_input( false )
 665           end
 666         rescue ArgumentError => e
 667           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 668         end
 669       end
 670       if ( cla.is_option_set?( OUTPUT_TYPE_OPTION ) )
 671         begin
 672           type = cla.get_option_value( OUTPUT_TYPE_OPTION )
 673           if ( type == "p" )
 674             set_pi_output( true )
 675             set_fasta_output( false )
 676             set_nexus_output( false )
 677           elsif ( type == "f" )
 678             set_pi_output( false )
 679             set_fasta_output( true )
 680             set_nexus_output( false )
 681           elsif ( type == "n" )
 682             set_pi_output( false )
 683             set_fasta_output( false )
 684             set_nexus_output( true )
 685           end
 686         rescue ArgumentError => e
 687           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 688         end
 689       end
 690       if ( cla.is_option_set?( MAXIMAL_NAME_LENGTH_OPTION ) )
 691         begin
 692           l = cla.get_option_value_as_int( MAXIMAL_NAME_LENGTH_OPTION )
 693           set_name_length( l )
 694         rescue ArgumentError => e
 695           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 696         end
 697       end
 698       if ( cla.is_option_set?( WIDTH_OPTION ) )
 699         begin
 700           w = cla.get_option_value_as_int( WIDTH_OPTION )
 701           set_width( w )
 702         rescue ArgumentError => e
 703           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 704         end
 705       end
 706       if ( cla.is_option_set?( CLEAN_UP_SEQ_OPTION ) )
 707         set_clean( true )
 708       end
 709       if ( cla.is_option_set?( REMOVE_GAP_COLUMNS_OPTION ) )
 710         set_remove_gap_columns( true )
 711       end
 712       if ( cla.is_option_set?( REM_RED_OPTION ) )
 713         set_rem_red( true )
 714       end
 715       if ( cla.is_option_set?( REMOVE_GAP_ONLY_COLUMNS ) )
 716         set_remove_gap_only_columns( true )
 717       end
 718       if ( cla.is_option_set?( REMOVE_ALL_GAP_CHARACTERS_OPTION ) )
 719         set_remove_gaps( true )
 720       end
 721       if ( cla.is_option_set?( REMOVE_COLUMNS_GAP_RATIO_OPTION ) )
 722         begin
 723           f = cla.get_option_value_as_float( REMOVE_COLUMNS_GAP_RATIO_OPTION )
 724           set_remove_gap_ratio( f )
 725         rescue ArgumentError => e
 726           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 727         end
 728       end
 729       if ( cla.is_option_set?( REMOVE_ALL_SEQUENCES_LISTED_OPTION ) )
 730         begin
 731           s = cla.get_option_value( REMOVE_ALL_SEQUENCES_LISTED_OPTION )
 732           set_remove_seqs( s )
 733         rescue ArgumentError => e
 734           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 735         end
 736       end
 737       if ( cla.is_option_set?( KEEP_ONLY_SEQUENCES_LISTED_OPTION ) )
 738         begin
 739           s = cla.get_option_value( KEEP_ONLY_SEQUENCES_LISTED_OPTION )
 740           set_keep_seqs( s )
 741         rescue ArgumentError => e
 742           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 743         end
 744       end
 745       if ( cla.is_option_set?( TRIM_OPTION ) )
 746         begin
 747           s = cla.get_option_value( TRIM_OPTION )
 748           if ( s =~ /(\d+)-(\d+)/ )
 749             set_trim( $1.to_i(), $2.to_i() )
 750           else
 751             puts( "illegal argument" )
 752             print_help
 753             exit( -1 )
 754           end
 755         rescue ArgumentError => e
 756           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 757         end
 758       end
 759       if ( cla.is_option_set?( REMOVE_SEQS_GAP_RATIO_OPTION ) )
 760         begin
 761           f = cla.get_option_value_as_float( REMOVE_SEQS_GAP_RATIO_OPTION )
 762           set_remove_seqs_gap_ratio( f )
 763         rescue ArgumentError => e
 764           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 765         end
 766       end
 767       if ( cla.is_option_set?( REMOVE_SEQS_NON_GAP_LENGTH_OPTION ) )
 768         begin
 769           f = cla.get_option_value_as_int( REMOVE_SEQS_NON_GAP_LENGTH_OPTION )
 770           set_remove_seqs_min_non_gap_length( f )
 771         rescue ArgumentError => e
 772           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 773         end
 774       end
 775       if ( cla.is_option_set?( SPLIT ) )
 776         begin
 777           s = cla.get_option_value_as_int( SPLIT )
 778           set_split( s )
 779         rescue ArgumentError => e
 780           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 781         end
 782
 783       end
 784       if ( cla.is_option_set?( REMOVE_MATCHING_SEQUENCES_OPTION ) )
 785         begin
 786           s = cla.get_option_value( REMOVE_MATCHING_SEQUENCES_OPTION )
 787           set_remove_matching( s )
 788         rescue ArgumentError => e
 789           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 790         end
 791       end
 792       if ( cla.is_option_set?( KEEP_MATCHING_SEQUENCES_OPTION ) )
 793         begin
 794           s = cla.get_option_value( KEEP_MATCHING_SEQUENCES_OPTION )
 795           set_keep_matching( s )
 796         rescue ArgumentError => e
 797           Util.fatal_error( PRG_NAME, "error: " + e.to_s, STDOUT )
 798         end
 799       end
 800
 801
 802     end
 803
 804     def print_help()
 805       puts()
 806       puts( "Usage:" )
 807       puts()
 808       puts( "  " + PRG_NAME + ".rb [options] <input alignment> <output>" )
 809       puts()
 810       puts( "  options: -" + INPUT_TYPE_OPTION + "=<input type>: f for fasta, p for phylip selex type" )
 811       puts( "           -" + OUTPUT_TYPE_OPTION + "=<output type>: f for fasta, n for nexus, p for phylip sequential (default)" )
 812       puts( "           -" + MAXIMAL_NAME_LENGTH_OPTION + "=<n>: n=maximal name length (default for phylip 10, for fasta: unlimited )" )
 813       puts( "           -" + WIDTH_OPTION + "=<n>: n=width (fasta output only, default is 60)" )
 814       puts( "           -" + CLEAN_UP_SEQ_OPTION + ": clean up sequences" )
 815       puts( "           -" + REMOVE_GAP_COLUMNS_OPTION + ": remove gap columns" )
 816       puts( "           -" + REMOVE_GAP_ONLY_COLUMNS + ": remove gap-only columns" )
 817       puts( "           -" + REMOVE_COLUMNS_GAP_RATIO_OPTION + "=<n>: remove columns for which ( seqs with gap / number of sequences > n )" )
 818       puts( "           -" + REMOVE_ALL_GAP_CHARACTERS_OPTION + ": remove all gap characters (destroys alignment, fasta output only)" )
 819       puts( "           -" + REMOVE_ALL_SEQUENCES_LISTED_OPTION + "=<file>: remove all sequences listed in file" )
 820       puts( "           -" + KEEP_ONLY_SEQUENCES_LISTED_OPTION + "=<file>: keep only sequences listed in file" )
 821       puts( "           -" + TRIM_OPTION + "=<first>-<last>: remove columns before first and after last" )
 822       puts( "           -" + REMOVE_SEQS_GAP_RATIO_OPTION + "=<n>: remove sequences for which the gap ratio > n (after column operations)" )
 823       puts( "           -" + REMOVE_SEQS_NON_GAP_LENGTH_OPTION + "=<n> remove sequences with less than n non-gap characters (after column operations)" )
 824       puts( "           -" + REMOVE_MATCHING_SEQUENCES_OPTION + "=<s> remove all sequences with names containing s" )
 825       puts( "           -" + KEEP_MATCHING_SEQUENCES_OPTION + "=<s> keep only sequences with names containing s" )
 826       puts( "           -" + SPLIT + "=<n> split a fasta file into n files of equal number of sequences (expect for " )
 827       puts( "            last one), cannot be used with other options" )
 828       puts( "           -" + REM_RED_OPTION + ": remove redundant sequences" )
 829       puts()
 830     end
 831
 832
 833
 834
 835
 836   end # class MsaProcessor
 837
 838
 839 end # module Evoruby