☰

다른 파일의 ID를 사용하여 gff3 파일에서 데이터를 추출하는 방법은 무엇입니까?

tag-icon • tag-icon perl
다른 파일의 ID를 사용하여 gff3 파일에서 데이터를 추출하는 방법은 무엇입니까?

여러 ID가 포함된 파일이 있습니다.

File 1:
g24007
g51692

그리고 gff3 파일은 다음과 같습니다

File2:
    # start gene g24007
    scaffold591 method  gene    3322458 3376057 0.41    -   .   ID=g24007
    scaffold591 method  transcript  3322458 3376057 0.41    -   .   ID=g24007.t1;Parent=g24007
    scaffold591 method  transcription_end_site  3322458 3322458 .   -   .   Parent=g24007.t1
    scaffold591 method  CDS 3323084 3323326 1   -   0   ID=g24007.t1.cds;Parent=g24007.t1
    # coding sequence = [atggaaaaagctaaagatggcgaagagagcccaagtgaggcatctcctccagcccaggtggggcttgaaaatatccctg
    # cgacggtgtctggggaggagggccagctgctgtatcacgaggagactatcgatcttggtggagacgagtttgggtctgaagagaatgaggaaccctca
    --
    # end gene g24007
    # start gene g20000
    scaffold591 method  gene    3322458 3376057 0.41    -   .   ID=g20000
    scaffold591 method  transcript  3322458 3376057 0.41    -   .   ID=g20000.t1;Parent=g20000
    ffold591    method  intron  3356166 3369049 1   -   .   Parent=g20000.t1
    scaffold591 method  CDS 3323084 3323326 1   -   0   ID=g20000.t1.cds;Parent=g20000.t1
    # coding sequence = [atggaaaaagctaaagatggcgaagagagcccaagtgaggcatctcctccagcccaggtggggcttgaaaatatccctg
    --
    # end gene g20000

여기서는 file1의 ID를 매핑하고 file2에서 해당 데이터, 즉 "Start Gene"과 "End Gene" 사이에 있는 데이터를 추출하려고 합니다. 또한 원하는 출력에서 ​​"Coding Sequence"를 제외하고 싶습니다.

Expected output:
# start gene g24007
scaffold591 method  gene    3322458 3376057 0.41    -   .   ID=g24007
scaffold591 method  transcript  3322458 3376057 0.41    -   .   ID=g24007.t1;Parent=g24007
scaffold591 method  transcription_end_site  3322458 3322458 .   -   .   Parent=g24007.t1
scaffold591 method  CDS 3323084 3323326 1   -   0   ID=g24007.t1.cds;Parent=g24007.t1
# end gene g24007

저는 펄을 사용해 보았습니다.

My code:
use strict;
use warnings;
use Data::Dumper;

my $file1 = 'IDs.txt';
open FILE1, "<", $file1 or die $!;
my $file2 = 'gff3.txt';

open FILE2, "<", $file2 or die $!;
my %id;
my @array;

while(<FILE1>)
{
    $id{$_} = 1;
}
#print Dumper \%id;

my $gene_id = 0;
while (<FILE2>)
{
    if($_ !~ /^#/)
    {
        @array = split(/\t/,$_);
        $array[8] =~ s/ID=//g;
        if($id{$_})
        {
            print $_, @array;
        }
    }
}
close FILE1;
close FILE2;

답변1

@Hari 예상되는 출력을 보지 않고 표준 gff3 파일을 사용해 보았습니다. 그러나 내 스크립트는 "#startgene" 및 "#endgene" 줄을 인쇄하지 않습니다. 이것이 당신에게 도움이 되기를 바랍니다

Code:
#!/usr/local/perl

use strict;
use warnings;

my $file1 = $ARGV[0];
my $file2 = $ARGV[1];
my $output_file = $ARGV[2];

my %id;
my $ctr = 0;
open(IN, $file1);
while(<IN>)
{

    $_ =~ s/\n|\r//g;
    $ctr++;
    $id{$_} = $ctr;
}
close IN;

open(IN, $file2);
open(OUT, ">".$output_file);
while(<IN>)
{
    $_ =~ s/\n|\r//g;

    if($_ !~ /^#/)
    {
        my @tmp = split(/\t/, $_);

        if($tmp[8] =~ /ID=g(\d+)/)
        {
            my $gene_id = "g".$1;
            if(exists $id{$gene_id})
            {
                print OUT $_."\n";
            }
        }
        elsif($tmp[8] =~ /Parent=g(\d+)\.t(\d+)/)
        {
            my $gene_id = "g".$1;
            if(exists $id{$gene_id})
            {
                print OUT $_."\n";
            }
        }
    }
}
close IN;
close OUT;

관련 정보