Compare two matrix

Question

[Amy]

I have two matrix. The first one is a genotype matrix in which:

• Rows represent locus
• Columns represent samples
• Each value represents a genotype which could be either P1/P1, P2/P2, P1/P2 or NA if the genotype is not determined.

The second matrix is a matrix of counts. As the first one:

• Rows represent locus
• Columns represent genotypes
• Each value represent of a count in each locus for each sample.

I would like to use the genotypic information to treat the second matrix. The aim is to replace the count value with NA when the genotype is not determined (i.e NA).

Here is an example of my two matrix:

-Genotypic Matric

CDS             BC1-III     BC1-IV      BC10-II     
LOC105031928    P1/P2       P1/P2      P1/P2    
LOC105031930    NA          NA         NA   
LOC105031931    P1/P1       P1/P1      P1/P1    
LOC105031933    P1/P1       P1/P1      P1/P1    
LOC105031934    NA          NA         NA   
LOC105031935    P1/P1       P1/P1      P1/P1    
LOC105031937    NA          NA         NA   
LOC105031938    P1/P1       P1/P1      P1/P1

-Matrix of Counts

CDS             BC1-III     BC1-IV      BC10-II     
LOC105031928    175         181.5       99
LOC105031930    10          50          0
LOC105031931    401         691         572
LOC105031933    17          69          15.75
LOC105031934    0           0           0
LOC105031935    6           0           17
LOC105031937    0           0           0
LOC105031938    408         520.1       165

What my script should give:

CDS             BC1-III     BC1-IV      BC10-II     
LOC105031928    175         181.5       99
LOC105031930    NA          NA          NA
LOC105031931    401         691         572
LOC105031933    17          69          15.75
LOC105031934    NA          NA          NA
LOC105031935    6           0           17
LOC105031937    NA          NA          NA
LOC105031938    408         520.1       165

I could read the genotypic matrix line by line and link the two matrix by their CDS as ID but i want to make sure that one value is specifi to its CDS and its sample. I am beginner to perl and by now i don't know yet how to extract from a matrix the header and row information and then assign them to one value. Thanks for your help.

PS: This is what I have done from now:

open(GENOTYPE, '<', "$matrix_geno") or die ("Cannot open $matrix_geno\n");
my %hash_Loc_line = ();
while (my $line = <GENOTYPE>)
{
    chomp $line;
    next if ($line =~ /^CDS/);
    my @columns = split (/\s+/, $line);
    my $nb_col = scalar(@columns)-1;
    my $locus = $columns[0];
    my @BC = @columns[1..$nb_col];
    foreach my $BC (@BC)
    {
        push @{$hash_Loc_line{$locus}}, $BC;
    }
}

[WouterDeCoster]

Hi Amy,

If an answer was helpful you should upvote it, if the answer resolved your question you should mark it as accepted.

[Michael]

Have you tried anything? From your description it seems you want to simply have NA in the same locations as Table 1 and values from Table 2 otherwise. That is very easy to solve in R for example.

[Amy]

Yes I'm currently working on it but still can't find the right way to solve the problem. Unfortunately i must write the script only in perl. I'am going to add in my post what I did from now. Thanks.

[Michael]

So this is an assignment? Because it would be a one-liner in R, that's a pitty.

[Amy]

Yes it is. Unfortunately :(

Answer 1

[Paul]

Hi, what about this bash solution:

paste -d '\n' matrix1 matrix2 | awk 'NR>1' | grep -v "\bP1\b" | awk '{if(! a[$1]){print; a[$1]++}}'
  

CDS             BC1-III     BC1-IV      BC10-II     
LOC105031928    175         181.5       99
LOC105031930    NA          NA         NA   
LOC105031931    401         691         572
LOC105031933    17          69          15.75
LOC105031934    NA          NA         NA   
LOC105031935    6           0           17
LOC105031937    NA          NA         NA   
LOC105031938    408         520.1       165

note: First you interleave your matrix, then you remove "double" header, then remove is there is an genotype and last you can remove duplicates lines (first allays should be "NA")

[Amy]

Thank you for your help. I tried it with my data and it worked well. I'll try to combine in with my perl script. :)

[Paul]

Your are welcome. Sorry for bash solution - I am not familiar with perl.

Answer 2

[Pierre Lindenbaum]

bash:

paste matrix1.txt matrix1.txt | awk '/^CDS/ {print $1"\t"$2"\t"$3"\t"$4; next;} {printf("%s",$1);for(i=2;i<=4;i++) {j=i+4;printf("\t%s",($i == "NA" ? "NA" : $j ));} printf("\n");}'



CDS BC1-III BC1-IV  BC10-II
LOC105031928    175 181.5   99
LOC105031930    NA  NA  NA
LOC105031931    401 691 572
LOC105031933    17  69  15.75
LOC105031934    NA  NA  NA
LOC105031935    6   0   17
LOC105031937    NA  NA  NA
LOC105031938    408 520.1   165

[Amy]

Thank you for your answer. I tried it on my axamples and it worked well. But actually the two matrix I gave as examples are small part of my whole files. In fact I have 51 samples and 27011 rows. So do you know how I can modify this code for general situations ?

Thank you.

[Pierre Lindenbaum]

In fact I have 51 samples and 27011 rows.

I'm sure you can try to modify the awk script by yourself.

[Amy]

Okay. Thanks anyway. :)