Re^3: Replacing substrings within hash values

[BillKSmith]

The main problem with your code is that you must read the entire IN file for each sequence. This would work if you rewind (seek IN, 0, 0) the IN file at the end of the sequence loop. This approach is very inefficient. Run-time would become unacceptable for larger files.

Your use of last is fine for sequence I. For all other sequences, it would exit the loop before it gets to the processing. Use next instead.

[K_Edw]

I see, thank you. The script is now working as desired after fixing the usage of substr, using next and the addition of seek(IN,0,0). Is there a more efficient way to accomplish this or is this an unavoidable consequence of the way the script is structured/written? My thought was that i'd only need to iterate over IN once as it's sorted and that once $F[0] no longer equalled the current key, it would begin reading IN where it left off for the next sequence before it terminated.

[poj]

Since you can access a sequence using it's key, no need to loop through them until you need to print.

#!/usr/bin/env perl
use strict;
use warnings;

my %sequences = (
 I   => 'CATCAGTATAAAATGACTAGTAGCTAGATACCACAGATACGATACAACA',
 II  => 'TACCACAGATACGATACAACACATCAGTATAAAATGACTAGTAGCAGAC',
);

while (<DATA>) {
  my @f = split(/\s+/, $_);
  substr ($sequences{$f[0]},$f[1]-1,1) = $f[3];
}

print "CGTTGGCATAAAATGACTAGTAGCTAGATACCACAGATACGATACAACA\n";
for my $key (keys %sequences){
  print $sequences{$key}."\n";
}
__DATA__
I    2    A    G
I    4    C    T
I    5    A    G
I    7    T    C
II    1    T    C
II    2    A    G
II    3    C    T
II    5    A    C
II    8    G    T
II    10    T    G

poj

[K_Edw]

Thanks! That is a lot more succinct and efficient :) I had not considered using $F0 to specify the sequence directly.

[K_Edw]

Now that i'm no longer iterating over keys%sequences, is there a way to reset a variable for each key? I now wish to add a second type of edit where letters can be deleted or inserted (rather than replaced). This will shift all remaining positions out of sync and so the changes would need to be monitored so that they can compensated for. My original idea to do this was to keep a cumulative total of the insertion/deletion sizes and adjust each remaining position accordingly - which would reset to 0 for each key.

Or would it be better to just create an array, with an entry for each key?

[BillKSmith]

I did not realize that you were trying to exploit the ordering of the edits. One way that your original code is wrong is that ignores an edit which does not belong to the current sequence rather than applying it to the next sequence. This would not be simple to fix.

I strongly recommend you change to BrowserUK's algorithm. It would be very easy to add a test to verify that the character at the position to be edited is what the edit expects. The additional effort would be paid back, the first time that it finds an example of inconsistent data.

[BrowserUk]

It would be very easy to add a test to verify that the character at the position to be edited is what the edit expects. The additional effort would be paid back, the first time that it finds an example of inconsistent data.

Whilst it would be easy to add; it is almost certainly unnecessary.

When you know how these edit lists are produced, you realise that the sequences being edited were (half of) the input to the processing that produced the edit list; thus with real data, if the sequence name/id -- which tend to look like uc002yje.1 chr21:13973492-13976330 or 32_Illumina_Multiplexing_PCR_Primer_1.01 or ceti albus: chrom 1 or SVN001-12|RMNH.ARA.14133|ANA0001|CP|M etc. -- is found in the hash, then the likeihood that the edit file will contain a different initial character at the specified position is very small indeed.

Perhaps the worst that could happen is that the post-edited sequence file could be (re)paired with the same edit file and re-run. The result would be that the entire file would be "edited", and the resulting output file would be identical to the input file. Ie. No harm done.

What my code did lack was a check for/handling of, the existence of the sequence (from the edit file) in the hash (from the sequence file), which would almost certainly indicate that the wrong edit file was being paired with the sequence file (or vice versa).

But then, my purpose was (as always) to provide the OP with the minimum demonstration that would explain the problem he was asking about -- in this case his misconception regarding 0-based and 1-based indexing -- and not production level, ready-to-run code.

That said; the addition of the consistency check would do no harm either :)