Using Bioperl To Match Primer Sites. Difficulty With Iupac.Pm
1
0
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.

Thanks.

This is where I am so far:

#!/usr/bin/perl
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";
        #$out->write_seq($trim53);
}

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

bioperl • 3.5k views
ADD COMMENT
4
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 IUPAC.pm has a great method for generating degenerate primer regexs. Came in very useful.

ADD COMMENT

Login before adding your answer.

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

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

Powered by the version 2.3.6