Overhaul SVDS input validation and test suite #20

mdhaber · 2021-06-27T01:40:54Z

Some work is not complete:

~~I haven't looked closely at the accuracy tests~~ Update: I confirmed that the basic test makes sense. I refactored it to cut lines, add cases, and strengthen the checks. I haven't looked carefully at test_svd_linop. I can't tell if it's an important addition or not, so rather than improving it or cutting it, I just left it alone. One thing that I think is missing is an accuracy test for a real-world, sparse matrix, but actually that is already in the svdp test suite, so I don't think I'll add it here.
~~I haven't looked closely at the edge cases~~ Update: I confirmed that the edge case tests make sense and refactored them to cut lines. Note that we're still not perfect here, especially with PROPACK.
~~There is no support for the PROPACK-specific options.~~ Update: save this for follow-up work.

There are a few kinks to work out:

The meaning of maxiter for PROPACK should be investigated more carefully. I think you had hooked it up to maxiter of svdp, but that is only used if irl_mode is True, so it seemed more appropriate to hook it up to kmax which is described as "Maximal number of iterations...". However, there is a requirement that kmax >= k that is more stringent than the usual requirements on maxiter. If we want to keep it this way, we should adjust the input validation test for solver='propack'. But I don't really know what's going on under the hood, so we should talk about this.
For arpack and lobpcg, v0 is supposed to be of length min(A.shape). This is not the case for PROPACK, and this is causing failure of one of the IV tests and also test_svd_linop. I think we can change the documentation and IV for PROPACK and this will fix the tests. Update: fixed.
The behavior of return_singular_value is really restrictive for PROPACK; it doesn't care about the shape of the matrix when you tell it to not calculate some of the eigenvectors. Currently I am enforcing the shape restrictions of ARPACK/LOBPCG, but that should be reconsidered.Update: PROPACK now respects the return_singular_value option regardless of matrix shape.
~~PROPACK is inaccurate when return_singular_vectors is FALSE in one test case. Weird.~~ Update: this was fixed by changing the random state.
PROPACK is failing one of the edge case tests. I think you noticed this before. I think I looked into it a bit and it might not be appropriate to raise an error here. If I remember correctly
- there is only one singular vector/value, and PROPACK finds it.
- When the remaining singular values are zero, it doesn't really attempt to find orthogonal singular vectors because they are not unique. I think the test checks for orthogonality of the singular vectors, so even when I commented out the error, PROPACK failed the test - but I think the test can be relaxed.
The documentation of return_singular_value did not reflect the behavior of ARPACK and LOBPCG after all. I think the inequalities on the shapes need to switch which is strict and which isn't to reflect the actual behavior.
~~AFAIK propack can handle k=min(A.shape); it does not need to be a strict inequality. This requirement of the IV should be relaxed.~~ Update: done.
~~I added random_state keyword support for PROPACK; I need to include it for the others.~~ Update: added for the others, but peculiarities of existing behavior needs to be documented (or we should drop the requirement that the outputs must be exactly the same as they were before.)
I've found a bunch of accuracy issues in LOBPCG in informal testing.
I should probably open an issue about LOBPCG's treatment of the v0 parameter.

I am continually tempted to abandon svds and just expose PROPACK in a function with the interface that PROPACK wants. But I figure I'll do the harder thing first - force PROPACK into SVDS - and then splitting it into its own function, if desired, would be very easy.

…xpansion

There was never a guarantee that solver=None was supported; solver='arpack' has always been the default value. No need to carry around this baggage.

…lar_vectors of svds

…r_vectors of svds

… svds

mdhaber · 2021-06-27T01:46:36Z

scipy/sparse/linalg/eigen/_svds_doc.py

-        If specified, must satistify ``k + 1 < ncv < N``; ``ncv > 2*k`` is
-        recommended.
+        If specified, must satistify ``k + 1 < ncv < min(M, N)``; ``ncv > 2*k``
+        is recommended.


This seemed to be wrong. I let's double check master before changing this.

mdhaber · 2021-06-27T01:48:07Z

scipy/sparse/linalg/eigen/_svds.py

-    if k <= 0 or k >= min(n, m):
-        raise ValueError("k must be between 1 and min(A.shape), k=%d" % k)
-
-    if isinstance(A, LinearOperator):


I just converted everything to a LinearOperator in the input validation so we can get rid of this branching.

