<p>To scan motifs in a genome (or database) I would use <a href="http://meme-suite.org/tools/fimo" rel="nofollow noreferrer">FIMO</a> which will give you the exact locations of these motifs in your genome.</p>

<p>Once you have the locations, you can use a <code>phastCons</code> bigiwig from UCSC to calculate the basewise conservation scores. However, please remember that <code>phastCons</code> scores are smoothed across windows and it might not be the best metric if you are trying to compare the conservation levels at your motif matching sites as compared to the sequences flanking them.</p>

<p>I wrote a <a href="https://github.com/saketkc/moca" rel="nofollow noreferrer">package</a> a while back to do this, including doing de-novo motif discovery. However, it might be an overkill for your use case. </p>

CC BY-SA 3.0
2017-05-31T01:26:37.957Z
I've used FIMO before actually, but for comparing motifs against PWMs (from JASPAR). Might be a dumb question, but is it pretty easy to compare against a genome rather than a database of PWM's?
Eric Brenner
CC BY-SA 3.0
2017-05-31T01:51:00.220Z
@EricBrenner You can simply upload a sequence file in the web version.
Saket Choudhary
CC BY-SA 3.0
2017-05-30T22:46:52.913Z
2017-05-31T01:51:14.620Z
2017-05-31T01:51:14.620Z
Saket Choudhary
<p><strong>In case you have only ACGT in your motifs</strong></p>

<p>The short motifs make it sound as if you are in the business of looking for a kmer counter. You can choose to use existing software or build your own.</p>

<ol>
<li>Using existing software might be your easiest path. An older post from 2014 will probably give you a first idea what's out there: <a href="http://homolog.us/blogs/blog/2014/04/07/kmer-counting-a-2014-recap/" rel="noreferrer">http://homolog.us/blogs/blog/2014/04/07/kmer-counting-a-2014-recap/</a> . Note that a couple of algorithms mentioned there got successors, so it is worthwhile digging a bit around. The small kmer size will make most of them usable for your needs.</li>
<li>As the maximum size of your kmers is comparatively small (12 nt need 24 bits, i.e., max 16.7 million entries in your kmer table), you should be able to easily roll your own kmer counting in about any language you like and on about any of nowadays computer. The pseudocode section on the <a href="https://en.wikipedia.org/wiki/K-mer" rel="noreferrer">Wikipedia</a> entry for kmers will give you first pointers for that. Might be a bit more work, but maybe more flexible depending on your needs.</li>
</ol>

<p><strong>In case you have IUPAC bases (N, W, etc.) in your motifs</strong></p>

<p>I don't know any any pre-existing software doing what you need. I could imagine that the short motifs make using <a href="https://en.wikipedia.org/wiki/Regular_expression" rel="noreferrer">regular expressions</a> doable for this kind of search, but I may be wrong. Testing this should be easy in a simple script as all major programming languages have modules or libraries for REs. Even if it should take a couple of hours to run on your data set, that would be good enough for a one-off calculation.</p>

CC BY-SA 3.0
2017-05-30T21:53:30.923Z
2017-05-30T21:53:30.923Z
BaCh
<p>Where transcription factor binding sites are concerned, in our papers we have used per-base <a href="http://compgen.cshl.edu/phast/faq.php" rel="nofollow noreferrer">phyloP</a> conservation data instead of smoothed phastCons scores. </p>

<p>We use BEDOPS <code>bedmap</code> to map scores over multiple, padded binding sites (BED-formatted) for a given motif model. The phyloP scores are WIG files obtained from UCSC goldenpath and converted to BED via <code>wig2bed</code>.</p>

<p>The resulting matrix of binding site intervals and their per-base scores can be turned into a ranked heatmap or aggregated to determine average per-base conservation for a motif model. </p>

<p>While TF binding sites innately show high information content and therefore high conservation, matrices can be further sorted by score-maps of ChIP-seq or DNaseI-seq tag density. </p>

<p>This can help filter for interesting low-noise, high-signal patterns and is useful when the matrix is windowed and the window includes low-information regions with a few high-information residues outside of the main motif — as is the case with CTCF, for instance.</p>

CC BY-SA 3.0
2018-03-07T00:10:38.133Z
2018-03-07T00:10:38.133Z
Alex Reynolds
<p>I have about 200 short nucleotide motifs (6-12 bp in length) from the human genome, and I'm trying to see how conserved they are across vertebrates. </p>

<p>I was thinking that I'd need to make a bed file for each motif that lists all of its occurrences in the human genome. From there, I could map the beds to a bigwig files of <a href="http://compgen.cshl.edu/phast/phastCons-HOWTO.html" rel="noreferrer">PhastCons</a> scores (essentially doing the reverse of what the PhastCons software was designed to do). Does that sound like the best approach? </p>

<p>I'm getting stuck at the step of going from motifs to bed files. I've tried using BLAST to find all occurrences of motifs, but their short length is causing issues.<br>
I've tried messing with the e-value threshold, word size, and filter parameters, but I still don't get any hits. </p>

<p>Is there a work-around for this issue, or should I just rethink my entire approach?</p>

