редкая регрессия наименьших квадратов

Question

Я пытаюсь установить линейную регрессию Ax = b где A является разреженной матрицей и b редким вектором. Я пробовал scipy.sparse.linalg.lsqr, но, видимо, b должен быть массивным (плотным) массивом. Действительно, если я бегу

A = [list(range(0,10)) for i in range(0,15)] 
A = scipy.sparse.coo_matrix(A) 
b = list(range(0,15)) 
b = scipy.sparse.coo_matrix(b) 
scipy.sparse.linalg.lsqr(A,b) 

Я в конечном итоге с:

AttributeError: squeeze not found

Хотя

scipy.sparse.linalg.lsqr(A,b.toarray()) 

, кажется, работает.

К сожалению, в моем случае b представляет собой вектор 1,5 миллиарда x 1, и я просто не могу использовать плотный массив. Кто-нибудь знает обходной путь или другие библиотеки для выполнения линейной регрессии с разреженной матрицей и вектором?

Answer 1

Кажется, что документация специально запрашивает массив numpy. Однако, учитывая масштаб вашей проблемы, возможно, ее проще использовать замкнутое решение линейных наименьших квадратов?

Учитывая, что вы хотите решить Ax = b, вы можете отличить нормальные уравнения и решить их. Другими словами, вы бы решили min ||Ax-b||.

Закрытое решение формы будет x = (A.T*A)^{-1} * A.T *b. Конечно, это закрытое решение формы имеет свои собственные требования (в частности, в ранге матрицы A).

Вы можете решить для x, используя spsolve, или если это слишком дорого, а затем используйте итеративный решатель (например, Conjugate Gradients), чтобы получить неточное решение.

Код будет:

A = scipy.sparse.rand(1500,1000,0.5) #Create a random instance 
b = scipy.sparse.rand(1500,1,0.5) 
x = scipy.sparse.linalg.spsolve(A.T*A,A.T*b) 
x_lsqr = scipy.sparse.linalg.lsqr(A,b.toarray()) #Just for comparison 
print scipy.linalg.norm(x_lsqr[0]-x) 

, который на несколько случайных случаях, последовательно дал мне значения меньше, чем 1E-7.

Answer 2

По-видимому, миллиарды наблюдений слишком много для моей машины. Я в конечном итоге:

1. Изменение алгоритма стохастического градиентного спуска (SGD): быстрее многих набл
2. Удаление полностью редкие примеры (т.е. функции и метки, равные нулю)

Действительно, обновление правило SGD с функцией наименьших квадратов потерь всегда равно нулю для obs в 2. Это уменьшенные наблюдения от миллиардов до миллионов, что оказалось возможным при SGD на моей машине.