Shell Programming and Scripting

Hi I have a multi-step project I am working on and have been finding it difficult to come up with the correct approach.
The data I have been given resembles:

Code:

 
Index      Chr       Genotype   Mutation Type
1           Chr1           TT            Intronic
2           Chr1           AA            Exonic
3           Chr1           AG            Exonic
4           Chr1           CC            Frameshift
5           Chr1           CA            Intronic
...         ...              ....

My goal here is to compare a large file such as this one to a set of references that are a single letter (T,A,C,G). I need to split the genotype in the given file into two individual letters and then compare this to the reference. If either matches the reference, then the program should move on to the next reference. If, however, neither letter is consistant, then the program should give an output of the index number and mutation type that corresponds to that data.
I am still fairly new to this so all help would be greatly appreciated.
Thanks!

Can you post sample data and desired output?

Yes sorry for not specifying more.
Say the given reference was:

Code:

 
Index 
1)T
2)A
3)C
4)G
5)C

and that the matix in my first post remains the same, indexes 3 and 4 would not have an allele in their genotype that matches that of the reference for that position. My desired output in this situation would be:

Code:

 
3   Chr1   Exonic
4   Chr1   Frameshift

I hope this explanation is more useful.
Thanks again

If I understand correctly, you will want something like this:

Code:

#!/usr/bin/perl


use strict;
use warnings;

my ($index, $chr, $geno, $mutation);
my $ref = "A";
open(FILE,"<","file.txt") or die $!;
while (<FILE>) {
	next unless $_ !~ /Index/;
	next unless $_ =~ /^(\d*)\s*([a-z]*\d)\s*([a-z]*)\s*([a-z]*)/i;
	$index = $1;
	$geno = $3;
	next unless $geno !~ /$ref/;
	print "Genotype '" . $geno . "' (Index: " . $index . ") does not match the reference type '" . $ref . "'.\n";
}

Which will result in this:

Code:

Genotype 'TT' (Index: 1) does not match the reference type 'A'.
Genotype 'CC' (Index: 4) does not match the reference type 'A'.

Obviously, I have printed my own formatted line, however you could simply use the '$_' to leave as is ( print $_ . "\n"; ).

Thank you ddreggors for your help. It is similar to that but the reference is not the same for each index. For example, for index 1, 'T' is the reference and thus the genotype, 'TT', for that index fits the requirement of having one or more alleles that matches the reference. In index 3 and 4, however, the reference alleles are 'C' and 'G' respectively, and and both of these cases, the coordinating genotype does not have an allele that matches the reference. It is in these cases that I wish for output.
Sorry if I had not specified that well

@drossy

I am not exactly sure I follow but that's ok. Fortunately for us, I do not have to understand genome reference indexes or alleles for that matter to understand logic.

Using the example I have given you, you should be able see that while I set $ref as a static value (A), you can follow the exact logic inside the while loop to pull a value or array of values from another file. If all needed values exist in "THIS" file, then I have already given you what you need.

Example:

Code:

next unless $_ =~ /^(\d*)\s*([a-z]*\d)\s*([a-z]*)\s*([a-z]*)/i;
$index = $1;
$geno = $3;

these lines do something very nice for you, namely it grabs all text and splits it into separate variables delimiting on white space (or multiple white space) characters.All text surrounded in parenthesis are "kept".

You can see that I have given friendly names to the index and genotype columns, but $2 would contain the "Chr" column and $4 would contain the "Mutation Type" column with that regex match. You can easily give them friendly names to reuse as well...


Example:

Code:

next unless $_ =~ /^(\d*)\s*([a-z]*\d)\s*([a-z]*)\s*([a-z]*)/i;
$index = $1;
$chr = $2;
$geno = $3;
$mutation = $4;

Going a bit further, if you need to split geno into 2 separate characters you can now take the $geno variable and do something like this:

Example:

Code:

my ($geno1,$geno2) = split(undef,$geno);

The framework is all here, for you to do everything you want now, but for me (or others) to give you a better solution it would require a more logical approach in explaining the problem I fear.


Maybe I am slow, but I do not see (based on your explanation) the correlation you are trying to make with these references/alleles.

You say:
What is considered a "coordinating genotype?
What are the reference alleles?
Last letter in the pair of 2 is reference?


More precisely it would be easier to phrase as:

The second letter in that column must match the first.

cheers

I still don't believe I am explaining it correctly:
The index number for the given data must match that of the reference.
Data:

Code:

Index      Chr       Genotype   Mutation Type
1           Chr1           TT            Intronic
2           Chr1           AA            Exonic
3           Chr1           AG            Exonic
4           Chr1           CC            Frameshift
5           Chr1           CA            Intronic

Reference:

Code:

1)T
2)A
3)C
4)G
5)C

For example: The genotype at position 1 is TT and the reference is T. Because atleast one of the letters in the genotype matches the reference, then the program should move on.
At position 3, the genotype is AG and the reference is C. Neither of the letters that make up the genotype match the reference, so the program should take note of this and display '3 Chr1 AG Exonic' in the output.

Thanks again

Moderator's Comments:

edit by bakunin: please use [CODE]..[/CODE]-tags ehen posting code or file content. It makes it easier to read and it stays formatted. Thank you.