ic23oluk
ic23oluk

Reputation: 145

perl check for valid DNA sequence with regex

I want to write a subroutine that takes a FASTA file as an argument and prints out the sequence (without the header). The subroutine should check if the sequence contains any other letters than DNA bases (A, T, G, C).

Here's my code:

scalar_sequence ("sequence.fa");

sub scalar_sequence {
    my $file = $_[0];
    my $sequence;
    open (READ, $file) || die "Cannot open $file: $!.\n";
    while (<READ>){
        if (/^>/){
            next;
        } 
        if (/^[ATCG]/){
            $sequence .= $_;
        } else {
            die "invalid sequence\n";
        }
    }
    print $sequence, "\n";
}

When I run this code, I get 'invalid sequence' as output. When I leave the 'else' out, it prints out the sequence even when the sequence contains another letter.

What's the problem?

Thanks in advance!

Upvotes: 1

Views: 1147

Answers (1)

mkHun
mkHun

Reputation: 5927

The problem is here /^[ATCG]/ this line should be /^[ATCG]+$/

Your code should be

chomp;  
next if (/^>/); # skip for header
next if(/^\s*$/);  #skip for empty line
if (/^[ATCG]+$/){
        $sequence .= $_;
    } else {
        die "invalid sequence\n";
    }

You are only consider the beginning of the line start wit A or T or G or C. You should expand the matches.

Upvotes: 1

Related Questions