从非标准格式的文本文件中提取核苷酸序列

Question

我的合作者在word文档中给了我一些DNA序列，我想将它们转换为一个文件中的一系列fasta序列。

我已经把它做成了一个文本文件，我想使用正则表达式来提取基因名称和序列：

use warnings;
use strict;

die "usage: make_fasta.pl <sequence file>" unless (@ARGV == 1);

my $seq_filename    = shift;
my $fasta_db_name   = $seq_filename . "_db.fa";

open(my $seq_file, '<', $seq_filename) 
          or die "can't open file $seq_filename, $!";
open(my $fasta_file, '>', $fasta_db_name) 
          or die "can't open file $fasta_db_name, $!";

while (my $line = <$seq_file>) {
    chomp $line;
    if ($line =~ /^[ATCG]+$/) { # if the line is entirely DNA seqence
       print $fasta_file "$line\n";            

    } elsif ($line =~ /Full-length (\w+) cDNA/) { # if the line has gene info
        print $fasta_file ">$1\n";

    } else {
        next;

    } 
}

但这给了我第一个基因的名字。显然，我对DNA正则表达式做了一些错误的事情，但我无论如何也解决不了它。在我看来，这与我在这个网站和其他网站上找到的其他建议的DNA测试完全相同。

我试图解析的文件是这样配置的：

Collaborators name

title of gene set

Full-length clock cDNA coding sequence 

ATGGTAGGATGTGTAATGCGTACGTGATCGT

Full-length per cDNA coding sequence

ATGCTAGCTACGTACGTAGCTACGTAGTACG

我希望输出是一个fasta文件，所以：

>clock
ATGGTAGGATGTGTAATGCGTACGTGATCGT
>per
ATGCTAGCTACGTACGTAGCTACGTAGTACG

实际输入文件的前几行是：

Dr Lin Zhang (Leicester University 10/2012) 

Canonical clock genes 

Full-length per cDNA coding seq (3693bp) 

ATGGACACAGGAACACCCCATGAAGATGTGCCCTCAGAGGACCACACCTTGGAAGAAGGGGACAGCAAGAACCCCTCGTGCCAGCAAGAGTCAGCCTACGGCTCCCTCGAGTCATCCTCCAATGGACAGTCTCAGAAAAGTTTCGGAGGAAGTGGAAGCAAAAGCTTAAATAGTGGTTCGAGTCACAGCAGCGGCTTTGGGGACCAAAATGATTTCAAGGGTATCCATCTTCACGAAGCGAAACACATAGCGTTGAAGAAGAAGAAAACTGGGAAAGGAGGTGAAAAGGTAGCAGAAATCCCCTTTCAAACTGCCTCTGAGGCAGAACTGTCCTCCAAAGGAAACGAAACAGAAAAGGAGAAAGAAACAAGCCTCGAGGAGTCTCCTGCTGCAAAAGAGGAAGCAATTATCGAAAAGGAGTCTCGTTACATCCACCCGAGGAACT

Answer 1

如果看不到实际输入文件的一部分，就很难回答这个问题。

您的示例输入和REGEX之间存在不匹配：

# looking for verbatim('Full-length') then <space> then one WORD_WITH_ALPHNUMERICS  then <space> and then verbatim 'cDNA'
$line =~ /Full-length (\w+) cDNA/;

您的示例输入行有‘全长’，没有破折号，基因名称有多个单词，而不是只有一个，并且末尾没有'cDNA‘。

如果你的输入行有'Full-length gene name with multiple words cDNA'，你的REGEX可以是：

$line=~/Full-length\s+(.*?)\s+cDNA/;

Answer 2

问题显然出在你的输入数据上。我修改了您发布的代码，生成了以下程序：

#!/usr/bin/env perl    

use warnings;
use strict;

while (my $line = <DATA>) {
    chomp $line;
    if ($line =~ /^[ATCG]+$/) { # if the line is entirely DNA seqence
       print "$line\n";            
    } elsif ($line =~ /Full-length (\w+) cDNA/) { # if the line has gene info
        print ">$1\n";
    } 
}


__DATA__
Collaborators name

title of gene set

Full-length clock cDNA coding sequence 

ATGGTAGGATGTGTAATGCGTACGTGATCGT

Full-length per cDNA coding sequence

ATGCTAGCTACGTACGTAGCTACGTAGTACG

并生成您指定的输出：

~$ src/tmp/cdna 
>clock
ATGGTAGGATGTGTAATGCGTACGTGATCGT
>per
ATGCTAGCTACGTACGTAGCTACGTAGTACG

我的修改只是为了使它成为自包含的，除了删除无用的else { next }子句之外，并没有改变任何流程控制或逻辑。

由于提供的虚拟数据似乎工作正常，您能否找到并发布几行为您失败的实际数据？