Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2025-09-18 08:48:30 

0001 // boost\math\distributions\non_central_f.hpp
0002 
0003 // Copyright John Maddock 2008.
0004 // Copyright Matt Borland 2024.
0005 // Use, modification and distribution are subject to the
0006 // Boost Software License, Version 1.0.
0007 // (See accompanying file LICENSE_1_0.txt
0008 // or copy at http://www.boost.org/LICENSE_1_0.txt)
0009 
0010 #ifndef BOOST_MATH_SPECIAL_NON_CENTRAL_F_HPP
0011 #define BOOST_MATH_SPECIAL_NON_CENTRAL_F_HPP
0012 
0013 #include <boost/math/tools/config.hpp>
0014 #include <boost/math/tools/tuple.hpp>
0015 #include <boost/math/tools/promotion.hpp>
0016 #include <boost/math/distributions/non_central_beta.hpp>
0017 #include <boost/math/distributions/detail/generic_mode.hpp>
0018 #include <boost/math/special_functions/pow.hpp>
0019 #include <boost/math/policies/policy.hpp>
0020 
0021 namespace boost
0022 {
0023    namespace math
0024    {
0025       template <class RealType = double, class Policy = policies::policy<> >
0026       class non_central_f_distribution
0027       {
0028       public:
0029          typedef RealType value_type;
0030          typedef Policy policy_type;
0031 
0032          BOOST_MATH_GPU_ENABLED non_central_f_distribution(RealType v1_, RealType v2_, RealType lambda) : v1(v1_), v2(v2_), ncp(lambda)
0033          {
0034             constexpr auto function = "boost::math::non_central_f_distribution<%1%>::non_central_f_distribution(%1%,%1%)";
0035             RealType r;
0036             detail::check_df(
0037                function,
0038                v1, &r, Policy());
0039             detail::check_df(
0040                function,
0041                v2, &r, Policy());
0042             detail::check_non_centrality(
0043                function,
0044                lambda,
0045                &r,
0046                Policy());
0047          } // non_central_f_distribution constructor.
0048 
0049          BOOST_MATH_GPU_ENABLED RealType degrees_of_freedom1()const
0050          {
0051             return v1;
0052          }
0053          BOOST_MATH_GPU_ENABLED RealType degrees_of_freedom2()const
0054          {
0055             return v2;
0056          }
0057          BOOST_MATH_GPU_ENABLED RealType non_centrality() const
0058          { // Private data getter function.
0059             return ncp;
0060          }
0061       private:
0062          // Data member, initialized by constructor.
0063          RealType v1;   // alpha.
0064          RealType v2;   // beta.
0065          RealType ncp; // non-centrality parameter
0066       }; // template <class RealType, class Policy> class non_central_f_distribution
0067 
0068       typedef non_central_f_distribution<double> non_central_f; // Reserved name of type double.
0069 
0070       #ifdef __cpp_deduction_guides
0071       template <class RealType>
0072       non_central_f_distribution(RealType,RealType,RealType)->non_central_f_distribution<typename boost::math::tools::promote_args<RealType>::type>;
0073       #endif
0074 
0075       // Non-member functions to give properties of the distribution.
0076 
0077       template <class RealType, class Policy>
0078       BOOST_MATH_GPU_ENABLED inline const boost::math::pair<RealType, RealType> range(const non_central_f_distribution<RealType, Policy>& /* dist */)
0079       { // Range of permissible values for random variable k.
0080          using boost::math::tools::max_value;
0081          return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>());
0082       }
0083 
0084       template <class RealType, class Policy>
0085       BOOST_MATH_GPU_ENABLED inline const boost::math::pair<RealType, RealType> support(const non_central_f_distribution<RealType, Policy>& /* dist */)
0086       { // Range of supported values for random variable k.
0087          // This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
0088          using boost::math::tools::max_value;
0089          return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>());
0090       }
0091 
0092       template <class RealType, class Policy>
0093       BOOST_MATH_GPU_ENABLED inline RealType mean(const non_central_f_distribution<RealType, Policy>& dist)
0094       {
0095          constexpr auto function = "mean(non_central_f_distribution<%1%> const&)";
0096          RealType v1 = dist.degrees_of_freedom1();
0097          RealType v2 = dist.degrees_of_freedom2();
0098          RealType l = dist.non_centrality();
0099          RealType r;
0100          if(!detail::check_df(
0101             function,
0102             v1, &r, Policy())
0103                ||
0104             !detail::check_df(
0105                function,
0106                v2, &r, Policy())
0107                ||
0108             !detail::check_non_centrality(
0109                function,
0110                l,
0111                &r,
0112                Policy()))
0113                   return r;
0114          if(v2 <= 2)
0115             return policies::raise_domain_error(
0116                function,
0117                "Second degrees of freedom parameter was %1%, but must be > 2 !",
0118                v2, Policy());
0119          return v2 * (v1 + l) / (v1 * (v2 - 2));
0120       } // mean
0121 
0122       template <class RealType, class Policy>
0123       BOOST_MATH_GPU_ENABLED inline RealType mode(const non_central_f_distribution<RealType, Policy>& dist)
0124       { // mode.
0125          constexpr auto function = "mode(non_central_chi_squared_distribution<%1%> const&)";
0126 
0127          RealType n = dist.degrees_of_freedom1();
0128          RealType m = dist.degrees_of_freedom2();
0129          RealType l = dist.non_centrality();
0130          RealType r;
0131          if(!detail::check_df(
0132             function,
0133             n, &r, Policy())
0134                ||
0135             !detail::check_df(
0136                function,
0137                m, &r, Policy())
0138                ||
0139             !detail::check_non_centrality(
0140                function,
0141                l,
0142                &r,
0143                Policy()))
0144                   return r;
0145          RealType guess = m > 2 ? RealType(m * (n + l) / (n * (m - 2))) : RealType(1);
0146          return detail::generic_find_mode(
0147             dist,
0148             guess,
0149             function);
0150       }
0151 
0152       template <class RealType, class Policy>
0153       BOOST_MATH_GPU_ENABLED inline RealType variance(const non_central_f_distribution<RealType, Policy>& dist)
0154       { // variance.
0155          constexpr auto function = "variance(non_central_f_distribution<%1%> const&)";
0156          RealType n = dist.degrees_of_freedom1();
0157          RealType m = dist.degrees_of_freedom2();
0158          RealType l = dist.non_centrality();
0159          RealType r;
0160          if(!detail::check_df(
0161             function,
0162             n, &r, Policy())
0163                ||
0164             !detail::check_df(
0165                function,
0166                m, &r, Policy())
0167                ||
0168             !detail::check_non_centrality(
0169                function,
0170                l,
0171                &r,
0172                Policy()))
0173                   return r;
0174          if(m <= 4)
0175             return policies::raise_domain_error(
0176                function,
0177                "Second degrees of freedom parameter was %1%, but must be > 4 !",
0178                m, Policy());
0179          RealType result = 2 * m * m * ((n + l) * (n + l)
0180             + (m - 2) * (n + 2 * l));
0181          result /= (m - 4) * (m - 2) * (m - 2) * n * n;
0182          return result;
0183       }
0184 
0185       // RealType standard_deviation(const non_central_f_distribution<RealType, Policy>& dist)
0186       // standard_deviation provided by derived accessors.
0187 
0188       template <class RealType, class Policy>
0189       BOOST_MATH_GPU_ENABLED inline RealType skewness(const non_central_f_distribution<RealType, Policy>& dist)
0190       { // skewness = sqrt(l).
0191          constexpr auto function = "skewness(non_central_f_distribution<%1%> const&)";
0192          BOOST_MATH_STD_USING
0193          RealType n = dist.degrees_of_freedom1();
0194          RealType m = dist.degrees_of_freedom2();
0195          RealType l = dist.non_centrality();
0196          RealType r;
0197          if(!detail::check_df(
0198             function,
0199             n, &r, Policy())
0200                ||
0201             !detail::check_df(
0202                function,
0203                m, &r, Policy())
0204                ||
0205             !detail::check_non_centrality(
0206                function,
0207                l,
0208                &r,
0209                Policy()))
0210                   return r;
0211          if(m <= 6)
0212             return policies::raise_domain_error(
0213                function,
0214                "Second degrees of freedom parameter was %1%, but must be > 6 !",
0215                m, Policy());
0216          RealType result = 2 * constants::root_two<RealType>();
0217          result *= sqrt(m - 4);
0218          result *= (n * (m + n - 2) *(m + 2 * n - 2)
0219             + 3 * (m + n - 2) * (m + 2 * n - 2) * l
0220             + 6 * (m + n - 2) * l * l + 2 * l * l * l);
0221          result /= (m - 6) * pow(n * (m + n - 2) + 2 * (m + n - 2) * l + l * l, RealType(1.5f));
0222          return result;
0223       }
0224 
0225       template <class RealType, class Policy>
0226       BOOST_MATH_GPU_ENABLED inline RealType kurtosis_excess(const non_central_f_distribution<RealType, Policy>& dist)
0227       {
0228          constexpr auto function = "kurtosis_excess(non_central_f_distribution<%1%> const&)";
0229          BOOST_MATH_STD_USING
0230          RealType n = dist.degrees_of_freedom1();
0231          RealType m = dist.degrees_of_freedom2();
0232          RealType l = dist.non_centrality();
0233          RealType r;
0234          if(!detail::check_df(
0235             function,
0236             n, &r, Policy())
0237                ||
0238             !detail::check_df(
0239                function,
0240                m, &r, Policy())
0241                ||
0242             !detail::check_non_centrality(
0243                function,
0244                l,
0245                &r,
0246                Policy()))
0247                   return r;
0248          if(m <= 8)
0249             return policies::raise_domain_error(
0250                function,
0251                "Second degrees of freedom parameter was %1%, but must be > 8 !",
0252                m, Policy());
0253          RealType l2 = l * l;
0254          RealType l3 = l2 * l;
0255          RealType l4 = l2 * l2;
0256          RealType result = (3 * (m - 4) * (n * (m + n - 2)
0257             * (4 * (m - 2) * (m - 2)
0258             + (m - 2) * (m + 10) * n
0259             + (10 + m) * n * n)
0260             + 4 * (m + n - 2) * (4 * (m - 2) * (m - 2)
0261             + (m - 2) * (10 + m) * n
0262             + (10 + m) * n * n) * l + 2 * (10 + m)
0263             * (m + n - 2) * (2 * m + 3 * n - 4) * l2
0264             + 4 * (10 + m) * (-2 + m + n) * l3
0265             + (10 + m) * l4))
0266             /
0267             ((-8 + m) * (-6 + m) * boost::math::pow<2>(n * (-2 + m + n)
0268             + 2 * (-2 + m + n) * l + l2));
0269             return result;
0270       } // kurtosis_excess
0271 
0272       template <class RealType, class Policy>
0273       BOOST_MATH_GPU_ENABLED inline RealType kurtosis(const non_central_f_distribution<RealType, Policy>& dist)
0274       {
0275          return kurtosis_excess(dist) + 3;
0276       }
0277 
0278       template <class RealType, class Policy>
0279       BOOST_MATH_GPU_ENABLED inline RealType pdf(const non_central_f_distribution<RealType, Policy>& dist, const RealType& x)
0280       { // Probability Density/Mass Function.
0281          typedef typename policies::evaluation<RealType, Policy>::type value_type;
0282          typedef typename policies::normalise<
0283             Policy,
0284             policies::promote_float<false>,
0285             policies::promote_double<false>,
0286             policies::discrete_quantile<>,
0287             policies::assert_undefined<> >::type forwarding_policy;
0288 
0289          value_type alpha = dist.degrees_of_freedom1() / 2;
0290          value_type beta = dist.degrees_of_freedom2() / 2;
0291          value_type y = x * alpha / beta;
0292          value_type r = pdf(boost::math::non_central_beta_distribution<value_type, forwarding_policy>(alpha, beta, dist.non_centrality()), y / (1 + y));
0293          return policies::checked_narrowing_cast<RealType, forwarding_policy>(
0294             r * (dist.degrees_of_freedom1() / dist.degrees_of_freedom2()) / ((1 + y) * (1 + y)),
0295             "pdf(non_central_f_distribution<%1%>, %1%)");
0296       } // pdf
0297 
0298       template <class RealType, class Policy>
0299       BOOST_MATH_GPU_ENABLED RealType cdf(const non_central_f_distribution<RealType, Policy>& dist, const RealType& x)
0300       {
0301          constexpr auto function = "cdf(const non_central_f_distribution<%1%>&, %1%)";
0302          RealType r;
0303          if(!detail::check_df(
0304             function,
0305             dist.degrees_of_freedom1(), &r, Policy())
0306                ||
0307             !detail::check_df(
0308                function,
0309                dist.degrees_of_freedom2(), &r, Policy())
0310                ||
0311             !detail::check_non_centrality(
0312                function,
0313                dist.non_centrality(),
0314                &r,
0315                Policy()))
0316                   return r;
0317 
0318          if((x < 0) || !(boost::math::isfinite)(x))
0319          {
0320             return policies::raise_domain_error<RealType>(
0321                function, "Random Variable parameter was %1%, but must be > 0 !", x, Policy());
0322          }
0323 
0324          RealType alpha = dist.degrees_of_freedom1() / 2;
0325          RealType beta = dist.degrees_of_freedom2() / 2;
0326          RealType y = x * alpha / beta;
0327          RealType c = y / (1 + y);
0328          RealType cp = 1 / (1 + y);
0329          //
0330          // To ensure accuracy, we pass both x and 1-x to the
0331          // non-central beta cdf routine, this ensures accuracy
0332          // even when we compute x to be ~ 1:
0333          //
0334          r = detail::non_central_beta_cdf(c, cp, alpha, beta,
0335             dist.non_centrality(), false, Policy());
0336          return r;
0337       } // cdf
0338 
0339       template <class RealType, class Policy>
0340       BOOST_MATH_GPU_ENABLED RealType cdf(const complemented2_type<non_central_f_distribution<RealType, Policy>, RealType>& c)
0341       { // Complemented Cumulative Distribution Function
0342          constexpr auto function = "cdf(complement(const non_central_f_distribution<%1%>&, %1%))";
0343          RealType r;
0344          if(!detail::check_df(
0345             function,
0346             c.dist.degrees_of_freedom1(), &r, Policy())
0347                ||
0348             !detail::check_df(
0349                function,
0350                c.dist.degrees_of_freedom2(), &r, Policy())
0351                ||
0352             !detail::check_non_centrality(
0353                function,
0354                c.dist.non_centrality(),
0355                &r,
0356                Policy()))
0357                   return r;
0358 
0359          if((c.param < 0) || !(boost::math::isfinite)(c.param))
0360          {
0361             return policies::raise_domain_error<RealType>(
0362                function, "Random Variable parameter was %1%, but must be > 0 !", c.param, Policy());
0363          }
0364 
0365          RealType alpha = c.dist.degrees_of_freedom1() / 2;
0366          RealType beta = c.dist.degrees_of_freedom2() / 2;
0367          RealType y = c.param * alpha / beta;
0368          RealType x = y / (1 + y);
0369          RealType cx = 1 / (1 + y);
0370          //
0371          // To ensure accuracy, we pass both x and 1-x to the
0372          // non-central beta cdf routine, this ensures accuracy
0373          // even when we compute x to be ~ 1:
0374          //
0375          r = detail::non_central_beta_cdf(x, cx, alpha, beta,
0376             c.dist.non_centrality(), true, Policy());
0377          return r;
0378       } // ccdf
0379 
0380       template <class RealType, class Policy>
0381       BOOST_MATH_GPU_ENABLED inline RealType quantile(const non_central_f_distribution<RealType, Policy>& dist, const RealType& p)
0382       { // Quantile (or Percent Point) function.
0383          RealType alpha = dist.degrees_of_freedom1() / 2;
0384          RealType beta = dist.degrees_of_freedom2() / 2;
0385          RealType x = quantile(boost::math::non_central_beta_distribution<RealType, Policy>(alpha, beta, dist.non_centrality()), p);
0386          if(x == 1)
0387             return policies::raise_overflow_error<RealType>(
0388                "quantile(const non_central_f_distribution<%1%>&, %1%)",
0389                "Result of non central F quantile is too large to represent.",
0390                Policy());
0391          return (x / (1 - x)) * (dist.degrees_of_freedom2() / dist.degrees_of_freedom1());
0392       } // quantile
0393 
0394       template <class RealType, class Policy>
0395       BOOST_MATH_GPU_ENABLED inline RealType quantile(const complemented2_type<non_central_f_distribution<RealType, Policy>, RealType>& c)
0396       { // Quantile (or Percent Point) function.
0397          RealType alpha = c.dist.degrees_of_freedom1() / 2;
0398          RealType beta = c.dist.degrees_of_freedom2() / 2;
0399          RealType x = quantile(complement(boost::math::non_central_beta_distribution<RealType, Policy>(alpha, beta, c.dist.non_centrality()), c.param));
0400          if(x == 1)
0401             return policies::raise_overflow_error<RealType>(
0402                "quantile(complement(const non_central_f_distribution<%1%>&, %1%))",
0403                "Result of non central F quantile is too large to represent.",
0404                Policy());
0405          return (x / (1 - x)) * (c.dist.degrees_of_freedom2() / c.dist.degrees_of_freedom1());
0406       } // quantile complement.
0407 
0408    } // namespace math
0409 } // namespace boost
0410 
0411 // This include must be at the end, *after* the accessors
0412 // for this distribution have been defined, in order to
0413 // keep compilers that support two-phase lookup happy.
0414 #include <boost/math/distributions/detail/derived_accessors.hpp>
0415 
0416 #endif // BOOST_MATH_SPECIAL_NON_CENTRAL_F_HPP
0417 
0418 
0419
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!