← Back to branch summary

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

« back to all changes in this revision

Viewing changes to scripts/utilities/bp_mrtrans.pl

  • Committer: Package Import Robot
  • Author(s): Charles Plessy
  • Date: 2013-09-22 13:39:48 UTC
  • mfrom: (3.1.11 sid)
  • Revision ID: package-import@ubuntu.com-20130922133948-c6z62zegjyp7ztou
Tags: 1.6.922-1
http://bugs.debian.org/722910
* New upstream release.
* Replaces and Breaks grinder (<< 0.5.3-3~) because of overlaping contents.
  Closes: #722910
* Stop Replacing and Breaking bioperl ( << 1.6.9 ): not needed anymore. 

Show diffs side-by-side

added added

removed removed

Lines of Context:
scripts/utilities/bp_mrtrans.pl
 
1
#!perl
 
2
use strict;
 
3
 
 
4
# Author:      Jason Stajich <jason-at-bioperl-dot-org>
 
5
# Description: Perl implementation of Bill Pearson's mrtrans
 
6
#              to project protein alignment back into cDNA coordinates
 
7
#
 
8
 
 
9
=head1 NAME
 
10
 
 
11
bp_mrtrans - implement a transformer of alignments from protein to mrna coordinates
 
12
 
 
13
=head1 SYNOPSIS
 
14
 
 
15
Usage:
 
16
  bp_mrtrans -i inputfile -o outputfile [-if input format] [-of output format] [-s cDNA sequence database]  [-sf cDNA sequence format] [-h]
 
17
 
 
18
=head1 DESCRIPTION
 
19
 
 
20
This script will convert a protein alignment back into a cDNA.  Loosely
 
21
based on Bill Pearson's mrtrans.
 
22
 
 
23
The options are:
 
24
 
 
25
   -o filename          - the output filename [default STDOUT]
 
26
   -of format           - output sequence format
 
27
                          (multiple sequence alignment)
 
28
                          [default phylip]
 
29
   -i filename          - the input filename [required]
 
30
   -if format           - input sequence format
 
31
                          (multiple sequence alignment)
 
32
                          [ default clustalw]
 
33
   -s --seqdb filename  - the cDNA sequence database file
 
34
   -sf --seqformat      - the cDNA seq db format (flatfile sequence format)
 
35
   -h                   - this help menu
 
36
 
 
37
=head1 AUTHOR
 
38
 
 
39
Jason Stajich, jason-at-bioperl-dot-org
 
40
 
 
41
=cut
 
42
 
 
43
use strict;
 
44
use warnings;
 
45
use Bio::Align::Utilities qw(aa_to_dna_aln);
 
46
use Bio::AlignIO;
 
47
use Bio::SeqIO;
 
48
use Getopt::Long;
 
49
 
 
50
# TODO - finish documentation,
 
51
#      - add support for extra options in output alignment formats
 
52
#        such as idnewline in phylip out to support Molphy input files
 
53
 
 
54
my ($iformat,$seqformat,$oformat,$seqdb,$input,$output) = ('clustalw','fasta',
 
55
                                                           'phylip');
 
56
my ($help,$usage);
 
57
 
 
58
$usage = "usage: bp_mrtrans.pl -i prot_alignment -if align_format -o out_dna_align -of output_format -s cDNA_seqdb -sf fasta\n".
 
59
"defaults: -if clustalw
 
60
          -of phylip
 
61
          -sf fasta\n";
 
62
 
 
63
GetOptions(
 
64
           'if|iformat:s'  => \$iformat,
 
65
           'i|input:s'     => \$input,
 
66
           'o|output:s'    => \$output,
 
67
           'of|outformat:s'=> \$oformat,
 
68
           's|seqdb|db:s'  => \$seqdb,
 
69
           'sf|seqformat:s'=> \$seqformat,
 
70
           'h|help'        => sub{ exec('perldoc',$0);
 
71
                                   exit(0)
 
72
                                   },
 
73
           );
 
74
 
 
75
$input ||= shift;
 
76
$seqdb ||= shift;
 
77
$output ||= shift;
 
78
if( ! defined $seqdb ) {
 
79
    die("cannot proceed without a valid seqdb\n$usage");
 
80
}
 
81
if( ! defined $input ) {
 
82
    die("cannot proceed without an input file\n$usage");
 
83
}
 
84
my $indb = new Bio::SeqIO(-file => $seqdb,
 
85
                          -format=>$seqformat);
 
86
my %seqs;
 
87
while( my $seq = $indb->next_seq ) {
 
88
    $seqs{$seq->id} = $seq;
 
89
}
 
90
 
 
91
my $in = new Bio::AlignIO(-format => $iformat,
 
92
                          -file   => $input);
 
93
my $out = new Bio::AlignIO(-format => $oformat,
 
94
                           -idlength => 22,
 
95
                           -interleaved => 0,
 
96
                           defined $output ? (-file   => ">$output") : () );
 
97
 
 
98
while( my $aln = $in->next_aln ) {
 
99
    my $dnaaln = aa_to_dna_aln($aln,\%seqs);
 
100
    $dnaaln->set_displayname_flat(1);
 
101
    $out->write_aln($dnaaln);
 
102
}
 
103
 
 
104
__END__