I have a system of equations that I am trying to solve. In matrix form, it's written as

$$x(I - S) = b.$$

I am solving for $x$, where $I$ is the identity matrix and $S$ is a matrix where each column sums up to one. $b$ is strictly positive.

Now, unfortunately, the determinant of $I - S$ is very small. Inverting it leads to a matrix that is huge everywhere, and even numpy.linalg.solve comes up with a solution for $x$ that is huge (every element around e+19). What is an alternative way of solving this? The matrix itself is not large, there are around 30-50 elements in $x$.

  • 1
    $\begingroup$ I don't know, but maybe some regularization method could be useful. $\endgroup$ – VoB Apr 7 '20 at 13:37
  • $\begingroup$ Usually SVD helps for such situations $\endgroup$ – Maxim Umansky Apr 8 '20 at 18:24
  • 1
    $\begingroup$ You've got a singular system of equations. There will either be no solutions or infinitely many solutions. What do you want to do if there are no solutions? If there are infinitely many solutions, which one do you want? $\endgroup$ – Brian Borchers Jan 3 at 21:40

You cannot solve what has no solutions

That matrix is singular, so the system has either zero or infinite solutions. In the case of your system, I think some Perron-Frobenius theory can be used to prove that there are zero: there is a non-negative vector such that $Sv=v$, hence $0=x(I-S)v=bv > 0$, contradiction. (Those $x$ and $b$ are row vectors, right?)

  • $\begingroup$ $x$ and $b$ are column vectors $\endgroup$ – FooBar Apr 7 '20 at 14:39
  • 3
    $\begingroup$ Note that you can't multiply a column vector $x$ by $(I-S)$. You could write this as $x^{T}(I-S)=b^{T}$ if $x$ and $b$ are to be column vectors. $\endgroup$ – Brian Borchers Apr 7 '20 at 15:38
  • $\begingroup$ Yeah sorry, my bad. I was tired when answering this. For the sake of this system they are row vectors. $\endgroup$ – FooBar Apr 7 '20 at 16:22
  • 2
    $\begingroup$ If you change the question I can edit the answer accordingly. Anyway, the algebra is slightly different, but all claims here still hold if $x,b$ are column vectors: your system has no solution. $\endgroup$ – Federico Poloni Apr 7 '20 at 17:15

Basic linear algebra states that $\det(I-S)$ must be non-zero so that a solution to your linear system exists. On the other hand, if your determinant is (numerically) zero, the basis vectors in your coefficient matrix $(I-S)$ do not span the whole vector space, so that you can't construct a general right-hand side $b$. Think of this picture:

enter image description here

In the example, the rhs $b$ is a 3D-vector with a $z$-component, but the vectors in the coefficient matrix only span the $x-y$-plane. So, the best you can do is to reproduce the vector in the $x-y$-plane, and forget about the $z$-direction. That is, you essentially solve for the projection $b^{||}$ of the vector $b$ onto the $x-y$-plane. By this you determine only two of the three parameters $x_1$ and $x_2$ in your solution vector. The third parameter $x_3$ can't be determined and is usually set to zero (such that $||x||$ is as small as possible).

To accomplish this, there are several regularization methods available. The two basic ones are the SVD and ridge regression:

  • The Singular Value Decomposisiton gives an easy method to form the pseudoinverse. It is used to decompose your matrix as $(I-S) = U D V^T$, where $D$ is the diagonal matrix of singular values. Now the pseudoinverse is formed by inverting all entries in $D$, and setting those $D_i$ larger than $1/\epsilon$ to zero (where $\epsilon$ is a small cutoff parameter).

  • In ridge regression, aka L2- or Tikhonov-regularization, you add a small term $\alpha$ to the diagonal of the moment matrix $(I-S)^T(I-S)$, which then becomes regular so that the linear system can be solved.

In the limit $\alpha,\epsilon \rightarrow 0$, both methods give the same result, namely the solution in the subspace that is spanned by the vectors in $(I-S)$.

By the way, to clarify: this is no solution, but an approximation to the solution. So the answer of @FredericiPoloni still holds. This here is just a workaround.

  • $\begingroup$ @FredericoPoloni: sure, thanks for pointing. $\endgroup$ – davidhigh Jan 3 at 19:56

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.