mdhaber · 2021-06-27T01:49:37Z

scipy/sparse/linalg/eigen/tests/test_svds.py

@@ -1,8 +1,10 @@
+import re


I ended up basically rewriting all this. There were just a lot of holes, and it didn't make sense to leave the existing tests and add new ones that were more comprehensive.

mdhaber · 2021-06-27T01:51:12Z

scipy/sparse/linalg/interface.py

@@ -823,4 +823,5 @@ def aslinearoperator(A):
                                  rmatmat=rmatmat, dtype=dtype)

        else:
-            raise TypeError('type not understood')
+            A = np.atleast_2d(np.asarray(A))
+            return MatrixLinearOperator(A)


See scipygh-14299

scipy/sparse/linalg/eigen/tests/test_svds.py

…opack

…vers

mckib2 · 2021-06-27T20:45:20Z

scipy/sparse/linalg/eigen/_svds.py

+
+    # input validation/standardization for `maxiter`
+    if maxiter is not None and (int(maxiter) != maxiter or maxiter <= 0):
+        message = "`maxiter` must be a non-negative integer."


"must be a positive integer"

mckib2 · 2021-06-27T20:46:54Z

scipy/sparse/linalg/eigen/_svds.py


+    largest = (which == 'LM')
    n, m = A.shape


Is this where svds doesn't match docs? I would have expected m, n = A.shape

No, surprisingly. Despite flipping n, m between code and documentation, everything is right except for the inequality strictness being swapped in the documentation.

…tors parameter

…ular_vectors

…strength

…tive'

mdhaber · 2021-06-29T01:02:17Z

I think I'm done with this PR. With these changes and any other cleanups you want to do, I think we're about ready to open a PR on the main repo. When we do, there are a few things I'd like to bring to the attention of reviewers there:

To keep things simple, I decided not to expose any PROPACK-specific options yet. We proposed to add a wrapper for PROPACK; we've done that and much more. There are additional things that would be nice to have, but I don't think that should stop merging a PR that is already huge.
PROPACK doesn't return orthogonal singular vectors corresponding with zero singular values. This is basically the same as svds returning nans for zero input matrix scipy/scipy#3452 reported for ARPACK. I see now why you questioned how in svds you initially just returned the results of svdp without looking for small singular values. At the time, I didn't know PROPACK would have this problem, so it seemed reasonable to just return the results given by PROPACK. But yes, it would be nice to re-use some of that code that ensures that zeros out sufficiently small eigenvalues and makes sure that corresponding entries of u and vh are orthogonal. I don't think this is essential for the initial PR, but a reviewer might disagree.
Reviewers should think carefully about how we call svdp from svds - that is, how the keywords passed to svds map to PROPACK arguments. We've had to make some decisions about that and how to do the input validation. I think it's reasonable, but it's all stuff that could use another set of eyes.
Reviewers should note how the random_state parameter is implemented. Did I need to be so careful to avoid changing the default behavior, or should the default be like other functions with a random_state keyword and default to using the np.RandomState singleton? Also, it was intentional to use random_state, supporting the old np.RandomState, which we technically don't need to do anymore. I understand that use of the singleton is not considered good practice when you care about randomness, but I don't think we do care about randomness, and I think it is more convenient to set np.random.seed once to get repeatable results than it is to pass a random seed into every single function call.

mdhaber added 19 commits June 24, 2021 23:13

TST: sparse.linalg: strength svds k parameter IV

6f80f76

TST: sparse.linalg: strengthen svds which parameter IV

9e46ae9

MAINT: sparse.linalg: refactor svds input validation to prepare for e…

b1da2ef

…xpansion

MAINT: sparse.linalg: add input validation for parameter A of svds

4e6d28a

MAINT: sparse.linalg: add input validation for parameter solver of svds

6278f87

MAINT: sparse.linalg: add input validation for parameter ncv of svds

736702b

MAINT: sparse.linalg: remove support for svds solver=None

18d7312

There was never a guarantee that solver=None was supported; solver='arpack' has always been the default value. No need to carry around this baggage.

MAINT: sparse.linalg: add input validation for parameter v0 of svds

6dd1f08

MAINT: sparse.linalg: add input validation for parameter tol of svds

94a2ace

MAINT: sparse.linalg: add input validation for parameter maxiter of svds

48dd1c9

