clean up XMFA parsing, allow spaces in parsing (ende++, from IRC)
[bioperl-live.git] / t / RestrictionAnalysis.t
blobcfee2d5df614709f6f87f2b412d83d598299d05a
1 # -*-Perl-*- Test Harness script for Bioperl
2 # $Id$
5 use strict;
7 BEGIN {
8     use lib 't/lib';
9     use BioperlTest;
10     
11     test_begin(-tests => 177);
12         
13     use_ok('Bio::Restriction::Enzyme');
14     use_ok('Bio::Restriction::Enzyme::MultiCut');
15     use_ok('Bio::Restriction::Enzyme::MultiSite');
16     use_ok('Bio::Restriction::EnzymeCollection');
17     use_ok('Bio::Restriction::Analysis');
18     use_ok('Bio::SeqIO');
22 # Bio::Restriction::Enzyme
25 my ($re, $seq, $iso, %meth, $microbe, $source, @vendors, @refs, $name);
26 ok $re=Bio::Restriction::Enzyme->new(-enzyme=>'EcoRI', -site=>'G^AATTC');
27 isa_ok($re, 'Bio::Restriction::EnzymeI');
28 is $re->cut, 1;
29 ok ! $re->cut(0);
30 is $re->complementary_cut, 6;
31 ok $re->cut(1);
33 is $re->complementary_cut,5;
34 is $re->site,'G^AATTC';
35 ok $seq = $re->seq;
36 isa_ok($seq, 'Bio::PrimarySeqI');
37 is $seq->seq, 'GAATTC';
38 is $re->string,'GAATTC';
39 is $re->revcom, 'GAATTC';
40 is $re->recognition_length, 6;
41 is $re->cutter, 6;
42 is $re->palindromic, 1;
43 is $re->overhang, "5'";
44 is $re->overhang_seq, 'AATT';
45 is $re->is_ambiguous, 0;
47 ok $re->compatible_ends($re);
49 ok $re->isoschizomers('BamHI', 'AvaI'); # not really true :)
51 is my @isos=$re->isoschizomers, 2;
52 is $isos[0],'BamHI';
53 ok $re->purge_isoschizomers;
54 is scalar($re->isoschizomers), 0;
55 ok $re->methylation_sites(2,5); # not really true :)
56 ok %meth = $re->methylation_sites;
57 is $meth{2}, 5;
58 ok $re->purge_methylation_sites;
59 is scalar($re->methylation_sites), 0;
62 ok $re->microbe('E. coli');
63 ok $microbe = $re->microbe;
64 is $microbe, "E. coli";
65 ok $re->source("Rob"); # not really true :)
67 ok $source = $re->source;
68 is $source, "Rob";
70 ok !$re->vendor;
71 ok $re->vendors('NEB'); # my favorite
72 ok @vendors = $re->vendors;
73 is $vendors[0], "NEB";
74 $re->purge_vendors;
75 is scalar($re->vendors), 0;
77 ok $re->references('Rob et al');
78 ok @refs = $re->references;
79 is $refs[0], "Rob et al";
80 $re->purge_references;
81 is scalar($re->references), 0;
83 ok $re->name('BamHI');
84 ok $name = $re->name;
85 is $name, "BamHI";
87 $re->verbose(2);
89 eval {$re->is_prototype};
90 ok($@);
91 like($@, qr/Can't unequivocally assign prototype based on input format alone/, 'bug 2179');
92 $re->verbose(2);
94 is $re->is_prototype(0), 0;
95 is $re->is_prototype, 0;
96 is $re->is_prototype(1), 1;
97 is $re->is_prototype, 1;
99 is $re->prototype_name, $re->name;
100 ok ! $re->is_prototype(0);
101 is $re->prototype_name('XxxI'), 'XxxI';
102 is $re->prototype_name, 'XxxI';
105 is $re->cutter, 6;
106 ok $re->seq->seq('RCATGY');
107 is $re->cutter, 5;
109 ok my $re2 = $re->clone;
110 isnt $re, $re2;
111 is $re->name, $re2->name;
113 ok $re = Bio::Restriction::Enzyme->new(-enzyme=>'AciI', 
114                                                                                 -site=>'C^CGC');
115 is $re->palindromic, 0;
116 is $re->is_palindromic, 0;
119 # Bio::Restriction::Enzyme::MultiSite
122 ok $re=Bio::Restriction::Enzyme::MultiSite->new(-enzyme=>'TaqII',
123                                               -site=>'GACCGA',
124                                               -cut=>17,
125                                               -complementary_cut=>15
126                                              );
127 ok $re2=Bio::Restriction::Enzyme::MultiSite->new(-enzyme=>'TaqII',
128                                                 -site=>'CACCCA',
129                                                 -cut=>17,
130                                                 -complementary_cut=>15
131                                                );
132 isa_ok( $re, 'Bio::Restriction::EnzymeI');
133 isa_ok( $re2, 'Bio::Restriction::EnzymeI');
134 ok $re->others($re2);
135 ok $re2->others($re);
137 is $re->others, 1;
138 is $re2->others, 1;
140 ok my $re3 = $re->clone;
141 isnt $re, $re3;
142 is $re->name , $re3->name; # wouldn't this be a circular reference???
143 #print Dumper $re, $re3;exit;
146 # Bio::Restriction::Enzyme::MultiCut
148 #Hin4I has four cut sites [(8/13)GAYNNNNNVTC(13/8)],
150 ok $re = Bio::Restriction::Enzyme::MultiCut->new(-enzyme=>'Hin4I',
151                                               -site=>'GAYNNNNNVTC',
152                                               -cut=>-8,
153                                               -complementary_cut=>-13
154                                              );
155 ok $re2 = Bio::Restriction::Enzyme::MultiCut->new(-enzyme=>'Hin4I',
156                                                -site=>'GAYNNNNNVTC',
157                                                -cut=>13,
158                                                -complementary_cut=>8
159                                               );
160 isa_ok($re, 'Bio::Restriction::EnzymeI');
161 isa_ok($re2, 'Bio::Restriction::EnzymeI');
162 ok $re->others($re2);
163 ok $re2->others($re);
165 ok $re3 = $re->clone;
166 isnt $re, $re3;
167 is $re->name, $re3->name;
168 #print Dumper $re, $re3;exit;
171 # Bio::Restriction::EnzymeCollection
174 my ($collection, $enz, $new_set);
176 ok $collection = Bio::Restriction::EnzymeCollection->new(-empty=>1);
177 is $collection->each_enzyme, 0;
178 # default set
179 $collection = Bio::Restriction::EnzymeCollection->new;
180 is $collection->each_enzyme, 532;
181 is $collection->each_enzyme, 532;
183 ok $enz = $collection->get_enzyme('AclI');
184 isa_ok($enz, 'Bio::Restriction::Enzyme');
185 is my @enzymes=$collection->available_list, 532;
187 ok $new_set = $collection->blunt_enzymes;
188 isa_ok($enz, 'Bio::Restriction::Enzyme');
189 is $new_set->each_enzyme, 114;
191 #map {print $_->name, ": ", $_->cutter, "\n"; } $collection->each_enzyme;
193 ok $new_set = $collection->cutters(8);
194 is $new_set->each_enzyme, 17;
196 ok $new_set=$collection->cutters(-start => 8, -end => 8);
197 is $new_set->each_enzyme, 17;
199 ok $new_set=$collection->cutters(-start => 6, -end => 8);
200 is $new_set->each_enzyme, 293;
202 ok $new_set=$collection->cutters(-start => 6, -end => 8,  -exclusive => 1);
203 is $new_set->each_enzyme, 10;
205 ok $new_set = $collection->cutters([4,8]);
206 is $new_set->each_enzyme, 129;
208 # bug 2128; enhancement request to pass array ref of sizes
211 # Restriction::Analysis
215 ok my $seqio=Bio::SeqIO->new(-file=>test_input_file('dna1.fa'),
216                          -format=>'fasta');
217 ok $seq=$seqio->next_seq;
219 ok my $ra = Bio::Restriction::Analysis->new(-seq=>$seq);
220 ok my $uniq = $ra->unique_cutters;
222 # test most objects
223 is $ra->unique_cutters->each_enzyme, 42, 'number of unique cutters';
224 is $ra->fragments('RsaI'), 2, 'number of RsaI fragments';
225 is $ra->max_cuts, 9, 'number of maximum cutters';
226 is $ra->zero_cutters->each_enzyme, 477, 'number of zero cutters';
227 is $ra->cutters->each_enzyme, 55, 'number of cutters';
228 is $ra->cutters(3)->each_enzyme, 8, 'number of 3x cutters';
229 is $ra->fragments('MseI'), 4, '4 MseI fragments';
230 is $ra->cuts_by_enzyme('MseI'), 3, '3 MseI cut sites';
232 #my $z = $ra->cutters(3);
233 #my $out=Bio::Restriction::IO->new;
234 #$out->write($z);
237 is $ra->fragments('PspEI'), 2, 'expected 2 PspEI fragments';
238 is $ra->cuts_by_enzyme('PspEI'), 1;
239 is $ra->cuts_by_enzyme('XxxI'), undef;
242 is my @ss = $ra->sizes('PspEI'), 2, 'expected 2 sizes for PspEI';
243 is $ss[0] + $ss[1], $seq->length;
245 is $ra->fragments('MwoI'), 1, 'not circular expected 1 fragments for MwoI as it doesnt cut';
247 # circularise the sequence, regenerate the cuts and test again
248 # note that there is one less fragment now!
249 ok $seq->is_circular(1);
251 # we need to regenerate all the cuts
252 ok $ra->cut;
254 is $ra->fragments('RsaI'), 1, 'number of RsaI fragments';
255 is $ra->fragments('MseI'), 3, '3 circular MseI fragments';
256 is $ra->cuts_by_enzyme('MseI'), 3, '3 circular MseI cut sites';
257 is $ra->fragments('AciI'), 1, 'number for AciI a non-palindromic enzyme';
259 is $ra->fragments('MwoI'), 1, '1 fragment for MwoI as it cuts across the circ point';
261 ok my @rb=($collection->get_enzyme("AluI"), $collection->get_enzyme("MseI"), $collection->get_enzyme("MaeIII"));
263 # test multiple digests
264 ok my $rbc=Bio::Restriction::EnzymeCollection->new(-empty=>1);
265 ok $rbc->enzymes(@rb);
266 ok $ra->cut('multiple', $rbc);
267 is $ra->fragments('multiple_digest'),7, '7 fragments in the multiple digest';
268 is my @pos=$ra->positions('multiple_digest'),7, '7 positions in the multiple digest';
269 is my @ssm = $ra->sizes('multiple_digest'),7, '7 sizes in the multiple digest';
270 my $check_len;
271 map {$check_len+=$_}@ssm;
272 is $check_len, $seq->length;
274 # now test the non-palindromic part
275 # HindI is a non palindromic enzyme that cuts 9 times
276 is $ra->positions('HindI'), 9, ' expected 9 cuts for HindI';
278 # now we need to test the fragment maps
279 # lets do this for HindI
280 is my @fm=$ra->fragment_maps('HindI'), 9, 'expect 9 fragment maps for HindI';
281 foreach my $fm (@fm) {
282  is exists $fm->{'seq'}, 1, "sequence for ".$fm->{'seq'};
283  is exists $fm->{'start'}, 1, "start at ".$fm->{'start'};
284  is exists $fm->{'end'}, 1, "end at ".$fm->{'end'};
287 # bug 2139
289 eval {$re = Bio::Restriction::Enzyme->new(
290         -name    => 'Invalid',
291         -site    => 'G^AATTE' );};
293 ok $@;
294 like($@, qr(Unrecognized characters in site), 'bug 2139');