CC BY-SA 3.0
2017-05-31T10:04:07.103Z
I did something similar, but with a list of well defined n-mers instead of "motifs", all with the same value of n. In case this can be of any help, the code I used to make the bed file is the following: https://bitbucket.org/blaiseli/conservation_mirna_targets/src/20b7a0894b957e69912f637ef2d0493779036029/trie.py?at=master&fileviewer=file-view-default
The code might not be super memory efficient.
bli
<k-mer><motifs>
Have DNA motifs 6-12bp long, trying to get conservation scores
# Have DNA motifs 6-12bp long, trying to get conservation scores

(Asked by: username_0 on 2017-05-30)

I have about 200 short nucleotide motifs (6-12 bp in length) from the human genome, and I'm trying to see how conserved they are across vertebrates.

I was thinking that I'd need to make a bed file for each motif that lists all of its occurrences in the human genome. From there, I could map the beds to a bigwig files of [PhastCons](http://compgen.cshl.edu/phast/phastCons-HOWTO.html) scores (essentially doing the reverse of what the PhastCons software was designed to do). Does that sound like the best approach?

I'm getting stuck at the step of going from motifs to bed files. I've tried using BLAST to find all occurrences of motifs, but their short length is causing issues.  
I've tried messing with the e-value threshold, word size, and filter parameters, but I still don't get any hits.

Is there a work-around for this issue, or should I just rethink my entire approach?

**Tags:** k-mer, motifs
> **Comments:**
> **username_1**: I did something similar, but with a list of well defined n-mers instead of "motifs", all with the same value of n. In case this can be of any help, the code I used to make the bed file is the following: https://bitbucket.org/blaiseli/conservation_mirna_targets/src/20b7a0894b957e69912f637ef2d0493779036029/trie.py?at=master&fileviewer=file-view-default The code might not be super memory efficient.  

---

## 3 answers

### Answer by username_3 on 2017-05-30. Score: 7

**In case you have only ACGT in your motifs**

The short motifs make it sound as if you are in the business of looking for a kmer counter. You can choose to use existing software or build your own.

1. Using existing software might be your easiest path. An older post from 2014 will probably give you a first idea what's out there: <http://homolog.us/blogs/blog/2014/04/07/kmer-counting-a-2014-recap/> . Note that a couple of algorithms mentioned there got successors, so it is worthwhile digging a bit around. The small kmer size will make most of them usable for your needs.
2. As the maximum size of your kmers is comparatively small (12 nt need 24 bits, i.e., max 16.7 million entries in your kmer table), you should be able to easily roll your own kmer counting in about any language you like and on about any of nowadays computer. The pseudocode section on the [Wikipedia](https://en.wikipedia.org/wiki/K-mer) entry for kmers will give you first pointers for that. Might be a bit more work, but maybe more flexible depending on your needs.

**In case you have IUPAC bases (N, W, etc.) in your motifs**

I don't know any any pre-existing software doing what you need. I could imagine that the short motifs make using [regular expressions](https://en.wikipedia.org/wiki/Regular_expression) doable for this kind of search, but I may be wrong. Testing this should be easy in a simple script as all major programming languages have modules or libraries for REs. Even if it should take a couple of hours to run on your data set, that would be good enough for a one-off calculation.


---

### Answer by username_2 on 2017-05-30. Score: 3

To scan motifs in a genome (or database) I would use [FIMO](http://meme-suite.org/tools/fimo) which will give you the exact locations of these motifs in your genome.

Once you have the locations, you can use a `phastCons` bigiwig from UCSC to calculate the basewise conservation scores. However, please remember that `phastCons` scores are smoothed across windows and it might not be the best metric if you are trying to compare the conservation levels at your motif matching sites as compared to the sequences flanking them.

I wrote a [package](https://github.com/saketkc/moca) a while back to do this, including doing de-novo motif discovery. However, it might be an overkill for your use case.

> **Comments:**
> **username_0**: I've used FIMO before actually, but for comparing motifs against PWMs (from JASPAR). Might be a dumb question, but is it pretty easy to compare against a genome rather than a database of PWM's?  
> **username_2**: @EricBrenner You can simply upload a sequence file in the web version.  

---

### Answer by username_4 on 2018-03-07. Score: 1

Where transcription factor binding sites are concerned, in our papers we have used per-base [phyloP](http://compgen.cshl.edu/phast/faq.php) conservation data instead of smoothed phastCons scores.

We use BEDOPS `bedmap` to map scores over multiple, padded binding sites (BED-formatted) for a given motif model. The phyloP scores are WIG files obtained from UCSC goldenpath and converted to BED via `wig2bed`.

The resulting matrix of binding site intervals and their per-base scores can be turned into a ranked heatmap or aggregated to determine average per-base conservation for a motif model.

While TF binding sites innately show high information content and therefore high conservation, matrices can be further sorted by score-maps of ChIP-seq or DNaseI-seq tag density.

This can help filter for interesting low-noise, high-signal patterns and is useful when the matrix is windowed and the window includes low-information regions with a few high-information residues outside of the main motif — as is the case with CTCF, for instance.


---

