t/seq_quality.t

   1 # -*-Perl-*- Test Harness script for Bioperl
   2 # $Id$
   3
   4 use strict;
   5
   6 BEGIN {
   7         use lib 't/lib';
   8     use BioperlTest;
   9
  10     test_begin(-tests => 53);
  11
  12         use_ok('Bio::Seq::Quality');
  13 }
  14
  15 my $DEBUG = test_debug();
  16
  17 # create some random sequence object with no id
  18 my $seqobj_broken = Bio::Seq::Quality->new( -seq => "ATCGATCGA",
  19                                           );
  20
  21 my $seqobj;
  22 lives_ok {
  23         $seqobj = Bio::Seq::Quality->new( -seq => "ATCGATCGA",
  24                                      -id  => 'QualityFragment-12',
  25                                      -accession_number => 'X78121',
  26                                    );
  27 };
  28
  29
  30 # create some random quality object with the same number of qualities and the same identifiers
  31 my $string_quals = "10 20 30 40 50 40 30 20 10";
  32 my $qualobj;
  33 lives_ok {
  34         $qualobj = Bio::Seq::Quality->new( -qual => $string_quals,
  35                                       -id  => 'QualityFragment-12',
  36                                       -accession_number => 'X78121',
  37                                         );
  38 };
  39
  40 # check to see what happens when you construct the Quality object
  41 ok my $swq1 = Bio::Seq::Quality->new( -seq => "ATCGATCGA",
  42                                       -id  => 'QualityFragment-12',
  43                                       -accession_number => 'X78121',
  44                                       -qual     =>      $string_quals);
  45
  46
  47 print("Testing various weird constructors...\n") if $DEBUG;
  48 print("\ta) No ids, Sequence object, no quality...\n") if $DEBUG;
  49 # w for weird
  50 my $wswq1;
  51 lives_ok {
  52         $wswq1 = Bio::Seq::Quality->new( -seq  =>  "ATCGATCGA",
  53                                          -qual  =>      "");
  54 };
  55 print $@ if $DEBUG;
  56
  57
  58 print("\tb) No ids, no sequence, quality object...\n") if $DEBUG;
  59         # note that you must provide a alphabet for this one.
  60 $wswq1 = Bio::Seq::Quality->new( -seq => "",
  61                                         -qual => $string_quals,
  62                                         -alphabet => 'dna'
  63 );
  64 print("\tc) Absolutely nothing. (HAHAHAHA)...\n") if $DEBUG;
  65 lives_ok {
  66         $wswq1 = Bio::Seq::Quality->new( -seq => "",
  67                                                 -qual => "",
  68                                                 -alphabet => 'dna'
  69         );
  70 };
  71
  72
  73 print("\td) Absolutely nothing but an ID\n") if $DEBUG;
  74 lives_ok {
  75     $wswq1 = Bio::Seq::Quality->new( -seq => "",
  76                                             -qual => "",
  77                                             -alphabet => 'dna',
  78                                             -id => 'an object with no sequence and no quality but with an id'
  79         );
  80 };
  81
  82 print("\td) No sequence, No quality, No ID...\n") if $DEBUG;
  83 warning_is { ok $wswq1 = Bio::Seq::Quality->new( -seq  =>       "",
  84                                     -qual =>    "",
  85                                     -verbose => 0
  86 ) } 'Got a sequence with no letters in it cannot guess alphabet';
  87
  88 print("Testing various methods and behaviors...\n") if $DEBUG;
  89
  90 print("1. Testing the seq() method...\n") if $DEBUG;
  91         print("\t1a) get\n") if $DEBUG;
  92         my $original_seq = $swq1->seq();
  93         is ($original_seq, "ATCGATCGA");
  94         print("\t1b) set\n") if $DEBUG;
  95         ok ($swq1->seq("AAAAAAAAAAAA"));
  96         print("\t1c) get (again, to make sure the set was done.)\n") if $DEBUG;
  97         is($swq1->seq(), "AAAAAAAAAAAA");
  98         print("\tSetting the sequence back to the original value...\n") if $DEBUG;
  99         $swq1->seq($original_seq);
 100
 101
 102 print("2. Testing the qual() method...\n") if $DEBUG;
 103         print("\t2a) get\n") if $DEBUG;
 104         my @qual = @{$swq1->qual()};
 105         my $str_qual = join(' ',@qual);
 106         is $str_qual, "10 20 30 40 50 40 30 20 10";
 107         print("\t2b) set\n") if $DEBUG;
 108         ok $swq1->qual("10 10 10 10 10");
 109         print("\t2c) get (again, to make sure the set was done.)\n") if $DEBUG;
 110         my @qual2 = @{$swq1->qual()};
 111         my $str_qual2 = join(' ',@qual2);
 112         is($str_qual2, "10 10 10 10 10 0 0 0 0"); ###!
 113         print("\tSetting the quality back to the original value...\n") if $DEBUG;
 114         $swq1->qual($str_qual);
 115
 116 print("3. Testing the length() method...\n") if $DEBUG;
 117         print("\t3a) When lengths are equal...\n") if $DEBUG;
 118         is($swq1->length(), 9);
 119         print("\t3b) When lengths are different\n") if $DEBUG;
 120         $swq1->qual("10 10 10 10 10");
 121         isnt ($swq1->length(), "DIFFERENT");
 122
 123
 124 print("6. Testing the subqual() method...\n") if $DEBUG;
 125      my $t_subqual = "10 20 30 40 50 60 70 80 90";
 126      $swq1->qual($t_subqual);
 127      print("\t6d) Testing the subqual at the start (border condition)\n") if $DEBUG;
 128           # ok ('1 2 3' eq join(' ',@{$swq1->subqual(1,3)}));
 129      print("\t6d) Testing the subqual at the end (border condition)\n") if $DEBUG;
 130           # ok ('7 8 9' eq join(' ',@{$swq1->subqual(7,9)}));
 131      print("\t6d) Testing the subqual in the middle\n") if $DEBUG;
 132           # ok ('4 5 6' eq join(' ',@{$swq1->subqual(4,6)}));
 133
 134 print("7. Testing cases where quality is zero...\n") if $DEBUG;
 135 $swq1 = Bio::Seq::Quality->new(-seq =>  'G',
 136                                -qual => '0',
 137                                      );
 138 my $swq2 = Bio::Seq::Quality->new(-seq =>  'G',
 139                                   -qual => '65',
 140                                      );
 141 is $swq1->length, $swq2->length;
 142
 143 $swq1 = Bio::Seq::Quality->new(-seq =>  'GC',
 144                                -qual => '0 0',
 145                                      );
 146 $swq2 = Bio::Seq::Quality->new(-seq =>  'GT',
 147                                -qual => '65 0',
 148                                      );
 149 is $swq1->length, $swq2->length;
 150
 151
 152 #
 153 # end of test inherited from seqwithquality.t
 154 #
 155 #################################################################
 156 #
 157 # testing new functionality
 158 #
 159
 160 my $qual = '0 1 2 3 4 5 6 7 8 9 11 12';
 161 my $trace = '0 5 10 15 20 25 30 35 40 45 50 55';
 162
 163 ok my $seq = Bio::Seq::Quality->new
 164     ( -qual => $qual,
 165       -trace_indices => $trace,
 166       -seq =>  'atcgatcgatcg',
 167       -id  => 'human_id',
 168       -accession_number => 'S000012',
 169       -verbose => $DEBUG >= 0 ? $DEBUG : 0
 170 );
 171
 172 is_deeply $seq->qual, [split / /, $qual];
 173 is_deeply $seq->trace, [split / /, $trace];
 174 is_deeply $seq->trace_indices, [split / /, $trace]; #deprecated
 175
 176 is $seq->qual_text, $qual;
 177 is $seq->trace_text, $trace;
 178
 179 is join (' ', @{$seq->subqual(2, 3)}), '1 2';
 180 is $seq->subqual_text(2, 3), '1 2';
 181 is join (' ', @{$seq->subqual(2, 3, "9 9")}), '9 9';
 182 is $seq->subqual_text(2, 3, "8 8"), '8 8';
 183
 184 is join (' ', @{$seq->subtrace(2, 3)}), '5 10';
 185 is $seq->subtrace_text(2, 3), '5 10';
 186 is join (' ', @{$seq->subtrace(2, 3, "9 9")}), '9 9';
 187 is $seq->subtrace_text(2, 3, "8 8"), '8 8';
 188
 189
 190 is $seq->trace_index_at(5), 20;
 191 is join(' ', @{$seq->sub_trace_index(5,6)}), "20 25";
 192
 193 is $seq->baseat(2), 't';
 194
 195
 196 #############################################
 197 #
 198 # same tests using Seq::Meta::Array methods follow ...
 199 #
 200
 201 my $meta = '0 1 2 3 4 5 6 7 8 9 11 12';
 202 $trace = '0 5 10 15 20 25 30 35 40 45 50 55';
 203 my @trace_array = qw(0 5 10 15 20 25 30 35 40 45 50 55);
 204
 205 ok $seq = Bio::Seq::Quality->new
 206     ( -meta => $meta,
 207       -seq =>  'atcgatcgatcg',
 208       -id  => 'human_id',
 209       -accession_number => 'S000012',
 210       -verbose => $DEBUG >= 0 ? $DEBUG : 0
 211 );
 212
 213 $seq->named_meta('trace', \@trace_array);
 214
 215 is_deeply $seq->meta, [split / /, $meta];
 216 is_deeply $seq->named_meta('trace'), [split / /, $trace];
 217
 218 is $seq->meta_text, $meta;
 219 is $seq->named_meta_text('trace'), $trace;
 220
 221 is join (' ', @{$seq->submeta(2, 3)}), '1 2';
 222 is $seq->submeta_text(2, 3), '1 2';
 223 is join (' ', @{$seq->submeta(2, 3, "9 9")}), '9 9';
 224 is $seq->submeta_text(2, 3, "8 8"), '8 8';
 225
 226 is join (' ', @{$seq->named_submeta('trace', 2, 3)}), '5 10';
 227 is $seq->named_submeta_text('trace', 2, 3), '5 10';
 228 is join (' ', @{$seq->named_submeta('trace', 2, 3, "9 9")}), '9 9';
 229 is $seq->named_submeta_text('trace', 2, 3, "8 8"), '8 8';
 230
 231
 232 ok $seq = Bio::Seq::Quality->new(
 233     -seq => "ATGGGGGTGGTGGTACCCTATGGGGGTGGTGGTACCCT",
 234     -qual => "10 59 12 75 63 76 84 36 42 10 35 97 81 50 81 53 93 13 38 10 59 12 75 63 76 84 36 42 10 35 97 81 50 81 53 93 13 38",
 235     -trace_indices => "1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38"
 236 );
 237
 238 my $rev;
 239 ok $rev = $seq->revcom;
 240 is $rev->seq, 'AGGGTACCACCACCCCCATAGGGTACCACCACCCCCAT';
 241 is $rev->qual_text, "38 13 93 53 81 50 81 97 35 10 42 36 84 76 63 75 12 59 10 38 13 93 53 81 50 81 97 35 10 42 36 84 76 63 75 12 59 10";