1
$\begingroup$

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

$\endgroup$
3
  • $\begingroup$ Hi @Lipsa, please provide a test example and your output together with any possible error messages. This appears certainly solvable. Is this method in Biopython? $\endgroup$
    – M__
    Nov 21 '20 at 9:52
  • $\begingroup$ I've had a guess at what I think the question should read $\endgroup$
    – M__
    Nov 21 '20 at 10:12
  • $\begingroup$ Thank you for your response Michael . However , this approach is getting very very slow when FASTA is of very large size and kmers are in millions. Any help is appreciated. $\endgroup$
    – Lipsa
    Dec 2 '20 at 7:27
1
$\begingroup$

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'))

Essentially count_overlap only counts the forward strand if this is the Biopython method. Not debugged the above.

$\endgroup$
1
$\begingroup$

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.

$\endgroup$
0
$\begingroup$

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
$\endgroup$
2
  • $\begingroup$ I am lost here. Please help! I am using python 3. One vaex dataframe has 1million kmers and I want to check the count of these kmers in the FASTA file. How will I use KMC for this? $\endgroup$
    – Lipsa
    Dec 2 '20 at 7:25
  • $\begingroup$ Yes, KMC can count kmer coverages directly from a fasta file. Just follow the code I posted there on how to use KMC to get kmers and their corresponding coverages. $\endgroup$
    – Kamil S Jaron
    Dec 2 '20 at 10:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.