Using Bioperl To Match Primer Sites. Difficulty With Iupac.Pm
Entering edit mode
8.3 years ago
Daniel ★ 3.8k

Hi All

I have been trying to write a script to take IUPAC compliment primers and pull out the corresponding subsequence from a full length (Aiming to select sequence regions from an amplicon).

The code works for direct matches but I am having trouble implementing the IUPAC module for generating a regexp on the fly. The code below only includes generating the forward primer regex at the moment. I also want to change the trimming portion to a substring argument but that's for when I've worked out how to find the character position of the match.

I have been mainly using Bio::Tools::IUPAC and Degenerate_primers but I have hit a brick wall, cant find any similar example scripts and would really appreciate any advice.


This is where I am so far:

use strict;
use warnings;

use Bio::Seq;
use Bio::SeqIO;
use Bio::Tools::IUPAC;

my $infile = $ARGV[0];
my $primer_f = $ARGV[1];
my $primer_r = $ARGV[2];

my $in = Bio::SeqIO->new(-file => "$infile", -format => 'fasta');
my $out = Bio::SeqIO->new(-file => ">out.fasta", -format => 'fasta');

##Generate forward matching regex 
my $ambiseq_f = Bio::Seq->new(-seq => "$primer_f", -alphabet => 'dna');

# Create all possible non-degenerate sequences
my $iupac = Bio::Tools::IUPAC->new(-seq => $ambiseq_f);
while (my $uniqueseq = $iupac->next_seq()) {
     # process the unique Bio::Seq object.    <---This is from the module documentation, Dont know what referring to.
        print $uniqueseq->seq . "\n";

my $regexp = $iupac->regexp();
#print "\n$regexp\n";

##Trim to matching region
while ( my $seq = $in->next_seq() ){
        my $trim5 = 0;
        my $trim53 = 0;

        if ( $seq->seq =~ m/^.*$primer_f(.*)/i ){
         $trim5 = $1;

        if ($trim5 =~ m/^(.*)$primer_r.*/i ){
                $trim53 = $1;
        #print ">" . $seq->id . "\n" . $trim53 . "\n";

EDIT: removed an array - pretty sure that's not right but forgot to remove it before.

bioperl • 3.5k views
Entering edit mode
8.3 years ago
Daniel ★ 3.8k

The issue was that I was reading the manual for the bioperl-live development package, then trying to use the stable version without realising.

For reference the development has a great method for generating degenerate primer regexs. Came in very useful.


Login before adding your answer.

Traffic: 1296 users visited in the last hour
Help About
Access RSS

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6