MAINT: sparse.linalg: add input validation for parameter return_singu…

5df6de3

…lar_vectors of svds

MAINT: sparse.linalg: check functionality of parameter k of svds

3eae3f4

MAINT: sparse.linalg: check functionality of parameter tol of svds

2b2b467

MAINT: sparse.linalg: can't check functionality of parameter ncv of svds

21126ad

MAINT: sparse.linalg: check functionality of parameter which of svds

ed8bdf7

MAINT: sparse.linalg: check functionality of parameter v0 of svds

57e10ca

MAINT: sparse.linalg: check functionality of parameter return_singula…

0e4a045

…r_vectors of svds

MAINT: sparse.linalg: check functionality of parameter maxiter of svds

fcd5129

TST: sparse.linalg: refactor test of functionality of parameter v0 of…

2e38f04

… svds

mdhaber commented Jun 27, 2021

View reviewed changes

scipy/sparse/linalg/eigen/tests/test_svds.py Show resolved Hide resolved

mdhaber commented Jun 27, 2021

View reviewed changes

scipy/sparse/linalg/eigen/tests/test_svds.py Outdated Show resolved Hide resolved

mdhaber changed the title ~~Overhaul PROPACK input validation and test suite~~ Overhaul SVDS input validation and test suite Jun 27, 2021

mdhaber added 3 commits June 26, 2021 22:24

MAINT: sparse.linalg: fix v0 length requirement for PROPACK

9d2864d

MAINT: sparse.linalg: add support for random_state keyword, solver pr…

762d423

…opack

MAINT: sparse.linalg: add support for random_state keyword, other sol…

2c3d61b

…vers

mckib2 reviewed Jun 27, 2021

View reviewed changes

mckib2 mentioned this pull request Jun 28, 2021

TST: complex example checks to PROPACK tests #21

Merged

1 task

mdhaber added 13 commits June 28, 2021 09:42

TST: sparse.linalg: adjust tests due to peculiarities of solvers

fb01269

TST: sparse.linalg: add svds random_state tests

2542ae7

TST: sparse.linalg: refactor svds tests for concision

da54a62

ENH: sparse.linalg: propack can do complete SVD

1ab3156

DOC: sparse.linalg: document that PROPACK can do complete SVD

7fe879e

DOC: sparse.linalg: improve documentation of svds return_singular_vec…

45e3856

…tors parameter

DOC: sparse.linalg: improve svds example

2b0c051

ENH: sparse.linalg: svds solver='propack' always respects return_sing…

72f5a6a

…ular_vectors

TST: sparse.linalg: xfail edge case test for svds with solver='propack'

f3781b1

TST: sparse.linalg: svds: refactor (@parametrize) edge case tests

148d93c

TST: sparse.linalg: svds: refactor test_svd_simple for concision and …

0a033a6

…strength

MAINT: sparse.linalg: svds: fix error message 'non-negative' -> 'posi…

263791c

…tive'

TST: sparse.linalg: svds: combine similar tests

3c23240

mckib2 merged commit 33d8c1b into mckib2:propack Jun 29, 2021

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Overhaul SVDS input validation and test suite #20

Overhaul SVDS input validation and test suite #20

mdhaber commented Jun 27, 2021 •

edited

Loading

mdhaber Jun 27, 2021

mdhaber Jun 27, 2021

mdhaber Jun 27, 2021

mdhaber Jun 27, 2021

mckib2 Jun 27, 2021

mckib2 Jun 27, 2021

mdhaber Jun 28, 2021

mdhaber commented Jun 29, 2021

Overhaul SVDS input validation and test suite #20

Overhaul SVDS input validation and test suite #20

Conversation

mdhaber commented Jun 27, 2021 • edited Loading

mdhaber Jun 27, 2021

Choose a reason for hiding this comment

mdhaber Jun 27, 2021

Choose a reason for hiding this comment

mdhaber Jun 27, 2021

Choose a reason for hiding this comment

mdhaber Jun 27, 2021

Choose a reason for hiding this comment

mckib2 Jun 27, 2021

Choose a reason for hiding this comment

mckib2 Jun 27, 2021

Choose a reason for hiding this comment

mdhaber Jun 28, 2021

Choose a reason for hiding this comment

mdhaber commented Jun 29, 2021

mdhaber commented Jun 27, 2021 •

edited

Loading