shape_of_distribution.hpp

Go to the documentation of this file.

 
#pragma once
 
#include <cmath>
#include <cstddef>
#include <functional>
#include <iterator>
#include <stdexcept>
#include <type_traits>
 
#include "statcpp/basic_statistics.hpp"
 
namespace statcpp {
 
// ============================================================================
// Skewness
// ============================================================================
 
template <typename Iterator>
double population_skewness(Iterator first, Iterator last)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_skewness: empty range");
    }
 
    double m = statcpp::mean(first, last);
    double sum_cubed = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(*it) - m;
        sum_sq += diff * diff;
        sum_cubed += diff * diff * diff;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_skewness: zero variance");
    }
 
    double stddev = std::sqrt(var);
    return (sum_cubed / static_cast<double>(n)) / (stddev * stddev * stddev);
}
 
template <typename Iterator>
double population_skewness(Iterator first, Iterator last, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_skewness: empty range");
    }
 
    double sum_cubed = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(*it) - precomputed_mean;
        sum_sq += diff * diff;
        sum_cubed += diff * diff * diff;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_skewness: zero variance");
    }
 
    double stddev = std::sqrt(var);
    return (sum_cubed / static_cast<double>(n)) / (stddev * stddev * stddev);
}
 
template <typename Iterator, typename Projection,
          typename = std::enable_if_t<
              std::is_invocable_v<Projection,
                  typename std::iterator_traits<Iterator>::value_type>>>
double population_skewness(Iterator first, Iterator last, Projection proj)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_skewness: empty range");
    }
 
    double m = statcpp::mean(first, last, proj);
    double sum_cubed = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(std::invoke(proj, *it)) - m;
        sum_sq += diff * diff;
        sum_cubed += diff * diff * diff;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_skewness: zero variance");
    }
 
    double stddev = std::sqrt(var);
    return (sum_cubed / static_cast<double>(n)) / (stddev * stddev * stddev);
}
 
template <typename Iterator, typename Projection>
double population_skewness(Iterator first, Iterator last, Projection proj, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_skewness: empty range");
    }
 
    double sum_cubed = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(std::invoke(proj, *it)) - precomputed_mean;
        sum_sq += diff * diff;
        sum_cubed += diff * diff * diff;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_skewness: zero variance");
    }
 
    double stddev = std::sqrt(var);
    return (sum_cubed / static_cast<double>(n)) / (stddev * stddev * stddev);
}
 
template <typename Iterator>
double sample_skewness(Iterator first, Iterator last)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 3) {
        throw std::invalid_argument("statcpp::sample_skewness: need at least 3 elements");
    }
 
    double g1 = population_skewness(first, last);
    double correction = std::sqrt(static_cast<double>(n * (n - 1))) / static_cast<double>(n - 2);
    return correction * g1;
}
 
template <typename Iterator>
double sample_skewness(Iterator first, Iterator last, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 3) {
        throw std::invalid_argument("statcpp::sample_skewness: need at least 3 elements");
    }
 
    double g1 = population_skewness(first, last, precomputed_mean);
    double correction = std::sqrt(static_cast<double>(n * (n - 1))) / static_cast<double>(n - 2);
    return correction * g1;
}
 
template <typename Iterator, typename Projection,
          typename = std::enable_if_t<
              std::is_invocable_v<Projection,
                  typename std::iterator_traits<Iterator>::value_type>>>
double sample_skewness(Iterator first, Iterator last, Projection proj)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 3) {
        throw std::invalid_argument("statcpp::sample_skewness: need at least 3 elements");
    }
 
    double g1 = population_skewness(first, last, proj);
    double correction = std::sqrt(static_cast<double>(n * (n - 1))) / static_cast<double>(n - 2);
    return correction * g1;
}
 
template <typename Iterator, typename Projection>
double sample_skewness(Iterator first, Iterator last, Projection proj, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 3) {
        throw std::invalid_argument("statcpp::sample_skewness: need at least 3 elements");
    }
 
    double g1 = population_skewness(first, last, proj, precomputed_mean);
    double correction = std::sqrt(static_cast<double>(n * (n - 1))) / static_cast<double>(n - 2);
    return correction * g1;
}
 
template <typename Iterator>
double skewness(Iterator first, Iterator last)
{
    return sample_skewness(first, last);
}
 
template <typename Iterator>
double skewness(Iterator first, Iterator last, double precomputed_mean)
{
    return sample_skewness(first, last, precomputed_mean);
}
 
template <typename Iterator, typename Projection,
          typename = std::enable_if_t<
              std::is_invocable_v<Projection,
                  typename std::iterator_traits<Iterator>::value_type>>>
double skewness(Iterator first, Iterator last, Projection proj)
{
    return sample_skewness(first, last, proj);
}
 
template <typename Iterator, typename Projection>
double skewness(Iterator first, Iterator last, Projection proj, double precomputed_mean)
{
    return sample_skewness(first, last, proj, precomputed_mean);
}
 
// ============================================================================
// Kurtosis
// ============================================================================
 
template <typename Iterator>
double population_kurtosis(Iterator first, Iterator last)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_kurtosis: empty range");
    }
 
    double m = statcpp::mean(first, last);
    double sum_fourth = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(*it) - m;
        double diff_sq = diff * diff;
        sum_sq += diff_sq;
        sum_fourth += diff_sq * diff_sq;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_kurtosis: zero variance");
    }
 
    double fourth_moment = sum_fourth / static_cast<double>(n);
    return (fourth_moment / (var * var)) - 3.0;
}
 
