diff --git a/lib/mlbackend/php/classes/processor.php b/lib/mlbackend/php/classes/processor.php new file mode 100644 index 00000000000..065f745c709 --- /dev/null +++ b/lib/mlbackend/php/classes/processor.php @@ -0,0 +1,340 @@ +. + +/** + * Php predictions processor + * + * @package mlbackend_php + * @copyright 2016 David Monllao {@link http://www.davidmonllao.com} + * @license http://www.gnu.org/copyleft/gpl.html GNU GPL v3 or later + */ + +namespace mlbackend_php; + +// TODO No support for 3rd party plugins psr4?? +spl_autoload_register(function($class) { + // Autoload Phpml classes. + $path = __DIR__ . '/../phpml/src/' . str_replace('\\', '/', $class) . '.php'; + if (file_exists($path)) { + require_once($path); + } +}); + +use Phpml\NeuralNetwork\Network\MultilayerPerceptron; +use Phpml\NeuralNetwork\Training\Backpropagation; +use Phpml\CrossValidation\RandomSplit; +use Phpml\Dataset\ArrayDataset; + +defined('MOODLE_INTERNAL') || die(); + +/** + * PHP predictions processor. + * + * @package mlbackend_php + * @copyright 2016 David Monllao {@link http://www.davidmonllao.com} + * @license http://www.gnu.org/copyleft/gpl.html GNU GPL v3 or later + */ +class processor implements \core_analytics\predictor { + + const BATCH_SIZE = 1000; + const TRAIN_ITERATIONS = 20; + const MODEL_FILENAME = 'model.ser'; + + protected $limitedsize = false; + + public function is_ready() { + return true; + } + + public function train($uniqueid, \stored_file $dataset, $outputdir) { + + // Output directory is already unique to the model. + $modelfilepath = $outputdir . DIRECTORY_SEPARATOR . self::MODEL_FILENAME; + + $modelmanager = new \Phpml\ModelManager(); + + if (file_exists($modelfilepath)) { + $classifier = $modelmanager->restoreFromFile($modelfilepath); + } else { + $classifier = new \Phpml\Classification\Linear\Perceptron(0.001, self::TRAIN_ITERATIONS, false); + } + + $fh = $dataset->get_content_file_handle(); + + // The first lines are var names and the second one values. + $metadata = $this->extract_metadata($fh); + + // Skip headers. + fgets($fh); + + $samples = array(); + $targets = array(); + while (($data = fgetcsv($fh)) !== false) { + $sampledata = array_map('floatval', $data); + $samples[] = array_slice($sampledata, 0, $metadata['nfeatures']); + $targets[] = intval($data[$metadata['nfeatures']]); + + if (count($samples) === self::BATCH_SIZE) { + // Training it batches to avoid running out of memory. + + $classifier->partialTrain($samples, $targets, array(0, 1)); + $samples = array(); + $targets = array(); + } + } + fclose($fh); + + // Train the remaining samples. + if ($samples) { + $classifier->partialTrain($samples, $targets, array(0, 1)); + } + + $resultobj = new \stdClass(); + $resultobj->status = \core_analytics\model::OK; + $resultobj->info = array(); + + // Store the trained model. + $modelmanager->saveToFile($classifier, $modelfilepath); + + return $resultobj; + } + + public function predict($uniqueid, \stored_file $dataset, $outputdir) { + + // Output directory is already unique to the model. + $modelfilepath = $outputdir . DIRECTORY_SEPARATOR . self::MODEL_FILENAME; + + if (!file_exists($modelfilepath)) { + throw new \moodle_exception('errorcantloadmodel', 'analytics', '', $modelfilepath); + } + + $modelmanager = new \Phpml\ModelManager(); + $classifier = $modelmanager->restoreFromFile($modelfilepath); + + $fh = $dataset->get_content_file_handle(); + + // The first lines are var names and the second one values. + $metadata = $this->extract_metadata($fh); + + // Skip headers. + fgets($fh); + + $sampleids = array(); + $samples = array(); + $predictions = array(); + while (($data = fgetcsv($fh)) !== false) { + $sampledata = array_map('floatval', $data); + $sampleids[] = $data[0]; + $samples[] = array_slice($sampledata, 1, $metadata['nfeatures']); + + if (count($samples) === self::BATCH_SIZE) { + // Prediction it batches to avoid running out of memory. + + // Append predictions incrementally, we want $sampleids keys in sync with $predictions keys. + $newpredictions = $classifier->predict($samples); + foreach ($newpredictions as $prediction) { + array_push($predictions, $prediction); + } + $samples = array(); + } + } + fclose($fh); + + // Finish the remaining predictions. + if ($samples) { + $predictions = $predictions + $classifier->predict($samples); + } + + $resultobj = new \stdClass(); + $resultobj->status = \core_analytics\model::OK; + $resultobj->info = array(); + + foreach ($predictions as $index => $prediction) { + $resultobj->predictions[$index] = array($sampleids[$index], $prediction); + } + + return $resultobj; + } + + /** + * Evaluates the provided dataset. + * + * During evaluation we need to shuffle the evaluation dataset samples to detect deviated results, + * if the dataset is massive we can not load everything into memory. We know that 2GB is the + * minimum memory limit we should have (\core_analytics\model::increase_memory), if we substract the memory + * that we already consumed and the memory that Phpml algorithms will need we should still have at + * least 500MB of memory, which should be enough to evaluate a model. In any case this is a robust + * solution that will work for all sites but it should minimize memory limit problems. Site admins + * can still set $CFG->mlbackend_php_no_evaluation_limits to true to skip this 500MB limit. + * + * @param string $uniqueid + * @param float $maxdeviation + * @param int $niterations + * @param \stored_file $dataset + * @param string $outputdir + * @return \stdClass + */ + public function evaluate($uniqueid, $maxdeviation, $niterations, \stored_file $dataset, $outputdir) { + $fh = $dataset->get_content_file_handle(); + + // The first lines are var names and the second one values. + $metadata = $this->extract_metadata($fh); + + // Skip headers. + fgets($fh); + + if (empty($CFG->mlbackend_php_no_evaluation_limits)) { + $samplessize = 0; + $limit = get_real_size('500MB'); + + // Just an approximation, will depend on PHP version, compile options... + // Double size + zval struct (6 bytes + 8 bytes + 16 bytes) + array bucket (96 bytes) + // https://nikic.github.io/2011/12/12/How-big-are-PHP-arrays-really-Hint-BIG.html + $floatsize = (PHP_INT_SIZE * 2) + 6 + 8 + 16 + 96; + } + + $samples = array(); + $targets = array(); + while (($data = fgetcsv($fh)) !== false) { + $sampledata = array_map('floatval', $data); + + $samples[] = array_slice($sampledata, 0, $metadata['nfeatures']); + $targets[] = array(intval($data[$metadata['nfeatures']])); + + if (empty($CFG->mlbackend_php_no_evaluation_limits)) { + // We allow admins to disable evaluation memory usage limits by modifying config.php. + + // We will have plenty of missing values in the dataset so it should be a conservative approximation: + $samplessize = $samplessize + (count($sampledata) * $floatsize); + + // Stop fetching more samples. + if ($samplessize >= $limit) { + $this->limitedsize = true; + break; + } + } + } + fclose($fh); + + $phis = array(); + + // Evaluate the model multiple times to confirm the results are not significantly random due to a short amount of data. + for ($i = 0; $i < $niterations; $i++) { + + //$classifier = new \Phpml\Classification\Linear\Perceptron(0.001, self::TRAIN_ITERATIONS, false); + $network = new MultilayerPerceptron([intval($metadata['nfeatures']), 2, 1]); + $training = new Backpropagation($network); + + // Split up the dataset in classifier and testing. + $data = new RandomSplit(new ArrayDataset($samples, $targets), 0.2); + + $training->train($data->getTrainSamples(), $data->getTrainLabels(), 0, 1); + + $predictedlabels = array(); + foreach ($data->getTestSamples() as $input) { + $output = $network->setInput($input)->getOutput(); + $predictedlabels[] = reset($output); + } + $phis[] = $this->get_phi($data->getTestLabels(), $predictedlabels); + } + + // Let's fill the results changing the returned status code depending on the phi-related calculated metrics. + return $this->get_evaluation_result_object($dataset, $phis, $maxdeviation); + } + + protected function get_evaluation_result_object(\stored_file $dataset, $phis, $maxdeviation) { + + if (count($phis) === 1) { + $avgphi = reset($phis); + } else { + $avgphi = \Phpml\Math\Statistic\Mean::arithmetic($phis); + } + + // Standard deviation should ideally be calculated against the area under the curve. + if (count($phis) === 1) { + $modeldev = 0; + } else { + $modeldev = \Phpml\Math\Statistic\StandardDeviation::population($phis); + } + + // Let's fill the results object. + $resultobj = new \stdClass(); + + // Zero is ok, now we add other bits if something is not right. + $resultobj->status = \core_analytics\model::OK; + $resultobj->info = array(); + + // Convert phi to a standard score (from -1 to 1 to a value between 0 and 1). + $resultobj->score = ($avgphi + 1) / 2; + + // If each iteration results varied too much we need more data to confirm that this is a valid model. + if ($modeldev > $maxdeviation) { + $resultobj->status = $resultobj->status + \core_analytics\model::EVALUATE_NOT_ENOUGH_DATA; + $a = new \stdClass(); + $a->deviation = $modeldev; + $a->accepteddeviation = $maxdeviation; + $resultobj->info[] = get_string('errornotenoughdata', 'mlbackend_php', $a); + } + + if ($resultobj->score < \core_analytics\model::MIN_SCORE) { + $resultobj->status = $resultobj->status + \core_analytics\model::EVALUATE_LOW_SCORE; + $a = new \stdClass(); + $a->score = $resultobj->score; + $a->minscore = \core_analytics\model::MIN_SCORE; + $resultobj->info[] = get_string('errorlowscore', 'mlbackend_php', $a); + } + + if ($this->limitedsize === true) { + $resultobj->info[] = get_string('datasetsizelimited', 'mlbackend_php', display_size($dataset->get_filesize())); + } + + return $resultobj; + } + + protected function get_phi($testlabels, $predictedlabels) { + + foreach ($testlabels as $key => $element) { + $value = reset($element); + $testlabels[$key] = $value; + } + + foreach ($predictedlabels as $key => $element) { + $predictedlabels[$key] = ($element > 0.5) ? 1 : 0; + } + + // Binary here only as well. + $matrix = \Phpml\Metric\ConfusionMatrix::compute($testlabels, $predictedlabels, array(0, 1)); + + $tptn = $matrix[0][0] * $matrix[1][1]; + $fpfn = $matrix[1][0] * $matrix[0][1]; + $tpfp = $matrix[0][0] + $matrix[1][0]; + $tpfn = $matrix[0][0] + $matrix[0][1]; + $tnfp = $matrix[1][1] + $matrix[1][0]; + $tnfn = $matrix[1][1] + $matrix[0][1]; + if ($tpfp === 0 || $tpfn === 0 || $tnfp === 0 || $tnfn === 0) { + $phi = 0; + } else { + $phi = ( $tptn - $fpfn ) / sqrt( $tpfp * $tpfn * $tnfp * $tnfn); + } + + return $phi; + } + + protected function extract_metadata($fh) { + $metadata = fgetcsv($fh); + return array_combine($metadata, fgetcsv($fh)); + } +} diff --git a/lib/mlbackend/php/lang/en/mlbackend_php.php b/lib/mlbackend/php/lang/en/mlbackend_php.php new file mode 100644 index 00000000000..6fed20fc534 --- /dev/null +++ b/lib/mlbackend/php/lang/en/mlbackend_php.php @@ -0,0 +1,7 @@ +deviation}, maximum recommended standard deviation = {$a->accepteddeviation}'; +$string['errorlowscore'] = 'The evaluated model prediction accuracy is not very high, some predictions may not be accurate. Model score = {$a->score}, minimum score = {$a->minscore}'; +$string['datasetsizelimited'] = 'Only a part of the evaluation dataset has been evaluated due to its size. Set $CFG->mlbackend_php_no_memory_limit if you are confident that your system can cope a {$a} dataset'; diff --git a/lib/mlbackend/php/phpml/LICENSE b/lib/mlbackend/php/phpml/LICENSE new file mode 100644 index 00000000000..bd5cb2fd6e4 --- /dev/null +++ b/lib/mlbackend/php/phpml/LICENSE @@ -0,0 +1,21 @@ +The MIT License (MIT) + +Copyright (c) 2016 Arkadiusz Kondas + +Permission is hereby granted, free of charge, to any person obtaining a copy +of this software and associated documentation files (the "Software"), to deal +in the Software without restriction, including without limitation the rights +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the Software is +furnished to do so, subject to the following conditions: + +The above copyright notice and this permission notice shall be included in all +copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE +SOFTWARE. diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-predict b/lib/mlbackend/php/phpml/bin/libsvm/svm-predict new file mode 100644 index 00000000000..e6df243d3c2 Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-predict differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-predict-osx b/lib/mlbackend/php/phpml/bin/libsvm/svm-predict-osx new file mode 100644 index 00000000000..480a60b7899 Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-predict-osx differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-predict.exe b/lib/mlbackend/php/phpml/bin/libsvm/svm-predict.exe new file mode 100644 index 00000000000..29760b09f7e Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-predict.exe differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-scale b/lib/mlbackend/php/phpml/bin/libsvm/svm-scale new file mode 100644 index 00000000000..edfb98a6bd7 Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-scale differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-scale-osx b/lib/mlbackend/php/phpml/bin/libsvm/svm-scale-osx new file mode 100644 index 00000000000..0ac83f662a6 Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-scale-osx differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-scale.exe b/lib/mlbackend/php/phpml/bin/libsvm/svm-scale.exe new file mode 100644 index 00000000000..75489ae97d2 Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-scale.exe differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-train b/lib/mlbackend/php/phpml/bin/libsvm/svm-train new file mode 100644 index 00000000000..9cb7f11166b Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-train differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-train-osx b/lib/mlbackend/php/phpml/bin/libsvm/svm-train-osx new file mode 100644 index 00000000000..a716cdabd5f Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-train-osx differ diff --git a/lib/mlbackend/php/phpml/bin/libsvm/svm-train.exe b/lib/mlbackend/php/phpml/bin/libsvm/svm-train.exe new file mode 100644 index 00000000000..23368b5b23b Binary files /dev/null and b/lib/mlbackend/php/phpml/bin/libsvm/svm-train.exe differ diff --git a/lib/mlbackend/php/phpml/bin/phpunit b/lib/mlbackend/php/phpml/bin/phpunit new file mode 120000 index 00000000000..4ba325648fb --- /dev/null +++ b/lib/mlbackend/php/phpml/bin/phpunit @@ -0,0 +1 @@ +../vendor/phpunit/phpunit/phpunit \ No newline at end of file diff --git a/lib/mlbackend/php/phpml/readme_moodle.txt b/lib/mlbackend/php/phpml/readme_moodle.txt new file mode 100644 index 00000000000..4ffb5debb66 --- /dev/null +++ b/lib/mlbackend/php/phpml/readme_moodle.txt @@ -0,0 +1,6 @@ +Current version is 12b8b11 + +# Download latest stable version from https://github.com/php-ai/php-ml +# Remove all files but: + * src/ + * LICENSE diff --git a/lib/mlbackend/php/phpml/src/Phpml/Association/Apriori.php b/lib/mlbackend/php/phpml/src/Phpml/Association/Apriori.php new file mode 100644 index 00000000000..362f25a6ce8 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Association/Apriori.php @@ -0,0 +1,345 @@ +support = $support; + $this->confidence = $confidence; + } + + /** + * Get all association rules which are generated for every k-length frequent item set. + * + * @return mixed[][] + */ + public function getRules() : array + { + if (!$this->large) { + $this->large = $this->apriori(); + } + + if ($this->rules) { + return $this->rules; + } + + $this->rules = []; + + $this->generateAllRules(); + + return $this->rules; + } + + /** + * Generates frequent item sets. + * + * @return mixed[][][] + */ + public function apriori() : array + { + $L = []; + $L[1] = $this->items(); + $L[1] = $this->frequent($L[1]); + + for ($k = 2; !empty($L[$k - 1]); ++$k) { + $L[$k] = $this->candidates($L[$k - 1]); + $L[$k] = $this->frequent($L[$k]); + } + + return $L; + } + + /** + * @param mixed[] $sample + * + * @return mixed[][] + */ + protected function predictSample(array $sample) : array + { + $predicts = array_values(array_filter($this->getRules(), function ($rule) use ($sample) { + return $this->equals($rule[self::ARRAY_KEY_ANTECEDENT], $sample); + })); + + return array_map(function ($rule) { + return $rule[self::ARRAY_KEY_CONSEQUENT]; + }, $predicts); + } + + /** + * Generate rules for each k-length frequent item set. + */ + private function generateAllRules() + { + for ($k = 2; !empty($this->large[$k]); ++$k) { + foreach ($this->large[$k] as $frequent) { + $this->generateRules($frequent); + } + } + } + + /** + * Generate confident rules for frequent item set. + * + * @param mixed[] $frequent + */ + private function generateRules(array $frequent) + { + foreach ($this->antecedents($frequent) as $antecedent) { + if ($this->confidence <= ($confidence = $this->confidence($frequent, $antecedent))) { + $consequent = array_values(array_diff($frequent, $antecedent)); + $this->rules[] = [ + self::ARRAY_KEY_ANTECEDENT => $antecedent, + self::ARRAY_KEY_CONSEQUENT => $consequent, + self::ARRAY_KEY_SUPPORT => $this->support($consequent), + self::ARRAY_KEY_CONFIDENCE => $confidence, + ]; + } + } + } + + /** + * Generates the power set for given item set $sample. + * + * @param mixed[] $sample + * + * @return mixed[][] + */ + private function powerSet(array $sample) : array + { + $results = [[]]; + foreach ($sample as $item) { + foreach ($results as $combination) { + $results[] = array_merge([$item], $combination); + } + } + + return $results; + } + + /** + * Generates all proper subsets for given set $sample without the empty set. + * + * @param mixed[] $sample + * + * @return mixed[][] + */ + private function antecedents(array $sample) : array + { + $cardinality = count($sample); + $antecedents = $this->powerSet($sample); + + return array_filter($antecedents, function ($antecedent) use ($cardinality) { + return (count($antecedent) != $cardinality) && ($antecedent != []); + }); + } + + /** + * Calculates frequent k = 1 item sets. + * + * @return mixed[][] + */ + private function items() : array + { + $items = []; + + foreach ($this->samples as $sample) { + foreach ($sample as $item) { + if (!in_array($item, $items, true)) { + $items[] = $item; + } + } + } + + return array_map(function ($entry) { + return [$entry]; + }, $items); + } + + /** + * Returns frequent item sets only. + * + * @param mixed[][] $samples + * + * @return mixed[][] + */ + private function frequent(array $samples) : array + { + return array_filter($samples, function ($entry) { + return $this->support($entry) >= $this->support; + }); + } + + /** + * Calculates frequent k item sets, where count($samples) == $k - 1. + * + * @param mixed[][] $samples + * + * @return mixed[][] + */ + private function candidates(array $samples) : array + { + $candidates = []; + + foreach ($samples as $p) { + foreach ($samples as $q) { + if (count(array_merge(array_diff($p, $q), array_diff($q, $p))) != 2) { + continue; + } + + $candidate = array_unique(array_merge($p, $q)); + + if ($this->contains($candidates, $candidate)) { + continue; + } + + foreach ((array) $this->samples as $sample) { + if ($this->subset($sample, $candidate)) { + $candidates[] = $candidate; + continue 2; + } + } + } + } + + return $candidates; + } + + /** + * Calculates confidence for $set. Confidence is the relative amount of sets containing $subset which also contain + * $set. + * + * @param mixed[] $set + * @param mixed[] $subset + * + * @return float + */ + private function confidence(array $set, array $subset) : float + { + return $this->support($set) / $this->support($subset); + } + + /** + * Calculates support for item set $sample. Support is the relative amount of sets containing $sample in the data + * pool. + * + * @see \Phpml\Association\Apriori::samples + * + * @param mixed[] $sample + * + * @return float + */ + private function support(array $sample) : float + { + return $this->frequency($sample) / count($this->samples); + } + + /** + * Counts occurrences of $sample as subset in data pool. + * + * @see \Phpml\Association\Apriori::samples + * + * @param mixed[] $sample + * + * @return int + */ + private function frequency(array $sample) : int + { + return count(array_filter($this->samples, function ($entry) use ($sample) { + return $this->subset($entry, $sample); + })); + } + + /** + * Returns true if set is an element of system. + * + * @see \Phpml\Association\Apriori::equals() + * + * @param mixed[][] $system + * @param mixed[] $set + * + * @return bool + */ + private function contains(array $system, array $set) : bool + { + return (bool) array_filter($system, function ($entry) use ($set) { + return $this->equals($entry, $set); + }); + } + + /** + * Returns true if subset is a (proper) subset of set by its items string representation. + * + * @param mixed[] $set + * @param mixed[] $subset + * + * @return bool + */ + private function subset(array $set, array $subset) : bool + { + return !array_diff($subset, array_intersect($subset, $set)); + } + + /** + * Returns true if string representation of items does not differ. + * + * @param mixed[] $set1 + * @param mixed[] $set2 + * + * @return bool + */ + private function equals(array $set1, array $set2) : bool + { + return array_diff($set1, $set2) == array_diff($set2, $set1); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Association/Associator.php b/lib/mlbackend/php/phpml/src/Phpml/Association/Associator.php new file mode 100644 index 00000000000..c339b5e65f3 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Association/Associator.php @@ -0,0 +1,11 @@ +maxDepth = $maxDepth; + } + + /** + * @param array $samples + * @param array $targets + */ + public function train(array $samples, array $targets) + { + $this->samples = array_merge($this->samples, $samples); + $this->targets = array_merge($this->targets, $targets); + + $this->featureCount = count($this->samples[0]); + $this->columnTypes = self::getColumnTypes($this->samples); + $this->labels = array_keys(array_count_values($this->targets)); + $this->tree = $this->getSplitLeaf(range(0, count($this->samples) - 1)); + + // Each time the tree is trained, feature importances are reset so that + // we will have to compute it again depending on the new data + $this->featureImportances = null; + + // If column names are given or computed before, then there is no + // need to init it and accidentally remove the previous given names + if ($this->columnNames === null) { + $this->columnNames = range(0, $this->featureCount - 1); + } elseif (count($this->columnNames) > $this->featureCount) { + $this->columnNames = array_slice($this->columnNames, 0, $this->featureCount); + } elseif (count($this->columnNames) < $this->featureCount) { + $this->columnNames = array_merge($this->columnNames, + range(count($this->columnNames), $this->featureCount - 1)); + } + } + + /** + * @param array $samples + * @return array + */ + public static function getColumnTypes(array $samples) : array + { + $types = []; + $featureCount = count($samples[0]); + for ($i=0; $i < $featureCount; $i++) { + $values = array_column($samples, $i); + $isCategorical = self::isCategoricalColumn($values); + $types[] = $isCategorical ? self::NOMINAL : self::CONTINUOUS; + } + + return $types; + } + + /** + * @param array $records + * @param int $depth + * @return DecisionTreeLeaf + */ + protected function getSplitLeaf(array $records, int $depth = 0) : DecisionTreeLeaf + { + $split = $this->getBestSplit($records); + $split->level = $depth; + if ($this->actualDepth < $depth) { + $this->actualDepth = $depth; + } + + // Traverse all records to see if all records belong to the same class, + // otherwise group the records so that we can classify the leaf + // in case maximum depth is reached + $leftRecords = []; + $rightRecords= []; + $remainingTargets = []; + $prevRecord = null; + $allSame = true; + + foreach ($records as $recordNo) { + // Check if the previous record is the same with the current one + $record = $this->samples[$recordNo]; + if ($prevRecord && $prevRecord != $record) { + $allSame = false; + } + $prevRecord = $record; + + // According to the split criteron, this record will + // belong to either left or the right side in the next split + if ($split->evaluate($record)) { + $leftRecords[] = $recordNo; + } else { + $rightRecords[]= $recordNo; + } + + // Group remaining targets + $target = $this->targets[$recordNo]; + if (! array_key_exists($target, $remainingTargets)) { + $remainingTargets[$target] = 1; + } else { + $remainingTargets[$target]++; + } + } + + if ($allSame || $depth >= $this->maxDepth || count($remainingTargets) === 1) { + $split->isTerminal = 1; + arsort($remainingTargets); + $split->classValue = key($remainingTargets); + } else { + if ($leftRecords) { + $split->leftLeaf = $this->getSplitLeaf($leftRecords, $depth + 1); + } + if ($rightRecords) { + $split->rightLeaf= $this->getSplitLeaf($rightRecords, $depth + 1); + } + } + + return $split; + } + + /** + * @param array $records + * @return DecisionTreeLeaf + */ + protected function getBestSplit(array $records) : DecisionTreeLeaf + { + $targets = array_intersect_key($this->targets, array_flip($records)); + $samples = array_intersect_key($this->samples, array_flip($records)); + $samples = array_combine($records, $this->preprocess($samples)); + $bestGiniVal = 1; + $bestSplit = null; + $features = $this->getSelectedFeatures(); + foreach ($features as $i) { + $colValues = []; + foreach ($samples as $index => $row) { + $colValues[$index] = $row[$i]; + } + $counts = array_count_values($colValues); + arsort($counts); + $baseValue = key($counts); + $gini = $this->getGiniIndex($baseValue, $colValues, $targets); + if ($bestSplit === null || $bestGiniVal > $gini) { + $split = new DecisionTreeLeaf(); + $split->value = $baseValue; + $split->giniIndex = $gini; + $split->columnIndex = $i; + $split->isContinuous = $this->columnTypes[$i] == self::CONTINUOUS; + $split->records = $records; + + // If a numeric column is to be selected, then + // the original numeric value and the selected operator + // will also be saved into the leaf for future access + if ($this->columnTypes[$i] == self::CONTINUOUS) { + $matches = []; + preg_match("/^([<>=]{1,2})\s*(.*)/", strval($split->value), $matches); + $split->operator = $matches[1]; + $split->numericValue = floatval($matches[2]); + } + + $bestSplit = $split; + $bestGiniVal = $gini; + } + } + + return $bestSplit; + } + + /** + * Returns available features/columns to the tree for the decision making + * process.
+ * + * If a number is given with setNumFeatures() method, then a random selection + * of features up to this number is returned.
+ * + * If some features are manually selected by use of setSelectedFeatures(), + * then only these features are returned
+ * + * If any of above methods were not called beforehand, then all features + * are returned by default. + * + * @return array + */ + protected function getSelectedFeatures() : array + { + $allFeatures = range(0, $this->featureCount - 1); + if ($this->numUsableFeatures === 0 && ! $this->selectedFeatures) { + return $allFeatures; + } + + if ($this->selectedFeatures) { + return $this->selectedFeatures; + } + + $numFeatures = $this->numUsableFeatures; + if ($numFeatures > $this->featureCount) { + $numFeatures = $this->featureCount; + } + shuffle($allFeatures); + $selectedFeatures = array_slice($allFeatures, 0, $numFeatures, false); + sort($selectedFeatures); + + return $selectedFeatures; + } + + /** + * @param $baseValue + * @param array $colValues + * @param array $targets + * @return float + */ + public function getGiniIndex($baseValue, array $colValues, array $targets) : float + { + $countMatrix = []; + foreach ($this->labels as $label) { + $countMatrix[$label] = [0, 0]; + } + foreach ($colValues as $index => $value) { + $label = $targets[$index]; + $rowIndex = $value === $baseValue ? 0 : 1; + $countMatrix[$label][$rowIndex]++; + } + $giniParts = [0, 0]; + for ($i=0; $i<=1; $i++) { + $part = 0; + $sum = array_sum(array_column($countMatrix, $i)); + if ($sum > 0) { + foreach ($this->labels as $label) { + $part += pow($countMatrix[$label][$i] / floatval($sum), 2); + } + } + $giniParts[$i] = (1 - $part) * $sum; + } + + return array_sum($giniParts) / count($colValues); + } + + /** + * @param array $samples + * @return array + */ + protected function preprocess(array $samples) : array + { + // Detect and convert continuous data column values into + // discrete values by using the median as a threshold value + $columns = []; + for ($i=0; $i<$this->featureCount; $i++) { + $values = array_column($samples, $i); + if ($this->columnTypes[$i] == self::CONTINUOUS) { + $median = Mean::median($values); + foreach ($values as &$value) { + if ($value <= $median) { + $value = "<= $median"; + } else { + $value = "> $median"; + } + } + } + $columns[] = $values; + } + // Below method is a strange yet very simple & efficient method + // to get the transpose of a 2D array + return array_map(null, ...$columns); + } + + /** + * @param array $columnValues + * @return bool + */ + protected static function isCategoricalColumn(array $columnValues) : bool + { + $count = count($columnValues); + + // There are two main indicators that *may* show whether a + // column is composed of discrete set of values: + // 1- Column may contain string values and non-float values + // 2- Number of unique values in the column is only a small fraction of + // all values in that column (Lower than or equal to %20 of all values) + $numericValues = array_filter($columnValues, 'is_numeric'); + $floatValues = array_filter($columnValues, 'is_float'); + if ($floatValues) { + return false; + } + if (count($numericValues) !== $count) { + return true; + } + + $distinctValues = array_count_values($columnValues); + + return count($distinctValues) <= $count / 5; + } + + /** + * This method is used to set number of columns to be used + * when deciding a split at an internal node of the tree.
+ * If the value is given 0, then all features are used (default behaviour), + * otherwise the given value will be used as a maximum for number of columns + * randomly selected for each split operation. + * + * @param int $numFeatures + * @return $this + * @throws InvalidArgumentException + */ + public function setNumFeatures(int $numFeatures) + { + if ($numFeatures < 0) { + throw new InvalidArgumentException('Selected column count should be greater or equal to zero'); + } + + $this->numUsableFeatures = $numFeatures; + + return $this; + } + + /** + * Used to set predefined features to consider while deciding which column to use for a split + * + * @param array $selectedFeatures + */ + protected function setSelectedFeatures(array $selectedFeatures) + { + $this->selectedFeatures = $selectedFeatures; + } + + /** + * A string array to represent columns. Useful when HTML output or + * column importances are desired to be inspected. + * + * @param array $names + * @return $this + * @throws InvalidArgumentException + */ + public function setColumnNames(array $names) + { + if ($this->featureCount !== 0 && count($names) !== $this->featureCount) { + throw new InvalidArgumentException(sprintf('Length of the given array should be equal to feature count %s', $this->featureCount)); + } + + $this->columnNames = $names; + + return $this; + } + + /** + * @return string + */ + public function getHtml() + { + return $this->tree->getHTML($this->columnNames); + } + + /** + * This will return an array including an importance value for + * each column in the given dataset. The importance values are + * normalized and their total makes 1.
+ * + * @return array + */ + public function getFeatureImportances() + { + if ($this->featureImportances !== null) { + return $this->featureImportances; + } + + $sampleCount = count($this->samples); + $this->featureImportances = []; + foreach ($this->columnNames as $column => $columnName) { + $nodes = $this->getSplitNodesByColumn($column, $this->tree); + + $importance = 0; + foreach ($nodes as $node) { + $importance += $node->getNodeImpurityDecrease($sampleCount); + } + + $this->featureImportances[$columnName] = $importance; + } + + // Normalize & sort the importances + $total = array_sum($this->featureImportances); + if ($total > 0) { + foreach ($this->featureImportances as &$importance) { + $importance /= $total; + } + arsort($this->featureImportances); + } + + return $this->featureImportances; + } + + /** + * Collects and returns an array of internal nodes that use the given + * column as a split criterion + * + * @param int $column + * @param DecisionTreeLeaf $node + * @return array + */ + protected function getSplitNodesByColumn(int $column, DecisionTreeLeaf $node) : array + { + if (!$node || $node->isTerminal) { + return []; + } + + $nodes = []; + if ($node->columnIndex === $column) { + $nodes[] = $node; + } + + $lNodes = []; + $rNodes = []; + if ($node->leftLeaf) { + $lNodes = $this->getSplitNodesByColumn($column, $node->leftLeaf); + } + if ($node->rightLeaf) { + $rNodes = $this->getSplitNodesByColumn($column, $node->rightLeaf); + } + $nodes = array_merge($nodes, $lNodes, $rNodes); + + return $nodes; + } + + /** + * @param array $sample + * @return mixed + */ + protected function predictSample(array $sample) + { + $node = $this->tree; + do { + if ($node->isTerminal) { + break; + } + if ($node->evaluate($sample)) { + $node = $node->leftLeaf; + } else { + $node = $node->rightLeaf; + } + } while ($node); + + return $node ? $node->classValue : $this->labels[0]; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/DecisionTree/DecisionTreeLeaf.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/DecisionTree/DecisionTreeLeaf.php new file mode 100644 index 00000000000..bbb3175112f --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/DecisionTree/DecisionTreeLeaf.php @@ -0,0 +1,171 @@ +columnIndex]; + + if ($this->isContinuous) { + $op = $this->operator; + $value= $this->numericValue; + $recordField = strval($recordField); + eval("\$result = $recordField $op $value;"); + return $result; + } + + return $recordField == $this->value; + } + + /** + * Returns Mean Decrease Impurity (MDI) in the node. + * For terminal nodes, this value is equal to 0 + * + * @return float + */ + public function getNodeImpurityDecrease(int $parentRecordCount) + { + if ($this->isTerminal) { + return 0.0; + } + + $nodeSampleCount = (float)count($this->records); + $iT = $this->giniIndex; + + if ($this->leftLeaf) { + $pL = count($this->leftLeaf->records)/$nodeSampleCount; + $iT -= $pL * $this->leftLeaf->giniIndex; + } + + if ($this->rightLeaf) { + $pR = count($this->rightLeaf->records)/$nodeSampleCount; + $iT -= $pR * $this->rightLeaf->giniIndex; + } + + return $iT * $nodeSampleCount / $parentRecordCount; + } + + /** + * Returns HTML representation of the node including children nodes + * + * @param $columnNames + * @return string + */ + public function getHTML($columnNames = null) + { + if ($this->isTerminal) { + $value = "$this->classValue"; + } else { + $value = $this->value; + if ($columnNames !== null) { + $col = $columnNames[$this->columnIndex]; + } else { + $col = "col_$this->columnIndex"; + } + if (! preg_match("/^[<>=]{1,2}/", $value)) { + $value = "=$value"; + } + $value = "$col $value
Gini: ". number_format($this->giniIndex, 2); + } + $str = ""; + if ($this->leftLeaf || $this->rightLeaf) { + $str .=''; + if ($this->leftLeaf) { + $str .=""; + } else { + $str .=''; + } + $str .=''; + if ($this->rightLeaf) { + $str .=""; + } else { + $str .=''; + } + $str .= ''; + } + $str .= '
+ $value
| Yes
" . $this->leftLeaf->getHTML($columnNames) . "
 No |
" . $this->rightLeaf->getHTML($columnNames) . "
'; + return $str; + } + + /** + * HTML representation of the tree without column names + * + * @return string + */ + public function __toString() + { + return $this->getHTML(); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/AdaBoost.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/AdaBoost.php new file mode 100644 index 00000000000..3d1e4187380 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/AdaBoost.php @@ -0,0 +1,265 @@ +maxIterations = $maxIterations; + } + + /** + * Sets the base classifier that will be used for boosting (default = DecisionStump) + * + * @param string $baseClassifier + * @param array $classifierOptions + */ + public function setBaseClassifier(string $baseClassifier = DecisionStump::class, array $classifierOptions = []) + { + $this->baseClassifier = $baseClassifier; + $this->classifierOptions = $classifierOptions; + } + + /** + * @param array $samples + * @param array $targets + */ + public function train(array $samples, array $targets) + { + // Initialize usual variables + $this->labels = array_keys(array_count_values($targets)); + if (count($this->labels) != 2) { + throw new \Exception("AdaBoost is a binary classifier and can classify between two classes only"); + } + + // Set all target values to either -1 or 1 + $this->labels = [1 => $this->labels[0], -1 => $this->labels[1]]; + foreach ($targets as $target) { + $this->targets[] = $target == $this->labels[1] ? 1 : -1; + } + + $this->samples = array_merge($this->samples, $samples); + $this->featureCount = count($samples[0]); + $this->sampleCount = count($this->samples); + + // Initialize AdaBoost parameters + $this->weights = array_fill(0, $this->sampleCount, 1.0 / $this->sampleCount); + $this->classifiers = []; + $this->alpha = []; + + // Execute the algorithm for a maximum number of iterations + $currIter = 0; + while ($this->maxIterations > $currIter++) { + + // Determine the best 'weak' classifier based on current weights + $classifier = $this->getBestClassifier(); + $errorRate = $this->evaluateClassifier($classifier); + + // Update alpha & weight values at each iteration + $alpha = $this->calculateAlpha($errorRate); + $this->updateWeights($classifier, $alpha); + + $this->classifiers[] = $classifier; + $this->alpha[] = $alpha; + } + } + + /** + * Returns the classifier with the lowest error rate with the + * consideration of current sample weights + * + * @return Classifier + */ + protected function getBestClassifier() + { + $ref = new \ReflectionClass($this->baseClassifier); + if ($this->classifierOptions) { + $classifier = $ref->newInstanceArgs($this->classifierOptions); + } else { + $classifier = $ref->newInstance(); + } + + if (is_subclass_of($classifier, WeightedClassifier::class)) { + $classifier->setSampleWeights($this->weights); + $classifier->train($this->samples, $this->targets); + } else { + list($samples, $targets) = $this->resample(); + $classifier->train($samples, $targets); + } + + return $classifier; + } + + /** + * Resamples the dataset in accordance with the weights and + * returns the new dataset + * + * @return array + */ + protected function resample() + { + $weights = $this->weights; + $std = StandardDeviation::population($weights); + $mean= Mean::arithmetic($weights); + $min = min($weights); + $minZ= (int)round(($min - $mean) / $std); + + $samples = []; + $targets = []; + foreach ($weights as $index => $weight) { + $z = (int)round(($weight - $mean) / $std) - $minZ + 1; + for ($i=0; $i < $z; $i++) { + if (rand(0, 1) == 0) { + continue; + } + $samples[] = $this->samples[$index]; + $targets[] = $this->targets[$index]; + } + } + + return [$samples, $targets]; + } + + /** + * Evaluates the classifier and returns the classification error rate + * + * @param Classifier $classifier + */ + protected function evaluateClassifier(Classifier $classifier) + { + $total = (float) array_sum($this->weights); + $wrong = 0; + foreach ($this->samples as $index => $sample) { + $predicted = $classifier->predict($sample); + if ($predicted != $this->targets[$index]) { + $wrong += $this->weights[$index]; + } + } + + return $wrong / $total; + } + + /** + * Calculates alpha of a classifier + * + * @param float $errorRate + * @return float + */ + protected function calculateAlpha(float $errorRate) + { + if ($errorRate == 0) { + $errorRate = 1e-10; + } + return 0.5 * log((1 - $errorRate) / $errorRate); + } + + /** + * Updates the sample weights + * + * @param Classifier $classifier + * @param float $alpha + */ + protected function updateWeights(Classifier $classifier, float $alpha) + { + $sumOfWeights = array_sum($this->weights); + $weightsT1 = []; + foreach ($this->weights as $index => $weight) { + $desired = $this->targets[$index]; + $output = $classifier->predict($this->samples[$index]); + + $weight *= exp(-$alpha * $desired * $output) / $sumOfWeights; + + $weightsT1[] = $weight; + } + + $this->weights = $weightsT1; + } + + /** + * @param array $sample + * @return mixed + */ + public function predictSample(array $sample) + { + $sum = 0; + foreach ($this->alpha as $index => $alpha) { + $h = $this->classifiers[$index]->predict($sample); + $sum += $h * $alpha; + } + + return $this->labels[ $sum > 0 ? 1 : -1]; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/Bagging.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/Bagging.php new file mode 100644 index 00000000000..1bb20273ec7 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/Bagging.php @@ -0,0 +1,194 @@ + 20]; + + /** + * @var array + */ + protected $classifiers; + + /** + * @var float + */ + protected $subsetRatio = 0.7; + + /** + * @var array + */ + private $samples = []; + + /** + * Creates an ensemble classifier with given number of base classifiers
+ * Default number of base classifiers is 100. + * The more number of base classifiers, the better performance but at the cost of procesing time + * + * @param int $numClassifier + */ + public function __construct($numClassifier = 50) + { + $this->numClassifier = $numClassifier; + } + + /** + * This method determines the ratio of samples used to create the 'bootstrap' subset, + * e.g., random samples drawn from the original dataset with replacement (allow repeats), + * to train each base classifier. + * + * @param float $ratio + * @return $this + * @throws Exception + */ + public function setSubsetRatio(float $ratio) + { + if ($ratio < 0.1 || $ratio > 1.0) { + throw new \Exception("Subset ratio should be between 0.1 and 1.0"); + } + $this->subsetRatio = $ratio; + return $this; + } + + /** + * This method is used to set the base classifier. Default value is + * DecisionTree::class, but any class that implements the Classifier + * can be used.
+ * While giving the parameters of the classifier, the values should be + * given in the order they are in the constructor of the classifier and parameter + * names are neglected. + * + * @param string $classifier + * @param array $classifierOptions + * @return $this + */ + public function setClassifer(string $classifier, array $classifierOptions = []) + { + $this->classifier = $classifier; + $this->classifierOptions = $classifierOptions; + return $this; + } + + /** + * @param array $samples + * @param array $targets + */ + public function train(array $samples, array $targets) + { + $this->samples = array_merge($this->samples, $samples); + $this->targets = array_merge($this->targets, $targets); + $this->featureCount = count($samples[0]); + $this->numSamples = count($this->samples); + + // Init classifiers and train them with bootstrap samples + $this->classifiers = $this->initClassifiers(); + $index = 0; + foreach ($this->classifiers as $classifier) { + list($samples, $targets) = $this->getRandomSubset($index); + $classifier->train($samples, $targets); + ++$index; + } + } + + /** + * @param int $index + * @return array + */ + protected function getRandomSubset(int $index) + { + $samples = []; + $targets = []; + srand($index); + $bootstrapSize = $this->subsetRatio * $this->numSamples; + for ($i=0; $i < $bootstrapSize; $i++) { + $rand = rand(0, $this->numSamples - 1); + $samples[] = $this->samples[$rand]; + $targets[] = $this->targets[$rand]; + } + return [$samples, $targets]; + } + + /** + * @return array + */ + protected function initClassifiers() + { + $classifiers = []; + for ($i=0; $i<$this->numClassifier; $i++) { + $ref = new \ReflectionClass($this->classifier); + if ($this->classifierOptions) { + $obj = $ref->newInstanceArgs($this->classifierOptions); + } else { + $obj = $ref->newInstance(); + } + $classifiers[] = $this->initSingleClassifier($obj, $i); + } + return $classifiers; + } + + /** + * @param Classifier $classifier + * @param int $index + * @return Classifier + */ + protected function initSingleClassifier($classifier, $index) + { + return $classifier; + } + + /** + * @param array $sample + * @return mixed + */ + protected function predictSample(array $sample) + { + $predictions = []; + foreach ($this->classifiers as $classifier) { + /* @var $classifier Classifier */ + $predictions[] = $classifier->predict($sample); + } + + $counts = array_count_values($predictions); + arsort($counts); + reset($counts); + return key($counts); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/RandomForest.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/RandomForest.php new file mode 100644 index 00000000000..273eb21aaa4 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Ensemble/RandomForest.php @@ -0,0 +1,153 @@ +setSubsetRatio(1.0); + } + + /** + * This method is used to determine how many of the original columns (features) + * will be used to construct subsets to train base classifiers.
+ * + * Allowed values: 'sqrt', 'log' or any float number between 0.1 and 1.0
+ * + * Default value for the ratio is 'log' which results in log(numFeatures, 2) + 1 + * features to be taken into consideration while selecting subspace of features + * + * @param mixed $ratio string or float should be given + * @return $this + * @throws Exception + */ + public function setFeatureSubsetRatio($ratio) + { + if (is_float($ratio) && ($ratio < 0.1 || $ratio > 1.0)) { + throw new \Exception("When a float given, feature subset ratio should be between 0.1 and 1.0"); + } + if (is_string($ratio) && $ratio != 'sqrt' && $ratio != 'log') { + throw new \Exception("When a string given, feature subset ratio can only be 'sqrt' or 'log' "); + } + $this->featureSubsetRatio = $ratio; + return $this; + } + + /** + * RandomForest algorithm is usable *only* with DecisionTree + * + * @param string $classifier + * @param array $classifierOptions + * @return $this + */ + public function setClassifer(string $classifier, array $classifierOptions = []) + { + if ($classifier != DecisionTree::class) { + throw new \Exception("RandomForest can only use DecisionTree as base classifier"); + } + + return parent::setClassifer($classifier, $classifierOptions); + } + + /** + * This will return an array including an importance value for + * each column in the given dataset. Importance values for a column + * is the average importance of that column in all trees in the forest + * + * @return array + */ + public function getFeatureImportances() + { + // Traverse each tree and sum importance of the columns + $sum = []; + foreach ($this->classifiers as $tree) { + /* @var $tree DecisionTree */ + $importances = $tree->getFeatureImportances(); + + foreach ($importances as $column => $importance) { + if (array_key_exists($column, $sum)) { + $sum[$column] += $importance; + } else { + $sum[$column] = $importance; + } + } + } + + // Normalize & sort the importance values + $total = array_sum($sum); + foreach ($sum as &$importance) { + $importance /= $total; + } + + arsort($sum); + + return $sum; + } + + /** + * A string array to represent the columns is given. They are useful + * when trying to print some information about the trees such as feature importances + * + * @param array $names + * @return $this + */ + public function setColumnNames(array $names) + { + $this->columnNames = $names; + + return $this; + } + + /** + * @param DecisionTree $classifier + * @param int $index + * @return DecisionTree + */ + protected function initSingleClassifier($classifier, $index) + { + if (is_float($this->featureSubsetRatio)) { + $featureCount = (int)($this->featureSubsetRatio * $this->featureCount); + } elseif ($this->featureCount == 'sqrt') { + $featureCount = (int)sqrt($this->featureCount) + 1; + } else { + $featureCount = (int)log($this->featureCount, 2) + 1; + } + + if ($featureCount >= $this->featureCount) { + $featureCount = $this->featureCount; + } + + if ($this->columnNames === null) { + $this->columnNames = range(0, $this->featureCount - 1); + } + + return $classifier + ->setColumnNames($this->columnNames) + ->setNumFeatures($featureCount); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/KNearestNeighbors.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/KNearestNeighbors.php new file mode 100644 index 00000000000..b52c95bd5be --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/KNearestNeighbors.php @@ -0,0 +1,82 @@ +k = $k; + $this->samples = []; + $this->targets = []; + $this->distanceMetric = $distanceMetric; + } + + /** + * @param array $sample + * + * @return mixed + */ + protected function predictSample(array $sample) + { + $distances = $this->kNeighborsDistances($sample); + + $predictions = array_combine(array_values($this->targets), array_fill(0, count($this->targets), 0)); + + foreach ($distances as $index => $distance) { + ++$predictions[$this->targets[$index]]; + } + + arsort($predictions); + reset($predictions); + + return key($predictions); + } + + /** + * @param array $sample + * + * @return array + * + * @throws \Phpml\Exception\InvalidArgumentException + */ + private function kNeighborsDistances(array $sample) + { + $distances = []; + + foreach ($this->samples as $index => $neighbor) { + $distances[$index] = $this->distanceMetric->distance($sample, $neighbor); + } + + asort($distances); + + return array_slice($distances, 0, $this->k, true); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/Adaline.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/Adaline.php new file mode 100644 index 00000000000..f34dc5c4086 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/Adaline.php @@ -0,0 +1,77 @@ + + * + * Learning rate should be a float value between 0.0(exclusive) and 1.0 (inclusive)
+ * Maximum number of iterations can be an integer value greater than 0
+ * If normalizeInputs is set to true, then every input given to the algorithm will be standardized + * by use of standard deviation and mean calculation + * + * @param int $learningRate + * @param int $maxIterations + */ + public function __construct(float $learningRate = 0.001, int $maxIterations = 1000, + bool $normalizeInputs = true, int $trainingType = self::BATCH_TRAINING) + { + if (! in_array($trainingType, [self::BATCH_TRAINING, self::ONLINE_TRAINING])) { + throw new \Exception("Adaline can only be trained with batch and online/stochastic gradient descent algorithm"); + } + + $this->trainingType = $trainingType; + + parent::__construct($learningRate, $maxIterations, $normalizeInputs); + } + + /** + * Adapts the weights with respect to given samples and targets + * by use of gradient descent learning rule + * + * @param array $samples + * @param array $targets + */ + protected function runTraining(array $samples, array $targets) + { + // The cost function is the sum of squares + $callback = function ($weights, $sample, $target) { + $this->weights = $weights; + + $output = $this->output($sample); + $gradient = $output - $target; + $error = $gradient ** 2; + + return [$error, $gradient]; + }; + + $isBatch = $this->trainingType == self::BATCH_TRAINING; + + return parent::runGradientDescent($samples, $targets, $callback, $isBatch); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/DecisionStump.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/DecisionStump.php new file mode 100644 index 00000000000..99f982ff11c --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/DecisionStump.php @@ -0,0 +1,362 @@ + + * + * If columnIndex is given, then the stump tries to produce a decision node + * on this column, otherwise in cases given the value of -1, the stump itself + * decides which column to take for the decision (Default DecisionTree behaviour) + * + * @param int $columnIndex + */ + public function __construct(int $columnIndex = self::AUTO_SELECT) + { + $this->givenColumnIndex = $columnIndex; + } + + /** + * @param array $samples + * @param array $targets + * @throws \Exception + */ + protected function trainBinary(array $samples, array $targets, array $labels) + { + $this->binaryLabels = $labels; + $this->featureCount = count($samples[0]); + + // If a column index is given, it should be among the existing columns + if ($this->givenColumnIndex > count($samples[0]) - 1) { + $this->givenColumnIndex = self::AUTO_SELECT; + } + + // Check the size of the weights given. + // If none given, then assign 1 as a weight to each sample + if ($this->weights) { + $numWeights = count($this->weights); + if ($numWeights != count($samples)) { + throw new \Exception("Number of sample weights does not match with number of samples"); + } + } else { + $this->weights = array_fill(0, count($samples), 1); + } + + // Determine type of each column as either "continuous" or "nominal" + $this->columnTypes = DecisionTree::getColumnTypes($samples); + + // Try to find the best split in the columns of the dataset + // by calculating error rate for each split point in each column + $columns = range(0, count($samples[0]) - 1); + if ($this->givenColumnIndex != self::AUTO_SELECT) { + $columns = [$this->givenColumnIndex]; + } + + $bestSplit = [ + 'value' => 0, 'operator' => '', + 'prob' => [], 'column' => 0, + 'trainingErrorRate' => 1.0]; + foreach ($columns as $col) { + if ($this->columnTypes[$col] == DecisionTree::CONTINUOUS) { + $split = $this->getBestNumericalSplit($samples, $targets, $col); + } else { + $split = $this->getBestNominalSplit($samples, $targets, $col); + } + + if ($split['trainingErrorRate'] < $bestSplit['trainingErrorRate']) { + $bestSplit = $split; + } + } + + // Assign determined best values to the stump + foreach ($bestSplit as $name => $value) { + $this->{$name} = $value; + } + } + + /** + * While finding best split point for a numerical valued column, + * DecisionStump looks for equally distanced values between minimum and maximum + * values in the column. Given $count value determines how many split + * points to be probed. The more split counts, the better performance but + * worse processing time (Default value is 10.0) + * + * @param float $count + */ + public function setNumericalSplitCount(float $count) + { + $this->numSplitCount = $count; + } + + /** + * Determines best split point for the given column + * + * @param array $samples + * @param array $targets + * @param int $col + * + * @return array + */ + protected function getBestNumericalSplit(array $samples, array $targets, int $col) + { + $values = array_column($samples, $col); + // Trying all possible points may be accomplished in two general ways: + // 1- Try all values in the $samples array ($values) + // 2- Artificially split the range of values into several parts and try them + // We choose the second one because it is faster in larger datasets + $minValue = min($values); + $maxValue = max($values); + $stepSize = ($maxValue - $minValue) / $this->numSplitCount; + + $split = null; + + foreach (['<=', '>'] as $operator) { + // Before trying all possible split points, let's first try + // the average value for the cut point + $threshold = array_sum($values) / (float) count($values); + list($errorRate, $prob) = $this->calculateErrorRate($targets, $threshold, $operator, $values); + if ($split == null || $errorRate < $split['trainingErrorRate']) { + $split = ['value' => $threshold, 'operator' => $operator, + 'prob' => $prob, 'column' => $col, + 'trainingErrorRate' => $errorRate]; + } + + // Try other possible points one by one + for ($step = $minValue; $step <= $maxValue; $step+= $stepSize) { + $threshold = (float)$step; + list($errorRate, $prob) = $this->calculateErrorRate($targets, $threshold, $operator, $values); + if ($errorRate < $split['trainingErrorRate']) { + $split = ['value' => $threshold, 'operator' => $operator, + 'prob' => $prob, 'column' => $col, + 'trainingErrorRate' => $errorRate]; + } + }// for + } + + return $split; + } + + /** + * @param array $samples + * @param array $targets + * @param int $col + * + * @return array + */ + protected function getBestNominalSplit(array $samples, array $targets, int $col) : array + { + $values = array_column($samples, $col); + $valueCounts = array_count_values($values); + $distinctVals= array_keys($valueCounts); + + $split = null; + + foreach (['=', '!='] as $operator) { + foreach ($distinctVals as $val) { + list($errorRate, $prob) = $this->calculateErrorRate($targets, $val, $operator, $values); + + if ($split == null || $split['trainingErrorRate'] < $errorRate) { + $split = ['value' => $val, 'operator' => $operator, + 'prob' => $prob, 'column' => $col, + 'trainingErrorRate' => $errorRate]; + } + } + } + + return $split; + } + + + /** + * + * @param type $leftValue + * @param type $operator + * @param type $rightValue + * + * @return boolean + */ + protected function evaluate($leftValue, $operator, $rightValue) + { + switch ($operator) { + case '>': return $leftValue > $rightValue; + case '>=': return $leftValue >= $rightValue; + case '<': return $leftValue < $rightValue; + case '<=': return $leftValue <= $rightValue; + case '=': return $leftValue === $rightValue; + case '!=': + case '<>': return $leftValue !== $rightValue; + } + + return false; + } + + /** + * Calculates the ratio of wrong predictions based on the new threshold + * value given as the parameter + * + * @param array $targets + * @param float $threshold + * @param string $operator + * @param array $values + * + * @return array + */ + protected function calculateErrorRate(array $targets, float $threshold, string $operator, array $values) : array + { + $wrong = 0.0; + $prob = []; + $leftLabel = $this->binaryLabels[0]; + $rightLabel= $this->binaryLabels[1]; + + foreach ($values as $index => $value) { + if ($this->evaluate($value, $operator, $threshold)) { + $predicted = $leftLabel; + } else { + $predicted = $rightLabel; + } + + $target = $targets[$index]; + if (strval($predicted) != strval($targets[$index])) { + $wrong += $this->weights[$index]; + } + + if (! isset($prob[$predicted][$target])) { + $prob[$predicted][$target] = 0; + } + $prob[$predicted][$target]++; + } + + // Calculate probabilities: Proportion of labels in each leaf + $dist = array_combine($this->binaryLabels, array_fill(0, 2, 0.0)); + foreach ($prob as $leaf => $counts) { + $leafTotal = (float)array_sum($prob[$leaf]); + foreach ($counts as $label => $count) { + if (strval($leaf) == strval($label)) { + $dist[$leaf] = $count / $leafTotal; + } + } + } + + return [$wrong / (float) array_sum($this->weights), $dist]; + } + + /** + * Returns the probability of the sample of belonging to the given label + * + * Probability of a sample is calculated as the proportion of the label + * within the labels of the training samples in the decision node + * + * @param array $sample + * @param mixed $label + * + * @return float + */ + protected function predictProbability(array $sample, $label) : float + { + $predicted = $this->predictSampleBinary($sample); + if (strval($predicted) == strval($label)) { + return $this->prob[$label]; + } + + return 0.0; + } + + /** + * @param array $sample + * + * @return mixed + */ + protected function predictSampleBinary(array $sample) + { + if ($this->evaluate($sample[$this->column], $this->operator, $this->value)) { + return $this->binaryLabels[0]; + } + + return $this->binaryLabels[1]; + } + + /** + * @return void + */ + protected function resetBinary() + { + } + + /** + * @return string + */ + public function __toString() + { + return "IF $this->column $this->operator $this->value " . + "THEN " . $this->binaryLabels[0] . " ". + "ELSE " . $this->binaryLabels[1]; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/LogisticRegression.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/LogisticRegression.php new file mode 100644 index 00000000000..bd56d347a50 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/LogisticRegression.php @@ -0,0 +1,281 @@ + + * - 'log' : log likelihood
+ * - 'sse' : sum of squared errors
+ * + * @var string + */ + protected $costFunction = 'sse'; + + /** + * Regularization term: only 'L2' is supported + * + * @var string + */ + protected $penalty = 'L2'; + + /** + * Lambda (λ) parameter of regularization term. If λ is set to 0, then + * regularization term is cancelled. + * + * @var float + */ + protected $lambda = 0.5; + + /** + * Initalize a Logistic Regression classifier with maximum number of iterations + * and learning rule to be applied
+ * + * Maximum number of iterations can be an integer value greater than 0
+ * If normalizeInputs is set to true, then every input given to the algorithm will be standardized + * by use of standard deviation and mean calculation
+ * + * Cost function can be 'log' for log-likelihood and 'sse' for sum of squared errors
+ * + * Penalty (Regularization term) can be 'L2' or empty string to cancel penalty term + * + * @param int $maxIterations + * @param bool $normalizeInputs + * @param int $trainingType + * @param string $cost + * @param string $penalty + * + * @throws \Exception + */ + public function __construct(int $maxIterations = 500, bool $normalizeInputs = true, + int $trainingType = self::CONJUGATE_GRAD_TRAINING, string $cost = 'sse', + string $penalty = 'L2') + { + $trainingTypes = range(self::BATCH_TRAINING, self::CONJUGATE_GRAD_TRAINING); + if (! in_array($trainingType, $trainingTypes)) { + throw new \Exception("Logistic regression can only be trained with " . + "batch (gradient descent), online (stochastic gradient descent) " . + "or conjugate batch (conjugate gradients) algorithms"); + } + + if (! in_array($cost, ['log', 'sse'])) { + throw new \Exception("Logistic regression cost function can be one of the following: \n" . + "'log' for log-likelihood and 'sse' for sum of squared errors"); + } + + if ($penalty != '' && strtoupper($penalty) !== 'L2') { + throw new \Exception("Logistic regression supports only 'L2' regularization"); + } + + $this->learningRate = 0.001; + + parent::__construct($this->learningRate, $maxIterations, $normalizeInputs); + + $this->trainingType = $trainingType; + $this->costFunction = $cost; + $this->penalty = $penalty; + } + + /** + * Sets the learning rate if gradient descent algorithm is + * selected for training + * + * @param float $learningRate + */ + public function setLearningRate(float $learningRate) + { + $this->learningRate = $learningRate; + } + + /** + * Lambda (λ) parameter of regularization term. If 0 is given, + * then the regularization term is cancelled + * + * @param float $lambda + */ + public function setLambda(float $lambda) + { + $this->lambda = $lambda; + } + + /** + * Adapts the weights with respect to given samples and targets + * by use of selected solver + * + * @param array $samples + * @param array $targets + */ + protected function runTraining(array $samples, array $targets) + { + $callback = $this->getCostFunction(); + + switch ($this->trainingType) { + case self::BATCH_TRAINING: + return $this->runGradientDescent($samples, $targets, $callback, true); + + case self::ONLINE_TRAINING: + return $this->runGradientDescent($samples, $targets, $callback, false); + + case self::CONJUGATE_GRAD_TRAINING: + return $this->runConjugateGradient($samples, $targets, $callback); + } + } + + /** + * Executes Conjugate Gradient method to optimize the + * weights of the LogReg model + */ + protected function runConjugateGradient(array $samples, array $targets, \Closure $gradientFunc) + { + if (empty($this->optimizer)) { + $this->optimizer = (new ConjugateGradient($this->featureCount)) + ->setMaxIterations($this->maxIterations); + } + + $this->weights = $this->optimizer->runOptimization($samples, $targets, $gradientFunc); + $this->costValues = $this->optimizer->getCostValues(); + } + + /** + * Returns the appropriate callback function for the selected cost function + * + * @return \Closure + */ + protected function getCostFunction() + { + $penalty = 0; + if ($this->penalty == 'L2') { + $penalty = $this->lambda; + } + + switch ($this->costFunction) { + case 'log': + /* + * Negative of Log-likelihood cost function to be minimized: + * J(x) = ∑( - y . log(h(x)) - (1 - y) . log(1 - h(x))) + * + * If regularization term is given, then it will be added to the cost: + * for L2 : J(x) = J(x) + λ/m . w + * + * The gradient of the cost function to be used with gradient descent: + * ∇J(x) = -(y - h(x)) = (h(x) - y) + */ + $callback = function ($weights, $sample, $y) use ($penalty) { + $this->weights = $weights; + $hX = $this->output($sample); + + // In cases where $hX = 1 or $hX = 0, the log-likelihood + // value will give a NaN, so we fix these values + if ($hX == 1) { + $hX = 1 - 1e-10; + } + if ($hX == 0) { + $hX = 1e-10; + } + $error = -$y * log($hX) - (1 - $y) * log(1 - $hX); + $gradient = $hX - $y; + + return [$error, $gradient, $penalty]; + }; + + return $callback; + + case 'sse': + /** + * Sum of squared errors or least squared errors cost function: + * J(x) = ∑ (y - h(x))^2 + * + * If regularization term is given, then it will be added to the cost: + * for L2 : J(x) = J(x) + λ/m . w + * + * The gradient of the cost function: + * ∇J(x) = -(h(x) - y) . h(x) . (1 - h(x)) + */ + $callback = function ($weights, $sample, $y) use ($penalty) { + $this->weights = $weights; + $hX = $this->output($sample); + + $error = ($y - $hX) ** 2; + $gradient = -($y - $hX) * $hX * (1 - $hX); + + return [$error, $gradient, $penalty]; + }; + + return $callback; + } + } + + /** + * Returns the output of the network, a float value between 0.0 and 1.0 + * + * @param array $sample + * + * @return float + */ + protected function output(array $sample) + { + $sum = parent::output($sample); + + return 1.0 / (1.0 + exp(-$sum)); + } + + /** + * Returns the class value (either -1 or 1) for the given input + * + * @param array $sample + * @return int + */ + protected function outputClass(array $sample) + { + $output = $this->output($sample); + + if (round($output) > 0.5) { + return 1; + } + + return -1; + } + + /** + * Returns the probability of the sample of belonging to the given label. + * + * The probability is simply taken as the distance of the sample + * to the decision plane. + * + * @param array $sample + * @param mixed $label + */ + protected function predictProbability(array $sample, $label) + { + $predicted = $this->predictSampleBinary($sample); + + if (strval($predicted) == strval($label)) { + $sample = $this->checkNormalizedSample($sample); + return abs($this->output($sample) - 0.5); + } + + return 0.0; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/Perceptron.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/Perceptron.php new file mode 100644 index 00000000000..91ffacf91a4 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/Linear/Perceptron.php @@ -0,0 +1,290 @@ + + * + * Learning rate should be a float value between 0.0(exclusive) and 1.0(inclusive)
+ * Maximum number of iterations can be an integer value greater than 0 + * @param int $learningRate + * @param int $maxIterations + */ + public function __construct(float $learningRate = 0.001, int $maxIterations = 1000, + bool $normalizeInputs = true) + { + if ($learningRate <= 0.0 || $learningRate > 1.0) { + throw new \Exception("Learning rate should be a float value between 0.0(exclusive) and 1.0(inclusive)"); + } + + if ($maxIterations <= 0) { + throw new \Exception("Maximum number of iterations should be an integer greater than 0"); + } + + if ($normalizeInputs) { + $this->normalizer = new Normalizer(Normalizer::NORM_STD); + } + + $this->learningRate = $learningRate; + $this->maxIterations = $maxIterations; + } + + /** + * @param array $samples + * @param array $targets + * @param array $labels + */ + public function partialTrain(array $samples, array $targets, array $labels = []) + { + return $this->trainByLabel($samples, $targets, $labels); + } + + /** + * @param array $samples + * @param array $targets + * @param array $labels + */ + public function trainBinary(array $samples, array $targets, array $labels) + { + if ($this->normalizer) { + $this->normalizer->transform($samples); + } + + // Set all target values to either -1 or 1 + $this->labels = [1 => $labels[0], -1 => $labels[1]]; + foreach ($targets as $key => $target) { + $targets[$key] = strval($target) == strval($this->labels[1]) ? 1 : -1; + } + + // Set samples and feature count vars + $this->featureCount = count($samples[0]); + + $this->runTraining($samples, $targets); + } + + protected function resetBinary() + { + $this->labels = []; + $this->optimizer = null; + $this->featureCount = 0; + $this->weights = null; + $this->costValues = []; + } + + /** + * Normally enabling early stopping for the optimization procedure may + * help saving processing time while in some cases it may result in + * premature convergence.
+ * + * If "false" is given, the optimization procedure will always be executed + * for $maxIterations times + * + * @param bool $enable + */ + public function setEarlyStop(bool $enable = true) + { + $this->enableEarlyStop = $enable; + + return $this; + } + + /** + * Returns the cost values obtained during the training. + * + * @return array + */ + public function getCostValues() + { + return $this->costValues; + } + + /** + * Trains the perceptron model with Stochastic Gradient Descent optimization + * to get the correct set of weights + * + * @param array $samples + * @param array $targets + */ + protected function runTraining(array $samples, array $targets) + { + // The cost function is the sum of squares + $callback = function ($weights, $sample, $target) { + $this->weights = $weights; + + $prediction = $this->outputClass($sample); + $gradient = $prediction - $target; + $error = $gradient**2; + + return [$error, $gradient]; + }; + + $this->runGradientDescent($samples, $targets, $callback); + } + + /** + * Executes a Gradient Descent algorithm for + * the given cost function + * + * @param array $samples + * @param array $targets + */ + protected function runGradientDescent(array $samples, array $targets, \Closure $gradientFunc, bool $isBatch = false) + { + $class = $isBatch ? GD::class : StochasticGD::class; + + if (empty($this->optimizer)) { + $this->optimizer = (new $class($this->featureCount)) + ->setLearningRate($this->learningRate) + ->setMaxIterations($this->maxIterations) + ->setChangeThreshold(1e-6) + ->setEarlyStop($this->enableEarlyStop); + } + + $this->weights = $this->optimizer->runOptimization($samples, $targets, $gradientFunc); + $this->costValues = $this->optimizer->getCostValues(); + } + + /** + * Checks if the sample should be normalized and if so, returns the + * normalized sample + * + * @param array $sample + * + * @return array + */ + protected function checkNormalizedSample(array $sample) + { + if ($this->normalizer) { + $samples = [$sample]; + $this->normalizer->transform($samples); + $sample = $samples[0]; + } + + return $sample; + } + + /** + * Calculates net output of the network as a float value for the given input + * + * @param array $sample + * @return int + */ + protected function output(array $sample) + { + $sum = 0; + foreach ($this->weights as $index => $w) { + if ($index == 0) { + $sum += $w; + } else { + $sum += $w * $sample[$index - 1]; + } + } + + return $sum; + } + + /** + * Returns the class value (either -1 or 1) for the given input + * + * @param array $sample + * @return int + */ + protected function outputClass(array $sample) + { + return $this->output($sample) > 0 ? 1 : -1; + } + + /** + * Returns the probability of the sample of belonging to the given label. + * + * The probability is simply taken as the distance of the sample + * to the decision plane. + * + * @param array $sample + * @param mixed $label + */ + protected function predictProbability(array $sample, $label) + { + $predicted = $this->predictSampleBinary($sample); + + if (strval($predicted) == strval($label)) { + $sample = $this->checkNormalizedSample($sample); + return abs($this->output($sample)); + } + + return 0.0; + } + + /** + * @param array $sample + * @return mixed + */ + protected function predictSampleBinary(array $sample) + { + $sample = $this->checkNormalizedSample($sample); + + $predictedClass = $this->outputClass($sample); + + return $this->labels[ $predictedClass ]; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/NaiveBayes.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/NaiveBayes.php new file mode 100644 index 00000000000..af81b00a086 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/NaiveBayes.php @@ -0,0 +1,184 @@ +samples = array_merge($this->samples, $samples); + $this->targets = array_merge($this->targets, $targets); + $this->sampleCount = count($this->samples); + $this->featureCount = count($this->samples[0]); + + $labelCounts = array_count_values($this->targets); + $this->labels = array_keys($labelCounts); + foreach ($this->labels as $label) { + $samples = $this->getSamplesByLabel($label); + $this->p[$label] = count($samples) / $this->sampleCount; + $this->calculateStatistics($label, $samples); + } + } + + /** + * Calculates vital statistics for each label & feature. Stores these + * values in private array in order to avoid repeated calculation + * @param string $label + * @param array $samples + */ + private function calculateStatistics($label, $samples) + { + $this->std[$label] = array_fill(0, $this->featureCount, 0); + $this->mean[$label]= array_fill(0, $this->featureCount, 0); + $this->dataType[$label] = array_fill(0, $this->featureCount, self::CONTINUOS); + $this->discreteProb[$label] = array_fill(0, $this->featureCount, self::CONTINUOS); + for ($i=0; $i<$this->featureCount; $i++) { + // Get the values of nth column in the samples array + // Mean::arithmetic is called twice, can be optimized + $values = array_column($samples, $i); + $numValues = count($values); + // if the values contain non-numeric data, + // then it should be treated as nominal/categorical/discrete column + if ($values !== array_filter($values, 'is_numeric')) { + $this->dataType[$label][$i] = self::NOMINAL; + $this->discreteProb[$label][$i] = array_count_values($values); + $db = &$this->discreteProb[$label][$i]; + $db = array_map(function ($el) use ($numValues) { + return $el / $numValues; + }, $db); + } else { + $this->mean[$label][$i] = Mean::arithmetic($values); + // Add epsilon in order to avoid zero stdev + $this->std[$label][$i] = 1e-10 + StandardDeviation::population($values, false); + } + } + } + + /** + * Calculates the probability P(label|sample_n) + * + * @param array $sample + * @param int $feature + * @param string $label + * @return float + */ + private function sampleProbability($sample, $feature, $label) + { + $value = $sample[$feature]; + if ($this->dataType[$label][$feature] == self::NOMINAL) { + if (! isset($this->discreteProb[$label][$feature][$value]) || + $this->discreteProb[$label][$feature][$value] == 0) { + return self::EPSILON; + } + return $this->discreteProb[$label][$feature][$value]; + } + $std = $this->std[$label][$feature] ; + $mean= $this->mean[$label][$feature]; + // Calculate the probability density by use of normal/Gaussian distribution + // Ref: https://en.wikipedia.org/wiki/Normal_distribution + // + // In order to avoid numerical errors because of small or zero values, + // some libraries adopt taking log of calculations such as + // scikit-learn did. + // (See : https://github.com/scikit-learn/scikit-learn/blob/master/sklearn/naive_bayes.py) + $pdf = -0.5 * log(2.0 * pi() * $std * $std); + $pdf -= 0.5 * pow($value - $mean, 2) / ($std * $std); + return $pdf; + } + + /** + * Return samples belonging to specific label + * @param string $label + * @return array + */ + private function getSamplesByLabel($label) + { + $samples = []; + for ($i=0; $i<$this->sampleCount; $i++) { + if ($this->targets[$i] == $label) { + $samples[] = $this->samples[$i]; + } + } + return $samples; + } + + /** + * @param array $sample + * @return mixed + */ + protected function predictSample(array $sample) + { + // Use NaiveBayes assumption for each label using: + // P(label|features) = P(label) * P(feature0|label) * P(feature1|label) .... P(featureN|label) + // Then compare probability for each class to determine which label is most likely + $predictions = []; + foreach ($this->labels as $label) { + $p = $this->p[$label]; + for ($i=0; $i<$this->featureCount; $i++) { + $Plf = $this->sampleProbability($sample, $i, $label); + $p += $Plf; + } + $predictions[$label] = $p; + } + arsort($predictions, SORT_NUMERIC); + reset($predictions); + return key($predictions); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Classification/SVC.php b/lib/mlbackend/php/phpml/src/Phpml/Classification/SVC.php new file mode 100644 index 00000000000..38ae9c45015 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Classification/SVC.php @@ -0,0 +1,31 @@ +weights = $weights; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Clustering/Clusterer.php b/lib/mlbackend/php/phpml/src/Phpml/Clustering/Clusterer.php new file mode 100644 index 00000000000..0c58b2e9bad --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Clustering/Clusterer.php @@ -0,0 +1,15 @@ +epsilon = $epsilon; + $this->minSamples = $minSamples; + $this->distanceMetric = $distanceMetric; + } + + /** + * @param array $samples + * + * @return array + */ + public function cluster(array $samples) + { + $clusters = []; + $visited = []; + + foreach ($samples as $index => $sample) { + if (isset($visited[$index])) { + continue; + } + $visited[$index] = true; + + $regionSamples = $this->getSamplesInRegion($sample, $samples); + if (count($regionSamples) >= $this->minSamples) { + $clusters[] = $this->expandCluster($regionSamples, $visited); + } + } + + return $clusters; + } + + /** + * @param array $localSample + * @param array $samples + * + * @return array + */ + private function getSamplesInRegion($localSample, $samples) + { + $region = []; + + foreach ($samples as $index => $sample) { + if ($this->distanceMetric->distance($localSample, $sample) < $this->epsilon) { + $region[$index] = $sample; + } + } + + return $region; + } + + /** + * @param array $samples + * @param array $visited + * + * @return array + */ + private function expandCluster($samples, &$visited) + { + $cluster = []; + + foreach ($samples as $index => $sample) { + if (!isset($visited[$index])) { + $visited[$index] = true; + $regionSamples = $this->getSamplesInRegion($sample, $samples); + if (count($regionSamples) > $this->minSamples) { + $cluster = array_merge($regionSamples, $cluster); + } + } + + $cluster[] = $sample; + } + + return $cluster; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Clustering/FuzzyCMeans.php b/lib/mlbackend/php/phpml/src/Phpml/Clustering/FuzzyCMeans.php new file mode 100644 index 00000000000..424f2f15094 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Clustering/FuzzyCMeans.php @@ -0,0 +1,243 @@ +clustersNumber = $clustersNumber; + $this->fuzziness = $fuzziness; + $this->epsilon = $epsilon; + $this->maxIterations = $maxIterations; + } + + protected function initClusters() + { + // Membership array is a matrix of cluster number by sample counts + // We initilize the membership array with random values + $dim = $this->space->getDimension(); + $this->generateRandomMembership($dim, $this->sampleCount); + $this->updateClusters(); + } + + /** + * @param int $rows + * @param int $cols + */ + protected function generateRandomMembership(int $rows, int $cols) + { + $this->membership = []; + for ($i=0; $i < $rows; $i++) { + $row = []; + $total = 0.0; + for ($k=0; $k < $cols; $k++) { + $val = rand(1, 5) / 10.0; + $row[] = $val; + $total += $val; + } + $this->membership[] = array_map(function ($val) use ($total) { + return $val / $total; + }, $row); + } + } + + protected function updateClusters() + { + $dim = $this->space->getDimension(); + if (! $this->clusters) { + $this->clusters = []; + for ($i=0; $i<$this->clustersNumber; $i++) { + $this->clusters[] = new Cluster($this->space, array_fill(0, $dim, 0.0)); + } + } + + for ($i=0; $i<$this->clustersNumber; $i++) { + $cluster = $this->clusters[$i]; + $center = $cluster->getCoordinates(); + for ($k=0; $k<$dim; $k++) { + $a = $this->getMembershipRowTotal($i, $k, true); + $b = $this->getMembershipRowTotal($i, $k, false); + $center[$k] = $a / $b; + } + $cluster->setCoordinates($center); + } + } + + protected function getMembershipRowTotal(int $row, int $col, bool $multiply) + { + $sum = 0.0; + for ($k = 0; $k < $this->sampleCount; $k++) { + $val = pow($this->membership[$row][$k], $this->fuzziness); + if ($multiply) { + $val *= $this->samples[$k][$col]; + } + $sum += $val; + } + return $sum; + } + + protected function updateMembershipMatrix() + { + for ($i = 0; $i < $this->clustersNumber; $i++) { + for ($k = 0; $k < $this->sampleCount; $k++) { + $distCalc = $this->getDistanceCalc($i, $k); + $this->membership[$i][$k] = 1.0 / $distCalc; + } + } + } + + /** + * + * @param int $row + * @param int $col + * @return float + */ + protected function getDistanceCalc(int $row, int $col) + { + $sum = 0.0; + $distance = new Euclidean(); + $dist1 = $distance->distance( + $this->clusters[$row]->getCoordinates(), + $this->samples[$col]); + for ($j = 0; $j < $this->clustersNumber; $j++) { + $dist2 = $distance->distance( + $this->clusters[$j]->getCoordinates(), + $this->samples[$col]); + $val = pow($dist1 / $dist2, 2.0 / ($this->fuzziness - 1)); + $sum += $val; + } + return $sum; + } + + /** + * The objective is to minimize the distance between all data points + * and all cluster centers. This method returns the summation of all + * these distances + */ + protected function getObjective() + { + $sum = 0.0; + $distance = new Euclidean(); + for ($i = 0; $i < $this->clustersNumber; $i++) { + $clust = $this->clusters[$i]->getCoordinates(); + for ($k = 0; $k < $this->sampleCount; $k++) { + $point = $this->samples[$k]; + $sum += $distance->distance($clust, $point); + } + } + return $sum; + } + + /** + * @return array + */ + public function getMembershipMatrix() + { + return $this->membership; + } + + /** + * @param array|Point[] $samples + * @return array + */ + public function cluster(array $samples) + { + // Initialize variables, clusters and membership matrix + $this->sampleCount = count($samples); + $this->samples =& $samples; + $this->space = new Space(count($samples[0])); + $this->initClusters(); + + // Our goal is minimizing the objective value while + // executing the clustering steps at a maximum number of iterations + $lastObjective = 0.0; + $difference = 0.0; + $iterations = 0; + do { + // Update the membership matrix and cluster centers, respectively + $this->updateMembershipMatrix(); + $this->updateClusters(); + + // Calculate the new value of the objective function + $objectiveVal = $this->getObjective(); + $difference = abs($lastObjective - $objectiveVal); + $lastObjective = $objectiveVal; + } while ($difference > $this->epsilon && $iterations++ <= $this->maxIterations); + + // Attach (hard cluster) each data point to the nearest cluster + for ($k=0; $k<$this->sampleCount; $k++) { + $column = array_column($this->membership, $k); + arsort($column); + reset($column); + $i = key($column); + $cluster = $this->clusters[$i]; + $cluster->attach(new Point($this->samples[$k])); + } + + // Return grouped samples + $grouped = []; + foreach ($this->clusters as $cluster) { + $grouped[] = $cluster->getPoints(); + } + return $grouped; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans.php b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans.php new file mode 100644 index 00000000000..a9e90833f7f --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans.php @@ -0,0 +1,60 @@ +clustersNumber = $clustersNumber; + $this->initialization = $initialization; + } + + /** + * @param array $samples + * + * @return array + */ + public function cluster(array $samples) + { + $space = new Space(count($samples[0])); + foreach ($samples as $sample) { + $space->addPoint($sample); + } + + $clusters = []; + foreach ($space->cluster($this->clustersNumber, $this->initialization) as $cluster) { + $clusters[] = $cluster->getPoints(); + } + + return $clusters; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Cluster.php b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Cluster.php new file mode 100644 index 00000000000..7cb9f126e3c --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Cluster.php @@ -0,0 +1,147 @@ +space = $space; + $this->points = new SplObjectStorage(); + } + + /** + * @return array + */ + public function getPoints() + { + $points = []; + foreach ($this->points as $point) { + $points[] = $point->toArray(); + } + + return $points; + } + + /** + * @return array + */ + public function toArray() + { + return [ + 'centroid' => parent::toArray(), + 'points' => $this->getPoints(), + ]; + } + + /** + * @param Point $point + * + * @return Point + * + * @throws \LogicException + */ + public function attach(Point $point) + { + if ($point instanceof self) { + throw new LogicException('cannot attach a cluster to another'); + } + + $this->points->attach($point); + + return $point; + } + + /** + * @param Point $point + * + * @return Point + */ + public function detach(Point $point) + { + $this->points->detach($point); + + return $point; + } + + /** + * @param SplObjectStorage $points + */ + public function attachAll(SplObjectStorage $points) + { + $this->points->addAll($points); + } + + /** + * @param SplObjectStorage $points + */ + public function detachAll(SplObjectStorage $points) + { + $this->points->removeAll($points); + } + + public function updateCentroid() + { + if (!$count = count($this->points)) { + return; + } + + $centroid = $this->space->newPoint(array_fill(0, $this->dimension, 0)); + + foreach ($this->points as $point) { + for ($n = 0; $n < $this->dimension; ++$n) { + $centroid->coordinates[$n] += $point->coordinates[$n]; + } + } + + for ($n = 0; $n < $this->dimension; ++$n) { + $this->coordinates[$n] = $centroid->coordinates[$n] / $count; + } + } + + /** + * @return Point[]|SplObjectStorage + */ + public function getIterator() + { + return $this->points; + } + + /** + * @return mixed + */ + public function count() + { + return count($this->points); + } + + /** + * @param array $newCoordinates + */ + public function setCoordinates(array $newCoordinates) + { + $this->coordinates = $newCoordinates; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Point.php b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Point.php new file mode 100644 index 00000000000..ce1c44ee5f9 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Point.php @@ -0,0 +1,124 @@ +dimension = count($coordinates); + $this->coordinates = $coordinates; + } + + /** + * @return array + */ + public function toArray() + { + return $this->coordinates; + } + + /** + * @param Point $point + * @param bool $precise + * + * @return int|mixed + */ + public function getDistanceWith(self $point, $precise = true) + { + $distance = 0; + for ($n = 0; $n < $this->dimension; ++$n) { + $difference = $this->coordinates[$n] - $point->coordinates[$n]; + $distance += $difference * $difference; + } + + return $precise ? sqrt((float) $distance) : $distance; + } + + /** + * @param array $points + * + * @return mixed + */ + public function getClosest(array $points) + { + foreach ($points as $point) { + $distance = $this->getDistanceWith($point, false); + + if (!isset($minDistance)) { + $minDistance = $distance; + $minPoint = $point; + continue; + } + + if ($distance < $minDistance) { + $minDistance = $distance; + $minPoint = $point; + } + } + + return $minPoint; + } + + /** + * @return array + */ + public function getCoordinates() + { + return $this->coordinates; + } + + /** + * @param mixed $offset + * + * @return bool + */ + public function offsetExists($offset) + { + return isset($this->coordinates[$offset]); + } + + /** + * @param mixed $offset + * + * @return mixed + */ + public function offsetGet($offset) + { + return $this->coordinates[$offset]; + } + + /** + * @param mixed $offset + * @param mixed $value + */ + public function offsetSet($offset, $value) + { + $this->coordinates[$offset] = $value; + } + + /** + * @param mixed $offset + */ + public function offsetUnset($offset) + { + unset($this->coordinates[$offset]); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Space.php b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Space.php new file mode 100644 index 00000000000..5a4d5305e73 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Clustering/KMeans/Space.php @@ -0,0 +1,259 @@ +dimension = $dimension; + } + + /** + * @return array + */ + public function toArray() + { + $points = []; + foreach ($this as $point) { + $points[] = $point->toArray(); + } + + return ['points' => $points]; + } + + /** + * @param array $coordinates + * + * @return Point + */ + public function newPoint(array $coordinates) + { + if (count($coordinates) != $this->dimension) { + throw new LogicException('('.implode(',', $coordinates).') is not a point of this space'); + } + + return new Point($coordinates); + } + + /** + * @param array $coordinates + * @param null $data + */ + public function addPoint(array $coordinates, $data = null) + { + $this->attach($this->newPoint($coordinates), $data); + } + + /** + * @param Point $point + * @param null $data + */ + public function attach($point, $data = null) + { + if (!$point instanceof Point) { + throw new InvalidArgumentException('can only attach points to spaces'); + } + + parent::attach($point, $data); + } + + /** + * @return int + */ + public function getDimension() + { + return $this->dimension; + } + + /** + * @return array|bool + */ + public function getBoundaries() + { + if (!count($this)) { + return false; + } + + $min = $this->newPoint(array_fill(0, $this->dimension, null)); + $max = $this->newPoint(array_fill(0, $this->dimension, null)); + + foreach ($this as $point) { + for ($n = 0; $n < $this->dimension; ++$n) { + ($min[$n] > $point[$n] || $min[$n] === null) && $min[$n] = $point[$n]; + ($max[$n] < $point[$n] || $max[$n] === null) && $max[$n] = $point[$n]; + } + } + + return [$min, $max]; + } + + /** + * @param Point $min + * @param Point $max + * + * @return Point + */ + public function getRandomPoint(Point $min, Point $max) + { + $point = $this->newPoint(array_fill(0, $this->dimension, null)); + + for ($n = 0; $n < $this->dimension; ++$n) { + $point[$n] = random_int($min[$n], $max[$n]); + } + + return $point; + } + + /** + * @param int $clustersNumber + * @param int $initMethod + * + * @return array|Cluster[] + */ + public function cluster(int $clustersNumber, int $initMethod = KMeans::INIT_RANDOM) + { + $clusters = $this->initializeClusters($clustersNumber, $initMethod); + + do { + } while (!$this->iterate($clusters)); + + return $clusters; + } + + /** + * @param $clustersNumber + * @param $initMethod + * + * @return array|Cluster[] + */ + protected function initializeClusters(int $clustersNumber, int $initMethod) + { + switch ($initMethod) { + case KMeans::INIT_RANDOM: + $clusters = $this->initializeRandomClusters($clustersNumber); + break; + + case KMeans::INIT_KMEANS_PLUS_PLUS: + $clusters = $this->initializeKMPPClusters($clustersNumber); + break; + } + $clusters[0]->attachAll($this); + + return $clusters; + } + + /** + * @param $clusters + * + * @return bool + */ + protected function iterate($clusters) + { + $convergence = true; + + $attach = new SplObjectStorage(); + $detach = new SplObjectStorage(); + + foreach ($clusters as $cluster) { + foreach ($cluster as $point) { + $closest = $point->getClosest($clusters); + + if ($closest !== $cluster) { + isset($attach[$closest]) || $attach[$closest] = new SplObjectStorage(); + isset($detach[$cluster]) || $detach[$cluster] = new SplObjectStorage(); + + $attach[$closest]->attach($point); + $detach[$cluster]->attach($point); + + $convergence = false; + } + } + } + + foreach ($attach as $cluster) { + $cluster->attachAll($attach[$cluster]); + } + + foreach ($detach as $cluster) { + $cluster->detachAll($detach[$cluster]); + } + + foreach ($clusters as $cluster) { + $cluster->updateCentroid(); + } + + return $convergence; + } + + /** + * @param int $clustersNumber + * + * @return array + */ + private function initializeRandomClusters(int $clustersNumber) + { + $clusters = []; + list($min, $max) = $this->getBoundaries(); + + for ($n = 0; $n < $clustersNumber; ++$n) { + $clusters[] = new Cluster($this, $this->getRandomPoint($min, $max)->getCoordinates()); + } + + return $clusters; + } + + /** + * @param int $clustersNumber + * + * @return array + */ + protected function initializeKMPPClusters(int $clustersNumber) + { + $clusters = []; + $this->rewind(); + + $clusters[] = new Cluster($this, $this->current()->getCoordinates()); + + $distances = new SplObjectStorage(); + + for ($i = 1; $i < $clustersNumber; ++$i) { + $sum = 0; + foreach ($this as $point) { + $distance = $point->getDistanceWith($point->getClosest($clusters)); + $sum += $distances[$point] = $distance; + } + + $sum = random_int(0, (int) $sum); + foreach ($this as $point) { + if (($sum -= $distances[$point]) > 0) { + continue; + } + + $clusters[] = new Cluster($this, $point->getCoordinates()); + break; + } + } + + return $clusters; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/RandomSplit.php b/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/RandomSplit.php new file mode 100644 index 00000000000..69c44c12249 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/RandomSplit.php @@ -0,0 +1,30 @@ +getSamples(); + $labels = $dataset->getTargets(); + $datasetSize = count($samples); + $testCount = count($this->testSamples); + + for ($i = $datasetSize; $i > 0; --$i) { + $key = mt_rand(0, $datasetSize - 1); + $setName = (count($this->testSamples) - $testCount) / $datasetSize >= $testSize ? 'train' : 'test'; + + $this->{$setName.'Samples'}[] = $samples[$key]; + $this->{$setName.'Labels'}[] = $labels[$key]; + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/Split.php b/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/Split.php new file mode 100644 index 00000000000..add181cb0be --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/Split.php @@ -0,0 +1,94 @@ += $testSize || 1 <= $testSize) { + throw InvalidArgumentException::percentNotInRange('testSize'); + } + $this->seedGenerator($seed); + + $this->splitDataset($dataset, $testSize); + } + + abstract protected function splitDataset(Dataset $dataset, float $testSize); + + /** + * @return array + */ + public function getTrainSamples() + { + return $this->trainSamples; + } + + /** + * @return array + */ + public function getTestSamples() + { + return $this->testSamples; + } + + /** + * @return array + */ + public function getTrainLabels() + { + return $this->trainLabels; + } + + /** + * @return array + */ + public function getTestLabels() + { + return $this->testLabels; + } + + /** + * @param int|null $seed + */ + protected function seedGenerator(int $seed = null) + { + if (null === $seed) { + mt_srand(); + } else { + mt_srand($seed); + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/StratifiedRandomSplit.php b/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/StratifiedRandomSplit.php new file mode 100644 index 00000000000..e6c80a2f3d5 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/CrossValidation/StratifiedRandomSplit.php @@ -0,0 +1,62 @@ +splitByTarget($dataset); + + foreach ($datasets as $targetSet) { + parent::splitDataset($targetSet, $testSize); + } + } + + /** + * @param Dataset $dataset + * + * @return Dataset[]|array + */ + private function splitByTarget(Dataset $dataset): array + { + $targets = $dataset->getTargets(); + $samples = $dataset->getSamples(); + + $uniqueTargets = array_unique($targets); + $split = array_combine($uniqueTargets, array_fill(0, count($uniqueTargets), [])); + + foreach ($samples as $key => $sample) { + $split[$targets[$key]][] = $sample; + } + + $datasets = $this->createDatasets($uniqueTargets, $split); + + return $datasets; + } + + /** + * @param array $uniqueTargets + * @param array $split + * + * @return array + */ + private function createDatasets(array $uniqueTargets, array $split): array + { + $datasets = []; + foreach ($uniqueTargets as $target) { + $datasets[$target] = new ArrayDataset($split[$target], array_fill(0, count($split[$target]), $target)); + } + + return $datasets; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Dataset/ArrayDataset.php b/lib/mlbackend/php/phpml/src/Phpml/Dataset/ArrayDataset.php new file mode 100644 index 00000000000..6f765fe5715 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Dataset/ArrayDataset.php @@ -0,0 +1,52 @@ +samples = $samples; + $this->targets = $targets; + } + + /** + * @return array + */ + public function getSamples(): array + { + return $this->samples; + } + + /** + * @return array + */ + public function getTargets(): array + { + return $this->targets; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Dataset/CsvDataset.php b/lib/mlbackend/php/phpml/src/Phpml/Dataset/CsvDataset.php new file mode 100644 index 00000000000..8bcd3c49034 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Dataset/CsvDataset.php @@ -0,0 +1,54 @@ +columnNames = array_slice($data, 0, $features); + } else { + $this->columnNames = range(0, $features - 1); + } + + while (($data = fgetcsv($handle, 1000, $delimiter)) !== false) { + $this->samples[] = array_slice($data, 0, $features); + $this->targets[] = $data[$features]; + } + fclose($handle); + } + + /** + * @return array + */ + public function getColumnNames() + { + return $this->columnNames; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Dataset/Dataset.php b/lib/mlbackend/php/phpml/src/Phpml/Dataset/Dataset.php new file mode 100644 index 00000000000..f851d8527df --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Dataset/Dataset.php @@ -0,0 +1,18 @@ +scanRootPath($rootPath); + } + + /** + * @param string $rootPath + */ + private function scanRootPath(string $rootPath) + { + foreach (glob($rootPath.DIRECTORY_SEPARATOR.'*', GLOB_ONLYDIR) as $dir) { + $this->scanDir($dir); + } + } + + /** + * @param string $dir + */ + private function scanDir(string $dir) + { + $target = basename($dir); + + foreach (array_filter(glob($dir.DIRECTORY_SEPARATOR.'*'), 'is_file') as $file) { + $this->samples[] = [file_get_contents($file)]; + $this->targets[] = $target; + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/DimensionReduction/KernelPCA.php b/lib/mlbackend/php/phpml/src/Phpml/DimensionReduction/KernelPCA.php new file mode 100644 index 00000000000..86070c72bbc --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/DimensionReduction/KernelPCA.php @@ -0,0 +1,246 @@ +
+ * Example: $kpca = new KernelPCA(KernelPCA::KERNEL_RBF, null, 2, 15.0); + * will initialize the algorithm with an RBF kernel having the gamma parameter as 15,0.
+ * This transformation will return the same number of rows with only 2 columns. + * + * @param int $kernel + * @param float $totalVariance Total variance to be preserved if numFeatures is not given + * @param int $numFeatures Number of columns to be returned + * @param float $gamma Gamma parameter is used with RBF and Sigmoid kernels + * + * @throws \Exception + */ + public function __construct(int $kernel = self::KERNEL_RBF, $totalVariance = null, $numFeatures = null, $gamma = null) + { + $availableKernels = [self::KERNEL_RBF, self::KERNEL_SIGMOID, self::KERNEL_LAPLACIAN, self::KERNEL_LINEAR]; + if (! in_array($kernel, $availableKernels)) { + throw new \Exception("KernelPCA can be initialized with the following kernels only: Linear, RBF, Sigmoid and Laplacian"); + } + + parent::__construct($totalVariance, $numFeatures); + + $this->kernel = $kernel; + $this->gamma = $gamma; + } + + /** + * Takes a data and returns a lower dimensional version + * of this data while preserving $totalVariance or $numFeatures.
+ * $data is an n-by-m matrix and returned array is + * n-by-k matrix where k <= m + * + * @param array $data + * + * @return array + */ + public function fit(array $data) + { + $numRows = count($data); + $this->data = $data; + + if ($this->gamma === null) { + $this->gamma = 1.0 / $numRows; + } + + $matrix = $this->calculateKernelMatrix($this->data, $numRows); + $matrix = $this->centerMatrix($matrix, $numRows); + + list($this->eigValues, $this->eigVectors) = $this->eigenDecomposition($matrix, $numRows); + + $this->fit = true; + + return Matrix::transposeArray($this->eigVectors); + } + + /** + * Calculates similarity matrix by use of selected kernel function
+ * An n-by-m matrix is given and an n-by-n matrix is returned + * + * @param array $data + * @param int $numRows + * + * @return array + */ + protected function calculateKernelMatrix(array $data, int $numRows) + { + $kernelFunc = $this->getKernel(); + + $matrix = []; + for ($i=0; $i < $numRows; $i++) { + for ($k=0; $k < $numRows; $k++) { + if ($i <= $k) { + $matrix[$i][$k] = $kernelFunc($data[$i], $data[$k]); + } else { + $matrix[$i][$k] = $matrix[$k][$i]; + } + } + } + + return $matrix; + } + + /** + * Kernel matrix is centered in its original space by using the following + * conversion: + * + * K′ = K − N.K − K.N + N.K.N where N is n-by-n matrix filled with 1/n + * + * @param array $matrix + * @param int $n + */ + protected function centerMatrix(array $matrix, int $n) + { + $N = array_fill(0, $n, array_fill(0, $n, 1.0/$n)); + $N = new Matrix($N, false); + $K = new Matrix($matrix, false); + + // K.N (This term is repeated so we cache it once) + $K_N = $K->multiply($N); + // N.K + $N_K = $N->multiply($K); + // N.K.N + $N_K_N = $N->multiply($K_N); + + return $K->subtract($N_K) + ->subtract($K_N) + ->add($N_K_N) + ->toArray(); + } + + /** + * Returns the callable kernel function + * + * @return \Closure + */ + protected function getKernel() + { + switch ($this->kernel) { + case self::KERNEL_LINEAR: + // k(x,y) = xT.y + return function ($x, $y) { + return Matrix::dot($x, $y)[0]; + }; + case self::KERNEL_RBF: + // k(x,y)=exp(-γ.|x-y|) where |..| is Euclidean distance + $dist = new Euclidean(); + return function ($x, $y) use ($dist) { + return exp(-$this->gamma * $dist->sqDistance($x, $y)); + }; + + case self::KERNEL_SIGMOID: + // k(x,y)=tanh(γ.xT.y+c0) where c0=1 + return function ($x, $y) { + $res = Matrix::dot($x, $y)[0] + 1.0; + return tanh($this->gamma * $res); + }; + + case self::KERNEL_LAPLACIAN: + // k(x,y)=exp(-γ.|x-y|) where |..| is Manhattan distance + $dist = new Manhattan(); + return function ($x, $y) use ($dist) { + return exp(-$this->gamma * $dist->distance($x, $y)); + }; + } + } + + /** + * @param array $sample + * + * @return array + */ + protected function getDistancePairs(array $sample) + { + $kernel = $this->getKernel(); + + $pairs = []; + foreach ($this->data as $row) { + $pairs[] = $kernel($row, $sample); + } + + return $pairs; + } + + /** + * @param array $pairs + * + * @return array + */ + protected function projectSample(array $pairs) + { + // Normalize eigenvectors by eig = eigVectors / eigValues + $func = function ($eigVal, $eigVect) { + $m = new Matrix($eigVect, false); + $a = $m->divideByScalar($eigVal)->toArray(); + + return $a[0]; + }; + $eig = array_map($func, $this->eigValues, $this->eigVectors); + + // return k.dot(eig) + return Matrix::dot($pairs, $eig); + } + + /** + * Transforms the given sample to a lower dimensional vector by using + * the variables obtained during the last run of fit. + * + * @param array $sample + * + * @return array + */ + public function transform(array $sample) + { + if (!$this->fit) { + throw new \Exception("KernelPCA has not been fitted with respect to original dataset, please run KernelPCA::fit() first"); + } + + if (is_array($sample[0])) { + throw new \Exception("KernelPCA::transform() accepts only one-dimensional arrays"); + } + + $pairs = $this->getDistancePairs($sample); + + return $this->projectSample($pairs); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/DimensionReduction/PCA.php b/lib/mlbackend/php/phpml/src/Phpml/DimensionReduction/PCA.php new file mode 100644 index 00000000000..422dae4d787 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/DimensionReduction/PCA.php @@ -0,0 +1,228 @@ + + * + * @param float $totalVariance Total explained variance to be preserved + * @param int $numFeatures Number of features to be preserved + * + * @throws \Exception + */ + public function __construct($totalVariance = null, $numFeatures = null) + { + if ($totalVariance !== null && ($totalVariance < 0.1 || $totalVariance > 0.99)) { + throw new \Exception("Total variance can be a value between 0.1 and 0.99"); + } + if ($numFeatures !== null && $numFeatures <= 0) { + throw new \Exception("Number of features to be preserved should be greater than 0"); + } + if ($totalVariance !== null && $numFeatures !== null) { + throw new \Exception("Either totalVariance or numFeatures should be specified in order to run the algorithm"); + } + + if ($numFeatures !== null) { + $this->numFeatures = $numFeatures; + } + if ($totalVariance !== null) { + $this->totalVariance = $totalVariance; + } + } + + /** + * Takes a data and returns a lower dimensional version + * of this data while preserving $totalVariance or $numFeatures.
+ * $data is an n-by-m matrix and returned array is + * n-by-k matrix where k <= m + * + * @param array $data + * + * @return array + */ + public function fit(array $data) + { + $n = count($data[0]); + + $data = $this->normalize($data, $n); + + $covMatrix = Covariance::covarianceMatrix($data, array_fill(0, $n, 0)); + + list($this->eigValues, $this->eigVectors) = $this->eigenDecomposition($covMatrix, $n); + + $this->fit = true; + + return $this->reduce($data); + } + + /** + * @param array $data + * @param int $n + */ + protected function calculateMeans(array $data, int $n) + { + // Calculate means for each dimension + $this->means = []; + for ($i=0; $i < $n; $i++) { + $column = array_column($data, $i); + $this->means[] = Mean::arithmetic($column); + } + } + + /** + * Normalization of the data includes subtracting mean from + * each dimension therefore dimensions will be centered to zero + * + * @param array $data + * @param int $n + * + * @return array + */ + protected function normalize(array $data, int $n) + { + if (empty($this->means)) { + $this->calculateMeans($data, $n); + } + + // Normalize data + foreach ($data as $i => $row) { + for ($k=0; $k < $n; $k++) { + $data[$i][$k] -= $this->means[$k]; + } + } + + return $data; + } + + /** + * Calculates eigenValues and eigenVectors of the given matrix. Returns + * top eigenVectors along with the largest eigenValues. The total explained variance + * of these eigenVectors will be no less than desired $totalVariance value + * + * @param array $matrix + * @param int $n + * + * @return array + */ + protected function eigenDecomposition(array $matrix, int $n) + { + $eig = new EigenvalueDecomposition($matrix); + $eigVals = $eig->getRealEigenvalues(); + $eigVects= $eig->getEigenvectors(); + + $totalEigVal = array_sum($eigVals); + // Sort eigenvalues in descending order + arsort($eigVals); + + $explainedVar = 0.0; + $vectors = []; + $values = []; + foreach ($eigVals as $i => $eigVal) { + $explainedVar += $eigVal / $totalEigVal; + $vectors[] = $eigVects[$i]; + $values[] = $eigVal; + + if ($this->numFeatures !== null) { + if (count($vectors) == $this->numFeatures) { + break; + } + } else { + if ($explainedVar >= $this->totalVariance) { + break; + } + } + } + + return [$values, $vectors]; + } + + /** + * Returns the reduced data + * + * @param array $data + * + * @return array + */ + protected function reduce(array $data) + { + $m1 = new Matrix($data); + $m2 = new Matrix($this->eigVectors); + + return $m1->multiply($m2->transpose())->toArray(); + } + + /** + * Transforms the given sample to a lower dimensional vector by using + * the eigenVectors obtained in the last run of fit. + * + * @param array $sample + * + * @return array + */ + public function transform(array $sample) + { + if (!$this->fit) { + throw new \Exception("PCA has not been fitted with respect to original dataset, please run PCA::fit() first"); + } + + if (! is_array($sample[0])) { + $sample = [$sample]; + } + + $sample = $this->normalize($sample, count($sample[0])); + + return $this->reduce($sample); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Estimator.php b/lib/mlbackend/php/phpml/src/Phpml/Estimator.php new file mode 100644 index 00000000000..8b98bb637c7 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Estimator.php @@ -0,0 +1,21 @@ +stopWords = array_fill_keys($stopWords, true); + } + + /** + * @param string $token + * + * @return bool + */ + public function isStopWord(string $token): bool + { + return isset($this->stopWords[$token]); + } + + /** + * @param string $language + * + * @return StopWords + * + * @throws InvalidArgumentException + */ + public static function factory($language = 'English'): StopWords + { + $className = __NAMESPACE__."\\StopWords\\$language"; + + if (!class_exists($className)) { + throw InvalidArgumentException::invalidStopWordsLanguage($language); + } + + return new $className(); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/StopWords/English.php b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/StopWords/English.php new file mode 100644 index 00000000000..fab079b1221 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/StopWords/English.php @@ -0,0 +1,33 @@ +stopWords); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/StopWords/Polish.php b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/StopWords/Polish.php new file mode 100644 index 00000000000..e452ebf2551 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/StopWords/Polish.php @@ -0,0 +1,30 @@ +stopWords); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/TfIdfTransformer.php b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/TfIdfTransformer.php new file mode 100644 index 00000000000..93357752444 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/TfIdfTransformer.php @@ -0,0 +1,66 @@ +fit($samples); + } + } + + /** + * @param array $samples + */ + public function fit(array $samples) + { + $this->countTokensFrequency($samples); + + $count = count($samples); + foreach ($this->idf as &$value) { + $value = log((float)($count / $value), 10.0); + } + } + + /** + * @param array $samples + */ + public function transform(array &$samples) + { + foreach ($samples as &$sample) { + foreach ($sample as $index => &$feature) { + $feature *= $this->idf[$index]; + } + } + } + + /** + * @param array $samples + */ + private function countTokensFrequency(array $samples) + { + $this->idf = array_fill_keys(array_keys($samples[0]), 0); + + foreach ($samples as $sample) { + foreach ($sample as $index => $count) { + if ($count > 0) { + ++$this->idf[$index]; + } + } + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/TokenCountVectorizer.php b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/TokenCountVectorizer.php new file mode 100644 index 00000000000..f5fab21c29f --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/FeatureExtraction/TokenCountVectorizer.php @@ -0,0 +1,214 @@ +tokenizer = $tokenizer; + $this->stopWords = $stopWords; + $this->minDF = $minDF; + + $this->vocabulary = []; + $this->frequencies = []; + } + + /** + * @param array $samples + */ + public function fit(array $samples) + { + $this->buildVocabulary($samples); + } + + /** + * @param array $samples + */ + public function transform(array &$samples) + { + foreach ($samples as &$sample) { + $this->transformSample($sample); + } + + $this->checkDocumentFrequency($samples); + } + + /** + * @return array + */ + public function getVocabulary() + { + return array_flip($this->vocabulary); + } + + /** + * @param array $samples + */ + private function buildVocabulary(array &$samples) + { + foreach ($samples as $index => $sample) { + $tokens = $this->tokenizer->tokenize($sample); + foreach ($tokens as $token) { + $this->addTokenToVocabulary($token); + } + } + } + + /** + * @param string $sample + */ + private function transformSample(string &$sample) + { + $counts = []; + $tokens = $this->tokenizer->tokenize($sample); + + foreach ($tokens as $token) { + $index = $this->getTokenIndex($token); + if (false !== $index) { + $this->updateFrequency($token); + if (!isset($counts[$index])) { + $counts[$index] = 0; + } + + ++$counts[$index]; + } + } + + foreach ($this->vocabulary as $index) { + if (!isset($counts[$index])) { + $counts[$index] = 0; + } + } + + ksort($counts); + + $sample = $counts; + } + + /** + * @param string $token + * + * @return int|bool + */ + private function getTokenIndex(string $token) + { + if ($this->isStopWord($token)) { + return false; + } + + return $this->vocabulary[$token] ?? false; + } + + /** + * @param string $token + */ + private function addTokenToVocabulary(string $token) + { + if ($this->isStopWord($token)) { + return; + } + + if (!isset($this->vocabulary[$token])) { + $this->vocabulary[$token] = count($this->vocabulary); + } + } + + /** + * @param string $token + * + * @return bool + */ + private function isStopWord(string $token): bool + { + return $this->stopWords && $this->stopWords->isStopWord($token); + } + + /** + * @param string $token + */ + private function updateFrequency(string $token) + { + if (!isset($this->frequencies[$token])) { + $this->frequencies[$token] = 0; + } + + ++$this->frequencies[$token]; + } + + /** + * @param array $samples + */ + private function checkDocumentFrequency(array &$samples) + { + if ($this->minDF > 0) { + $beyondMinimum = $this->getBeyondMinimumIndexes(count($samples)); + foreach ($samples as &$sample) { + $this->resetBeyondMinimum($sample, $beyondMinimum); + } + } + } + + /** + * @param array $sample + * @param array $beyondMinimum + */ + private function resetBeyondMinimum(array &$sample, array $beyondMinimum) + { + foreach ($beyondMinimum as $index) { + $sample[$index] = 0; + } + } + + /** + * @param int $samplesCount + * + * @return array + */ + private function getBeyondMinimumIndexes(int $samplesCount) + { + $indexes = []; + foreach ($this->frequencies as $token => $frequency) { + if (($frequency / $samplesCount) < $this->minDF) { + $indexes[] = $this->getTokenIndex($token); + } + } + + return $indexes; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/OneVsRest.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/OneVsRest.php new file mode 100644 index 00000000000..e207c46b5df --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/OneVsRest.php @@ -0,0 +1,197 @@ +reset(); + + return $this->trainBylabel($samples, $targets); + } + + /** + * @param array $samples + * @param array $targets + * @param array $allLabels All training set labels + * @return void + */ + protected function trainByLabel(array $samples, array $targets, array $allLabels = []) + { + + // Overwrites the current value if it exist. $allLabels must be provided for each partialTrain run. + if (!empty($allLabels)) { + $this->allLabels = $allLabels; + } else { + $this->allLabels = array_keys(array_count_values($targets)); + } + sort($this->allLabels, SORT_STRING); + + // If there are only two targets, then there is no need to perform OvR + if (count($this->allLabels) == 2) { + + // Init classifier if required. + if (empty($this->classifiers)) { + $this->classifiers[0] = $this->getClassifierCopy(); + } + + $this->classifiers[0]->trainBinary($samples, $targets, $this->allLabels); + } else { + // Train a separate classifier for each label and memorize them + + foreach ($this->allLabels as $label) { + + // Init classifier if required. + if (empty($this->classifiers[$label])) { + $this->classifiers[$label] = $this->getClassifierCopy(); + } + + list($binarizedTargets, $classifierLabels) = $this->binarizeTargets($targets, $label); + $this->classifiers[$label]->trainBinary($samples, $binarizedTargets, $classifierLabels); + } + } + + // If the underlying classifier is capable of giving the cost values + // during the training, then assign it to the relevant variable + // Adding just the first classifier cost values to avoid complex average calculations. + $classifierref = reset($this->classifiers); + if (method_exists($classifierref, 'getCostValues')) { + $this->costValues = $classifierref->getCostValues(); + } + } + + /** + * Resets the classifier and the vars internally used by OneVsRest to create multiple classifiers. + */ + public function reset() + { + $this->classifiers = []; + $this->allLabels = []; + $this->costValues = []; + + $this->resetBinary(); + } + + /** + * Returns an instance of the current class after cleaning up OneVsRest stuff. + * + * @return \Phpml\Estimator + */ + protected function getClassifierCopy() + { + + // Clone the current classifier, so that + // we don't mess up its variables while training + // multiple instances of this classifier + $classifier = clone $this; + $classifier->reset(); + return $classifier; + } + + /** + * Groups all targets into two groups: Targets equal to + * the given label and the others + * + * $targets is not passed by reference nor contains objects so this method + * changes will not affect the caller $targets array. + * + * @param array $targets + * @param mixed $label + * @return array Binarized targets and target's labels + */ + private function binarizeTargets($targets, $label) + { + $notLabel = "not_$label"; + foreach ($targets as $key => $target) { + $targets[$key] = $target == $label ? $label : $notLabel; + } + + $labels = [$label, $notLabel]; + return [$targets, $labels]; + } + + + /** + * @param array $sample + * + * @return mixed + */ + protected function predictSample(array $sample) + { + if (count($this->allLabels) == 2) { + return $this->classifiers[0]->predictSampleBinary($sample); + } + + $probs = []; + + foreach ($this->classifiers as $label => $predictor) { + $probs[$label] = $predictor->predictProbability($sample, $label); + } + + arsort($probs, SORT_NUMERIC); + return key($probs); + } + + /** + * Each classifier should implement this method instead of train(samples, targets) + * + * @param array $samples + * @param array $targets + * @param array $labels + */ + abstract protected function trainBinary(array $samples, array $targets, array $labels); + + /** + * To be overwritten by OneVsRest classifiers. + * + * @return void + */ + abstract protected function resetBinary(); + + /** + * Each classifier that make use of OvR approach should be able to + * return a probability for a sample to belong to the given label. + * + * @param array $sample + * + * @return mixed + */ + abstract protected function predictProbability(array $sample, string $label); + + /** + * Each classifier should implement this method instead of predictSample() + * + * @param array $sample + * + * @return mixed + */ + abstract protected function predictSampleBinary(array $sample); +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/ConjugateGradient.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/ConjugateGradient.php new file mode 100644 index 00000000000..18ae89a09e8 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/ConjugateGradient.php @@ -0,0 +1,361 @@ +samples = $samples; + $this->targets = $targets; + $this->gradientCb = $gradientCb; + $this->sampleCount = count($samples); + $this->costValues = []; + + $d = mp::muls($this->gradient($this->theta), -1); + + for ($i=0; $i < $this->maxIterations; $i++) { + // Obtain α that minimizes f(θ + α.d) + $alpha = $this->getAlpha(array_sum($d)); + + // θ(k+1) = θ(k) + α.d + $thetaNew = $this->getNewTheta($alpha, $d); + + // β = ||∇f(x(k+1))||² ∕ ||∇f(x(k))||² + $beta = $this->getBeta($thetaNew); + + // d(k+1) =–∇f(x(k+1)) + β(k).d(k) + $d = $this->getNewDirection($thetaNew, $beta, $d); + + // Save values for the next iteration + $oldTheta = $this->theta; + $this->costValues[] = $this->cost($thetaNew); + + $this->theta = $thetaNew; + if ($this->enableEarlyStop && $this->earlyStop($oldTheta)) { + break; + } + } + + $this->clear(); + + return $this->theta; + } + + /** + * Executes the callback function for the problem and returns + * sum of the gradient for all samples & targets. + * + * @param array $theta + * + * @return float + */ + protected function gradient(array $theta) + { + list($_, $gradient, $_) = parent::gradient($theta); + + return $gradient; + } + + /** + * Returns the value of f(x) for given solution + * + * @param array $theta + * + * @return float + */ + protected function cost(array $theta) + { + list($cost, $_, $_) = parent::gradient($theta); + + return array_sum($cost) / $this->sampleCount; + } + + /** + * Calculates alpha that minimizes the function f(θ + α.d) + * by performing a line search that does not rely upon the derivation. + * + * There are several alternatives for this function. For now, we + * prefer a method inspired from the bisection method for its simplicity. + * This algorithm attempts to find an optimum alpha value between 0.0001 and 0.01 + * + * Algorithm as follows: + * a) Probe a small alpha (0.0001) and calculate cost function + * b) Probe a larger alpha (0.01) and calculate cost function + * b-1) If cost function decreases, continue enlarging alpha + * b-2) If cost function increases, take the midpoint and try again + * + * @param float $d + * + * @return array + */ + protected function getAlpha(float $d) + { + $small = 0.0001 * $d; + $large = 0.01 * $d; + + // Obtain θ + α.d for two initial values, x0 and x1 + $x0 = mp::adds($this->theta, $small); + $x1 = mp::adds($this->theta, $large); + + $epsilon = 0.0001; + $iteration = 0; + do { + $fx1 = $this->cost($x1); + $fx0 = $this->cost($x0); + + // If the difference between two values is small enough + // then break the loop + if (abs($fx1 - $fx0) <= $epsilon) { + break; + } + + if ($fx1 < $fx0) { + $x0 = $x1; + $x1 = mp::adds($x1, 0.01); // Enlarge second + } else { + $x1 = mp::divs(mp::add($x1, $x0), 2.0); + } // Get to the midpoint + + $error = $fx1 / $this->dimensions; + } while ($error <= $epsilon || $iteration++ < 10); + + // Return α = θ / d + if ($d == 0) { + return $x1[0] - $this->theta[0]; + } + + return ($x1[0] - $this->theta[0]) / $d; + } + + /** + * Calculates new set of solutions with given alpha (for each θ(k)) and + * gradient direction. + * + * θ(k+1) = θ(k) + α.d + * + * @param float $alpha + * @param array $d + * + * return array + */ + protected function getNewTheta(float $alpha, array $d) + { + $theta = $this->theta; + + for ($i=0; $i < $this->dimensions + 1; $i++) { + if ($i == 0) { + $theta[$i] += $alpha * array_sum($d); + } else { + $sum = 0.0; + foreach ($this->samples as $si => $sample) { + $sum += $sample[$i - 1] * $d[$si] * $alpha; + } + + $theta[$i] += $sum; + } + } + + return $theta; + } + + /** + * Calculates new beta (β) for given set of solutions by using + * Fletcher–Reeves method. + * + * β = ||f(x(k+1))||² ∕ ||f(x(k))||² + * + * See: + * R. Fletcher and C. M. Reeves, "Function minimization by conjugate gradients", Comput. J. 7 (1964), 149–154. + * + * @param array $newTheta + * + * @return float + */ + protected function getBeta(array $newTheta) + { + $dNew = array_sum($this->gradient($newTheta)); + $dOld = array_sum($this->gradient($this->theta)) + 1e-100; + + return $dNew ** 2 / $dOld ** 2; + } + + /** + * Calculates the new conjugate direction + * + * d(k+1) =–∇f(x(k+1)) + β(k).d(k) + * + * @param array $theta + * @param float $beta + * @param array $d + * + * @return array + */ + protected function getNewDirection(array $theta, float $beta, array $d) + { + $grad = $this->gradient($theta); + + return mp::add(mp::muls($grad, -1), mp::muls($d, $beta)); + } +} + +/** + * Handles element-wise vector operations between vector-vector + * and vector-scalar variables + */ +class mp +{ + /** + * Element-wise multiplication of two vectors of the same size + * + * @param array $m1 + * @param array $m2 + * + * @return array + */ + public static function mul(array $m1, array $m2) + { + $res = []; + foreach ($m1 as $i => $val) { + $res[] = $val * $m2[$i]; + } + + return $res; + } + + /** + * Element-wise division of two vectors of the same size + * + * @param array $m1 + * @param array $m2 + * + * @return array + */ + public static function div(array $m1, array $m2) + { + $res = []; + foreach ($m1 as $i => $val) { + $res[] = $val / $m2[$i]; + } + + return $res; + } + + /** + * Element-wise addition of two vectors of the same size + * + * @param array $m1 + * @param array $m2 + * + * @return array + */ + public static function add(array $m1, array $m2, $mag = 1) + { + $res = []; + foreach ($m1 as $i => $val) { + $res[] = $val + $mag * $m2[$i]; + } + + return $res; + } + + /** + * Element-wise subtraction of two vectors of the same size + * + * @param array $m1 + * @param array $m2 + * + * @return array + */ + public static function sub(array $m1, array $m2) + { + return self::add($m1, $m2, -1); + } + + /** + * Element-wise multiplication of a vector with a scalar + * + * @param array $m1 + * @param float $m2 + * + * @return array + */ + public static function muls(array $m1, float $m2) + { + $res = []; + foreach ($m1 as $val) { + $res[] = $val * $m2; + } + + return $res; + } + + /** + * Element-wise division of a vector with a scalar + * + * @param array $m1 + * @param float $m2 + * + * @return array + */ + public static function divs(array $m1, float $m2) + { + $res = []; + foreach ($m1 as $val) { + $res[] = $val / ($m2 + 1e-32); + } + + return $res; + } + + /** + * Element-wise addition of a vector with a scalar + * + * @param array $m1 + * @param float $m2 + * + * @return array + */ + public static function adds(array $m1, float $m2, $mag = 1) + { + $res = []; + foreach ($m1 as $val) { + $res[] = $val + $mag * $m2; + } + + return $res; + } + + /** + * Element-wise subtraction of a vector with a scalar + * + * @param array $m1 + * @param float $m2 + * + * @return array + */ + public static function subs(array $m1, array $m2) + { + return self::adds($m1, $m2, -1); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/GD.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/GD.php new file mode 100644 index 00000000000..8974c8e769c --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/GD.php @@ -0,0 +1,121 @@ +samples = $samples; + $this->targets = $targets; + $this->gradientCb = $gradientCb; + $this->sampleCount = count($this->samples); + + // Batch learning is executed: + $currIter = 0; + $this->costValues = []; + while ($this->maxIterations > $currIter++) { + $theta = $this->theta; + + // Calculate update terms for each sample + list($errors, $updates, $totalPenalty) = $this->gradient($theta); + + $this->updateWeightsWithUpdates($updates, $totalPenalty); + + $this->costValues[] = array_sum($errors)/$this->sampleCount; + + if ($this->earlyStop($theta)) { + break; + } + } + + $this->clear(); + + return $this->theta; + } + + /** + * Calculates gradient, cost function and penalty term for each sample + * then returns them as an array of values + * + * @param array $theta + * + * @return array + */ + protected function gradient(array $theta) + { + $costs = []; + $gradient= []; + $totalPenalty = 0; + + foreach ($this->samples as $index => $sample) { + $target = $this->targets[$index]; + + $result = ($this->gradientCb)($theta, $sample, $target); + list($cost, $grad, $penalty) = array_pad($result, 3, 0); + + $costs[] = $cost; + $gradient[]= $grad; + $totalPenalty += $penalty; + } + + $totalPenalty /= $this->sampleCount; + + return [$costs, $gradient, $totalPenalty]; + } + + /** + * @param array $updates + * @param float $penalty + */ + protected function updateWeightsWithUpdates(array $updates, float $penalty) + { + // Updates all weights at once + for ($i=0; $i <= $this->dimensions; $i++) { + if ($i == 0) { + $this->theta[0] -= $this->learningRate * array_sum($updates); + } else { + $col = array_column($this->samples, $i - 1); + + $error = 0; + foreach ($col as $index => $val) { + $error += $val * $updates[$index]; + } + + $this->theta[$i] -= $this->learningRate * + ($error + $penalty * $this->theta[$i]); + } + } + } + + /** + * Clears the optimizer internal vars after the optimization process. + * + * @return void + */ + protected function clear() + { + $this->sampleCount = null; + parent::clear(); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/Optimizer.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/Optimizer.php new file mode 100644 index 00000000000..9ef4c4d0949 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/Optimizer.php @@ -0,0 +1,61 @@ +dimensions = $dimensions; + + // Inits the weights randomly + $this->theta = []; + for ($i=0; $i < $this->dimensions; $i++) { + $this->theta[] = rand() / (float) getrandmax(); + } + } + + /** + * Sets the weights manually + * + * @param array $theta + */ + public function setInitialTheta(array $theta) + { + if (count($theta) != $this->dimensions) { + throw new \Exception("Number of values in the weights array should be $this->dimensions"); + } + + $this->theta = $theta; + + return $this; + } + + /** + * Executes the optimization with the given samples & targets + * and returns the weights + * + */ + abstract protected function runOptimization(array $samples, array $targets, \Closure $gradientCb); +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/StochasticGD.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/StochasticGD.php new file mode 100644 index 00000000000..e9e318a8a5f --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/Optimizer/StochasticGD.php @@ -0,0 +1,285 @@ + + * + * Larger values of lr may overshoot the optimum or even cause divergence + * while small values slows down the convergence and increases the time + * required for the training + * + * @var float + */ + protected $learningRate = 0.001; + + /** + * Minimum amount of change in the weights and error values + * between iterations that needs to be obtained to continue the training + * + * @var float + */ + protected $threshold = 1e-4; + + /** + * Enable/Disable early stopping by checking the weight & cost values + * to see whether they changed large enough to continue the optimization + * + * @var bool + */ + protected $enableEarlyStop = true; + /** + * List of values obtained by evaluating the cost function at each iteration + * of the algorithm + * + * @var array + */ + protected $costValues= []; + + /** + * Initializes the SGD optimizer for the given number of dimensions + * + * @param int $dimensions + */ + public function __construct(int $dimensions) + { + // Add one more dimension for the bias + parent::__construct($dimensions + 1); + + $this->dimensions = $dimensions; + } + + /** + * Sets minimum value for the change in the theta values + * between iterations to continue the iterations.
+ * + * If change in the theta is less than given value then the + * algorithm will stop training + * + * @param float $threshold + * + * @return $this + */ + public function setChangeThreshold(float $threshold = 1e-5) + { + $this->threshold = $threshold; + + return $this; + } + + /** + * Enable/Disable early stopping by checking at each iteration + * whether changes in theta or cost value are not large enough + * + * @param bool $enable + * + * @return $this + */ + public function setEarlyStop(bool $enable = true) + { + $this->enableEarlyStop = $enable; + + return $this; + } + + /** + * @param float $learningRate + * + * @return $this + */ + public function setLearningRate(float $learningRate) + { + $this->learningRate = $learningRate; + + return $this; + } + + /** + * @param int $maxIterations + * + * @return $this + */ + public function setMaxIterations(int $maxIterations) + { + $this->maxIterations = $maxIterations; + + return $this; + } + + /** + * Optimization procedure finds the unknow variables for the equation A.Ï´ = y + * for the given samples (A) and targets (y).
+ * + * The cost function to minimize and the gradient of the function are to be + * handled by the callback function provided as the third parameter of the method. + * + * @param array $samples + * @param array $targets + * @param \Closure $gradientCb + * + * @return array + */ + public function runOptimization(array $samples, array $targets, \Closure $gradientCb) + { + $this->samples = $samples; + $this->targets = $targets; + $this->gradientCb = $gradientCb; + + $currIter = 0; + $bestTheta = null; + $bestScore = 0.0; + $bestWeightIter = 0; + $this->costValues = []; + + while ($this->maxIterations > $currIter++) { + $theta = $this->theta; + + // Update the guess + $cost = $this->updateTheta(); + + // Save the best theta in the "pocket" so that + // any future set of theta worse than this will be disregarded + if ($bestTheta == null || $cost <= $bestScore) { + $bestTheta = $theta; + $bestScore = $cost; + $bestWeightIter = $currIter; + } + + // Add the cost value for this iteration to the list + $this->costValues[] = $cost; + + // Check for early stop + if ($this->enableEarlyStop && $this->earlyStop($theta)) { + break; + } + } + + $this->clear(); + + // Solution in the pocket is better than or equal to the last state + // so, we use this solution + return $this->theta = $bestTheta; + } + + /** + * @return float + */ + protected function updateTheta() + { + $jValue = 0.0; + $theta = $this->theta; + + foreach ($this->samples as $index => $sample) { + $target = $this->targets[$index]; + + $result = ($this->gradientCb)($theta, $sample, $target); + + list($error, $gradient, $penalty) = array_pad($result, 3, 0); + + // Update bias + $this->theta[0] -= $this->learningRate * $gradient; + + // Update other values + for ($i=1; $i <= $this->dimensions; $i++) { + $this->theta[$i] -= $this->learningRate * + ($gradient * $sample[$i - 1] + $penalty * $this->theta[$i]); + } + + // Sum error rate + $jValue += $error; + } + + return $jValue / count($this->samples); + } + + /** + * Checks if the optimization is not effective enough and can be stopped + * in case large enough changes in the solution do not happen + * + * @param array $oldTheta + * + * @return boolean + */ + protected function earlyStop($oldTheta) + { + // Check for early stop: No change larger than threshold (default 1e-5) + $diff = array_map( + function ($w1, $w2) { + return abs($w1 - $w2) > $this->threshold ? 1 : 0; + }, + $oldTheta, $this->theta); + + if (array_sum($diff) == 0) { + return true; + } + + // Check if the last two cost values are almost the same + $costs = array_slice($this->costValues, -2); + if (count($costs) == 2 && abs($costs[1] - $costs[0]) < $this->threshold) { + return true; + } + + return false; + } + + /** + * Returns the list of cost values for each iteration executed in + * last run of the optimization + * + * @return array + */ + public function getCostValues() + { + return $this->costValues; + } + + /** + * Clears the optimizer internal vars after the optimization process. + * + * @return void + */ + protected function clear() + { + $this->samples = []; + $this->targets = []; + $this->gradientCb = null; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/Predictable.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/Predictable.php new file mode 100644 index 00000000000..097edaabeab --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/Predictable.php @@ -0,0 +1,34 @@ +predictSample($samples); + } else { + $predicted = []; + foreach ($samples as $index => $sample) { + $predicted[$index] = $this->predictSample($sample); + } + } + + return $predicted; + } + + /** + * @param array $sample + * + * @return mixed + */ + abstract protected function predictSample(array $sample); +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Helper/Trainable.php b/lib/mlbackend/php/phpml/src/Phpml/Helper/Trainable.php new file mode 100644 index 00000000000..3d011ac47f8 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Helper/Trainable.php @@ -0,0 +1,28 @@ +samples = array_merge($this->samples, $samples); + $this->targets = array_merge($this->targets, $targets); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/IncrementalEstimator.php b/lib/mlbackend/php/phpml/src/Phpml/IncrementalEstimator.php new file mode 100644 index 00000000000..fc6912d1109 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/IncrementalEstimator.php @@ -0,0 +1,16 @@ + $val) { + $distance += ($val - $b[$i]) ** 2; + } + + return sqrt((float) $distance); + } + + /** + * Square of Euclidean distance + * + * @param array $a + * @param array $b + * + * @return float + */ + public function sqDistance(array $a, array $b): float + { + return $this->distance($a, $b) ** 2; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Distance/Manhattan.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Distance/Manhattan.php new file mode 100644 index 00000000000..b6f6eb8bdd2 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Distance/Manhattan.php @@ -0,0 +1,35 @@ +lambda = $lambda; + } + + /** + * @param array $a + * @param array $b + * + * @return float + * + * @throws InvalidArgumentException + */ + public function distance(array $a, array $b): float + { + if (count($a) !== count($b)) { + throw InvalidArgumentException::arraySizeNotMatch(); + } + + $distance = 0; + $count = count($a); + + for ($i = 0; $i < $count; ++$i) { + $distance += pow(abs($a[$i] - $b[$i]), $this->lambda); + } + + return (float)pow($distance, 1 / $this->lambda); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Kernel.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Kernel.php new file mode 100644 index 00000000000..6d1461fde0c --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Kernel.php @@ -0,0 +1,16 @@ +gamma = $gamma; + } + + /** + * @param float $a + * @param float $b + * + * @return float + */ + public function compute($a, $b) + { + $score = 2 * Product::scalar($a, $b); + $squares = Product::scalar($a, $a) + Product::scalar($b, $b); + $result = exp(-$this->gamma * ($squares - $score)); + + return $result; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/LinearAlgebra/EigenvalueDecomposition.php b/lib/mlbackend/php/phpml/src/Phpml/Math/LinearAlgebra/EigenvalueDecomposition.php new file mode 100644 index 00000000000..27557bbd83a --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/LinearAlgebra/EigenvalueDecomposition.php @@ -0,0 +1,890 @@ +d = $this->V[$this->n-1]; + // Householder reduction to tridiagonal form. + for ($i = $this->n-1; $i > 0; --$i) { + $i_ = $i -1; + // Scale to avoid under/overflow. + $h = $scale = 0.0; + $scale += array_sum(array_map('abs', $this->d)); + if ($scale == 0.0) { + $this->e[$i] = $this->d[$i_]; + $this->d = array_slice($this->V[$i_], 0, $i_); + for ($j = 0; $j < $i; ++$j) { + $this->V[$j][$i] = $this->V[$i][$j] = 0.0; + } + } else { + // Generate Householder vector. + for ($k = 0; $k < $i; ++$k) { + $this->d[$k] /= $scale; + $h += pow($this->d[$k], 2); + } + $f = $this->d[$i_]; + $g = sqrt($h); + if ($f > 0) { + $g = -$g; + } + $this->e[$i] = $scale * $g; + $h = $h - $f * $g; + $this->d[$i_] = $f - $g; + for ($j = 0; $j < $i; ++$j) { + $this->e[$j] = 0.0; + } + // Apply similarity transformation to remaining columns. + for ($j = 0; $j < $i; ++$j) { + $f = $this->d[$j]; + $this->V[$j][$i] = $f; + $g = $this->e[$j] + $this->V[$j][$j] * $f; + for ($k = $j+1; $k <= $i_; ++$k) { + $g += $this->V[$k][$j] * $this->d[$k]; + $this->e[$k] += $this->V[$k][$j] * $f; + } + $this->e[$j] = $g; + } + $f = 0.0; + for ($j = 0; $j < $i; ++$j) { + if ($h === 0) { + $h = 1e-20; + } + $this->e[$j] /= $h; + $f += $this->e[$j] * $this->d[$j]; + } + $hh = $f / (2 * $h); + for ($j=0; $j < $i; ++$j) { + $this->e[$j] -= $hh * $this->d[$j]; + } + for ($j = 0; $j < $i; ++$j) { + $f = $this->d[$j]; + $g = $this->e[$j]; + for ($k = $j; $k <= $i_; ++$k) { + $this->V[$k][$j] -= ($f * $this->e[$k] + $g * $this->d[$k]); + } + $this->d[$j] = $this->V[$i-1][$j]; + $this->V[$i][$j] = 0.0; + } + } + $this->d[$i] = $h; + } + + // Accumulate transformations. + for ($i = 0; $i < $this->n-1; ++$i) { + $this->V[$this->n-1][$i] = $this->V[$i][$i]; + $this->V[$i][$i] = 1.0; + $h = $this->d[$i+1]; + if ($h != 0.0) { + for ($k = 0; $k <= $i; ++$k) { + $this->d[$k] = $this->V[$k][$i+1] / $h; + } + for ($j = 0; $j <= $i; ++$j) { + $g = 0.0; + for ($k = 0; $k <= $i; ++$k) { + $g += $this->V[$k][$i+1] * $this->V[$k][$j]; + } + for ($k = 0; $k <= $i; ++$k) { + $this->V[$k][$j] -= $g * $this->d[$k]; + } + } + } + for ($k = 0; $k <= $i; ++$k) { + $this->V[$k][$i+1] = 0.0; + } + } + + $this->d = $this->V[$this->n-1]; + $this->V[$this->n-1] = array_fill(0, $j, 0.0); + $this->V[$this->n-1][$this->n-1] = 1.0; + $this->e[0] = 0.0; + } + + + /** + * Symmetric tridiagonal QL algorithm. + * + * This is derived from the Algol procedures tql2, by + * Bowdler, Martin, Reinsch, and Wilkinson, Handbook for + * Auto. Comp., Vol.ii-Linear Algebra, and the corresponding + * Fortran subroutine in EISPACK. + */ + private function tql2() + { + for ($i = 1; $i < $this->n; ++$i) { + $this->e[$i-1] = $this->e[$i]; + } + $this->e[$this->n-1] = 0.0; + $f = 0.0; + $tst1 = 0.0; + $eps = pow(2.0, -52.0); + + for ($l = 0; $l < $this->n; ++$l) { + // Find small subdiagonal element + $tst1 = max($tst1, abs($this->d[$l]) + abs($this->e[$l])); + $m = $l; + while ($m < $this->n) { + if (abs($this->e[$m]) <= $eps * $tst1) { + break; + } + ++$m; + } + // If m == l, $this->d[l] is an eigenvalue, + // otherwise, iterate. + if ($m > $l) { + $iter = 0; + do { + // Could check iteration count here. + $iter += 1; + // Compute implicit shift + $g = $this->d[$l]; + $p = ($this->d[$l+1] - $g) / (2.0 * $this->e[$l]); + $r = hypot($p, 1.0); + if ($p < 0) { + $r *= -1; + } + $this->d[$l] = $this->e[$l] / ($p + $r); + $this->d[$l+1] = $this->e[$l] * ($p + $r); + $dl1 = $this->d[$l+1]; + $h = $g - $this->d[$l]; + for ($i = $l + 2; $i < $this->n; ++$i) { + $this->d[$i] -= $h; + } + $f += $h; + // Implicit QL transformation. + $p = $this->d[$m]; + $c = 1.0; + $c2 = $c3 = $c; + $el1 = $this->e[$l + 1]; + $s = $s2 = 0.0; + for ($i = $m-1; $i >= $l; --$i) { + $c3 = $c2; + $c2 = $c; + $s2 = $s; + $g = $c * $this->e[$i]; + $h = $c * $p; + $r = hypot($p, $this->e[$i]); + $this->e[$i+1] = $s * $r; + $s = $this->e[$i] / $r; + $c = $p / $r; + $p = $c * $this->d[$i] - $s * $g; + $this->d[$i+1] = $h + $s * ($c * $g + $s * $this->d[$i]); + // Accumulate transformation. + for ($k = 0; $k < $this->n; ++$k) { + $h = $this->V[$k][$i+1]; + $this->V[$k][$i+1] = $s * $this->V[$k][$i] + $c * $h; + $this->V[$k][$i] = $c * $this->V[$k][$i] - $s * $h; + } + } + $p = -$s * $s2 * $c3 * $el1 * $this->e[$l] / $dl1; + $this->e[$l] = $s * $p; + $this->d[$l] = $c * $p; + // Check for convergence. + } while (abs($this->e[$l]) > $eps * $tst1); + } + $this->d[$l] = $this->d[$l] + $f; + $this->e[$l] = 0.0; + } + + // Sort eigenvalues and corresponding vectors. + for ($i = 0; $i < $this->n - 1; ++$i) { + $k = $i; + $p = $this->d[$i]; + for ($j = $i+1; $j < $this->n; ++$j) { + if ($this->d[$j] < $p) { + $k = $j; + $p = $this->d[$j]; + } + } + if ($k != $i) { + $this->d[$k] = $this->d[$i]; + $this->d[$i] = $p; + for ($j = 0; $j < $this->n; ++$j) { + $p = $this->V[$j][$i]; + $this->V[$j][$i] = $this->V[$j][$k]; + $this->V[$j][$k] = $p; + } + } + } + } + + + /** + * Nonsymmetric reduction to Hessenberg form. + * + * This is derived from the Algol procedures orthes and ortran, + * by Martin and Wilkinson, Handbook for Auto. Comp., + * Vol.ii-Linear Algebra, and the corresponding + * Fortran subroutines in EISPACK. + */ + private function orthes() + { + $low = 0; + $high = $this->n-1; + + for ($m = $low+1; $m <= $high-1; ++$m) { + // Scale column. + $scale = 0.0; + for ($i = $m; $i <= $high; ++$i) { + $scale = $scale + abs($this->H[$i][$m-1]); + } + if ($scale != 0.0) { + // Compute Householder transformation. + $h = 0.0; + for ($i = $high; $i >= $m; --$i) { + $this->ort[$i] = $this->H[$i][$m-1] / $scale; + $h += $this->ort[$i] * $this->ort[$i]; + } + $g = sqrt($h); + if ($this->ort[$m] > 0) { + $g *= -1; + } + $h -= $this->ort[$m] * $g; + $this->ort[$m] -= $g; + // Apply Householder similarity transformation + // H = (I -u * u' / h) * H * (I -u * u') / h) + for ($j = $m; $j < $this->n; ++$j) { + $f = 0.0; + for ($i = $high; $i >= $m; --$i) { + $f += $this->ort[$i] * $this->H[$i][$j]; + } + $f /= $h; + for ($i = $m; $i <= $high; ++$i) { + $this->H[$i][$j] -= $f * $this->ort[$i]; + } + } + for ($i = 0; $i <= $high; ++$i) { + $f = 0.0; + for ($j = $high; $j >= $m; --$j) { + $f += $this->ort[$j] * $this->H[$i][$j]; + } + $f = $f / $h; + for ($j = $m; $j <= $high; ++$j) { + $this->H[$i][$j] -= $f * $this->ort[$j]; + } + } + $this->ort[$m] = $scale * $this->ort[$m]; + $this->H[$m][$m-1] = $scale * $g; + } + } + + // Accumulate transformations (Algol's ortran). + for ($i = 0; $i < $this->n; ++$i) { + for ($j = 0; $j < $this->n; ++$j) { + $this->V[$i][$j] = ($i == $j ? 1.0 : 0.0); + } + } + for ($m = $high-1; $m >= $low+1; --$m) { + if ($this->H[$m][$m-1] != 0.0) { + for ($i = $m+1; $i <= $high; ++$i) { + $this->ort[$i] = $this->H[$i][$m-1]; + } + for ($j = $m; $j <= $high; ++$j) { + $g = 0.0; + for ($i = $m; $i <= $high; ++$i) { + $g += $this->ort[$i] * $this->V[$i][$j]; + } + // Double division avoids possible underflow + $g = ($g / $this->ort[$m]) / $this->H[$m][$m-1]; + for ($i = $m; $i <= $high; ++$i) { + $this->V[$i][$j] += $g * $this->ort[$i]; + } + } + } + } + } + + + /** + * Performs complex division. + */ + private function cdiv($xr, $xi, $yr, $yi) + { + if (abs($yr) > abs($yi)) { + $r = $yi / $yr; + $d = $yr + $r * $yi; + $this->cdivr = ($xr + $r * $xi) / $d; + $this->cdivi = ($xi - $r * $xr) / $d; + } else { + $r = $yr / $yi; + $d = $yi + $r * $yr; + $this->cdivr = ($r * $xr + $xi) / $d; + $this->cdivi = ($r * $xi - $xr) / $d; + } + } + + + /** + * Nonsymmetric reduction from Hessenberg to real Schur form. + * + * Code is derived from the Algol procedure hqr2, + * by Martin and Wilkinson, Handbook for Auto. Comp., + * Vol.ii-Linear Algebra, and the corresponding + * Fortran subroutine in EISPACK. + */ + private function hqr2() + { + // Initialize + $nn = $this->n; + $n = $nn - 1; + $low = 0; + $high = $nn - 1; + $eps = pow(2.0, -52.0); + $exshift = 0.0; + $p = $q = $r = $s = $z = 0; + // Store roots isolated by balanc and compute matrix norm + $norm = 0.0; + + for ($i = 0; $i < $nn; ++$i) { + if (($i < $low) or ($i > $high)) { + $this->d[$i] = $this->H[$i][$i]; + $this->e[$i] = 0.0; + } + for ($j = max($i-1, 0); $j < $nn; ++$j) { + $norm = $norm + abs($this->H[$i][$j]); + } + } + + // Outer loop over eigenvalue index + $iter = 0; + while ($n >= $low) { + // Look for single small sub-diagonal element + $l = $n; + while ($l > $low) { + $s = abs($this->H[$l-1][$l-1]) + abs($this->H[$l][$l]); + if ($s == 0.0) { + $s = $norm; + } + if (abs($this->H[$l][$l-1]) < $eps * $s) { + break; + } + --$l; + } + // Check for convergence + // One root found + if ($l == $n) { + $this->H[$n][$n] = $this->H[$n][$n] + $exshift; + $this->d[$n] = $this->H[$n][$n]; + $this->e[$n] = 0.0; + --$n; + $iter = 0; + // Two roots found + } elseif ($l == $n-1) { + $w = $this->H[$n][$n-1] * $this->H[$n-1][$n]; + $p = ($this->H[$n-1][$n-1] - $this->H[$n][$n]) / 2.0; + $q = $p * $p + $w; + $z = sqrt(abs($q)); + $this->H[$n][$n] = $this->H[$n][$n] + $exshift; + $this->H[$n-1][$n-1] = $this->H[$n-1][$n-1] + $exshift; + $x = $this->H[$n][$n]; + // Real pair + if ($q >= 0) { + if ($p >= 0) { + $z = $p + $z; + } else { + $z = $p - $z; + } + $this->d[$n-1] = $x + $z; + $this->d[$n] = $this->d[$n-1]; + if ($z != 0.0) { + $this->d[$n] = $x - $w / $z; + } + $this->e[$n-1] = 0.0; + $this->e[$n] = 0.0; + $x = $this->H[$n][$n-1]; + $s = abs($x) + abs($z); + $p = $x / $s; + $q = $z / $s; + $r = sqrt($p * $p + $q * $q); + $p = $p / $r; + $q = $q / $r; + // Row modification + for ($j = $n-1; $j < $nn; ++$j) { + $z = $this->H[$n-1][$j]; + $this->H[$n-1][$j] = $q * $z + $p * $this->H[$n][$j]; + $this->H[$n][$j] = $q * $this->H[$n][$j] - $p * $z; + } + // Column modification + for ($i = 0; $i <= $n; ++$i) { + $z = $this->H[$i][$n-1]; + $this->H[$i][$n-1] = $q * $z + $p * $this->H[$i][$n]; + $this->H[$i][$n] = $q * $this->H[$i][$n] - $p * $z; + } + // Accumulate transformations + for ($i = $low; $i <= $high; ++$i) { + $z = $this->V[$i][$n-1]; + $this->V[$i][$n-1] = $q * $z + $p * $this->V[$i][$n]; + $this->V[$i][$n] = $q * $this->V[$i][$n] - $p * $z; + } + // Complex pair + } else { + $this->d[$n-1] = $x + $p; + $this->d[$n] = $x + $p; + $this->e[$n-1] = $z; + $this->e[$n] = -$z; + } + $n = $n - 2; + $iter = 0; + // No convergence yet + } else { + // Form shift + $x = $this->H[$n][$n]; + $y = 0.0; + $w = 0.0; + if ($l < $n) { + $y = $this->H[$n-1][$n-1]; + $w = $this->H[$n][$n-1] * $this->H[$n-1][$n]; + } + // Wilkinson's original ad hoc shift + if ($iter == 10) { + $exshift += $x; + for ($i = $low; $i <= $n; ++$i) { + $this->H[$i][$i] -= $x; + } + $s = abs($this->H[$n][$n-1]) + abs($this->H[$n-1][$n-2]); + $x = $y = 0.75 * $s; + $w = -0.4375 * $s * $s; + } + // MATLAB's new ad hoc shift + if ($iter == 30) { + $s = ($y - $x) / 2.0; + $s = $s * $s + $w; + if ($s > 0) { + $s = sqrt($s); + if ($y < $x) { + $s = -$s; + } + $s = $x - $w / (($y - $x) / 2.0 + $s); + for ($i = $low; $i <= $n; ++$i) { + $this->H[$i][$i] -= $s; + } + $exshift += $s; + $x = $y = $w = 0.964; + } + } + // Could check iteration count here. + $iter = $iter + 1; + // Look for two consecutive small sub-diagonal elements + $m = $n - 2; + while ($m >= $l) { + $z = $this->H[$m][$m]; + $r = $x - $z; + $s = $y - $z; + $p = ($r * $s - $w) / $this->H[$m+1][$m] + $this->H[$m][$m+1]; + $q = $this->H[$m+1][$m+1] - $z - $r - $s; + $r = $this->H[$m+2][$m+1]; + $s = abs($p) + abs($q) + abs($r); + $p = $p / $s; + $q = $q / $s; + $r = $r / $s; + if ($m == $l) { + break; + } + if (abs($this->H[$m][$m-1]) * (abs($q) + abs($r)) < + $eps * (abs($p) * (abs($this->H[$m-1][$m-1]) + abs($z) + abs($this->H[$m+1][$m+1])))) { + break; + } + --$m; + } + for ($i = $m + 2; $i <= $n; ++$i) { + $this->H[$i][$i-2] = 0.0; + if ($i > $m+2) { + $this->H[$i][$i-3] = 0.0; + } + } + // Double QR step involving rows l:n and columns m:n + for ($k = $m; $k <= $n-1; ++$k) { + $notlast = ($k != $n-1); + if ($k != $m) { + $p = $this->H[$k][$k-1]; + $q = $this->H[$k+1][$k-1]; + $r = ($notlast ? $this->H[$k+2][$k-1] : 0.0); + $x = abs($p) + abs($q) + abs($r); + if ($x != 0.0) { + $p = $p / $x; + $q = $q / $x; + $r = $r / $x; + } + } + if ($x == 0.0) { + break; + } + $s = sqrt($p * $p + $q * $q + $r * $r); + if ($p < 0) { + $s = -$s; + } + if ($s != 0) { + if ($k != $m) { + $this->H[$k][$k-1] = -$s * $x; + } elseif ($l != $m) { + $this->H[$k][$k-1] = -$this->H[$k][$k-1]; + } + $p = $p + $s; + $x = $p / $s; + $y = $q / $s; + $z = $r / $s; + $q = $q / $p; + $r = $r / $p; + // Row modification + for ($j = $k; $j < $nn; ++$j) { + $p = $this->H[$k][$j] + $q * $this->H[$k+1][$j]; + if ($notlast) { + $p = $p + $r * $this->H[$k+2][$j]; + $this->H[$k+2][$j] = $this->H[$k+2][$j] - $p * $z; + } + $this->H[$k][$j] = $this->H[$k][$j] - $p * $x; + $this->H[$k+1][$j] = $this->H[$k+1][$j] - $p * $y; + } + // Column modification + for ($i = 0; $i <= min($n, $k+3); ++$i) { + $p = $x * $this->H[$i][$k] + $y * $this->H[$i][$k+1]; + if ($notlast) { + $p = $p + $z * $this->H[$i][$k+2]; + $this->H[$i][$k+2] = $this->H[$i][$k+2] - $p * $r; + } + $this->H[$i][$k] = $this->H[$i][$k] - $p; + $this->H[$i][$k+1] = $this->H[$i][$k+1] - $p * $q; + } + // Accumulate transformations + for ($i = $low; $i <= $high; ++$i) { + $p = $x * $this->V[$i][$k] + $y * $this->V[$i][$k+1]; + if ($notlast) { + $p = $p + $z * $this->V[$i][$k+2]; + $this->V[$i][$k+2] = $this->V[$i][$k+2] - $p * $r; + } + $this->V[$i][$k] = $this->V[$i][$k] - $p; + $this->V[$i][$k+1] = $this->V[$i][$k+1] - $p * $q; + } + } // ($s != 0) + } // k loop + } // check convergence + } // while ($n >= $low) + + // Backsubstitute to find vectors of upper triangular form + if ($norm == 0.0) { + return; + } + + for ($n = $nn-1; $n >= 0; --$n) { + $p = $this->d[$n]; + $q = $this->e[$n]; + // Real vector + if ($q == 0) { + $l = $n; + $this->H[$n][$n] = 1.0; + for ($i = $n-1; $i >= 0; --$i) { + $w = $this->H[$i][$i] - $p; + $r = 0.0; + for ($j = $l; $j <= $n; ++$j) { + $r = $r + $this->H[$i][$j] * $this->H[$j][$n]; + } + if ($this->e[$i] < 0.0) { + $z = $w; + $s = $r; + } else { + $l = $i; + if ($this->e[$i] == 0.0) { + if ($w != 0.0) { + $this->H[$i][$n] = -$r / $w; + } else { + $this->H[$i][$n] = -$r / ($eps * $norm); + } + // Solve real equations + } else { + $x = $this->H[$i][$i+1]; + $y = $this->H[$i+1][$i]; + $q = ($this->d[$i] - $p) * ($this->d[$i] - $p) + $this->e[$i] * $this->e[$i]; + $t = ($x * $s - $z * $r) / $q; + $this->H[$i][$n] = $t; + if (abs($x) > abs($z)) { + $this->H[$i+1][$n] = (-$r - $w * $t) / $x; + } else { + $this->H[$i+1][$n] = (-$s - $y * $t) / $z; + } + } + // Overflow control + $t = abs($this->H[$i][$n]); + if (($eps * $t) * $t > 1) { + for ($j = $i; $j <= $n; ++$j) { + $this->H[$j][$n] = $this->H[$j][$n] / $t; + } + } + } + } + // Complex vector + } elseif ($q < 0) { + $l = $n-1; + // Last vector component imaginary so matrix is triangular + if (abs($this->H[$n][$n-1]) > abs($this->H[$n-1][$n])) { + $this->H[$n-1][$n-1] = $q / $this->H[$n][$n-1]; + $this->H[$n-1][$n] = -($this->H[$n][$n] - $p) / $this->H[$n][$n-1]; + } else { + $this->cdiv(0.0, -$this->H[$n-1][$n], $this->H[$n-1][$n-1] - $p, $q); + $this->H[$n-1][$n-1] = $this->cdivr; + $this->H[$n-1][$n] = $this->cdivi; + } + $this->H[$n][$n-1] = 0.0; + $this->H[$n][$n] = 1.0; + for ($i = $n-2; $i >= 0; --$i) { + // double ra,sa,vr,vi; + $ra = 0.0; + $sa = 0.0; + for ($j = $l; $j <= $n; ++$j) { + $ra = $ra + $this->H[$i][$j] * $this->H[$j][$n-1]; + $sa = $sa + $this->H[$i][$j] * $this->H[$j][$n]; + } + $w = $this->H[$i][$i] - $p; + if ($this->e[$i] < 0.0) { + $z = $w; + $r = $ra; + $s = $sa; + } else { + $l = $i; + if ($this->e[$i] == 0) { + $this->cdiv(-$ra, -$sa, $w, $q); + $this->H[$i][$n-1] = $this->cdivr; + $this->H[$i][$n] = $this->cdivi; + } else { + // Solve complex equations + $x = $this->H[$i][$i+1]; + $y = $this->H[$i+1][$i]; + $vr = ($this->d[$i] - $p) * ($this->d[$i] - $p) + $this->e[$i] * $this->e[$i] - $q * $q; + $vi = ($this->d[$i] - $p) * 2.0 * $q; + if ($vr == 0.0 & $vi == 0.0) { + $vr = $eps * $norm * (abs($w) + abs($q) + abs($x) + abs($y) + abs($z)); + } + $this->cdiv($x * $r - $z * $ra + $q * $sa, $x * $s - $z * $sa - $q * $ra, $vr, $vi); + $this->H[$i][$n-1] = $this->cdivr; + $this->H[$i][$n] = $this->cdivi; + if (abs($x) > (abs($z) + abs($q))) { + $this->H[$i+1][$n-1] = (-$ra - $w * $this->H[$i][$n-1] + $q * $this->H[$i][$n]) / $x; + $this->H[$i+1][$n] = (-$sa - $w * $this->H[$i][$n] - $q * $this->H[$i][$n-1]) / $x; + } else { + $this->cdiv(-$r - $y * $this->H[$i][$n-1], -$s - $y * $this->H[$i][$n], $z, $q); + $this->H[$i+1][$n-1] = $this->cdivr; + $this->H[$i+1][$n] = $this->cdivi; + } + } + // Overflow control + $t = max(abs($this->H[$i][$n-1]), abs($this->H[$i][$n])); + if (($eps * $t) * $t > 1) { + for ($j = $i; $j <= $n; ++$j) { + $this->H[$j][$n-1] = $this->H[$j][$n-1] / $t; + $this->H[$j][$n] = $this->H[$j][$n] / $t; + } + } + } // end else + } // end for + } // end else for complex case + } // end for + + // Vectors of isolated roots + for ($i = 0; $i < $nn; ++$i) { + if ($i < $low | $i > $high) { + for ($j = $i; $j < $nn; ++$j) { + $this->V[$i][$j] = $this->H[$i][$j]; + } + } + } + + // Back transformation to get eigenvectors of original matrix + for ($j = $nn-1; $j >= $low; --$j) { + for ($i = $low; $i <= $high; ++$i) { + $z = 0.0; + for ($k = $low; $k <= min($j, $high); ++$k) { + $z = $z + $this->V[$i][$k] * $this->H[$k][$j]; + } + $this->V[$i][$j] = $z; + } + } + } // end hqr2 + + + /** + * Constructor: Check for symmetry, then construct the eigenvalue decomposition + * + * @param array $Arg + */ + public function __construct(array $Arg) + { + $this->A = $Arg; + $this->n = count($Arg[0]); + + $issymmetric = true; + for ($j = 0; ($j < $this->n) & $issymmetric; ++$j) { + for ($i = 0; ($i < $this->n) & $issymmetric; ++$i) { + $issymmetric = ($this->A[$i][$j] == $this->A[$j][$i]); + } + } + + if ($issymmetric) { + $this->V = $this->A; + // Tridiagonalize. + $this->tred2(); + // Diagonalize. + $this->tql2(); + } else { + $this->H = $this->A; + $this->ort = []; + // Reduce to Hessenberg form. + $this->orthes(); + // Reduce Hessenberg to real Schur form. + $this->hqr2(); + } + } + + /** + * Return the eigenvector matrix + * + * @access public + * @return array + */ + public function getEigenvectors() + { + $vectors = $this->V; + + // Always return the eigenvectors of length 1.0 + $vectors = new Matrix($vectors); + $vectors = array_map(function ($vect) { + $sum = 0; + for ($i=0; $itranspose()->toArray()); + + return $vectors; + } + + + /** + * Return the real parts of the eigenvalues
+ * d = real(diag(D)); + * + * @return array + */ + public function getRealEigenvalues() + { + return $this->d; + } + + + /** + * Return the imaginary parts of the eigenvalues
+ * d = imag(diag(D)) + * + * @return array + */ + public function getImagEigenvalues() + { + return $this->e; + } + + + /** + * Return the block diagonal eigenvalue matrix + * + * @return array + */ + public function getDiagonalEigenvalues() + { + for ($i = 0; $i < $this->n; ++$i) { + $D[$i] = array_fill(0, $this->n, 0.0); + $D[$i][$i] = $this->d[$i]; + if ($this->e[$i] == 0) { + continue; + } + $o = ($this->e[$i] > 0) ? $i + 1 : $i - 1; + $D[$i][$o] = $this->e[$i]; + } + return $D; + } +} // class EigenvalueDecomposition diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Matrix.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Matrix.php new file mode 100644 index 00000000000..25101f3f4ab --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Matrix.php @@ -0,0 +1,385 @@ +rows = 1; + $this->columns = count($matrix); + $matrix = [$matrix]; + } else { + $this->rows = count($matrix); + $this->columns = count($matrix[0]); + } + + if ($validate) { + for ($i = 0; $i < $this->rows; ++$i) { + if (count($matrix[$i]) !== $this->columns) { + throw InvalidArgumentException::matrixDimensionsDidNotMatch(); + } + } + } + + $this->matrix = $matrix; + } + + /** + * @param array $array + * + * @return Matrix + */ + public static function fromFlatArray(array $array) + { + $matrix = []; + foreach ($array as $value) { + $matrix[] = [$value]; + } + + return new self($matrix); + } + + /** + * @return array + */ + public function toArray() + { + return $this->matrix; + } + + /** + * @return float + */ + public function toScalar() + { + return $this->matrix[0][0]; + } + + /** + * @return int + */ + public function getRows() + { + return $this->rows; + } + + /** + * @return int + */ + public function getColumns() + { + return $this->columns; + } + + /** + * @param $column + * + * @return array + * + * @throws MatrixException + */ + public function getColumnValues($column) + { + if ($column >= $this->columns) { + throw MatrixException::columnOutOfRange(); + } + + return array_column($this->matrix, $column); + } + + + /** + * @return float|int + * + * @throws MatrixException + */ + public function getDeterminant() + { + if ($this->determinant) { + return $this->determinant; + } + + if (!$this->isSquare()) { + throw MatrixException::notSquareMatrix(); + } + + return $this->determinant = $this->calculateDeterminant(); + } + + /** + * @return float|int + * + * @throws MatrixException + */ + private function calculateDeterminant() + { + $determinant = 0; + if ($this->rows == 1 && $this->columns == 1) { + $determinant = $this->matrix[0][0]; + } elseif ($this->rows == 2 && $this->columns == 2) { + $determinant = + $this->matrix[0][0] * $this->matrix[1][1] - + $this->matrix[0][1] * $this->matrix[1][0]; + } else { + for ($j = 0; $j < $this->columns; ++$j) { + $subMatrix = $this->crossOut(0, $j); + $minor = $this->matrix[0][$j] * $subMatrix->getDeterminant(); + $determinant += fmod((float) $j, 2.0) == 0 ? $minor : -$minor; + } + } + + return $determinant; + } + + /** + * @return bool + */ + public function isSquare() + { + return $this->columns === $this->rows; + } + + /** + * @return Matrix + */ + public function transpose() + { + if ($this->rows == 1) { + $matrix = array_map(function ($el) { + return [$el]; + }, $this->matrix[0]); + } else { + $matrix = array_map(null, ...$this->matrix); + } + + return new self($matrix, false); + } + + /** + * @param Matrix $matrix + * + * @return Matrix + * + * @throws InvalidArgumentException + */ + public function multiply(Matrix $matrix) + { + if ($this->columns != $matrix->getRows()) { + throw InvalidArgumentException::inconsistentMatrixSupplied(); + } + + $product = []; + $multiplier = $matrix->toArray(); + for ($i = 0; $i < $this->rows; ++$i) { + $columns = $matrix->getColumns(); + for ($j = 0; $j < $columns; ++$j) { + $product[$i][$j] = 0; + for ($k = 0; $k < $this->columns; ++$k) { + $product[$i][$j] += $this->matrix[$i][$k] * $multiplier[$k][$j]; + } + } + } + + return new self($product, false); + } + + /** + * @param $value + * + * @return Matrix + */ + public function divideByScalar($value) + { + $newMatrix = []; + for ($i = 0; $i < $this->rows; ++$i) { + for ($j = 0; $j < $this->columns; ++$j) { + $newMatrix[$i][$j] = $this->matrix[$i][$j] / $value; + } + } + + return new self($newMatrix, false); + } + + /** + * @param $value + * + * @return Matrix + */ + public function multiplyByScalar($value) + { + $newMatrix = []; + for ($i = 0; $i < $this->rows; ++$i) { + for ($j = 0; $j < $this->columns; ++$j) { + $newMatrix[$i][$j] = $this->matrix[$i][$j] * $value; + } + } + + return new self($newMatrix, false); + } + + /** + * Element-wise addition of the matrix with another one + * + * @param Matrix $other + */ + public function add(Matrix $other) + { + return $this->_add($other); + } + + /** + * Element-wise subtracting of another matrix from this one + * + * @param Matrix $other + */ + public function subtract(Matrix $other) + { + return $this->_add($other, -1); + } + + /** + * Element-wise addition or substraction depending on the given sign parameter + * + * @param Matrix $other + * @param type $sign + */ + protected function _add(Matrix $other, $sign = 1) + { + $a1 = $this->toArray(); + $a2 = $other->toArray(); + + $newMatrix = []; + for ($i=0; $i < $this->rows; $i++) { + for ($k=0; $k < $this->columns; $k++) { + $newMatrix[$i][$k] = $a1[$i][$k] + $sign * $a2[$i][$k]; + } + } + + return new Matrix($newMatrix, false); + } + + /** + * @return Matrix + * + * @throws MatrixException + */ + public function inverse() + { + if (!$this->isSquare()) { + throw MatrixException::notSquareMatrix(); + } + + if ($this->isSingular()) { + throw MatrixException::singularMatrix(); + } + + $newMatrix = []; + for ($i = 0; $i < $this->rows; ++$i) { + for ($j = 0; $j < $this->columns; ++$j) { + $minor = $this->crossOut($i, $j)->getDeterminant(); + $newMatrix[$i][$j] = fmod((float) ($i + $j), 2.0) == 0 ? $minor : -$minor; + } + } + + $cofactorMatrix = new self($newMatrix, false); + + return $cofactorMatrix->transpose()->divideByScalar($this->getDeterminant()); + } + + /** + * @param int $row + * @param int $column + * + * @return Matrix + */ + public function crossOut(int $row, int $column) + { + $newMatrix = []; + $r = 0; + for ($i = 0; $i < $this->rows; ++$i) { + $c = 0; + if ($row != $i) { + for ($j = 0; $j < $this->columns; ++$j) { + if ($column != $j) { + $newMatrix[$r][$c] = $this->matrix[$i][$j]; + ++$c; + } + } + ++$r; + } + } + + return new self($newMatrix, false); + } + + /** + * @return bool + */ + public function isSingular() : bool + { + return 0 == $this->getDeterminant(); + } + + /** + * Returns the transpose of given array + * + * @param array $array + * + * @return array + */ + public static function transposeArray(array $array) + { + return (new Matrix($array, false))->transpose()->toArray(); + } + + /** + * Returns the dot product of two arrays
+ * Matrix::dot(x, y) ==> x.y' + * + * @param array $array1 + * @param array $array2 + * + * @return array + */ + public static function dot(array $array1, array $array2) + { + $m1 = new Matrix($array1, false); + $m2 = new Matrix($array2, false); + + return $m1->multiply($m2->transpose())->toArray()[0]; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Product.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Product.php new file mode 100644 index 00000000000..35ef79cbfaa --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Product.php @@ -0,0 +1,26 @@ + $value) { + if (is_numeric($value) && is_numeric($b[$index])) { + $product += $value * $b[$index]; + } + } + + return $product; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Set.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Set.php new file mode 100644 index 00000000000..20fc78099e6 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Set.php @@ -0,0 +1,211 @@ +elements = self::sanitize($elements); + } + + /** + * Creates the union of A and B. + * + * @param Set $a + * @param Set $b + * + * @return Set + */ + public static function union(Set $a, Set $b) : Set + { + return new self(array_merge($a->toArray(), $b->toArray())); + } + + /** + * Creates the intersection of A and B. + * + * @param Set $a + * @param Set $b + * + * @return Set + */ + public static function intersection(Set $a, Set $b) : Set + { + return new self(array_intersect($a->toArray(), $b->toArray())); + } + + /** + * Creates the difference of A and B. + * + * @param Set $a + * @param Set $b + * + * @return Set + */ + public static function difference(Set $a, Set $b) : Set + { + return new self(array_diff($a->toArray(), $b->toArray())); + } + + /** + * Creates the Cartesian product of A and B. + * + * @param Set $a + * @param Set $b + * + * @return Set[] + */ + public static function cartesian(Set $a, Set $b) : array + { + $cartesian = []; + + foreach ($a as $multiplier) { + foreach ($b as $multiplicand) { + $cartesian[] = new self(array_merge([$multiplicand], [$multiplier])); + } + } + + return $cartesian; + } + + /** + * Creates the power set of A. + * + * @param Set $a + * + * @return Set[] + */ + public static function power(Set $a) : array + { + $power = [new self()]; + + foreach ($a as $multiplicand) { + foreach ($power as $multiplier) { + $power[] = new self(array_merge([$multiplicand], $multiplier->toArray())); + } + } + + return $power; + } + + /** + * Removes duplicates and rewrites index. + * + * @param string[]|int[]|float[] $elements + * + * @return string[]|int[]|float[] + */ + private static function sanitize(array $elements) : array + { + sort($elements, SORT_ASC); + + return array_values(array_unique($elements, SORT_ASC)); + } + + /** + * @param string|int|float $element + * + * @return Set + */ + public function add($element) : Set + { + return $this->addAll([$element]); + } + + /** + * @param string[]|int[]|float[] $elements + * + * @return Set + */ + public function addAll(array $elements) : Set + { + $this->elements = self::sanitize(array_merge($this->elements, $elements)); + + return $this; + } + + /** + * @param string|int|float $element + * + * @return Set + */ + public function remove($element) : Set + { + return $this->removeAll([$element]); + } + + /** + * @param string[]|int[]|float[] $elements + * + * @return Set + */ + public function removeAll(array $elements) : Set + { + $this->elements = self::sanitize(array_diff($this->elements, $elements)); + + return $this; + } + + /** + * @param string|int|float $element + * + * @return bool + */ + public function contains($element) : bool + { + return $this->containsAll([$element]); + } + + /** + * @param string[]|int[]|float[] $elements + * + * @return bool + */ + public function containsAll(array $elements) : bool + { + return !array_diff($elements, $this->elements); + } + + /** + * @return string[]|int[]|float[] + */ + public function toArray() : array + { + return $this->elements; + } + + /** + * @return \ArrayIterator + */ + public function getIterator() : \ArrayIterator + { + return new \ArrayIterator($this->elements); + } + + /** + * @return bool + */ + public function isEmpty() : bool + { + return $this->cardinality() == 0; + } + + /** + * @return int + */ + public function cardinality() : int + { + return count($this->elements); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Correlation.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Correlation.php new file mode 100644 index 00000000000..0f60223fc01 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Correlation.php @@ -0,0 +1,45 @@ + $xi) { + $yi = $y[$index]; + $sum += ($xi - $meanX) * ($yi - $meanY); + } + + if ($sample) { + --$n; + } + + return $sum / $n; + } + + /** + * Calculates covariance of two dimensions, i and k in the given data. + * + * @param array $data + * @param int $i + * @param int $k + * @param type $sample + * @param int $n + * @param float $meanX + * @param float $meanY + */ + public static function fromDataset(array $data, int $i, int $k, $sample = true, float $meanX = null, float $meanY = null) + { + if (empty($data)) { + throw InvalidArgumentException::arrayCantBeEmpty(); + } + + $n = count($data); + if ($sample && $n === 1) { + throw InvalidArgumentException::arraySizeToSmall(2); + } + + if ($i < 0 || $k < 0 || $i >= $n || $k >= $n) { + throw new \Exception("Given indices i and k do not match with the dimensionality of data"); + } + + if ($meanX === null || $meanY === null) { + $x = array_column($data, $i); + $y = array_column($data, $k); + + $meanX = Mean::arithmetic($x); + $meanY = Mean::arithmetic($y); + $sum = 0.0; + foreach ($x as $index => $xi) { + $yi = $y[$index]; + $sum += ($xi - $meanX) * ($yi - $meanY); + } + } else { + // In the case, whole dataset given along with dimension indices, i and k, + // we would like to avoid getting column data with array_column and operate + // over this extra copy of column data for memory efficiency purposes. + // + // Instead we traverse through the whole data and get what we actually need + // without copying the data. This way, memory use will be reduced + // with a slight cost of CPU utilization. + $sum = 0.0; + foreach ($data as $row) { + $val = []; + foreach ($row as $index => $col) { + if ($index == $i) { + $val[0] = $col - $meanX; + } + if ($index == $k) { + $val[1] = $col - $meanY; + } + } + $sum += $val[0] * $val[1]; + } + } + + if ($sample) { + --$n; + } + + return $sum / $n; + } + + /** + * Returns the covariance matrix of n-dimensional data + * + * @param array $data + * + * @return array + */ + public static function covarianceMatrix(array $data, array $means = null) + { + $n = count($data[0]); + + if ($means === null) { + $means = []; + for ($i=0; $i < $n; $i++) { + $means[] = Mean::arithmetic(array_column($data, $i)); + } + } + + $cov = []; + for ($i=0; $i < $n; $i++) { + for ($k=0; $k < $n; $k++) { + if ($i > $k) { + $cov[$i][$k] = $cov[$k][$i]; + } else { + $cov[$i][$k] = Covariance::fromDataset( + $data, $i, $k, true, $means[$i], $means[$k]); + } + } + } + + return $cov; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Gaussian.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Gaussian.php new file mode 100644 index 00000000000..df27f076dc6 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Gaussian.php @@ -0,0 +1,60 @@ +mean = $mean; + $this->std = $std; + } + + /** + * Returns probability density of the given $value + * + * @param float $value + * + * @return type + */ + public function pdf(float $value) + { + // Calculate the probability density by use of normal/Gaussian distribution + // Ref: https://en.wikipedia.org/wiki/Normal_distribution + $std2 = $this->std ** 2; + $mean = $this->mean; + return exp(- (($value - $mean) ** 2) / (2 * $std2)) / sqrt(2 * $std2 * pi()); + } + + /** + * Returns probability density value of the given $value based on + * given standard deviation and the mean + * + * @param float $mean + * @param float $std + * @param float $value + * + * @return float + */ + public static function distributionPdf(float $mean, float $std, float $value) + { + $normal = new self($mean, $std); + return $normal->pdf($value); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Mean.php b/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Mean.php new file mode 100644 index 00000000000..581a1225903 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Math/Statistic/Mean.php @@ -0,0 +1,75 @@ + $label) { + if ($label == $predictedLabels[$index]) { + ++$score; + } + } + + if ($normalize) { + $score /= count($actualLabels); + } + + return $score; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Metric/ClassificationReport.php b/lib/mlbackend/php/phpml/src/Phpml/Metric/ClassificationReport.php new file mode 100644 index 00000000000..c7cc147898f --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Metric/ClassificationReport.php @@ -0,0 +1,183 @@ +support = self::getLabelIndexedArray($actualLabels, $predictedLabels); + + foreach ($actualLabels as $index => $actual) { + $predicted = $predictedLabels[$index]; + ++$this->support[$actual]; + + if ($actual === $predicted) { + ++$truePositive[$actual]; + } else { + ++$falsePositive[$predicted]; + ++$falseNegative[$actual]; + } + } + + $this->computeMetrics($truePositive, $falsePositive, $falseNegative); + $this->computeAverage(); + } + + /** + * @return array + */ + public function getPrecision() + { + return $this->precision; + } + + /** + * @return array + */ + public function getRecall() + { + return $this->recall; + } + + /** + * @return array + */ + public function getF1score() + { + return $this->f1score; + } + + /** + * @return array + */ + public function getSupport() + { + return $this->support; + } + + /** + * @return array + */ + public function getAverage() + { + return $this->average; + } + + /** + * @param array $truePositive + * @param array $falsePositive + * @param array $falseNegative + */ + private function computeMetrics(array $truePositive, array $falsePositive, array $falseNegative) + { + foreach ($truePositive as $label => $tp) { + $this->precision[$label] = $this->computePrecision($tp, $falsePositive[$label]); + $this->recall[$label] = $this->computeRecall($tp, $falseNegative[$label]); + $this->f1score[$label] = $this->computeF1Score((float) $this->precision[$label], (float) $this->recall[$label]); + } + } + + private function computeAverage() + { + foreach (['precision', 'recall', 'f1score'] as $metric) { + $values = array_filter($this->$metric); + if (0 == count($values)) { + $this->average[$metric] = 0.0; + continue; + } + $this->average[$metric] = array_sum($values) / count($values); + } + } + + /** + * @param int $truePositive + * @param int $falsePositive + * + * @return float|string + */ + private function computePrecision(int $truePositive, int $falsePositive) + { + if (0 == ($divider = $truePositive + $falsePositive)) { + return 0.0; + } + + return $truePositive / $divider; + } + + /** + * @param int $truePositive + * @param int $falseNegative + * + * @return float|string + */ + private function computeRecall(int $truePositive, int $falseNegative) + { + if (0 == ($divider = $truePositive + $falseNegative)) { + return 0.0; + } + + return $truePositive / $divider; + } + + /** + * @param float $precision + * @param float $recall + * + * @return float + */ + private function computeF1Score(float $precision, float $recall): float + { + if (0 == ($divider = $precision + $recall)) { + return 0.0; + } + + return 2.0 * (($precision * $recall) / $divider); + } + + /** + * @param array $actualLabels + * @param array $predictedLabels + * + * @return array + */ + private static function getLabelIndexedArray(array $actualLabels, array $predictedLabels): array + { + $labels = array_values(array_unique(array_merge($actualLabels, $predictedLabels))); + sort($labels); + $labels = array_combine($labels, array_fill(0, count($labels), 0)); + + return $labels; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Metric/ConfusionMatrix.php b/lib/mlbackend/php/phpml/src/Phpml/Metric/ConfusionMatrix.php new file mode 100644 index 00000000000..664e355e945 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Metric/ConfusionMatrix.php @@ -0,0 +1,71 @@ + $actual) { + $predicted = $predictedLabels[$index]; + + if (!isset($labels[$actual]) || !isset($labels[$predicted])) { + continue; + } + + if ($predicted === $actual) { + $row = $column = $labels[$actual]; + } else { + $row = $labels[$actual]; + $column = $labels[$predicted]; + } + + $matrix[$row][$column] += 1; + } + + return $matrix; + } + + /** + * @param array $labels + * + * @return array + */ + private static function generateMatrixWithZeros(array $labels): array + { + $count = count($labels); + $matrix = []; + + for ($i = 0; $i < $count; ++$i) { + $matrix[$i] = array_fill(0, $count, 0); + } + + return $matrix; + } + + /** + * @param array $labels + * + * @return array + */ + private static function getUniqueLabels(array $labels): array + { + $labels = array_values(array_unique($labels)); + sort($labels); + $labels = array_flip($labels); + + return $labels; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/ModelManager.php b/lib/mlbackend/php/phpml/src/Phpml/ModelManager.php new file mode 100644 index 00000000000..c03d0ed25c7 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/ModelManager.php @@ -0,0 +1,54 @@ += 0 ? 1.0 : 0.0; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/ActivationFunction/Gaussian.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/ActivationFunction/Gaussian.php new file mode 100644 index 00000000000..0e3e848d0f1 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/ActivationFunction/Gaussian.php @@ -0,0 +1,20 @@ +beta = $beta; + } + + /** + * @param float|int $value + * + * @return float + */ + public function compute($value): float + { + return tanh($this->beta * $value); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/ActivationFunction/Sigmoid.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/ActivationFunction/Sigmoid.php new file mode 100644 index 00000000000..23ac7ce12f8 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/ActivationFunction/Sigmoid.php @@ -0,0 +1,33 @@ +beta = $beta; + } + + /** + * @param float|int $value + * + * @return float + */ + public function compute($value): float + { + return 1 / (1 + exp(-$this->beta * $value)); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Layer.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Layer.php new file mode 100644 index 00000000000..632bc009efd --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Layer.php @@ -0,0 +1,65 @@ +nodes[] = $this->createNode($nodeClass, $activationFunction); + } + } + + /** + * @param string $nodeClass + * @param ActivationFunction|null $activationFunction + * + * @return Neuron + */ + private function createNode(string $nodeClass, ActivationFunction $activationFunction = null) + { + if (Neuron::class == $nodeClass) { + return new Neuron($activationFunction); + } + + return new $nodeClass(); + } + + /** + * @param Node $node + */ + public function addNode(Node $node) + { + $this->nodes[] = $node; + } + + /** + * @return Node[] + */ + public function getNodes() + { + return $this->nodes; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Network.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Network.php new file mode 100644 index 00000000000..c6c25af2e99 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Network.php @@ -0,0 +1,30 @@ +layers[] = $layer; + } + + /** + * @return Layer[] + */ + public function getLayers(): array + { + return $this->layers; + } + + /** + * @return Layer + */ + public function getOutputLayer(): Layer + { + return $this->layers[count($this->layers) - 1]; + } + + /** + * @return array + */ + public function getOutput(): array + { + $result = []; + foreach ($this->getOutputLayer()->getNodes() as $neuron) { + $result[] = $neuron->getOutput(); + } + + return $result; + } + + /** + * @param mixed $input + * + * @return $this + */ + public function setInput($input) + { + $firstLayer = $this->layers[0]; + + foreach ($firstLayer->getNodes() as $key => $neuron) { + if ($neuron instanceof Input) { + $neuron->setInput($input[$key]); + } + } + + foreach ($this->getLayers() as $layer) { + foreach ($layer->getNodes() as $node) { + if ($node instanceof Neuron) { + $node->refresh(); + } + } + } + + return $this; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Network/MultilayerPerceptron.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Network/MultilayerPerceptron.php new file mode 100644 index 00000000000..04664f9c7db --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Network/MultilayerPerceptron.php @@ -0,0 +1,95 @@ +addInputLayer(array_shift($layers)); + $this->addNeuronLayers($layers, $activationFunction); + $this->addBiasNodes(); + $this->generateSynapses(); + } + + /** + * @param int $nodes + */ + private function addInputLayer(int $nodes) + { + $this->addLayer(new Layer($nodes, Input::class)); + } + + /** + * @param array $layers + * @param ActivationFunction|null $activationFunction + */ + private function addNeuronLayers(array $layers, ActivationFunction $activationFunction = null) + { + foreach ($layers as $neurons) { + $this->addLayer(new Layer($neurons, Neuron::class, $activationFunction)); + } + } + + private function generateSynapses() + { + $layersNumber = count($this->layers) - 1; + for ($i = 0; $i < $layersNumber; ++$i) { + $currentLayer = $this->layers[$i]; + $nextLayer = $this->layers[$i + 1]; + $this->generateLayerSynapses($nextLayer, $currentLayer); + } + } + + private function addBiasNodes() + { + $biasLayers = count($this->layers) - 1; + for ($i = 0; $i < $biasLayers; ++$i) { + $this->layers[$i]->addNode(new Bias()); + } + } + + /** + * @param Layer $nextLayer + * @param Layer $currentLayer + */ + private function generateLayerSynapses(Layer $nextLayer, Layer $currentLayer) + { + foreach ($nextLayer->getNodes() as $nextNeuron) { + if ($nextNeuron instanceof Neuron) { + $this->generateNeuronSynapses($currentLayer, $nextNeuron); + } + } + } + + /** + * @param Layer $currentLayer + * @param Neuron $nextNeuron + */ + private function generateNeuronSynapses(Layer $currentLayer, Neuron $nextNeuron) + { + foreach ($currentLayer->getNodes() as $currentNeuron) { + $nextNeuron->addSynapse(new Synapse($currentNeuron)); + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node.php new file mode 100644 index 00000000000..65d5cdcdcb4 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node.php @@ -0,0 +1,13 @@ +input = $input; + } + + /** + * @return float + */ + public function getOutput(): float + { + return $this->input; + } + + /** + * @param float $input + */ + public function setInput(float $input) + { + $this->input = $input; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node/Neuron.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node/Neuron.php new file mode 100644 index 00000000000..519443844a5 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node/Neuron.php @@ -0,0 +1,75 @@ +activationFunction = $activationFunction ?: new ActivationFunction\Sigmoid(); + $this->synapses = []; + $this->output = 0; + } + + /** + * @param Synapse $synapse + */ + public function addSynapse(Synapse $synapse) + { + $this->synapses[] = $synapse; + } + + /** + * @return Synapse[] + */ + public function getSynapses() + { + return $this->synapses; + } + + /** + * @return float + */ + public function getOutput(): float + { + if (0 === $this->output) { + $sum = 0; + foreach ($this->synapses as $synapse) { + $sum += $synapse->getOutput(); + } + + $this->output = $this->activationFunction->compute($sum); + } + + return $this->output; + } + + public function refresh() + { + $this->output = 0; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node/Neuron/Synapse.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node/Neuron/Synapse.php new file mode 100644 index 00000000000..b9c036fb211 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Node/Neuron/Synapse.php @@ -0,0 +1,70 @@ +node = $node; + $this->weight = $weight ?: $this->generateRandomWeight(); + } + + /** + * @return float + */ + protected function generateRandomWeight(): float + { + return 1 / random_int(5, 25) * (random_int(0, 1) ? -1 : 1); + } + + /** + * @return float + */ + public function getOutput(): float + { + return $this->weight * $this->node->getOutput(); + } + + /** + * @param float $delta + */ + public function changeWeight($delta) + { + $this->weight += $delta; + } + + /** + * @return float + */ + public function getWeight() + { + return $this->weight; + } + + /** + * @return Node + */ + public function getNode() + { + return $this->node; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Training.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Training.php new file mode 100644 index 00000000000..d876af2e50f --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Training.php @@ -0,0 +1,16 @@ +network = $network; + $this->theta = $theta; + } + + /** + * @param array $samples + * @param array $targets + * @param float $desiredError + * @param int $maxIterations + */ + public function train(array $samples, array $targets, float $desiredError = 0.001, int $maxIterations = 10000) + { + for ($i = 0; $i < $maxIterations; ++$i) { + $resultsWithinError = $this->trainSamples($samples, $targets, $desiredError); + + if ($resultsWithinError == count($samples)) { + break; + } + } + } + + /** + * @param array $samples + * @param array $targets + * @param float $desiredError + * + * @return int + */ + private function trainSamples(array $samples, array $targets, float $desiredError): int + { + $resultsWithinError = 0; + foreach ($targets as $key => $target) { + $result = $this->network->setInput($samples[$key])->getOutput(); + + if ($this->isResultWithinError($result, $target, $desiredError)) { + ++$resultsWithinError; + } else { + $this->trainSample($samples[$key], $target); + } + } + + return $resultsWithinError; + } + + /** + * @param array $sample + * @param array $target + */ + private function trainSample(array $sample, array $target) + { + $this->network->setInput($sample)->getOutput(); + $this->sigmas = []; + + $layers = $this->network->getLayers(); + $layersNumber = count($layers); + + for ($i = $layersNumber; $i > 1; --$i) { + foreach ($layers[$i - 1]->getNodes() as $key => $neuron) { + if ($neuron instanceof Neuron) { + $sigma = $this->getSigma($neuron, $target, $key, $i == $layersNumber); + foreach ($neuron->getSynapses() as $synapse) { + $synapse->changeWeight($this->theta * $sigma * $synapse->getNode()->getOutput()); + } + } + } + } + } + + /** + * @param Neuron $neuron + * @param array $target + * @param int $key + * @param bool $lastLayer + * + * @return float + */ + private function getSigma(Neuron $neuron, array $target, int $key, bool $lastLayer): float + { + $neuronOutput = $neuron->getOutput(); + $sigma = $neuronOutput * (1 - $neuronOutput); + + if ($lastLayer) { + $sigma *= ($target[$key] - $neuronOutput); + } else { + $sigma *= $this->getPrevSigma($neuron); + } + + $this->sigmas[] = new Sigma($neuron, $sigma); + + return $sigma; + } + + /** + * @param Neuron $neuron + * + * @return float + */ + private function getPrevSigma(Neuron $neuron): float + { + $sigma = 0.0; + + foreach ($this->sigmas as $neuronSigma) { + $sigma += $neuronSigma->getSigmaForNeuron($neuron); + } + + return $sigma; + } + + /** + * @param array $result + * @param array $target + * @param float $desiredError + * + * @return bool + */ + private function isResultWithinError(array $result, array $target, float $desiredError) + { + foreach ($target as $key => $value) { + if ($result[$key] > $value + $desiredError || $result[$key] < $value - $desiredError) { + return false; + } + } + + return true; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Training/Backpropagation/Sigma.php b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Training/Backpropagation/Sigma.php new file mode 100644 index 00000000000..62e2e7a5641 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/NeuralNetwork/Training/Backpropagation/Sigma.php @@ -0,0 +1,64 @@ +neuron = $neuron; + $this->sigma = $sigma; + } + + /** + * @return Neuron + */ + public function getNeuron() + { + return $this->neuron; + } + + /** + * @return float + */ + public function getSigma() + { + return $this->sigma; + } + + /** + * @param Neuron $neuron + * + * @return float + */ + public function getSigmaForNeuron(Neuron $neuron): float + { + $sigma = 0.0; + + foreach ($this->neuron->getSynapses() as $synapse) { + if ($synapse->getNode() == $neuron) { + $sigma += $synapse->getWeight() * $this->getSigma(); + } + } + + return $sigma; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Pipeline.php b/lib/mlbackend/php/phpml/src/Phpml/Pipeline.php new file mode 100644 index 00000000000..a6b3d562673 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Pipeline.php @@ -0,0 +1,106 @@ +addTransformer($transformer); + } + + $this->estimator = $estimator; + } + + /** + * @param Transformer $transformer + */ + public function addTransformer(Transformer $transformer) + { + $this->transformers[] = $transformer; + } + + /** + * @param Estimator $estimator + */ + public function setEstimator(Estimator $estimator) + { + $this->estimator = $estimator; + } + + /** + * @return array|Transformer[] + */ + public function getTransformers() + { + return $this->transformers; + } + + /** + * @return Estimator + */ + public function getEstimator() + { + return $this->estimator; + } + + /** + * @param array $samples + * @param array $targets + */ + public function train(array $samples, array $targets) + { + $this->fitTransformers($samples); + $this->transformSamples($samples); + $this->estimator->train($samples, $targets); + } + + /** + * @param array $samples + * + * @return mixed + */ + public function predict(array $samples) + { + $this->transformSamples($samples); + + return $this->estimator->predict($samples); + } + + /** + * @param array $samples + */ + private function fitTransformers(array &$samples) + { + foreach ($this->transformers as $transformer) { + $transformer->fit($samples); + } + } + + /** + * @param array $samples + */ + private function transformSamples(array &$samples) + { + foreach ($this->transformers as $transformer) { + $transformer->transform($samples); + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Imputer.php b/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Imputer.php new file mode 100644 index 00000000000..805d3f62096 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Imputer.php @@ -0,0 +1,99 @@ +missingValue = $missingValue; + $this->strategy = $strategy; + $this->axis = $axis; + $this->samples = $samples; + } + + /** + * @param array $samples + */ + public function fit(array $samples) + { + $this->samples = $samples; + } + + /** + * @param array $samples + */ + public function transform(array &$samples) + { + foreach ($samples as &$sample) { + $this->preprocessSample($sample); + } + } + + /** + * @param array $sample + */ + private function preprocessSample(array &$sample) + { + foreach ($sample as $column => &$value) { + if ($value === $this->missingValue) { + $value = $this->strategy->replaceValue($this->getAxis($column, $sample)); + } + } + } + + /** + * @param int $column + * @param array $currentSample + * + * @return array + */ + private function getAxis(int $column, array $currentSample): array + { + if (self::AXIS_ROW === $this->axis) { + return array_diff($currentSample, [$this->missingValue]); + } + + $axis = []; + foreach ($this->samples as $sample) { + if ($sample[$column] !== $this->missingValue) { + $axis[] = $sample[$column]; + } + } + + return $axis; + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Imputer/Strategy.php b/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Imputer/Strategy.php new file mode 100644 index 00000000000..9125e06fb5c --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Imputer/Strategy.php @@ -0,0 +1,15 @@ +norm = $norm; + } + + /** + * @param array $samples + */ + public function fit(array $samples) + { + if ($this->fitted) { + return; + } + + if ($this->norm == self::NORM_STD) { + $features = range(0, count($samples[0]) - 1); + foreach ($features as $i) { + $values = array_column($samples, $i); + $this->std[$i] = StandardDeviation::population($values); + $this->mean[$i] = Mean::arithmetic($values); + } + } + + $this->fitted = true; + } + + /** + * @param array $samples + */ + public function transform(array &$samples) + { + $methods = [ + self::NORM_L1 => 'normalizeL1', + self::NORM_L2 => 'normalizeL2', + self::NORM_STD=> 'normalizeSTD' + ]; + $method = $methods[$this->norm]; + + $this->fit($samples); + + foreach ($samples as &$sample) { + $this->$method($sample); + } + } + + /** + * @param array $sample + */ + private function normalizeL1(array &$sample) + { + $norm1 = 0; + foreach ($sample as $feature) { + $norm1 += abs($feature); + } + + if (0 == $norm1) { + $count = count($sample); + $sample = array_fill(0, $count, 1.0 / $count); + } else { + foreach ($sample as &$feature) { + $feature /= $norm1; + } + } + } + + /** + * @param array $sample + */ + private function normalizeL2(array &$sample) + { + $norm2 = 0; + foreach ($sample as $feature) { + $norm2 += $feature * $feature; + } + $norm2 = sqrt((float)$norm2); + + if (0 == $norm2) { + $sample = array_fill(0, count($sample), 1); + } else { + foreach ($sample as &$feature) { + $feature /= $norm2; + } + } + } + + /** + * @param array $sample + */ + private function normalizeSTD(array &$sample) + { + foreach ($sample as $i => $val) { + if ($this->std[$i] != 0) { + $sample[$i] = ($sample[$i] - $this->mean[$i]) / $this->std[$i]; + } else { + // Same value for all samples. + $sample[$i] = 0; + } + } + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Preprocessor.php b/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Preprocessor.php new file mode 100644 index 00000000000..3ec1566b88c --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Preprocessing/Preprocessor.php @@ -0,0 +1,11 @@ +samples = array_merge($this->samples, $samples); + $this->targets = array_merge($this->targets, $targets); + + $this->computeCoefficients(); + } + + /** + * @param array $sample + * + * @return mixed + */ + public function predictSample(array $sample) + { + $result = $this->intercept; + foreach ($this->coefficients as $index => $coefficient) { + $result += $coefficient * $sample[$index]; + } + + return $result; + } + + /** + * @return array + */ + public function getCoefficients() + { + return $this->coefficients; + } + + /** + * @return float + */ + public function getIntercept() + { + return $this->intercept; + } + + /** + * coefficient(b) = (X'X)-1X'Y. + */ + private function computeCoefficients() + { + $samplesMatrix = $this->getSamplesMatrix(); + $targetsMatrix = $this->getTargetsMatrix(); + + $ts = $samplesMatrix->transpose()->multiply($samplesMatrix)->inverse(); + $tf = $samplesMatrix->transpose()->multiply($targetsMatrix); + + $this->coefficients = $ts->multiply($tf)->getColumnValues(0); + $this->intercept = array_shift($this->coefficients); + } + + /** + * Add one dimension for intercept calculation. + * + * @return Matrix + */ + private function getSamplesMatrix() + { + $samples = []; + foreach ($this->samples as $sample) { + array_unshift($sample, 1); + $samples[] = $sample; + } + + return new Matrix($samples); + } + + /** + * @return Matrix + */ + private function getTargetsMatrix() + { + if (is_array($this->targets[0])) { + return new Matrix($this->targets); + } + + return Matrix::fromFlatArray($this->targets); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Regression/MLPRegressor.php b/lib/mlbackend/php/phpml/src/Phpml/Regression/MLPRegressor.php new file mode 100644 index 00000000000..72e6a81e52b --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Regression/MLPRegressor.php @@ -0,0 +1,80 @@ +hiddenLayers = $hiddenLayers; + $this->desiredError = $desiredError; + $this->maxIterations = $maxIterations; + $this->activationFunction = $activationFunction; + } + + /** + * @param array $samples + * @param array $targets + */ + public function train(array $samples, array $targets) + { + $layers = $this->hiddenLayers; + array_unshift($layers, count($samples[0])); + $layers[] = count($targets[0]); + + $this->perceptron = new MultilayerPerceptron($layers, $this->activationFunction); + + $trainer = new Backpropagation($this->perceptron); + $trainer->train($samples, $targets, $this->desiredError, $this->maxIterations); + } + + /** + * @param array $sample + * + * @return array + */ + protected function predictSample(array $sample) + { + return $this->perceptron->setInput($sample)->getOutput(); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/Regression/Regression.php b/lib/mlbackend/php/phpml/src/Phpml/Regression/Regression.php new file mode 100644 index 00000000000..542685ccb28 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/Regression/Regression.php @@ -0,0 +1,11 @@ + $label) { + $set .= sprintf('%s %s %s', ($targets ? $label : $numericLabels[$label]), self::sampleRow($samples[$index]), PHP_EOL); + } + + return $set; + } + + /** + * @param array $samples + * + * @return string + */ + public static function testSet(array $samples): string + { + if (!is_array($samples[0])) { + $samples = [$samples]; + } + + $set = ''; + foreach ($samples as $sample) { + $set .= sprintf('0 %s %s', self::sampleRow($sample), PHP_EOL); + } + + return $set; + } + + /** + * @param string $rawPredictions + * @param array $labels + * + * @return array + */ + public static function predictions(string $rawPredictions, array $labels): array + { + $numericLabels = self::numericLabels($labels); + $results = []; + foreach (explode(PHP_EOL, $rawPredictions) as $result) { + if (strlen($result) > 0) { + $results[] = array_search($result, $numericLabels); + } + } + + return $results; + } + + /** + * @param array $labels + * + * @return array + */ + public static function numericLabels(array $labels): array + { + $numericLabels = []; + foreach ($labels as $label) { + if (isset($numericLabels[$label])) { + continue; + } + + $numericLabels[$label] = count($numericLabels); + } + + return $numericLabels; + } + + /** + * @param array $sample + * + * @return string + */ + private static function sampleRow(array $sample): string + { + $row = []; + foreach ($sample as $index => $feature) { + $row[] = sprintf('%s:%s', $index + 1, $feature); + } + + return implode(' ', $row); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/SupportVectorMachine/Kernel.php b/lib/mlbackend/php/phpml/src/Phpml/SupportVectorMachine/Kernel.php new file mode 100644 index 00000000000..9918a3fc252 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/SupportVectorMachine/Kernel.php @@ -0,0 +1,28 @@ +type = $type; + $this->kernel = $kernel; + $this->cost = $cost; + $this->nu = $nu; + $this->degree = $degree; + $this->gamma = $gamma; + $this->coef0 = $coef0; + $this->epsilon = $epsilon; + $this->tolerance = $tolerance; + $this->cacheSize = $cacheSize; + $this->shrinking = $shrinking; + $this->probabilityEstimates = $probabilityEstimates; + + $rootPath = realpath(implode(DIRECTORY_SEPARATOR, [__DIR__, '..', '..', '..'])).DIRECTORY_SEPARATOR; + + $this->binPath = $rootPath.'bin'.DIRECTORY_SEPARATOR.'libsvm'.DIRECTORY_SEPARATOR; + $this->varPath = $rootPath.'var'.DIRECTORY_SEPARATOR; + } + + /** + * @param array $samples + * @param array $targets + */ + public function train(array $samples, array $targets) + { + $this->samples = array_merge($this->samples, $samples); + $this->targets = array_merge($this->targets, $targets); + + $trainingSet = DataTransformer::trainingSet($this->samples, $this->targets, in_array($this->type, [Type::EPSILON_SVR, Type::NU_SVR])); + file_put_contents($trainingSetFileName = $this->varPath.uniqid('phpml', true), $trainingSet); + $modelFileName = $trainingSetFileName.'-model'; + + $command = $this->buildTrainCommand($trainingSetFileName, $modelFileName); + $output = ''; + exec(escapeshellcmd($command), $output); + + $this->model = file_get_contents($modelFileName); + + unlink($trainingSetFileName); + unlink($modelFileName); + } + + /** + * @return string + */ + public function getModel() + { + return $this->model; + } + + /** + * @param array $samples + * + * @return array + */ + public function predict(array $samples) + { + $testSet = DataTransformer::testSet($samples); + file_put_contents($testSetFileName = $this->varPath.uniqid('phpml', true), $testSet); + file_put_contents($modelFileName = $testSetFileName.'-model', $this->model); + $outputFileName = $testSetFileName.'-output'; + + $command = sprintf('%ssvm-predict%s %s %s %s', $this->binPath, $this->getOSExtension(), $testSetFileName, $modelFileName, $outputFileName); + $output = ''; + exec(escapeshellcmd($command), $output); + + $predictions = file_get_contents($outputFileName); + + unlink($testSetFileName); + unlink($modelFileName); + unlink($outputFileName); + + if (in_array($this->type, [Type::C_SVC, Type::NU_SVC])) { + $predictions = DataTransformer::predictions($predictions, $this->targets); + } else { + $predictions = explode(PHP_EOL, trim($predictions)); + } + + if (!is_array($samples[0])) { + return $predictions[0]; + } + + return $predictions; + } + + /** + * @return string + */ + private function getOSExtension() + { + $os = strtoupper(substr(PHP_OS, 0, 3)); + if ($os === 'WIN') { + return '.exe'; + } elseif ($os === 'DAR') { + return '-osx'; + } + + return ''; + } + + /** + * @param $trainingSetFileName + * @param $modelFileName + * + * @return string + */ + private function buildTrainCommand(string $trainingSetFileName, string $modelFileName): string + { + return sprintf('%ssvm-train%s -s %s -t %s -c %s -n %s -d %s%s -r %s -p %s -m %s -e %s -h %d -b %d %s %s', + $this->binPath, + $this->getOSExtension(), + $this->type, + $this->kernel, + $this->cost, + $this->nu, + $this->degree, + $this->gamma !== null ? ' -g '.$this->gamma : '', + $this->coef0, + $this->epsilon, + $this->cacheSize, + $this->tolerance, + $this->shrinking, + $this->probabilityEstimates, + escapeshellarg($trainingSetFileName), + escapeshellarg($modelFileName) + ); + } +} diff --git a/lib/mlbackend/php/phpml/src/Phpml/SupportVectorMachine/Type.php b/lib/mlbackend/php/phpml/src/Phpml/SupportVectorMachine/Type.php new file mode 100644 index 00000000000..1b454a5d5f5 --- /dev/null +++ b/lib/mlbackend/php/phpml/src/Phpml/SupportVectorMachine/Type.php @@ -0,0 +1,33 @@ + + + + phpml + PHP-ML + MIT + 0.4.1+ + + + diff --git a/lib/mlbackend/php/version.php b/lib/mlbackend/php/version.php new file mode 100644 index 00000000000..592a6c76047 --- /dev/null +++ b/lib/mlbackend/php/version.php @@ -0,0 +1,29 @@ +. + +/** + * Version details. + * + * @package mlbackend_php + * @copyright 2017 David Monllao {@link http://www.davidmonllao.com/} + * @license http://www.gnu.org/copyleft/gpl.html GNU GPL v3 or later + */ + +defined('MOODLE_INTERNAL') || die(); + +$plugin->version = 2017051500; // The current plugin version (Date: YYYYMMDDXX). +$plugin->requires = 2017050500; // Requires this Moodle version. +$plugin->component = 'mlbackend_php'; // Full name of the plugin (used for diagnostics).