"megalodon_extras modified_bases create_motif_bed" does not return all of the positions

[ziczhang]

Hi,

I used the following command to produce a bed file of the location of CpGs in hg38.

megalodon_extras modified_bases create_motif_bed --motif CG 0 --out-filename motif_sites.bed UCSC.hg38.fa

and the motif_sites.bed started from

chr1	11460	11461	.	.	+
chr1	11461	11462	.	.	-
chr1	11479	11480	.	.	+
chr1	11480	11481	.	.	-

However, the actual CpGs in hg38 start at

chr1	10468	10469	.	.	+
chr1	10469	10470	.	.	-
chr1	10470	10471	.	.	+
chr1	10471	10472	.	.	-

I don't know if this is the intended behavior,
but this is because "megalodon_extras modified_bases create_motif_bed" ignores lowercase sequences in the fasta file.

This can be fixed by modifying

def compile_motif_pat(raw_motif):
    ambig_pat_str = "".join(
        "[{}]".format(SINGLE_LETTER_CODE[letter]) for letter in raw_motif
    )
    # add lookahead group to serach for overlapping motif hits
-    return re.compile("(?=({}))".format(ambig_pat_str))
+    return re.compile("(?=({}))".format(ambig_pat_str),flags=re.IGNORECASE)

def compile_rev_comp_motif_pat(raw_motif):
    ambig_pat_str = "".join(
        "[{}]".format("".join(comp(b) for b in SINGLE_LETTER_CODE[letter]))
        for letter in raw_motif[::-1]
    )
    # add lookahead group to serach for overlapping motif hits
-    return re.compile("(?=({}))".format(ambig_pat_str))
+    return re.compile("(?=({}))".format(ambig_pat_str),flags=re.IGNORECASE)

in megalodon_helper.py.

Best,

Zicong