Code Coverage
 
Classes and Traits
Functions and Methods
Lines
Total
0.00% covered (danger)
0.00%
0 / 1
93.33% covered (success)
93.33%
14 / 15
CRAP
94.12% covered (success)
94.12%
96 / 102
Stats
0.00% covered (danger)
0.00%
0 / 1
93.33% covered (success)
93.33%
14 / 15
31.20
94.12% covered (success)
94.12%
96 / 102
 __construct
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
6 / 6
 pushArray
100.00% covered (success)
100.00%
1 / 1
2
100.00% covered (success)
100.00%
3 / 3
 push
0.00% covered (danger)
0.00%
0 / 1
10.18
87.76% covered (success)
87.76%
43 / 49
 count
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
1 / 1
 min
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
1 / 1
 avg
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
1 / 1
 max
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
1 / 1
 sum
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
1 / 1
 variance
100.00% covered (success)
100.00%
1 / 1
2
100.00% covered (success)
100.00%
3 / 3
 stdDev
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
1 / 1
 skewness
100.00% covered (success)
100.00%
1 / 1
2
100.00% covered (success)
100.00%
7 / 7
 kurtosis
100.00% covered (success)
100.00%
1 / 1
2
100.00% covered (success)
100.00%
9 / 9
 cv
100.00% covered (success)
100.00%
1 / 1
2
100.00% covered (success)
100.00%
3 / 3
 stats
100.00% covered (success)
100.00%
1 / 1
1
100.00% covered (success)
100.00%
10 / 10
 bcmulArray
100.00% covered (success)
100.00%
1 / 1
3
100.00% covered (success)
100.00%
6 / 6
<?php
namespace ia\Math;
/**
 * Calculate Running basic statistics, using bcMath: min, avg, max, sum, var, stdDev, skewness, kurtosis, coefficient of variation (cv)
 *
 *  based on https://www.johndcook.com/blog/skewness_kurtosis, thanks
 *  recomended: include iaBcMathShim
 *
 * @author Raul Santos
 * @license MIT
 * @copyright 2017
 * @version 1.1.3
 */
