/* *
* Title:
* Description:
* Copyright: Copyright (c)
* Company:
* @author
* @version 1.0
*/
package umontreal.iro.lecuyer.stochprocess;
import umontreal.iro.lecuyer.rng.*;
import umontreal.iro.lecuyer.randvar.*;
import umontreal.iro.lecuyer.probdist.GammaDist;
import umontreal.iro.lecuyer.util.Num;
/**
* This class represents a gamma process
*
* {S(t) = G(t;μ, ν) : t >= 0} with mean parameter μ and
* variance parameter ν.
* It is a continuous-time process with stationary,
* independent gamma increments such that for any
* Δt > 0,
*
*
* where X is a random variate from the gamma distribution
* Gamma
* (μ2Δt/ν, μ/ν).
*
*
* In this class, the gamma process is sampled sequentially using
* equation.
*
*/
public class GammaProcess extends StochasticProcess {
// Make sure that process is strictly increasing: when two successive steps
// are equal, reset the new one to old*(1 + EPS)
protected static final double EPS = 1.0e-15;
protected boolean usesAnti = false;
protected RandomStream stream; // used in the generatePath(), so must be kept
protected GammaGen Ggen;
protected double mu,
nu;
protected double muOverNu,
mu2OverNu;
// For precomputations for standard GP
protected double[] mu2dtOverNu;
protected void setLarger (double[] path, int left, int mid, int right) {
/* make sure that path[left] < path[mid] < path[right] by adding a
small epsilon to the last two */
double myEps;
if (path[left] >= 0)
myEps = EPS;
else
myEps = -EPS;
path[mid] = path[left] * (1.0 + myEps);
if (path[mid] >= path[right])
path[right] = path[mid] * (1.0 + myEps);
}
protected double setLarger (double[] path, int left, int right) {
/* make sure that path[left] < s < path[right] by adding a
small epsilon to the last two; return s. */
double myEps;
if (path[left] >= 0)
myEps = EPS;
else
myEps = -EPS;
double s = path[left] * (1.0 + myEps);
if (s >= path[right])
path[right] = s * (1.0 + myEps);
return s;
}
protected double setLarger (double v) {
/* return a number that is slightly larger than v */
if (v >= 0)
return v*(1.0 + EPS) + Num.DBL_EPSILON;
else
return v*(1.0 - EPS);
}
/**
* Constructs a new GammaProcess with parameters
* μ = mu ,
* ν = nu and initial value
* S(t0) = s0 .
* The gamma variates X in are generated by inversion
* using stream.
*
*/
public GammaProcess (double s0, double mu, double nu,
RandomStream stream) {
this (s0, mu, nu, new GammaGen (stream, 1.0));
// Note that the parameters (the 1.0) supplied to
// the GammaGen constructors are meaningless.
// The 'real' parameters are supplied at time of generation
}
/**
* Constructs a new GammaProcess with parameters
* μ = mu ,
* ν = nu and initial value
* S(t0) = s0 .
* The gamma variates X in are supplied by the gamma random
* variate generator Ggen. Note that the parameters of the
* {@link umontreal.iro.lecuyer.randvar.GammaGen GammaGen} object
* Ggen are not important since the implementation forces the generator
* to use the correct parameters (as defined above).
*
*/
public GammaProcess (double s0, double mu, double nu, GammaGen Ggen) {
this.mu = mu;
this.nu = nu;
this.x0 = s0;
this.Ggen = Ggen;
stream = Ggen.getStream();
}
public double nextObservation() {
double s = path[observationIndex];
double v = s;
s += Ggen.nextDouble (stream, mu2dtOverNu[observationIndex], muOverNu);
if (s <= v)
s = setLarger (v);
observationIndex++;
observationCounter = observationIndex;
path[observationIndex] = s;
return s;
}
/**
* Generates and returns the next observation at time
* tj+1 = nextTime ,
* using the previous observation time tj defined earlier
* (either by this method or by setObservationTimes),
* as well as the value of the previous observation X(tj).
* Warning: This method will reset the observations time tj+1
* for this process to nextT. The user must make sure that
* the tj+1 supplied is
* >= tj.
*
*/
public double nextObservation (double nextT) {
// This method is useful for generating variance gamma processes
double s = path[observationIndex];
double v = s;
double previousT = t[observationIndex];
observationIndex++;
observationCounter = observationIndex;
t[observationIndex] = nextT;
double dt = nextT - previousT;
s += Ggen.nextDouble (stream, mu2OverNu*dt, muOverNu);
if (s <= v)
s = setLarger (v);
path[observationIndex] = s;
return s;
}
/**
* Generates, returns and saves the path
*
* {X(t0), X(t1),…, X(td)}.
* The gamma variates X in are generated using the
* {@link umontreal.iro.lecuyer.rng.RandomStream RandomStream} stream
* or the {@link umontreal.iro.lecuyer.rng.RandomStream RandomStream} included
* in the {@link umontreal.iro.lecuyer.randvar.GammaGen GammaGen} Ggen.
*
*/
public double[] generatePath() {
double s = x0;
double v;
for (int i = 0; i < d; i++) {
v = s;
s += Ggen.nextDouble(stream, mu2dtOverNu[i], muOverNu);
if (v <= s)
s = setLarger (v);
path[i+1] = s;
}
observationIndex = d;
observationCounter = d;
return path;
}
/**
* Generates, returns and saves the path
* {X(t0), X(t1),…, X(td)}. This method does not use the
* {@link umontreal.iro.lecuyer.rng.RandomStream RandomStream} stream nor the
* {@link umontreal.iro.lecuyer.randvar.GammaGen GammaGen} Ggen. It
* uses the vector of uniform random numbers U(0, 1) provided by the user
* and generates the path by inversion. The vector uniform01 must be of
* dimension d.
*
*/
public double[] generatePath (double[] uniform01) {
double s = x0;
double v;
for (int i = 0; i < d; i++) {
v = s;
s += GammaDist.inverseF(mu2dtOverNu[i], muOverNu, 10, uniform01[i]);
if (v <= s)
s = setLarger (v);
path[i+1] = s;
}
observationIndex = d;
observationCounter = d;
return path;
}
/**
* Sets the parameters
* S(t0) = s0 ,
*
* μ = mu and
* ν = nu of the process.
* Warning: This method will recompute some quantities stored internally,
* which may be slow if called repeatedly.
*
*/
public void setParams (double s0, double mu, double nu) {
this.x0 = s0;
this.mu = mu;
this.nu = nu;
if (observationTimesSet) init(); // Otherwise no need to.
}
/**
* Returns the value of the parameter μ.
*
*/
public double getMu() { return mu; }
/**
* Returns the value of the parameter ν.
*
*/
public double getNu() { return nu; }
/**
* Resets the {@link umontreal.iro.lecuyer.rng.RandomStream RandomStream}
* of the {@link umontreal.iro.lecuyer.randvar.GammaGen GammaGen}
* to stream.
*
*/
public void setStream (RandomStream stream) {
this.stream = stream;
Ggen.setStream (stream);
}
/**
* Returns the {@link umontreal.iro.lecuyer.rng.RandomStream RandomStream} stream.
*
*/
public RandomStream getStream() { return stream; }
protected void init() {
// init of StochasticProcess only sets path[0] = x0 if observation times are set.
super.init();
muOverNu = mu / nu;
mu2OverNu = mu * mu / nu;
mu2dtOverNu = new double[d];
if (observationTimesSet) {
for (int i = 0; i < d; i++)
mu2dtOverNu[i] = mu2OverNu * (t[i+1] - t[i]);
}
}
}
Ads help maintain this website.