***** Acces problem to record *****\
Identifieur interne :
000F93 ( Pmc/Corpus );
précédent :
000F929;
suivant :
000F940
***** probable Xml problem with record *****
Links to Exploration step
Le document en format XML
<record><TEI><teiHeader><fileDesc><titleStmt><title xml:lang="en">Enhanced Regulatory Sequence Prediction Using Gapped <italic>k</italic>
-mer Features</title><author><name sortKey="Ghandi, Mahmoud" sort="Ghandi, Mahmoud" uniqKey="Ghandi M" first="Mahmoud" last="Ghandi">Mahmoud Ghandi</name>
<affiliation><nlm:aff id="aff1"><addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><author><name sortKey="Lee, Dongwon" sort="Lee, Dongwon" uniqKey="Lee D" first="Dongwon" last="Lee">Dongwon Lee</name>
<affiliation><nlm:aff id="aff1"><addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><author><name sortKey="Mohammad Noori, Morteza" sort="Mohammad Noori, Morteza" uniqKey="Mohammad Noori M" first="Morteza" last="Mohammad-Noori">Morteza Mohammad-Noori</name>
<affiliation><nlm:aff id="aff2"><addr-line>School of Mathematics, Statistics and Computer Science, University of Tehran, Tehran, Iran</addr-line>
</nlm:aff><affiliation><nlm:aff id="aff3"><addr-line>School of Computer Science, Institute for Research in Fundamental Sciences (IPM), Tehran, Iran</addr-line>
</nlm:aff><author><name sortKey="Beer, Michael A" sort="Beer, Michael A" uniqKey="Beer M" first="Michael A." last="Beer">Michael A. Beer</name>
<affiliation><nlm:aff id="aff1"><addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><affiliation><nlm:aff id="aff4"><addr-line>McKusick-Nathans Institute of Genetic Medicine, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><publicationStmt><idno type="wicri:source">PMC</idno>
<idno type="pmid">25033408</idno>
<idno type="pmc">4102394</idno>
<idno type="url">http://www.ncbi.nlm.nih.gov/pmc/articles/PMC4102394</idno>
<idno type="RBID">PMC:4102394</idno>
<idno type="doi">10.1371/journal.pcbi.1003711</idno>
<date when="2014">2014</date>
<idno type="wicri:Area/Pmc/Corpus">000F93</idno>
<idno type="wicri:explorRef" wicri:stream="Pmc" wicri:step="Corpus" wicri:corpus="PMC">000F93</idno>
</publicationStmt><sourceDesc><biblStruct><analytic><title xml:lang="en" level="a" type="main">Enhanced Regulatory Sequence Prediction Using Gapped <italic>k</italic>
-mer Features</title><author><name sortKey="Ghandi, Mahmoud" sort="Ghandi, Mahmoud" uniqKey="Ghandi M" first="Mahmoud" last="Ghandi">Mahmoud Ghandi</name>
<affiliation><nlm:aff id="aff1"><addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><author><name sortKey="Lee, Dongwon" sort="Lee, Dongwon" uniqKey="Lee D" first="Dongwon" last="Lee">Dongwon Lee</name>
<affiliation><nlm:aff id="aff1"><addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><author><name sortKey="Mohammad Noori, Morteza" sort="Mohammad Noori, Morteza" uniqKey="Mohammad Noori M" first="Morteza" last="Mohammad-Noori">Morteza Mohammad-Noori</name>
<affiliation><nlm:aff id="aff2"><addr-line>School of Mathematics, Statistics and Computer Science, University of Tehran, Tehran, Iran</addr-line>
</nlm:aff><affiliation><nlm:aff id="aff3"><addr-line>School of Computer Science, Institute for Research in Fundamental Sciences (IPM), Tehran, Iran</addr-line>
</nlm:aff><author><name sortKey="Beer, Michael A" sort="Beer, Michael A" uniqKey="Beer M" first="Michael A." last="Beer">Michael A. Beer</name>
<affiliation><nlm:aff id="aff1"><addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><affiliation><nlm:aff id="aff4"><addr-line>McKusick-Nathans Institute of Genetic Medicine, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</nlm:aff><series><title level="j">PLoS Computational Biology</title>
<idno type="ISSN">1553-734X</idno>
<idno type="eISSN">1553-7358</idno>
<imprint><date when="2014">2014</date>
</imprint><profileDesc><textClass></textClass>
</profileDesc><front><div type="abstract" xml:lang="en"><title>Abstract</title>
<p>Oligomers of length <italic>k</italic>
, or <italic>k</italic>
-mers, are convenient and widely used features for modeling the properties and functions of DNA and protein sequences. However, <italic>k</italic>
-mers suffer from the inherent limitation that if the parameter <italic>k</italic>
is increased to resolve longer features, the probability of observing any specific <italic>k</italic>
-mer becomes very small, and <italic>k-</italic>
mer counts approach a binary variable, with most <italic>k</italic>
-mers absent and a few present once. Thus, any statistical learning approach using <italic>k</italic>
-mers as features becomes susceptible to noisy training set <italic>k</italic>
-mer frequencies once <italic>k</italic>
becomes large. To address this problem, we introduce alternative feature sets using gapped <italic>k</italic>
-mers, a new classifier, gkm-SVM, and a general method for robust estimation of <italic>k</italic>
-mer frequencies. To make the method applicable to large-scale genome wide applications, we develop an efficient tree data structure for computing the kernel matrix. We show that compared to our original kmer-SVM and alternative approaches, our gkm-SVM predicts functional genomic regulatory elements and tissue specific enhancers with significantly improved accuracy, increasing the precision by up to a factor of two. We then show that gkm-SVM consistently outperforms kmer-SVM on human ENCODE ChIP-seq datasets, and further demonstrate the general utility of our method using a Naïve-Bayes classifier. Although developed for regulatory sequence analysis, these methods can be applied to any sequence classification problem.</p><back><div1 type="bibliography"><listBibl><biblStruct><analytic><author><name sortKey="Manolio, Ta" uniqKey="Manolio T">TA Manolio</name>
</author><biblStruct><analytic><author><name sortKey="Maurano, Mt" uniqKey="Maurano M">MT Maurano</name>
</author><author><name sortKey="Humbert, R" uniqKey="Humbert R">R Humbert</name>
</author><author><name sortKey="Rynes, E" uniqKey="Rynes E">E Rynes</name>
</author><author><name sortKey="Thurman, Re" uniqKey="Thurman R">RE Thurman</name>
</author><author><name sortKey="Haugen, E" uniqKey="Haugen E">E Haugen</name>
</author><biblStruct><analytic><author><name sortKey="Lee, D" uniqKey="Lee D">D Lee</name>
</author><author><name sortKey="Karchin, R" uniqKey="Karchin R">R Karchin</name>
</author><author><name sortKey="Beer, Ma" uniqKey="Beer M">MA Beer</name>
</author><biblStruct><analytic><author><name sortKey="Stormo, Gd" uniqKey="Stormo G">GD Stormo</name>
</author><biblStruct><analytic><author><name sortKey="Beer, Ma" uniqKey="Beer M">MA Beer</name>
</author><author><name sortKey="Tavazoie, S" uniqKey="Tavazoie S">S Tavazoie</name>
</author><biblStruct><analytic><author><name sortKey="Fletez Brant, C" uniqKey="Fletez Brant C">C Fletez-Brant</name>
</author><author><name sortKey="Lee, D" uniqKey="Lee D">D Lee</name>
</author><author><name sortKey="Mccallion, As" uniqKey="Mccallion A">AS McCallion</name>
</author><author><name sortKey="Beer, Ma" uniqKey="Beer M">MA Beer</name>
</author><biblStruct></biblStruct>
<biblStruct><analytic><author><name sortKey="Leslie, C" uniqKey="Leslie C">C Leslie</name>
</author><author><name sortKey="Eskin, E" uniqKey="Eskin E">E Eskin</name>
</author><author><name sortKey="Cohen, A" uniqKey="Cohen A">A Cohen</name>
</author><author><name sortKey="Weston, J" uniqKey="Weston J">J Weston</name>
</author><author><name sortKey="Noble, Ws" uniqKey="Noble W">WS Noble</name>
</author><biblStruct><analytic><author><name sortKey="Leslie, C" uniqKey="Leslie C">C Leslie</name>
</author><author><name sortKey="Kuang, R" uniqKey="Kuang R">R Kuang</name>
</author><biblStruct><analytic><author><name sortKey="Agius, P" uniqKey="Agius P">P Agius</name>
</author><author><name sortKey="Arvey, A" uniqKey="Arvey A">A Arvey</name>
</author><author><name sortKey="Chang, W" uniqKey="Chang W">W Chang</name>
</author><author><name sortKey="Noble, Ws" uniqKey="Noble W">WS Noble</name>
</author><author><name sortKey="Leslie, C" uniqKey="Leslie C">C Leslie</name>
</author><biblStruct><analytic><author><name sortKey="Arvey, A" uniqKey="Arvey A">A Arvey</name>
</author><author><name sortKey="Agius, P" uniqKey="Agius P">P Agius</name>
</author><author><name sortKey="Noble, Ws" uniqKey="Noble W">WS Noble</name>
</author><author><name sortKey="Leslie, C" uniqKey="Leslie C">C Leslie</name>
</author><biblStruct><analytic><author><name sortKey="Amanchy, R" uniqKey="Amanchy R">R Amanchy</name>
</author><author><name sortKey="Kandasamy, K" uniqKey="Kandasamy K">K Kandasamy</name>
</author><author><name sortKey="Mathivanan, S" uniqKey="Mathivanan S">S Mathivanan</name>
</author><author><name sortKey="Periaswamy, B" uniqKey="Periaswamy B">B Periaswamy</name>
</author><author><name sortKey="Reddy, R" uniqKey="Reddy R">R Reddy</name>
</author><biblStruct><analytic><author><name sortKey="Mcdaniell, R" uniqKey="Mcdaniell R">R McDaniell</name>
</author><author><name sortKey="Lee, B K" uniqKey="Lee B">B-K Lee</name>
</author><author><name sortKey="Song, L" uniqKey="Song L">L Song</name>
</author><author><name sortKey="Liu, Z" uniqKey="Liu Z">Z Liu</name>
</author><author><name sortKey="Boyle, Ap" uniqKey="Boyle A">AP Boyle</name>
</author><biblStruct><analytic><author><name sortKey="Kim, Th" uniqKey="Kim T">TH Kim</name>
</author><author><name sortKey="Abdullaev, Zk" uniqKey="Abdullaev Z">ZK Abdullaev</name>
</author><author><name sortKey="Smith, Ad" uniqKey="Smith A">AD Smith</name>
</author><author><name sortKey="Ching, Ka" uniqKey="Ching K">KA Ching</name>
</author><author><name sortKey="Loukinov, Di" uniqKey="Loukinov D">DI Loukinov</name>
</author><biblStruct><analytic><author><name sortKey="Visel, A" uniqKey="Visel A">A Visel</name>
</author><author><name sortKey="Blow, Mj" uniqKey="Blow M">MJ Blow</name>
</author><author><name sortKey="Zhang, T" uniqKey="Zhang T">T Zhang</name>
</author><author><name sortKey="Akiyama, Ja" uniqKey="Akiyama J">JA Akiyama</name>
</author><author><name sortKey="Holt, A" uniqKey="Holt A">A Holt</name>
</author><biblStruct><analytic><author><name sortKey="Zhang, Y" uniqKey="Zhang Y">Y Zhang</name>
</author><author><name sortKey="Liu, T" uniqKey="Liu T">T Liu</name>
</author><author><name sortKey="Meyer, Ca" uniqKey="Meyer C">CA Meyer</name>
</author><author><name sortKey="Eeckhoute, J" uniqKey="Eeckhoute J">J Eeckhoute</name>
</author><author><name sortKey="Johnson, Ds" uniqKey="Johnson D">DS Johnson</name>
</author><biblStruct><analytic><author><name sortKey="Gerstein, Mb" uniqKey="Gerstein M">MB Gerstein</name>
</author><author><name sortKey="Kundaje, A" uniqKey="Kundaje A">A Kundaje</name>
</author><author><name sortKey="Hariharan, M" uniqKey="Hariharan M">M Hariharan</name>
</author><author><name sortKey="Landt, Sg" uniqKey="Landt S">SG Landt</name>
</author><author><name sortKey="Yan, K K" uniqKey="Yan K">K-K Yan</name>
</author><biblStruct><analytic><author><name sortKey="Wang, J" uniqKey="Wang J">J Wang</name>
</author><author><name sortKey="Zhuang, J" uniqKey="Zhuang J">J Zhuang</name>
</author><author><name sortKey="Iyer, S" uniqKey="Iyer S">S Iyer</name>
</author><author><name sortKey="Lin, X" uniqKey="Lin X">X Lin</name>
</author><author><name sortKey="Whitfield, Tw" uniqKey="Whitfield T">TW Whitfield</name>
</author><biblStruct><analytic><author><name sortKey="Parelho, V" uniqKey="Parelho V">V Parelho</name>
</author><author><name sortKey="Hadjur, S" uniqKey="Hadjur S">S Hadjur</name>
</author><author><name sortKey="Spivakov, M" uniqKey="Spivakov M">M Spivakov</name>
</author><author><name sortKey="Leleu, M" uniqKey="Leleu M">M Leleu</name>
</author><author><name sortKey="Sauer, S" uniqKey="Sauer S">S Sauer</name>
</author><biblStruct><analytic><author><name sortKey="Gorkin, Du" uniqKey="Gorkin D">DU Gorkin</name>
</author><author><name sortKey="Lee, D" uniqKey="Lee D">D Lee</name>
</author><author><name sortKey="Reed, X" uniqKey="Reed X">X Reed</name>
</author><author><name sortKey="Fletez Brant, C" uniqKey="Fletez Brant C">C Fletez-Brant</name>
</author><author><name sortKey="Bessling, Sl" uniqKey="Bessling S">SL Bessling</name>
</author><biblStruct><analytic><author><name sortKey="Machanick, P" uniqKey="Machanick P">P Machanick</name>
</author><author><name sortKey="Bailey, Tl" uniqKey="Bailey T">TL Bailey</name>
</author><biblStruct><analytic><author><name sortKey="Pevny, L" uniqKey="Pevny L">L Pevny</name>
</author><author><name sortKey="Simon, Mc" uniqKey="Simon M">MC Simon</name>
</author><author><name sortKey="Robertson, E" uniqKey="Robertson E">E Robertson</name>
</author><author><name sortKey="Klein, Wh" uniqKey="Klein W">WH Klein</name>
</author><author><name sortKey="Tsai, S F" uniqKey="Tsai S">S-F Tsai</name>
</author><biblStruct><analytic><author><name sortKey="Lelievre, E" uniqKey="Lelievre E">E Lelièvre</name>
</author><author><name sortKey="Lionneton, F" uniqKey="Lionneton F">F Lionneton</name>
</author><author><name sortKey="Soncin, F" uniqKey="Soncin F">F Soncin</name>
</author><author><name sortKey="Vandenbunder, B" uniqKey="Vandenbunder B">B Vandenbunder</name>
</author><biblStruct><analytic><author><name sortKey="Leslie, C" uniqKey="Leslie C">C Leslie</name>
</author><author><name sortKey="Eskin, E" uniqKey="Eskin E">E Eskin</name>
</author><author><name sortKey="Noble, Ws" uniqKey="Noble W">WS Noble</name>
</author><biblStruct><analytic><author><name sortKey="Berger, Mf" uniqKey="Berger M">MF Berger</name>
</author><author><name sortKey="Philippakis, Aa" uniqKey="Philippakis A">AA Philippakis</name>
</author><author><name sortKey="Qureshi, Am" uniqKey="Qureshi A">AM Qureshi</name>
</author><author><name sortKey="He, Fs" uniqKey="He F">FS He</name>
</author><author><name sortKey="Estep, Pw" uniqKey="Estep P">PW Estep</name>
</author><biblStruct><analytic><author><name sortKey="Weirauch, Mt" uniqKey="Weirauch M">MT Weirauch</name>
</author><author><name sortKey="Cote, A" uniqKey="Cote A">A Cote</name>
</author><author><name sortKey="Norel, R" uniqKey="Norel R">R Norel</name>
</author><author><name sortKey="Annala, M" uniqKey="Annala M">M Annala</name>
</author><author><name sortKey="Zhao, Y" uniqKey="Zhao Y">Y Zhao</name>
</author><biblStruct></biblStruct>
<biblStruct><analytic><author><name sortKey="Sandberg, R" uniqKey="Sandberg R">R Sandberg</name>
</author><author><name sortKey="Winberg, G" uniqKey="Winberg G">G Winberg</name>
</author><author><name sortKey="Br Nden, C I" uniqKey="Br Nden C">C-I Bränden</name>
</author><author><name sortKey="Kaske, A" uniqKey="Kaske A">A Kaske</name>
</author><author><name sortKey="Ernberg, I" uniqKey="Ernberg I">I Ernberg</name>
</author><biblStruct></biblStruct>
<biblStruct></biblStruct>
<biblStruct></biblStruct>
<biblStruct><analytic><author><name sortKey="Sonnenburg, S" uniqKey="Sonnenburg S">S Sonnenburg</name>
</author><author><name sortKey="R Tsch, G" uniqKey="R Tsch G">G Rätsch</name>
</author><author><name sortKey="Henschel, S" uniqKey="Henschel S">S Henschel</name>
</author><author><name sortKey="Widmer, C" uniqKey="Widmer C">C Widmer</name>
</author><author><name sortKey="Behr, J" uniqKey="Behr J">J Behr</name>
</author><biblStruct></biblStruct>
<biblStruct><analytic><author><name sortKey="Jaakkola, T" uniqKey="Jaakkola T">T Jaakkola</name>
</author><author><name sortKey="Diekhans, M" uniqKey="Diekhans M">M Diekhans</name>
</author><author><name sortKey="Haussler, D" uniqKey="Haussler D">D Haussler</name>
</author><biblStruct></biblStruct>
<biblStruct><analytic><author><name sortKey="Gupta, S" uniqKey="Gupta S">S Gupta</name>
</author><author><name sortKey="Stamatoyannopoulos, Ja" uniqKey="Stamatoyannopoulos J">JA Stamatoyannopoulos</name>
</author><author><name sortKey="Bailey, Tl" uniqKey="Bailey T">TL Bailey</name>
</author><author><name sortKey="Noble, Ws" uniqKey="Noble W">WS Noble</name>
</author><biblStruct><analytic><author><name sortKey="Sing, T" uniqKey="Sing T">T Sing</name>
</author><author><name sortKey="Sander, O" uniqKey="Sander O">O Sander</name>
</author><author><name sortKey="Beerenwinkel, N" uniqKey="Beerenwinkel N">N Beerenwinkel</name>
</author><author><name sortKey="Lengauer, T" uniqKey="Lengauer T">T Lengauer</name>
</author><biblStruct><analytic><author><name sortKey="Bryne, Jc" uniqKey="Bryne J">JC Bryne</name>
</author><author><name sortKey="Valen, E" uniqKey="Valen E">E Valen</name>
</author><author><name sortKey="Tang, M He" uniqKey="Tang M">M-HE Tang</name>
</author><author><name sortKey="Marstrand, T" uniqKey="Marstrand T">T Marstrand</name>
</author><author><name sortKey="Winther, O" uniqKey="Winther O">O Winther</name>
</author><pmc article-type="research-article"><pmc-dir>properties open_access</pmc-dir>
<front><journal-meta><journal-id journal-id-type="nlm-ta">PLoS Comput Biol</journal-id>
<journal-id journal-id-type="iso-abbrev">PLoS Comput. Biol</journal-id>
<journal-id journal-id-type="publisher-id">plos</journal-id>
<journal-id journal-id-type="pmc">ploscomp</journal-id>
<journal-title-group><journal-title>PLoS Computational Biology</journal-title>
</journal-title-group><issn pub-type="ppub">1553-734X</issn>
<issn pub-type="epub">1553-7358</issn>
<publisher><publisher-name>Public Library of Science</publisher-name>
<publisher-loc>San Francisco, USA</publisher-loc>
</publisher><article-meta><article-id pub-id-type="pmid">25033408</article-id>
<article-id pub-id-type="pmc">4102394</article-id>
<article-id pub-id-type="publisher-id">PCOMPBIOL-D-13-01373</article-id>
<article-id pub-id-type="doi">10.1371/journal.pcbi.1003711</article-id>
<article-categories><subj-group subj-group-type="heading"><subject>Research Article</subject>
</subj-group><subj-group subj-group-type="Discipline-v2"><subject>Biology and Life Sciences</subject>
<subj-group><subject>Computational Biology</subject>
</subj-group><subj-group><subject>Genetics</subject>
<subj-group><subject>Genomics</subject>
<subj-group><subject>Functional Genomics</subject>
</subj-group><subj-group><subject>Molecular Biology</subject>
<subj-group><subject>Molecular Biology Techniques</subject>
<subj-group><subject>Sequencing Techniques</subject>
<subj-group><subject>Sequence Analysis</subject>
</subj-group><subj-group subj-group-type="Discipline-v2"><subject>Computer and Information Sciences</subject>
<subj-group><subject>Network Analysis</subject>
<subj-group><subject>Regulatory Networks</subject>
</subj-group><title-group><article-title>Enhanced Regulatory Sequence Prediction Using Gapped <italic>k</italic>
-mer Features</article-title><alt-title alt-title-type="running-head">Gapped <italic>k</italic>
-mer SVM for Regulatory Elements</alt-title><contrib-group><contrib contrib-type="author" equal-contrib="yes"><name><surname>Ghandi</surname>
<given-names>Mahmoud</given-names>
</name><xref ref-type="aff" rid="aff1"><sup>1</sup>
</xref><xref ref-type="author-notes" rid="fn1"><sup>¤</sup>
</xref><contrib contrib-type="author" equal-contrib="yes"><name><surname>Lee</surname>
<given-names>Dongwon</given-names>
</name><xref ref-type="aff" rid="aff1"><sup>1</sup>
</xref><contrib contrib-type="author"><name><surname>Mohammad-Noori</surname>
<given-names>Morteza</given-names>
</name><xref ref-type="aff" rid="aff2"><sup>2</sup>
</xref><xref ref-type="aff" rid="aff3"><sup>3</sup>
</xref><contrib contrib-type="author"><name><surname>Beer</surname>
<given-names>Michael A.</given-names>
</name><xref ref-type="aff" rid="aff1"><sup>1</sup>
</xref><xref ref-type="aff" rid="aff4"><sup>4</sup>
</xref><xref ref-type="corresp" rid="cor1"><sup>*</sup>
</xref><aff id="aff1"><label>1</label>
<addr-line>Department of Biomedical Engineering, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</aff><aff id="aff2"><label>2</label>
<addr-line>School of Mathematics, Statistics and Computer Science, University of Tehran, Tehran, Iran</addr-line>
</aff><aff id="aff3"><label>3</label>
<addr-line>School of Computer Science, Institute for Research in Fundamental Sciences (IPM), Tehran, Iran</addr-line>
</aff><aff id="aff4"><label>4</label>
<addr-line>McKusick-Nathans Institute of Genetic Medicine, Johns Hopkins University, Baltimore, Maryland, United States of America</addr-line>
</aff><contrib-group><contrib contrib-type="editor"><name><surname>Morris</surname>
<given-names>Quaid</given-names>
</name><role>Editor</role>
<xref ref-type="aff" rid="edit1"></xref>
</contrib><aff id="edit1"><addr-line>University of Toronto, Canada</addr-line>
</aff><author-notes><corresp id="cor1">* E-mail: <email>mbeer@jhu.edu</email>
</corresp><fn fn-type="COI-statement"><p>The authors have declared that no competing interests exist.</p>
</fn><fn fn-type="con"><p>Conceived and designed the experiments: MG DL MAB. Performed the experiments: MG DL. Analyzed the data: MG DL. Contributed reagents/ materials/analysis tools: MG DL MMN. Wrote the paper: MG DL MAB.</p>
</fn><fn id="fn1" fn-type="current-aff"><label>¤</label>
<p>Current address: Broad Institute of MIT and Harvard, Cambridge, Massachusetts, United States of America</p>
</fn><pub-date pub-type="collection"><month>7</month>
<year>2014</year>
</pub-date><pub-date pub-type="epub"><day>17</day>
<month>7</month>
<year>2014</year>
</pub-date><volume>10</volume>
<issue>7</issue>
<elocation-id>e1003711</elocation-id>
<history><date date-type="received"><day>1</day>
<month>8</month>
<year>2013</year>
</date><date date-type="accepted"><day>28</day>
<month>5</month>
<year>2014</year>
</date><permissions><copyright-statement>© 2014 Ghandi et al</copyright-statement>
<copyright-year>2014</copyright-year>
<copyright-holder>Ghandi et al</copyright-holder>
<license xlink:href="http://creativecommons.org/licenses/by/4.0/"><license-p>This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are properly credited.</license-p>
</license><abstract><title>Abstract</title>
<p>Oligomers of length <italic>k</italic>
, or <italic>k</italic>
-mers, are convenient and widely used features for modeling the properties and functions of DNA and protein sequences. However, <italic>k</italic>
-mers suffer from the inherent limitation that if the parameter <italic>k</italic>
is increased to resolve longer features, the probability of observing any specific <italic>k</italic>
-mer becomes very small, and <italic>k-</italic>
mer counts approach a binary variable, with most <italic>k</italic>
-mers absent and a few present once. Thus, any statistical learning approach using <italic>k</italic>
-mers as features becomes susceptible to noisy training set <italic>k</italic>
-mer frequencies once <italic>k</italic>
becomes large. To address this problem, we introduce alternative feature sets using gapped <italic>k</italic>
-mers, a new classifier, gkm-SVM, and a general method for robust estimation of <italic>k</italic>
-mer frequencies. To make the method applicable to large-scale genome wide applications, we develop an efficient tree data structure for computing the kernel matrix. We show that compared to our original kmer-SVM and alternative approaches, our gkm-SVM predicts functional genomic regulatory elements and tissue specific enhancers with significantly improved accuracy, increasing the precision by up to a factor of two. We then show that gkm-SVM consistently outperforms kmer-SVM on human ENCODE ChIP-seq datasets, and further demonstrate the general utility of our method using a Naïve-Bayes classifier. Although developed for regulatory sequence analysis, these methods can be applied to any sequence classification problem.</p><abstract abstract-type="summary"><title>Author Summary</title>
<p>Genomic regulatory elements (enhancers, promoters, and insulators) control the expression of their target genes and are widely believed to play a key role in human development and disease by altering protein concentrations. A fundamental step in understanding enhancers is the development of DNA sequence-based models to predict the tissue specific activity of regulatory elements. Such models facilitate both the identification of the molecular pathways which impinge on enhancer activity through direct transcription factor binding, and the direct evaluation of the impact of specific common or rare genetic variants on enhancer function. We have previously developed a successful sequence-based model for enhancer prediction using a <italic>k</italic>
-mer support vector machine (kmer-SVM). Here, we address a significant limitation of the kmer-SVM approach and present an alternative method using gapped k-mers (gkm-SVM) which exhibits dramatically improved accuracy in all test cases. While we focus on enhancers and transcription factor binding, our method can be applied to improve a much broader class of sequence analysis problems, including proteins and RNA. In addition, we expect that most <italic>k</italic>
-mer based methods can be significantly improved by simply using the generalized <italic>k</italic>
-mer count method that we present in this paper. We believe this improved model will enable significant contributions to our understanding of the human regulatory system.</p><funding-group><funding-statement>Searle Scholars Program NIH R01 HG007348. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.</funding-statement>
</funding-group><counts><page-count count="15"></page-count>
</counts><body><sec id="s1"><title>Introduction</title>
<p>Predicting the function of regulatory elements from primary DNA sequence still remains a major problem in computational biology. These elements typically contain combinations of several binding sites for regulatory factors whose activity together specifies the developmental times, cell-types, or environmental signals in which the element will be active. Genetic variation in regulatory elements is increasingly thought to play a significant role in the etiology and heritability of common diseases, and surveys of Genome Wide Association Studies have highlighted the preponderance of significant variants in regulatory DNA <xref rid="pcbi.1003711-Manolio1" ref-type="bibr">[1]</xref>
,<xref rid="pcbi.1003711-Maurano1" ref-type="bibr">[2]</xref>
. An accurate computational model to predict regulatory elements can 1) help identify and link core sets of regulatory factors with specific diseases, and 2) predict the functional consequences of variation or mutations in specific sites within regulatory elements.</p><p>We have recently introduced a successful method for regulatory DNA sequence prediction, kmer-SVM, which uses combinations of short (6–8 bp) <italic>k</italic>
-mer frequencies to predict the activity of larger functional genomic sequence elements, typically ranging from 500 to 2000 bp in length <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
. An advantage of <italic>k</italic>
-mer based approaches relative to the alternative position weight matrix (PWM) approach is that PWMs can require large amounts of data to optimize and determine appropriate scoring thresholds <xref rid="pcbi.1003711-Stormo1" ref-type="bibr">[4]</xref>
,<xref rid="pcbi.1003711-Beer1" ref-type="bibr">[5]</xref>
, while <italic>k</italic>
-mers are simple features which are either present or absent. However, in our previous implementation of the kmer-SVM <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
, the choice to use a single <italic>k</italic>
, and which <italic>k</italic>
, is somewhat arbitrary and based on performance on a limited selection of datasets. A major contribution of the present work is an extension of this single <italic>k</italic>
approach to include longer and much more general sequence features. The function of these DNA regulatory elements is generally thought to be specified at the molecular level by the binding of combinations of Transcription Factors (TFs) or other DNA binding regulatory factors, and many of these binding sites are short and fall within the range of <italic>k</italic>
(6–8) where our kmer-SVM approach was successful. However, Transcription Factor Binding Sites (TFBS) can vary from 6–20 bp, so some are much longer (such as ABF1, CTCF, etc.), and thus cannot be completely represented by the short <italic>k</italic>
-mers. Alternatively, TFBS can be defined by a set of sequences with some gaps (non-informative positions) as each given DNA sequence has some binding affinity for the TF. Although the kmer-SVM method can model TFBS longer than <italic>k</italic>
by tiling across TFBS with overlapping <italic>k</italic>
-mers, this loses some spatial information in the binding site, and overall classification accuracy can be significantly impaired when long TFBS are important predictive features <xref rid="pcbi.1003711-FletezBrant1" ref-type="bibr">[6]</xref>
.</p><p>Naively one could address this issue by using longer <italic>k</italic>
's or combinations of <italic>k-</italic>
mers spanning the expected size range of TFBS, but a major limitation of this approach is that longer <italic>k</italic>
-mers generate extremely sparse feature vectors (i.e. most <italic>k</italic>
-mers simply do not appear in a training sequence and thus receive zero counts, or appear only once), which causes a severe overfitting problem even at quite moderate <italic>k</italic>
. Therefore, the original kmer-SVM was limited in practice to <italic>k</italic>
-mer lengths from 6 to 10, with performance already degrading at <italic>k</italic>
= 9 or 10, depending on the dataset. Thus in practice, the parameter <italic>k</italic>
was chosen by a tradeoff between resolving longer features and robust estimation of their frequencies.</p><p>We recently introduced gapped <italic>k-</italic>
mers as a way to resolve this fundamental limitation with <italic>k</italic>
-mer features and showed that they can be used to more robustly estimate <italic>k-</italic>
mer frequencies in real biological sequences <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
. In this paper, we present a simple and efficient method for calculation of the robust <italic>k</italic>
-mer count estimates. We also expand our kmer-SVM method <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
to use gapped <italic>k</italic>
-mers or robust <italic>k</italic>
-mer count estimates as feature sets and present efficient methods to compute these new kernels. We show that our new method, gkm-SVM, consistently and significantly outperforms a kmer-SVM using both CTCF and EP300 genomic bound regions over a wide range of varying feature lengths. Furthermore, we show that, while kmer-SVM suffers significantly from overfitting as <italic>k</italic>
is increased, gkm-SVM performance is only very modestly affected by changes in the chosen feature length parameters. Next, we systematically compare the two approaches on the complete human ENCODE ChIP-seq data sets, and show that gkm-SVM either significantly outperforms or is comparable to kmer-SVM in all cases. Of biological interest, on the ENCODE ChIP-seq data sets, we also show that gkm-SVM outperforms the best known single PWM by detecting necessary co-factors. We also systematically compare gkm-SVM to similar earlier SVM approaches <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
–<xref rid="pcbi.1003711-Arvey1" ref-type="bibr">[11]</xref>
, and show that they perform comparably for optimal parameters in terms of accuracy, but that gkm-SVM is less sensitive to parameter choice and is computationally more efficient. To further demonstrate the more general utility of the <italic>k</italic>
-mer count estimates, we apply them in a simple Naïve-Bayes classifier, and show that using <italic>k</italic>
-mer count estimates instead of <italic>k</italic>
-mer counts consistently improves classification accuracy. Since our proposed method is general, we anticipate that many other sequence classification problems will also benefit from using these features. For example, word based methods can also be used to detect functional motifs in protein sequences, where the length of the functional domain is unknown <xref rid="pcbi.1003711-Amanchy1" ref-type="bibr">[12]</xref>
.</p><sec id="s2"><title>Results</title>
<sec id="s2a"><title>Calculation of sequence similarity score using gapped <italic>k</italic>
-mers</title><p>To overcome the limitations associated with using <italic>k</italic>
-mers as features described above, we introduce a new method called gkm-SVM, which uses as features a full set of <italic>k</italic>
-mers <italic>with gaps.</italic>
At the heart of most classification methods is a distance or similarity score, often called a kernel function in the SVM context, which calculates the similarity between any two elements in the chosen feature space. Therefore, in this section, we first describe the feature set and how to efficiently calculate the similarity score. This new feature set, called <italic>gapped k-mers</italic>
, is characterized by two parameters; (1) <italic>l</italic>
, the whole word length including gaps, and (2) <italic>k</italic>
, the number of informative, or non-gapped, positions in each word. The number of gaps is thus <italic>l – k</italic>
.</p><p>We first define a feature vector for a given sequence <italic>S</italic>
to be <inline-formula><inline-graphic xlink:href="pcbi.1003711.e001.jpg"></inline-graphic>
</inline-formula><italic>M</italic>
is the number of all gapped <italic>k</italic>
-mers (i.e. for DNA sequences, <inline-formula><inline-graphic xlink:href="pcbi.1003711.e002.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e003.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
-mers appeared in the sequence <italic>S</italic>
. We then define a similarity score, or a kernel function, between two sequences, <italic>S</italic>
<sub>1</sub>
and <italic>S</italic>
<sub>2</sub>
, as the normalized inner product of the corresponding feature vectors as follows:<disp-formula><graphic xlink:href="pcbi.1003711.e004.jpg" position="float" orientation="portrait"></graphic>
<label>(1)</label>
</disp-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e005.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e006.jpg"></inline-graphic>
</inline-formula><italic>K</italic>
(<italic>S</italic>
<sub>1</sub>
, <italic>S</italic>
<sub>2</sub>
), is always between 0 and 1, and <italic>K</italic>
(<italic>S</italic>
, <italic>S</italic>
) is equal to 1. We will refer to Equation (1) as the gkm-kernel. It is similar to the wildcard kernel introduced in Ref. <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
, but our approach differs in that we do not sum over the number of wild-cards, or gaps, as formulated in Ref. <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
.</p><p>Since the number of all possible gapped <italic>k</italic>
-mers grows extremely rapidly as <italic>k</italic>
increases, direct calculation of Equation (1) quickly becomes intractable. To implement gapped <italic>k</italic>
-mers as features, it is necessary to overcome this serious issue, by deriving a new equation for <italic>K</italic>
(<italic>S</italic>
<sub>1</sub>
, <italic>S</italic>
<sub>2</sub>
) that does not involve the computation of all possible gapped <italic>k</italic>
-mer counts. The key idea is that only the full <italic>l</italic>
-mers present in the two sequences can contribute to the similarity score via all gapped <italic>k</italic>
-mers derived from them. Thus the inner product in Equation (1), which involves a sum over all gapped <italic>k</italic>
-mers, can be computed by a much more compact sum, which involves only a double sum over the sequential <italic>l</italic>
-mers present in each of the two sequences:<disp-formula><graphic xlink:href="pcbi.1003711.e007.jpg" position="float" orientation="portrait"></graphic>
<label>(2)</label>
</disp-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e008.jpg"></inline-graphic>
</inline-formula><italic>i</italic>
'th <italic>l</italic>
-mer appearing in <italic>S</italic>
<sub>1</sub>
, <inline-formula><inline-graphic xlink:href="pcbi.1003711.e009.jpg"></inline-graphic>
</inline-formula><italic>j</italic>
'th <italic>l</italic>
-mer appearing in <italic>S</italic>
<sub>2</sub>
, and <italic>n</italic>
<sub>1</sub>
and <italic>n</italic>
<sub>2</sub>
are the numbers of full <italic>l</italic>
-mers in <italic>S</italic>
<sub>1</sub>
and <italic>S</italic>
<sub>2</sub>
respectively, i.e. <italic>n</italic>
<sub>1</sub>
= length(<italic>S</italic>
<sub>1</sub>
)−<italic>l</italic>
+1 and <italic>n</italic>
<sub>2</sub>
= length(<italic>S</italic>
<sub>2</sub>
)−<italic>l</italic>
+1. Evaluation of Equation (2) is much more efficient than Equation (1) because almost always, <inline-formula><inline-graphic xlink:href="pcbi.1003711.e010.jpg"></inline-graphic>
</inline-formula><italic>h<sub>lk</sub>
</italic><italic>u</italic>
<sub>1</sub>
, <italic>u</italic>
<sub>2</sub>
) only depends on the number of mismatches, <italic>m</italic>
, between the two full <italic>l</italic>
-mers, <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
<sub>2</sub>
, i.e. <italic>h<sub>lk</sub>
</italic><italic>u</italic>
<sub>1</sub>
, <italic>u</italic>
<sub>2</sub>
) = <italic>h<sub>lk</sub>
</italic><italic>m</italic>
). Therefore, we can rewrite Equation (2) by grouping all the <italic>l</italic>
-mer pairs of the same number of mismatches together as follows:<disp-formula><graphic xlink:href="pcbi.1003711.e011.jpg" position="float" orientation="portrait"></graphic>
<label>(3)</label>
</disp-formula><italic>N<sub>m</sub>
</italic><italic>S</italic>
<sub>1</sub>
, <italic>S</italic>
<sub>2</sub>
) is the number of pairs of <italic>l</italic>
-mers with <italic>m</italic>
mismatches, and <italic>h<sub>lk</sub>
</italic><italic>m</italic>
) is the corresponding coefficient. We refer to <italic>N<sub>m</sub>
</italic><italic>S</italic>
<sub>1</sub>
, <italic>S</italic>
<sub>2</sub>
) as the <italic>mismatch profile of S</italic>
<sub>1</sub>
<italic>and S</italic>
<sub>2</sub>
. Since each <italic>l</italic>
-mer pair with <italic>m</italic>
mismatches contributes to <inline-formula><inline-graphic xlink:href="pcbi.1003711.e012.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
-mers, the coefficient <italic>h<sub>lk</sub>
</italic><italic>m</italic>
), denoted in short by <italic>h<sub>m</sub>
</italic><disp-formula><graphic xlink:href="pcbi.1003711.e013.jpg" position="float" orientation="portrait"></graphic>
<label>(4)</label>
</disp-formula><p>Determining a mismatch profile in Equation (3) is still computationally challenging since the numbers of mismatches between all possible <italic>l</italic>
-mer pairs has yet to be determined. To address this issue, we have developed two different algorithms. We first considered direct evaluation of the mismatch profiles between all pairs of training sequences. To minimize the cost of counting mismatches between two words, we develop an efficient mismatch counting algorithm that practically runs in constant time, independent of <italic>k</italic>
and <italic>l</italic>
parameters (see <xref ref-type="sec" rid="s5">Methods</xref>
). We then use Equation (3) to obtain the inner products for every pair of sequences.</p><p>The direct and sequential evaluation of the kernel function between all training sequences becomes less practical as the number of training sequences gets larger, since it requires O(<italic>N</italic>
<sup>2</sup>
<italic>L</italic>
<sup>2</sup>
) operations of mismatch counting between <italic>l</italic>
-mer pairs, where <italic>N</italic>
is the number of training sequences and <italic>L</italic>
is the average sequence length. Because of this unfavorable scaling, we implemented an alternative method using a <italic>k</italic>
-mer tree data structure, similar to one previously introduced in Ref. <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
, but with some modifications (see <xref ref-type="sec" rid="s5">Methods</xref>
). This method simultaneously calculates the mismatch profile for all the sequence pairs, and, therefore, can significantly reduce the computation time especially when the number of gaps is relatively small, typically when <italic>l – k</italic>
< = 4. We can further improve the efficiency if we truncate the sum in Equation (3) to only consider up to a maximum number of mismatches, <italic>m<sub>max</sub>
</italic><xref ref-type="sec" rid="s5">Methods</xref>
). This approximate method is especially favorable when the number of gaps is large, but the efficiency comes at the cost of exact evaluation of the kernel and classification accuracy, which we will discuss in greater detail below. Therefore, we used one of the two algorithms depending on the size of data sets and the number of gaps we choose for analysis.</p><sec id="s2b"><title>Gapped <italic>k</italic>
-mer SVM classifier outperforms <italic>k</italic>
-mer SVM classifier</title><p>Because of the difficulty of reliably estimating long <italic>k</italic>
-mer counts, we hypothesized that gkm-SVM would perform better than kmer-SVM, and that gapped <italic>k</italic>
-mers would be most advantageous as features, when long TFBSs are important sequence elements in a given data set. To directly test this idea, we compared the classification performance of gkm-SVM to kmer-SVM in predicting the binding sites of CTCF <xref rid="pcbi.1003711-McDaniell1" ref-type="bibr">[13]</xref>
in the human genome, a TF whose binding specificity has been well-characterized <xref rid="pcbi.1003711-Kim1" ref-type="bibr">[14]</xref>
. As shown in <xref ref-type="supplementary-material" rid="pcbi.1003711.s001">Figure S1</xref>
, CTCF recognizes very long DNA sequences (the full PWM is 19 bp), and the genomic CTCF bound regions are almost perfectly predicted by matches to the CTCF PWM (<xref ref-type="supplementary-material" rid="pcbi.1003711.s002">Figure S2</xref>
): in the PWM analysis, we used as a predictor the best matching log-odd score to the PWM model in the region, and achieved area under the ROC curve (AUC) of 0.983. It is very rare for a single PWM to perform this well, and in our experience CTCF is unique in this regard. The CTCF dataset therefore provides an excellent opportunity to test our gapped <italic>k-</italic>
mer classifier. We used the top 2,500 CTCF ChIP-seq signal enriched regions in the GM12878 cell line available at Gene Expression Omnibus (GSE19622) <xref rid="pcbi.1003711-McDaniell1" ref-type="bibr">[13]</xref>
as a positive dataset, and equal numbers of random genomic sequences (1×) as a negative dataset. We generated these negative sequences by matching length, GC and repeat fraction of the positive set <xref rid="pcbi.1003711-FletezBrant1" ref-type="bibr">[6]</xref>
.</p><p>We compared the performance of gkm-SVM and kmer-SVM on the CTCF data set for a range of oligomer lengths by varying either <italic>k</italic>
(for kmer-SVM) or <italic>l</italic>
(for gkm-SVM) from 6 to 20. We fixed the parameter <italic>k = </italic>
6 for gkm-SVM. We then quantified the classification performance of each by calculating test-set AUC with standard five-fold cross validation (CV) (see <xref ref-type="sec" rid="s5">Methods</xref>
). <xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1A</xref>
shows a summary of the comparisons. As anticipated, gkm-SVM performs consistently better than kmer-SVM for all lengths. More significantly, while kmer-SVM suffers severely from overfitting when <italic>k</italic>
is greater than 10, gkm-SVM is virtually unaffected by <italic>l</italic>
. In fact, gkm-SVM achieves the best result (AUC = 0.967) when <italic>l</italic>
= 14 and <italic>k</italic>
= 6, which is significantly better than the kmer-SVM (AUC = 0.912 when <italic>k</italic>
= 10); the best ROC curve is shown in <xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1C</xref>
. It should be noted, however, that the PWM classification result (<xref ref-type="supplementary-material" rid="pcbi.1003711.s002">Figure S2</xref>
) is still the best (AUC = 0.983) among the three methods we tested in this analysis. A complicating factor is that while both kmer-SVM and gkm-SVM use entire sequences (average length is 316 bp) to calculate the prediction scores, the PWM scores are from the best matching 19 bp sub-sequence in the region. It may be that the extra ∼300 bp sequences contribute noise in the SVM prediction scores, which slightly impairs the overall classification accuracy. In any event, the gkm-SVM is a significant improvement in accuracy over the kmer-SVM, and both gkm-SVM and the PWM are excellent predictors on this dataset.</p><fig id="pcbi-1003711-g001" orientation="portrait" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1003711.g001</object-id>
<label>Figure 1</label>
<caption><title>gkm-SVM outperforms kmer-SVM over a wide range of <italic>k</italic>
-mer length.</title><p>Both gkm-SVM and kmer-SVM were trained on (A) CTCF bound and (B) EP300 bound genomic regions using different word lengths (<italic>k</italic>
for kmer-SVM and <italic>l</italic>
for gkm-SVM). The parameter <italic>k</italic>
for gkm-SVM was fixed at 6. While AUCs of the kmer-SVMs show significant overfitting in both cases as <italic>k</italic>
gets larger (dotted), gkm-SVMs accuracy is higher for a broad range of larger <italic>l</italic>
(solid). Results using the truncated Gkm-SVM with <italic>m</italic>
<sub>max</sub>
= 3 are shown as dashed lines and AUCs of these faster approximations are comparable when the difference between m<sub>max</sub>
and <italic>l</italic>
– <italic>k</italic>
are relatively small. ROC for the optimal <italic>k</italic>
or <italic>l</italic>
for each case are shown in (C) and (D). Gkm-SVMs (solid) consistently outperform kmer-SVMs (dashed) on both data sets. Error bars here and below represent 5-fold CV standard deviation.</p><graphic xlink:href="pcbi.1003711.g001"></graphic>
</fig><p>Interestingly, gkm-SVM shows consistently better performance than kmer-SVM even if <italic>l</italic>
is relatively small (<italic>l</italic>
<10) (<xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1A</xref>
). This suggests that gkm-SVM may also be better at modeling diverse combinations of TFBSs than kmer-SVM. To test this hypothesis, we analyzed a mouse enhancer dataset of more varied sequence composition: genomic EP300 bound regions in embryonic mouse forebrain <xref rid="pcbi.1003711-Visel1" ref-type="bibr">[15]</xref>
. We have previously shown that our original kmer-SVM classifiers can accurately predict EP300 binding when mediated by sets of active TFBSs <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
. This EP300 data set thus provides a direct test of the effectiveness of using gapped <italic>k</italic>
-mer features to detect more complex regulatory features. For this analysis, we defined a new set of the 1,693 400 bp sites that maximize the EP300 ChIP-seq signal within each of the peaks determined by MACS <xref rid="pcbi.1003711-Zhang1" ref-type="bibr">[16]</xref>
after removing any regions which were more than 70% repeats. We repeated the <italic>k</italic>
and <italic>l</italic>
scaling with the EP300 data set and a 1× negative set, and again found that gkm-SVM consistently outperforms kmer-SVM for all feature lengths (<xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1B</xref>
). Analogous to the observations modeling CTCF binding, gkm-SVM AUC is high and does not degrade with large <italic>l</italic>
. In contrast, the kmer-SVM accuracy drops rapidly as <italic>k</italic>
increases. Moreover, although the difference in performance is smaller than found for the CTCF data set, the gkm-SVM achieves the best AUC (0.947) with <italic>l = </italic>
9 and <italic>k = </italic>
6, while the kmer-SVM achieves 0.932 with <italic>k = </italic>
7, suggesting that longer <italic>k</italic>
-mers with some flexibility do contain more complete information about TF binding (<xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1D</xref>
). At the same time, the gapped <italic>k-</italic>
mer features are more robustly estimated (having more counts) and for this reason make more reliable predictors. The consequences of these improvements in AUC are significant when considering the genome-scale precision of the improved gkm-SVM classifiers. The rate of false positive predictions is dominated by the large neutral fraction of the genome, so the precision of a genome-scale classifier is best assessed by a Precision-Recall curve in combination with a much larger negative set, as discussed in Ref. <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
. The Precision-Recall curves for the gkm-SVM and kmer-SVM classifiers on a 100× negative set are shown in <xref ref-type="supplementary-material" rid="pcbi.1003711.s003">Figure S3</xref>
. For CTCF, at a recall of 50%, the precision increases from 36% to 59%. These ranges of precision and recall are in the relevant range of experiments aiming to discover and test novel enhancers, and we therefore expect that predictions based on gkm-SVM will have up to a two-fold higher successful validation rate.</p><p>One further modification can substantially reduce the computational cost of using gapped <italic>k</italic>
-mers with little degradation in performance. The algorithm using the <italic>k</italic>
-mer tree data structure produces identical results to the direct evaluation of Equation (3), but typically is much faster when the number of mismatches, <italic>l – k</italic>
, is smaller than four, and the number of training sequences is large. The <italic>k</italic>
-mer tree algorithm can be made even more computationally efficient, if we prune the traversal of the tree, by ignoring any <italic>k</italic>
-mer pairs that have more mismatches than a predetermined parameter, <italic>m<sub>max</sub>
</italic><italic>h<sub>m</sub>
</italic><italic>m</italic>
. This approximation significantly reduces the total number of calculations and allows the user to control the running time of the algorithm by setting the parameter <italic>m<sub>max</sub>
</italic><italic>l</italic>
feasible for any given <italic>k</italic>
. To systematically investigate the classification performance of this approximation, we applied the same analysis above using both CTCF and EP300 data sets (<xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1A, B</xref>
), and found that AUCs from the approximate method are virtually identical to the exact method when the difference between <italic>m<sub>max</sub>
</italic><italic>l</italic>
– <italic>k</italic>
are small. Interestingly, the approximation method achieved even higher AUC with CTCF data set in some cases.</p><p>Encouraged by the analyses of CTCF and EP300 data sets above, we systematically compared gkm-SVM to kmer-SVM using a very broad range of human data sets generated by the ENCODE project <xref rid="pcbi.1003711-Gerstein1" ref-type="bibr">[17]</xref>
,<xref rid="pcbi.1003711-Wang1" ref-type="bibr">[18]</xref>
. We used 467 sets of ChIP-seq peaks produced by the ENCODE uniform processing pipeline containing at least 500 regions (see <xref ref-type="sec" rid="s5">Methods</xref>
). We truncated any data set with greater than 5,000 regions by random sampling. We then trained both kmer-SVM and gkm-SVM on each set against an equal size (1×) negative set of random genomic regions and calculated AUCs with five-fold cross validation. We used <italic>k</italic>
= 6 for kmer-SVM, and <italic>l</italic>
= 10 and <italic>k</italic>
= 6 for gkm-SVM, but as shown in <xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1</xref>
the improvements are generally insensitive to these parameter choices. Strikingly, we find that gkm-SVM almost always outperforms kmer-SVM (<xref ref-type="fig" rid="pcbi-1003711-g002">Figure 2A</xref>
). We also find that variances of AUCs from test CV sets are generally reduced, suggesting that gkm-SVM is more robust than kmer-SVM (<xref ref-type="supplementary-material" rid="pcbi.1003711.s004">Figure S4</xref>
). More significantly, gkm-SVM performs much better especially for TFs with long binding sites. In this dataset, most of these long binding sites arise in ChIP-seq data sets for CTCF and members of the cohesin complex (RAD21, SMC3) known to be physically associated with CTCF <xref rid="pcbi.1003711-Parelho1" ref-type="bibr">[19]</xref>
. On these CTCF associated factors gkm-SVM exhibits much higher AUC than kmer-SVM, as highlighted by the cluster of purple circles in <xref ref-type="fig" rid="pcbi-1003711-g002">Figure 2A</xref>
. We have also compared gkm-SVM to the best single PWM AUC as shown in Ref. <xref rid="pcbi.1003711-FletezBrant1" ref-type="bibr">[6]</xref>
(<xref ref-type="fig" rid="pcbi-1003711-g002">Figure 2B</xref>
). As expected, gkm-SVM outperforms all datasets except CTCF, for which gkm-SVM performance is only marginally reduced. For a consistent analysis of this dataset, we used <italic>l</italic>
= 10 and <italic>k</italic>
= 6, although for CTCF the gkm-SVM performance is optimal at larger <italic>l</italic>
, as seen in <xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1A</xref>
.</p><fig id="pcbi-1003711-g002" orientation="portrait" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1003711.g002</object-id>
<label>Figure 2</label>
<caption><title>gkm-SVM consistently outperforms kmer-SVM and the best known PWM on human ENCODE ChIP-seq data sets.</title>
<p>(A) We trained gkm-SVM and kmer-SVM on the complete set of 467 ENCODE ChIP-seq data sets (with <italic>k</italic>
= 6 for kmer-SVM, and <italic>l</italic>
= 10 and <italic>k</italic>
= 6 for gkm-SVM). gkm-SVM AUC is consistently higher than kmer-SVM with only a few very minor exceptions. The gkm-SVM method specially outperforms the kmer-SVM for the data sets bound by members of the CTCF complex, highlighted as purple circles. (B) We also compared gkm-SVM and the best known PWM on the same data sets, and gkm-SVM AUCs are significantly higher than the PWM AUC in almost all cases. (C) The ENCODE data sets were divided into four groups: (1) no PWM, (2) only one PWM, (3) two PWMs, and (4) three or more PWMs identified by Wang <italic>et al</italic>
. Then, for each group except the first one, we calculated the number of PWMs recovered by our method. At least one PWM was recovered for more than ∼90% of the data sets.</p><graphic xlink:href="pcbi.1003711.g002"></graphic>
</fig><sec id="s2c"><title>Motif analysis of the ENCODE ChIP-seq data sets</title>
<p>The predictive sequence features that allow gkm-SVM to outperform the single best PWM imply that cooperative binding is the underlying molecular mechanisms that targets TFs to these regulatory regions. Previously we have typically focused on a handful of the highest SVM weight <italic>k</italic>
-mers (say top ten positive and top ten negative weight <italic>k</italic>
-mers) to interpret the classification results <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
,<xref rid="pcbi.1003711-FletezBrant1" ref-type="bibr">[6]</xref>
,<xref rid="pcbi.1003711-Gorkin1" ref-type="bibr">[20]</xref>
. This simple method becomes unwieldy when applied to the gkm-SVM results because of the large number of very similar significant features (when <italic>l</italic>
and/or <italic>k</italic>
are large). Although the <italic>k</italic>
-mers at the extreme top and bottom tails of the <italic>k</italic>
-mer weight distribution are still important and biologically meaningful, those <italic>k</italic>
-mers usually cover only a fraction of the significant feature set, and many more important features are included in the larger tails of the <italic>k</italic>
-mer weight distribution. Therefore, more sophisticated algorithms are needed to extract the biologically relevant features from the classification results.</p><p>To directly address this issue, we developed a new method to combine multiple similar <italic>k</italic>
-mers into more compact and interpretable PWMs and analyzed the 467 ENCODE data sets <xref rid="pcbi.1003711-Wang1" ref-type="bibr">[18]</xref>
. In this approach, we used a larger number of predictive <italic>k</italic>
-mers to build <italic>de novo</italic>
PWMs (see <xref ref-type="sec" rid="s5">Methods</xref>
). We used the top 1% of 10-mers from each of the gkm-SVMs trained on the ENCODE data sets and identified up to three distinct PWMs (<xref ref-type="supplementary-material" rid="pcbi.1003711.s005">Figure S5</xref>
) from <italic>k</italic>
-mers in this set. We then compared our results with the previous PWMs found in the same data sets using a conventional tool (MEME-ChIP) <xref rid="pcbi.1003711-Wang1" ref-type="bibr">[18]</xref>
<xref rid="pcbi.1003711-Machanick1" ref-type="bibr">[21]</xref>
. Similar to our approach, Wang <italic>et al</italic>
. analyzed 457 ENCODE ChIP-seq data sets (440 sets are in common with those we analyzed above) and identified five PWMs from each data set. Collectively, Wang <italic>et al.</italic>
found 79 distinct PWMs enriched, of which our method recovered 74. Comparing each ChIP-seq data set individually, we recovered most of the PWMs reported by <italic>Wang et al.</italic>
using our method (<xref ref-type="fig" rid="pcbi-1003711-g002">Figure 2C</xref>
). Interestingly, while Wang <italic>et al.</italic>
largely failed to identify biologically meaningful PWMs from most of the POL2 ChIP-seq data sets (47 out of 58 sets returned no meaningful PWMs), our methods frequently identified cell-specific TFs as well as promoter specific TFs (<xref ref-type="supplementary-material" rid="pcbi.1003711.s005">Figure S5</xref>
). For example, the GATA1 TF identified from POL2 ChIP-seq in the erythroleukemic cell line K562 is known to play central roles in erythroid differentiation <xref rid="pcbi.1003711-Pevny1" ref-type="bibr">[22]</xref>
. The ETS1 TF from HUVEC is another extensively studied TF, known to be important for angiogenesis <xref rid="pcbi.1003711-Lelivre1" ref-type="bibr">[23]</xref>
. A major difference between the two methods is the number of training sequences. While the previous study was limited to the top 500 of ChIP-seq peaks (ranked by ChIP-seq signal), we were able to use 10× larger numbers of ChIP-seq peaks (5,000 regions), and the large training sizes enabled us to robustly identify diverse combinatorial sequence features.</p><sec id="s2d"><title>Comparison to previous kernels</title>
<p>Since the early development of <italic>k</italic>
-mer based supervised machine learning techniques <xref rid="pcbi.1003711-Leslie3" ref-type="bibr">[24]</xref>
, there have been a number of improvements. Some of these extend the feature set to include imperfect matches, similar in spirit to our gkm-SVM. The mismatch string kernel <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
is one such method, originally motivated by the fact that homologous protein sequences are not usually identical and have many frequently mutated positions. The mismatch kernel also uses <italic>k</italic>
-mers as features, but allows some mismatches when counting <italic>k</italic>
-mers and building feature vectors. The wildcard kernel <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
is another variant of the original string kernel, which introduces a wildcard character that matches any single letter in the given alphabet. More recently, an alternative di-mismatch kernel <xref rid="pcbi.1003711-Agius1" ref-type="bibr">[10]</xref>
has been proposed to directly model TFBSs, and has been successfully applied to protein binding microarray (PBM) data sets <xref rid="pcbi.1003711-Berger1" ref-type="bibr">[25]</xref>
and several other ChIP-seq data sets <xref rid="pcbi.1003711-Agius1" ref-type="bibr">[10]</xref>
,<xref rid="pcbi.1003711-Arvey1" ref-type="bibr">[11]</xref>
. The di-mismatch method tries to overcome the limitation of the mismatch kernel by favoring <italic>k</italic>
-mers with consecutive mismatches. However, in a recent comparison of methods for modeling transcription factor sequence specificity, full <italic>k-</italic>
mer methods outperformed the di-nucleotide approaches when applied to PBM data <xref rid="pcbi.1003711-Weirauch1" ref-type="bibr">[26]</xref>
.</p><p>To further evaluate our proposed method, we directly compared the gkm-kernel with the aforementioned three alternative methods, Mismatch kernel <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
, Wildcard kernel <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
, and Di-mismatch kernel <xref rid="pcbi.1003711-Agius1" ref-type="bibr">[10]</xref>
,<xref rid="pcbi.1003711-Arvey1" ref-type="bibr">[11]</xref>
, using the mouse forebrain EP300 data set. As shown in <xref ref-type="fig" rid="pcbi-1003711-g003">Figure 3</xref>
, gkm-kernel outperforms the other three existing methods both in terms of the classification accuracy and running time. The best AUC we achieved for gkm-kernel is 0.947 as compared to 0.937, 0.935, and 0.944 for the wildcard kernel, mismatch kernel, and di-mismatch kernel, respectively (<xref ref-type="fig" rid="pcbi-1003711-g003">Figure 3A</xref>
). Although the wildcard kernel and gkm-kernel are quite similar, the systematic improvement in gkm-kernel AUCs is primarily due to the incorporation of reverse complement sequences. We directly tested this by adding reverse complement sequences to the feature set for the previously published methods, and indeed found that with this modification, these methods were also able to achieve comparable AUCs (<xref ref-type="supplementary-material" rid="pcbi.1003711.s006">Figure S6</xref>
).</p><p>To further evaluate our proposed method, we directly compared the gkm-kernel with the aforementioned three alternative methods, Mismatch kernel <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
, Wildcard kernel <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
, and Di-mismatch kernel <xref rid="pcbi.1003711-Agius1" ref-type="bibr">[10]</xref>
,<xref rid="pcbi.1003711-Arvey1" ref-type="bibr">[11]</xref>
, using the mouse forebrain EP300 data set. As shown in <xref ref-type="fig" rid="pcbi-1003711-g003">Figure 3</xref>
, gkm-kernel outperforms the other three existing methods both in terms of the classification accuracy and running time. The best AUC we achieved for gkm-kernel is 0.947 as compared to 0.937, 0.935, and 0.944 for the wildcard kernel, mismatch kernel, and di-mismatch kernel, respectively (<xref ref-type="fig" rid="pcbi-1003711-g003">Figure 3A</xref>
). Although the wildcard kernel and gkm-kernel are quite similar, the systematic improvement in gkm-kernel AUCs is primarily due to the incorporation of reverse complement sequences. We directly tested this by adding reverse complement sequences to the feature set for the previously published methods, and indeed found that with this modification, these methods were also able to achieve comparable AUCs (<xref ref-type="supplementary-material" rid="pcbi.1003711.s006">Figure S6</xref>
).</p><fig id="pcbi-1003711-g003" orientation="portrait" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1003711.g003</object-id>
<label>Figure 3</label>
<caption><title>Comparison of gkm-SVM and existing methods on the mouse forebrain EP300 data set.</title>
<p>(A) For each method, averages of 5-CV AUCs are shown as a function of the word length with the optimal number of mismatches, <italic>m</italic>
, held fixed. Also shown are gkm-SVM results using fixed <italic>k</italic>
= 6 and varying <italic>m</italic>
<sub>max</sub>
. (B) Running time for each of the kernel computations shown in (A). Gkm-kernels show better classification performance and significantly more efficient computation at peak AUC.</p><graphic xlink:href="pcbi.1003711.g003"></graphic>
</fig><p>More significantly, when we compare running times at parameters which maximize AUC for each method, our gkm-SVM implementation (<italic>l</italic>
= 9,<italic>l-k</italic>
= 3) is roughly two orders of magnitude faster than di-mismatch (10,3), and slightly more than one order of magnitude faster than mismatch (<italic>l</italic>
= 10,<italic>m</italic>
= 2) and wildcard (<italic>l</italic>
= 8,<italic>m</italic>
= 3) on the EP300 data set (<xref ref-type="fig" rid="pcbi-1003711-g003">Figure 3B</xref>
and <xref ref-type="supplementary-material" rid="pcbi.1003711.s007">Figure S7</xref>
). Also, by fixing <italic>k</italic>
= 6 and the parameter <italic>m<sub>max</sub>
</italic><italic>l</italic>
, compared to a scan at fixed <italic>m</italic>
, varying <italic>k</italic>
(<xref ref-type="fig" rid="pcbi-1003711-g003">Figure 3A</xref>
). Direct running time comparisons using our tree structure in the mismatch and wildcard kernels (described below) are shown in <xref ref-type="supplementary-material" rid="pcbi.1003711.s007">Figure S7B and S7C</xref>
. We should note that we were only able to test the di-mismatch kernel up to <italic>l</italic>
= 10, because it required more than 128 GB of memory and did not finish within 2000 minutes when using <italic>l</italic>
= 11.</p><p>Interestingly, we also note that both Mismatch kernel and Wildcard kernel are special cases of the more general class of kernels, defined by Equation (3). This unification allows direct application of the methods we developed for mismatch profile computation and therefore gives more efficient methods for computation of these existing methods (see <xref ref-type="sec" rid="s5">Methods</xref>
).</p><sec id="s2e"><title>Calculation and performance of estimated <italic>l</italic>
-mer frequencies for gkm-SVM</title><p>As an alternative to the gapped <italic>k</italic>
-mer feature set, we also developed an alternative kernel by replacing the <italic>k</italic>
-mer counts with robust <italic>l</italic>
-mer count estimates <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
in our original kmer-SVM framework. We have developed efficient methods to compute this new kernel (see <xref ref-type="sec" rid="s5">Methods</xref>
). In Ref. <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
, we considered the mapping from <italic>l</italic>
-mers to gapped <italic>k</italic>
-mers. Among all possible sets of <italic>l</italic>
-mer frequencies that could produce the same gapped <italic>k</italic>
-mer frequency distribution, we developed a method to estimate the “most likely” <italic>l</italic>
-mer frequency set. Full details of this method are described in the Ref. <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
. In brief, we first define a gapped <italic>k</italic>
-mer count vector <inline-formula><inline-graphic xlink:href="pcbi.1003711.e014.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
-mer feature vector for gkm-SVM as shown above. Then, the count estimate, <inline-formula><inline-graphic xlink:href="pcbi.1003711.e015.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
-mer <italic>u</italic>
is given by<disp-formula><graphic xlink:href="pcbi.1003711.e016.jpg" position="float" orientation="portrait"></graphic>
<label>(5)</label>
</disp-formula><p>The weight <inline-formula><inline-graphic xlink:href="pcbi.1003711.e017.jpg"></inline-graphic>
</inline-formula><italic>m</italic>
, between the gapped <italic>k</italic>
-mer corresponding to <inline-formula><inline-graphic xlink:href="pcbi.1003711.e018.jpg"></inline-graphic>
</inline-formula><italic>u</italic>
, and takes the following form:<disp-formula><graphic xlink:href="pcbi.1003711.e019.jpg" position="float" orientation="portrait"></graphic>
<label>(6)</label>
</disp-formula><italic>b</italic>
is the alphabet size and is equal to four in case of DNA sequences (A,C,G and T). Since the above equation is applied to every <italic>l</italic>
-mer, it would provide a non-zero frequency even for an <italic>l</italic>
-mer that does not have any exact match appearing in any training set sequence.</p><p>Direct calculation of Equation (5), however, requires actual counting of all of the <italic>M</italic>
gapped <italic>k</italic>
-mers, which becomes computationally intractable for large <italic>l</italic>
and <italic>k</italic>
in a way similar to Equation (1). Besides, summing up a large set of floating point numbers may result in poor numerical precision. To overcome these issues, we developed a simple method, referred to as the <italic>gkm-filter</italic>
, to more efficiently calculate the robust <italic>l</italic>
-mer count estimates, <inline-formula><inline-graphic xlink:href="pcbi.1003711.e020.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
-mer counts (see <xref ref-type="sec" rid="s5">Methods</xref>
). In summary, in the calculation of the robust <italic>l</italic>
-mer count estimates, we give a non-zero weight to <italic>l</italic>
-mers with few numbers of mismatches. The <italic>k</italic>
-mer frequency estimation method is not constrained to produce non-negative frequencies and may occasionally generate negative count estimates. To obtain strictly positive frequencies, we used a revised version of the gkm-filter method, which we call the <italic>truncated gkm-filter</italic>
. Finally, we developed a method to directly calculate the kernels using these feature sets (see <xref ref-type="sec" rid="s5">Methods</xref>
). An important result here is that the evaluation of the gkm-kernel (the inner product of the <italic>l</italic>
-mer count estimates vectors) is still given by Equation (3), but with a new set of weights <italic>c<sub>lk</sub>
</italic><italic>m</italic>
) given by Equation (14), below, replacing <italic>h<sub>lk</sub>
</italic><italic>m</italic>
). Therefore, efficient algorithms for pairwise mismatch profiles that we developed for the gkm-kernel can be directly used for this new feature set without any modification. Because of this symmetry, we also refer to this method as gkm-kernel <italic>with (full or truncated) filter</italic>
. A numerical example using count estimates on two short sequences is provided in <xref ref-type="supplementary-material" rid="pcbi.1003711.s014">Text S1</xref>
.</p><p>To systematically compare the classification performance of these new methods with the original gapped <italic>k</italic>
-mers, we repeated the previous analysis with the ENCODE ChIP-seq data sets. Using the truncated gkm-filter yields results highly comparable to the original gkm-SVM for most datasets with modestly but consistently better relative performance when AUC is greater than 0.9 (shown as purple circles in <xref ref-type="supplementary-material" rid="pcbi.1003711.s008">Figure S8A</xref>
). Any improvement in the range of high AUC (>0.9) typically strongly reduces the classifier's False Prediction Rate <xref rid="pcbi.1003711-Lee2" ref-type="bibr">[27]</xref>
, therefore, we generally recommend the truncated filter method as the method of choice for most analyses. Compared to the original gkm-SVM, using the gkm-SVM with full filter yields lower AUCs (<xref ref-type="supplementary-material" rid="pcbi.1003711.s008">Figure S8B</xref>
) although it is still significantly higher compared to the kmer-SVM method.</p><sec id="s2f"><title>Application of the robust <italic>l</italic>
-mer count estimates for Naïve-Bayes classifier</title><p>So far, we have focused on using gapped <italic>k</italic>
-mer based methods for improving sequence kernel methods. We have shown that, by direct use of gapped <italic>k</italic>
-mers as features or by using the robust <italic>l</italic>
-mer count estimates, we can significantly overcome the long <italic>k</italic>
-mers' sparse count problem for these methods. We further demonstrate the general utility of the robust <italic>l</italic>
-mer count estimates in sequence classification problems by applying it to a simple Naïve-Bayes (NB) classifier similar to the one previously introduced in Ref. <xref rid="pcbi.1003711-Sandberg1" ref-type="bibr">[28]</xref>
and show that by using robust count estimates instead of conventional <italic>k</italic>
-mer counts we can significantly boost the performance of the Naïve-Bayes classifier for long <italic>k</italic>
-mers.</p><p>Here, we used the log-likelihood ratio of the estimated <italic>l</italic>
-mer frequencies in the positive and negative sets as a predictor, using the NB assumption of feature independence. The prediction score of any given sequence of length <italic>n</italic>
, denoted by <bold><italic>S</italic>
</bold><italic>s</italic>
<sub>0</sub>
<italic>s</italic>
<sub>1</sub>
…<italic>s<sub>n</sub>
</italic><sub>–1</sub>
, is then given by:<disp-formula><graphic xlink:href="pcbi.1003711.e021.jpg" position="float" orientation="portrait"></graphic>
<label>(7)</label>
</disp-formula><italic>N<sub>P</sub>
</italic><italic>N<sub>N</sub>
</italic><italic>l</italic>
-mers, <italic>s<sub>i</sub>
s<sub>i</sub>
</italic><sub>+1</sub>
…<italic>s<sub>i+l</sub>
</italic><sub>−1</sub>
, in the positive and negative training set, and are given by Equation (11) below. We used the <italic>truncated</italic>
gkm-filter method adding pseudo-count (half of the smallest positive coefficient of the truncated gkm-filter) to each of the estimated frequencies to obtain strictly positive frequencies for log-likelihood ratio. As a comparison, we also implemented the NB classifier without the gkm-filter, using actual <italic>l</italic>
-mer counts with a pseudo-count (0.5) for <italic>N<sub>P</sub>
</italic><italic>N<sub>N</sub>
</italic><xref ref-type="supplementary-material" rid="pcbi.1003711.s002">Figure S2</xref>
), and the local sequence features around the CTCF binding motif do not seem to significantly contribute to the prediction. Thus, to precisely detect the CTCF binding motif and achieve the best classification performance, we scored every substring of length <italic>n</italic>
= 15+<italic>l</italic>
−1 for each sequence and assigned the maximum as the final score for the sequence. The window size of 15 was chosen to optimize the detection of the CTCF site within a small window of flanking sequence, which maximizes the performance of the NB classifier without the gkm-filter. For the EP300 genomic bound regions, in contrast, we used the full sequence in both classifiers. We compare the performance of these NB classifiers on both data sets in <xref ref-type="fig" rid="pcbi-1003711-g004">Figure 4</xref>
for a range of feature length (6–20 bp). Similar to the previous analysis using gkm-SVM and kmer-SVM (<xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1</xref>
), using robust count estimates (gkm-filter) significantly improves the classification accuracy especially for longer k-mers (<xref ref-type="fig" rid="pcbi-1003711-g004">Figure 4</xref>
). On the CTCF data set, the NB classifier using the gkm-filter achieves best performance with <italic>l</italic>
= 20 (AUC = 0.99), which is even better than that of the CTCF PWM (red dotted line, AUC = 0.983) (<xref ref-type="fig" rid="pcbi-1003711-g004">Figure 4A</xref>
). Also on the EP300 dataset, the gkm-filter significantly improves the overall performance of NB classifier (<xref ref-type="fig" rid="pcbi-1003711-g004">Figure 4B</xref>
). The superior classification performance using gapped <italic>k</italic>
-mer based features is thus consistent for both SVM and NB classifiers, and strongly suggests that the robust <italic>l</italic>
-mer count estimates provide a more complete and robust set of sequence features than simple <italic>k</italic>
-mers in most sequence classification problems, as hinted at in our preliminary analysis of <italic>k-</italic>
mer frequency spectra in Ref. <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
.</p><fig id="pcbi-1003711-g004" orientation="portrait" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1003711.g004</object-id>
<label>Figure 4</label>
<caption><title>Gapped <italic>k</italic>
-mer features also improve performance of Naïve-Bayes classifiers.</title><p>Naïve-Bayes classifiers were trained on (A) CTCF bound and (B) EP300 bound genomic regions using different word lengths, <italic>k</italic>
, using both actual <italic>k</italic>
-mer counts (dashed), and estimated <italic>k</italic>
-mer counts from the gkm-filter (solid). As shown above for SVM, the Naïve-Bayes accuracy as measured by AUC is systematically higher using gapped <italic>k</italic>
-mer estimated frequencies instead of actual <italic>k</italic>
-mer counts, further supporting the utility of gapped <italic>k</italic>
-mer based features. For CTCF the Naïve-Bayes AUC is comparable to the best SVM (dotted red lines), but for EP300 the SVM outperforms the Naïve-Bayes classifier.</p><graphic xlink:href="pcbi.1003711.g004"></graphic>
</fig><sec id="s3"><title>Discussion</title>
<p>In this paper, we presented a significantly improved method for sequence prediction using gapped <italic>k</italic>
-mers as features, gkm-SVM. We introduced a new set of algorithms to efficiently calculate the kernel matrix, and demonstrated that by using these new methods we can significantly overcome the sparse <italic>k</italic>
-mer count problem for long <italic>k</italic>
-mers and hence significantly improve the classification accuracy especially for long TFBSs. Detailed comparisons of our proposed method with some existing methods show that our gkm-SVM outperforms existing methods in terms of classification accuracy on benchmark data and is also typically orders of magnitude faster. We also introduced the concept of gkm-filters for efficient calculation of the robust <italic>k</italic>
-mer count estimates and derived optimal weights for penalizing different number of mismatches. We showed that one could successfully replace <italic>k</italic>
-mers with robust <italic>k</italic>
-mer count estimates to avoid long <italic>k</italic>
-mer sparse count problem, and demonstrated the effectiveness of this method by showing examples in SVM and Naïve-Bayes classifiers. We thus expect that most <italic>k</italic>
-mer based methods can be significantly improved by simply using this generalized <italic>k</italic>
-mer count.</p><p>The main biological relevance of the computational method we present in this paper is that gkm-SVM is capable of accurately predicting a wide range of specific classes of functional regulatory elements based on DNA sequence features in those elements alone. This in itself is interesting and implies that the epigenomic state of a DNA regulatory element primarily is specified by its sequence. In addition, our predictions facilitate direct investigation of how these elements function, either by targeted mutation of the predictive elements within the larger regulatory region, or by modulating the activity of the TFs which bind the predictive sequence elements. We are currently using changes in the gkm-SVM score to systematically evaluate the predicted impact of human regulatory variation (single nucleotide polymorphisms (SNPs) or indels) to interpret significant SNPs identified in genome wide association studies. We demonstrated that gkm-SVM is better at predicting all ENCODE ChIP-seq data than the best single PWM found from the ChIP-seq regions, or previously known PWMs. The gkm-SVM is able to do so by integrating cofactor sequences which may not be directly bound by the ChIP-ed TF but facilitate its occupancy. To predict this ChIP-seq set accurately required the improved accuracy of the gkm-SVM and its ability to describe longer binding sites such as CTCF, which were very difficult for our earlier kmer-SVM approach. We recovered most of the cofactors found by traditional PWM discovery methods, but we further show that these combinations of cofactors are predictive in the sense that they are sufficient to define the experimentally bound regions.</p>
<p>There are some further issues that need to be considered in the application of these methods. First, one will typically be interested in finding an optimal set of the parameters (<italic>l</italic>
and <italic>k</italic>
) to achieve the best classification performance. A significant advantage of gapped <italic>k</italic>
-mer methods over <italic>k</italic>
-mer methods is that they are more robust and are less sensitive to the particular choices of <italic>l</italic>
or <italic>k</italic>
compared to kmer-SVM or NB classifiers, as shown in <xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1</xref>
and <xref ref-type="fig" rid="pcbi-1003711-g004">Figure 4</xref>
. Nevertheless, these parameters can still be optimized to maximize cross validation AUC. As a general rule, we have found that when choosing the parameter <italic>k</italic>
, which determines how different numbers of mismatches are weighted, given a whole word length <italic>l</italic>
, smaller values of <italic>k</italic>
(typically less than 8) are usually better when important sequence elements are believed to be more degenerate or when only small amount of training data is available. Although the choice of <italic>k</italic>
directly affects the feature set, our analysis of several datasets shows that the overall performance of the classifier is not very sensitive to changes in <italic>k</italic>
. The parameter <italic>l</italic>
is directly related to TFBS lengths and should be comparable to or slightly larger than the longest important feature, as demonstrated by our analysis of the CTCF and EP300 data sets in <xref ref-type="fig" rid="pcbi-1003711-g001">Figure 1</xref>
and <xref ref-type="fig" rid="pcbi-1003711-g004">Figure 4</xref>
.</p><p>Our approach also avoids an issue that would arise if one chose instead to directly use Equation (5) for computing count estimates. This would involve a large number of floating point operations, and accumulated round-off error could become significant in the large summations. There are some algorithms, such as Kahan compensated summation <xref rid="pcbi.1003711-Higham1" ref-type="bibr">[29]</xref>
, which can significantly reduce this error, however, we explicitly avoided evaluating this sum by first computing the mismatch profiles between sequences, which involves only integer calculations. Then, we calculate the weighted sum of the number of mismatches using Equation (11), which involves a much smaller number of floating point operations.</p><p>Two issues which are left for future investigation are different treatment of end vs. internal gaps, and allowing imperfect mismatches. We currently do not make special consideration for gaps which occur at the end of a <italic>k</italic>
-mer instead of internal gaps. Also, our implementation of a mismatch treats all nucleotides equally, but often TF binding sites can prefer an A or T in a given position, or a purine vs. pyrimidine pair. Our approach recovers these preferences by assigning different weights to <italic>k</italic>
-mers which do not have gaps at these positions, but including a wider alphabet including (W,S,Y,R) for (AT,GC,AG,CT) may have some advantages.</p><p>Throughout this paper, we have focused on using DNA sequences as features for classifying the molecular or biological function of a genomic region. However, in principle, our method can be applied to any classification or prediction problem involving a large feature set. In general, when the number of features used by a classifier increases, the number of samples in the training set for each point in the feature space becomes smaller, and small sample count issues occur (which we have resolved using gapped <italic>k</italic>
-mers). One approach to the large feature space is <italic>feature selection</italic>
, which selects a subset of features and builds a classifier only using those features, ignoring all the other features. However, usually a limited subset of features cannot explain all the variation in the predicted quantity. While hypothetical at this point, our analysis suggests that an alternative approach might be of general value. Analogous to the way we have used gapped <italic>k</italic>
-mers to more robustly estimate <italic>k</italic>
-mer feature frequencies, we speculate that there may be a general approach which uses subsets of a larger feature set to combine observed feature counts with weights reflecting the similarity to some generalized feature. These estimated feature frequencies will be less susceptible to statistical noise by construction, and thus may provide consistently better classification performance, as we have shown for gapped <italic>k</italic>
-mers.</p><sec sec-type="methods" id="s4"><title>Methods</title>
<sec id="s4a"><title>Support Vector Machine</title>
<p>The Support Vector Machine (SVM) <xref rid="pcbi.1003711-Boser1" ref-type="bibr">[301]</xref>
,<xref rid="pcbi.1003711-Vapnik1" ref-type="bibr">[31]</xref>
is one of the most successful binary classifiers and has been widely used in many classification problems. We have previously developed an SVM based framework, or “kmer-SVM”, for enhancer prediction and have successfully applied to embryonic mouse enhancers <xref rid="pcbi.1003711-Lee1" ref-type="bibr">[3]</xref>
and many other regulatory datasets <xref rid="pcbi.1003711-FletezBrant1" ref-type="bibr">[6]</xref>
,<xref rid="pcbi.1003711-Gorkin1" ref-type="bibr">[20]</xref>
. Briefly, our kmer-SVM method finds a decision boundary that maximally discriminates a set of regulatory sequences from random genomic non-regulatory sequences in the <italic>k</italic>
-mer frequency feature vector space. Here, we developed new kernel functions using gapped <italic>k</italic>
-mers and <italic>l</italic>
-mer count estimates as features, and software that calculates the kernel matrix. For SVM training, we developed a custom Python script that takes the kernel matrix as input and learns support vectors. We used Shogun Machine Learning Toolbox <xref rid="pcbi.1003711-Sonnenburg1" ref-type="bibr">[32]</xref>
and SVM-light <xref rid="pcbi.1003711-Joachims1" ref-type="bibr">[33]</xref>
for the SVM training script. As an alternative method, we also implemented an SVM classifier based on the iterative algorithm described in Ref. <xref rid="pcbi.1003711-Jaakkola1" ref-type="bibr">[34]</xref>
.</p><sec id="s4b"><title>Direct computation of Gkm-kernel</title>
<p>For direct computation of the gkm-SVM kernel matrix, we represent each training sequence with a list of <italic>l</italic>
-mers and corresponding count for each <italic>l</italic>
-mer. Then for each pair of sequences, we compute the number of mismatches for all pairs of <italic>l</italic>
-mers and use the corresponding coefficient <italic>h<sub>m</sub>
</italic><italic>l</italic>
-mers in each sequence is <italic>L</italic>
and the number of sequences is <italic>N</italic>
, this algorithm would require <italic>O</italic>
(<italic>N<sup>2</sup>
L</italic><sup>2</sup>
) comparisons. In addition, a naive algorithm for counting the number of mismatches between two <italic>l</italic>
-mers (i.e. the hamming distance) would be <italic>O</italic>
(<italic>l</italic>
). Our implementation employs bitwise operators, providing a constant-factor speedup. Briefly, using two bits to represent each base (A,C,G and T), we used an integer variable to represent non-overlapping substrings of <italic>t</italic>
base pairs of the <italic>l</italic>
-mer, therefore using total <inline-formula><inline-graphic xlink:href="pcbi.1003711.e022.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
-mer, where <inline-formula><inline-graphic xlink:href="pcbi.1003711.e023.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
-mers and use a precomputed look-up table to obtain the total number of mismatches using the XOR result. This method requires a look-up table of size 2<sup>2<italic>t</italic>
</sup><italic>t</italic>
depends on the processor architecture and amount of cache memory. We used <italic>t = </italic>
6 for our analysis.</p><sec id="s4c"><title>Gkm-kernel with <italic>k</italic>
-mer tree data structure</title><p>As depicted in <xref ref-type="fig" rid="pcbi-1003711-g005">Figure 5</xref>
, we use a <italic>k</italic>
-mer tree to hold all the <italic>l</italic>
-mers in the collection of all of the sequences. We construct the tree by adding a path for every <italic>l</italic>
-mer observed in a training sequence. Each node <italic>t<sub>i</sub>
</italic><italic>d</italic>
represents a sub-sequence of length <italic>d</italic>
, denoted by <italic>s</italic>
(<italic>t<sub>i</sub>
</italic><italic>t<sub>i</sub>
</italic><italic>l</italic>
-mer, and holds the list of training sequence labels in which that <italic>l</italic>
-mer appeared and the number of times that <italic>l</italic>
-mer appeared in each sequence. As an example, <xref ref-type="fig" rid="pcbi-1003711-g005">Figure 5</xref>
shows the tree that stores all the substrings of length <italic>l</italic>
= 3 in three sequences <italic>S<sub>1</sub>
</italic><italic>S<sub>2</sub>
</italic><italic>S<sub>3</sub>
</italic><xref rid="pcbi.1003711-Cormen1" ref-type="bibr">[35]</xref>
order. In contrast to the mismatch tree used in Ref. <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
, here for each node <italic>t<sub>i</sub>
</italic><italic>d</italic>
, we store the list of pointers to all the nodes <italic>t<sub>j</sub>
</italic><italic>d</italic>
for which <italic>s</italic>
(<italic>t<sub>i</sub>
</italic><italic>s</italic>
(<italic>t<sub>j</sub>
</italic><italic>l</italic>
– <italic>k</italic>
number of mismatches. We also store the number of mismatches between <italic>s</italic>
(<italic>t<sub>i</sub>
</italic><italic>s</italic>
(<italic>t<sub>j</sub>
</italic><xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
, we do not need to store these values for all the nodes in the tree, but we compute them recursively as we traverse the tree. When reaching a leaf node, we increment the corresponding mismatch profile <italic>N<sub>m</sub>
</italic><italic>S<sub>i</sub>
</italic><italic>S<sub>j</sub>
</italic><italic>S<sub>i</sub>
</italic><italic>S<sub>j</sub>
</italic><fig id="pcbi-1003711-g005" orientation="portrait" position="float"><object-id pub-id-type="doi">10.1371/journal.pcbi.1003711.g005</object-id>
<label>Figure 5</label>
<caption><title>Fast computation of mismatch profiles using <italic>k</italic>
-mer tree structure.</title><p>As an example, we use <italic>l</italic>
= 3 and three sequences <italic>S<sub>1</sub>
</italic><italic>S<sub>2</sub>
</italic><italic>S<sub>3</sub>
</italic><italic>k</italic>
-mer tree. The leaves (nodes at depth <italic>d</italic>
= <italic>l</italic>
= 3) correspond to 3-mers AAA, AAC, ACC, and CCC. The sequence ID and the number of times each 3-mer appeared in each sequence are stored for each leaf. Each node <italic>t<sub>i</sub>
</italic><italic>d</italic>
represents a sequence of length <italic>d</italic>
, denoted by <italic>s</italic>
(<italic>t<sub>i</sub>
</italic><italic>t<sub>i</sub>
</italic><italic>s</italic>
(<italic>t<sub>2</sub>
</italic><italic>s</italic>
(<italic>t<sub>4</sub>
</italic><italic>t<sub>0</sub>
</italic><italic>t<sub>i</sub>
</italic><italic>d</italic>
, we compute the list of all the nodes <italic>t<sub>j</sub>
</italic><italic>d</italic>
for which <italic>s</italic>
(<italic>t<sub>i</sub>
</italic><italic>s</italic>
(<italic>t<sub>j</sub>
</italic><italic>m<sub>max</sub>
</italic><italic>s</italic>
(<italic>t<sub>i</sub>
</italic><italic>s</italic>
(<italic>t<sub>j</sub>
</italic><italic>N<sub>m</sub>
</italic><italic>S<sub>i</sub>
</italic><italic>S<sub>j</sub>
</italic><italic>S<sub>i</sub>
</italic><italic>S<sub>j</sub>
</italic><graphic xlink:href="pcbi.1003711.g005"></graphic>
</fig><p>To increase the speed further, we also introduce an optional parameter <italic>m<sub>max</sub>
</italic><italic>m<sub>max</sub>
</italic><italic>l</italic>
– <italic>k</italic>
, we only consider <italic>l</italic>
-mer pairs that have at most <italic>m<sub>max</sub>
</italic><italic>l</italic>
-mer pairs which potentially contribute less to the overall similarity scores. This method provides fast and efficient approximations of the exact solution. In addition, we only compute the lower triangle of the matrix because of the symmetry in the kernel matrix. Hence, at each node <italic>t<sub>i</sub>
</italic><italic>t<sub>j</sub>
</italic><italic>t<sub>i</sub>
</italic><italic>t<sub>j</sub>
</italic><italic>t<sub>i</sub>
</italic><italic>t<sub>j</sub>
</italic><italic>t<sub>i</sub>
</italic><italic>t<sub>j</sub>
</italic><sec id="s4d"><title>Analysis of <italic>de novo</italic>
PWMs from gkm-SVM</title><p>We developed a new method for building <italic>de novo</italic>
PWMs by systematically merging the most predictive <italic>k</italic>
-mers from a trained gkm-SVM. We first determined a set of predictive <italic>k</italic>
-mers by scoring all possible 10-mers and selecting the top 1% of the high-scoring 10-mers. We then found a set of distinct PWM models from these predictive 10-mers using a heuristic iterated greedy algorithm. Specifically, we first built an initial PWM model from the highest scoring 10-mer. Then, for each of the remaining predictive 10-mers, we calculated the log-odd ratios of all possible alignments of the 10-mer to the PWM model, and identified the best alignment (i.e. the position and the orientation that give rise to the highest log-odd ratio value). Since multiple distinct classes of TFBSs are expected to be identified in most cases, we only considered 10-mers with good alignments (i.e. we used threshold of 5.0 for log-odd ratio scores relative to a genomic GC = 0.42 background). After each of the 10-mers was aligned, we updated the PWM model only with successfully aligned 10-mers. To further refine the PWM, we repeated this by iterating through all of the top 1% 10-mers until no changes were made. When updating the PWM model, we assumed that the contribution of each <italic>k</italic>
-mer is exponentially weighted proportional to its SVM score, using exp(α <italic>w<sub>i</sub>
</italic><sup>st</sup>
PWM were then removed from the list, and the process was repeated on the remaining predictive <italic>k</italic>
-mers, to find up to three PWMs. Lastly, we matched our PWMs to the previously identified PWMs <xref rid="pcbi.1003711-Wang1" ref-type="bibr">[18]</xref>
using TOMTOM <xref rid="pcbi.1003711-Gupta1" ref-type="bibr">[36]</xref>
software. Each of the PWMs identified by our method were associated with Ref. <xref rid="pcbi.1003711-Wang1" ref-type="bibr">[18]</xref>
PWMs if the <italic>q</italic>
-value (false discovery rate) <0.05.</p><sec id="s4e"><title>Implementation of mismatch and wildcard kernels using gkm-kernel framework</title>
<p>In the gkm-kernel, we define the feature vector to consist of the frequency of all the <italic>l</italic>
-mers with exactly <italic>k</italic>
known bases and <italic>l</italic>
– <italic>k</italic>
gaps. In contrast, the wildcard kernel <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
also includes all the <italic>l</italic>
-mers with <italic>l</italic>
– <italic>k</italic>
wildcards, where <italic>l</italic>
– <italic>k</italic>
ranges from 0 to the maximum number of wildcards allowed, <italic>M</italic>
. Thus in the wildcard kernel, the parameter <italic>M</italic>
replaces <italic>k</italic>
in our gkm-kernel. In the sum, these are weighted by λ<italic><sup>l - k</sup>
</italic><xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
. We derived an equation to directly compute the inner products from the mismatch profiles without the need to calculate the actual gapped <italic>k</italic>
-mer counts. Here we show that a similar approach can be used to calculate the wildcard kernel. We derive a new set of coefficients <inline-formula><inline-graphic xlink:href="pcbi.1003711.e024.jpg"></inline-graphic>
</inline-formula><italic>h<sub>m</sub>
</italic><inline-formula><inline-graphic xlink:href="pcbi.1003711.e025.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
-mers with <italic>m</italic>
mismatches in the inner product of the corresponding feature vectors of the two sequences. Equation (8) gives those weights:<disp-formula><graphic xlink:href="pcbi.1003711.e026.jpg" position="float" orientation="portrait"></graphic>
<label>(8)</label>
</disp-formula><p>Using the above form allows us to directly use the fast algorithms we have developed for calculation of the mismatch profiles to calculate the wildcard kernels. Although there are similarities between our tree algorithm and the tree algorithm described in Ref. <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
, there are some key differences. In the Ref. <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
, the algorithm literally transverses all the possible gapped <italic>l</italic>
-mers (with maximum <italic>M</italic>
number of gaps) while our algorithm takes advantage of the fact that the final inner product will only depend on the number of pairwise mismatches and hence only traverses all the <italic>l</italic>
-mers that are present in the input data. Another difference is that Ref. <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
uses a list of all partially matching <italic>l</italic>
-mers at each node of the tree, while we use a list of pointers to tree nodes instead. So, for example, at the beginning of the algorithm (at depth <italic>d</italic>
= 0) they start with a large list consisting of all the possible <italic>l</italic>
-mers in the input data, while in our algorithm the list at depth <italic>d</italic>
= 0 consists of only one node (the root of the tree). Using this representation of all the partially matching <italic>l</italic>
-mers, we can more efficiently perform the comparisons at each step of the algorithm when the tree is dense.</p><p>In the mismatch string kernel described in Ref. <xref rid="pcbi.1003711-Leslie1" ref-type="bibr">[8]</xref>
and <xref rid="pcbi.1003711-Leslie2" ref-type="bibr">[9]</xref>
, the feature vectors consist of the counts for all the <italic>l</italic>
-mers with maximum distance <italic>M</italic>
from the <italic>l</italic>
-mers in the sequence. Here we show that the approach above can be used to implement the mismatch kernel. Again, the only difference is in the set of weights used in Equation (3). To calculate the mismatch string kernel value for two sequences we replace <italic>h<sub>lk</sub>
</italic><italic>m</italic>
) in Equation (3) by <inline-formula><inline-graphic xlink:href="pcbi.1003711.e027.jpg"></inline-graphic>
</inline-formula><disp-formula><graphic xlink:href="pcbi.1003711.e028.jpg" position="float" orientation="portrait"></graphic>
<label>(9)</label>
</disp-formula><italic>b</italic>
is the alphabet size (<italic>b</italic>
= 4 for DNA sequences) and <italic>r</italic>
= <italic>m</italic>
<sub>1</sub>
+<italic>m</italic>
<sub>2</sub>
−<italic>m</italic>
−2<italic>t</italic>
. Given two <italic>l</italic>
-mers <italic>x</italic>
<sub>1</sub>
and <italic>x</italic>
<sub>2</sub>
where <italic>x</italic>
<sub>1</sub>
and <italic>x</italic>
<sub>2</sub>
differ in exactly <italic>m</italic>
places, the term inside the summations counts the number of all possible <italic>l</italic>
-mers that exactly differ <italic>x</italic>
<sub>1</sub>
in <italic>m</italic>
<sub>1</sub>
places and <italic>x</italic>
<sub>2</sub>
in <italic>m</italic>
<sub>2</sub>
places <italic>t</italic>
of which fall in the common <italic>l</italic>
-<italic>m</italic>
bases of <italic>x</italic>
<sub>1</sub>
and <italic>x</italic>
<sub>2</sub>
. (See <xref ref-type="supplementary-material" rid="pcbi.1003711.s009">Figure S9</xref>
). So the result of the summation is the number of all <italic>l</italic>
-mers that differ <italic>x</italic>
<sub>1</sub>
and <italic>x</italic>
<sub>2</sub>
in at most <italic>M</italic>
places. This form for the mismatch string kernel has the advantage that we can directly use equation (3) to compute the kernels by only having the mismatch profiles that can be computed more efficiently.</p><sec id="s4f"><title>Gkm-filters for computation of the robust <italic>l</italic>
-mer count estimates</title><p>To compute the <italic>l</italic>
-mer count estimates by using Equation (5), one should first calculate the gapped <italic>k</italic>
-mer counts, <italic>y<sub>i</sub>
</italic><italic>y<sub>i</sub>
</italic><xref ref-type="supplementary-material" rid="pcbi.1003711.s010">Figure S10</xref>
. The mapping from observed <italic>l</italic>
-mer counts to gapped <italic>k</italic>
-mer counts is performed by the matrix <italic>A</italic>
, whose elements are <italic>a<sub>ij</sub>
</italic><italic>k</italic>
-mer <italic>v<sub>i</sub>
</italic><italic>l</italic>
-mer <italic>u<sub>j</sub>
</italic><italic>a<sub>ij</sub>
</italic><italic>a<sub>ij</sub>
</italic><italic>W</italic>
, which performs the mapping from gapped <italic>k</italic>
-mer counts to <italic>estimated l</italic>
-mer counts, and whose elements are <italic>w<sub>ij</sub>
</italic><italic>W</italic>
is the Penrose-Moore pseudo-inverse of <italic>A</italic>
<xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
. The element <italic>w<sub>ij</sub>
</italic><italic>l</italic>
-mer <italic>u<sub>i</sub>
</italic><italic>k</italic>
-mer <italic>v<sub>j</sub>
</italic><italic>A</italic>
and <italic>W</italic>
, and directly calculate the minimum norm <italic>l</italic>
-mer count estimates from actual <italic>l</italic>
-mer counts in a sequence. We refer to this combined mapping as the <italic>gkm-filter</italic>
. The combined mapping matrix <italic>G</italic>
= <italic>WA</italic>
, has elements <italic>g<sub>ij</sub>
</italic><xref ref-type="supplementary-material" rid="pcbi.1003711.s010">Figure S10</xref>
. As shown below, <italic>g<sub>ij</sub>
</italic><italic>m</italic>
, between the <italic>l</italic>
-mers <italic>u<sub>i</sub>
</italic><italic>u<sub>j</sub>
</italic><italic>g<sub>lk</sub>
</italic><italic>m</italic>
) and refer to this as the <italic>gkm-filter</italic>
since the domain and range of this mapping is the same.</p><p>To obtain the element <italic>g<sub>lk</sub>
</italic><italic>m</italic>
), that gives the weight for the contribution of an observed <italic>l</italic>
-mer <italic>u<sub>i</sub>
</italic><italic>l</italic>
-mer count estimate <italic>u<sub>j</sub>
</italic><italic>m</italic>
mismatches with <italic>u<sub>i</sub>
</italic><italic>k</italic>
-mers <italic>v<sub>τ</sub>
</italic><italic>u<sub>i</sub>
</italic><italic>a<sub>ij</sub>
</italic><italic>k</italic>
-mers. There exist <inline-formula><inline-graphic xlink:href="pcbi.1003711.e029.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
-mers that match <italic>u<sub>i</sub>
</italic><italic>m</italic>
mismatches with <italic>u<sub>j</sub>
</italic><xref ref-type="supplementary-material" rid="pcbi.1003711.s011">Figure S11</xref>
shows how we enumerate all these gapped <italic>k</italic>
-mers. The black solid circles denote the <italic>m</italic>
mismatch positions of <italic>u<sub>i</sub>
</italic><italic>u<sub>j</sub>
</italic><italic>l – m</italic>
match positions and the empty dotted circles denote the <italic>l – k</italic>
gap positions. For a gapped <italic>k</italic>
-mer to have exactly <italic>t</italic>
mismatches with <italic>u<sub>j</sub>
</italic><inline-formula><inline-graphic xlink:href="pcbi.1003711.e030.jpg"></inline-graphic>
</inline-formula><italic>t</italic>
mismatch positions and <inline-formula><inline-graphic xlink:href="pcbi.1003711.e031.jpg"></inline-graphic>
</inline-formula><italic>k – t</italic>
match positions. Now considering the weight <italic>w</italic>
(<italic>t</italic>
) for the gapped <italic>k</italic>
-mers with <italic>t</italic>
mismatches, the gapped <italic>k</italic>
-mer filter elements, <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) can be obtained as follows:<disp-formula><graphic xlink:href="pcbi.1003711.e032.jpg" position="float" orientation="portrait"></graphic>
<label>(10)</label>
</disp-formula><p>In other words, there are <inline-formula><inline-graphic xlink:href="pcbi.1003711.e033.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
-mer that matches <italic>u<sub>i</sub>
</italic><italic>t</italic>
mismatches with <italic>u<sub>j</sub>
</italic><italic>t</italic>
positions from the <italic>m</italic>
mismatch positions and <italic>k</italic>
– <italic>t</italic>
positions from the <italic>l</italic>
– <italic>m</italic>
match positions as explained above (<xref ref-type="supplementary-material" rid="pcbi.1003711.s011">Figure S11</xref>
). It can be easily shown that <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) is a polynomial of degree <italic>k</italic>
in <italic>m</italic>
. Now using the weights given in Equation (10), for any given <italic>l</italic>
-mer, <italic>u</italic>
we finally obtain the minimum norm <italic>l</italic>
-mer count estimate as follows:<disp-formula><graphic xlink:href="pcbi.1003711.e034.jpg" position="float" orientation="portrait"></graphic>
<label>(11)</label>
</disp-formula><italic>N</italic>
<sub>tr</sub>
(<italic>u</italic>
, <italic>m</italic>
) is the number of <italic>l</italic>
-mers with exactly <italic>m</italic>
mismatches with <italic>u</italic>
in the training set. For large values of <italic>l</italic>
and <italic>k</italic>
, the number of all possible gapped <italic>k</italic>
-mers gets exponentially large and since this method avoids evaluating the gapped <italic>k</italic>
-mer counts, it significantly reduces the cost of calculating the <italic>l</italic>
-mer count estimates compared to the original method we developed in Ref. <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
.</p><p>In summary, we defined a generalized <italic>k</italic>
-mer count (referred to as the robust <italic>l</italic>
-mer count estimates) by giving a non-zero weight to <italic>l</italic>
-mers with few number of mismatches (In the conventional <italic>k</italic>
-mer count only perfectly matching <italic>k</italic>
-mers are counted). These weights are given by <italic>g<sub>lk</sub>
</italic><italic>m</italic>
). <xref ref-type="supplementary-material" rid="pcbi.1003711.s012">Figure S12A</xref>
shows the plots for <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) for <italic>l</italic>
= 20 and various values of <italic>k</italic>
. Each plot is normalized so that weight corresponding to zero mismatches is equal to one. The case with <italic>l</italic>
= <italic>k</italic>
is equivalent to the conventional <italic>k</italic>
-mer count. Also <xref ref-type="supplementary-material" rid="pcbi.1003711.s012">Figure S12B</xref>
shows <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) for <inline-formula><inline-graphic xlink:href="pcbi.1003711.e035.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
. With a fixed length <italic>l</italic>
, higher values of <italic>k</italic>
result in smaller coefficients for larger mismatches, and therefore less smoothing of the estimated counts (<xref ref-type="supplementary-material" rid="pcbi.1003711.s012">Figure S12</xref>
). Moreover, <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) can become slightly negative for large numbers of mismatches. This is because in our estimation of the frequencies we did not restrict the frequencies to be positive, and doing so would yield a more complicated expression. The assumed Gaussian distribution allows non-physical negative frequencies to have non-zero probability. A beta-distribution would not have this problem but would introduce offsetting complications. In cases where the estimated counts are required to be strictly positive, such as when we need to calculate the logarithm or ratios of the estimated frequencies, we truncate the gkm-filter <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) by setting <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) = 0 for every <italic>m</italic>
≥<italic>m</italic>
<sub>0</sub>
, where <italic>m</italic>
<sub>0</sub>
is the smallest number of mismatches for which <italic>g<sub>lk</sub>
</italic><italic>m</italic>
<sub>0</sub>
)<0. This will give an approximation to the value of <inline-formula><inline-graphic xlink:href="pcbi.1003711.e036.jpg"></inline-graphic>
</inline-formula><sec id="s4g"><title>Gkm-kernel with <italic>l</italic>
-mer count estimates</title><p>Given a sequence <italic>S</italic>
, we define an <italic>l</italic>
-mer count estimate vector <inline-formula><inline-graphic xlink:href="pcbi.1003711.e037.jpg"></inline-graphic>
</inline-formula><italic>N</italic>
is the number of all <italic>l</italic>
-mers (4<italic><sup>l</sup>
</italic><inline-formula><inline-graphic xlink:href="pcbi.1003711.e038.jpg"></inline-graphic>
</inline-formula><italic>i</italic>
<sup>th</sup>
<italic>l</italic>
-mer appearing in sequence <italic>S</italic>
using Equation (11). Then, we can calculate a standard linear kernel simply by using this vector in Equation (1). Similar to the gkm-kernel method, we can further simplify this equation using the same technique introduced in Equation (2) which does not involve the computation of individual <italic>l</italic>
-mer estimates. We show that the inner product of the two <italic>l</italic>
-mer count estimate vectors can be obtained as follows:<disp-formula><graphic xlink:href="pcbi.1003711.e039.jpg" position="float" orientation="portrait"></graphic>
<label>(12)</label>
</disp-formula><italic>n</italic>
<sub>1</sub>
and <italic>n</italic>
<sub>2</sub>
are the number of <italic>l</italic>
-mers in <italic>S</italic>
<sub>1</sub>
and <italic>S</italic>
<sub>2</sub>
, and <inline-formula><inline-graphic xlink:href="pcbi.1003711.e040.jpg"></inline-graphic>
</inline-formula><italic>i</italic>
'th <italic>l</italic>
-mer in <italic>S</italic>
<sub>1</sub>
and <inline-formula><inline-graphic xlink:href="pcbi.1003711.e041.jpg"></inline-graphic>
</inline-formula><italic>j</italic>
'th <italic>l</italic>
-mer in S<sub>2</sub>
. If <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
<sub>2</sub>
have exactly <italic>m</italic>
mismatches then <italic>c</italic>
(<italic>u</italic>
<sub>1</sub>
, <italic>u</italic>
<sub>2</sub>
) = <italic>c<sub>m</sub>
</italic><italic>l</italic>
-mer pairs with <italic>m</italic>
mismatches, we can rewrite Equation (12) as follows:<disp-formula><graphic xlink:href="pcbi.1003711.e042.jpg" position="float" orientation="portrait"></graphic>
<label>(13)</label>
</disp-formula><italic>N<sub>m</sub>
</italic><italic>S</italic>
<sub>1</sub>
, <italic>S</italic>
<sub>2</sub>
) is the mismatch profile of <italic>S</italic>
<sub>1</sub>
and <italic>S</italic>
<sub>2</sub>
as previously defined in Equation (3). We show that the weight <italic>c<sub>lk</sub>
</italic><italic>m</italic>
), denoted in short by <italic>c<sub>m</sub>
</italic><disp-formula><graphic xlink:href="pcbi.1003711.e043.jpg" position="float" orientation="portrait"></graphic>
<label>(14)</label>
</disp-formula><italic>r</italic>
= <italic>m</italic>
<sub>1</sub>
+<italic>m</italic>
<sub>2</sub>
−2<italic>t</italic>
−<italic>m</italic>
, <italic>b</italic>
is the alphabet size. The summations are taken over the range 0 to <italic>l</italic>
. <xref ref-type="supplementary-material" rid="pcbi.1003711.s013">Figure S13</xref>
shows how we obtained the equation for <italic>c<sub>m</sub>
</italic><xref ref-type="supplementary-material" rid="pcbi.1003711.s011">Figure S11</xref>
. Given two <italic>l</italic>
-mers <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
<sub>2</sub>
, with <inline-formula><inline-graphic xlink:href="pcbi.1003711.e044.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
– <italic>m</italic>
matched positions, we want to enumerate the number of all possible <italic>l</italic>
-mers, <italic>u</italic>
, that have <italic>m</italic>
<sub>1</sub>
mismatches with <italic>u</italic>
<sub>1</sub>
and <italic>m</italic>
<sub>2</sub>
mismatches with <italic>u</italic>
<sub>2</sub>
. For this, we assume that <italic>t</italic>
of the <italic>m</italic>
<sub>1</sub>
mismatches are among the <italic>l</italic>
– <italic>m</italic>
match positions and <italic>m</italic>
<sub>1</sub>
– <italic>t</italic>
of them are among the <italic>m</italic>
mismatch positions. There are <inline-formula><inline-graphic xlink:href="pcbi.1003711.e045.jpg"></inline-graphic>
</inline-formula><italic>m</italic>
<sub>1</sub>
positions and (<italic>b</italic>
– 1)<italic><sup>t</sup>
</italic><italic>t</italic>
mismatches. These <italic>t</italic>
mismatches plus the <italic>m</italic>
−(<italic>m</italic>
<sub>1</sub>
−<italic>t</italic>
) unselected mismatch positions also do not match <italic>u</italic>
<sub>2</sub>
. Then, for the remaining <italic>r</italic>
= <italic>m</italic>
<sub>2</sub>
−(<italic>t</italic>
+<italic>m</italic>
−(<italic>m</italic>
<sub>1</sub>
−<italic>t</italic>
)) mismatches for <italic>u</italic>
<sub>2</sub>
there are <inline-formula><inline-graphic xlink:href="pcbi.1003711.e046.jpg"></inline-graphic>
</inline-formula><italic>b</italic>
– 2)<italic><sup>r</sup>
</italic><italic>l</italic>
-mers, <italic>u</italic>
with <italic>m</italic>
<sub>1</sub>
mismatches with <italic>u</italic>
<sub>1</sub>
and <italic>m</italic>
<sub>2</sub>
mismatches with <italic>u</italic>
<sub>2</sub>
, where <italic>t</italic>
of the mismatches of <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
are among the <italic>l</italic>
– <italic>m</italic>
match positions of <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
<sub>2</sub>
is given by <inline-formula><inline-graphic xlink:href="pcbi.1003711.e047.jpg"></inline-graphic>
</inline-formula><p>Using matrix notation, we can further show that <italic>c<sub>m</sub>
</italic><italic>g<sub>m</sub>
</italic><italic>g<sub>lk</sub>
</italic><italic>m</italic>
). To see this, note that <inline-formula><inline-graphic xlink:href="pcbi.1003711.e048.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e049.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e050.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
-mer count vectors for <italic>S</italic>
<sub>1</sub>
and <italic>S</italic>
<sub>2</sub>
. Given <italic>G</italic>
= <italic>WA</italic>
, we have <inline-formula><inline-graphic xlink:href="pcbi.1003711.e051.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e052.jpg"></inline-graphic>
</inline-formula><italic>A</italic>
is the binary incidence matrix that maps <italic>l</italic>
-mer counts to gapped <italic>k</italic>
-mer counts as defined in Ref. <xref rid="pcbi.1003711-Ghandi1" ref-type="bibr">[7]</xref>
and <italic>W</italic>
is the Moore-Penrose pseudo-inverse of <italic>A</italic>
. Note that this result does not hold for the truncated filter <italic>g<sub>m</sub>
</italic><italic>c<sub>m</sub>
</italic><sec id="s4h"><title>ROC curves</title>
<p>To compare the performance of different classification methods, we calculated the area under the receiver operating characteristic (ROC) curve for each classifier. To plot the ROC curves and calculate area under the curves (AUCs) we used the ROCR package <xref rid="pcbi.1003711-Sing1" ref-type="bibr">[37]</xref>
in R.</p><sec id="s4i"><title>Cross validation</title>
<p>Following standard five-fold cross validation procedures, we divided the positive and negative sets into five segments, left one segment out as the test set and used the other four segments for training. We repeated for all of the five segments and calculated the mean and standard error of the prediction accuracy on the test set elements.</p>
</sec><sec id="s4j"><title>ENCODE ChIP-seq datasets</title>
<p>The ENCODE ChIP-seq datasets were downloaded from <ext-link ext-link-type="ftp" xlink:href="ftp://ftp.ebi.ac.uk/pub/databases/ensembl/encode/integration_data_jan2011/byDataType/peaks/jan2011/spp/optimal/hub/">ftp://ftp.ebi.ac.uk/pub/databases/ensembl/encode/integration_data_jan2011/byDataType/peaks/jan2011/spp/optimal/hub/</ext-link>
.</p><sec id="s4k"><title>Implementation and source code</title>
<p>We have implemented these algorithms in C++, and the source code and executable files are available on our website at <ext-link ext-link-type="uri" xlink:href="http://www.beerlab.org/gkmsvm/">http://www.beerlab.org/gkmsvm/</ext-link>
.</p><sec sec-type="supplementary-material" id="s5"><title>Supporting Information</title>
<supplementary-material content-type="local-data" id="pcbi.1003711.s001"><label>Figure S1</label>
<caption><p><bold>PWM model for CTCF binding sites.</bold>
CTCF specifically binds to a set of very long sequences via its eleven zinc finger domains, which can be effectively modeled by a PWM. This CTCF logo was directly obtained from the JASPAR database <xref rid="pcbi.1003711-Bryne1" ref-type="bibr">[38]</xref>
(available at <ext-link ext-link-type="uri" xlink:href="http://http://jaspar.cgb.ki.se/">http://http://jaspar.cgb.ki.se/</ext-link>
).</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s001.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s002"><label>Figure S2</label>
<caption><p><bold>Classification results for CTCF binding using the CTCF PWM.</bold>
The CTCF bound regions and the corresponding negative regions were scored by the CTCF PWM and the best log-odd score for each sequence was then used to calculate the ROC curve. Extremely high AUC was achieved, indicating that CTCF binding is well-modeled by the PWM.</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s002.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s003"><label>Figure S3</label>
<caption><p><bold>Precision of gkm-SVM is significantly higher than kmer-SVM.</bold>
To extrapolate to larger negative sets, we re-trained both gkm-SVM and kmer-SVM on each of the positive data sets (CTCF and EP300) against a 10× larger negative set. We independently selected the parameter <italic>l</italic>
and <italic>k</italic>
which exhibited the best performance when trained on the 1× negative set as shown in <xref ref-type="fig" rid="pcbi-1003711-g002">Figure 2C and D</xref>
. We additionally applied <italic>m<sub>max</sub>
</italic><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s003.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s004"><label>Figure S4</label>
<caption><p><bold>gkm-SVM is generally more robust than kmer-SVM.</bold>
We calculated standard deviation (SD) of the AUCs from the test CV sets for both gkm-SVM and kmer-SVM. In most cases, gkm-SVM AUC SD is significantly smaller than kmer-SVM AUC SD.</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s004.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s005"><label>Figure S5</label>
<caption><p><bold><italic>de novo</italic>
</bold><bold> PWMs from gkm-SVMs trained on the 467 ENCODE ChIP-seq data sets.</bold>
The name of the best matching known PWM (from Wang <italic>et al.</italic>
) was assigned to each of the PWMs.</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s005.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s006"><label>Figure S6</label>
<caption><p><bold>Classification results on the mouse forebrain EP300 data set with various methods.</bold>
For each of the methods, we examined combinations of the parameters and measured AUCs with 5-CVs. Red dash-dotted line in each plot denotes the best AUC achieved by gkm-kernel with <italic>l</italic>
= 9 and <italic>k</italic>
= 6. (A) Mismatch kernel for <italic>k</italic>
= 6∼12 and <italic>m</italic>
= 1∼2 using original implementation. Note that we obtained the kernel of (<italic>k</italic>
, <italic>m</italic>
) = (6, 2) using our implementation due to the prohibitive computing time of the original method. (B) The previous experiments were repeated using our new implementation with the “adding reverse complement sequences option” enabled. (C) Wildcard kernel for <italic>k</italic>
= 6∼12 and <italic>m</italic>
= 1∼3 with <italic>λ</italic>
= 1. Note that we obtained the kernels with <italic>k</italic>
= 6, <italic>m</italic>
= 2∼3 using our implementation. (D) The previous experiments were repeated using our new implementation with the reverse complement sequence option enabled. (E) di-mismatch kernels for <italic>k</italic>
= 6∼10 and <italic>m</italic>
= 1∼3. (F) The previous experiments using only the top 1000 most discriminative features as recommended in the original study.</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s006.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s007"><label>Figure S7</label>
<caption><p><bold>Comparisons of running times between different methods.</bold>
(A) gkm-kernel, (B) mismatch-kernel, (C) wildcard-kernel, and (D) di-mismatch kernel. For mismatch and wildcard we also show results using our tree structure (dashed). For consistency, we used a single machine equipped with Intel Core i5-2410M (2.30 GHz) processor and 6 GB RAM, except di-mismatch kernels. Due to the prohibitive memory requirement of the di-mismatch kernels for large <italic>k</italic>
, we separately measured the running times on different machines.</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s007.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s008"><label>Figure S8</label>
<caption><p><bold>Comparison of different filters in gkm-SVM on human ENCODE ChIP-seq data sets.</bold>
We compared the performance of the gkm-SVM using the gapped kmers as features (gkm-SVM) to gkm-SVM using <italic>l</italic>
-mer count estimates (A) with truncated gkm-filter and (B) with full gkm-filter. We used <italic>l</italic>
= 10 and <italic>k</italic>
= 6 for all the methods. Those data sets where using the truncated-filter gives higher AUCs are marked as purple circles. The truncated filter method is marginally but systematically better when AUC is greater than 0.9.</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s008.pdf"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s009"><label>Figure S9</label>
<caption><p><bold>Deriving the weights for calculation of the mismatch kernel using </bold>
Equation (3)<bold>.</bold>
x<sub>1</sub>
and x<sub>2</sub>
differ in m places. u differs x<sub>1</sub>
in m<sub>1</sub>
places and x<sub>2</sub>
in m<sub>2</sub>
places. t of the u mismatch places are among the l – m ,x<sub>1</sub>
, x<sub>2</sub>
common places. There are <inline-formula><inline-graphic xlink:href="pcbi.1003711.e053.jpg"></inline-graphic>
</inline-formula><sub>1</sub>
, m<sub>2</sub>
, t≤M.</p><p>(PNG)</p>
</caption><media xlink:href="pcbi.1003711.s009.png"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s010"><label>Figure S10</label>
<caption><p><bold>Block diagrams of the proposed method for the gkm-filter.</bold>
(Top) Gapped <italic>k</italic>
-mer counts are obtained from <italic>l</italic>
-mer counts in the training set. Then minimum norm <italic>l</italic>
-mer count estimates are obtained from the gapped <italic>k</italic>
-mer counts. The <italic>a<sub>ij</sub>
</italic><italic>A</italic>
, that maps the <italic>l</italic>
-mer counts in the training set to the gapped <italic>k</italic>
-mer counts. <italic>a<sub>ij</sub>
</italic><italic>k</italic>
-mer <italic>v<sub>i</sub>
</italic><italic>l</italic>
-mer <italic>u<sub>j</sub>
</italic><italic>w<sub>ij</sub>
</italic><italic>W</italic>
(the pseudo-inverse of <italic>A</italic>
) mapping gapped <italic>k</italic>
-mer frequencies to estimated <italic>l</italic>
-mer frequencies. (Bottom) We combine the two mapping matrices <italic>A</italic>
and <italic>W</italic>
to directly calculate the minimum norm <italic>l</italic>
-mer count estimates from the <italic>l</italic>
-mer counts in the training set. <italic>g<sub>ij</sub>
</italic><italic>G</italic>
mapping the <italic>l</italic>
-mer counts in the training set to the minimum norm <italic>l</italic>
-mer count estimates.</p><p>(PNG)</p>
</caption><media xlink:href="pcbi.1003711.s010.png"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s011"><label>Figure S11</label>
<caption><p><bold>Enumeration of gapped </bold>
<bold><italic>k</italic>
</bold><bold>-mers with exactly </bold>
<bold><italic>t</italic>
</bold><bold> mismatches.</bold>
Given the <italic>l</italic>
-mers <italic>u<sub>i</sub>
</italic><italic>u<sub>j</sub>
</italic><italic>k</italic>
-mer that matches <italic>u<sub>i</sub>
</italic><italic>t</italic>
mismatches with <italic>u<sub>j</sub>
</italic><inline-formula><inline-graphic xlink:href="pcbi.1003711.e054.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e055.jpg"></inline-graphic>
</inline-formula><italic>t</italic>
mismatch positions and <inline-formula><inline-graphic xlink:href="pcbi.1003711.e056.jpg"></inline-graphic>
</inline-formula><italic>k</italic>
– <italic>t</italic>
match positions. The black solid circles denote the <italic>m</italic>
mismatch positions of <italic>u<sub>i</sub>
</italic><italic>u<sub>j</sub>
</italic><italic>l</italic>
– <italic>m</italic>
match positions, and the empty dotted circles denote the <italic>l</italic>
– <italic>k</italic>
unselected (gap) positions.</p><p>(PNG)</p>
</caption><media xlink:href="pcbi.1003711.s011.png"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s012"><label>Figure S12</label>
<caption><p><bold>Plot of </bold>
<bold><italic>g<sub>lk</sub>
</italic><bold>(</bold>
<bold><italic>m</italic>
</bold><bold>).</bold>
Plot of the normalized filter function <italic>g<sub>lk</sub>
</italic><italic>m</italic>
) for (A) <italic>l</italic>
= 20 and various values of <italic>k</italic>
and (B) <italic>k</italic>
= 6 and various values of <italic>l</italic>
.</p><p>(PNG)</p>
</caption><media xlink:href="pcbi.1003711.s012.png"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s013"><label>Figure S13</label>
<caption><p><bold>Enumeration of </bold>
<bold><italic>l</italic>
</bold><bold>-mers with </bold>
<bold><italic>m<sub>1</sub>
</italic><bold> and </bold>
<bold><italic>m<sub>2</sub>
</italic><bold> mismatches.</bold>
Given two <italic>l</italic>
-mers <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
<sub>2</sub>
, with <italic>m</italic>
mismatches and <italic>l</italic>
– <italic>m</italic>
matched positions, we want to enumerate the number of all possible <italic>l</italic>
-mers, <italic>u′</italic>
, that have <italic>m</italic>
<sub>1</sub>
mismatches with <italic>u</italic>
<sub>1</sub>
and <italic>m</italic>
<sub>2</sub>
mismatches with <italic>u</italic>
<sub>2</sub>
. For this, we assume that <italic>t</italic>
of the <italic>m</italic>
<sub>1</sub>
mismatches are among the <italic>l</italic>
– <italic>m</italic>
match positions and <italic>m</italic>
<sub>1</sub>
– <italic>t</italic>
of them are among the <italic>m</italic>
mismatch positions. There are <inline-formula><inline-graphic xlink:href="pcbi.1003711.e057.jpg"></inline-graphic>
</inline-formula><italic>m</italic>
<sub>1</sub>
positions and <inline-formula><inline-graphic xlink:href="pcbi.1003711.e058.jpg"></inline-graphic>
</inline-formula><italic>t</italic>
mismatches. These <italic>t</italic>
mismatches plus the (<italic>m</italic>
−(<italic>m</italic>
<sub>1</sub>
−<italic>t</italic>
)) unselected mismatch positions are also mismatches for <italic>u</italic>
<sub>2</sub>
. For the remaining <italic>r</italic>
= <italic>m</italic>
<sub>2</sub>
−(<italic>t</italic>
+<italic>m</italic>
−(<italic>m</italic>
<sub>1</sub>
−<italic>t</italic>
)) mismatches for <italic>u</italic>
<sub>2</sub>
there are <inline-formula><inline-graphic xlink:href="pcbi.1003711.e059.jpg"></inline-graphic>
</inline-formula><inline-formula><inline-graphic xlink:href="pcbi.1003711.e060.jpg"></inline-graphic>
</inline-formula><italic>l</italic>
-mers, <italic>u′</italic>
, with <italic>m</italic>
<sub>1</sub>
mismatches with <italic>u</italic>
<sub>1</sub>
and <italic>m</italic>
<sub>2</sub>
mismatches with <italic>u</italic>
<sub>2</sub>
, where <italic>t</italic>
of the mismatches of <italic>u</italic>
<sub>1</sub>
and <italic>u′</italic>
are among the (<italic>l</italic>
– <italic>m</italic>
) match positions of <italic>u</italic>
<sub>1</sub>
and <italic>u</italic>
<sub>2</sub>
is given by <inline-formula><inline-graphic xlink:href="pcbi.1003711.e061.jpg"></inline-graphic>
</inline-formula><p>(PNG)</p>
</caption><media xlink:href="pcbi.1003711.s013.png"><caption><p>Click here for additional data file.</p>
</caption><supplementary-material content-type="local-data" id="pcbi.1003711.s014"><label>Text S1</label>
<caption><p><bold>Numerical example.</bold>
</p><p>(PDF)</p>
</caption><media xlink:href="pcbi.1003711.s014.pdf"><caption><p>Click here for additional data file.</p>
</caption><back><ack><p>We gratefully acknowledge contribution of source code for the wildcard and di-mismatch approaches from Rui Kuang and Christina Leslie for the direct comparisons of these methods.</p>
</ack><ref-list><title>Reference</title>
<ref id="pcbi.1003711-Manolio1"><label>1</label>
<mixed-citation publication-type="journal"><name><surname>Manolio</surname>
<given-names>TA</given-names>
</name><year>2010</year>
) <article-title>Genomewide Association Studies and Assessment of the Risk of Disease</article-title>
. <source>N Engl J Med</source>
<volume>363</volume>
: <fpage>166</fpage>
–<lpage>176</lpage>
<pub-id pub-id-type="doi">10.1056/NEJMra0905980</pub-id>
<pub-id pub-id-type="pmid">20647212</pub-id>
</mixed-citation><ref id="pcbi.1003711-Maurano1"><label>2</label>
<mixed-citation publication-type="journal"><name><surname>Maurano</surname>
<given-names>MT</given-names>
</name><name><surname>Humbert</surname>
<given-names>R</given-names>
</name><name><surname>Rynes</surname>
<given-names>E</given-names>
</name><name><surname>Thurman</surname>
<given-names>RE</given-names>
</name><name><surname>Haugen</surname>
<given-names>E</given-names>
</name><etal>et al</etal>
(<year>2012</year>
) <article-title>Systematic Localization of Common Disease-Associated Variation in Regulatory DNA</article-title>
. <source>Science</source>
<volume>337</volume>
: <fpage>1190</fpage>
–<lpage>1195</lpage>
<pub-id pub-id-type="doi">10.1126/science.1222794</pub-id>
<pub-id pub-id-type="pmid">22955828</pub-id>
</mixed-citation><ref id="pcbi.1003711-Lee1"><label>3</label>
<mixed-citation publication-type="journal"><name><surname>Lee</surname>
<given-names>D</given-names>
</name><name><surname>Karchin</surname>
<given-names>R</given-names>
</name><name><surname>Beer</surname>
<given-names>MA</given-names>
</name><year>2011</year>
) <article-title>Discriminative prediction of mammalian enhancers from DNA sequence</article-title>
. <source>Genome Res</source>
<volume>21</volume>
: <fpage>2167</fpage>
–<lpage>2180</lpage>
<pub-id pub-id-type="doi">10.1101/gr.121905.111</pub-id>
<pub-id pub-id-type="pmid">21875935</pub-id>
</mixed-citation><ref id="pcbi.1003711-Stormo1"><label>4</label>
<mixed-citation publication-type="journal"><name><surname>Stormo</surname>
<given-names>GD</given-names>
</name><year>2000</year>
) <article-title>DNA binding sites: representation and discovery</article-title>
. <source>Bioinformatics</source>
<volume>16</volume>
: <fpage>16</fpage>
–<lpage>23</lpage>
<pub-id pub-id-type="doi">10.1093/bioinformatics/16.1.16</pub-id>
<pub-id pub-id-type="pmid">10812473</pub-id>
</mixed-citation><ref id="pcbi.1003711-Beer1"><label>5</label>
<mixed-citation publication-type="journal"><name><surname>Beer</surname>
<given-names>MA</given-names>
</name><name><surname>Tavazoie</surname>
<given-names>S</given-names>
</name><year>2004</year>
) <article-title>Predicting Gene Expression from Sequence</article-title>
. <source>Cell</source>
<volume>117</volume>
: <fpage>185</fpage>
–<lpage>198</lpage>
.<pub-id pub-id-type="pmid">15084257</pub-id>
</mixed-citation><ref id="pcbi.1003711-FletezBrant1"><label>6</label>
<mixed-citation publication-type="journal"><name><surname>Fletez-Brant</surname>
<given-names>C</given-names>
</name><name><surname>Lee</surname>
<given-names>D</given-names>
</name><name><surname>McCallion</surname>
<given-names>AS</given-names>
</name><name><surname>Beer</surname>
<given-names>MA</given-names>
</name><year>2013</year>
) <article-title>kmer-SVM: a web server for identifying predictive regulatory sequence features in genomic data sets</article-title>
. <source>Nucleic Acids Res</source>
<volume>41</volume>
: <fpage>W544</fpage>
–<lpage>W556</lpage>
<pub-id pub-id-type="doi">10.1093/nar/gkt519</pub-id>
<pub-id pub-id-type="pmid">23771147</pub-id>
</mixed-citation><ref id="pcbi.1003711-Ghandi1"><label>7</label>
<mixed-citation publication-type="other">Ghandi M, Mohammad-Noori M, Beer MA (2013) Robust k-mer frequency estimation using gapped k-mers. J Math Biol: 1–32. doi:10.1007/s00285-013-0705-3.</mixed-citation>
</ref><ref id="pcbi.1003711-Leslie1"><label>8</label>
<mixed-citation publication-type="journal"><name><surname>Leslie</surname>
<given-names>C</given-names>
</name><name><surname>Eskin</surname>
<given-names>E</given-names>
</name><name><surname>Cohen</surname>
<given-names>A</given-names>
</name><name><surname>Weston</surname>
<given-names>J</given-names>
</name><name><surname>Noble</surname>
<given-names>WS</given-names>
</name><year>2004</year>
) <article-title>Mismatch string kernels for discriminative protein classification</article-title>
. <source>Bioinformatics</source>
<volume>20</volume>
: <fpage>467</fpage>
–<lpage>476</lpage>
<pub-id pub-id-type="doi">10.1093/bioinformatics/btg431</pub-id>
<pub-id pub-id-type="pmid">14990442</pub-id>
</mixed-citation><ref id="pcbi.1003711-Leslie2"><label>9</label>
<mixed-citation publication-type="journal"><name><surname>Leslie</surname>
<given-names>C</given-names>
</name><name><surname>Kuang</surname>
<given-names>R</given-names>
</name><year>2004</year>
) <article-title>Fast String Kernels using Inexact Matching for Protein Sequences</article-title>
. <source>J Mach Learn Res</source>
<volume>5</volume>
: <fpage>1435</fpage>
–<lpage>1455</lpage>
.</mixed-citation><ref id="pcbi.1003711-Agius1"><label>10</label>
<mixed-citation publication-type="journal"><name><surname>Agius</surname>
<given-names>P</given-names>
</name><name><surname>Arvey</surname>
<given-names>A</given-names>
</name><name><surname>Chang</surname>
<given-names>W</given-names>
</name><name><surname>Noble</surname>
<given-names>WS</given-names>
</name><name><surname>Leslie</surname>
<given-names>C</given-names>
</name><year>2010</year>
) <article-title>High Resolution Models of Transcription Factor-DNA Affinities Improve In Vitro and In Vivo Binding Predictions</article-title>
. <source>PLoS Comput Biol</source>
<volume>6</volume>
: <fpage>e1000916</fpage>
<pub-id pub-id-type="doi">10.1371/journal.pcbi.1000916</pub-id>
<pub-id pub-id-type="pmid">20838582</pub-id>
</mixed-citation><ref id="pcbi.1003711-Arvey1"><label>11</label>
<mixed-citation publication-type="journal"><name><surname>Arvey</surname>
<given-names>A</given-names>
</name><name><surname>Agius</surname>
<given-names>P</given-names>
</name><name><surname>Noble</surname>
<given-names>WS</given-names>
</name><name><surname>Leslie</surname>
<given-names>C</given-names>
</name><year>2012</year>
) <article-title>Sequence and chromatin determinants of cell-type–specific transcription factor binding</article-title>
. <source>Genome Res</source>
<volume>22</volume>
: <fpage>1723</fpage>
–<lpage>1734</lpage>
<pub-id pub-id-type="doi">10.1101/gr.127712.111</pub-id>
<pub-id pub-id-type="pmid">22955984</pub-id>
</mixed-citation><ref id="pcbi.1003711-Amanchy1"><label>12</label>
<mixed-citation publication-type="journal"><name><surname>Amanchy</surname>
<given-names>R</given-names>
</name><name><surname>Kandasamy</surname>
<given-names>K</given-names>
</name><name><surname>Mathivanan</surname>
<given-names>S</given-names>
</name><name><surname>Periaswamy</surname>
<given-names>B</given-names>
</name><name><surname>Reddy</surname>
<given-names>R</given-names>
</name><etal>et al</etal>
(<year>2011</year>
) <article-title>Identification of Novel Phosphorylation Motifs Through an Integrative Computational and Experimental Analysis of the Human Phosphoproteome</article-title>
. <source>J Proteomics Bioinform</source>
<volume>4</volume>
: <fpage>22</fpage>
–<lpage>35</lpage>
<pub-id pub-id-type="doi">10.4172/jpb.1000163</pub-id>
<pub-id pub-id-type="pmid">21720494</pub-id>
</mixed-citation><ref id="pcbi.1003711-McDaniell1"><label>13</label>
<mixed-citation publication-type="journal"><name><surname>McDaniell</surname>
<given-names>R</given-names>
</name><name><surname>Lee</surname>
<given-names>B-K</given-names>
</name><name><surname>Song</surname>
<given-names>L</given-names>
</name><name><surname>Liu</surname>
<given-names>Z</given-names>
</name><name><surname>Boyle</surname>
<given-names>AP</given-names>
</name><etal>et al</etal>
(<year>2010</year>
) <article-title>Heritable Individual-Specific and Allele-Specific Chromatin Signatures in Humans</article-title>
. <source>Science</source>
<volume>328</volume>
: <fpage>235</fpage>
–<lpage>239</lpage>
<pub-id pub-id-type="doi">10.1126/science.1184655</pub-id>
<pub-id pub-id-type="pmid">20299549</pub-id>
</mixed-citation><ref id="pcbi.1003711-Kim1"><label>14</label>
<mixed-citation publication-type="journal"><name><surname>Kim</surname>
<given-names>TH</given-names>
</name><name><surname>Abdullaev</surname>
<given-names>ZK</given-names>
</name><name><surname>Smith</surname>
<given-names>AD</given-names>
</name><name><surname>Ching</surname>
<given-names>KA</given-names>
</name><name><surname>Loukinov</surname>
<given-names>DI</given-names>
</name><etal>et al</etal>
(<year>2007</year>
) <article-title>Analysis of the Vertebrate Insulator Protein CTCF-Binding Sites in the Human Genome</article-title>
. <source>Cell</source>
<volume>128</volume>
: <fpage>1231</fpage>
–<lpage>1245</lpage>
<pub-id pub-id-type="doi">10.1016/j.cell.2006.12.048</pub-id>
<pub-id pub-id-type="pmid">17382889</pub-id>
</mixed-citation><ref id="pcbi.1003711-Visel1"><label>15</label>
<mixed-citation publication-type="journal"><name><surname>Visel</surname>
<given-names>A</given-names>
</name><name><surname>Blow</surname>
<given-names>MJ</given-names>
</name><name><surname>Zhang</surname>
<given-names>T</given-names>
</name><name><surname>Akiyama</surname>
<given-names>JA</given-names>
</name><name><surname>Holt</surname>
<given-names>A</given-names>
</name><etal>et al</etal>
(<year>2009</year>
) <article-title>ChIP-seq accurately predicts tissue-specific activity of enhancers</article-title>
. <source>Nature</source>
<volume>457</volume>
: <fpage>854</fpage>
–<lpage>858</lpage>
<pub-id pub-id-type="doi">10.1038/nature07730</pub-id>
<pub-id pub-id-type="pmid">19212405</pub-id>
</mixed-citation><ref id="pcbi.1003711-Zhang1"><label>16</label>
<mixed-citation publication-type="journal"><name><surname>Zhang</surname>
<given-names>Y</given-names>
</name><name><surname>Liu</surname>
<given-names>T</given-names>
</name><name><surname>Meyer</surname>
<given-names>CA</given-names>
</name><name><surname>Eeckhoute</surname>
<given-names>J</given-names>
</name><name><surname>Johnson</surname>
<given-names>DS</given-names>
</name><etal>et al</etal>
(<year>2008</year>
) <article-title>Model-based Analysis of ChIP-Seq (MACS)</article-title>
. <source>Genome Biol</source>
<volume>9</volume>
: <fpage>R137</fpage>
<pub-id pub-id-type="doi">10.1186/gb-2008-9-9-r137</pub-id>
<pub-id pub-id-type="pmid">18798982</pub-id>
</mixed-citation><ref id="pcbi.1003711-Gerstein1"><label>17</label>
<mixed-citation publication-type="journal"><name><surname>Gerstein</surname>
<given-names>MB</given-names>
</name><name><surname>Kundaje</surname>
<given-names>A</given-names>
</name><name><surname>Hariharan</surname>
<given-names>M</given-names>
</name><name><surname>Landt</surname>
<given-names>SG</given-names>
</name><name><surname>Yan</surname>
<given-names>K-K</given-names>
</name><etal>et al</etal>
(<year>2012</year>
) <article-title>Architecture of the human regulatory network derived from ENCODE data</article-title>
. <source>Nature</source>
<volume>489</volume>
: <fpage>91</fpage>
–<lpage>100</lpage>
<pub-id pub-id-type="doi">10.1038/nature11245</pub-id>
<pub-id pub-id-type="pmid">22955619</pub-id>
</mixed-citation><ref id="pcbi.1003711-Wang1"><label>18</label>
<mixed-citation publication-type="journal"><name><surname>Wang</surname>
<given-names>J</given-names>
</name><name><surname>Zhuang</surname>
<given-names>J</given-names>
</name><name><surname>Iyer</surname>
<given-names>S</given-names>
</name><name><surname>Lin</surname>
<given-names>X</given-names>
</name><name><surname>Whitfield</surname>
<given-names>TW</given-names>
</name><etal>et al</etal>
(<year>2012</year>
) <article-title>Sequence features and chromatin structure around the genomic regions bound by 119 human transcription factors</article-title>
. <source>Genome Res</source>
<volume>22</volume>
: <fpage>1798</fpage>
–<lpage>1812</lpage>
<pub-id pub-id-type="doi">10.1101/gr.139105.112</pub-id>
<pub-id pub-id-type="pmid">22955990</pub-id>
</mixed-citation><ref id="pcbi.1003711-Parelho1"><label>19</label>
<mixed-citation publication-type="journal"><name><surname>Parelho</surname>
<given-names>V</given-names>
</name><name><surname>Hadjur</surname>
<given-names>S</given-names>
</name><name><surname>Spivakov</surname>
<given-names>M</given-names>
</name><name><surname>Leleu</surname>
<given-names>M</given-names>
</name><name><surname>Sauer</surname>
<given-names>S</given-names>
</name><etal>et al</etal>
(<year>2008</year>
) <article-title>Cohesins Functionally Associate with CTCF on Mammalian Chromosome Arms</article-title>
. <source>Cell</source>
<volume>132</volume>
: <fpage>422</fpage>
–<lpage>433</lpage>
<pub-id pub-id-type="doi">10.1016/j.cell.2008.01.011</pub-id>
<pub-id pub-id-type="pmid">18237772</pub-id>
</mixed-citation><ref id="pcbi.1003711-Gorkin1"><label>20</label>
<mixed-citation publication-type="journal"><name><surname>Gorkin</surname>
<given-names>DU</given-names>
</name><name><surname>Lee</surname>
<given-names>D</given-names>
</name><name><surname>Reed</surname>
<given-names>X</given-names>
</name><name><surname>Fletez-Brant</surname>
<given-names>C</given-names>
</name><name><surname>Bessling</surname>
<given-names>SL</given-names>
</name><etal>et al</etal>
(<year>2012</year>
) <article-title>Integration of ChIP-seq and machine learning reveals enhancers and a predictive regulatory sequence vocabulary in melanocytes</article-title>
. <source>Genome Res</source>
<volume>22</volume>
: <fpage>2290</fpage>
–<lpage>2301</lpage>
<pub-id pub-id-type="doi">10.1101/gr.139360.112</pub-id>
<pub-id pub-id-type="pmid">23019145</pub-id>
</mixed-citation><ref id="pcbi.1003711-Machanick1"><label>21</label>
<mixed-citation publication-type="journal"><name><surname>Machanick</surname>
<given-names>P</given-names>
</name><name><surname>Bailey</surname>
<given-names>TL</given-names>
</name><year>2011</year>
) <article-title>MEME-ChIP: motif analysis of large DNA datasets</article-title>
. <source>Bioinformatics</source>
<volume>27</volume>
: <fpage>1696</fpage>
–<lpage>1697</lpage>
<pub-id pub-id-type="doi">10.1093/bioinformatics/btr189</pub-id>
<pub-id pub-id-type="pmid">21486936</pub-id>
</mixed-citation><ref id="pcbi.1003711-Pevny1"><label>22</label>
<mixed-citation publication-type="journal"><name><surname>Pevny</surname>
<given-names>L</given-names>
</name><name><surname>Simon</surname>
<given-names>MC</given-names>
</name><name><surname>Robertson</surname>
<given-names>E</given-names>
</name><name><surname>Klein</surname>
<given-names>WH</given-names>
</name><name><surname>Tsai</surname>
<given-names>S-F</given-names>
</name><etal>et al</etal>
(<year>1991</year>
) <article-title>Erythroid differentiation in chimaeric mice blocked by a targeted mutation in the gene for transcription factor GATA-1</article-title>
. <source>Nature</source>
<volume>349</volume>
: <fpage>257</fpage>
–<lpage>260</lpage>
<pub-id pub-id-type="doi">10.1038/349257a0</pub-id>
<pub-id pub-id-type="pmid">1987478</pub-id>
</mixed-citation><ref id="pcbi.1003711-Lelivre1"><label>23</label>
<mixed-citation publication-type="journal"><name><surname>Lelièvre</surname>
<given-names>E</given-names>
</name><name><surname>Lionneton</surname>
<given-names>F</given-names>
</name><name><surname>Soncin</surname>
<given-names>F</given-names>
</name><name><surname>Vandenbunder</surname>
<given-names>B</given-names>
</name><year>2001</year>
) <article-title>The Ets family contains transcriptional activators and repressors involved in angiogenesis</article-title>
. <source>Int J Biochem Cell Biol</source>
<volume>33</volume>
: <fpage>391</fpage>
–<lpage>407</lpage>
<pub-id pub-id-type="doi">10.1016/S1357-2725(01)00025-5</pub-id>
<pub-id pub-id-type="pmid">11312108</pub-id>
</mixed-citation><ref id="pcbi.1003711-Leslie3"><label>24</label>
<mixed-citation publication-type="journal"><name><surname>Leslie</surname>
<given-names>C</given-names>
</name><name><surname>Eskin</surname>
<given-names>E</given-names>
</name><name><surname>Noble</surname>
<given-names>WS</given-names>
</name><year>2002</year>
) <article-title>The spectrum kernel: a string kernel for SVM protein classification</article-title>
. <source>Pac Symp Biocomput</source>
<volume>2002</volume>
: <fpage>564</fpage>
–<lpage>575</lpage>
.</mixed-citation><ref id="pcbi.1003711-Berger1"><label>25</label>
<mixed-citation publication-type="journal"><name><surname>Berger</surname>
<given-names>MF</given-names>
</name><name><surname>Philippakis</surname>
<given-names>AA</given-names>
</name><name><surname>Qureshi</surname>
<given-names>AM</given-names>
</name><name><surname>He</surname>
<given-names>FS</given-names>
</name><name><surname>Estep</surname>
<given-names>PW</given-names>
</name><etal>et al</etal>
(<year>2006</year>
) <article-title>Compact, universal DNA microarrays to comprehensively determine transcription-factor binding site specificities</article-title>
. <source>Nat Biotechnol</source>
<volume>24</volume>
: <fpage>1429</fpage>
–<lpage>1435</lpage>
<pub-id pub-id-type="doi">10.1038/nbt1246</pub-id>
<pub-id pub-id-type="pmid">16998473</pub-id>
</mixed-citation><ref id="pcbi.1003711-Weirauch1"><label>26</label>
<mixed-citation publication-type="journal"><name><surname>Weirauch</surname>
<given-names>MT</given-names>
</name><name><surname>Cote</surname>
<given-names>A</given-names>
</name><name><surname>Norel</surname>
<given-names>R</given-names>
</name><name><surname>Annala</surname>
<given-names>M</given-names>
</name><name><surname>Zhao</surname>
<given-names>Y</given-names>
</name><etal>et al</etal>
(<year>2013</year>
) <article-title>Evaluation of methods for modeling transcription factor sequence specificity</article-title>
. <source>Nat Biotechnol</source>
<volume>31</volume>
: <fpage>126</fpage>
–<lpage>134</lpage>
<pub-id pub-id-type="doi">10.1038/nbt.2486</pub-id>
<pub-id pub-id-type="pmid">23354101</pub-id>
</mixed-citation><ref id="pcbi.1003711-Lee2"><label>27</label>
<mixed-citation publication-type="book">Lee D, Beer MA (2014) Mammalian Enhancer Prediction. Genome Analysis: Current Procedures and Applications. Caister Academic Press.</mixed-citation>
</ref><ref id="pcbi.1003711-Sandberg1"><label>28</label>
<mixed-citation publication-type="journal"><name><surname>Sandberg</surname>
<given-names>R</given-names>
</name><name><surname>Winberg</surname>
<given-names>G</given-names>
</name><name><surname>Bränden</surname>
<given-names>C-I</given-names>
</name><name><surname>Kaske</surname>
<given-names>A</given-names>
</name><name><surname>Ernberg</surname>
<given-names>I</given-names>
</name><etal>et al</etal>
(<year>2001</year>
) <article-title>Capturing Whole-Genome Characteristics in Short Sequences Using a Naïve Bayesian Classifier</article-title>
. <source>Genome Res</source>
<volume>11</volume>
: <fpage>1404</fpage>
–<lpage>1409</lpage>
<pub-id pub-id-type="doi">10.1101/gr.186401</pub-id>
<pub-id pub-id-type="pmid">11483581</pub-id>
</mixed-citation><ref id="pcbi.1003711-Higham1"><label>29</label>
<mixed-citation publication-type="book">Higham NJ (2002) Accuracy and Stability of Numerical Algorithms. 2nd ed. Philadelphia, PA, USA: Society for Industrial and Applied Mathematics.</mixed-citation>
</ref><ref id="pcbi.1003711-Boser1"><label>30</label>
<mixed-citation publication-type="other">Boser BE, Guyon IM, Vapnik VN (1992) A training algorithm for optimal margin classifiers. Proceedings of the fifth annual workshop on Computational learning theory. COLT '92. New York, NY: ACM. pp. 144–152. doi:10.1145/130385.130401.</mixed-citation>
</ref><ref id="pcbi.1003711-Vapnik1"><label>31</label>
<mixed-citation publication-type="book">Vapnik VN (1995) The Nature of Statistical Learning Theory. New York, NY: Springer.</mixed-citation>
</ref><ref id="pcbi.1003711-Sonnenburg1"><label>32</label>
<mixed-citation publication-type="journal"><name><surname>Sonnenburg</surname>
<given-names>S</given-names>
</name><name><surname>Rätsch</surname>
<given-names>G</given-names>
</name><name><surname>Henschel</surname>
<given-names>S</given-names>
</name><name><surname>Widmer</surname>
<given-names>C</given-names>
</name><name><surname>Behr</surname>
<given-names>J</given-names>
</name><etal>et al</etal>
(<year>2010</year>
) <article-title>The SHOGUN Machine Learning Toolbox</article-title>
. <source>J Mach Learn Res</source>
<volume>11</volume>
: <fpage>1799</fpage>
–<lpage>1802</lpage>
.</mixed-citation><ref id="pcbi.1003711-Joachims1"><label>33</label>
<mixed-citation publication-type="book">Joachims T (1999) Making large-scale support vector machine learning practical. Advances in Kernel Methods. Cambridge, MA: MIT Press. pp. 169–184.</mixed-citation>
</ref><ref id="pcbi.1003711-Jaakkola1"><label>34</label>
<mixed-citation publication-type="journal"><name><surname>Jaakkola</surname>
<given-names>T</given-names>
</name><name><surname>Diekhans</surname>
<given-names>M</given-names>
</name><name><surname>Haussler</surname>
<given-names>D</given-names>
</name><year>2000</year>
) <article-title>A Discriminative Framework for Detecting Remote Protein Homologies</article-title>
. <source>J Comput Biol</source>
<volume>7</volume>
: <fpage>95</fpage>
–<lpage>114</lpage>
<pub-id pub-id-type="doi">10.1089/10665270050081405</pub-id>
<pub-id pub-id-type="pmid">10890390</pub-id>
</mixed-citation><ref id="pcbi.1003711-Cormen1"><label>35</label>
<mixed-citation publication-type="book">Cormen TH, Leiserson CE, Rivest RL, Stein C (2001) Introduction to algorithms. Cambridge, MA: MIT press.</mixed-citation>
</ref><ref id="pcbi.1003711-Gupta1"><label>36</label>
<mixed-citation publication-type="journal"><name><surname>Gupta</surname>
<given-names>S</given-names>
</name><name><surname>Stamatoyannopoulos</surname>
<given-names>JA</given-names>
</name><name><surname>Bailey</surname>
<given-names>TL</given-names>
</name><name><surname>Noble</surname>
<given-names>WS</given-names>
</name><year>2007</year>
) <article-title>Quantifying similarity between motifs</article-title>
. <source>Genome Biol</source>
<volume>8</volume>
: <fpage>R24</fpage>
<pub-id pub-id-type="doi">10.1186/gb-2007-8-2-r24</pub-id>
<pub-id pub-id-type="pmid">17324271</pub-id>
</mixed-citation><ref id="pcbi.1003711-Sing1"><label>37</label>
<mixed-citation publication-type="journal"><name><surname>Sing</surname>
<given-names>T</given-names>
</name><name><surname>Sander</surname>
<given-names>O</given-names>
</name><name><surname>Beerenwinkel</surname>
<given-names>N</given-names>
</name><name><surname>Lengauer</surname>
<given-names>T</given-names>
</name><year>2005</year>
) <article-title>ROCR: visualizing classifier performance in R</article-title>
. <source>Bioinformatics</source>
<volume>21</volume>
: <fpage>3940</fpage>
–<lpage>3941</lpage>
<pub-id pub-id-type="doi">10.1093/bioinformatics/bti623</pub-id>
<pub-id pub-id-type="pmid">16096348</pub-id>
</mixed-citation><ref id="pcbi.1003711-Bryne1"><label>38</label>
<mixed-citation publication-type="journal"><name><surname>Bryne</surname>
<given-names>JC</given-names>
</name><name><surname>Valen</surname>
<given-names>E</given-names>
</name><name><surname>Tang</surname>
<given-names>M-HE</given-names>
</name><name><surname>Marstrand</surname>
<given-names>T</given-names>
</name><name><surname>Winther</surname>
<given-names>O</given-names>
</name><etal>et al</etal>
(<year>2008</year>
) <article-title>JASPAR, the open access database of transcription factor-binding profiles: new content and tools in the 2008 update</article-title>
. <source>Nucleic Acids Res</source>
<volume>36</volume>
: <fpage>D102</fpage>
–<lpage>106</lpage>
<pub-id pub-id-type="doi">10.1093/nar/gkm955</pub-id>
<pub-id pub-id-type="pmid">18006571</pub-id>
</mixed-citation>
