Return to Answer

UPDATE: Upon further reflection, the more complicated example you listed actually makes sense. In the reference we have the sequence AACAC, and the alternate alleles represent two variations on this: deletion of the final two bp (in the first case), and a point mutation of the middle C to T (in both cases). Usually, only a single bp precedes the definition of an indel, so I would've encoded this complex variant as ref=ACAC alt=AT,ATAC.

So the SNV is "implied by"/"encoded in"/"redundant with" the complex variant, but it's not strictly a duplicate. I'm curious whether the VCF validator complains about these cases as well?

The script below will convert these funky SNVs to the canonical representation, and will also work on MNVs. You'd still need a tool to remove the duplicates, though.

The script below will convert these funky SNVs to the canonical representation, and will also work on MNVs. Standard tools should then work to remove the duplicates.

#!/usr/bin/env python3


def canonicalize(instream):
    for line in instream:
        if not line.startswith('#'):
            values = line.split('\t')
            pos = int(values[1])
            ref, alt = values[3:5]
            if len(ref) > 1 and len(ref) == len(alt):
                offset = None
                revoffset = None
                newref, newalt = '', ''
                for n, (r, a) in enumerate(zip(ref[::-1], alt[::-1])):
                    if r != a:
                        revoffset = -1 * (n)
                        break
                for n, (r, a) in enumerate(zip(ref, alt)):
                    if r != a:
                        offset = n
                        values[1] = str(pos + offset)
                        values[3] = ref[offset:revoffset]
                        values[4] = alt[offset:revoffset]
                        break
                line = '\t'.join(values)
        yield line


if __name__ == '__main__':
    import sys
    for line in canonicalize(sys.stdin):
        print(line, end='')

#!/usr/bin/env python3


def canonicalize(instream):
    for line in instream:
        if not line.startswith('#'):
            values = line.split('\t')
            pos = int(values[1])
            ref, alt = values[3:5]
            if len(ref) > 1 and len(ref) == len(alt):
                # How many bp to trim off the end
                for n, (r, a) in enumerate(zip(ref[::-1], alt[::-1])):
                    if r != a:
                        revoffset = -1 * n
                        break

                # How many bp to trim off the front
                for n, (r, a) in enumerate(zip(ref, alt)):
                    if r != a:
                        offset = n
                        values[1] = str(pos + offset)
                        values[3] = ref[offset:revoffset]
                        values[4] = alt[offset:revoffset]
                        break
                line = '\t'.join(values)
        yield line


if __name__ == '__main__':
    import sys
    for line in canonicalize(sys.stdin):
        print(line, end='')