math/distributions/laplace.hpp

0001 //  Copyright Thijs van den Berg, 2008.
0002 //  Copyright John Maddock 2008.
0003 //  Copyright Paul A. Bristow 2008, 2014.
0004
0005 //  Use, modification and distribution are subject to the
0006 //  Boost Software License, Version 1.0. (See accompanying file
0007 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
0008
0009 // This module implements the Laplace distribution.
0010 // Weisstein, Eric W. "Laplace Distribution." From MathWorld--A Wolfram Web Resource.
0011 // http://mathworld.wolfram.com/LaplaceDistribution.html
0012 // http://en.wikipedia.org/wiki/Laplace_distribution
0013 //
0014 // Abramowitz and Stegun 1972, p 930
0015 // http://www.math.sfu.ca/~cbm/aands/page_930.htm
0016
0017 #ifndef BOOST_STATS_LAPLACE_HPP
0018 #define BOOST_STATS_LAPLACE_HPP
0019
0020 #include <boost/math/special_functions/log1p.hpp>
0021 #include <boost/math/distributions/detail/common_error_handling.hpp>
0022 #include <boost/math/distributions/complement.hpp>
0023 #include <boost/math/constants/constants.hpp>
0024 #include <limits>
0025
0026 namespace boost{ namespace math{
0027
0028 #ifdef _MSC_VER
0029 #  pragma warning(push)
0030 #  pragma warning(disable:4127) // conditional expression is constant
0031 #endif
0032
0033 template <class RealType = double, class Policy = policies::policy<> >
0034 class laplace_distribution
0035 {
0036 public:
0037    // ----------------------------------
0038    // public Types
0039    // ----------------------------------
0040    using value_type = RealType;
0041    using policy_type = Policy;
0042
0043    // ----------------------------------
0044    // Constructor(s)
0045    // ----------------------------------
0046    explicit laplace_distribution(RealType l_location = 0, RealType l_scale = 1)
0047       : m_location(l_location), m_scale(l_scale)
0048    {
0049       RealType result;
0050       check_parameters("boost::math::laplace_distribution<%1%>::laplace_distribution()", &result);
0051    }
0052
0053
0054    // ----------------------------------
0055    // Public functions
0056    // ----------------------------------
0057
0058    RealType location() const
0059    {
0060       return m_location;
0061    }
0062
0063    RealType scale() const
0064    {
0065       return m_scale;
0066    }
0067
0068    bool check_parameters(const char* function, RealType* result) const
0069    {
0070          if(false == detail::check_scale(function, m_scale, result, Policy())) return false;
0071          if(false == detail::check_location(function, m_location, result, Policy())) return false;
0072          return true;
0073    }
0074
0075 private:
0076    RealType m_location;
0077    RealType m_scale;
0078 }; // class laplace_distribution
0079
0080 //
0081 // Convenient type synonym for double.
0082 using laplace = laplace_distribution<double>;
0083
0084 #ifdef __cpp_deduction_guides
0085 template <class RealType>
0086 laplace_distribution(RealType)->laplace_distribution<typename boost::math::tools::promote_args<RealType>::type>;
0087 template <class RealType>
0088 laplace_distribution(RealType,RealType)->laplace_distribution<typename boost::math::tools::promote_args<RealType>::type>;
0089 #endif
0090
0091 //
0092 // Non-member functions.
0093 template <class RealType, class Policy>
0094 inline std::pair<RealType, RealType> range(const laplace_distribution<RealType, Policy>&)
0095 {
0096   if (std::numeric_limits<RealType>::has_infinity)
0097   {  // Can use infinity.
0098      return std::pair<RealType, RealType>(-std::numeric_limits<RealType>::infinity(), std::numeric_limits<RealType>::infinity()); // - to + infinity.
0099   }
0100   else
0101   { // Can only use max_value.
0102     using boost::math::tools::max_value;
0103     return std::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.
0104   }
0105
0106 }
0107
0108 template <class RealType, class Policy>
0109 inline std::pair<RealType, RealType> support(const laplace_distribution<RealType, Policy>&)
0110 {
0111   if (std::numeric_limits<RealType>::has_infinity)
0112   { // Can Use infinity.
0113      return std::pair<RealType, RealType>(-std::numeric_limits<RealType>::infinity(), std::numeric_limits<RealType>::infinity()); // - to + infinity.
0114   }
0115   else
0116   { // Can only use max_value.
0117     using boost::math::tools::max_value;
0118     return std::pair<RealType, RealType>(-max_value<RealType>(), max_value<RealType>()); // - to + max value.
0119   }
0120 }
0121
0122 template <class RealType, class Policy>
0123 inline RealType pdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
0124 {
0125    BOOST_MATH_STD_USING // for ADL of std functions
0126
0127    // Checking function argument
0128    RealType result = 0;
0129    const char* function = "boost::math::pdf(const laplace_distribution<%1%>&, %1%))";
0130
0131    // Check scale and location.
0132    if (false == dist.check_parameters(function, &result)) return result;
0133    // Special pdf values.
0134    if((boost::math::isinf)(x))
0135    {
0136       return 0; // pdf + and - infinity is zero.
0137    }
0138    if (false == detail::check_x(function, x, &result, Policy())) return result;
0139
0140    // General case
0141    RealType scale( dist.scale() );
0142    RealType location( dist.location() );
0143
0144    RealType exponent = x - location;
0145    if (exponent>0) exponent = -exponent;
0146    exponent /= scale;
0147
0148    result = exp(exponent);
0149    result /= 2 * scale;
0150
0151    return result;
0152 } // pdf
0153
0154 template <class RealType, class Policy>
0155 inline RealType logpdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
0156 {
0157    BOOST_MATH_STD_USING // for ADL of std functions
0158
0159    // Checking function argument
0160    RealType result = -std::numeric_limits<RealType>::infinity();
0161    const char* function = "boost::math::logpdf(const laplace_distribution<%1%>&, %1%))";
0162
0163    // Check scale and location.
0164    if (false == dist.check_parameters(function, &result))
0165    {
0166        return result;
0167    }
0168    // Special pdf values.
0169    if((boost::math::isinf)(x))
0170    {
0171       return result; // pdf + and - infinity is zero so logpdf is -INF
0172    }
0173    if (false == detail::check_x(function, x, &result, Policy()))
0174    {
0175        return result;
0176    }
0177
0178    const RealType mu = dist.scale();
0179    const RealType b = dist.location();
0180
0181    // if b is 0 avoid divde by 0 error
0182    if(abs(b) < std::numeric_limits<RealType>::epsilon())
0183    {
0184       result = log(pdf(dist, x));
0185    }
0186    else
0187    {
0188       // General case
0189       const RealType log2 = boost::math::constants::ln_two<RealType>();
0190       result = -abs(x-mu)/b - log(b) - log2;
0191    }
0192
0193    return result;
0194 } // logpdf
0195
0196 template <class RealType, class Policy>
0197 inline RealType cdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
0198 {
0199    BOOST_MATH_STD_USING  // For ADL of std functions.
0200
0201    RealType result = 0;
0202    // Checking function argument.
0203    const char* function = "boost::math::cdf(const laplace_distribution<%1%>&, %1%)";
0204    // Check scale and location.
0205    if (false == dist.check_parameters(function, &result)) return result;
0206
0207    // Special cdf values:
0208    if((boost::math::isinf)(x))
0209    {
0210      if(x < 0) return 0; // -infinity.
0211      return 1; // + infinity.
0212    }
0213    if (false == detail::check_x(function, x, &result, Policy())) return result;
0214
0215    // General cdf  values
0216    RealType scale( dist.scale() );
0217    RealType location( dist.location() );
0218
0219    if (x < location)
0220    {
0221       result = exp( (x-location)/scale )/2;
0222    }
0223    else
0224    {
0225       result = 1 - exp( (location-x)/scale )/2;
0226    }
0227    return result;
0228 } // cdf
0229
0230 template <class RealType, class Policy>
0231 inline RealType logcdf(const laplace_distribution<RealType, Policy>& dist, const RealType& x)
0232 {
0233    BOOST_MATH_STD_USING  // For ADL of std functions.
0234
0235    RealType result = 0;
0236    // Checking function argument.
0237    const char* function = "boost::math::logcdf(const laplace_distribution<%1%>&, %1%)";
0238    // Check scale and location.
0239    if (false == dist.check_parameters(function, &result))
0240    {
0241       return result;
0242    }
0243
0244    // Special cdf values:
0245    if((boost::math::isinf)(x))
0246    {
0247       if(x < 0)
0248       {
0249          return 0; // -infinity.
0250       }
0251       return 1; // + infinity.
0252    }
0253
0254    if (false == detail::check_x(function, x, &result, Policy()))
0255    {
0256       return result;
0257    }
0258
0259    // General cdf  values
0260    RealType scale( dist.scale() );
0261    RealType location( dist.location() );
0262
0263    if (x < location)
0264    {
0265       result = ((x - location) / scale) - boost::math::constants::ln_two<RealType>();
0266    }
0267    else
0268    {
0269       result = log1p(-exp((location - x) / scale) / 2);
0270    }
0271
0272    return result;
0273 } // logcdf
0274
0275 template <class RealType, class Policy>
0276 inline RealType quantile(const laplace_distribution<RealType, Policy>& dist, const RealType& p)
0277 {
0278    BOOST_MATH_STD_USING // for ADL of std functions.
0279
0280    // Checking function argument
0281    RealType result = 0;
0282    const char* function = "boost::math::quantile(const laplace_distribution<%1%>&, %1%)";
0283    if (false == dist.check_parameters(function, &result)) return result;
0284    if(false == detail::check_probability(function, p, &result, Policy())) return result;
0285
0286    // Extreme values of p:
0287    if(p == 0)
0288    {
0289       result = policies::raise_overflow_error<RealType>(function,
0290         "probability parameter is 0, but must be > 0!", Policy());
0291       return -result; // -inf
0292    }
0293
0294    if(p == 1)
0295    {
0296       result = policies::raise_overflow_error<RealType>(function,
0297         "probability parameter is 1, but must be < 1!", Policy());
0298       return result; // inf
0299    }
0300    // Calculate Quantile
0301    RealType scale( dist.scale() );
0302    RealType location( dist.location() );
0303
0304    if (p - 0.5 < 0.0)
0305       result = location + scale*log( static_cast<RealType>(p*2) );
0306    else
0307       result = location - scale*log( static_cast<RealType>(-p*2 + 2) );
0308
0309    return result;
0310 } // quantile
0311
0312
0313 template <class RealType, class Policy>
0314 inline RealType cdf(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
0315 {
0316    // Calculate complement of cdf.
0317    BOOST_MATH_STD_USING // for ADL of std functions
0318
0319    RealType scale = c.dist.scale();
0320    RealType location = c.dist.location();
0321    RealType x = c.param;
0322    RealType result = 0;
0323
0324    // Checking function argument.
0325    const char* function = "boost::math::cdf(const complemented2_type<laplace_distribution<%1%>, %1%>&)";
0326
0327    // Check scale and location.
0328     if (false == c.dist.check_parameters(function, &result)) return result;
0329
0330    // Special cdf values.
0331    if((boost::math::isinf)(x))
0332    {
0333      if(x < 0) return 1; // cdf complement -infinity is unity.
0334      return 0; // cdf complement +infinity is zero.
0335    }
0336    if(false == detail::check_x(function, x, &result, Policy()))return result;
0337
0338    // Cdf interval value.
0339    if (-x < -location)
0340    {
0341       result = exp( (-x+location)/scale )/2;
0342    }
0343    else
0344    {
0345       result = 1 - exp( (-location+x)/scale )/2;
0346    }
0347    return result;
0348 } // cdf complement
0349
0350 template <class RealType, class Policy>
0351 inline RealType logcdf(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
0352 {
0353    // Calculate complement of logcdf.
0354    BOOST_MATH_STD_USING // for ADL of std functions
0355
0356    RealType scale = c.dist.scale();
0357    RealType location = c.dist.location();
0358    RealType x = c.param;
0359    RealType result = 0;
0360
0361    // Checking function argument.
0362    const char* function = "boost::math::logcdf(const complemented2_type<laplace_distribution<%1%>, %1%>&)";
0363
0364    // Check scale and location.
0365     if (false == c.dist.check_parameters(function, &result)) return result;
0366
0367    // Special cdf values.
0368    if((boost::math::isinf)(x))
0369    {
0370      if(x < 0)
0371      {
0372        return 1; // cdf complement -infinity is unity.
0373      }
0374
0375      return 0; // cdf complement +infinity is zero.
0376    }
0377    if(false == detail::check_x(function, x, &result, Policy()))return result;
0378
0379    // Cdf interval value.
0380    if (-x < -location)
0381    {
0382       result = (-x+location)/scale - boost::math::constants::ln_two<RealType>();
0383    }
0384    else
0385    {
0386       result = log1p(-exp( (-location+x)/scale )/2, Policy());
0387    }
0388    return result;
0389 } // cdf complement
0390
0391 template <class RealType, class Policy>
0392 inline RealType quantile(const complemented2_type<laplace_distribution<RealType, Policy>, RealType>& c)
0393 {
0394    BOOST_MATH_STD_USING // for ADL of std functions.
0395
0396    // Calculate quantile.
0397    RealType scale = c.dist.scale();
0398    RealType location = c.dist.location();
0399    RealType q = c.param;
0400    RealType result = 0;
0401
0402    // Checking function argument.
0403    const char* function = "quantile(const complemented2_type<laplace_distribution<%1%>, %1%>&)";
0404    if (false == c.dist.check_parameters(function, &result)) return result;
0405
0406    // Extreme values.
0407    if(q == 0)
0408    {
0409        return std::numeric_limits<RealType>::infinity();
0410    }
0411    if(q == 1)
0412    {
0413        return -std::numeric_limits<RealType>::infinity();
0414    }
0415    if(false == detail::check_probability(function, q, &result, Policy())) return result;
0416
0417    if (0.5 - q < 0.0)
0418       result = location + scale*log( static_cast<RealType>(-q*2 + 2) );
0419    else
0420       result = location - scale*log( static_cast<RealType>(q*2) );
0421
0422
0423    return result;
0424 } // quantile
0425
0426 template <class RealType, class Policy>
0427 inline RealType mean(const laplace_distribution<RealType, Policy>& dist)
0428 {
0429    return dist.location();
0430 }
0431
0432 template <class RealType, class Policy>
0433 inline RealType standard_deviation(const laplace_distribution<RealType, Policy>& dist)
0434 {
0435    return constants::root_two<RealType>() * dist.scale();
0436 }
0437
0438 template <class RealType, class Policy>
0439 inline RealType mode(const laplace_distribution<RealType, Policy>& dist)
0440 {
0441    return dist.location();
0442 }
0443
0444 template <class RealType, class Policy>
0445 inline RealType median(const laplace_distribution<RealType, Policy>& dist)
0446 {
0447    return dist.location();
0448 }
0449
0450 template <class RealType, class Policy>
0451 inline RealType skewness(const laplace_distribution<RealType, Policy>& /*dist*/)
0452 {
0453    return 0;
0454 }
0455
0456 template <class RealType, class Policy>
0457 inline RealType kurtosis(const laplace_distribution<RealType, Policy>& /*dist*/)
0458 {
0459    return 6;
0460 }
0461
0462 template <class RealType, class Policy>
0463 inline RealType kurtosis_excess(const laplace_distribution<RealType, Policy>& /*dist*/)
0464 {
0465    return 3;
0466 }
0467
0468 template <class RealType, class Policy>
0469 inline RealType entropy(const laplace_distribution<RealType, Policy> & dist)
0470 {
0471    using std::log;
0472    return log(2*dist.scale()*constants::e<RealType>());
0473 }
0474
0475 #ifdef _MSC_VER
0476 #  pragma warning(pop)
0477 #endif
0478
0479 } // namespace math
0480 } // namespace boost
0481
0482 // This include must be at the end, *after* the accessors
0483 // for this distribution have been defined, in order to
0484 // keep compilers that support two-phase lookup happy.
0485 #include <boost/math/distributions/detail/derived_accessors.hpp>
0486
0487 #endif // BOOST_STATS_LAPLACE_HPP
0488
0489