3

Following is a question that did not receive attention at math.SE at all. I am aware that is would be better suited there, but given its Economic background, perhaps it will get more attention here. If not, I'm happy to call this a failed attempt and delete it.

I need to solve $Ag=b$ for $b$, where $A$ is not symmetric.

$$A g = b \\ A = \left(\begin{matrix} A_1 \\ E_{IJ} \otimes \mathbb 1 \end{matrix}\right)$$

For some positive integers $I, J$: $A$ has size $(IJ\cdot J)\times IJ$, $g$ has size $IJ \times 1$, and $b$ has size $IJ \times 1$. $A$ is vertically stacked $A_1, A_2$, where $A_1$ is $IJ\times IJ$. $\mathbb 1$ is $J\times 1$ vector of ones, and $E$ is the identity matrix.

I need to solve this for $g$ on a computer for the case where $J$ small and $I$ large. I am fairly certain (through the way I get $A_1$ and $A_2$) that with $A_1$ being singular, my system of equations is exactly identifying $g$.

I'm used to the case where $A$ is symmetric and invertible. How do I proceed here?

FooBar
  • 10,712
  • 1
  • 29
  • 60
  • I don't understand. If $E_{IJ}$ is the (IJ x IJ) identity and $1$ is $1 X J$, then it seems that the kronecker product $E_{IJ} \otimes 1$ should be $IJ \times IJ^2$. If $A_1$ is $IJ \times IJ$, then there seems to be a conformity problem. – jmbejara Mar 06 '15 at 22:02
  • I don't think that's correct. You probably mean that $1 \otimes E_{IJ}$ should be $J* IJ \times IJ$. If A is an m × n matrix and B is a p × q matrix, then the Kronecker product $A \otimes B$ is the mp × nq block matrix. – jmbejara Mar 06 '15 at 22:35
  • @jmbejara yes, I agree - but the matrices are still conform. – FooBar Mar 06 '15 at 22:44
  • Yes, now that the notation is fixed. So, what exactly is the question? If the system is exactly identified and you just want to solve this on a computer, A\b in Matlab will solve this using an appropriate method. Do you want something more? – jmbejara Mar 06 '15 at 22:53
  • @jmbejara A is not symmetric and not invertible. I am not using Matlab. What would be an "appropriate method" here? Standard "linsolve" methods that I use in Python require invertibility. – FooBar Mar 06 '15 at 22:55
  • Note: I'm pretty sure if the system is singular there won't necessarily be a unique solution (Edit: It will have either infinite or possibly none -- see this question. Also recommended was least squares). That being said, if there is a method to solve it then it will be listed here – cc7768 Mar 06 '15 at 23:23
  • @cc7768 the system is not singular. $A_1$ is singular, $A$ is not. I can expand on my question on why I am pretty certain about that. I'll have a look at that link once I have sufficient time. – FooBar Mar 06 '15 at 23:27
  • I could be way off (I sometimes make bold dumb claims when it comes to linear algebra), but if $A$ isn't invertible then isn't it singular? – cc7768 Mar 06 '15 at 23:36
  • Notation is still unclear. As it is written now, $A_2$ will be $IJ^2 \times IJ$, and thus $A$ will be $IJ(J+1) \times IJ$. More importantly, $b$ must have as many elements as $A$ has rows. The system is overdetermined, so I'd expect typically no exact solution (unless $b$ has also some special structure?). – ivansml Mar 07 '15 at 02:50

1 Answers1

2

Maybe this helps. One method to do this is to use least squares (as was mentioned by @cc7768). Matlab does this by default when using the \ operator. This is also useful because it's sensible even if the system is overidentified as in the following example:

A = [1 22; 2 100; 7 6]
b = [1 2 5]'
A\b

produces

ans =

    0.7108
    0.0061

. In Python, this can be done with

import numpy
from numpy import linalg
A = numpy.matrix("1 22; 2 100; 7 6")
b = numpy.matrix("1 2 5").T
linalg.lstsq(A,b)[0]

which gives

matrix([[ 0.71084144],
        [ 0.00611577]])

.

jmbejara
  • 9,345
  • 5
  • 30
  • 80
  • Instead of Matlab you could stick with open source and use Julia ;) Python also has this functionality through scipy's solve method scipy.linalg.solve. – cc7768 Mar 06 '15 at 23:37