From Code to Community: Sponsoring The Perl and Raku Conference 2025 Learn more

SENDU / bioperl-1.5.2_004-RCf / t / SeqStats.t 

            

              
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
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121

              
# -*-Perl-*-
## Bioperl Test Harness Script for Modules
##
use strict;
BEGIN {
    # to handle systems with no installed Test module
    # we include the t dir (where a copy of Test.pm is located)
    # as a fallback
    eval { require Test; };
    if( $@ ) { use lib 't'; }
    use Test;
    plan tests => 36;
}
use Bio::PrimarySeq;
use Bio::Tools::SeqStats;
use vars ('$DEBUG');
my ($seqobj, $count, $seqobj_stats, $wt);
$seqobj = Bio::PrimarySeq->new(-seq=>'ACTGTGGCGTCAACTG',
                               -alphabet=>'dna', -id=>'test');
$seqobj_stats  =  Bio::Tools::SeqStats->new(-seq=>$seqobj);
ok defined($seqobj_stats) && ref($seqobj_stats) &&
    $seqobj_stats->isa('Bio::Tools::SeqStats');
$count = $seqobj_stats->count_monomers();  # for DNA sequence
ok $count->{'A'}, 3;
ok $count->{'C'}, 4;
ok $count->{'G'}, 5;
ok $count->{'T'}, 4;
$count = $seqobj_stats->count_codons();
ok $count->{'ACT'}, 2;
ok $count->{'GTG'}, 1;
ok $count->{'GCG'}, 1;
ok $count->{'TCA'}, 1;
$seqobj = Bio::PrimarySeq->new(-seq=>'ACTACTTCA', -alphabet=>'dna',
                               -id=>'test');
$seqobj_stats  =  Bio::Tools::SeqStats->new('-seq' => $seqobj);
$wt = $seqobj_stats->get_mol_wt();  # for DNA sequence
ok &round($$wt[0]), 2738 ;
$seqobj = Bio::PrimarySeq->new(-seq=>'ACXACNNCA',
                               -alphabet=>'dna', -id=>'test');
$wt = Bio::Tools::SeqStats->get_mol_wt($seqobj);
ok &round($$wt[0]), 2693;
ok &round($$wt[1]), 2813;
$seqobj = Bio::PrimarySeq->new(-seq=>'ACTGTGGCGTCAACTG',
                               -alphabet=>'dna', -id=>'test');
$count = Bio::Tools::SeqStats->count_monomers($seqobj);  # for DNA sequence
ok $count->{'A'}, 3;
ok $count->{'C'}, 4;
ok $count->{'G'}, 5;
ok $count->{'T'}, 4;
$seqobj = Bio::PrimarySeq->new(-seq=>'MQSERGITIDISLWKFETSKYYVT',
                               -alphabet=>'protein', -id=>'test');
$seqobj_stats  =  Bio::Tools::SeqStats->new('-seq' => $seqobj);
$count = $seqobj_stats->count_monomers();  # for amino sequence
ok $$count{'M'}, 1;
ok $$count{'I'}, 3;
ok $$count{'Y'}, 2;
ok $$count{'T'}, 3;
$wt = Bio::Tools::SeqStats->get_mol_wt($seqobj);
ok int $$wt[0], 2896;
ok int $$wt[1], 2896;
$seqobj = Bio::PrimarySeq->new(-seq=>'UYXUYNNYU', -alphabet=>'rna');
$wt = Bio::Tools::SeqStats->get_mol_wt($seqobj);
ok &round($$wt[0]), 2768;
ok &round($$wt[1]), 2891;
ok $seqobj = Bio::PrimarySeq->new(-seq=>'TGCCGTGTGTGCTGCTGCT', -alphabet=>'rna');
$wt = Bio::Tools::SeqStats->get_mol_wt($seqobj);
ok &round($$wt[0]), 6104 ;
# selenocysteine
ok $seqobj = Bio::PrimarySeq->new(-seq=>'MQSERGITIDISLWKFETSKYYVT',
                                  -alphabet=>'protein');
$wt = Bio::Tools::SeqStats->get_mol_wt($seqobj);
ok &round($$wt[0]), 2896 ;
#
# hydropathicity aka "gravy" score
#
# normal seq (should succeed)
ok $seqobj = Bio::PrimarySeq->new(-seq=>'MSFVLVAPDMLATAAADVVQIGSAVSAGS',
                                  -alphabet=>'protein');
my $gravy = Bio::Tools::SeqStats->hydropathicity($seqobj);
ok int($gravy*1000) == 1224;  # check to nearest 0.1%
# ambiguous sequence (should fail)
ok $seqobj = Bio::PrimarySeq->new(-seq=>'XXXB**BS', -alphabet=>'protein');
eval { Bio::Tools::SeqStats->hydropathicity($seqobj) };
ok $@ =~ /ambiguous amino acids/i;
# empty sequence (should fail)
ok $seqobj = Bio::PrimarySeq->new(-seq=>'', -alphabet=>'protein');
eval { Bio::Tools::SeqStats->hydropathicity($seqobj) };
ok $@ =~ /hydropathicity not defined/i;
# DNA sequence (should fail)
ok $seqobj = Bio::PrimarySeq->new(-seq=>'GATTACA', -alphabet=>'dna');
eval { Bio::Tools::SeqStats->hydropathicity($seqobj) };
ok $@ =~ /only meaningful for protein/;
#
# Extra functions
#
# perl does not have an explicit rounding function
sub round { return int ((shift @_) + 0.5 ) }