| Server IP : 79.136.114.73 / Your IP : 216.73.216.55 Web Server : Apache/2.4.7 (Ubuntu) PHP/5.5.9-1ubuntu4.29 OpenSSL/1.0.1f System : Linux b8009 3.13.0-170-generic #220-Ubuntu SMP Thu May 9 12:40:49 UTC 2019 x86_64 User : www-data ( 33) PHP Version : 5.5.9-1ubuntu4.29 Disable Function : pcntl_alarm,pcntl_fork,pcntl_waitpid,pcntl_wait,pcntl_wifexited,pcntl_wifstopped,pcntl_wifsignaled,pcntl_wexitstatus,pcntl_wtermsig,pcntl_wstopsig,pcntl_signal,pcntl_signal_dispatch,pcntl_get_last_error,pcntl_strerror,pcntl_sigprocmask,pcntl_sigwaitinfo,pcntl_sigtimedwait,pcntl_exec,pcntl_getpriority,pcntl_setpriority, MySQL : ON | cURL : ON | WGET : ON | Perl : ON | Python : ON | Sudo : ON | Pkexec : ON Directory : /var/www/appsrv.astacus.se/energicontrol/PolynomialRegression/ |
Upload File : |
<?php
/*=========================================================================*/
/* Name: PolynomialRegression.php */
/* Uses: Calculates and returns coefficients for polynomial regression. */
/* Date: 06/01/2009 */
/* Author: Andrew Que (http://www.DrQue.net/) */
/* Revisions: */
/* 0.8 - 06/01/2009- QUE - Creation. */
/* 0.9 - 06/14/2012- QUE - */
/* + Bug fix: removed notice causes by uninitialized variable. */
/* + Converted naming convention. */
/* + Fix spelling errors (or the ones I found). */
/* + Changed to row-echelon method for solving matrix which is much */
/* faster than the determinant method. */
/* 0.91 - 05/17/2013- QUE - */
/* = Changed name to Polynonial regression as this is more fitting to */
/* to the function. */
/* 0.92 - 12/28/2013- QUE - */
/* + Added forced offset. */
/* 1.00 - 12/29/2013 - QUE - */
/* + Forced offset changed to allow any term to be forced. */
/* Unit complete. Correlation coefficient (r-squared) has been */
/* implemented externally in the demos. */
/* 1.1 - 2015/05/05 - QUE - */
/* + 'interpolate' is now static as it does not need an instance to */
/* operate. Useful if coefficients have been calculated elsewhere. */
/* - Deprecated 'setDegree' function. This is the wrong terminology for */
/* what the function does. It actually sets the number of */
/* coefficients for the polynomial. The degree of the polynomial is */
/* the number of coefficients less one. Made the identical function */
/* 'setNumberOfCoefficient' to replace it. */
/* + Added getter functions for anything that has a set function. */
/* 1.2 - 2015/02/09 - QUE - */
/* + Support for weighting. */
/* + Slight improvement in data accumulation. */
/* 1.2.1 - 2015/02/17 - QUE - */
/* + Bug fix to LinearWeighting class. */
/* */
/* This project is maintained at: */
/* http://PolynomialRegression.drque.net/ */
/* */
/* ----------------------------------------------------------------------- */
/* */
/* Polynomial regression class. */
/* Copyright (C) 2009, 2012-2015 Andrew Que */
/* */
/* This program is free software: you can redistribute it and/or modify */
/* it under the terms of the GNU General Public License as published by */
/* the Free Software Foundation, either version 3 of the License, or */
/* (at your option) any later version. */
/* */
/* This program is distributed in the hope that it will be useful, */
/* but WITHOUT ANY WARRANTY; without even the implied warranty of */
/* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
/* GNU General Public License for more details. */
/* */
/* You should have received a copy of the GNU General Public License */
/* along with this program. If not, see <http://www.gnu.org/licenses/>. */
/* */
/* ----------------------------------------------------------------------- */
/* */
/* (C) Copyright 2009, 2012-2015 */
/* Andrew Que */
/*=========================================================================*/
/**
* Polynomial regression.
*
* <p>
* Used for calculating polynomial regression coefficients. Useful for
* linear and non-linear regression, and polynomial curve fitting.
*
* @package PolynomialRegression
* @author Andrew Que ({@link http://www.DrQue.net/})
* @link http://PolynomialRegression.drque.net/ Project home page.
* @copyright Copyright (c) 2009, 2012-2015, Andrew Que
* @license http://opensource.org/licenses/gpl-license.php GNU Public License
* @version 1.2.1
*/
/**
* Interface for weighting the regression. Implementations of this interface can be used to
* apply a custom weighting value that emphasizes some values more than others.
*
* Note: BC math in the weighting interface is optional as long as the number are not expected
* to become too small or large.
*
* @since Version 1.2.1
*/
interface WeightingInterface
{
/**
* Return the weighting term for the given index.
* Note that the index begins with 1, and counts up for each data point added. Thus when
* data is weighted, the order it is entered is important.
*
* @param int $index Current index for which to return weighting term.
*/
public function getWeight( $index );
}
/**
* Used for calculating polynomial regression coefficients and interpolation using
* those coefficients. Useful for linear and non-linear regression, and polynomial
* curve fitting.
*
* Note: Requires BC math to be compiled into PHP. Higher-degree polynomials end up
* with very large/small numbers, requiring an arbitrary precision arithmetic. Make sure
* to set "bcscale" as coefficients will likely have decimal values.
*
* Quick example of using this unit to calculate linear regression (1st degree polynomial):
*
* <pre>
* $regression = new PolynomialRegression( 2 );
* // ...
* $regression->addData( $x, $y );
* // ...
* $coefficients = $regression->getCoefficients();
* // ...
* $y = $regression->interpolate( $coefficients, $x );
* </pre>
*
*/
class PolynomialRegression
{
/** @internal Array of running sum of x^n. */
private $xPowers;
/** @internal Array of running sum of x*y powers. */
private $xyPowers;
/** @internal Number of coefficients. */
private $numberOfCoefficient;
/** @internal Array of forcing terms. */
private $forcedValue;
/** @internal The index of the current element. Basiclly a count of data added. */
private $index = 0;
/** @internal The weighting interface. NULL is unused. */
private $weightingInterface = NULL;
/**
* Constructor
*
* Create new class.
* @param $numberOfCoefficient Number of coefficients in polynomial (degree
* of polynomial + 1).
*/
public function __construct( $numberOfCoefficient = 3 )
{
$this->numberOfCoefficient = $numberOfCoefficient;
$this->reset();
} // __construct
/**
* Reset data.
*
* Clear all internal data and prepare for new calculation.
* Must be called *after* setNumberOfCoefficient if number of coefficients has
* changed.
*/
public function reset()
{
$this->forcedValue = array();
$this->xPowers = array();
$this->xyPowers = array();
$this->index = 0;
$this->weightingInterface = NULL;
$squares = ( $this->numberOfCoefficient - 1 ) * 2;
// Initialize power arrays.
for ( $index = 0; $index <= $squares; ++$index )
{
$this->xPowers[ $index ] = 0;
$this->xyPowers[ $index ] = 0;
}
} // reset
/**
* Set degree (Deprecated).
*
* This is the maximum number of coefficients the polynomial function that
* will be calculate. Note that the request for coefficients can be lower
* then this value. If number is higher, data must be reset and
* added again.
* @param int $numberOfCoefficient Number of coefficients.
* @deprecated Deprecated in version 1.1 because the name doesn't reflect what
* the function actually does. Use 'setNumberOfCoefficient' instead--operations
* are identical.
*/
public function setDegree( $numberOfCoefficient )
{
$this->numberOfCoefficient = $numberOfCoefficient;
} // setDegree
/**
* Set the number of coefficients to calculate.
*
* This is the maximum number of coefficients the polynomial function that
* will be calculate. Note that the request for coefficients can be lower
* then this value. If number is higher, data must be reset and
* added again.
* @param int $numberOfCoefficient Number of coefficients.
* @since Version 1.1
*/
public function setNumberOfCoefficient( $numberOfCoefficient )
{
$this->numberOfCoefficient = $numberOfCoefficient;
} // setNumberOfCoefficient
/**
* Get the number of coefficients to calculate.
*
* Returns the number of coefficients calculated.
* @return int Number of coefficients.
* @since Version 1.1
*/
public function getNumberOfCoefficient( $numberOfCoefficient )
{
return $this->numberOfCoefficient;
} // getnumberOfCoefficient
/**
* Set a forced coefficient.
*
* Force a coefficient to be assumed a specific value for the calculation.
* Most often used to force an offset of zero, but can be used to set any
* known coefficient for a set of data.
* @param int $coefficient Which coefficient to force.
* @param float $value Value to force this coefficient.
* @since Version 1.0
*/
public function setForcedCoefficient( $coefficient, $value )
{
$this->forcedValue[ $coefficient ] = $value;
} // setForcedCoefficient
/**
* Get a forced coefficient.
*
* Get a previously set forced coefficient.
* @param int $coefficient Which coefficient.
* @return float Value of this force this coefficient. Null if the
* coefficient isn't being forced.
* @since Version 1.1
*/
public function getForcedCoefficient( $coefficient, $value )
{
$result = null;
if ( isset( $this->forcedValue[ $coefficient ] ) )
$result = $this->forcedValue[ $coefficient ];
return $result;
} // getForcedCoefficient
/**
* Set a weighting interface.
*
* The regression can be weighted on a per-index basis using a weighting
* interface. An instance of this interface can be set here.
* @param WeightingInterface $weightingInterface Instance of weighting system
* to be used.
* @since Version 1.2
*/
public function setWeighting( WeightingInterface $weightingInterface )
{
$this->weightingInterface = $weightingInterface;
}
/**
* Get the weighting interface.
*
* Return the current weighting interface being used. Returns NULL if no
* interface is used.
* @param WeightingInterface $weightingInterface Instance of weighting system
* to be used.
* @since Version 1.2
*/
public function getWeighting()
{
return $this->weightingInterface;
}
/**
* Add data.
*
* Add a data point to calculation.
* @param float $x Some real value.
* @param float $y Some real value corresponding to $x.
*/
public function addData( $x, $y )
{
$squares = ( $this->numberOfCoefficient - 1 ) * 2;
// Get weighting term for this index.
$this->index += 1;
$weight = NULL;
if ( NULL !== $this->weightingInterface )
$weight = $this->weightingInterface->getWeight( $this->index );
// Remove the effect of the forced coefficient from this value.
foreach ( $this->forcedValue as $coefficient => $value )
{
$sub = bcpow( $x, $coefficient );
$sub = bcmul( $sub, $value );
$y = bcsub( $y, $sub );
}
// Accumulation of $x raised to the power of the loop iteration.
// $xSum = pow( $x, $index ) starting with pow( $x, 0 ) which is 1.
$xSum = 1;
// Accumulate new data to power sums.
for ( $index = 0; $index <= $squares; ++$index )
{
$accumulator = $xSum;
// Add weighting term (if applicable).
if ( NULL !== $weight )
$accumulator = bcmul( $accumulator, $weight );
$this->xPowers[ $index ] =
bcadd( $this->xPowers[ $index ], $accumulator );
$this->xyPowers[ $index ] =
bcadd( $this->xyPowers[ $index ], bcmul( $y, $accumulator ) );
$xSum = bcmul( $xSum, $x );
}
} // addData
/**
* Get coefficients.
*
* Calculate and return coefficients based on current data.
* @param int $numberOfCoefficient Integer value of the degree polynomial desired. Default
* is -1 which is the max number of coefficients set by class.
* @return array Array of coefficients (as BC strings).
*/
public function getCoefficients( $numberOfCoefficient = -1 )
{
// If no number of coefficients specified, use standard.
if ( $numberOfCoefficient == -1 )
$numberOfCoefficient = $this->numberOfCoefficient;
// Build a matrix.
// The matrix is made up of the sum of powers. So if the number represents the power,
// the matrix will look like this for a 4th degree polynomial:
// [ 0 1 2 3 4 ]
// [ 1 2 3 4 5 ]
// [ 2 3 4 5 6 ]
// [ 3 4 5 6 7 ]
// [ 4 5 6 7 8 ]
//
$matrix = array();
for ( $row = 0; $row < $numberOfCoefficient; ++$row )
{
$matrix[ $row ] = array();
for ( $column = 0; $column < $numberOfCoefficient; ++$column )
$matrix[ $row ][ $column ] =
$this->xPowers[ $row + $column ];
}
// Create augmented matrix by adding X*Y powers.
for ( $row = 0; $row < $numberOfCoefficient; ++$row )
$matrix[ $row ][ $numberOfCoefficient ] = $this->xyPowers[ $row ];
// Add in the forced coefficients. This is done by nulling the row and column
// for each forced coefficient. For example, a 3th degree polynomial
// matrix with have the 2nd coefficient set to F:
// [ a b c d w ] [ a 0 c d w ]
// [ b c d e x ] -> [ 0 1 0 0 F ]
// [ c d e f y ] [ c 0 e f y ]
// [ d e f g z ] [ d 0 f g z ]
foreach ( $this->forcedValue as $coefficient => $value )
{
for ( $index = 0; $index < $numberOfCoefficient; ++$index )
{
$matrix[ $index ][ $coefficient ] = "0";
$matrix[ $coefficient ][ $index ] = "0";
}
$matrix[ $coefficient ][ $coefficient ] = "1";
$matrix[ $coefficient ][ $numberOfCoefficient ] = $value;
}
// Determine number of rows in matrix.
$rows = count( $matrix );
// Initialize done.
$isDone = array();
for ( $column = 0; $column < $rows; ++$column )
$isDone[ $column ] = false;
// This loop will result in an upper-triangle matrix with the
// diagonals all 1--the first part of row-reduction--using 2
// elementary row operations: multiplying a row by a scalar, and
// subtracting a row by a multiple of an other row.
// NOTE: This loop can be done out-of-order. That is, the first
// row may not begin with the first term. Order is tracked in the
// "order" array.
$order = array();
for ( $column = 0; $column < $rows; ++$column )
{
// Find a row to work with.
// A row that has a term in this column, and has not yet been
// reduced.
$activeRow = 0;
while ( ( ( 0 == $matrix[ $activeRow ][ $column ] )
|| ( $isDone[ $activeRow ] ) )
&& ( $activeRow < $rows ) )
{
++$activeRow;
}
// Do we have a term in this row?
if ( $activeRow < $rows )
{
// Remember the order.
$order[ $column ] = $activeRow;
// Normalize row--results in the first term being 1.
$firstTerm = $matrix[ $activeRow ][ $column ];
for ( $subColumn = $column; $subColumn <= $rows; ++$subColumn )
$matrix[ $activeRow ][ $subColumn ] =
bcdiv( $matrix[ $activeRow ][ $subColumn ], $firstTerm );
// This row is finished.
$isDone[ $activeRow ] = true;
// Subtract the active row from all rows that are not finished.
for ( $row = 0; $row < $rows; ++$row )
if ( ( ! $isDone[ $row ] )
&& ( 0 != $matrix[ $row ][ $column ] ) )
{
// Get first term in row.
$firstTerm = $matrix[ $row ][ $column ];
for ( $subColumn = $column; $subColumn <= $rows; ++$subColumn )
{
$accumulator = bcmul( $firstTerm, $matrix[ $activeRow ][ $subColumn ] );
$matrix[ $row ][ $subColumn ] =
bcsub( $matrix[ $row ][ $subColumn ], $accumulator );
}
}
}
}
// Reset done.
for ( $row = 0; $row < $rows; ++$row )
$isDone[ $row ] = false;
$coefficients = array();
// Back-substitution.
// This will solve the matrix completely, resulting in the identity
// matrix in the x-locations, and the coefficients in the last column.
// | 1 0 0 ... 0 c0 |
// | 0 1 0 ... 0 c1 |
// | . . . . . |
// | . . . . . |
// | 0 0 0 ... 1 cn |
for ( $column = ( $rows - 1 ); $column >= 0; --$column )
{
// The active row is based on order.
$activeRow = $order[ $column ];
// The active row is now finished.
$isDone[ $activeRow ] = true;
// For all rows not finished...
for ( $row = 0; $row < $rows; ++$row )
if ( ! $isDone[ $row ] )
{
$firstTerm = $matrix[ $row ][ $column ];
// Back substitution.
for ( $subColumn = $column; $subColumn <= $rows; ++$subColumn )
{
$accumulator =
bcmul( $firstTerm, $matrix[ $activeRow ][ $subColumn ] );
$matrix[ $row ][ $subColumn ] =
bcsub( $matrix[ $row ][ $subColumn ], $accumulator );
}
}
// Save this coefficient for the return.
$coefficients[ $column ] = $matrix[ $activeRow ][ $rows ];
}
// Coefficients are stored backward, so sort them.
ksort( $coefficients );
// Return the coefficients.
return $coefficients;
} // getCoefficients
/**
* Interpolate
*
* Return y point for given x and coefficient set. Function is static as it
* does not require any instance data to operate.
* @param array $coefficients Coefficients as calculated by 'getCoefficients'.
* @param float $x X-coordinate from which to calculate Y.
* @return float Y-coordinate (as floating-point).
*/
static public function interpolate( $coefficients, $x )
{
$numberOfCoefficient = count( $coefficients );
$y = 0;
for ( $coefficentIndex = 0; $coefficentIndex < $numberOfCoefficient; ++$coefficentIndex )
{
// y += coefficients[ coefficentIndex ] * x^coefficentIndex
$y =
bcadd
(
$y,
bcmul
(
$coefficients[ $coefficentIndex ],
bcpow( $x, $coefficentIndex )
)
);
}
return floatval( $y );
} // interpolate
} // Class
// "Science is like an inoculation against charlatans who would have you believe
// whatever it is they tell you." -- Neil deGrasse Tyson
?>