[Numpy-discussion] matrix default to column vector?
Charles R Harris
Sun Jun 14 12:02:53 CDT 2009
On Sun, Jun 14, 2009 at 10:20 AM, Tom K.<email@example.com> wrote:
> jseabold wrote:
>> On Mon, Jun 8, 2009 at 3:33 PM, Robert Kern<firstname.lastname@example.org> wrote:
>>> On Mon, Jun 8, 2009 at 14:10, Alan G Isaac<email@example.com> wrote:
>>>>>> Going back to Alan Isaac's example:
>>>>>> 1) beta = (X.T*X).I * X.T * Y
>>>>>> 2) beta = np.dot(np.dot(la.inv(np.dot(X.T,X)),X.T),Y)
>>>> Robert Kern wrote:
>>>>> 4) beta = la.lstsq(X, Y)
>>>>> I really hate that example.
> I propose the following alternative for discussion:
> U, s, Vh = np.linalg.svd(X, full_matrices=False)
> 1) beta = Vh.T*np.asmatrix(np.diag(1/s))*U.T*Y
> 2) beta = np.dot(Vh.T, np.dot(np.diag(1/s), np.dot(U.T, Y)))
> 1) is with X and Y starting out as matrices, 2) is with .
The problem is that I left the diagonal returned by svd as an array
rather than a matrix for backward compatibility. Diag returns the
diagonal when presented with a 2d matrix (array). I went back and
forth on that, but not doing so would have required everyone to change
their code to use diagflat instead of diag. I do note that diag
doesn't preserve matrices and that looks like a bug to me.
This is also an argument against over-overloaded functions such as
diag. Such functions are one of the blemishes of MatLab, IMHO. OTOH,
there is something of a shortage of short, meaningful names
More information about the Numpy-discussion