class Stats {
    const MISSING_ZERO = 3;
    const MISSING_SKIP = 2;
    const NA_ZERO = 5;
    const NA_SKIP = 4;
    const STAT_MIN = 'min';
    const STAT_AVG = 'avg';
    const STAT_MAX = 'max';
    const STAT_COUNT = 'count';
    const STAT_SUM = 'sum';
    const STAT_VARIANCE = 'variance';
    const STAT_STDDEV = 'stdDev';
    const STAT_CV = 'cv';
    const STAT_SKEWNESS = 'skewness';
    const STAT_KURTOSIS = 'kurtosis';
    protected $exceptions; // Stats::MISSING_*, Stats::NA_*
    protected $scale = 8;
    protected $n = 0;
    protected $min = null;
    protected $max = null;
    protected $sum = 0.00;
    protected $M1 = 0.00;
    protected $M2 = 0.00;
    protected $M3 = 0.00;
    protected $M4 = 0.00;
    /**
     * Sets decimals to use and wheather to use bcmath if avaliable
     *
     * @param integer $scale number of decimals to calculate, for bcmath
     * @return void
     */
    /**
     * Constructor
     *
     * @param integer $scale number of decimals
     * @param integer $exceptions  How to manage missing, null, N/A & not numeric. Stats::MISSING_ Stats::NA_ default Asume 0.00
     * @return void
     */
    public function __construct($scale = 8, $exceptions = 8) {
        $this->scale = $scale;
        $this->exceptions = $exceptions;
        $this->n = 0;
        $this->sum = $this->M1 = $this->M2 = $this->M3 = $this->M4 = 0.0;
        $this->min = $this->max = null;
    }
    /**
     * Caclulate stats for values in array
     *
     * @param array $valuesArray [val1, val2, ...]
     * @return Stats object this
     */
    public function pushArray($valuesArray) {
        foreach($valuesArray as $x) {
            $this->push($x);
        }
        return $this->stats();
    }
    /**
     * Add a value to caclulate
     *
     * @param string|int $x
     * @return void
     */
    public function push($x) {
        if($x === null || $x === '') {
            if(($this->exceptions & Stats::MISSING_SKIP) === Stats::MISSING_SKIP) {
                return;
            }
            $x = 0.00;
        }
        elseif(!is_numeric($x)) {
            if(($this->exceptions & Stats::NA_SKIP) === Stats::NA_SKIP) {
                return;
            }
            $x = 0.00;
        }
        $this->sum = $this->sum + $x;
        if($x < $this->min || $this->min === null) {
            $this->min = $x;
        }
        if($x > $this->max || $this->max === null) {
            $this->max = $x;
        }
        $nl = $this->n;
        $this->n = bcadd($this->n, "1", 0);
        $n = $this->n;
        $delta = bcsub($x, $this->M1, $this->scale);
        $delta_n = bcdiv($delta, $n, $this->scale);
        $delta_n2 = bcmul($delta_n, $delta_n, $this->scale);
        $term1 = bcmul(bcmul($delta, $delta_n, $this->scale), $nl, $this->scale);
        $this->M1 = bcadd($this->M1, $delta_n, $this->scale);
        // M4 += term1 * delta_n2 * (n*n - 3*n + 3) + 6 * delta_n2 * M2 - 4 * delta_n * M3;
        $m4_1 = $this->bcmulArray([
            $term1,
            $delta_n2,
            bcadd( bcsub( bcmul($n, $n, 0), bcmul($n, "3", 0), 0 ), "3", 0)
        ]);
        $m4_2 = $this->bcmulArray([
            "6",
            $delta_n2,
            $this->M2
        ]);
        $m4_3 = $this->bcmulArray([
            "4",
            $delta_n,
            $this->M3
        ]);
        $this->M4 = bcadd($this->M4, bcsub( bcadd($m4_1, $m4_2, $this->scale), $m4_3, $this->scale), $this->scale);
       // M3 += term1 * delta_n * (n - 2) - 3 * delta_n * M2;
        $this->M3 = bcadd(
                $this->M3,
                bcsub(
                    $this->bcmulArray([
                            $term1,
                            $delta_n,
                            bcsub($n, "2", 0)
                    ]),
                    $this->bcmulArray([
                            "3",
                            $delta_n,
                            $this->M2
                    ]),
                    $this->scale
                )
            , $this->scale);
        $this->M2 = bcadd($term1, $this->M2, $this->scale);
    }
    /**
     * Return count or n elements processed
     *
     * @return string|int number of items calculated
     */
    public function count() {
        return $this->n;
    }
    /**
     * Minimum value encountered
     *
     * @return string bcmath number
     */
    public function min() {
        return $this->min;
    }
    /**
     * Average, mean, of calculated values
     *
     * @return string bcmath number
     */
    public function avg() {
        return $this->M1;
    }
    /**
     * Maximum value encountered
     *
     * @return string bcmath number
     */
    public function max() {
        return $this->max;
    }
    /**
     * Sum or total of cacluated values
     *
     * @return string bcmath number
     */
    public function sum() {
        return $this->sum;
    }
    /**
     * Varaince
     *
     * @return string bcmath number
     */
    public function variance() {
        if(bccomp($this->n, 2, 0) < 0) {
            return "0";
        }
        return bcdiv($this->M2, bcsub($this->n, "1", 0), $this->scale);
    }
    /**
     * Standard deviation
     *
     * @return string bcmath number
     */
    public function stdDev() {
        return bcsqrt($this->variance(), $this->scale);
    }
    /**
     * Skewness
     *
     * Characterizes the degree of asymmetry of a distribution around its mean.
     * [-0.5, 0.5] aprox symetric, [0.5,1],[-0.5,-1] moderatly, else highly
     * Positive skewness indicates a distribution with an asymmetric tail extending toward more positive values.
     * Negative skewness indicates a distribution with an asymmetric tail extending toward more negative values.
     *
     * Skewness, has several estimators see:
     * https://www.researchgate.net/post/Is_there_any_difference_in_formula_when_calculating_Skewness_manually_and_using_SPSS_software
     *
     */
    public function skewness() {
        // sqrt(double(n)) * M3/ pow(M2, 1.5)
        $denom = bcmul($this->M2, bcsqrt($this->M2, $this->scale), $this->scale);
        if(bccomp($denom, "0.00", $this->scale) === 0) {
            return "0";
        }
        return
            bcdiv(
                bcmul(bcsqrt($this->n, $this->scale), $this->M3, $this->scale),
                $denom,
                $this->scale
            );
    }
    /**
     * Kurtosis
     *
     * Characterizes the relative peakedness or flatness of a distribution compared with the normal distribution.
     * Positive kurtosis indicates a relatively peaked distribution.
     * Negative kurtosis indicates a relatively flat distribution.
     *
     * A normal distribution has kurtosis exactly 3 (excess kurtosis exactly 0). Any distribution with kurtosis ≈3 (excess ≈0) is called mesokurtic.
     * A distribution with kurtosis <3 (excess kurtosis <0) is called platykurtic. Compared to a normal distribution, its tails are shorter and thinner, and often its central peak is lower and broader.
     * A distribution with kurtosis >3 (excess kurtosis >0) is called leptokurtic. Compared to a normal distribution, its tails are longer and fatter, and often its central peak is higher and sharper.
     *
     * @return string bcmath number
     */
    public function kurtosis() {
        //return double(n)*M4 / (M2*M2) - 3.0;
        if(bccomp($this->M2, "0.00", $this->scale) === 0) {
            return "0";
        }
        return  bcsub(
                    bcdiv(
                        bcmul($this->n, $this->M4, $this->scale),
                        bcmul($this->M2, $this->M2, $this->scale),
                    $this->scale),
                    "3",
                $this->scale);
    }
    /**
     * Coefficient of Variation (cv) stdDev/Avg * 100
     *
     * @return string bcmath number
     */
    public function cv() {
        if(bccomp("0", $this->avg(), $this->scale) == "0.00") {
            return "0";
        }
        return bcmul( bcdiv($this->stdDev(), $this->avg(), $this->scale), "100", $this->scale);
    }
    /**
     * Returns basci statistics on calculated values
     *
     * @return array ['min'=>1,'avg'=>1,'max'=>1,'count'=>1,'sum'=1,'variance'=>1,'stdDev'=>1,'skewness'=>1,'kurtosis'=>1]
     */
    public function stats() {
        return [
            Stats::STAT_MIN => $this->min,
            Stats::STAT_AVG => $this->M1,
            Stats::STAT_MAX => $this->max,
            Stats::STAT_COUNT => $this->n,
            Stats::STAT_SUM => $this->sum,
            Stats::STAT_VARIANCE => $this->variance(),
            Stats::STAT_STDDEV => $this->stdDev(),
            Stats::STAT_CV => $this->cv(),
            Stats::STAT_SKEWNESS => $this->skewness(),
            Stats::STAT_KURTOSIS => $this->kurtosis(),
        ];
    }
    /**
     * Multiply all values in array
     *
     * @param array $arr array of numbers to multiply [2,3,4]
     * @return @return string bcmath number
     */
    protected function bcmulArray($arr) {
        $result = null;
        foreach($arr as $n) {
            if($result === null) {
                $result = $n;
            } else {
                $result = bcmul($result, $n, $this->scale);
            }
        }
        return $result;
    }
}