I have a small bioinformatics problem that I think should be easy to solve. Related to "genotype phasing". But I'm not sure how to tackle it. In the extract below, the first column are identifiers, the subsequent columns are binary genotypes labelled with "a" or "b". "-" means missing value.
Si_gnF.scaffold10533.53688bp_tag414456 b a a b b a b a a a b a b b a b a a b b a a b b
Si_gnF.scaffold10533.76297bp_tag414484 a b b a a b a b b b a b a a b a b b a - b b a a
Si_gnF.scaffold10533.98416bp_tag414526 a b b a a b a b b b a b a a b a b b a a b b a a
Si_gnF.scaffold10534.48805bp_tag414546 b a a b a b a b b b b b b a a a a b a b b b b a
Si_gnF.scaffold10535.1091787bp_tag414684 a a a b b a a a b a b a a a a b b b a a b b a a
Si_gnF.scaffold10535.1151107bp_tag414765 b b b a a b b b a b a - b b b a a a b b a a b b
Si_gnF.scaffold10535.1220879bp_tag414877 a a a b b a a a b a b a a a a b b b a a b b a a
Si_gnF.scaffold10535.1304464bp_tag414988 b b b a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1347462bp_tag415047 b b b a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1379804bp_tag415090 b b b a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1540335bp_tag415345 a a a b b a a a b a b a a a a b b b a a b b a a
Si_gnF.scaffold10535.1585442bp_tag415410 a a a b b a a a b a b a a a a b b b a a b b a a
Si_gnF.scaffold10535.1609908bp_tag415431 b b b a a b a b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1711158bp_tag415567 b b b a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1744394bp_tag415609 b b b a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1751886bp_tag415620 a a a b b a a a b a b a a a a b b b a a b b a a
Si_gnF.scaffold10535.1752774bp_tag415622 a a a b b a a a b a b a a a a b b b a a b b a a
Si_gnF.scaffold10535.1789478bp_tag415675 b b - a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1800135bp_tag415687 b b b a a b b b a b a b b b b a a a b b a a b b
Si_gnF.scaffold10535.1885424bp_tag415814 a a a b b a a a b a b a a a a b b b a a b b a a
Basically, I want to minimize the number of differences between lines. (I cannot edit individual columns, but can flip the labels on whole lines). The result for the first four lines would be this:
Si_gnF.scaffold10533.53688bp_tag414456 b a a b b a b a a a b a b b a b a a b b a a b b
Si_gnF.scaffold10533.76297bp_tag414484 b a a b b a b a a a b a b b a b a a b - a a b b <-- this one flipped
Si_gnF.scaffold10533.98416bp_tag414526 b a a b b a b a a a b a b b a b a a b b a a b b <-- this one flipped
Si_gnF.scaffold10533.53688bp_tag414456 b a a b b a b a a a b a b b a b a a b b a a b b
As a first step I'll need to make pairwise comparisons. But what is a good way of quantifying the differences, so that I know for which lines labels must be flipped? (2 consecutive lines rarely match 100%; there can be multiple (even many) mismatches as well as missing values).
(ideally in ruby or R)
You can use the Levenshtein algorithm to quantify the difference between two strings. One way to do it:
require 'text' # See http://rubygems.org/gems/text
lines # => a array with each line
def compare(line1, line2)
Text::Levenshtein.distance(line1.sub(/.*\s/, '').sort,
line2.sub(/.*\s/, '').sort)
end
compare(lines[0], lines[1]) # => 1 (one value different)
(If "a b a a" is not equal to "a a a b", remove the sort from the method.)