    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    }