.\" Automatically generated by Pod::Man 2.27 (Pod::Simple 3.28) .\" .\" Standard preamble: .\" ======================================================================== .de Sp \" Vertical space (when we can't use .PP) .if t .sp .5v .if n .sp .. .de Vb \" Begin verbatim text .ft CW .nf .ne \\$1 .. .de Ve \" End verbatim text .ft R .fi .. .\" Set up some character translations and predefined strings. \*(-- will .\" give an unbreakable dash, \*(PI will give pi, \*(L" will give a left .\" double quote, and \*(R" will give a right double quote. \*(C+ will .\" give a nicer C++. Capital omega is used to do unbreakable dashes and .\" therefore won't be available. \*(C` and \*(C' expand to `' in nroff, .\" nothing in troff, for use with C<>. .tr \(*W- .ds C+ C\v'-.1v'\h'-1p'\s-2+\h'-1p'+\s0\v'.1v'\h'-1p' .ie n \{\ . ds -- \(*W- . ds PI pi . if (\n(.H=4u)&(1m=24u) .ds -- \(*W\h'-12u'\(*W\h'-12u'-\" diablo 10 pitch . if (\n(.H=4u)&(1m=20u) .ds -- \(*W\h'-12u'\(*W\h'-8u'-\" diablo 12 pitch . ds L" "" . ds R" "" . ds C` "" . ds C' "" 'br\} .el\{\ . ds -- \|\(em\| . ds PI \(*p . ds L" `` . ds R" '' . ds C` . ds C' 'br\} .\" .\" Escape single quotes in literal strings from groff's Unicode transform. .ie \n(.g .ds Aq \(aq .el .ds Aq ' .\" .\" If the F register is turned on, we'll generate index entries on stderr for .\" titles (.TH), headers (.SH), subsections (.SS), items (.Ip), and index .\" entries marked with X<> in POD. Of course, you'll have to process the .\" output yourself in some meaningful fashion. .\" .\" Avoid warning from groff about undefined register 'F'. .de IX .. .nr rF 0 .if \n(.g .if rF .nr rF 1 .if (\n(rF:(\n(.g==0)) \{ . if \nF \{ . de IX . tm Index:\\$1\t\\n%\t"\\$2" .. . if !\nF==2 \{ . nr % 0 . nr F 2 . \} . \} .\} .rr rF .\" .\" Accent mark definitions (@(#)ms.acc 1.5 88/02/08 SMI; from UCB 4.2). .\" Fear. Run. Save yourself. No user-serviceable parts. . \" fudge factors for nroff and troff .if n \{\ . ds #H 0 . ds #V .8m . ds #F .3m . ds #[ \f1 . ds #] \fP .\} .if t \{\ . ds #H ((1u-(\\\\n(.fu%2u))*.13m) . ds #V .6m . ds #F 0 . ds #[ \& . ds #] \& .\} . \" simple accents for nroff and troff .if n \{\ . ds ' \& . ds ` \& . ds ^ \& . ds , \& . ds ~ ~ . ds / .\} .if t \{\ . ds ' \\k:\h'-(\\n(.wu*8/10-\*(#H)'\'\h"|\\n:u" . ds ` \\k:\h'-(\\n(.wu*8/10-\*(#H)'\`\h'|\\n:u' . ds ^ \\k:\h'-(\\n(.wu*10/11-\*(#H)'^\h'|\\n:u' . ds , \\k:\h'-(\\n(.wu*8/10)',\h'|\\n:u' . ds ~ \\k:\h'-(\\n(.wu-\*(#H-.1m)'~\h'|\\n:u' . ds / \\k:\h'-(\\n(.wu*8/10-\*(#H)'\z\(sl\h'|\\n:u' .\} . \" troff and (daisy-wheel) nroff accents .ds : \\k:\h'-(\\n(.wu*8/10-\*(#H+.1m+\*(#F)'\v'-\*(#V'\z.\h'.2m+\*(#F'.\h'|\\n:u'\v'\*(#V' .ds 8 \h'\*(#H'\(*b\h'-\*(#H' .ds o \\k:\h'-(\\n(.wu+\w'\(de'u-\*(#H)/2u'\v'-.3n'\*(#[\z\(de\v'.3n'\h'|\\n:u'\*(#] .ds d- \h'\*(#H'\(pd\h'-\w'~'u'\v'-.25m'\f2\(hy\fP\v'.25m'\h'-\*(#H' .ds D- D\\k:\h'-\w'D'u'\v'-.11m'\z\(hy\v'.11m'\h'|\\n:u' .ds th \*(#[\v'.3m'\s+1I\s-1\v'-.3m'\h'-(\w'I'u*2/3)'\s-1o\s+1\*(#] .ds Th \*(#[\s+2I\s-2\h'-\w'I'u*3/5'\v'-.3m'o\v'.3m'\*(#] .ds ae a\h'-(\w'a'u*4/10)'e .ds Ae A\h'-(\w'A'u*4/10)'E . \" corrections for vroff .if v .ds ~ \\k:\h'-(\\n(.wu*9/10-\*(#H)'\s-2\u~\d\s+2\h'|\\n:u' .if v .ds ^ \\k:\h'-(\\n(.wu*10/11-\*(#H)'\v'-.4m'^\v'.4m'\h'|\\n:u' . \" for low resolution devices (crt and lpr) .if \n(.H>23 .if \n(.V>19 \ \{\ . ds : e . ds 8 ss . ds o a . ds d- d\h'-1'\(ga . ds D- D\h'-1'\(hy . ds th \o'bp' . ds Th \o'LP' . ds ae ae . ds Ae AE .\} .rm #[ #] #H #V #F C .\" ======================================================================== .\" .IX Title "Bio::Tools::dpAlign 3pm" .TH Bio::Tools::dpAlign 3pm "2014-07-13" "perl v5.18.2" "User Contributed Perl Documentation" .\" For nroff, turn off justification. Always turn off hyphenation; it makes .\" way too many mistakes in technical documents. .if n .ad l .nh .SH "NAME" Bio::Tools::dpAlign \- Perl extension to do pairwise dynamic programming sequence alignment .SH "SYNOPSIS" .IX Header "SYNOPSIS" .Vb 5 \& use Bio::Tools::dpAlign; \& use Bio::SeqIO; \& use Bio::SimpleAlign; \& use Bio::AlignIO; \& use Bio::Matrix::IO; \& \& $seq1 = Bio::SeqIO\->new(\-file => $ARGV[0], \-format => \*(Aqfasta\*(Aq); \& $seq2 = Bio::SeqIO\->new(\-file => $ARGV[1], \-format => \*(Aqfasta\*(Aq); \& \& # create a dpAlign object \& # to do global alignment, specify DPALIGN_GLOBAL_MILLER_MYERS \& # to do ends\-free alignment, specify DPALIGN_ENDSFREE_MILLER_MYERS \& $factory = new dpAlign(\-match => 3, \& \-mismatch => \-1, \& \-gap => 3, \& \-ext => 1, \& \-alg => Bio::Tools::dpAlign::DPALIGN_LOCAL_MILLER_MYERS); \& \& # actually do the alignment \& $out = $factory\->pairwise_alignment($seq1\->next_seq, $seq2\->next_seq); \& $alnout = Bio::AlignIO\->new(\-format => \*(Aqpfam\*(Aq, \-fh => \e*STDOUT); \& $alnout\->write_aln($out); \& \& # To do protein alignment, set the sequence type to protein \& # By default all protein alignments are using BLOSUM62 matrix \& # the gap opening cost is 7 and gap extension is 1. These \& # values are from ssearch. To use your own custom substitution \& # matrix, you can create a Bio::Matrix::MatrixI object. \& \& $parser = Bio::Matrix::IO\->new(\-format => \*(Aqscoring\*(Aq, \-file => \*(Aqblosum50.mat\*(Aq); \& $matrix = $parser\->next_matrix; \& $factory = Bio::Tools::dpAlign\->new(\-matrix => $matrix, \-alg => Bio::Tools::dpAlign::DPALIGN_LOCAL_MILLERMYERS); \& $seq1\->alphabet(\*(Aqprotein\*(Aq); \& $seq2\->alphabet(\*(Aqprotein\*(Aq); \& $out = $factory\->pairwise_alignment($seq1\->next_seq, $seq2\->next_seq); \& $alnout\->write_aln($out); \& \& # use the factory to make some output \& \& $factory\->align_and_show($seq1, $seq2, STDOUT); \& \& # use Phil Green\*(Aqs algorithm to calculate the optimal local \& # alignment score between two sequences quickly. It is very \& # useful when you are searching a query sequence in a database \& # of sequences. Since finding a alignment is more costly \& # than just calculating scores, you can save time if you only \& # align sequences that have a high alignment score. \& \& # To use this feature, first you call the sequence_profile function \& # to obtain the profile of the query sequence. \& $profile = $factory\->sequence_profile($query); \& \& %scores = (); \& # Then use a loop to run a database of sequences against the \& # profile to obtain a table of alignment scores \& $dbseq = Bio::SeqIO(\-file => \*(Aqdbseq.fa\*(Aq, \-format => \*(Aqfasta\*(Aq); \& while (defined($seq = $dbseq\->next_seq)) { \& $scores{$seq\->id} = $factory\->pairwise_alignment_score($profile, $seq); \& } .Ve .SH "DESCRIPTION" .IX Header "DESCRIPTION" Dynamic Programming approach is considered to be the most sensitive way to align two biological sequences. There are currently three major types of dynamic programming algorithms: Global Alignment, Local Alignment and Ends-free Alignment. .PP Global Alignment compares two sequences in their entirety. By inserting gaps in the two sequences, it aligns two sequences to minimize the edit distance as defined by the gap cost function and the substitution matrix. Global Alignment is generally applied to two sequences that are very similar in length and content. .PP Local Alignment instead attempts to find out the subsequences that has the minimal edit distance among all possible subsequences. It is good for sequences that has a stretch of subsequences that are similar to each other. .PP Ends-free Alignment is a special case of Global Alignment. There are no gap penalty imposed for the gaps that extended from the end points of two sequences. Therefore it will be a good application when you think one sequence is contained by the other or when you think two sequences overlap each other. .PP Dynamic Programming was first introduced by Needleman-Wunsch (1970) to globally align two sequences. The idea of local alignment was later introduced by Smith-Waterman (1981). Gotoh (1982) improved both algorithms by introducing auxillary arrays that reduced the time complexity of the algorithms to O(m*n). Miller-Myers (1988) exploits the divide-and-conquer idea introduced by Hirschberg (1975) to solve the affine gap cost dynamic programming using only linear space. At the time of this writing, it is accepted that Miller-Myers is the fastest single \s-1CPU\s0 implementation and using the least memory that is truly equivalent to original algorithm introduced by Needleman-Wunsch. According to Aaron Mackey, Phil Green's \s-1SWAT\s0 implementation introduced a heuristic that does not consider paths through the matrix where the score would be less than the gap opening penalty, yielding a 1.5\-2X speedup on most comparisons. to skip the calculation of some cells. However, his approach is only good for calculating the minimum edit distance and find out the corresponding subsequences (aka search phase). Bill Pearson's popular dynamic programming alignment program \s-1SSEARCH\s0 uses Phil Green's algorithm to find the subsequences and then Miller-Myers's algorithm to find the actual alignment. (aka alignment phase) .PP The current implementation supports local alignment of either \s-1DNA\s0 sequences or protein sequences. It allows you to specify either the Miller-Myers Global Alignment (\s-1DPALIGN_GLOBAL_MILLER_MYERS\s0) or Miller-Myers Local Alignment (\s-1DPALIGN_LOCAL_MILLER_MYERS\s0). For \s-1DNA\s0 alignment, you can specify the scores for match, mismatch, gap opening cost and gap extension cost. For protein alignment, it is using \&\s-1BLOSUM62\s0 by default. Currently the substitution matrix is not configurable. .PP Note: If you supply LocatableSeq objects to pairwise_alignment, pairwise_alignment_score, align_and_show or sequence_profile and the sequence supplied contains gaps, these functions will treat these sequences as if they are without gaps. .SH "DEPENDENCIES" .IX Header "DEPENDENCIES" This package comes with the main bioperl distribution. You also need to install the lastest bioperl-ext package which contains the \s-1XS\s0 code that implements the algorithms. This package won't work if you haven't compiled the bioperl-ext package. .SH "TO-DO" .IX Header "TO-DO" .IP "1." 3 Basic support for \s-1IUPAC\s0 code for \s-1DNA\s0 sequence is now implemented. X will mismatch any character. T will match U. For others, whenever there is a possibility for match, it is considered a full match, for example, W will match B. .IP "2." 3 Allow custom substitution matrix for \s-1DNA.\s0 Note that for proteins, you can now use your own subsitution matirx. .SH "FEEDBACK" .IX Header "FEEDBACK" .SS "Mailing Lists" .IX Subsection "Mailing Lists" User feedback is an integral part of the evolution of this and other Bioperl modules. Send your comments and suggestions preferably to one of the Bioperl mailing lists. Your participation is much appreciated. .PP .Vb 2 \& bioperl\-l@bioperl.org \- General discussion \& http://bioperl.org/wiki/Mailing_lists \- About the mailing lists .Ve .SS "Support" .IX Subsection "Support" Please direct usage questions or support issues to the mailing list: .PP \&\fIbioperl\-l@bioperl.org\fR .PP rather than to the module maintainer directly. Many experienced and reponsive experts will be able look at the problem and quickly address it. Please include a thorough description of the problem with code and data examples if at all possible. .SS "Reporting Bugs" .IX Subsection "Reporting Bugs" Report bugs to the Bioperl bug tracking system to help us keep track the bugs and their resolution. Bug reports can be submitted via the web: .PP .Vb 1 \& https://github.com/bioperl/bioperl\-live/issues .Ve .SH "AUTHOR" .IX Header "AUTHOR" .Vb 7 \& This implementation was written by Yee Man Chan (ymc@yahoo.com). \& Copyright (c) 2003 Yee Man Chan. All rights reserved. This program \& is free software; you can redistribute it and/or modify it under \& the same terms as Perl itself. Special thanks to Aaron Mackey \& and WIlliam Pearson for the helpful discussions. [The portion \& of code inside pgreen subdirectory was borrowed from ssearch. It \& should be distributed in the same terms as ssearch.] .Ve .SS "sequence_profile" .IX Subsection "sequence_profile" .Vb 6 \& Title : sequence_profile \& Usage : $prof = $factory\->sequence_profile($seq1) \& Function: Makes a dpAlign_SequenceProfile object from one sequence \& Returns : A dpAlign_SequenceProfile object \& Args : The lone argument is a Bio::PrimarySeqI that we want to \& build a profile for. Usually, this would be the Query sequence .Ve .SS "pairwise_alignment_score" .IX Subsection "pairwise_alignment_score" .Vb 11 \& Title : pairwise_alignment_score \& Usage : $score = $factory\->pairwise_alignment_score($prof,$seq2) \& Function: Makes a SimpleAlign object from two sequences \& Returns : An integer that is the score of the optimal alignment. \& Args : The first argument is the sequence profile obtained from a \& call to the sequence_profile function. The second argument \& is a Bio::PrimarySeqI object to be aligned. The second argument \& is usually a sequence in the database sequence. Note \& that this function only uses Phil Green\*(Aqs algorithm and \& therefore theoretically may not always give you the optimal \& score. .Ve .SS "pairwise_alignment" .IX Subsection "pairwise_alignment" .Vb 7 \& Title : pairwise_alignment \& Usage : $aln = $factory\->pairwise_alignment($seq1,$seq2) \& Function: Makes a SimpleAlign object from two sequences \& Returns : A SimpleAlign object if there is an alignment with positive \& score. Otherwise, return undef. \& Args : The first and second arguments are both Bio::PrimarySeqI \& objects that are to be aligned. .Ve .SS "align_and_show" .IX Subsection "align_and_show" .Vb 2 \& Title : align_and_show \& Usage : $factory\->align_and_show($seq1,$seq2,STDOUT) .Ve .SS "match" .IX Subsection "match" .Vb 7 \& Title : match \& Usage : $match = $factory\->match() #get \& : $factory\->match($value) #set \& Function : the set get for the match score \& Example : \& Returns : match value \& Arguments : new value .Ve .SS "mismatch" .IX Subsection "mismatch" .Vb 7 \& Title : mismatch \& Usage : $mismatch = $factory\->mismatch() #get \& : $factory\->mismatch($value) #set \& Function : the set get for the mismatch penalty \& Example : \& Returns : mismatch value \& Arguments : new value .Ve .SS "gap" .IX Subsection "gap" .Vb 7 \& Title : gap \& Usage : $gap = $factory\->gap() #get \& : $factory\->gap($value) #set \& Function : the set get for the gap penalty \& Example : \& Returns : gap value \& Arguments : new value .Ve .SS "ext" .IX Subsection "ext" .Vb 7 \& Title : ext \& Usage : $ext = $factory\->ext() #get \& : $factory\->ext($value) #set \& Function : the set get for the ext penalty \& Example : \& Returns : ext value \& Arguments : new value .Ve .SS "alg" .IX Subsection "alg" .Vb 7 \& Title : alg \& Usage : $alg = $factory\->alg() #get \& : $factory\->alg($value) #set \& Function : the set get for the algorithm \& Example : \& Returns : alg value \& Arguments : new value .Ve