<?php
require_once dirname(__DIR__) . '/base/hypothesis.php';
require_once dirname(dirname(__DIR__)) . '/zendasmath/basic/describe.php';
require_once dirname(dirname(__DIR__)) . '/zendasmath/distribution/studentT.php';
require_once dirname(dirname(__DIR__)) . '/zendasmath/distribution/f.php';
require_once dirname(dirname(__DIR__)) . '/zendasmath/linearAlgebra/matrix.php';

/**
 * oneWayAnova
 */
class oneWayAnova extends hypothesis
{
    public $datas;

    /**
     * constructor
     *
     * @param object $dataframe
     * @param array $datas
     * @param string $hypo
     * @param string $hypoVal
     * @param string $alterHypo
     * @access public
     * @return void
     */
    public function __construct($dataframe, $datas)
    {
        $this->datas = $datas;

        parent::__construct($dataframe);
    }

    /**
     * ANOVA confidence interval method.
     *
     * @param array $cols
     * @param string $type
     * @param bool $round
     * @access public
     * @return array
     */
    public function ci($round = true)
    {
        $count = $this->count();

        $α = 0.05;

        $nT = $count->nT; // all data count.
        $r  = $count->r;
        $t  = 0;
        if($nT > $r)
        {
            $studentT = new StudentT($nT - $r);
            $t        = $studentT->inverse(1 - $α / 2);
        }

        $S = $this->s(false);

        $result = array();
        foreach($this->datas as $index => $data)
        {
            $n    = count($data);
            $mean = Describe::mean($data, false);
            if(!self::is_num($t) || !self::is_num($S) || $t == 0)
            {
                $result[$index] = array('-', '-');
                continue;
            }
            $y    = $t * $S / sqrt($n);

            $lower = $mean - $y;
            $upper = $mean + $y;
            $lower = $this->roundDigit($lower, $round);
            $upper = $this->roundDigit($upper, $round);
            $result[$index] = array($lower, $upper);
        }

        return $result;
    }

    /**
     * Mean table result.
     *
     * @param bool $round
     * @access public
     * @return array
     */
    public function meanTable($round = true)
    {
        $result = array();
        $ci = $this->ci($round);
        foreach($this->datas as $index => $data)
        {
            $mean = Describe::mean($data);
            $σ    = Describe::standard($data);
            $obj  = new stdclass();
            $obj->n    = count($data);
            $obj->mean = $this->roundDigit($mean, $round);
            $obj->σ    = $this->roundDigit($σ, $round);
            $obj->ci   = $ci[$index];

            $result[$index] = $obj;
        }

        return $result;
    }

    /**
     * ANOVA S result.
     *
     * @param array $cols
     * @param string $type
     * @param bool $round
     * @access public
     * @return array
     */
    public function s($round = true)
    {
        $count = $this->count();

        $nT = $count->nT; // all data count.
        $r  = $count->r;

        $S² = 0;
        foreach($this->datas as $index => $data)
        {
            $n  = count($data);
            if($n < 2) continue;

            $σ  = Describe::standard($data, false);
            $σ² = Describe::variance($data, false);


            $S² += ($n - 1) * $σ² / ($nT - $r);
        }
        $S = $S² == 0 ? '-' : sqrt($S²);

        return $this->roundDigit($S, $round);
    }

    /**
     * ANOVA table result.
     *
     * @param array $cols
     * @param string $type
     * @param bool $round
     * @access public
     * @return array
     */
    public function table($round = true)
    {
        $count = $this->count();

        $nT = $count->nT; // all data count.
        $r  = $count->r;

        $allSum = 0;
        foreach($this->datas as $data)
        {
            $allSum += Describe::sum($data, false);
        }
        $ybarAll = $allSum / $nT;

        $result = new stdclass();
        $result->factorDF = $r - 1;
        $result->errorDF  = $nT - $r;
        $result->totalDF  = $nT - 1;

        $result->factorSS = 0;
        $result->errorSS  = 0;

        foreach($this->datas as $data)
        {
            $n = count($data);
            if($n != 0)
            {
                $mean = Describe::mean($data, false);
                $result->factorSS += $n * pow($mean - $ybarAll, 2);
            }

            if($n < 2) continue;

            foreach($data as $value)
            {
                $result->errorSS += pow($value - $mean, 2);
            }
        }

        $result->totalSS  = $result->factorSS + $result->errorSS;

        $result->factorMS = $result->factorSS / $result->factorDF;
        $result->errorMS  = $result->errorSS != 0 ? $result->errorSS / $result->errorDF : '-';

        $result->F = is_numeric($result->errorMS) ? $result->factorMS / $result->errorMS : '-';

        $result->p = '-';
        if($result->factorDF > 0 && $result->errorDF > 0)
        {
            $fDist = new F($result->factorDF, $result->errorDF);
            $result->p = is_numeric($result->F) ? 1 - $fDist->cdf($result->F) : '-';
        }

        $result->factorSS = $this->roundDigit($result->factorSS, $round);
        $result->errorSS  = $this->roundDigit($result->errorSS, $round);
        $result->totalSS  = $this->roundDigit($result->totalSS, $round);
        $result->factorMS = $this->roundDigit($result->factorMS, $round);
        $result->errorMS  = $this->roundDigit($result->errorMS, $round);
        $result->F        = $this->roundDigit($result->F, $round);
        $result->p        = $this->roundDigit($result->p, $round);

        return $result;
    }

