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 */ 017package org.apache.commons.math3.distribution; 018 019import org.apache.commons.math3.exception.NotStrictlyPositiveException; 020import org.apache.commons.math3.exception.OutOfRangeException; 021import org.apache.commons.math3.exception.util.LocalizedFormats; 022import org.apache.commons.math3.random.RandomGenerator; 023import org.apache.commons.math3.random.Well19937c; 024import org.apache.commons.math3.util.FastMath; 025 026/** 027 * Implementation of the Cauchy distribution. 028 * 029 * @see <a href="http://en.wikipedia.org/wiki/Cauchy_distribution">Cauchy distribution (Wikipedia)</a> 030 * @see <a href="http://mathworld.wolfram.com/CauchyDistribution.html">Cauchy Distribution (MathWorld)</a> 031 * @since 1.1 (changed to concrete class in 3.0) 032 */ 033public class CauchyDistribution extends AbstractRealDistribution { 034 /** 035 * Default inverse cumulative probability accuracy. 036 * @since 2.1 037 */ 038 public static final double DEFAULT_INVERSE_ABSOLUTE_ACCURACY = 1e-9; 039 /** Serializable version identifier */ 040 private static final long serialVersionUID = 8589540077390120676L; 041 /** The median of this distribution. */ 042 private final double median; 043 /** The scale of this distribution. */ 044 private final double scale; 045 /** Inverse cumulative probability accuracy */ 046 private final double solverAbsoluteAccuracy; 047 048 /** 049 * Creates a Cauchy distribution with the median equal to zero and scale 050 * equal to one. 051 */ 052 public CauchyDistribution() { 053 this(0, 1); 054 } 055 056 /** 057 * Creates a Cauchy distribution using the given median and scale. 058 * <p> 059 * <b>Note:</b> this constructor will implicitly create an instance of 060 * {@link Well19937c} as random generator to be used for sampling only (see 061 * {@link #sample()} and {@link #sample(int)}). In case no sampling is 062 * needed for the created distribution, it is advised to pass {@code null} 063 * as random generator via the appropriate constructors to avoid the 064 * additional initialisation overhead. 065 * 066 * @param median Median for this distribution. 067 * @param scale Scale parameter for this distribution. 068 */ 069 public CauchyDistribution(double median, double scale) { 070 this(median, scale, DEFAULT_INVERSE_ABSOLUTE_ACCURACY); 071 } 072 073 /** 074 * Creates a Cauchy distribution using the given median and scale. 075 * <p> 076 * <b>Note:</b> this constructor will implicitly create an instance of 077 * {@link Well19937c} as random generator to be used for sampling only (see 078 * {@link #sample()} and {@link #sample(int)}). In case no sampling is 079 * needed for the created distribution, it is advised to pass {@code null} 080 * as random generator via the appropriate constructors to avoid the 081 * additional initialisation overhead. 082 * 083 * @param median Median for this distribution. 084 * @param scale Scale parameter for this distribution. 085 * @param inverseCumAccuracy Maximum absolute error in inverse 086 * cumulative probability estimates 087 * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}). 088 * @throws NotStrictlyPositiveException if {@code scale <= 0}. 089 * @since 2.1 090 */ 091 public CauchyDistribution(double median, double scale, 092 double inverseCumAccuracy) { 093 this(new Well19937c(), median, scale, inverseCumAccuracy); 094 } 095 096 /** 097 * Creates a Cauchy distribution. 098 * 099 * @param rng Random number generator. 100 * @param median Median for this distribution. 101 * @param scale Scale parameter for this distribution. 102 * @throws NotStrictlyPositiveException if {@code scale <= 0}. 103 * @since 3.3 104 */ 105 public CauchyDistribution(RandomGenerator rng, double median, double scale) { 106 this(rng, median, scale, DEFAULT_INVERSE_ABSOLUTE_ACCURACY); 107 } 108 109 /** 110 * Creates a Cauchy distribution. 111 * 112 * @param rng Random number generator. 113 * @param median Median for this distribution. 114 * @param scale Scale parameter for this distribution. 115 * @param inverseCumAccuracy Maximum absolute error in inverse 116 * cumulative probability estimates 117 * (defaults to {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}). 118 * @throws NotStrictlyPositiveException if {@code scale <= 0}. 119 * @since 3.1 120 */ 121 public CauchyDistribution(RandomGenerator rng, 122 double median, 123 double scale, 124 double inverseCumAccuracy) { 125 super(rng); 126 if (scale <= 0) { 127 throw new NotStrictlyPositiveException(LocalizedFormats.SCALE, scale); 128 } 129 this.scale = scale; 130 this.median = median; 131 solverAbsoluteAccuracy = inverseCumAccuracy; 132 } 133 134 /** {@inheritDoc} */ 135 public double cumulativeProbability(double x) { 136 return 0.5 + (FastMath.atan((x - median) / scale) / FastMath.PI); 137 } 138 139 /** 140 * Access the median. 141 * 142 * @return the median for this distribution. 143 */ 144 public double getMedian() { 145 return median; 146 } 147 148 /** 149 * Access the scale parameter. 150 * 151 * @return the scale parameter for this distribution. 152 */ 153 public double getScale() { 154 return scale; 155 } 156 157 /** {@inheritDoc} */ 158 public double density(double x) { 159 final double dev = x - median; 160 return (1 / FastMath.PI) * (scale / (dev * dev + scale * scale)); 161 } 162 163 /** 164 * {@inheritDoc} 165 * 166 * Returns {@code Double.NEGATIVE_INFINITY} when {@code p == 0} 167 * and {@code Double.POSITIVE_INFINITY} when {@code p == 1}. 168 */ 169 @Override 170 public double inverseCumulativeProbability(double p) throws OutOfRangeException { 171 double ret; 172 if (p < 0 || p > 1) { 173 throw new OutOfRangeException(p, 0, 1); 174 } else if (p == 0) { 175 ret = Double.NEGATIVE_INFINITY; 176 } else if (p == 1) { 177 ret = Double.POSITIVE_INFINITY; 178 } else { 179 ret = median + scale * FastMath.tan(FastMath.PI * (p - .5)); 180 } 181 return ret; 182 } 183 184 /** {@inheritDoc} */ 185 @Override 186 protected double getSolverAbsoluteAccuracy() { 187 return solverAbsoluteAccuracy; 188 } 189 190 /** 191 * {@inheritDoc} 192 * 193 * The mean is always undefined no matter the parameters. 194 * 195 * @return mean (always Double.NaN) 196 */ 197 public double getNumericalMean() { 198 return Double.NaN; 199 } 200 201 /** 202 * {@inheritDoc} 203 * 204 * The variance is always undefined no matter the parameters. 205 * 206 * @return variance (always Double.NaN) 207 */ 208 public double getNumericalVariance() { 209 return Double.NaN; 210 } 211 212 /** 213 * {@inheritDoc} 214 * 215 * The lower bound of the support is always negative infinity no matter 216 * the parameters. 217 * 218 * @return lower bound of the support (always Double.NEGATIVE_INFINITY) 219 */ 220 public double getSupportLowerBound() { 221 return Double.NEGATIVE_INFINITY; 222 } 223 224 /** 225 * {@inheritDoc} 226 * 227 * The upper bound of the support is always positive infinity no matter 228 * the parameters. 229 * 230 * @return upper bound of the support (always Double.POSITIVE_INFINITY) 231 */ 232 public double getSupportUpperBound() { 233 return Double.POSITIVE_INFINITY; 234 } 235 236 /** {@inheritDoc} */ 237 public boolean isSupportLowerBoundInclusive() { 238 return false; 239 } 240 241 /** {@inheritDoc} */ 242 public boolean isSupportUpperBoundInclusive() { 243 return false; 244 } 245 246 /** 247 * {@inheritDoc} 248 * 249 * The support of this distribution is connected. 250 * 251 * @return {@code true} 252 */ 253 public boolean isSupportConnected() { 254 return true; 255 } 256}