001    /*
002     * Licensed to the Apache Software Foundation (ASF) under one or more
003     * contributor license agreements.  See the NOTICE file distributed with
004     * this work for additional information regarding copyright ownership.
005     * The ASF licenses this file to You under the Apache License, Version 2.0
006     * (the "License"); you may not use this file except in compliance with
007     * the License.  You may obtain a copy of the License at
008     *
009     *      http://www.apache.org/licenses/LICENSE-2.0
010     *
011     * Unless required by applicable law or agreed to in writing, software
012     * distributed under the License is distributed on an "AS IS" BASIS,
013     * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014     * See the License for the specific language governing permissions and
015     * limitations under the License.
016     */
017    
018    package org.apache.commons.math.distribution;
019    
020    import java.io.Serializable;
021    
022    import org.apache.commons.math.MathRuntimeException;
023    import org.apache.commons.math.exception.util.LocalizedFormats;
024    import org.apache.commons.math.util.FastMath;
025    
026    /**
027     * Default implementation of
028     * {@link org.apache.commons.math.distribution.CauchyDistribution}.
029     *
030     * @since 1.1
031     * @version $Revision: 1054524 $ $Date: 2011-01-03 05:59:18 +0100 (lun. 03 janv. 2011) $
032     */
033    public class CauchyDistributionImpl extends AbstractContinuousDistribution
034            implements CauchyDistribution, Serializable {
035    
036        /**
037         * Default inverse cumulative probability accuracy
038         * @since 2.1
039         */
040        public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9;
041    
042        /** Serializable version identifier */
043        private static final long serialVersionUID = 8589540077390120676L;
044    
045        /** The median of this distribution. */
046        private double median = 0;
047    
048        /** The scale of this distribution. */
049        private double scale = 1;
050    
051        /** Inverse cumulative probability accuracy */
052        private final double solverAbsoluteAccuracy;
053    
054        /**
055         * Creates cauchy distribution with the medain equal to zero and scale
056         * equal to one.
057         */
058        public CauchyDistributionImpl(){
059            this(0.0, 1.0);
060        }
061    
062        /**
063         * Create a cauchy distribution using the given median and scale.
064         * @param median median for this distribution
065         * @param s scale parameter for this distribution
066         */
067        public CauchyDistributionImpl(double median, double s){
068            this(median, s, DEFAULT_INVERSE_ABSOLUTE_ACCURACY);
069        }
070    
071        /**
072         * Create a cauchy distribution using the given median and scale.
073         * @param median median for this distribution
074         * @param s scale parameter for this distribution
075         * @param inverseCumAccuracy the maximum absolute error in inverse cumulative probability estimates
076         * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY})
077         * @since 2.1
078         */
079        public CauchyDistributionImpl(double median, double s, double inverseCumAccuracy) {
080            super();
081            setMedianInternal(median);
082            setScaleInternal(s);
083            solverAbsoluteAccuracy = inverseCumAccuracy;
084        }
085    
086        /**
087         * For this distribution, X, this method returns P(X &lt; <code>x</code>).
088         * @param x the value at which the CDF is evaluated.
089         * @return CDF evaluated at <code>x</code>.
090         */
091        public double cumulativeProbability(double x) {
092            return 0.5 + (FastMath.atan((x - median) / scale) / FastMath.PI);
093        }
094    
095        /**
096         * Access the median.
097         * @return median for this distribution
098         */
099        public double getMedian() {
100            return median;
101        }
102    
103        /**
104         * Access the scale parameter.
105         * @return scale parameter for this distribution
106         */
107        public double getScale() {
108            return scale;
109        }
110    
111        /**
112         * Returns the probability density for a particular point.
113         *
114         * @param x The point at which the density should be computed.
115         * @return The pdf at point x.
116         * @since 2.1
117         */
118        @Override
119        public double density(double x) {
120            final double dev = x - median;
121            return (1 / FastMath.PI) * (scale / (dev * dev + scale * scale));
122        }
123    
124        /**
125         * For this distribution, X, this method returns the critical point x, such
126         * that P(X &lt; x) = <code>p</code>.
127         * <p>
128         * Returns <code>Double.NEGATIVE_INFINITY</code> for p=0 and
129         * <code>Double.POSITIVE_INFINITY</code> for p=1.</p>
130         *
131         * @param p the desired probability
132         * @return x, such that P(X &lt; x) = <code>p</code>
133         * @throws IllegalArgumentException if <code>p</code> is not a valid
134         *         probability.
135         */
136        @Override
137        public double inverseCumulativeProbability(double p) {
138            double ret;
139            if (p < 0.0 || p > 1.0) {
140                throw MathRuntimeException.createIllegalArgumentException(
141                      LocalizedFormats.OUT_OF_RANGE_SIMPLE, p, 0.0, 1.0);
142            } else if (p == 0) {
143                ret = Double.NEGATIVE_INFINITY;
144            } else  if (p == 1) {
145                ret = Double.POSITIVE_INFINITY;
146            } else {
147                ret = median + scale * FastMath.tan(FastMath.PI * (p - .5));
148            }
149            return ret;
150        }
151    
152        /**
153         * Modify the median.
154         * @param median for this distribution
155         * @deprecated as of 2.1 (class will become immutable in 3.0)
156         */
157        @Deprecated
158        public void setMedian(double median) {
159            setMedianInternal(median);
160        }
161    
162        /**
163         * Modify the median.
164         * @param newMedian for this distribution
165         */
166        private void setMedianInternal(double newMedian) {
167            this.median = newMedian;
168        }
169    
170        /**
171         * Modify the scale parameter.
172         * @param s scale parameter for this distribution
173         * @throws IllegalArgumentException if <code>sd</code> is not positive.
174         * @deprecated as of 2.1 (class will become immutable in 3.0)
175         */
176        @Deprecated
177        public void setScale(double s) {
178            setScaleInternal(s);
179        }
180    
181        /**
182         * Modify the scale parameter.
183         * @param s scale parameter for this distribution
184         * @throws IllegalArgumentException if <code>sd</code> is not positive.
185         */
186        private void setScaleInternal(double s) {
187            if (s <= 0.0) {
188                throw MathRuntimeException.createIllegalArgumentException(
189                      LocalizedFormats.NOT_POSITIVE_SCALE, s);
190            }
191            scale = s;
192        }
193    
194        /**
195         * Access the domain value lower bound, based on <code>p</code>, used to
196         * bracket a CDF root.  This method is used by
197         * {@link #inverseCumulativeProbability(double)} to find critical values.
198         *
199         * @param p the desired probability for the critical value
200         * @return domain value lower bound, i.e.
201         *         P(X &lt; <i>lower bound</i>) &lt; <code>p</code>
202         */
203        @Override
204        protected double getDomainLowerBound(double p) {
205            double ret;
206    
207            if (p < .5) {
208                ret = -Double.MAX_VALUE;
209            } else {
210                ret = median;
211            }
212    
213            return ret;
214        }
215    
216        /**
217         * Access the domain value upper bound, based on <code>p</code>, used to
218         * bracket a CDF root.  This method is used by
219         * {@link #inverseCumulativeProbability(double)} to find critical values.
220         *
221         * @param p the desired probability for the critical value
222         * @return domain value upper bound, i.e.
223         *         P(X &lt; <i>upper bound</i>) &gt; <code>p</code>
224         */
225        @Override
226        protected double getDomainUpperBound(double p) {
227            double ret;
228    
229            if (p < .5) {
230                ret = median;
231            } else {
232                ret = Double.MAX_VALUE;
233            }
234    
235            return ret;
236        }
237    
238        /**
239         * Access the initial domain value, based on <code>p</code>, used to
240         * bracket a CDF root.  This method is used by
241         * {@link #inverseCumulativeProbability(double)} to find critical values.
242         *
243         * @param p the desired probability for the critical value
244         * @return initial domain value
245         */
246        @Override
247        protected double getInitialDomain(double p) {
248            double ret;
249    
250            if (p < .5) {
251                ret = median - scale;
252            } else if (p > .5) {
253                ret = median + scale;
254            } else {
255                ret = median;
256            }
257    
258            return ret;
259        }
260    
261        /**
262         * Return the absolute accuracy setting of the solver used to estimate
263         * inverse cumulative probabilities.
264         *
265         * @return the solver absolute accuracy
266         * @since 2.1
267         */
268        @Override
269        protected double getSolverAbsoluteAccuracy() {
270            return solverAbsoluteAccuracy;
271        }
272    
273        /**
274         * Returns the lower bound of the support for this distribution.
275         * The lower bound of the support of the Cauchy distribution is always
276         * negative infinity, regardless of the parameters.
277         *
278         * @return lower bound of the support (always Double.NEGATIVE_INFINITY)
279         * @since 2.2
280         */
281        public double getSupportLowerBound() {
282            return Double.NEGATIVE_INFINITY;
283        }
284    
285        /**
286         * Returns the upper bound of the support for this distribution.
287         * The upper bound of the support of the Cauchy distribution is always
288         * positive infinity, regardless of the parameters.
289         *
290         * @return upper bound of the support (always Double.POSITIVE_INFINITY)
291         * @since 2.2
292         */
293        public double getSupportUpperBound() {
294            return Double.POSITIVE_INFINITY;
295        }
296    
297        /**
298         * Returns the mean.
299         *
300         * The mean is always undefined, regardless of the parameters.
301         *
302         * @return mean (always Double.NaN)
303         * @since 2.2
304         */
305        public double getNumericalMean() {
306            return Double.NaN;
307        }
308    
309        /**
310         * Returns the variance.
311         *
312         * The variance is always undefined, regardless of the parameters.
313         *
314         * @return variance (always Double.NaN)
315         * @since 2.2
316         */
317        public double getNumericalVariance() {
318            return Double.NaN;
319        }
320    }