Closed
Bug 1206357
Opened 6 years ago
Closed 6 years ago
MFBT needs a fast way to conduct Bernoulli trials
Categories
(Core :: MFBT, defect)
Core
MFBT
Tracking
()
RESOLVED
FIXED
mozilla44
Tracking  Status  

firefox44    fixed 
People
(Reporter: jimb, Assigned: jimb)
References
Details
Attachments
(1 file, 2 obsolete files)
The comments in the patch explain everything. At length. Perhaps too much length.
Attachment #8663252 
Flags: review?(jwalden+bmo)
Comment 1•6 years ago


Comment on attachment 8663252 [details] [diff] [review] Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. Review of attachment 8663252 [details] [diff] [review]:  ::: mfbt/FastBernoulliTrial.h @@ +165,5 @@ > + > + public: > + // Construct a fast Bernoulli trial generator. Calls to trial() return true > + // with probability aProbability. Use aState0 and aState1 to seed the > + // random number generator; both may not be zero. MOZ_ASSERT this?
Comment 2•6 years ago


(In reply to Nick Fitzgerald [:fitzgen][:nf] from comment #1) > Comment on attachment 8663252 [details] [diff] [review] > Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. > > Review of attachment 8663252 [details] [diff] [review]: >  > > ::: mfbt/FastBernoulliTrial.h > @@ +165,5 @@ > > + > > + public: > > + // Construct a fast Bernoulli trial generator. Calls to trial() return true > > + // with probability aProbability. Use aState0 and aState1 to seed the > > + // random number generator; both may not be zero. > > MOZ_ASSERT this? Nevermind, I see it is asserted in the patch in bug 1206356.
Assignee  
Comment 3•6 years ago


Comment on attachment 8663252 [details] [diff] [review] Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. Review of attachment 8663252 [details] [diff] [review]:  ::: mfbt/FastBernoulliTrial.h @@ +165,5 @@ > + > + public: > + // Construct a fast Bernoulli trial generator. Calls to trial() return true > + // with probability aProbability. Use aState0 and aState1 to seed the > + // random number generator; both may not be zero. The constructor for XorShift128PlusRNG asserts this.
Assignee  
Updated•6 years ago

Attachment #8663252 
Flags: review?(jwalden+bmo) → review?(nfroyd)
Comment 4•6 years ago


Comment on attachment 8663252 [details] [diff] [review] Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. Review of attachment 8663252 [details] [diff] [review]:  Cool beans, neat idea, fun comments to read, mirabile lectu. (I'm gonna pretend that pseudoLatin is real :) Might well be!) The undercounting issue with trial(size_t), however, means this definitely isn't ready to go yet, unless I'm entirely misunderstanding what should happen for the Ntrial case. Please do correct me if I'm wrong. ::: mfbt/FastBernoulliTrial.h @@ +2,5 @@ > +/* vim: set ts=8 sts=2 et sw=2 tw=80: */ > +/* This Source Code Form is subject to the terms of the Mozilla Public > + * License, v. 2.0. If a copy of the MPL was not distributed with this > + * file, You can obtain one at http://mozilla.org/MPL/2.0/. */ > + Add an overview comment like /* An algorithm for sampling N trials with a desired probability. */ or something. @@ +15,5 @@ > + > +namespace mozilla { > + > +// class FastBernoulliTrial > +// Much comments! I agree it's a bit overkill, but (like you) I kind of don't want to take a knife to it anywhere, and it's nicely readable and all. So good enough asis. @@ +98,5 @@ > + // get "rolling the dice" afresh for every trial. > + // > + // In short, skip counts in Bernoulli trials of probability P obey a geometric > + // distribution. If a random variable X is uniformly distributed from [0..1), > + // then floor(log X base (1P)) has the appropriate geometric distribution for I would descend into LaTeX mode for this, I think, to be completely clear what's said. @@ +170,5 @@ > + FastBernoulliTrial(double aProbability, uint64_t aState0, uint64_t aState1) > + : mGenerator(aState0, aState1), > + mSkipCount(0) > + { > + setProbability(aProbability); I like how this also initializes the skipCount  clever. @@ +174,5 @@ > + setProbability(aProbability); > + } > + > + // Return true with probability mProbability. The lower mProbability is, > + // the faster this function runs. Use /** */ comments for all these userfacing documentation comments. @@ +186,5 @@ > + } > + > + // Return true with the same likelihood that at least one of the next aCount > + // calls to trial() would. This takes roughly as much time as a single call to > + // trial(). So, this is interesting. Won't this potentially *undersample* events, because it can't return true for a case where multiple trials succeed within the given count of trials? Like, say we have FastBernoulliTrial fbt(0.5, ..., ...); if (fbt.trial(42)) ...; // thing happens once for (size_t i = 0; i < 42; i++) { if (fbt.trial()) ...; // thing happens 21ish times } Seems to me like this particular method  while a nice idea  needs a different signature and a different algorithm, somehow. Maybe returning the number of times a trial succeeded? At that point, you probably could unify these into a single method: size_t trial(size_t aCount = 1) { ... } because the oncecase would be a degenerate form of the manycase. @@ +218,5 @@ > + // 1.0. Otherwise, this is not used. > + double mLogNotProbability; > + > + // Our random number generator. > + mozilla::XorShift128PlusRNG mGenerator; For what it's worth, this doesn't need mozilla:: as you're in the mozilla namespace. @@ +237,5 @@ > + > + // If the probabilility is zero, trial never returns true. Don't bother us > + // for a while. > + if (mProbability == 0.0) { > + mSkipCount = SIZE_MAX; lol ::: mfbt/tests/TestFastBernoulliTrial.cpp @@ +119,5 @@ > + > + { > + size_t count = 0; > + for (size_t i = 0; i < 10000; i++) > + count += bernoulli.trial(1000); Given the undercounting flub in trial(size_t), I think we probably should have a test that verifies that N calls to trial(1) has equivalent behavior to 1 call to trial(N) in terms of observations, somehow. @@ +133,5 @@ > + 0x74f298193fe1c5b1ULL); > + > + // Establish a very high skip count. > + bernoulli.setProbability(0.0); > + Worth a !trial() assert here?
Attachment #8663252 
Flags: review?(nfroyd) → feedback+
Assignee  
Comment 5•6 years ago


(In reply to Jeff Walden [:Waldo] (remove +bmo to email) from comment #4) > > + // In short, skip counts in Bernoulli trials of probability P obey a geometric > > + // distribution. If a random variable X is uniformly distributed from [0..1), > > + // then floor(log X base (1P)) has the appropriate geometric distribution for > > I would descend into LaTeX mode for this, I think, to be completely clear > what's said. Okay. You mean, like what Wikipedia has? \lfloor \ln(U) / \ln(1p)\rfloor As someone who hasn't used TeX in a very long time, that's not something I'd find helpful. I could write it out in real C++: std::floor(std::log(X) / std::log(1P)), or JS. Would that serve? > > + // Return true with the same likelihood that at least one of the next aCount > > + // calls to trial() would. This takes roughly as much time as a single call to > > + // trial(). > > So, this is interesting. Won't this potentially *undersample* events, > because it can't return true for a case where multiple trials succeed within > the given count of trials? Like, say we have > > FastBernoulliTrial fbt(0.5, ..., ...); > if (fbt.trial(42)) > ...; // thing happens once > for (size_t i = 0; i < 42; i++) { > if (fbt.trial()) > ...; // thing happens 21ish times > } The way you'd use it is actually more like this: FastBernoulliTrial fbt(0.5, ..., ...); if (fbt.trial(42)) note_event_with_size(42); for (size_t i = 0; i < 42; i++) { if (fbt.trial()) note_event_with_size(1); } and the events are going to be displayed visually in some way that presents the larger ones more prominently. This is the technique used by Google's tcmalloc library: https://github.com/gperftools/gperftools/blob/master/src/sampler.h#L73L100 The comments above the class for trial(n) are meant to be a restatement of the comments in the tcmalloc header. So I don't consider it a "flub"; it was designed rather carefully, and with respect for standard practices. > ::: mfbt/tests/TestFastBernoulliTrial.cpp > @@ +119,5 @@ > > + > > + { > > + size_t count = 0; > > + for (size_t i = 0; i < 10000; i++) > > + count += bernoulli.trial(1000); > > Given the undercounting flub in trial(size_t), I think we probably should > have a test that verifies that N calls to trial(1) has equivalent behavior > to 1 call to trial(N) in terms of observations, somehow. That's exactly what this test does. Those counts weren't just chosen to match the observed behavior; I computed the "expected" values mathematically, and checked that each count was plausibly close. I've added comments like these; does that help? { size_t count = 0; for (size_t i = 0; i < 10000; i++) count += bernoulli.trial(3); // Expected value: (1  (1  0.01) ** 3) == 0.0297, // 0.0297 * 10000 == 297 MOZ_RELEASE_ASSERT(count == 304); } { size_t count = 0; for (size_t i = 0; i < 10000; i++) count += bernoulli.trial(10); // Expected value: (1  (1  0.01) ** 10) == 0.0956, // 0.0956 * 10000 == 956 MOZ_RELEASE_ASSERT(count == 936); } { size_t count = 0; for (size_t i = 0; i < 10000; i++) count += bernoulli.trial(100); // Expected value: (1  (1  0.01) ** 100) == 0.6339 // 0.6339 * 10000 == 6339 MOZ_RELEASE_ASSERT(count == 6372); } > > + // Establish a very high skip count. > > + bernoulli.setProbability(0.0); > > + > > Worth a !trial() assert here? This one was a regression test; I need the state to be exactly as it is after a call to setProbability in order to reproduce the conditions of the original bug. Calling 'trial' changes the state.
Comment 6•6 years ago


(In reply to Jim Blandy :jimb from comment #5) > So I don't consider it a "flub"; it was designed rather carefully, and with > respect for standard practices. Easy, easy! As I said, > unless I'm entirely misunderstanding what should happen for the Ntrial case. > Please do correct me if I'm wrong. If your comment is accurate, then precisely this would be what had happened. Anyway. Looking more at your comment, and thinking more, I think I'm leaning toward concluding that I wasn't thinking hard enough, or was misunderstanding something, or something along those lines. But I need to think some more on this before I can convince myself that's the case. Please give me a little more time to convince myself one way or the other.
Assignee  
Comment 7•6 years ago


I'm sorry for being prickly. The other comments all sound good; I'll take care of them now.
Assignee  
Comment 8•6 years ago


(In reply to Jeff Walden [:Waldo] (remove +bmo to email) from comment #4) > @@ +170,5 @@ > > + FastBernoulliTrial(double aProbability, uint64_t aState0, uint64_t aState1) > > + : mGenerator(aState0, aState1), > > + mSkipCount(0) > > + { > > + setProbability(aProbability); > > I like how this also initializes the skipCount  clever. Actually, since you pointed it out, I think this became unnecessary after the fix for the bug that TestChangeProbability looks for. The constructor always calls setProbability, which always calls chooseSkipCount, which always sets mSkipCount. So I've taken this out.
Comment 9•6 years ago


Comment on attachment 8663252 [details] [diff] [review] Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. Review of attachment 8663252 [details] [diff] [review]:  ::: mfbt/FastBernoulliTrial.h @@ +204,5 @@ > + void setProbability(double aProbability) { > + MOZ_ASSERT(0.0 <= aProbability && aProbability <= 1.0); > + mProbability = aProbability; > + if (0.0 < mProbability && mProbability < 1.0) > + mLogNotProbability = std::log(1.0  mProbability); In IEEE floating point semantics, x < 1.0 does indeed imply that 1.0  x != 0.0. Getting this right in all cases depends on denormal numbers. However, not all hardware supports denormals, and sometimes they are turned off for performance. (Look up the FTZ and DAZ flags for Intel CPUs). There could be values of aProbability that would cause mLogNotProbability to become infinity. It looks like everything still works as expected in that case. @@ +215,5 @@ > + double mProbability; > + > + // The value of std::log(1  mProbability), if mProbability is neither 0.0 nor > + // 1.0. Otherwise, this is not used. > + double mLogNotProbability; Since this value is only used as a divisor, it may be faster to store the reciprocal value: 1/(1P). Floating point multiplication is faster than division. @@ +242,5 @@ > + return false; > + } > + > + mSkipCount = std::floor(std::log(mGenerator.nextDouble()) > + / mLogNotProbability); What happens if nextDouble() returns exactly 0.0? log(0.0) returns infinity. Similarly, there could be rare cases where floor() returns something finite, but larger than SIZE_MAX. You could clamp the value before converting it size_t: double x = std::floor(std::log(nextDouble()) * mRecipLogNotProbability); mSkipCount = x < SIZE_MAX ? x : SIZE_MAX;
Assignee  
Comment 10•6 years ago


(In reply to Jakob Stoklund Olesen [:jolesen] from comment #9) > Comment on attachment 8663252 [details] [diff] [review] > Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. > > Review of attachment 8663252 [details] [diff] [review]: >  > > ::: mfbt/FastBernoulliTrial.h > @@ +204,5 @@ > > + void setProbability(double aProbability) { > > + MOZ_ASSERT(0.0 <= aProbability && aProbability <= 1.0); > > + mProbability = aProbability; > > + if (0.0 < mProbability && mProbability < 1.0) > > + mLogNotProbability = std::log(1.0  mProbability); > > In IEEE floating point semantics, x < 1.0 does indeed imply that 1.0  x != > 0.0. Getting this right in all cases depends on denormal numbers. > > However, not all hardware supports denormals, and sometimes they are turned > off for performance. (Look up the FTZ and DAZ flags for Intel CPUs). > > There could be values of aProbability that would cause mLogNotProbability to > become infinity. It looks like everything still works as expected in that > case. Based on our inperson conversations I've revised some of this code and added comments (will post after this reply), but I don't think mLogNotProbability can be infinity, and I don't think avoiding infinity requires denormals. Assume all numbers are exact unless marked otherwise with "~".  In the 'then' branch, mProbability is (0,1). The gaps below 1 are 2**53, so that range is really (0, 12**53]. Since the gaps near zero are much smaller, I don't think it has any problem with this.  Gaps near 1 are larger than gaps near 0, so for small ε > 0, 1ε == 1. If denormals are not permitted, this is still true. So 1mProbability is in [2**53, 1].  log(2**53) is ~37. So log(1mProbability) is in [37, 0]. In chooseSkipCount, we would like to assume that mLogNotProbability is negative, so we need to catch zero in setProbability. C++11 has a 'std::nextafter' function to advance doubles by the smallest possible increment in one direction or another, which makes it possible to actually test these cases, which is nice. > @@ +215,5 @@ > > + double mProbability; > > + > > + // The value of std::log(1  mProbability), if mProbability is neither 0.0 nor > > + // 1.0. Otherwise, this is not used. > > + double mLogNotProbability; > > Since this value is only used as a divisor, it may be faster to store the > reciprocal value: 1/(1P). Floating point multiplication is faster than > division. Interesting. I'll give it a shot. > @@ +242,5 @@ > > + return false; > > + } > > + > > + mSkipCount = std::floor(std::log(mGenerator.nextDouble()) > > + / mLogNotProbability); > > What happens if nextDouble() returns exactly 0.0? log(0.0) returns > infinity. Similarly, there could be rare cases where floor() returns > something finite, but larger than SIZE_MAX. Yes, this can occur. We need to clamp the value. I've added code for this, and addressed the SIZE_MAX > double conversion issues we discussed. > You could clamp the value before converting it size_t: > > double x = std::floor(std::log(nextDouble()) * mRecipLogNotProbability); > mSkipCount = x < SIZE_MAX ? x : SIZE_MAX; I think this won't work: the 'then' and 'else' expressions of a ?: operator are promoted to the 'wider' type, which would be 'double' in this case; if SIZE_MAX is 2**641, it will be rounded to 2**64, and then converted to size_t, which is undefined behavior. The comments in the code explain how I dealt with this.
Comment 11•6 years ago


jimb: I don't suppose you want to modify DMD to use this new class to control its sampling...? :)
Assignee  
Comment 12•6 years ago


Here's the revised patch. I've addressed Waldo's review comments and the numeric issues Jakub raised, and added a bunch more comments to cover the latter; the fixes are reasonable but not obvious.
Attachment #8663252 
Attachment is obsolete: true
Attachment #8666499 
Flags: review?(jwalden+bmo)
Assignee  
Comment 13•6 years ago


(In reply to Nicholas Nethercote [:njn] from comment #11) > jimb: I don't suppose you want to modify DMD to use this new class to > control its sampling...? :) Not in this quarter! :D But yeah, that's hopefully the kind of thing this would serve well.
Assignee  
Comment 14•6 years ago


Try push, including cppunit tests, which includes TestFastBernoulliTrial: https://treeherder.mozilla.org/#/jobs?repo=try&revision=75789f10fd1d
Comment 15•6 years ago


Comment on attachment 8666499 [details] [diff] [review] Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. Review of attachment 8666499 [details] [diff] [review]:  (In reply to Jim Blandy :jimb from comment #5) > As someone who hasn't used TeX in a very long time, that's not something I'd > find helpful. I could write it out in real C++: std::floor(std::log(X) / > std::log(1P)), or JS. Would that serve? That'd work, yes. Mostly I think it's that the wording to write out a log of nonstandard base is not the easiest wording to parse, IMO. If LaTeX doesn't work for you, I'm happy avoiding it. C++ that parses unambiguously works, too. > The way you'd use it is actually more like this: > > FastBernoulliTrial fbt(0.5, ..., ...); > if (fbt.trial(42)) > note_event_with_size(42); > for (size_t i = 0; i < 42; i++) { > if (fbt.trial()) > note_event_with_size(1); > } > > and the events are going to be displayed visually in some way that presents > the larger ones more prominently. This is the technique used by Google's > tcmalloc library: I see, so it's the responsibility of the "noter" to weight accordingly. And because P is constant for each event, and we're assuming many events, any discounting for the number of trials that succeed in a multitrial doesn't require knowledge of the decremented number in the FBT  if desired, you could just use the average as a reasonable enough guess that'd even out over time. My mind still isn't grokking this scheme, only mostly understanding it, but I guess it's close enough. > > ::: mfbt/tests/TestFastBernoulliTrial.cpp > > @@ +119,5 @@ > > > + > > > + { > > > + size_t count = 0; > > > + for (size_t i = 0; i < 10000; i++) > > > + count += bernoulli.trial(1000); > > > > Given the undercounting flub in trial(size_t), I think we probably should > > have a test that verifies that N calls to trial(1) has equivalent behavior > > to 1 call to trial(N) in terms of observations, somehow. > > That's exactly what this test does. Those counts weren't just chosen to > match the observed behavior; I computed the "expected" values > mathematically, and checked that each count was plausibly close. I've added > comments like these; does that help? The changes are fine, but mostly I think it's just following up on the responsibilityofthenoter bit above, not an independently different issue needing much addressing. > > Worth a !trial() assert here? > > This one was a regression test; I need the state to be exactly as it is > after a call to setProbability in order to reproduce the conditions of the > original bug. Calling 'trial' changes the state. Fine by me. ::: mfbt/FastBernoulliTrial.h @@ +14,5 @@ > +#include <stdint.h> > + > +namespace mozilla { > + > +/* Supernitpick, but the idea is you use /** */ comments  two stars after the opening  for comments that are user docs, so this should have two stars. The implementation comments  everything but the comment on FastBernoulliTrial(double aProbability, uint64_t aState0, uint64_t aState1), bool trial(), bool trial(size_t aCount), void setRandomState(uint64_t aState0, uint64_t aState1), and void setProbability(double aProbability)  can be whatever. So just one tweak to make now, but the idea's worth belaboring slightly for next time. ::: mfbt/tests/TestFastBernoulliTrial.cpp @@ +6,5 @@ > + > +#include "mozilla/Assertions.h" > +#include "mozilla/FastBernoulliTrial.h" > + > +#include <math.h> If you're using the std:: version of nextafter, use <cmath> here. Technically FBT.h implies this, but IWYU and all.
Attachment #8666499 
Flags: review?(jwalden+bmo) → review+
Assignee  
Comment 16•6 years ago


https://treeherder.mozilla.org/#/jobs?repo=try&revision=6587556569f5
Assignee  
Comment 17•6 years ago


Thanks for the review, Waldo. I think this patch implements the changes you requested. I ran into difficulties earlier misusing cmath, so I've done a buildonly try push: https://treeherder.mozilla.org/#/jobs?repo=try&revision=4d421fde0f42
Attachment #8666499 
Attachment is obsolete: true
Assignee  
Comment 18•6 years ago


https://treeherder.mozilla.org/#/jobs?repo=try&revision=831166404548
Assignee  
Comment 19•6 years ago


(In reply to Jim Blandy :jimb from comment #17) > I ran into difficulties earlier misusing cmath, so I've done a buildonly > try push: > > https://treeherder.mozilla.org/#/jobs?repo=try&revision=4d421fde0f42 Using <cmath> and std::nextafter is less portable than <math.h> and unqualified nextafter.
Assignee  
Comment 20•6 years ago


Okay, that push looks clean.
Assignee  
Comment 21•6 years ago


https://hg.mozilla.org/integration/mozillainbound/rev/0140d126b76b4fe692b22d880f897496370ef480 Bug 1206357: Add mfbt/FastBernoulliTrial.h, implementing efficient random sampling. r=waldo
Assignee  
Updated•6 years ago

Flags: intestsuite+
Target Milestone:  → mozilla44
You need to log in
before you can comment on or make changes to this bug.
Description
•