Bio/Assembly/Singlet.pm

   1 #
   2 # BioPerl module for Bio::Assembly::Singlet
   3 #
   4 # Please direct questions and support issues to <bioperl-l@bioperl.org>
   5 #
   6 # Cared for by Chad Matsalla <bioinformatics1 at dieselwurks.com>
   7 #
   8 # Copyright Chad Matsalla
   9 #
  10 # You may distribute this module under the same terms as perl itself
  11
  12 # POD documentation - main docs before the code
  13
  14 =head1 NAME
  15
  16 Bio::Assembly::Singlet - Perl module to hold and manipulate
  17                      singlets from sequence assembly contigs.
  18
  19 =head1 SYNOPSIS
  20
  21     # Module loading
  22     use Bio::Assembly::IO;
  23
  24     # Assembly loading methods
  25     $aio = Bio::Assembly::IO->new( -file   => 'test.ace.1',
  26                                    -format => 'phrap'      );
  27
  28     $assembly = $aio->next_assembly;
  29     foreach $singlet ($assembly->all_singlets) {
  30       # do something
  31     }
  32
  33     # OR, if you want to build the singlet yourself,
  34
  35     use Bio::Assembly::Singlet;
  36     $singlet = Bio::Assembly::Singlet->new(
  37         -id     => 'Singlet1',
  38         -seqref => $seq
  39     );
  40
  41 =head1 DESCRIPTION
  42
  43 A singlet is a sequence that phrap was unable to align to any other sequences.
  44
  45 =head1 FEEDBACK
  46
  47 =head2 Mailing Lists
  48
  49 User feedback is an integral part of the evolution of this and other
  50 Bioperl modules. Send your comments and suggestions preferably to the
  51 Bioperl mailing lists  Your participation is much appreciated.
  52
  53   bioperl-l@bioperl.org                  - General discussion
  54   http://bioperl.org/wiki/Mailing_lists  - About the mailing lists
  55
  56 =head2 Support
  57
  58 Please direct usage questions or support issues to the mailing list:
  59
  60 I<bioperl-l@bioperl.org>
  61
  62 rather than to the module maintainer directly. Many experienced and
  63 reponsive experts will be able look at the problem and quickly
  64 address it. Please include a thorough description of the problem
  65 with code and data examples if at all possible.
  66
  67 =head2 Reporting Bugs
  68
  69 Report bugs to the Bioperl bug tracking system to help us keep track
  70 the bugs and their resolution.  Bug reports can be submitted via the
  71 web:
  72
  73   https://redmine.open-bio.org/projects/bioperl/
  74
  75 =head1 AUTHOR - Chad S. Matsalla
  76
  77 bioinformatics1 at dieselwurks.com
  78
  79 =head1 APPENDIX
  80
  81 The rest of the documentation details each of the object
  82 methods. Internal methods are usually preceded with a _
  83
  84 =cut
  85
  86 #'
  87 package Bio::Assembly::Singlet;
  88
  89 use strict;
  90
  91 use Bio::SeqFeature::Collection;
  92 use Bio::LocatableSeq;
  93 use Bio::Seq::PrimaryQual;
  94
  95 use base qw(Bio::Assembly::Contig Bio::Root::Root Bio::Align::AlignI);
  96
  97 =head2 new
  98
  99     Title   : new
 100     Usage   : $singlet = $io->new( -seqref => $seq )
 101     Function: Create a new singlet object
 102     Returns : A Bio::Assembly::Singlet object
 103     Args    : -seqref => Bio::Seq-compliant sequence object for the singlet
 104
 105 =cut
 106
 107 sub new {
 108     my ($class, @args) = @_;
 109     my $self = $class->SUPER::new(@args);
 110     my ($seqref) = $self->_rearrange([qw(SEQREF)], @args);
 111     $self->{'_seqref'} = undef;
 112     if (defined $seqref) {
 113         $self->seqref($seqref);
 114     }
 115     return $self;
 116 }
 117
 118 =head2 seqref
 119
 120     Title   : seqref
 121     Usage   : $seqref = $singlet->seqref($seq);
 122     Function: Get/set the sequence to which this singlet refers
 123     Returns : A Bio::Seq-compliant object
 124     Args    : A Bio::Seq-compliant or Bio::Seq::Quality object
 125
 126 =cut
 127
 128 sub seqref {
 129     my ($self,$seq) = @_;
 130     if (defined $seq) { $self->_seq_to_singlet($seq) };
 131     return $self->{'_seqref'};
 132 }
 133
 134 =head2 _seq_to_singlet
 135
 136     Title   : _seq_to_singlet
 137     Usage   : $singlet->seqref($seq)
 138     Function: Transform a sequence into a singlet
 139     Returns : 1 for sucess
 140     Args    : A Bio::Seq-compliant object
 141
 142 =cut
 143
 144 sub _seq_to_singlet {
 145     my ($self, $seq) = @_;
 146     # Object type checking
 147     $self->throw("Unable to process non Bio::Seq-compliant object [".ref($seq)."]")
 148       unless ( defined $seq && ($seq->isa('Bio::PrimarySeqI') || $seq->isa('Bio::Seq::Quality')) );
 149     # Sanity check
 150     $self->throw("Unable to have more than one sequence reference in a singlet")
 151       if (defined $self->{'_seqref'});
 152     # From sequence to locatable sequence
 153     my $lseq = Bio::LocatableSeq->new(
 154         -id     => $seq->id,
 155         -seq    => $seq->seq,
 156         -strand => $seq->isa('Bio::LocatableSeq') ? $seq->strand : 1,
 157         -start  => 1,
 158         #-end   => we let Bio::LocatableSeq calculate it (Seq and LocatableSeq)
 159     );
 160     # Add new sequence and its coordinates to the contig
 161     my $lcoord = Bio::SeqFeature::Generic->new( -start => $lseq->start,
 162                                                 -end   => $lseq->end    );
 163     $self->set_seq_coord( $lcoord, $lseq );
 164     $self->{'_seqref'} = $lseq;
 165     # Creating consensus
 166     $self->set_consensus_sequence($lseq);
 167     if ($seq->isa("Bio::Seq::Quality")) {
 168         my $qual = Bio::Seq::PrimaryQual->new( -id   => $seq->id,
 169                                                -qual => $seq->qual );
 170         $self->set_consensus_quality($qual);
 171     }
 172     return 1;
 173 }
 174
 175 1;