Skip to content

Commit 2702915

Browse files
committed
Added random generators for exponential, gamma and Nakagami distributions.
Fixes #421
1 parent b1b63a2 commit 2702915

File tree

12 files changed

+1186
-0
lines changed

12 files changed

+1186
-0
lines changed

hipparchus-core/src/changes/changes.xml

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -50,6 +50,12 @@ If the output is not quite correct, check for invisible trailing spaces!
5050
</properties>
5151
<body>
5252
<release version="4.1" date="TBD" description="TBD">
53+
<action dev="luc" type="add" issue="issues/421">
54+
Added random generators for exponential, gamma and Nakagami distributions.
55+
</action>
56+
<action dev="luc" type="add" due-to="Xavier Vasseur" issue="issues/419">
57+
Added half Cauchy, half normal and inverse gamma distributions.
58+
</action>
5359
<action dev="luc" type="add" issue="issues/406">
5460
Added AxisChecker helper class to check grid data consistency.
5561
</action>
Lines changed: 56 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,56 @@
1+
/*
2+
* Licensed to the Hipparchus project under one or more
3+
* contributor license agreements. See the NOTICE file distributed with
4+
* this work for additional information regarding copyright ownership.
5+
* The Hipparchus project licenses this file to You under the Apache License, Version 2.0
6+
* (the "License"); you may not use this file except in compliance with
7+
* the License. You may obtain a copy of the License at
8+
*
9+
* https://www.apache.org/licenses/LICENSE-2.0
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*/
17+
package org.hipparchus.random;
18+
19+
/**
20+
* Base class for generating non-uniform random variates from an underlying uniform generator.
21+
* @author Luc Maisonobe
22+
*/
23+
public abstract class AbstractNonUniformGenerator implements NonUniformGenerator
24+
{
25+
26+
/** Underlying random generator. */
27+
private final RandomGenerator random;
28+
29+
/**
30+
* Simple constructor.
31+
* @param random underlying random generator
32+
*/
33+
public AbstractNonUniformGenerator(final RandomGenerator random)
34+
{
35+
this.random = random;
36+
}
37+
38+
/**
39+
* Generate next uniform variate.
40+
* @return next uniform variate between 0 and 1
41+
*/
42+
protected double nextUniform()
43+
{
44+
return random.nextDouble();
45+
}
46+
47+
/**
48+
* Generate next normal variate.
49+
* @return next normal variate witm mean 0 and standard deviation 1
50+
*/
51+
protected double nextNormal()
52+
{
53+
return random.nextGaussian();
54+
}
55+
56+
}
Lines changed: 49 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,49 @@
1+
/*
2+
* Licensed to the Hipparchus project under one or more
3+
* contributor license agreements. See the NOTICE file distributed with
4+
* this work for additional information regarding copyright ownership.
5+
* The Hipparchus project licenses this file to You under the Apache License, Version 2.0
6+
* (the "License"); you may not use this file except in compliance with
7+
* the License. You may obtain a copy of the License at
8+
*
9+
* https://www.apache.org/licenses/LICENSE-2.0
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*/
17+
package org.hipparchus.random;
18+
19+
import org.hipparchus.util.FastMath;
20+
21+
/**
22+
* Procedure to generate variates for exponential distribution.
23+
* @author Luc Maisonobe
24+
*/
25+
public class ExponentialGenerator extends AbstractNonUniformGenerator
26+
{
27+
28+
/** Rate parameter */
29+
private final double lambda;
30+
31+
/**
32+
* Simple constructor.
33+
* @param random underlying random generator
34+
* @param lambda rate parameter
35+
*/
36+
public ExponentialGenerator(final RandomGenerator random, final double lambda)
37+
{
38+
super(random);
39+
this.lambda = lambda;
40+
}
41+
42+
/** {@inheritDoc} */
43+
@Override
44+
public double nextVariate()
45+
{
46+
return -FastMath.log(nextUniform()) / lambda;
47+
}
48+
49+
}
Lines changed: 138 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,138 @@
1+
/*
2+
* Licensed to the Hipparchus project under one or more
3+
* contributor license agreements. See the NOTICE file distributed with
4+
* this work for additional information regarding copyright ownership.
5+
* The Hipparchus project licenses this file to You under the Apache License, Version 2.0
6+
* (the "License"); you may not use this file except in compliance with
7+
* the License. You may obtain a copy of the License at
8+
*
9+
* https://www.apache.org/licenses/LICENSE-2.0
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*/
17+
package org.hipparchus.random;
18+
19+
import org.hipparchus.exception.LocalizedCoreFormats;
20+
import org.hipparchus.exception.MathIllegalStateException;
21+
import org.hipparchus.util.FastMath;
22+
23+
/**
24+
* Procedure to generate variates for gamma distribution.
25+
* <p>
26+
* The method used for generation is based on Marsaglia and Tsang
27+
* 2000 paper <a href="https://doi.org/10.1145/358407.358414">A
28+
* Simple Method for Generating Gamma Variables</a>
29+
* </p>
30+
* @author Luc Maisonobe
31+
*/
32+
public class GammaGenerator extends AbstractNonUniformGenerator
33+
{
34+
35+
/** Maximum number of rejections allowed before failing. */
36+
private static final int MAX_REJECTION = 100;
37+
38+
/** Factor used in squeeze method. */
39+
private static final double SQUEEZE_FACTOR = 0.0331;
40+
41+
/** Indicator for shape factor smaller than 1. */
42+
private final boolean smallShape;
43+
44+
/** D factor from Marsaglia and Tsang paper. */
45+
private final double d;
46+
47+
/** C factor from Marsaglia and Tsang paper. */
48+
private final double c;
49+
50+
/** Shape parameter */
51+
private final double alpha;
52+
53+
/** Scale parameter */
54+
private final double theta;
55+
56+
/**
57+
* Simple constructor.
58+
* @param random underlying random generator
59+
* @param alpha shape parameter
60+
* @param theta scale parameter (this is the inverse of the rate parameter λ)
61+
*/
62+
public GammaGenerator(final RandomGenerator random, final double alpha, final double theta)
63+
{
64+
super(random);
65+
66+
// fix small shape factor to ensure we use at least 1
67+
final double fixedShape;
68+
if (alpha < 1.0)
69+
{
70+
fixedShape = alpha + 1.0;
71+
smallShape = true;
72+
}
73+
else
74+
{
75+
fixedShape = alpha;
76+
smallShape = false;
77+
}
78+
this.d = (3.0 * fixedShape - 1.0) / 3.0;
79+
this.c = 1.0 / FastMath.sqrt(9.0 * fixedShape - 3.0);
80+
81+
this.alpha = alpha;
82+
this.theta = theta;
83+
84+
}
85+
86+
/** {@inheritDoc} */
87+
@Override
88+
public double nextVariate()
89+
{
90+
final double canonical = canonicalVariate();
91+
return theta * (smallShape ? canonical * FastMath.pow(nextUniform(), 1.0 / alpha): canonical);
92+
}
93+
94+
/**
95+
* Generate canonical variate for shape at least 1 and unit scale
96+
* @return canonical variate
97+
*/
98+
private double canonicalVariate()
99+
{
100+
for (int i = 0; i < MAX_REJECTION; ++i)
101+
{
102+
103+
// generate the normal variate
104+
double x = 0;
105+
double n = -1;
106+
for (int j = 0; n <= 0 && j < MAX_REJECTION; ++j)
107+
{
108+
x = nextNormal();
109+
n = 1 + c * x;
110+
}
111+
112+
// intermediate variate v
113+
final double v = n * n * n;
114+
final double dv = d * v;
115+
116+
// generate uniform variate
117+
final double u = nextUniform();
118+
final double x2 = x * x;
119+
120+
// the first alternative (with SQUEEZE_FACTOR) is evaluated first;
121+
// as Java short-circuits the || boolean operator, we often can
122+
// avoid the costly computation of logarithms in the second alternative
123+
if (u < 1 - SQUEEZE_FACTOR * x2 * x2 ||
124+
FastMath.log(u) < 0.5 * x * x + d - dv + d * FastMath.log(v))
125+
{
126+
return dv;
127+
}
128+
129+
// reject this attempt, we need to try again with new random variables
130+
131+
}
132+
133+
// this should never happen
134+
throw new MathIllegalStateException(LocalizedCoreFormats.INTERNAL_ERROR);
135+
136+
}
137+
138+
}
Lines changed: 50 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,50 @@
1+
/*
2+
* Licensed to the Hipparchus project under one or more
3+
* contributor license agreements. See the NOTICE file distributed with
4+
* this work for additional information regarding copyright ownership.
5+
* The Hipparchus project licenses this file to You under the Apache License, Version 2.0
6+
* (the "License"); you may not use this file except in compliance with
7+
* the License. You may obtain a copy of the License at
8+
*
9+
* https://www.apache.org/licenses/LICENSE-2.0
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*/
17+
package org.hipparchus.random;
18+
19+
import org.hipparchus.util.FastMath;
20+
21+
/**
22+
* Procedure to generate variates for Nakagami distribution.
23+
* @author Luc Maisonobe
24+
*/
25+
public class NakagamiGenerator extends AbstractNonUniformGenerator
26+
{
27+
28+
/** Associated gamma distribution generator. */
29+
private final GammaGenerator gamma;
30+
31+
/**
32+
* Simple constructor.
33+
* @param random underlying random generator
34+
* @param m shape parameter
35+
* @param omega scale parameter
36+
*/
37+
public NakagamiGenerator(final RandomGenerator random, final double m, final double omega)
38+
{
39+
super(random); // really unused
40+
this.gamma = new GammaGenerator(random, m, omega / m);
41+
}
42+
43+
/** {@inheritDoc} */
44+
@Override
45+
public double nextVariate()
46+
{
47+
return FastMath.sqrt(gamma.nextVariate());
48+
}
49+
50+
}
Lines changed: 32 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,32 @@
1+
/*
2+
* Licensed to the Hipparchus project under one or more
3+
* contributor license agreements. See the NOTICE file distributed with
4+
* this work for additional information regarding copyright ownership.
5+
* The Hipparchus project licenses this file to You under the Apache License, Version 2.0
6+
* (the "License"); you may not use this file except in compliance with
7+
* the License. You may obtain a copy of the License at
8+
*
9+
* https://www.apache.org/licenses/LICENSE-2.0
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*/
17+
package org.hipparchus.random;
18+
19+
/**
20+
* Interface for generating non-uniform random variates.
21+
* @author Luc Maisonobe
22+
*/
23+
public interface NonUniformGenerator
24+
{
25+
26+
/**
27+
* Generate next variate.
28+
* @return next variate
29+
*/
30+
double nextVariate();
31+
32+
}

0 commit comments

Comments
 (0)