Re: bioperl module to extract specific nucleotides given chromosome and exact location of that nucleotide

"any ideas if there is a BioPerl module that can do that, or how to approach this problem ?"

I'm not familiar with the BioPerl modules; however, based on your problem description, a script to achieve this is very straightforward:

#!/usr/bin/env perl -l

use strict;
use warnings;
use autodie;

use Text::CSV;

my $ref_file = 'pm_1079751_reference.csv';
my @gen_files = qw{pm_1079751_source1.gen pm_1079751_source2.gen};
my $allele_sep = '-';

my %ref_data;

my $ref_csv = Text::CSV::->new;
open my $ref_fh, '<', $ref_file;
while (my $row = $ref_csv->getline($ref_fh)) {
    $ref_data{$row->[0]} = join $allele_sep, @$row[1,2];
} 
close $ref_fh;

for my $gen_file (@gen_files) {
    my $gen_csv = Text::CSV::->new({allow_whitespace => 1});
    open my $gen_fh, '<', $gen_file;
    while (my $row = $gen_csv->getline($gen_fh)) {
        next unless exists $ref_data{$row->[2]};
        if ($ref_data{$row->[2]} eq join $allele_sep, @$row[3,4]) {
            print join(',' => @$row[2,3,4]),
                    "\n\tfound in '$gen_file':\n\t\t", join ', ' => @$
+row;
        }
    }
    close $gen_fh;
}
[download]

You stated the file you've already created is in CSV format, so I've used the Text::CSV module to parse it.

[Side note: CSV files have, by convention, a .csv extension, not a .txt extension (as you indicated you used).]

Your description of the .gen files seems to indicate a CSV format so I've used Text::CSV for these as well. If that's a false assumption on my part, you'll need to make some code changes; however, the logic I've used shouldn't need reworking.

Using these test files:

$ cat pm_1079751_reference.csv
34787638,A,C
34787800,A,G
[download]

$ cat pm_1079751_source1.gen 
C1, M1, 99999999, A, C, ... source1 info ...
C1, M2, 34787638, A, G, ... source1 info ...
C1, M3, 34787800, A, C, ... source1 info ...
C1, M4, 34787800, A, G, ... source1 info ...
[download]

$ cat pm_1079751_source2.gen 
C1, M1, 99999999, A, C, ... source2 info ...
C1, M2, 34787638, A, G, ... source2 info ...
C1, M4, 34787800, A, G, ... source2 info ...
C1, M3, 34787638, A, C, ... source2 info ...
[download]

I get this output:

34787800,A,G
    found in 'pm_1079751_source1.gen':
        C1, M4, 34787800, A, G, ... source1 info ...
34787800,A,G
    found in 'pm_1079751_source2.gen':
        C1, M4, 34787800, A, G, ... source2 info ...
34787638,A,C
    found in 'pm_1079751_source2.gen':
        C1, M3, 34787638, A, C, ... source2 info ...
[download]

-- Ken

Comment on Re: bioperl module to extract specific nucleotides given chromosome and exact location of that nucleotide Select or Download Code

Replies are listed 'Best First'.
Re^2: bioperl module to extract specific nucleotides given chromosome and exact location of that nucleotide by xxArwaxx (Novice) on Mar 27, 2014 at 19:36 UTC
I appreciate your response, but I believe that it should be much more simpler than this, check the FINAL OUT portion in the post that I've just added, it explains it much more. Thanks!	[reply]
Re^3: bioperl module to extract specific nucleotides given chromosome and exact location of that nucleotide by kcott (Archbishop) on Mar 27, 2014 at 19:58 UTC
"I appreciate your response, but I believe that it should be much more simpler than this, check the FINAL OUT portion in the post that I've just added, it explains it much more. Thanks!" So, two days after you posted your question you decide to tell us the exact output format you want. You then have a problem with my sample output which, in the absence of any other information, I had to guess at. If you're incapable of modifying `print join(',' => @$row[2,3,4]), "\n\tfound in '$gen_file':\n\t\t", join ', ' => @$row;` [download] which outputs `34787800,A,G found in 'pm_1079751_source1.gen': C1, M4, 34787800, A, G, ... source1 info ...` [download] to code which just outputs the first line of that, i.e. `34787800,A,G` [download] I suggest you read "Perl introduction for beginners". -- Ken	[reply] [d/l] [select]


Welcome to the Monastery
	PerlMonks