2012-11-06 18 views
13

Praca nad algebrą macierzową tutaj. Czasami potrzebuję odwrócić matrycę, która może być pojedyncza lub źle uwarunkowana. Rozumiem, że jest to pythonic, aby po prostu to zrobić:Efektywne i pythonicowe sprawdzanie pojedynczej matrycy

try: 
    i = linalg.inv(x) 
except LinAlgErr as err: 
    #handle it 

, ale nie jestem pewien, jak skuteczne jest to. Czy nie byłoby lepiej?

if linalg.cond(x) < 1/sys.float_info.epsilon: 
    i = linalg.inv(x) 
else: 
    #handle it 

Czy numpy.linalg po prostu wykonuje z góry test, który zakazałem?

+0

Próbowanie czegoś i reagowanie na warunki błędu po ich wystąpieniu jest zwykle Pythonicznym sposobem robienia rzeczy, ale sposób Pythonica nie zawsze jest najbardziej efektywnym sposobem. Zakłada się, że jeśli zależy ci na robieniu rzeczy "w sensie Pythonicznym", efektywność jest priorytetem drugorzędnym. –

+0

Nie mów, że efektywność jest priorytetem dodatkowym; zamiast tego powiem, że chcę przeprowadzić biwariacyjną optymalizację: pythonicity + efficiency (stąd tytuł postu). –

+2

W obliczeniach numerycznych zwykle uważane jest za niewłaściwe ćwiczenie, aby wyraźnie obliczyć odwrotność. W większości przypadków znacznie lepiej jest obliczyć dekompozycję LU za pomocą scipy.linalg.lu_factor, a później można go szybko rozwiązać dla wielu wektorów za pomocą scipy.linalg.lu_solve. – DaveP

Odpowiedz

7

Więc na podstawie tutaj wejść, jestem znakowania mój oryginalny blok kodu z wyraźnym testu jako rozwiązanie:

if linalg.cond(x) < 1/sys.float_info.epsilon: 
    i = linalg.inv(x) 
else: 
    #handle it 

Niespodziewanie, funkcja numpy.linalg.inv nie wykonuje ten test. Sprawdziłem kod i stwierdziłem, że przechodzi on przez wszystkie jego machinacje, a następnie po prostu wywołuje rutynową procedurę - wydaje się dość nieefektywny. Poza tym dodałbym drugi punkt DaveP: że odwrotność macierzy nie powinna być obliczana, chyba że jest to wyraźnie potrzebne.

+0

Zamiast 'sys.float_info.epsilon' wolę' numpy.spacing', ponieważ może on mieć wartość 1.0 jako float16, float32, float64, ... niezależnie od typu 'x'. –

+0

Powinieneś użyć epsilon dla dtype samej tablicy, np. 'np.finfo (x.dtype) .eps', które może różnić się od' sys.float_info.epsilon' –

3

Należy obliczyć wartość condition number macierzy, aby sprawdzić, czy jest odwracalna.

import numpy.linalg 

if numpy.isfinite(numpy.linalg.cond(A)): 
    B = numpy.linalg.inv(A) 
else: 
    # handle it 
+2

Matryca nie będzie kondycjonowana, gdy liczba warunków zbliży się do 1/epsilon. Tak więc drugie rozwiązanie, które zostało zaproponowane w pytaniu, jest lepsze niż rozwiązanie, które obejmuje tylko wyjątkowo dużą liczbę warunków. – DaveP

11

Twoje pierwsze rozwiązanie łapie przypadek gdzie matryca jest tak osobliwy, że numpy nie radzą sobie w ogóle - potencjalnie dość skrajny przypadek. Twoje drugie rozwiązanie jest lepsze, ponieważ łapie przypadek, w którym numpy daje odpowiedź, ale ta odpowiedź jest potencjalnie zepsuta błędem zaokrąglania - wydaje się to o wiele bardziej sensowne.

Jeśli próbujesz odwrócić źle przygotowane macierze, powinieneś rozważyć użycie singular value decomposition. Jeśli użyje się go ostrożnie, może dać rozsądną odpowiedź w przypadku niepowodzenia innych procedur.

Jeśli nie chcesz SVD, zobacz także mój komentarz na temat używania lu_factor zamiast inv.

+0

Pierwszy to, jak sądzę (i wydaje się potwierdzony), najbardziej pythonic, ale nie nazwałbym go "moim rozwiązaniem" :-). Po złapaniu złej kondycji, znam już wiele metod regulacji, których mogę użyć, aby to naprawić. –

+0

Domyślam się, że moim głównym przesłaniem jest to, że zanim się martwisz o bycie "najbardziej pytonistą", najpierw upewnij się, że jesteś najbardziej poprawny pod względem numerycznym. Drugie proponowane przez ciebie rozwiązanie jest lepsze z tego punktu widzenia. – DaveP

+0

Całkowicie się z Tobą zgadzam DaveP. –

0

Dlaczego nie sprawdzić, czy wyznacznik jest różny od zera?

det = numpy.linalg.det(A) 
if det != 0: 
    #proceed 
+0

To działałoby, gdyby komputery mogły reprezentować floats * dokładnie *. Jednak nie mogą tego zrobić, więc przez większość czasu otrzymasz odpowiedź typu 'det = 3.123e-14' i nigdy nie będziesz mieć pewności, czy jest to 0, czy tylko mała liczba. Ale oczywiście, gdybyś mógł obliczyć dokładne liczby, sprawdzanie niezerowej determinacji jest drogą do zrobienia. – tamasgal