Skip to content

Introduce consume_seqfile_banding method - #1571

Closed
standage wants to merge 35 commits into
masterfrom
feature/consume_bitsplit
Closed

standage wants to merge 35 commits into
masterfrom
feature/consume_bitsplit

Conversation

@standage

@standage standage commented Jan 3, 2017 •

Copy link
Copy Markdown
Member

WIP that, when finished, closes #1570. Implements a variant of consume_fasta (and all of the functions it calls) that consumes a subset of the observed k-mers. Given num_batches = 2^n and batch \in range(num_batches), the function only stores k-mers whose n least significant bits encode batch.

Remaining tasks:


  • Is it mergeable?
  • make test Did it pass the tests?
  • make clean diff-cover If it introduces new functionality in
    scripts/ is it tested?
  • make format diff_pylint_report cppcheck doc pydocstyle Is it well
    formatted?
  • Did it change the command-line interface? Only backwards-compatible
    additions are allowed without a major version increment. Changing file
    formats also requires a major version number increment.
  • For substantial changes or changes to the command-line interface, is it
    documented in CHANGELOG.md? See keepachangelog
    for more details.
  • Was a spellchecker run on the source code and documentation after
    changes were made?
  • Do the changes respect streaming IO? (Are they
    tested for streaming IO?)

@standage

standage commented Jan 3, 2017

Copy link
Copy Markdown
Member Author

This PR more or less duplicates the contents of several functions. If it turns out to be generally useful and something we want to keep long term, we'll need to consider how to refactor the code to minimize duplication and maintenance burden. Note, we've already acknowledged the need for some cleanup in #1378.

@codecov-io

codecov-io commented Jan 3, 2017 •

Copy link
Copy Markdown

Codecov Report

Merging #1571 into master will decrease coverage by 0.02%.
The diff coverage is 57.4%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master    #1571      +/-   ##
==========================================
- Coverage   64.16%   64.13%   -0.03%     
==========================================
  Files          79       79              
  Lines       12189    12244      +55     
  Branches     3072     3096      +24     
==========================================
+ Hits         7821     7853      +32     
- Misses       2670     2674       +4     
- Partials     1698     1717      +19
Impacted Files Coverage Δ
lib/hashtable.hh 68.42% <ø> (ø) ⬆️
lib/hashtable.cc 57.07% <56.75%> (-0.07%) ⬇️
khmer/_khmer.cc 57.49% <58.82%> (+0.01%) ⬆️
lib/hashgraph.cc 47.56% <0%> (+0.12%) ⬆️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update d40a81a...0dbe256. Read the comment docs.

@standage

standage commented Jan 3, 2017 •

Copy link
Copy Markdown
Member Author

Ready for initial review @ctb @luizirber @camillescott @betatim. This may be easier to merge before #1511. If so I'll clean it up for final review.

@standage standage mentioned this pull request Jan 4, 2017
12 of 13 tasks
Comment thread lib/hashtable.cc Outdated
delete parser;
}

