Решить линейную систему, которая имеет бесконечно много решений с обратной косой чертой (\) и lsqminnorm. Сравните результаты, используя 2-нормы решений.

Когда существуют бесконечные решения для $\mathrm{Ax}\mathit{=}$b, каждое из них $минимизирует\mathrm{\Vert }\mathit{}$Ax-b ‖. Команда обратной косой черты (\) вычисляет одно такое решение, но это решение обычно не минимизирует $\Vert \mathit{}x$‖. Решение, вычисленное поlsqminnorm минимизирует не только norm(A*x-b), но также norm(x).

Рассмотрим простую линейную систему с одним уравнением и двумя неизвестными, $\mathrm{}{\mathit{}}_{\mathrm{2x1}}\mathrm{+}{\mathit{}}_{\mathrm{}}\mathrm{3x2}$= 8. Эта система недооценивается, поскольку уравнений меньше, чем неизвестных. Решите уравнение, используя обратную косую черту иlsqminnorm.

A = [2 3];
b = 8;
x_a = A\b

x_a = 2×1

         0
    2.6667

x_b = lsqminnorm(A,b)

x_b = 2×1

    1.2308
    1.8462

Два метода получают разные решения, потому что обратная косая черта направлена только на минимизацию norm(A*x-b), тогда как lsqminnorm также стремится минимизировать norm(x). Рассчитайте эти нормы и поместите результаты в таблицу для простого сравнения.

s1 = {'Backslash'; 'lsqminnorm'};
s2 = {'norm_Ax_minus_b','norm_x'};
T = table([norm(A*x_a-b); norm(A*x_b-b)],[norm(x_a); norm(x_b)],'RowNames',s1,'VariableNames',s2)

T=2×2 table
                  norm_Ax_minus_b    norm_x
                  _______________    ______

    Backslash                0       2.6667
    lsqminnorm      8.8818e-16       2.2188

Этот рисунок иллюстрирует ситуацию и показывает, какие решения возвращаются к каждому из методов. Синяя линия представляет бесконечное число решений уравнения ${\mathit{}}_{\mathrm{x2}}=\frac{\mathrm{}}{\mathrm{-}}{\mathit{}}_{\mathrm{}}\frac{\mathrm{23x1}}{\mathrm{}}$+ 83. Оранжевая окружность представляет минимальное расстояние от начала координат до линии решений, и решение возвращаетсяlsqminnorm лежит точно в касательной точке между линией и окружностью, указывая, что именно решение наиболее близко к началу координат.

Показать, как задать допуск для вычисления ранга в lsqminnorm может помочь определить масштаб проблемы, чтобы случайный шум не повредил решение.

Создайте матрицу низкого ранга ранга 5 и вектор правой стороны b.

rng default % for reproducibility
U = randn(200,5);
V = randn(100,5);
A = U*V';
b = U*randn(5,1) + 1e-4*randn(200,1);

Решить линейную систему $\mathrm{Ax}\mathit{=}$ b с помощьюlsqminnorm. Вычислить нормы A*x-b и x для проверки качества решения.

x = lsqminnorm(A,b);
norm(A*x-b)

ans = 0.0014

norm(x)

ans = 0.1741

Теперь добавьте небольшое количество шума в матрицу A и снова решить линейную систему. Шум влияет на вектор решения x линейной системы непропорционально.

Anoise = A + 1e-12*randn(200,100);
xnoise = lsqminnorm(Anoise,b);
norm(Anoise*xnoise - b)

ans = 0.0010

norm(xnoise)

ans = 1.1216e+08

Причиной большой разницы в решениях является то, что шум влияет на приближение низкого ранга A. Другими словами, lsqminnorm обрабатывает небольшие значения по диагонали R матрица в QR разложении A как более важные, чем они есть. В идеале эти небольшие значения на диагонали R следует рассматривать как нули.

Постройте график диагональных элементов R матрица в QR разложении Anoise. Большое количество диагональных элементов находится порядка 1e-10.

[Q,R,p] = qr(Anoise,0);
semilogy(abs(diag(R)),'o')

Решение этой проблемы заключается в увеличении допуска, используемого lsqminnorm чтобы низкоранговое приближение Anoise при вычислении используется ошибка менее 1e-8. Это делает результат гораздо менее восприимчивым к шуму. Решение, использующее допуск, очень близко к исходному решению x.

xnoise = lsqminnorm(Anoise, b, 1e-8);
norm(Anoise*xnoise - b)

ans = 0.0014

norm(xnoise)

ans = 0.1741

norm(x - xnoise)

ans = 1.0804e-14

Решите линейную систему, включающую матрицу коэффициентов низкого ранга с включенными предупреждениями.

Создайте матрицу 3 на 3 ранга 2. В этой матрице можно получить третий столбец путем сложения первых двух столбцов.

A = [1 2 3; 4 5 9; 6 7 13]

A = 3×3

     1     2     3
     4     5     9
     6     7    13

Найдите минимальное нормальное решение для наименьших квадратов для задачи $\mathrm{Ax}\mathit{=}$b, $\mathit{где}$ b равно второму столбцу $\mathit{в}$A. Укажите 'warn' флаг для lsqminnorm для отображения предупреждения, если он обнаруживает, что A имеет низкий ранг.

b = A(:,2);
x = lsqminnorm(A,b,'warn')

Warning: Rank deficient, rank = 2, tol =  1.072041e-14.

x = 3×1

   -0.3333
    0.6667
    0.3333