F-Distribution

Table of contents


Density Function

The density function of the F distribution:

\[f(x; d_1, d_2) = \frac{1}{\mathcal{B} \left( \frac{d_1}{2}, \frac{d_2}{2} \right)} \left( \frac{d_1}{d_2} \right)^{\frac{d_1}{2}} x^{d_1/2 - 1} \left(1 + \frac{d_1}{d_2} x \right)^{- \frac{d_1+d_2}{2}} \times \mathbf{1} [x \geq 0]\]

where \(\mathcal{B}(a,b)\) denotes the Beta function.

Methods for scalar input, as well as for vector/matrix input, are listed below.

Scalar Input

template<typename T1, typename T2, typename T3>
constexpr common_return_t<T1, T2, T3> df(const T1 x, const T2 df1_par, const T3 df2_par, const bool log_form = false) noexcept

Density function of the F-distribution.

Example:

stats::df(1.5,10.0,12.0,false); 

Parameters
  • x – a real-valued input.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-density or the true form.

Returns

the density function evaluated at x.

Vector/Matrix Input

STL Containers

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>>
inline std::vector<rT> df(const std::vector<eT> &x, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Density function of the F-distribution.

Example:

std::vector<double> x = {0.3, 0.5, 0.9};
stats::df(x,3.0,2.0,false);

Parameters
  • x – a standard vector.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-density or the true form.

Returns

a vector of density function values corresponding to the elements of x.

Armadillo

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>>
inline ArmaMat<rT> df(const ArmaMat<eT> &X, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Density function of the F-distribution.

Example:

arma::mat X = { {0.2,  0.7,  0.1},
                {0.9, -0.3,  1.3} };
stats::df(X,3.0,2.0,false);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-density or the true form.

Returns

a matrix of density function values corresponding to the elements of X.

Blaze

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>, bool To = blaze::columnMajor>
inline BlazeMat<rT, To> df(const BlazeMat<eT, To> &X, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Density function of the F-distribution.

Example:

stats::df(X,3.0,2.0,false);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-density or the true form.

Returns

a matrix of density function values corresponding to the elements of X.

Eigen

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>, int iTr = Eigen::Dynamic, int iTc = Eigen::Dynamic>
inline EigenMat<rT, iTr, iTc> df(const EigenMat<eT, iTr, iTc> &X, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Density function of the F-distribution.

Example:

stats::df(X,3.0,2.0,false);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-density or the true form.

Returns

a matrix of density function values corresponding to the elements of X.


Cumulative Distribution Function

The cumulative distribution function of the F distribution:

\[F(x; d_1, d_2) = \int_0^x f(z; d_1, d_2) dz = I_{\frac{d_1 x}{d_2 + d_1 x}} (d_1 / 2, d_2 / 2)\]

where \(I_x (a,b)\) denotes the regularized incomplete Beta function.

Methods for scalar input, as well as for vector/matrix input, are listed below.

Scalar Input

template<typename T1, typename T2, typename T3>
constexpr common_return_t<T1, T2, T3> pf(const T1 x, const T2 df1_par, const T3 df2_par, const bool log_form = false) noexcept

Distribution function of the F-distribution.

Example:

stats::pf(1.5,10.0,12.0,false); 

Parameters
  • x – a real-valued input.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-probability or the true form.

Returns

the cumulative distribution function evaluated at x.

Vector/Matrix Input

STL Containers

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>>
inline std::vector<rT> pf(const std::vector<eT> &x, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Distribution function of the Beta distribution.

Example:

std::vector<double> x = {0.3, 0.5, 0.9};
stats::pf(x,3.0,2.0,false);

Parameters
  • x – a standard vector.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-probability or the true form.

Returns

a vector of CDF values corresponding to the elements of x.

Armadillo

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>>
inline ArmaMat<rT> pf(const ArmaMat<eT> &X, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Distribution function of the Beta distribution.

Example:

arma::mat X = { {0.2,  0.7,  0.1},
                {0.9, -0.3,  1.3} };
stats::pf(X,3.0,2.0,false);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-probability or the true form.

Returns

a matrix of CDF values corresponding to the elements of X.

Blaze

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>, bool To = blaze::columnMajor>
inline BlazeMat<rT, To> pf(const BlazeMat<eT, To> &X, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Distribution function of the Beta distribution.

Example:

stats::pf(X,3.0,2.0,false);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-probability or the true form.

Returns

a matrix of CDF values corresponding to the elements of X.

Eigen

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>, int iTr = Eigen::Dynamic, int iTc = Eigen::Dynamic>
inline EigenMat<rT, iTr, iTc> pf(const EigenMat<eT, iTr, iTc> &X, const T1 df1_par, const T2 df2_par, const bool log_form = false)

Distribution function of the Beta distribution.

Example:

stats::pf(X,3.0,2.0,false);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • log_form – return the log-probability or the true form.

Returns

a matrix of CDF values corresponding to the elements of X.


Quantile Function

The quantile function of the F distribution:

\[q(p; a,b) = \inf \left\{ x : p \leq I_{\frac{d_1 x}{d_2 + d_1 x}} (d_1 / 2, d_2 / 2) \right\}\]

Methods for scalar input, as well as for vector/matrix input, are listed below.

Scalar Input

template<typename T1, typename T2, typename T3>
constexpr common_return_t<T1, T2, T3> qf(const T1 p, const T2 df1_par, const T3 df2_par) noexcept

Quantile function of the F-distribution.

Example:

stats::qf(0.5,10.0,12.0); 

Parameters
  • p – a real-valued input.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

Returns

the quantile function evaluated at p.

Vector/Matrix Input

STL Containers

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>>
inline std::vector<rT> qf(const std::vector<eT> &x, const T1 df1_par, const T2 df2_par)

Quantile function of the F-distribution.

Example:

std::vector<double> x = {0.3, 0.5, 0.9};
stats::qf(x,3.0,2.0);

Parameters
  • x – a standard vector.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

Returns

a vector of quantile values corresponding to the elements of x.

Armadillo

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>>
inline ArmaMat<rT> qf(const ArmaMat<eT> &X, const T1 df1_par, const T2 df2_par)

Quantile function of the F-distribution.

Example:

arma::mat X = { {0.2,  0.7,  0.1},
                {0.9,  0.3,  0.87} };
stats::qf(X,3.0,2.0);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

Returns

a matrix of quantile values corresponding to the elements of X.

Blaze

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>, bool To = blaze::columnMajor>
inline BlazeMat<rT, To> qf(const BlazeMat<eT, To> &X, const T1 df1_par, const T2 df2_par)

Quantile function of the F-distribution.

Example:

stats::qf(X,3.0,2.0);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

Returns

a matrix of quantile values corresponding to the elements of X.

Eigen

template<typename eT, typename T1, typename T2, typename rT = common_return_t<eT, T1, T2>, int iTr = Eigen::Dynamic, int iTc = Eigen::Dynamic>
inline EigenMat<rT, iTr, iTc> qf(const EigenMat<eT, iTr, iTc> &X, const T1 df1_par, const T2 df2_par)

Quantile function of the F-distribution.

Example:

stats::qf(X,3.0,2.0);

Parameters
  • X – a matrix of input values.

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

Returns

a matrix of quantile values corresponding to the elements of X.


Random Sampling

Random sampling for the Beta distribution is achieved by simulating two independent \(\chi^2\)-distributed random variables, \(X \sim \chi^2(d_1), Y \sim \chi^2(d_2)\), then returning:

\[Z = \frac{d_1}{d_2} \frac{X}{Y}\]

Scalar Output

  1. Random number engines

template<typename T1, typename T2>
inline common_return_t<T1, T2> rf(const T1 df1_par, const T2 df2_par, rand_engine_t &engine)

Random sampling function for the F-distribution.

Example:

stats::rand_engine_t engine(1776);
stats::rf(3.0,2.0,engine);

Parameters
  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • engine – a random engine, passed by reference.

Returns

a pseudo-random draw from the F-distribution.

  1. Seed values

template<typename T1, typename T2>
inline common_return_t<T1, T2> rf(const T1 df1_par, const T2 df2_par, const ullint_t seed_val = std::random_device{}())

Random sampling function for the F-distribution.

Example:

stats::rf(3.0,2.0,1776);

Parameters
  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • seed_val – initialize the random engine with a non-negative integral-valued seed.

Returns

a pseudo-random draw from the F-distribution.

Vector/Matrix Output

  1. Random number engines

template<typename mT, typename T1, typename T2>
inline mT rf(const ullint_t n, const ullint_t k, const T1 df1_par, const T2 df2_par, rand_engine_t &engine)

Random matrix sampling function for the F-distribution.

Example:

stats::rand_engine_t engine(1776);
// std::vector
stats::rf<std::vector<double>>(5,4,3.0,2.0,engine);
// Armadillo matrix
stats::rf<arma::mat>(5,4,3.0,2.0,engine);
// Blaze dynamic matrix
stats::rf<blaze::DynamicMatrix<double,blaze::columnMajor>>(5,4,3.0,2.0,engine);
// Eigen dynamic matrix
stats::rf<Eigen::MatrixXd>(5,4,3.0,2.0,engine);

Note

This function requires template instantiation; acceptable output types include: std::vector, with element type float, double, etc., as well as Armadillo, Blaze, and Eigen dense matrices.

Parameters
  • n – the number of output rows

  • k – the number of output columns

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • engine – a random engine, passed by reference.

Returns

a matrix of pseudo-random draws from the F-distribution.

  1. Seed values

template<typename mT, typename T1, typename T2>
inline mT rf(const ullint_t n, const ullint_t k, const T1 df1_par, const T2 df2_par, const ullint_t seed_val = std::random_device{}())

Random matrix sampling function for the F-distribution.

Example:

// std::vector
stats::rf<std::vector<double>>(5,4,3.0,2.0);
// Armadillo matrix
stats::rf<arma::mat>(5,4,3.0,2.0);
// Blaze dynamic matrix
stats::rf<blaze::DynamicMatrix<double,blaze::columnMajor>>(5,4,3.0,2.0);
// Eigen dynamic matrix
stats::rf<Eigen::MatrixXd>(5,4,3.0,2.0);

Note

This function requires template instantiation; acceptable output types include: std::vector, with element type float, double, etc., as well as Armadillo, Blaze, and Eigen dense matrices.

Parameters
  • n – the number of output rows

  • k – the number of output columns

  • df1_par – a degrees of freedom parameter, a real-valued input.

  • df2_par – a degrees of freedom parameter, a real-valued input.

  • seed_val – initialize the random engine with a non-negative integral-valued seed.

Returns

a matrix of pseudo-random draws from the F-distribution.