void Hashtable::consume_fasta_bitsplit(

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

bitsplit doesn't mean anything to me :-/ Coudl we call it consume_fasta_sampling or some such?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I agree the name is wonky. consume_fasta_sampling is better but honestly I've thought of half a dozen names and I'm not in love with any of them.

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

my suggestion for the name is "banding".

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah, I quite like banding. And I think num_batches and batch are good replacements for num_bins and tag.

Comment thread lib/hashtable.cc Outdated
{
bool powerof2 = !(num_bins == 0) && !(num_bins & (num_bins - 1));
if (!powerof2 || tag >= num_bins) {
return;

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

why not raise a ValueError?

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed.

@betatim

betatim commented Jan 4, 2017 •

Copy link
Copy Markdown
Member

As mentioned here #1570 (comment) this needs a good hashfunction. Should make it either impossible to call with 2bit encoding or check at run time that the user knows what they are doing.

tables use murmurhash

@standage

standage commented Jan 4, 2017

Copy link
Copy Markdown
Member Author

tables use murmurhash

And to be clear, the concern with the hash function isn't correctness per se, but getting a roughly even split. Unless the sequence is truly random (simulated and biologically uninteresting), some k-mer suffixes are going to be overrepresented and skew the split. At least for my purposes, this wouldn't be a problem for correctness as much as for performance. I want to be able to split the work into evenly sized bins.

@standage

standage commented Jan 4, 2017

Copy link
Copy Markdown
Member Author

I did some basic profiling. For very small memory sizes, there is little difference in runtime between consume_fasta and consume_fasta_bitsplit on ≈1/8 of the k-mers. But as memory increases, the runtime of consume_fasta grows linearly-ish with memory while the runtime of consume_fasta_bitsplit increases only marginally. What's great is that consume_fasta_bitsplit can operate at a reasonable FPR using less memory, presumably taking less of a performance hit from NUMA.

@standage

standage commented Jan 4, 2017

Copy link
Copy Markdown
Member Author

But regardless, increasing memory only has a marginal effect on runtime. Perhaps it becomes more I/O bound?

@betatim

betatim commented Jan 5, 2017

Copy link
Copy Markdown
Member

When you say increasing memory you mean making the hashtable bigger?

So for small hashtable sizes we spend most of our time doing I/O, hence throwing away kmers after reading them in doesn't speed things up. For large hashtables we spend most of our effort on entering stuff into the hashtable, so reading in kmers and then chucking them makes things faster?! Interesting.

@ctb

ctb commented Jan 5, 2017

Copy link
Copy Markdown
Member

excellent. @betatim I've lost track (vacation brain) of your results on cache misses; this is probably not due to cache missing, right?

@ctb

ctb commented Jan 5, 2017

Copy link
Copy Markdown
Member

@standage, you should be able to separate out the IO effects from the memory size effects by an affine fit, yes?

@ctb

ctb commented Jan 5, 2017

Copy link
Copy Markdown
Member

What's old is new: #70

@ctb

ctb commented Jan 5, 2017

Copy link
Copy Markdown
Member

Commit 786a3be

Comment thread tests/test_banding.py
(21, 5e7, 2, 16),
])
def test_banding_in_memory(ksize, memory, epsilon, numbands):
infile = utils.get_test_data('test-reads.fa')

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure this function meets @ctb's simplicity criteria for automated tests. I'm not sure there's a simpler, clearer way to iterate over the data.

Comment thread tests/test_banding.py
def test_banding_to_disk(ksize, memory, numbands):
infile = utils.get_test_data('test-reads.fa')
path1 = utils.get_temp_filename('normal.ct')
path2 = utils.get_temp_filename('banding.ct')

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This function is simple enough, and it provides a stronger guarantee. Perhaps it's all we need?

@ctb

ctb commented Feb 10, 2017 via email

Copy link
Copy Markdown
Member

@standage

Copy link
Copy Markdown
Member Author

Procedural complexity == branching?

@standage

Copy link
Copy Markdown
Member Author

My recent tests have put my concerns regarding the hash functions and banding to rest. Now that #1511 is merged, this is now ready for review and merge. @ctb @betatim @camillescott @luizirber

Comment thread lib/hashtable.cc
continue;
}

bool is_valid = true;

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can I persuade you to use read.set_clean_seq(); plus read.cleaned_seq != read.sequence instead of doing your own processing/normalising?

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe what we need is a method that does the toupper() for us to avoid having "loop-over-seq-and-uppercase-then-check" here.

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Two concerns / observations.

  • Comparing cleaned seq vs original seq will disqualify reads with lowercase ACGT, which I think is a lousy default with unexpected consequences.
  • The current consume_fasta function in Hashtable doesn't use this utility yet either.

As @ctb mentioned, there's some work to be done with removing redundant read cleaning code (see #1590), and that's not all going to be resolved in this PR. Ultimately this will be handled by the parser instead of the consume function, right?

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah yes, humm, so we definitely need a method that can do the comparison that is centralised. And with cython'isation things will move again. Punt?

Comment thread tests/test_tabletype.py
code...
"""

from __future__ import absolute_import

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Was this missing before? Puzzled how travis didn't find it.

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah, weird.

byallaccountsitdoesntmakesense

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Now we need it to make the tests pass or do they work without as well?

@betatim

betatim commented Feb 16, 2017

Copy link
Copy Markdown
Member

Do you want to add one more test that computes the hash of two kmers, one that should be added and one that shouldn't, and explicitly checks that they ended up in the CT? A bit like https://github.com/dib-lab/khmer/pull/1571/files#diff-7352ff13b9e5cdbeacb24c9e7bcb49b1R1136 except here I am guessing that this is what you wanted because I don't know the hash values off the top of my head.

@ctb

ctb commented Feb 16, 2017 via email

Copy link
Copy Markdown
Member

@standage

Copy link
Copy Markdown
Member Author

@betatim Yeah, that's exactly what I was going for with that test. Are you saying to test the hash value instead of the k-mer itself? Please clarify.

@standage standage changed the title Introduce consume_fasta_banding method Introduce consume_seqfile_banding method Feb 16, 2017
@betatim

betatim commented Feb 17, 2017

Copy link
Copy Markdown
Member

I was going for a test like this:

ct = Counttable(...)
h = ct.hash(some_kmer)
assert h % 8 == 3
ct.consume_seqfile_banding(file_that_contains_some_kmer, 8, 3)
assert ct.get(some_kmer) == 1

we could even add the assert for the hash value to the current test. The motivation is to leave an explicit clue as to why you expect a certain kmer to be present and which ones not.

@betatim betatim left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM.

@ctb

ctb commented Feb 25, 2017 via email

Copy link
Copy Markdown
Member

@standage

standage commented Jun 2, 2017

Copy link
Copy Markdown
Member Author

Closing in favor of #1720.

@standage standage closed this Jun 2, 2017
@standage
standage deleted the feature/consume_bitsplit branch June 3, 2017 04:05
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Implement a new 'consume' method that only stores a fraction of the observed k-mers

4 participants