File Coverage

blib/lib/Bio/Tools/Run/Cap3.pm
Criterion Covered Total %
statement 17 45 37.7
branch 1 12 8.3
condition 0 8 0.0
subroutine 4 5 80.0
pod 1 1 100.0
total 23 71 32.3


line stmt bran cond sub pod time code
1             # You may distribute this module under the same terms as perl itself
2             #
3             # POD documentation - main docs before the code
4              
5             =head1 NAME
6              
7             Bio::Tools::Run::Cap3 - wrapper for CAP3
8              
9             =head1 SYNOPSIS
10              
11             use Bio::Tools::Run::Cap3;
12             # Run Cap3 using an input FASTA file
13             my $factory = Bio::Tools::Run::Cap3->new( -clipping_range => 150 );
14             my $asm_obj = $factory->run($fasta_file, $qual_file);
15             # An assembly object is returned by default
16             for my $contig ($assembly->all_contigs) {
17             ... do something ...
18             }
19              
20             # Read some sequences
21             use Bio::SeqIO;
22             my $sio = Bio::SeqIO->new(-file => $fasta_file, -format => 'fasta');
23             my @seqs;
24             while (my $seq = $sio->next_seq()) {
25             push @seqs,$seq;
26             }
27              
28             # Run Cap3 using input sequence objects and returning an assembly file
29             my $asm_file = 'results.ace';
30             $factory->out_type($asm_file);
31             $factory->run(\@seqs);
32              
33              
34             =head1 DESCRIPTION
35              
36             Wrapper module for CAP3 program
37              
38             =head1 FEEDBACK
39              
40             =head2 Mailing Lists
41              
42             User feedback is an integral part of the evolution of this and other
43             Bioperl modules. Send your comments and suggestions preferably to one
44             of the Bioperl mailing lists. Your participation is much appreciated.
45              
46             bioperl-l@bioperl.org - General discussion
47             http://bioperl.org/wiki/Mailing_lists - About the mailing lists
48              
49             =head2 Support
50              
51             Please direct usage questions or support issues to the mailing list:
52              
53             I
54              
55             rather than to the module maintainer directly. Many experienced and
56             reponsive experts will be able look at the problem and quickly
57             address it. Please include a thorough description of the problem
58             with code and data examples if at all possible.
59              
60             =head2 Reporting Bugs
61              
62             Report bugs to the Bioperl bug tracking system to help us keep track
63             the bugs and their resolution. Bug reports can be submitted via the
64             web:
65              
66             http://redmine.open-bio.org/projects/bioperl/
67              
68             =head1 AUTHORS
69              
70             Marc Logghe
71              
72             =head1 APPENDIX
73              
74             The rest of the documentation details each of the object
75             methods. Internal methods are usually preceded with a _
76              
77             =cut
78              
79             package Bio::Tools::Run::Cap3;
80              
81 1     1   110009 use strict;
  1         1  
  1         24  
82 1     1   415 use File::Copy;
  1         1581  
  1         45  
83              
84 1     1   4 use base qw(Bio::Root::Root Bio::Tools::Run::AssemblerBase);
  1         1  
  1         459  
