@Override public boolean apply(@Nullable Double input) { assert input != null; return rng.nextDouble() <= lambda.apply(input) / lambdaMax; } }
/** * Creates a non-homogenous Poisson process of the specified length. The * intensity is specified by the {@link IntensityFunction}. The non-homogenous * Poisson process is implemented using the thinning method as described in * [1]. * <p> * <b>References</b> * <ol> * <li>Lewis, P.A.W. and Shedler, G.S. Simulation of nonhomogenous Poisson * processes by thinning. Naval Research Logistic Quarterly 26, (1979), * 403–414.</li> * </ol> * @param length The length of Poisson process, all generated times will be in * the interval [0,length). * @param function The intensity function. * @return A newly constructed non-homogenous Poisson process * {@link TimeSeriesGenerator}. */ public static TimeSeriesGenerator nonHomogenousPoisson(double length, IntensityFunction function) { checkArgument(length > 0d); checkArgument(function.getMax() > 0d); return new NonHomogenous(length, function); }
/** * Tests whether the predicate approximates the intensity function when a * large number of runs is done. */ @Test public void intensityApproximationPredicateTest() { final RandomGenerator rng = new MersenneTwister(123); final Predicate<Double> pred = new NHPredicate(rng, intensityFunction); final ImmutableList.Builder<Double> b = ImmutableList.builder(); for (int i = 0; i < 100; i++) { b.add(new Double(i)); } final ImmutableList<Double> doubles = b.build(); final Multiset<Double> ms = TreeMultiset.create(); final int repetitions = 10000; for (int i = 0; i < repetitions; i++) { ms.addAll(Collections2.filter(doubles, pred)); } for (final Multiset.Entry<Double> entry : ms.entrySet()) { final double prob = intensityFunction.apply(entry.getElement()) / intensityFunction.getMax(); final double observation = entry.getCount() / (double) repetitions; assertEquals(prob, observation, 0.015); } }
@SafeVarargs private static void nonZeroCheck(Range<Double> globalRange, IntensityFunction intFunc, Range<Double>... nonZeroRanges) { final RangeSet<Double> nonZeroRangeSet = asSet(nonZeroRanges); for (double d = globalRange.lowerEndpoint(); d <= globalRange .upperEndpoint(); d = Math.round((d + .01) * 100d) / 100d) { if (nonZeroRangeSet.contains(d)) { assertTrue(intFunc.apply(d) > 0); } else { assertEquals(0d, intFunc.apply(d), 0.000001); } } } }
@Override public double value(double x) { return function.apply(x); } }
NHPredicate(RandomGenerator r, IntensityFunction l) { rng = r; lambda = l; lambdaMax = lambda.getMax(); }
NonHomogenous(double l, IntensityFunction func) { super(l, func.getMax()); lambd = func; }