~ubuntu-branches/ubuntu/raring/bioperl/raring

« back to all changes in this revision

Viewing changes to Bio/PopGen/IO/hapmap.pm

  • Committer: Bazaar Package Importer
  • Author(s): Charles Plessy
  • Date: 2008-03-18 14:44:57 UTC
  • mfrom: (4 hardy)
  • mto: This revision was merged to the branch mainline in revision 6.
  • Revision ID: james.westby@ubuntu.com-20080318144457-1jjoztrvqwf0gruk
* debian/control:
  - Removed MIA Matt Hope (dopey) from the Uploaders field.
    Thank you for your work, Matt. I hope you are doing well.
  - Downgraded some recommended package to the 'Suggests' priority,
    according to the following discussion on Upstream's mail list.
    http://bioperl.org/pipermail/bioperl-l/2008-March/027379.html
    (Closes: #448890)
* debian/copyright converted to machine-readable format.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
# $Id: hapmap.pm,v 1.8.4.1 2006/10/02 23:10:23 sendu Exp $
 
2
#
 
3
# BioPerl module for Bio::PopGen::IO::hapmap
 
4
#
 
5
# Cared for by Rich Dobson <r.j.dobson-at-qmul.ac.uk>
 
6
#
 
7
# Copyright Rich Dobson
 
8
#
 
9
# You may distribute this module under the same terms as perl itself
 
10
 
 
11
# POD documentation - main docs before the code
 
12
 
 
13
=head1 NAME
 
14
 
 
15
Bio::PopGen::IO::hapmap - A parser for HapMap output data
 
16
 
 
17
=head1 SYNOPSIS
 
18
 
 
19
  # Do not use directly, use through the Bio::PopGen::IO driver
 
20
 
 
21
  use Bio::PopGen::IO;
 
22
  my $io = new Bio::PopGen::IO(-format => 'hapmap',
 
23
                               -file   => 'data.hapmap');
 
24
 
 
25
  # Some IO might support reading in a population at a time
 
26
 
 
27
  my @population;
 
28
  while( my $ind = $io->next_individual ) {
 
29
      push @population, $ind;
 
30
  }
 
31
 
 
32
=head1 DESCRIPTION
 
33
 
 
34
A driver module for Bio::PopGen::IO for parsing hapmap data.
 
35
 
 
36
=head1 FEEDBACK
 
37
 
 
38
=head2 Mailing Lists
 
39
 
 
40
User feedback is an integral part of the evolution of this and other
 
41
Bioperl modules. Send your comments and suggestions preferably to
 
42
the Bioperl mailing list.  Your participation is much appreciated.
 
43
 
 
44
  bioperl-l@bioperl.org                  - General discussion
 
45
  http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
 
46
 
 
47
=head2 Reporting Bugs
 
48
 
 
49
Report bugs to the Bioperl bug tracking system to help us keep track
 
50
of the bugs and their resolution. Bug reports can be submitted via
 
51
the web:
 
52
 
 
53
  http://bugzilla.open-bio.org/
 
54
 
 
55
=head1 AUTHOR - Rich Dobson
 
56
 
 
57
Email r.j.dobson-at-qmul.ac.uk
 
58
 
 
59
=head1 CONTRIBUTORS
 
60
 
 
61
Jason Stajich, jason-at-bioperl.org
 
62
 
 
63
=head1 APPENDIX
 
64
 
 
65
The rest of the documentation details each of the object methods.
 
66
Internal methods are usually preceded with a _
 
67
 
 
68
=cut
 
69
 
 
70
 
 
71
# Let the code begin...
 
72
 
 
73
package Bio::PopGen::IO::hapmap;
 
74
use vars qw($FieldDelim $AlleleDelim $NoHeader $StartingCol);
 
75
use strict;
 
76
 
 
77
($FieldDelim,$AlleleDelim,$NoHeader,$StartingCol) =( '\s+','',0,11);
 
78
 
 
79
use Bio::PopGen::Individual;
 
80
use Bio::PopGen::Population;
 
81
use Bio::PopGen::Genotype;
 
82
 
 
83
use base qw(Bio::PopGen::IO);
 
84
 
 
85
 
 
86
=head2 new
 
87
 
 
88
 Title   : new
 
89
 Usage   : my $obj = new Bio::PopGen::IO::hapmap();
 
90
 Function: Builds a new Bio::PopGen::IO::hapmap object 
 
91
 Returns : an instance of Bio::PopGen::IO::hapmap
 
92
 Args    : [optional, these are the current defaults] 
 
93
           -field_delimiter => ','
 
94
           -allele_delimiter=> '\s+'
 
95
           -no_header       => 0,
 
96
           -starting_column => 11
 
97
 
 
98
=cut
 
99
 
 
100
 
 
101
sub _initialize  {
 
102
 
 
103
    my($self, @args) = @_;
 
104
 
 
105
    $Bio::PopGen::Genotype::BlankAlleles='';
 
106
 
 
107
    my ($fieldsep,$all_sep, 
 
108
        $noheader, $start_col) = $self->_rearrange([qw(FIELD_DELIMITER
 
109
                                                       ALLELE_DELIMITER
 
110
                                                       NO_HEADER
 
111
                                                       STARTING_COLUMN)],
 
112
                                                   @args);
 
113
 
 
114
    $self->flag('no_header', defined $noheader ? $noheader : $NoHeader);
 
115
    $self->flag('field_delimiter',defined $fieldsep ? $fieldsep : $FieldDelim);
 
116
    $self->flag('allele_delimiter',defined $all_sep ? $all_sep : $AlleleDelim);
 
117
    $self->starting_column(defined $start_col ? $start_col : $StartingCol );
 
118
 
 
119
    $self->{'_header'} = undef;
 
120
    return 1;
 
121
 
 
122
}
 
123
 
 
124
=head2 flag
 
125
 
 
126
 Title   : flag
 
127
 Usage   : $obj->flag($flagname,$newval)
 
128
 Function: Get/Set the flag value
 
129
 Returns : value of a flag (a boolean)
 
130
 Args    : A flag name, currently we expect 
 
131
           'no_header', 'field_delimiter', or 'allele_delimiter' 
 
132
           on set, new value (a boolean or undef, optional)
 
133
 
 
134
=cut
 
135
 
 
136
sub flag  {
 
137
 
 
138
    my $self = shift;
 
139
    my $fieldname = shift;
 
140
    return unless defined $fieldname;
 
141
    return $self->{'_flag'}->{$fieldname} = shift if @_;
 
142
    return $self->{'_flag'}->{$fieldname};
 
143
 
 
144
}
 
145
 
 
146
sub _pivot {
 
147
    my ($self) = @_;
 
148
 
 
149
    my (@cols,@rows,@idheader);
 
150
    while ($_ = $self->_readline){
 
151
        chomp($_);
 
152
        next if( /^\s*\#/ || /^\s+$/ || ! length($_) );
 
153
        if( /^rs\#\s+alleles\s+chrom\s+pos\s+strand/ ) {
 
154
            @idheader = split $self->flag('field_delimiter');
 
155
        } else { 
 
156
            push @cols, [split $self->flag('field_delimiter')];
 
157
        }
 
158
    }
 
159
    my $startingcol = $self->starting_column;
 
160
 
 
161
    $self->{'_header'} = [ map { $_->[0] } @cols];
 
162
    for my $n ($startingcol.. $#{ $cols[ 0 ]}) { 
 
163
        my $column = [ $idheader[$n],
 
164
                       map{ $_->[ $n ] } @cols ];       
 
165
        push (@rows, $column); 
 
166
    }
 
167
    $self->{'_pivot'} = [@rows];
 
168
    $self->{'_i'} = 0;
 
169
}
 
170
 
 
171
 
 
172
=head2 next_individual
 
173
 
 
174
 Title   : next_individual
 
175
 Usage   : my $ind = $popgenio->next_individual;
 
176
 Function: Retrieve the next individual from a dataset
 
177
 Returns : A Bio::PopGen::IndividualI object
 
178
 Args    : none
 
179
 
 
180
See L<Bio::PopGen::IndividualI>
 
181
 
 
182
=cut
 
183
 
 
184
sub next_individual  {
 
185
    my ($self) = @_;
 
186
    unless($self->{'_pivot'}){
 
187
        #if it's the first time then pivot the table and store.
 
188
        #Lines will now be read from the stored pivot version of the input file
 
189
        $self->_pivot;
 
190
    }
 
191
 
 
192
    $_ = $self->{'_pivot'}->[$self->{'_i'}++];
 
193
 
 
194
    return unless defined $_;
 
195
 
 
196
    # Store all the marker related info. Now that the pivot has taken
 
197
    # place this is in the first few lines of the file Maybe this
 
198
    # should be put in a marker object. Doesn't seem to fit too well
 
199
    # though
 
200
 
 
201
    my ($samp,@marker_results) = @$_;
 
202
 
 
203
    # at some point use all this info
 
204
    my $i = 1;
 
205
    foreach my $m ( @marker_results ) {
 
206
        $m =~ s/^\s+//;
 
207
        $m =~ s/\s+$//;
 
208
        my $markername;
 
209
        if( defined $self->{'_header'} ) {
 
210
            $markername = $self->{'_header'}->[$i-1];
 
211
        } else { 
 
212
            $markername = "Marker$i";
 
213
        }
 
214
 
 
215
        my @alleles = split($self->flag('allele_delimiter'), $m);
 
216
        if( @alleles != 2 ) { 
 
217
            $self->warn("$m for $samp\n");
 
218
        } else { 
 
219
            $m = new Bio::PopGen::Genotype(-alleles      => \@alleles,
 
220
                                           -marker_name  => $markername,
 
221
                                           -individual_id=> $samp);
 
222
        }
 
223
        $i++; 
 
224
    }
 
225
 
 
226
    return new Bio::PopGen::Individual(-unique_id => $samp,
 
227
                                       -genotypes => \@marker_results);
 
228
 
 
229
}
 
230
 
 
231
=head2 next_population
 
232
 
 
233
 Title   : next_population
 
234
 Usage   : my $ind = $popgenio->next_population;
 
235
 Function: Retrieve the next population from a dataset
 
236
 Returns : Bio::PopGen::PopulationI object
 
237
 Args    : none
 
238
 Note    : Many implementation will not implement this
 
239
 
 
240
See L<Bio::PopGen::PopulationI>
 
241
 
 
242
=cut
 
243
 
 
244
sub next_population {
 
245
    my ($self) = @_;
 
246
    my @inds;
 
247
    while( my $ind = $self->next_individual ) {
 
248
        push @inds, $ind;
 
249
    }
 
250
    Bio::PopGen::Population->new(-individuals => \@inds);
 
251
}
 
252
 
 
253
=head2 write_individual
 
254
 
 
255
 Title   : write_individual
 
256
 Usage   : $popgenio->write_individual($ind);
 
257
 Function: Write an individual out in the file format
 
258
           NOT SUPPORTED  BY hapmap format
 
259
 Returns : none
 
260
 Args    : Bio::PopGen::PopulationI object(s)
 
261
 
 
262
See L<Bio::PopGen::PopulationI>
 
263
 
 
264
=cut
 
265
 
 
266
sub write_individual {
 
267
    my ($self,@inds) = @_;
 
268
 
 
269
    # data from hapmap is output, not input, so 
 
270
    # we don't need a method for writing and input file
 
271
 
 
272
    $self->throw_not_implemented();
 
273
}
 
274
 
 
275
=head2 write_population
 
276
 
 
277
 Title   : write_population
 
278
 Usage   : $popgenio->write_population($pop);
 
279
 Function: Write a population out in the file format
 
280
           NOT SUPPORTED  BY hapmap format
 
281
 Returns : none
 
282
 Args    : Bio::PopGen::PopulationI object(s)
 
283
 Note    : Many implementation will not implement this
 
284
 
 
285
See L<Bio::PopGen::PopulationI>
 
286
 
 
287
=cut
 
288
 
 
289
sub write_population {
 
290
    my ($self,@inds) = @_;
 
291
    $self->throw_not_implemented();
 
292
}
 
293
 
 
294
 
 
295
=head2 starting_column
 
296
 
 
297
 Title   : starting_column
 
298
 Usage   : $obj->starting_column($newval)
 
299
 Function: Column where data starts
 
300
 Example : 
 
301
 Returns : value of starting_column (a scalar)
 
302
 Args    : on set, new value (a scalar or undef, optional)
 
303
 
 
304
=cut
 
305
 
 
306
sub starting_column{
 
307
    my $self = shift;
 
308
 
 
309
    return $self->{'starting_column'} = shift if @_;
 
310
    return $self->{'starting_column'};
 
311
}
 
312
 
 
313
1;