85              
86             our $program_name = 'cap3';
87             our @program_params = (qw( band_expansion_size differences_quality_cutoff
88             clipping_quality_cutoff max_qscore_sum extra_nof_differences max_gap_length
89             gap_penalty_factor max_overhang_percent match_score_factor mismatch_score_factor
90             overlap_length_cutoff overlap_identity_cutoff reverse_orientation_value
91             overlap_score_cutoff max_word_occurrences min_correction_constraints
92             min_linking_constraints clipping_info_file output_prefix_string clipping_range
93             min_clip_good_reads ));
94             our @program_switches;
95             our %param_translation = (
96             'band_expansion_size' => 'a',
97             'differences_quality_cutoff' => 'b',
98             'clipping_quality_cutoff' => 'c',
99             'max_qscore_sum' => 'd',
100             'extra_nof_differences' => 'e',
101             'max_gap_length' => 'f',
102             'gap_penalty_factor' => 'g',
103             'max_overhang_percent' => 'h',
104             'match_score_factor' => 'm',
105             'mismatch_score_factor' => 'n',
106             'overlap_length_cutoff' => 'o',
107             'overlap_identity_cutoff' => 'p',
108             'reverse_orientation_value' => 'r',
109             'overlap_score_cutoff' => 's',
110             'max_word_occurrences' => 't',
111             'min_correction_constraints' => 'u',
112             'min_linking_constraints' => 'v',
113             'clipping_info_file' => 'w',
114             'output_prefix_string' => 'x',
115             'clipping_range' => 'y',
116             'min_clip_good_reads' => 'z'
117             );
118             our $qual_param;
119             our $use_dash = 1;
120             our $join = ' ';
121             our $asm_format = 'ace';
122              
123             =head2 new
124              
125             Title : new
126             Usage : $factory->new(
127             -overlap_length_cutoff => 35,
128             -overlap_identity_cutoff => 98 # %
129             }
130             Function: Create a new Cap3 factory
131             Returns : A Bio::Tools::Run::Cap3 object
132             Args : Cap3 options available in this module:
133             band_expansion_size specify band expansion size N > 10 (20)
134             differences_quality_cutoff specify base quality cutoff for differences N > 15 (20)
135             clipping_quality_cutoff specify base quality cutoff for clipping N > 5 (12)
136             max_qscore_sum specify max qscore sum at differences N > 20 (200)
137             extra_nof_differences specify clearance between no. of diff N > 10 (30)
138             max_gap_length specify max gap length in any overlap N > 1 (20)
139             gap_penalty_factor specify gap penalty factor N > 0 (6)
140             max_overhang_percent specify max overhang percent length N > 2 (20)
141             match_score_factor specify match score factor N > 0 (2)
142             mismatch_score_factor specify mismatch score factor N < 0 (-5)
143             overlap_length_cutoff / minimum_overlap_length
144             specify overlap length cutoff > 20 (40)
145             overlap_identity_cutoff / minimum_overlap_similarity
146             specify overlap percent identity cutoff N > 65 (80)
147             reverse_orientation_value specify reverse orientation value N >= 0 (1)
148             overlap_score_cutoff specify overlap similarity score cutoff N > 400 (900)
149             max_word_occurrences specify max number of word matches N > 30 (300)
150             min_correction_constraints specify min number of constraints for correction N > 0 (3)
151             min_linking_constraints specify min number of constraints for linking N > 0 (2)
152             clipping_info_file specify file name for clipping information (none)
153             output_prefix_string specify prefix string for output file names (cap)
154             clipping_range specify clipping range N > 5 (250)
155             min_clip_good_reads specify min no. of good reads at clip pos N > 0 (3)
156              
157             =cut
158              
159             sub new {
160 1     1 1 78 my ($class,@args) = @_;
161 1         11 my $self = $class->SUPER::new(@args);
162 1         20 $self->_set_program_options(\@args, \@program_params, \@program_switches, \%param_translation,
163             $qual_param, $use_dash, $join);
164 1         3 *minimum_overlap_length = \&overlap_length_cutoff;
165 1         2 *minimum_overlap_similarity = \&overlap_identity_cutoff;
166 1 50       4 $self->program_name($program_name) if not defined $self->program_name();
167 1         6 $self->_assembly_format($asm_format);
168 1         5 return $self;
169             }
170              
171              
172             =head2 out_type
173              
174             Title : out_type
175             Usage : $assembler->out_type('Bio::Assembly::ScaffoldI')
176             Function: Get/set the desired type of output
177             Returns : The type of results to return
178             Args : Desired type of results to return (optional):
179             'Bio::Assembly::IO' object
180             'Bio::Assembly::ScaffoldI' object (default)
181             The name of a file to save the results in
182              
183             =cut
184              
185              
186             =head2 run
187              
188             Title : run
189             Usage : $asm = $factory->run($fasta_file);
190             Function: Run CAP3
191             Returns : Assembly results (file, IO object or assembly object)
192             Args : - sequence input (FASTA file or sequence object arrayref)
193             - optional quality score input (QUAL file or quality score object
194             arrayref)
195             =cut
196              
197              
198             =head2 _run
199              
200             Title : _run
201             Usage : $factory->_run()
202             Function: Make a system call and run Cap3
203             Returns : An assembly file
204             Args : - FASTA file
205             - optional QUAL file
206              
207             =cut
208              
209             sub _run {
210 0     0     my ($self, $fasta_file, $qual_file) = @_;
211              
212             # Move quality file to proper place
213 0           my $tmp_qual_file = "$fasta_file.qual";
214 0 0 0       if ($qual_file && not $qual_file eq $tmp_qual_file) {
215 0           $tmp_qual_file = "$fasta_file.qual"; # by Cap3 convention
216 0 0 0       link ($qual_file, $tmp_qual_file) or copy ($qual_file, $tmp_qual_file) or
217             $self->throw("Could not copy file '$qual_file' to '$tmp_qual_file': $!");
218             }
219              
220             # Setup needed files and filehandles
221 0           my ($output_fh, $output_file) = $self->_prepare_output_file( );
222              
223             # Get program executable
224 0           my $exe = $self->executable;
225              
226             # Get command-line options
227 0           my $options = join ' ', @{$self->_translate_params()};
  0            
228              
229             # Usage: cap3 File_of_reads [options]
230 0           my $commandstring = "$exe $fasta_file $options";
231              
232 0 0         if ($self->verbose() >= 0) {
233 0           $self->debug( "$exe command = $commandstring\n" );
234             }
235              
236 0 0         open(CAP3, "$commandstring |") ||
237             $self->throw(sprintf("%s call crashed: %s %s\n", $self->program_name, $!, $commandstring));
238 0           local $/ = undef;
239             #my ($result) = ; # standard output of the program
240 0           ;
241 0           close CAP3;
242 0           close $output_fh;
243              
244             # Result files
245 0   0       my $prefix = $self->output_prefix_string() || 'cap';
246 0           my $ace_file = "$fasta_file.$prefix.ace";
247 0           my $contigs_file = "$fasta_file.$prefix.contigs";
248 0           $qual_file = "$fasta_file.$prefix.contigs.links";
249 0           my $links_file = "$fasta_file.$prefix.contigs.qual";
250 0           my $info_file = "$fasta_file.$prefix.info";
251 0           my $singlet_file = "$fasta_file.$prefix.singlets";
252              
253             # Remove all files except for the ACE file
254 0           for my $file ($contigs_file, $qual_file, $links_file, $info_file, $singlet_file, $tmp_qual_file) {
255 0           unlink $file;
256             }
257              
258             # Move the ACE file to its final destination
259 0 0         move ($ace_file, $output_file) or $self->throw("Could not move file '$ace_file' to '$output_file': $!");
260              
261 0           return $output_file;
262             }
263              
264             1;