~ubuntu-branches/ubuntu/trusty/bioperl/trusty-proposed

« back to all changes in this revision

Viewing changes to t/SeqFeature/Unflattener2.t

  • Committer: Bazaar Package Importer
  • Author(s): Charles Plessy
  • Date: 2009-03-10 07:19:11 UTC
  • mfrom: (1.2.3 upstream)
  • Revision ID: james.westby@ubuntu.com-20090310071911-fukqzw54pyb1f0bd
Tags: 1.6.0-2
* Removed patch system (not used):
  - removed instuctions in debian/rules;
  - removed quilt from Build-Depends in debian/control.
* Re-enabled tests:
  - uncommented test command in debian/rules;
  - uncommented previously missing build-dependencies in debian/control.
  - Re-enabled tests and uncommented build-dependencies accordingly.
* Removed libmodule-build-perl and libtest-harness-perl from
  Build-Depends-Indep (provided by perl-modules).
* Better cleaning of empty directories using find -type d -empty -delete
  instead of rmdir in debian/rules (LP: #324001).

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
# -*-Perl-*- Test Harness script for Bioperl
 
2
# $Id: Unflattener2.t 15112 2008-12-08 18:12:38Z sendu $
 
3
 
 
4
use strict;
 
5
 
 
6
BEGIN { 
 
7
    use lib '.';
 
8
    use Bio::Root::Test;
 
9
    
 
10
    test_begin(-tests => 12);
 
11
        
 
12
        use_ok('Bio::SeqIO');
 
13
        use_ok('Bio::SeqFeature::Tools::Unflattener');
 
14
}
 
15
 
 
16
my $verbosity = test_debug();
 
17
 
 
18
my ($seq, @sfs);
 
19
my $unflattener = Bio::SeqFeature::Tools::Unflattener->new;
 
20
$unflattener->verbose($verbosity);
 
21
 
 
22
if (1) {
 
23
    
 
24
    # this is an arabidopsise gbk record. it has no mRNA features.
 
25
    # it has explicit exon/intron records
 
26
 
 
27
    my @path = ("ATF14F8.gbk");
 
28
    $seq = getseq(@path);
 
29
    
 
30
    is ($seq->accession_number, 'AL391144');
 
31
    my @topsfs = $seq->get_SeqFeatures;
 
32
    my @cdss = grep {$_->primary_tag eq 'CDS'} @topsfs;
 
33
    my $n = scalar(@topsfs);
 
34
    if( $verbosity > 0 ) {
 
35
        warn sprintf "TOP:%d\n", scalar(@topsfs);
 
36
        write_hier(@topsfs);
 
37
    }
 
38
    # UNFLATTEN
 
39
    @sfs = $unflattener->unflatten_seq(-seq=>$seq,
 
40
                                       -use_magic=>1,
 
41
                                      );
 
42
    @sfs = $seq->get_SeqFeatures;
 
43
    if( $verbosity > 0 ) {
 
44
        warn "\n\nPOST PROCESSING:\n";
 
45
        write_hier(@sfs);
 
46
        warn sprintf "PROCESSED/TOP:%d\n", scalar(@sfs);
 
47
    }
 
48
    is(@sfs,28);
 
49
    my @allsfs = $seq->get_all_SeqFeatures;
 
50
    is(@allsfs,202);
 
51
    my @mrnas = grep {$_->primary_tag eq 'mRNA'} @allsfs;
 
52
    if( $verbosity > 0 ) {
 
53
        warn sprintf "ALL:%d\n", scalar(@allsfs);
 
54
        warn sprintf "mRNAs:%d\n", scalar(@mrnas);
 
55
    }
 
56
 
 
57
    # relationship between mRNA and CDS should be one-one
 
58
    is(@mrnas,@cdss);
 
59
}
 
60
 
 
61
if (1) {
 
62
    
 
63
    # this is a record from FlyBase
 
64
    # it has mRNA features, and explicit exon/intron records
 
65
 
 
66
    my @path = ("AnnIX-v003.gbk");
 
67
    $seq = getseq(@path);
 
68
    
 
69
    my @topsfs = $seq->get_SeqFeatures;
 
70
    if( $verbosity > 0 ) {
 
71
        warn sprintf "TOP:%d\n", scalar(@topsfs);
 
72
        write_hier(@topsfs);
 
73
    }
 
74
    # UNFLATTEN
 
75
    @sfs = $unflattener->unflatten_seq(-seq=>$seq,
 
76
                                       -use_magic=>1,
 
77
                                      );
 
78
    @sfs = $seq->get_SeqFeatures;
 
79
    if( $verbosity > 0 ) {
 
80
        warn "\n\nPOST PROCESSING:\n";
 
81
        write_hier(@sfs);
 
82
        warn sprintf "PROCESSED/TOP:%d\n", scalar(@sfs);
 
83
    }
 
84
    is scalar(@sfs), 1;
 
85
    my @exons = grep {$_->primary_tag eq 'exon'} $seq->get_all_SeqFeatures;
 
86
    is scalar(@exons), 6;    # total number of exons per splice
 
87
    my %numberh = map {$_->get_tag_values("number") => 1} @exons;
 
88
    my @numbers = keys %numberh;
 
89
    if( $verbosity > 0 ) {
 
90
        warn sprintf "DISTINCT EXONS: %d [@numbers]\n", scalar(@numbers);
 
91
    }
 
92
    is scalar(@numbers), 6;  # distinct exons
 
93
}
 
94
 
 
95
if (1) {
 
96
    
 
97
    # example of a BAD genbank entry
 
98
 
 
99
    my @path = ("dmel_2Lchunk.gb");
 
100
    $seq = getseq(@path);
 
101
    
 
102
    my @topsfs = $seq->get_SeqFeatures;
 
103
    if( $verbosity > 0 ) {
 
104
        warn sprintf "TOP:%d\n", scalar(@topsfs);
 
105
        write_hier(@topsfs);
 
106
    }
 
107
    # UNFLATTEN
 
108
    #
 
109
    # we EXPECT problems with this erroneous record
 
110
    $unflattener->error_threshold(2);
 
111
    @sfs = $unflattener->unflatten_seq(-seq=>$seq,
 
112
                                       -use_magic=>1,
 
113
                                      );
 
114
    my @probs = $unflattener->get_problems;
 
115
    $unflattener->report_problems(\*STDERR) if $verbosity > 0;
 
116
    $unflattener->clear_problems;
 
117
    @sfs = $seq->get_SeqFeatures;
 
118
    if( $verbosity > 0 ) {
 
119
        warn "\n\nPOST PROCESSING:\n";
 
120
        write_hier(@sfs);
 
121
        warn sprintf "PROCESSED/TOP:%d\n", scalar(@sfs);
 
122
    }
 
123
    is scalar(@sfs), 2;
 
124
    my @exons = grep {$_->primary_tag eq 'exon'} $seq->get_all_SeqFeatures;
 
125
    is scalar(@exons), 2;    # total number of exons per splice
 
126
    if( $verbosity > 0 ) {
 
127
        warn sprintf "PROBLEMS ENCOUNTERED: %d (EXPECTED: 6)\n", scalar(@probs);
 
128
    }
 
129
    is scalar(@probs), 6;
 
130
}
 
131
 
 
132
 
 
133
sub write_hier {
 
134
    my @sfs = @_;
 
135
    _write_hier(0, @sfs);
 
136
}
 
137
 
 
138
sub _write_hier {
 
139
    my $indent = shift;
 
140
    my @sfs = @_;
 
141
    foreach my $sf (@sfs) {
 
142
        my $label = '?';
 
143
        if ($sf->has_tag('gene')) {
 
144
            ($label) = $sf->get_tag_values('gene');
 
145
        }
 
146
        if ($sf->has_tag('product')) {
 
147
            ($label) = $sf->get_tag_values('product');
 
148
        }
 
149
        if ($sf->has_tag('number')) {
 
150
            $label = join("; ", $sf->get_tag_values('number'));
 
151
        }
 
152
        printf "%s%s $label\n", '  ' x $indent, $sf->primary_tag;
 
153
        my @sub_sfs = $sf->sub_SeqFeature;
 
154
        _write_hier($indent+1, @sub_sfs);
 
155
    }
 
156
}
 
157
 
 
158
sub getseq {
 
159
    my @path = @_;
 
160
    my $seqio =
 
161
      Bio::SeqIO->new('-file'=> test_input_file(@path), 
 
162
                      '-format' => 'GenBank');
 
163
    $seqio->verbose($verbosity);
 
164
 
 
165
    my $seq = $seqio->next_seq();
 
166
    return $seq;
 
167
}