com.github.rinde.rinsim.scenario.generator.TimeSeries Maven / Gradle / Ivy
/*
* Copyright (C) 2011-2015 Rinde van Lon, iMinds-DistriNet, KU Leuven
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package com.github.rinde.rinsim.scenario.generator;
import static com.google.common.base.Preconditions.checkArgument;
import java.util.Iterator;
import java.util.List;
import javax.annotation.Nullable;
import org.apache.commons.math3.distribution.ExponentialDistribution;
import org.apache.commons.math3.distribution.NormalDistribution;
import org.apache.commons.math3.distribution.RealDistribution;
import org.apache.commons.math3.distribution.UniformRealDistribution;
import org.apache.commons.math3.random.MersenneTwister;
import org.apache.commons.math3.random.RandomGenerator;
import com.github.rinde.rinsim.scenario.generator.IntensityFunctions.IntensityFunction;
import com.github.rinde.rinsim.util.StochasticSupplier;
import com.github.rinde.rinsim.util.StochasticSuppliers;
import com.google.common.base.Predicate;
import com.google.common.collect.AbstractSequentialIterator;
import com.google.common.collect.ImmutableList;
import com.google.common.collect.Iterators;
import com.google.common.collect.Range;
/**
* Utilities for generating time series.
* @author Rinde van Lon
*/
public final class TimeSeries {
private TimeSeries() {}
/**
* Creates a homogenous Poisson process of the specified length. The intensity
* is calculated as numEvents / length.
* @param length The length of Poisson process, all generated times will be in
* the interval [0,length).
* @param numEvents The number of events that will be generated (on average).
* @return A newly constructed Poisson process {@link TimeSeriesGenerator}.
*/
public static TimeSeriesGenerator homogenousPoisson(double length,
int numEvents) {
checkArgument(length > 0d);
checkArgument(numEvents > 0);
return new PoissonProcess(length, numEvents / length);
}
/**
* 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].
*
* References
*
* - Lewis, P.A.W. and Shedler, G.S. Simulation of nonhomogenous Poisson
* processes by thinning. Naval Research Logistic Quarterly 26, (1979),
* 403–414.
*
* @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);
}
/**
* Creates a non-homogenous Poisson process of the specified length. The
* intensity is specified by the {@link StochasticSupplier}. Each time
* {@link TimeSeriesGenerator#generate(long)} is called, a new
* {@link IntensityFunction} is requested from the {@link StochasticSupplier}.
* The non-homogenous Poisson process is implemented using the thinning method
* as described in [1].
*
* References
*
* - Lewis, P.A.W. and Shedler, G.S. Simulation of nonhomogenous Poisson
* processes by thinning. Naval Research Logistic Quarterly 26, (1979),
* 403–414.
*
* @param length The length of Poisson process, all generated times will be in
* the interval [0,length).
* @param functionSupplier The intensity function supplier.
* @return A newly constructed non-homogenous Poisson process
* {@link TimeSeriesGenerator}.
*/
public static TimeSeriesGenerator nonHomogenousPoisson(double length,
StochasticSupplier functionSupplier) {
checkArgument(length > 0d);
return new SuppliedNonHomogenous(length, functionSupplier);
}
/**
* Creates a {@link TimeSeriesGenerator} using a uniform distribution for the
* inter arrival times of events.
* @param length The length of the time series, all generated times will be in
* the interval [0,length).
* @param numEvents The total number of events in the time series (on
* average).
* @param maxDeviation The maximum deviation from the mean for the uniform
* distribution.
* @return A {@link TimeSeriesGenerator} based on a uniform distribution.
*/
public static TimeSeriesGenerator uniform(double length, int numEvents,
double maxDeviation) {
checkArgument(length > 0d);
checkArgument(numEvents > 0);
checkArgument(maxDeviation > 0d);
final double average = length / numEvents;
return new UniformTimeSeries(length, average,
StochasticSuppliers.constant(maxDeviation));
}
/**
* Creates a {@link TimeSeriesGenerator} using a uniform distribution for the
* inter arrival times of events. The spread of the uniform distribution is
* defined by the maxDeviation {@link StochasticSupplier}, this
* means that each time a time series is generated a different max deviation
* settings is used.
* @param length The length of the time series, all generated times will be in
* the interval [0,length).
* @param numEvents The total number of events in the time series (on
* average).
* @param maxDeviation A supplier that is used for max deviation values.
* @return A {@link TimeSeriesGenerator} based on a uniform distribution, each
* time series that is generated has a different max deviation drawn
* from the supplier.
*/
public static TimeSeriesGenerator uniform(double length, int numEvents,
StochasticSupplier maxDeviation) {
checkArgument(length > 0d);
checkArgument(numEvents > 0);
final double average = length / numEvents;
return new UniformTimeSeries(length, average, StochasticSuppliers.checked(
maxDeviation, Range.atLeast(0d)));
}
/**
* Creates a {@link TimeSeriesGenerator} that uses a truncated normal
* distribution for the inter arrival times of events. The normal distribution
* is truncated at a lower bound of 0 since it is not allowed (it
* makes no sense) to have negative inter arrival times. For more information
* about the normal distribution see {@link StochasticSuppliers#normal()}.
* @param length The length of the time series, all generated times will be in
* the interval [0,length).
* @param numEvents The total number of events in the time series (on
* average).
* @param sd The standard deviation of the normal distribution.
* @return A {@link TimeSeriesGenerator} based on a normal distribution.
*/
public static TimeSeriesGenerator normal(double length, int numEvents,
double sd) {
checkArgument(length > 0d);
checkArgument(numEvents > 0);
final double average = length / numEvents;
return toTimeSeries(length, StochasticSuppliers.normal()
.mean(average)
.std(sd)
.lowerBound(0d)
.redrawWhenOutOfBounds()
.scaleMean()
.buildDouble());
}
/**
* Converts a {@link StochasticSupplier} of {@link Double}s into a
* {@link TimeSeriesGenerator}. Each time in the time series is created by
* t[n] = t[n-1] + ss.get(..), here ss is the
* {@link Double} supplier.
* @param length The length of the time series, all generated times will be in
* the interval [0,length).
* @param interArrivalTimesSupplier The supplier to use for computing event
* inter arrival times.
* @return A new {@link TimeSeriesGenerator} based on the supplier.
*/
public static TimeSeriesGenerator toTimeSeries(double length,
StochasticSupplier interArrivalTimesSupplier) {
return new SupplierTimeSeries(length, interArrivalTimesSupplier);
}
/**
* Decorates the specified {@link TimeSeriesGenerator} such that it only
* generates time series which conform to the specified {@link Predicate}.
* Predicates can be combined by using the methods provided by
* {@link com.google.common.base.Predicates}. Note that when an impossible
* {@link Predicate} is specified, such as
* {@link com.google.common.base.Predicates#alwaysFalse()} the resulting
* {@link TimeSeriesGenerator} will enter an infinite loop.
* @param tsg The {@link TimeSeriesGenerator} to filter.
* @param predicate All returned {@link TimeSeriesGenerator}s will conform to
* this predicate.
* @return A filtered generator.
*/
public static TimeSeriesGenerator filter(TimeSeriesGenerator tsg,
Predicate> predicate) {
return new FilteredTSG(tsg, predicate);
}
/**
* Creates a {@link Predicate} for a specified number of events. This
* predicate only accepts time series with exactly num events.
* @param num The number of events a time series should have.
* @return A newly created predicate.
*/
public static Predicate> numEventsPredicate(final int num) {
return new Predicate>() {
@Override
public boolean apply(@Nullable List input) {
assert input != null;
return input.size() == num;
}
};
}
/**
* Generator of a time series.
* @author Rinde van Lon
*/
public interface TimeSeriesGenerator {
/**
* Should generate a time series.
* @param seed The random seed to use.
* @return An immutable list of times in ascending order, may contain
* duplicates.
*/
ImmutableList generate(long seed);
}
static class FilteredTSG implements TimeSeriesGenerator {
private final TimeSeriesGenerator delegate;
private final Predicate> predicate;
private final RandomGenerator rng;
FilteredTSG(TimeSeriesGenerator tsg, Predicate> pred) {
delegate = tsg;
predicate = pred;
rng = new MersenneTwister();
}
@Override
public ImmutableList generate(long seed) {
rng.setSeed(seed);
while (true) {
final ImmutableList timeSeries = delegate.generate(rng
.nextLong());
if (predicate.apply(timeSeries)) {
return timeSeries;
}
}
}
}
static class PoissonProcess implements TimeSeriesGenerator {
/**
* Random generator used for drawing random numbers.
*/
protected final RandomGenerator rng;
final double length;
final double intensity;
PoissonProcess(double len, double intens) {
length = len;
intensity = intens;
rng = new MersenneTwister();
}
/**
* All times will be in the interval [0,length)
* @return The upper bound of the interval.
*/
public double getLength() {
return length;
}
// internal use only!
Iterator iterator() {
return new TimeSeriesIterator(new ExponentialDistribution(rng,
1d / intensity,
ExponentialDistribution.DEFAULT_INVERSE_ABSOLUTE_ACCURACY),
length);
}
@Override
public ImmutableList generate(long seed) {
rng.setSeed(seed);
return ImmutableList.copyOf(iterator());
}
}
static class NonHomogenous extends PoissonProcess {
final IntensityFunction lambd;
NonHomogenous(double l, IntensityFunction func) {
super(l, func.getMax());
lambd = func;
}
@Override
public Iterator iterator() {
return Iterators.filter(super.iterator(), new NHPredicate(rng, lambd));
}
}
static class SuppliedNonHomogenous implements TimeSeriesGenerator {
final double length;
final StochasticSupplier lambdSup;
final RandomGenerator rng;
SuppliedNonHomogenous(double l,
StochasticSupplier funcSup) {
length = l;
lambdSup = funcSup;
rng = new MersenneTwister();
}
@Override
public ImmutableList generate(long seed) {
rng.setSeed(seed);
final TimeSeriesGenerator tsg = new NonHomogenous(length,
lambdSup.get(rng.nextLong()));
return tsg.generate(rng.nextLong());
}
}
static class SupplierTimeSeries implements TimeSeriesGenerator {
private final double length;
private final StochasticSupplier supplier;
SupplierTimeSeries(double len, StochasticSupplier sup) {
length = len;
supplier = sup;
}
@Override
public ImmutableList generate(long seed) {
return ImmutableList.copyOf(new SupplierIterator(length, supplier,
new MersenneTwister(seed)));
}
}
static class SupplierIterator extends AbstractSequentialIterator {
private final double length;
private final StochasticSupplier supplier;
private final RandomGenerator randomNumberGenerator;
SupplierIterator(double len, StochasticSupplier sup,
RandomGenerator rng) {
super(next(0d, len, sup, rng));
length = len;
supplier = sup;
randomNumberGenerator = rng;
}
@SuppressWarnings("null")
@Nullable
@Override
protected Double computeNext(Double previous) {
return next(previous, length, supplier, randomNumberGenerator);
}
@Nullable
static Double next(Double prev, double len,
StochasticSupplier supplier, RandomGenerator rng) {
final double nextVal = prev + getValue(supplier, rng);
if (nextVal < len) {
return nextVal;
}
return null;
}
static double getValue(StochasticSupplier ed, RandomGenerator rng) {
final double sample = ed.get(rng.nextLong());
checkArgument(
sample >= 0d,
"A StochasticSupplier used in a TimeSeries may not return negative values, was: %s.",
sample);
return sample;
}
}
static class UniformTimeSeries implements TimeSeriesGenerator {
static final double SMALLEST_DEVIATION = .0000001;
private final RandomGenerator rng;
private final double length;
private final double average;
private final StochasticSupplier deviationSupplier;
UniformTimeSeries(double len, double avg, StochasticSupplier dev) {
rng = new MersenneTwister();
length = len;
average = avg;
deviationSupplier = dev;
}
@Override
public ImmutableList generate(long seed) {
rng.setSeed(seed);
double deviation = deviationSupplier.get(rng.nextLong());
deviation = Math.min(average, deviation);
final double lowerBound = average - deviation;
final double upperBound = average + deviation;
if (deviation < SMALLEST_DEVIATION) {
return ImmutableList.copyOf(new FixedTimeSeriesIterator(rng, length,
average));
}
return ImmutableList.copyOf(new TimeSeriesIterator(
new UniformRealDistribution(rng, lowerBound, upperBound), length));
}
}
static class NormalTimeSeries implements TimeSeriesGenerator {
private final double length;
private final RealDistribution distribution;
NormalTimeSeries(double len, double avg, double sd) {
length = len;
distribution = new NormalDistribution(avg, sd);
}
@Override
public ImmutableList generate(long seed) {
distribution.reseedRandomGenerator(seed);
return ImmutableList.copyOf(new TimeSeriesIterator(
distribution, length));
}
}
static class FixedTimeSeriesIterator extends
AbstractSequentialIterator {
private final double length;
private final double average;
protected FixedTimeSeriesIterator(RandomGenerator rng, double len,
double avg) {
super(rng.nextDouble() * avg);
length = len;
average = avg;
}
@SuppressWarnings("null")
@Nullable
@Override
protected Double computeNext(Double prev) {
final double nextVal = prev + average;
if (nextVal < length) {
return nextVal;
}
return null;
}
}
static class TimeSeriesIterator extends AbstractSequentialIterator {
private final RealDistribution ed;
private final double length;
TimeSeriesIterator(RealDistribution distr, double len) {
super(next(0d, len, distr));
length = len;
ed = distr;
}
@SuppressWarnings("null")
@Nullable
@Override
protected Double computeNext(Double previous) {
return next(previous, length, ed);
}
@Nullable
static Double next(Double prev, double len, RealDistribution ed) {
final double nextVal = prev + getFirstPositive(ed);
if (nextVal < len) {
return nextVal;
}
return null;
}
static double getFirstPositive(RealDistribution ed) {
double sample = ed.sample();
while (sample < 0) {
sample = ed.sample();
}
return sample;
}
}
static class NHPredicate implements Predicate {
private final RandomGenerator rng;
private final IntensityFunction lambda;
private final double lambdaMax;
NHPredicate(RandomGenerator r, IntensityFunction l) {
rng = r;
lambda = l;
lambdaMax = lambda.getMax();
}
@Override
public boolean apply(@Nullable Double input) {
assert input != null;
return rng.nextDouble() <= lambda.apply(input) / lambdaMax;
}
}
}