problinsolve¶
- probnum.linalg.problinsolve(A, b, A0=None, Ainv0=None, x0=None, assume_A='sympos', maxiter=None, atol=1e-06, rtol=1e-06, callback=None, **kwargs)[source]¶
Infer a solution to the linear system \(A x = b\) in a Bayesian framework.
Probabilistic linear solvers infer solutions to problems of the form
\[Ax=b,\]where \(A \in \mathbb{R}^{n \times n}\) and \(b \in \mathbb{R}^{n}\). They return a probability measure which quantifies uncertainty in the output arising from finite computational resources. This solver can take prior information either on the linear operator \(A\) or its inverse \(H=A^{-1}\) in the form of a random variable
A0
orAinv0
and outputs a posterior belief over \(A\) or \(H\). This code implements the method described in Wenger et al. 1 based on the work in Hennig et al. 2.- Parameters
A (
Union
[ndarray
,spmatrix
,LinearOperator
,RandomVariable
]) – shape=(n, n) – A square linear operator (or matrix). Only matrix-vector products \(v \mapsto Av\) are used internally.b (
Union
[ndarray
,RandomVariable
]) – shape=(n, ) or (n, nrhs) – Right-hand side vector, matrix or random variable in \(A x = b\). For multiple right hand sides,nrhs
problems are solved sequentially with the posteriors over the matrices acting as priors for subsequent solves. If the right-hand-side is assumed to be noisy, every iteration of the solver samples a realization fromb
.A0 (
Union
[ndarray
,spmatrix
,LinearOperator
,RandomVariable
,None
]) – shape=(n, n) – A square matrix, linear operator or random variable representing the prior belief over the linear operator \(A\). If an array or linear operator is given, a prior distribution is chosen automatically.Ainv0 (
Union
[ndarray
,spmatrix
,LinearOperator
,RandomVariable
,None
]) – shape=(n, n) – A square matrix, linear operator or random variable representing the prior belief over the inverse \(H=A^{-1}\). This can be viewed as taking the form of a pre-conditioner. If an array or linear operator is given, a prior distribution is chosen automatically.x0 (
Union
[ndarray
,RandomVariable
,None
]) – shape=(n, ) or (n, nrhs) – Prior belief for the solution of the linear system. Will be ignored ifAinv0
is given.assume_A (
str
) –Assumptions on the linear operator which can influence solver choice and behavior. The available options are (combinations of)
generic matrix
gen
symmetric
sym
positive definite
pos
(additive) noise
noise
maxiter (
Optional
[int
]) – Maximum number of iterations. Defaults to \(10n\), where \(n\) is the dimension of \(A\).atol (
float
) – Absolute convergence tolerance.rtol (
float
) – Relative convergence tolerance.callback (
Optional
[Callable
]) – User-supplied function called after each iteration of the linear solver. It is called ascallback(xk, Ak, Ainvk, sk, yk, alphak, resid, **kwargs)
and can be used to return quantities from the iteration. Note that depending on the function supplied, this can slow down the solver considerably.kwargs (optional) – Optional keyword arguments passed onto the solver iteration.
- Return type
- Returns
x – Approximate solution \(x\) to the linear system. Shape of the return matches the shape of
b
.A – Posterior belief over the linear operator.
Ainv – Posterior belief over the linear operator inverse \(H=A^{-1}\).
info – Information on convergence of the solver.
- Raises
ValueError – If size mismatches detected or input matrices are not square.
LinAlgError – If the matrix
A
is singular.LinAlgWarning – If an ill-conditioned input
A
is detected.
Notes
For a specific class of priors the posterior mean of \(x_k=Hb\) coincides with the iterates of the conjugate gradient method. The matrix-based view taken here recovers the solution-based inference of
bayescg()
3.References
- 1
Wenger, J. and Hennig, P., Probabilistic Linear Solvers for Machine Learning, Advances in Neural Information Processing Systems (NeurIPS), 2020
- 2
Hennig, P., Probabilistic Interpretation of Linear Solvers, SIAM Journal on Optimization, 2015, 25, 234-260
- 3
Bartels, S. et al., Probabilistic Linear Solvers: A Unifying View, Statistics and Computing, 2019
See also
bayescg
Solve linear systems with prior information on the solution.
Examples
>>> import numpy as np >>> np.random.seed(1) >>> n = 20 >>> A = np.random.rand(n, n) >>> A = 0.5 * (A + A.T) + 5 * np.eye(n) >>> b = np.random.rand(n) >>> x, A, Ainv, info = problinsolve(A=A, b=b) >>> print(info["iter"]) 9