Как узнать, является ли матрица сингулярной?
Как мы можем узнать, что матрица 4x4 является единственной или нет?
Можем ли мы узнать это, не дополняя данную матрицу единичной матрицей, а затем выполняя операции со строками?
4 ответа
Так как же определить, является ли матрица действительно единственной? (В MATLAB, без использования бумаги и карандаша, или символьных вычислений, или рукописных операций со строками?) Учебники иногда говорят ученикам использовать определитель, поэтому я начну с этого.
Теоретически, можно просто проверить, равен ли ноль определителю вашей матрицы. таким образом
M = magic(4)
M =
16 2 3 13
5 11 10 8
9 7 6 12
4 14 15 1
det(M)
ans =
-1.4495e-12
Как оказалось, эта матрица действительно единственная, поэтому есть способ записать строку M как линейную комбинацию других строк (также верно для столбцов.) Но мы получили значение, которое не было точно нулевым от det, Это действительно ноль, а MATLAB только что запутался? Это почему? Символический определитель действительно равен нулю.
det(sym(M))
ans =
0
Как выясняется, вычисление детерминанта - ужасно неэффективная вещь для больших массивов. Таким образом, хорошей альтернативой является использование произведения диагональных элементов конкретной матричной факторизации нашего квадратного массива. Фактически, это то, что MATLAB делает внутри себя для несимвольных входных данных.
[L,U,P] = lu(M)
L =
1 0 0 0
0.25 1 0 0
0.3125 0.76852 1 0
0.5625 0.43519 1 1
U =
16 2 3 13
0 13.5 14.25 -2.25
0 0 -1.8889 5.6667
0 0 0 3.5527e-15
P =
1 0 0 0
0 0 0 1
0 1 0 0
0 0 1 0
Посмотрите, что диагональные элементы L являются единицами, но U имеет ненулевые диагональные элементы. И есть хорошие свойства в отношении определителя произведения матриц, а также определителя верхних или нижних треугольных матриц.
prod(diag(U))
ans =
-1.4495e-12
Видите, мы получили тот же ответ. Поскольку факторизация LU достаточно быстрая даже для больших массивов, это хороший способ вычислить определитель. Проблема в том, что он использует арифметику с плавающей запятой. Таким образом, эти диагональные элементы U являются действительными числами с плавающей точкой. Когда мы берем продукт, мы получаем результат, который точно не равен нулю. Этот последний элемент был немного ненулевым.
Есть другие проблемы с дет. Например, мы знаем, что матричный глазок (100) чрезвычайно хорошо обусловлен. В конце концов, это матрица идентичности.
det(eye(100))
ans =
1
Теперь, если мы умножим матрицу на константу, это НЕ изменит статус матрицы как единичной. Но что происходит с определителем?
det(.1*eye(100))
ans =
1e-100
Так это матрица единственного числа? В конце концов, если det(magic(4)) дает нам 1e-12, то 1e-100 должно соответствовать особой матрице! Но это не так. Еще хуже,
det(.0001*eye(100))
ans =
0
Фактически, детерминант масштабируется на 0,0001^100, что в Matlab будет 1e-400. По крайней мере, это было бы верно, если бы matlab мог представлять малое число, используя double. Это не может сделать это. Это число будет занижено. Или так же легко, мы можем сделать это переполнением.
det(10000*eye(100))
ans =
Inf
Ясно, что все эти масштабированные матрицы тождественности одинаково неособы, но можно сделать вывод, чтобы дать нам любой ответ, который мы хотим видеть! Поэтому мы должны сделать вывод, что вычисление определителя - это ужасная вещь, которую нужно сделать с матрицей. Мне все равно, что ваш учебник сказал вам давно, или что сказал вам ваш начальник или учитель. Если кто-то сказал вам вычислить определитель для этой цели, используя компьютер, это был ужасный совет. Период. Детерминанты просто имеют слишком много проблем.
Мы можем сделать другие вещи, чтобы проверить единственность. Лучший инструмент - использовать ранг. Таким образом, если ранг матрицы NxM меньше min(N,M), то матрица является особой. Вот пара тестов:
rank(M)
ans =
3
rank(.0001*eye(100))
ans =
100
Таким образом, ранг может сказать нам, что магический квадрат 4х4 является единичным, но наша масштабированная единичная матрица не является единственной. (Приятно то, что ранг может проверять сингулярность неквадратной матрицы.)
Мы также можем использовать cond для проверки числовой особенности. Наименьшее возможное число условий равно 1,0, что соответствует матрице с очень хорошим поведением. Большие числа условий плохие. В двойной точности это означает, что когда номер условия больше, чем приблизительно 1e15, ваша матрица очень проблематична.
cond(M)
ans =
8.148e+16
cond(.0001*eye(100))
ans =
1
Фактически, cond понимает, что масштабированная единичная матрица все-таки хорошо обусловлена. Большие числа условий плохие. Для матрицы двойной точности число условия, которое находится где-то около 1e15 или около того, указывает матрицу, которая, вероятно, численно сингулярна. Итак, мы видим, что М явно в единственном числе. Опять же, cond может работать с неквадратными матрицами.
Или мы могли бы использовать rcond, инструмент, который оценивает обратную величину числа условий. Это хороший инструмент для действительно больших массивов. Когда rcond дает число, близкое к eps, СМОТРЕТЬ!
rcond(M)
ans =
1.3061e-17
rcond(.0001*eye(100))
ans =
1
Наконец, для математических головок передач (как я), мы могли бы вытащить SVD. В конце концов, svd - это инструмент, на котором основаны cond и rank. Когда одно или несколько сингулярных значений матрицы крошечные по сравнению с наибольшим сингулярным значением, мы снова имеем сингулярность.
svd(M)
ans =
34
17.889
4.4721
4.1728e-16
Здесь мы рассмотрим, когда сингулярное значение мало по сравнению с наибольшим сингулярным значением матрицы. Приятно то, что svd может сказать нам, насколько близка матрица к сингулярности, и, если имеется более одного небольшого сингулярного значения, если дает нам информацию о ранге матрицы.
Приятно то, что ни один из показанных мною инструментов не требует от пользователя выполнения элементарных операций со строками или чего-либо необычного.
Но, пожалуйста, не используйте DET! Да, это проявляется в учебниках. Да, может быть, ваш инструктор или ваш начальник сказал вам использовать его. Они были просто неправы, потому что подобные инструменты не работают, когда они применяются на компьютере, который использует арифметику с плавающей запятой. И вы просто не хотите вычислять символический определитель, который будет ужасно неэффективным.
Извините, если это было долгое чтение. Я сейчас сойду с мыла.
Рассчитайте ранг и сравните с измерением. Если ранг меньше размерности, то матрица является сингулярной.
Наиболее надежный подход - выполнить разложение по сингулярным числам на матрицу. Отношение наибольшего к наименьшему единственному значению должно находиться в пределах некоторого разумного допуска. Это соотношение является условным числом матрицы. При значениях двойной точности все становится очень сложным при значениях двойной точности, когда число условий превышает миллион или более, и это довольно высокий предел. Обратите внимание: если у вас есть SVD, это полезно для многих других вещей, помимо простого вычисления номера условия.
Разложение по сингулярным числам - это швейцарская армейская бензопила численного анализа; это может быть немного тяжелым инструментом, если вы знаете, что матрица не единственная / плохо обусловленная. Но если вы не знаете, это хороший инструмент, чтобы знать. Особенно так с Matlab, так как это встроенный инструмент.
Я хотел бы использовать cond, Это дает численную оценку того, насколько близка к особой матрица (где Inf - особая матрица).
Например:
m = randn(4);
cond(m) %Well conditioned, usually in the 10's
m = diag([1e-6 1 2 1e6]);
cond(m) %Less well conditioned, 1e12
m = diag([0 1 2 3]);
cond(m) %Singular: Inf