I am using count_overlap() for counting the kmers from Biopython .
Does it take the reverse complement of kmer into account?
I need to count reverse complement too
Bioinformatics Stack Exchange is a question and answer site for researchers, developers, students, teachers, and end users interested in bioinformatics. It only takes a minute to sign up.Sign up to join this community
There's a lot detail missing here, but at a guess could you try,
from Bio import SeqIO from Bio.Seq import Seq for record in SeqIO.parse("ls_orchid.gbk", "genbank"): reversec = record.seq.reverse_complement()) revcount = reversec.count_overlap('AT') print (revcount + count_overlap('AT'))
count_overlap only counts the forward strand if this is the Biopython method. Not debugged the above.
Here's a C++ kmer counter I wrote for my own purposes. It can count canonical and reverse complement kmers from FASTA files up to length 31. Run
kmer-counter --help or read the README for count options.
For longer kmers, you'd want to look at an approximation tool like Jellyfish.
I would NOT recommend using biopython for a task like this. First, it's slow, but also it's really not necessary. There are several other existing tools optimised for exactly this purpose.
My favourite is KMC. It's very fast, and by default, it reports canonical kmers (Parameter
-b would turn this off). So, assuming you have KMC installed, this is how you get all the kmers from a range of occurances (1 - 10000).
# just making a directory for temprary files mkdir tmp # building a kmer database with k = 21, # 16 cores and # up to 64GB of memory. # Counting occurances between 1 to 10000, # `-fa` specifies that the input is fasta, you can also input other formats # the database name will be kmer_counts # & saving temporary data to tmp directory; kmc -k21 -t16 -m64 -ci1 -cs10000 -fa my_seq.fasta kmer_counts tmp # run kmc # extract kmers from the database in a text file kmc_dump -ci1 -cx10000 kmer_counts kmer_k21.dump