Многовариантная множественная линейная регрессия с использованием только numpy

Question

Я хочу построить многомерную модель множественной линейной регрессии с N зависимыми переменными и M независимыми переменными. Я смотрел вокруг и, похоже, не нашел реализации. Я провел небольшое исследование и нашел здесь несколько примечаний: http://users.stat.umn.edu/~helwig/notes/mvlr-Notes.pdf на слайде 51. Это кажется очень простым в реализации:

import numpy as np


M = 10
N = 3
p = 15

Y = np.random.rand(p,N)

X = np.random.rand(p,M)

A = np.dot(np.transpose(X),X)
B = np.dot(np.transpose(X),Y)

sol = np.linalg.solve(A,B)

где sol выводит матрицу коэффициентов. В конечном итоге я буду масштабировать это до очень больших наборов данных. Меня больше всего беспокоит точность этого метода. Если честно, все кажется слишком простым. Может ли кто-нибудь взвесить, достаточно ли этого в многомерной, множественной регрессии или есть какой-то пакет или что-то еще, что я могу использовать, что лучше?

Спасибо

Answer 1

За кулисами происходит гораздо больше, чем вы думаете. То, что вы видите, - это решение с наименьшими квадратами системы чрезмерно ограниченной линейной системы уравнений. Позвольте мне объяснить это так. У вас есть p уравнения и q неизвестно.

Случай 1: p < q Бесконечно много решений. Матрица A размера p x q сингулярна, и поэтому существует бесконечно много решений. Их можно найти, найдя конкретное решение и основу для нулевого пространства. np.linalg.solve не может использоваться для решения такой системы, поскольку он ожидает только квадратные матрицы полного ранга. Вместо этого вы можете использовать np.linalg.lstsq.

Случай 2: p = q Уникальное решение. Это означает, что матрица p x q обратима, и вы можете решить систему Ax = b, используя x = A^(-1) b. Фактически, когда вы звоните np.linalg.solve, именно это и происходит.

Случай 3: p > q Решения не существует. Но мы можем сделать аппроксимацию, проецируя вектор b ортогонально на пространство столбцов матрицы A. Это означает, что мы хотим найти вектор b_hat такой, что b_hat + w = b, где w - некоторый произвольный вектор, перпендикулярный b_hat, а b_hat лежит в пространстве столбцов матрицы A. Следовательно, есть некоторые x_hat так, что A * x_hat = b и A^(T) w = 0 (потому что w лежит в левом пустом пространстве). Как w = b - b_hat, мы имеем A^(T) w = A^(T) * b - A^(T) * ( A * x_hat ). Следовательно, мы приходим к решению, x_hat = (A^(T) * A)^(-1) * A^(T) * b. Это x_hat является наилучшим приближенным решением системы линейных уравнений с p > q. Это можно доказать математически, но я думаю, что это правдоподобное объяснение. Опять же, np.linalg.solve нельзя использовать для решения этого типа линейной системы, но есть еще одна процедура np.linalg.lstsq, которую можно использовать.