123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405 |
- <?php
- namespace MathPHP\Tests\LinearAlgebra\Decomposition;
- use MathPHP\Functions\Support;
- use MathPHP\LinearAlgebra\MatrixFactory;
- use MathPHP\Exception;
- use MathPHP\LinearAlgebra\NumericMatrix;
- use MathPHP\LinearAlgebra\Vector;
- use MathPHP\Tests\LinearAlgebra\Fixture\MatrixDataProvider;
- class SVDTest extends \PHPUnit\Framework\TestCase
- {
- use MatrixDataProvider;
- /**
- * @test SVD returns the expected array of U, S, and Vt factorized matrices
- * @dataProvider dataProviderForSVD
- * @dataProvider dataProviderForLesserRankSVD
- * @param array $A
- * @param array $expected
- * @throws \Exception
- */
- public function testSVD(array $A, array $expected)
- {
- // Given
- $A = MatrixFactory::createNumeric($A);
- $expected_S = MatrixFactory::createNumeric($expected['S']);
- // When
- $svd = $A->svd();
- // And
- $U = $svd->U;
- $S = $svd->S;
- $V = $svd->V;
- // Then A = USVᵀ
- $this->assertEqualsWithDelta($A->getMatrix(), $U->multiply($S)->multiply($V->transpose())->getMatrix(), 0.00001);
- // And S is expected solution to SVD
- $this->assertEqualsWithDelta($expected_S->getMatrix(), $S->getMatrix(), 0.00001);
- }
- /**
- * Test data created with:
- * R: svd(A)
- * Python: scipy.linalg.svd(A)
- * @return array
- */
- public function dataProviderForSVD(): array
- {
- return [
- [
- [
- [1, 0, 0, 0, 2],
- [0, 0, 3, 0, 0],
- [0, 0, 0, 0, 0],
- [0, 2, 0, 0, 0],
- ],
- [ // Technically, the order of the diagonal elements can be in any order
- 'S' => [
- [3, 0, 0, 0, 0],
- [0, sqrt(5), 0, 0, 0],
- [0, 0, 2, 0, 0],
- [0, 0, 0, 0, 0],
- ],
- ],
- ],
- [
- [
- [8, -6, 2],
- [-6, 7, -4],
- [2, -4, -3],
- ],
- [
- 'S' => [
- [14.528807, 0, 0],
- [0, 4.404176, 0],
- [0, 0, 1.875369],
- ],
- ],
- ],
- [
- [
- [1, 2],
- [3, 4],
- [5, 6],
- ],
- [
- 'S' => [
- [9.52551809, 0],
- [0, 0.51430058],
- [0, 0],
- ],
- ],
- ],
- [
- [[3]],
- [
- 'S' => [[3]],
- ],
- ],
- [
- [[0]],
- [
- 'S' => [[0]],
- ],
- ],
- [
- [[1]],
- [
- 'S' => [[1]],
- ],
- ],
- [
- [
- [1, 2, 3],
- [4, 5, 6],
- [7, 8, 9],
- ],
- [
- 'S' => [
- [1.684810e+01, 0, 0],
- [0, 1.068370e+00, 0],
- [0, 0, 4.418425e-16],
- ],
- ],
- ],
- [
- [
- [2, 2, 2],
- [2, 2, 2],
- [2, 2, 2],
- ],
- [
- 'S' => [
- [6, 0, 0],
- [0, 0, 0],
- [0, 0, 0],
- ],
- ],
- ],
- [
- [
- [-2, -2, -2],
- [-2, -2, -2],
- [-2, -2, -2],
- ],
- [
- 'S' => [
- [6, 0, 0],
- [0, 0, 0],
- [0, 0, 0],
- ],
- ],
- ],
- [
- [
- [1, 2, 3],
- [0, 4, 5],
- [0, 0, 6],
- ],
- [
- 'S' => [
- [9.0125424, 0, 0],
- [0, 2.9974695, 0],
- [0, 0, 0.8884012],
- ],
- ],
- ],
- [
- [
- [1, 0, 0],
- [2, 3, 0],
- [4, 5, 6],
- ],
- [
- 'S' => [
- [9.2000960, 0, 0],
- [0, 2.3843001, 0],
- [0, 0, 0.8205768],
- ],
- ],
- ],
- // Singular
- [
- [
- [1, 0],
- [0, 0],
- ],
- [
- 'S' => [
- [1, 0],
- [0, 0],
- ],
- ],
- ],
- // Singular
- [
- [
- [1, 0, 1],
- [0, 1, -1],
- [0, 0, 0],
- ],
- [
- 'S' => [
- [1.732051, 0, 0],
- [0, 1.000000, 0],
- [0, 0, 0.0],
- ],
- ],
- ],
- // Idempotent
- [
- [
- [1, 0, 0],
- [0, 0, 0],
- [0, 0, 1],
- ],
- [
- 'S' => [
- [1, 0, 0],
- [0, 1, 0],
- [0, 0, 0],
- ],
- ],
- ],
- // Idempotent
- [
- [
- [2, -2, -4],
- [-1, 3, 4],
- [1, -2, -3],
- ],
- [
- 'S' => [
- [7.937254, 0, 0],
- [0, 1, 0],
- [0, 0, 2.198569e-17],
- ],
- ],
- ],
- // Floating point
- [
- [
- [2.5, 6.3, 9.1],
- [-1.4, 3.0, 4.45],
- [1.01, 8.5, -3.334],
- ],
- [
- 'S' => [
- [12.786005, 0, 0],
- [0, 8.663327, 0],
- [0, 0, 2.315812],
- ],
- ],
- ],
- ];
- }
- /**
- * @return array
- */
- public function dataProviderForLesserRankSVD(): array
- {
- return [
- [
- [
- [1, 1, 1, 1, 1],
- [1, 1, 1, 1, 1],
- [1, 1, 1, 1, 1],
- ],
- [
- 'S' => [
- [3.872983, 0, 0, 0, 0],
- [0, 1.812987e-16, 0, 0, 0],
- [0, 0, 1.509615e-32, 0, 0],
- ],
- ],
- ]
- ];
- }
- /**
- * @test SVD properties
- * @dataProvider dataProviderForSVD
- * @param array $A
- * @throws \Exception
- */
- public function testSVDProperties(array $A)
- {
- // Given
- $A = MatrixFactory::createNumeric($A);
- // When
- $svd = $A->svd();
- // And
- $U = $svd->U;
- $S = $svd->S;
- $V = $svd->V;
- $D = $svd->D;
- // Then U and V are orthogonal
- $this->assertTrue($svd->getU()->isOrthogonal());
- $this->assertTrue($svd->getV()->isOrthogonal());
- // And S is rectangular diagonal with non-negative real numbers on the diagonal
- $this->assertTrue($S->isRectangularDiagonal());
- foreach ($S->getDiagonalElements() as $diagonalElement) {
- $this->assertTrue($diagonalElement >= 0);
- }
- // And D contains the diagonal elements of S
- $this->assertEqualsWithDelta($D->getVector(), $S->getDiagonalElements(), 0.00001, '');
- // And the number of non-zero singular values is equal to the rank of M
- $nonZeroSingularValues = array_filter(
- $D->getVector(),
- function ($singularValue) {
- return Support::isNotZero($singularValue);
- }
- );
- $this->assertEquals($A->rank(), count($nonZeroSingularValues));
- // And UUᵀ = I
- $this->assertEqualsWithDelta(MatrixFactory::identity($U->getM())->getMatrix(), $U->multiply($U->transpose())->getMatrix(), 0.00001);
- // And VVᵀ = I
- $this->assertEqualsWithDelta(MatrixFactory::identity($V->getM())->getMatrix(), $V->multiply($V->transpose())->getMatrix(), 0.00001);
- }
- /**
- * @test SVD properties of less than full rank matrices
- * @dataProvider dataProviderForLesserRankSVD
- * @param array $A
- * @throws \Exception
- */
- public function testLesserRankSVDProperties(array $A)
- {
- // Given
- $A = MatrixFactory::createNumeric($A);
- // When
- $svd = $A->svd();
- // Then
- $this->assertTrue($svd->getU()->isOrthogonal());
- $this->assertTrue($svd->getS()->isRectangularDiagonal());
- $this->assertEqualsWithDelta($svd->D->getVector(), $svd->getS()->getDiagonalElements(), 0.00001, '');
- }
- /**
- * @test SVD get properties
- */
- public function testSVDGetProperties()
- {
- // Given
- $A = MatrixFactory::createNumeric([
- [4, 1, -1],
- [1, 2, 1],
- [-1, 1, 2],
- ]);
- $svd = $A->svd();
- // When
- $S = $svd->S;
- $V = $svd->V;
- $D = $svd->D;
- $U = $svd->U;
- // Then
- $this->assertInstanceOf(NumericMatrix::class, $S);
- $this->assertInstanceOf(NumericMatrix::class, $V);
- $this->assertInstanceOf(NumericMatrix::class, $U);
- $this->assertInstanceOf(Vector::class, $D);
- // And
- $this->assertEquals($svd->getS(), $S);
- $this->assertEquals($svd->getV(), $V);
- $this->assertEquals($svd->getD(), $D);
- $this->assertEquals($svd->getU(), $U);
- }
- /**
- * @test SVD invalid property
- */
- public function testSVDInvalidProperty()
- {
- // Given
- $A = MatrixFactory::createNumeric([
- [4, 1, -1],
- [1, 2, 1],
- [-1, 1, 2],
- ]);
- $svd = $A->svd();
- // Then
- $this->expectException(Exception\MathException::class);
- // When
- $doesNotExist = $svd->doesNotExist;
- }
- }
|