template <typename Iterator>
double population_kurtosis(Iterator first, Iterator last, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_kurtosis: empty range");
    }
 
    double sum_fourth = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(*it) - precomputed_mean;
        double diff_sq = diff * diff;
        sum_sq += diff_sq;
        sum_fourth += diff_sq * diff_sq;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_kurtosis: zero variance");
    }
 
    double fourth_moment = sum_fourth / static_cast<double>(n);
    return (fourth_moment / (var * var)) - 3.0;
}
 
template <typename Iterator, typename Projection,
          typename = std::enable_if_t<
              std::is_invocable_v<Projection,
                  typename std::iterator_traits<Iterator>::value_type>>>
double population_kurtosis(Iterator first, Iterator last, Projection proj)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_kurtosis: empty range");
    }
 
    double m = statcpp::mean(first, last, proj);
    double sum_fourth = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(std::invoke(proj, *it)) - m;
        double diff_sq = diff * diff;
        sum_sq += diff_sq;
        sum_fourth += diff_sq * diff_sq;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_kurtosis: zero variance");
    }
 
    double fourth_moment = sum_fourth / static_cast<double>(n);
    return (fourth_moment / (var * var)) - 3.0;
}
 
template <typename Iterator, typename Projection>
double population_kurtosis(Iterator first, Iterator last, Projection proj, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n == 0) {
        throw std::invalid_argument("statcpp::population_kurtosis: empty range");
    }
 
    double sum_fourth = 0.0;
    double sum_sq = 0.0;
 
    for (auto it = first; it != last; ++it) {
        double diff = static_cast<double>(std::invoke(proj, *it)) - precomputed_mean;
        double diff_sq = diff * diff;
        sum_sq += diff_sq;
        sum_fourth += diff_sq * diff_sq;
    }
 
    double var = sum_sq / static_cast<double>(n);
    if (var == 0.0) {
        throw std::invalid_argument("statcpp::population_kurtosis: zero variance");
    }
 
    double fourth_moment = sum_fourth / static_cast<double>(n);
    return (fourth_moment / (var * var)) - 3.0;
}
 
template <typename Iterator>
double sample_kurtosis(Iterator first, Iterator last)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 4) {
        throw std::invalid_argument("statcpp::sample_kurtosis: need at least 4 elements");
    }
 
    double g2 = population_kurtosis(first, last);
    double nd = static_cast<double>(n);
    double correction = ((nd + 1.0) * g2 + 6.0) * (nd - 1.0) / ((nd - 2.0) * (nd - 3.0));
    return correction;
}
 
template <typename Iterator>
double sample_kurtosis(Iterator first, Iterator last, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 4) {
        throw std::invalid_argument("statcpp::sample_kurtosis: need at least 4 elements");
    }
 
    double g2 = population_kurtosis(first, last, precomputed_mean);
    double nd = static_cast<double>(n);
    double correction = ((nd + 1.0) * g2 + 6.0) * (nd - 1.0) / ((nd - 2.0) * (nd - 3.0));
    return correction;
}
 
template <typename Iterator, typename Projection,
          typename = std::enable_if_t<
              std::is_invocable_v<Projection,
                  typename std::iterator_traits<Iterator>::value_type>>>
double sample_kurtosis(Iterator first, Iterator last, Projection proj)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 4) {
        throw std::invalid_argument("statcpp::sample_kurtosis: need at least 4 elements");
    }
 
    double g2 = population_kurtosis(first, last, proj);
    double nd = static_cast<double>(n);
    double correction = ((nd + 1.0) * g2 + 6.0) * (nd - 1.0) / ((nd - 2.0) * (nd - 3.0));
    return correction;
}
 
template <typename Iterator, typename Projection>
double sample_kurtosis(Iterator first, Iterator last, Projection proj, double precomputed_mean)
{
    auto n = static_cast<std::size_t>(std::distance(first, last));
    if (n < 4) {
        throw std::invalid_argument("statcpp::sample_kurtosis: need at least 4 elements");
    }
 
    double g2 = population_kurtosis(first, last, proj, precomputed_mean);
    double nd = static_cast<double>(n);
    double correction = ((nd + 1.0) * g2 + 6.0) * (nd - 1.0) / ((nd - 2.0) * (nd - 3.0));
    return correction;
}
 
template <typename Iterator>
double kurtosis(Iterator first, Iterator last)
{
    return sample_kurtosis(first, last);
}
 
template <typename Iterator>
double kurtosis(Iterator first, Iterator last, double precomputed_mean)
{
    return sample_kurtosis(first, last, precomputed_mean);
}
 
template <typename Iterator, typename Projection,
          typename = std::enable_if_t<
              std::is_invocable_v<Projection,
                  typename std::iterator_traits<Iterator>::value_type>>>
double kurtosis(Iterator first, Iterator last, Projection proj)
{
    return sample_kurtosis(first, last, proj);
}
 
template <typename Iterator, typename Projection>
double kurtosis(Iterator first, Iterator last, Projection proj, double precomputed_mean)
{
    return sample_kurtosis(first, last, proj, precomputed_mean);
}
 
} // namespace statcpp