    /**
     * ANOVA model sum reuslt.
     *
     * @param array $cols
     * @param string $type
     * @param bool $round
     * @access public
     * @return array
     */
    public function model($round = true)
    {
        $table = $this->table(false);

        $result         = new stdclass();
        $result->S      = $this->s($round);
        $result->Rsq    = is_numeric($table->errorSS) && is_numeric($table->totalSS) ? (1 - $table->errorSS / $table->totalSS) * 100 : 100;
        $result->RsqAdj = is_numeric($table->errorMS) && is_numeric($table->totalSS) ? (1 - $table->errorMS / ($table->totalSS / $table->totalDF)) * 100 : '-';
        $result->RsqPre = '-';

        if(is_numeric($result->RsqAdj) and $result->RsqAdj < 0) $result->RsqAdj = 0;

        $minColCount = 2;
        foreach($this->datas as $data)
        {
            $n = count($data);
            if($n < $minColCount) $minColCount = $n;
        }

        if(is_numeric($result->RsqAdj) && $minColCount > 1)
        {
            $matrix = $this->designMatrix();
            $X      = zMatrixFactory::create($matrix);
            $Xt     = $X->transpose(); // Xᵀ
            $B      = $Xt->multiply($X)->inverse();// ⟮XᵀX⟯⁻¹
            $H      = $X->multiply($B)->multiply($Xt); //  X⟮XᵀX⟯⁻¹Xᵀ

            $e        = array();
            $allData  = array();

            foreach($this->datas as $index => $data)
            {
                $mean = Describe::mean($data);
                foreach($data as $value)
                {
                    $e[]    = $value - $mean;
                    $allData[] = $value;
                }
            }

            $meanAll = Describe::mean($allData);
            $sum1 = 0;
            $sum2 = 0;
            foreach($allData as $index => $value)
            {
                $sum1 += pow($e[$index] / (1 - $H[$index][$index]), 2);
                $sum2 += pow($value - $meanAll, 2);
            }

            $result->RsqPre = (1 - ($sum1 / $sum2)) * 100;
            $result->RsqPre = $result->RsqPre < 0 ? 0 : $result->RsqPre;
        }

        $result->S      = $this->roundDigit($result->S, $round);
        $result->Rsq    = $this->roundDigit($result->Rsq, $round);
        $result->RsqAdj = $this->roundDigit($result->RsqAdj, $round);
        $result->RsqPre = $this->roundDigit($result->RsqPre, $round);

        return $result;
    }

    /**
     * ANOVA nT and r.
     *
     * @param array $cols
     * @param string $type
     * @param bool $round
     * @access public
     * @return array
     */
    public function count($round = true)
    {
        $nT = 0; // all data count.
        $r  = 0;

        foreach($this->datas as $data)
        {
            $n = count($data);
            if($n > 0)
            {
                $nT += $n;
                $r += 1;
            }
        }

        $result = new stdclass();
        $result->nT = $nT;
        $result->r  = $r;

        return $result;
    }

    /**
     * ANOVA design matrix.
     * https://en.wikipedia.org/wiki/Design_matrix
     *
     * @param array $cols
     * @param string $type
     * @param bool $round
     * @access public
     * @return array
     */
    public function designMatrix()
    {
        $factorNum = count($this->datas);
        $matrix    = array();

        foreach(array_values($this->datas) as $index => $data)
        {
            $n   = count($data);
            $row = array_fill(0, $factorNum, 0);

            $row[$index] = 1;
            $subMatrix   = array_fill(0, $n, $row);

            $matrix = array_merge($matrix, $subMatrix);
        }

        return $matrix;
    }
}
