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

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

Рассмотрите простую линейную систему с одним уравнением и двумя неизвестными, $2{\mathit{x}}_1+3{\mathit{x}}_2=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{x}}_2=-\frac{2}{3}{\mathit{x}}_1+\frac{8}{3}$. Оранжевый круг представляет минимальное расстояние от источника до строки решений, и решение, возвращенное 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