1 /*
2 * Licensed to the Apache Software Foundation (ASF) 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 ASF 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
18 /*
19 * This is not the original file distributed by the Apache Software Foundation
20 * It has been modified by the Hipparchus project
21 */
22 package org.hipparchus.distribution.continuous;
23
24 import org.hipparchus.exception.MathIllegalArgumentException;
25 import org.hipparchus.special.Erf;
26 import org.hipparchus.util.FastMath;
27 import org.hipparchus.util.MathUtils;
28
29 /**
30 * This class implements the <a href="http://en.wikipedia.org/wiki/L%C3%A9vy_distribution">
31 * Lévy distribution</a>.
32 */
33 public class LevyDistribution extends AbstractRealDistribution {
34
35 /** Serializable UID. */
36 private static final long serialVersionUID = 20130314L;
37
38 /** Location parameter. */
39 private final double mu;
40
41 /** Scale parameter. */
42 private final double c; // Setting this to 1 returns a cumProb of 1.0
43
44 /** Half of c (for calculations). */
45 private final double halfC;
46
47 /**
48 * Build a new instance.
49 *
50 * @param mu location parameter
51 * @param c scale parameter
52 */
53 public LevyDistribution(final double mu, final double c) {
54 super();
55 this.mu = mu;
56 this.c = c;
57 this.halfC = 0.5 * c;
58 }
59
60
61 /** {@inheritDoc}
62 * <p>
63 * From Wikipedia: The probability density function of the Lévy distribution
64 * over the domain is
65 * </p>
66 * \[
67 * f(x; \mu, c) = \sqrt{\frac{c}{2\pi}} \frac{e^{\frac{-c}{2 (x - \mu)}}}{(x - \mu)^\frac{3}{2}}
68 * \]
69 * <p>
70 * For this distribution, {@code X}, this method returns {@code P(X < x)}.
71 * If {@code x} is less than location parameter μ, {@code Double.NaN} is
72 * returned, as in these cases the distribution is not defined.
73 * </p>
74 */
75 @Override
76 public double density(final double x) {
77 if (x < mu) {
78 return Double.NaN;
79 }
80
81 final double delta = x - mu;
82 final double f = halfC / delta;
83 return FastMath.sqrt(f / FastMath.PI) * FastMath.exp(-f) /delta;
84 }
85
86 /** {@inheritDoc}
87 *
88 * See documentation of {@link #density(double)} for computation details.
89 */
90 @Override
91 public double logDensity(double x) {
92 if (x < mu) {
93 return Double.NaN;
94 }
95
96 final double delta = x - mu;
97 final double f = halfC / delta;
98 return 0.5 * FastMath.log(f / FastMath.PI) - f - FastMath.log(delta);
99 }
100
101 /** {@inheritDoc}
102 * <p>
103 * From Wikipedia: the cumulative distribution function is
104 * </p>
105 * <pre>
106 * f(x; u, c) = erfc (√ (c / 2 (x - u )))
107 * </pre>
108 */
109 @Override
110 public double cumulativeProbability(final double x) {
111 if (x < mu) {
112 return Double.NaN;
113 }
114 return Erf.erfc(FastMath.sqrt(halfC / (x - mu)));
115 }
116
117 /** {@inheritDoc} */
118 @Override
119 public double inverseCumulativeProbability(final double p) throws MathIllegalArgumentException {
120 MathUtils.checkRangeInclusive(p, 0, 1);
121 final double t = Erf.erfcInv(p);
122 return mu + halfC / (t * t);
123 }
124
125 /** Get the scale parameter of the distribution.
126 * @return scale parameter of the distribution
127 */
128 public double getScale() {
129 return c;
130 }
131
132 /** Get the location parameter of the distribution.
133 * @return location parameter of the distribution
134 */
135 public double getLocation() {
136 return mu;
137 }
138
139 /** {@inheritDoc} */
140 @Override
141 public double getNumericalMean() {
142 return Double.POSITIVE_INFINITY;
143 }
144
145 /** {@inheritDoc} */
146 @Override
147 public double getNumericalVariance() {
148 return Double.POSITIVE_INFINITY;
149 }
150
151 /** {@inheritDoc} */
152 @Override
153 public double getSupportLowerBound() {
154 return mu;
155 }
156
157 /** {@inheritDoc} */
158 @Override
159 public double getSupportUpperBound() {
160 return Double.POSITIVE_INFINITY;
161 }
162
163 /** {@inheritDoc} */
164 @Override
165 public boolean isSupportConnected() {
166 return true;
167 }
168
169 }