Несколько наборов коэффициентов для наименьших квадратов, подходящих в numpy / scipy - PullRequest
Новая
       78

Несколько наборов коэффициентов для наименьших квадратов, подходящих в numpy / scipy

0 голосов
/ 24 октября 2018

Есть ли способ выполнить несколько одновременных (но не связанных) подгонок наименьших квадратов с различными матрицами коэффициентов в numpy.linalg.lstsq или scipy.linalg.lstsq?Например, вот тривиальное линейное соответствие, которое я хотел бы иметь возможность делать с разными значениями x, но одинаковыми значениями y.В настоящее время я должен написать цикл: 

x = np.arange(12.0).reshape(4, 3)
y = np.arange(12.0, step=3.0)
m = np.stack((x, np.broadcast_to(1, x.shape)), axis=0)

fit = np.stack(tuple(np.linalg.lstsq(w, y, rcond=-1)[0] for w in m), axis=-1)

Это приводит к набору подгонок с одинаковым наклоном и разными перехватами, так что fit[n] соответствует коэффициентам m[n].

Линейный метод наименьших квадратов не является хорошим примером, поскольку он обратим, и у обеих функций есть опция для нескольких значений y.Тем не менее, это служит иллюстрацией моей точки зрения.

В идеале я хотел бы расширить это на любую «транслируемую» комбинацию a и b, где точно a.shape[-2] == b.shape[0], а последние измерения должнылибо совпадать, либо быть одним (или отсутствующим).Я не особо зациклен на том, какое измерение a является представлением различных матриц: просто было удобно сделать его первым, чтобы сократить цикл.

Есть ли встроенный метод в numpy?или scipy, чтобы избежать цикла Python?Я очень заинтересован в использовании lstsq вместо того, чтобы вручную транспонировать, умножать и инвертировать матрицы.

1 Ответ

0 голосов
/ 24 октября 2018

Вы можете использовать scipy.sparse.linalg.lsqr вместе с scipy.sparse.block_diag.Я просто не уверен, что это будет быстрее.

Пример: 

>>> import numpy as np
>>> from scipy.sparse import block_diag
>>> from scipy.sparse import linalg as sprsla
>>> 
>>> x = np.random.random((3,5,4))
>>> y = np.random.random((3,5))
>>> 
>>> for A, b in zip(x, y):
...     print(np.linalg.lstsq(A, b))
... 
(array([-0.11536962,  0.22575441,  0.03597646,  0.52014899]), array([0.22232195]), 4, array([2.27188101, 0.69355384, 0.63567141, 0.21700743]))
(array([-2.36307163,  2.27693405, -1.85653264,  3.63307554]), array([0.04810252]), 4, array([2.61853881, 0.74251282, 0.38701194, 0.06751288]))
(array([-0.6817038 , -0.02537582,  0.75882223,  0.03190649]), array([0.09892803]), 4, array([2.5094637 , 0.55673403, 0.39252624, 0.18598489]))
>>> 
>>> sprsla.lsqr(block_diag(x), y.ravel())
(array([-0.11536962,  0.22575441,  0.03597646,  0.52014899, -2.36307163,
        2.27693405, -1.85653264,  3.63307554, -0.6817038 , -0.02537582,
        0.75882223,  0.03190649]), 2, 15, 0.6077437777160813, 0.6077437777160813, 6.226368324510392, 106.63227777368986, 1.3277892240815807e-14, 5.36589277249043, array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]))
...