Мат. моделирование / ДЗЗ
Лекция 5

Связь математической теории и прикладных алгоритмов. Тематическая классификация

О чём эта тема

Как снимок ДЗЗ проходит путь от сырых данных до тематической карты: коррекция искажений, повышение контраста и фильтрация, классификация по эталонам и объективная оценка точности результата.

Аннотация

Тема связывает воедино темы предыдущих лекций и практик в единый рабочий процесс обработки снимка ДЗЗ. Разбираются три вида предварительной коррекции — радиометрическая (аппаратные искажения яркости, LUT), атмосферная (эмпирическая коррекция линий, вычитание тёмных объектов, модели переноса излучения) и геометрическая с ортотрансформированием, — а также физика взаимодействия излучения с поверхностью (зеркальное и диффузное отражение). Далее рассматриваются спектральные преобразования контраста (минимаксное растяжение, нормализация и эквализация гистограммы) и пространственная фильтрация как частный случай свёртки из Лекции 4. Основная часть посвящена классификации по эталонным сигнатурам — методам минимального расстояния, расстояния Махаланобиса и спектрального угла — и связанной с ней проблеме линейной разделимости классов, а также месту свёрточных нейронных сетей в этой схеме. Завершает тему оценка точности классификации: матрица ошибок, accuracy, precision и recall, AUC-ROC/AUC-PR и коэффициент каппа. После изучения темы студент сможет назвать и объяснить назначение каждого этапа обработки снимка, выбрать метод классификации по эталонам исходя из формы обучающей выборки и подобрать метрику точности, устойчивую к несбалансированным классам.

Пререквизиты
  • Лекция 4 — оператор свёртки, фильтрация изображений
  • Практика 4 — геопривязка, передискретизация
  • Практика 7–8 — обучающая выборка, спектральные сигнатуры, матрица ошибок
Мотивация

Каждый метод в этой лекции уже встречался в курсе по отдельности — свёртка, геопривязка, SVM, матрица ошибок из Практики 8. Здесь они впервые собираются в один сквозной процесс: без коррекции искажений результат классификации ненадёжен независимо от алгоритма, а без правильно подобранной метрики точности красивая карта классов может скрывать классификатор, который не различает объекты вовсе.

Лекция 5. Связь математической теории и прикладных алгоритмов. Тематическая классификация

Презентация доступна здесь

Лекция связывает разрозненные темы курса в единый рабочий процесс обработки снимка ДЗЗ: от предварительной коррекции искажений через улучшение изображения к классификации и, наконец, к объективной оценке того, насколько классификации вообще можно доверять.

Последовательность обработки данных ДЗЗ

Путь от космического снимка до тематической карты состоит из трёх крупных блоков:

Космическое изображение
Предварительная обработка
Радиометрическая коррекция
Атмосферная коррекция
Геометрическая коррекция, ортотрансформирование
Дешифрирование (тематическая обработка)
Классификация с обучением
Классификация без обучения
Тематическая карта, отчёт

Рис. Л5.1 — последовательность обработки данных ДЗЗ от снимка до тематической карты

Конкретный пример такого перехода — снимок с КА «Ресурс-О1» (сканер МСУ-Э) и построенная по нему лесная карта:

исходный снимок с КА «Ресурс-О1» (сканер МСУ-Э)
Рис. Л5.1 — исходный снимок с КА «Ресурс-О1» (сканер МСУ-Э)
карта лесов, построенная по этому снимку: типы леса, болота, пойма, облачность
Рис. Л5.2 — карта лесов, построенная по этому снимку: типы леса, болота, пойма, облачность

Сам этап дешифрирования (тематической классификации) традиционно раскладывают ещё на три шага:

  1. Выделение признаков. Преобразование многоспектрального снимка пространственным или спектральным преобразованием в характерные признаки изображения. Примеры — выбор подмножества спектральных каналов, PCT (метод главных компонент) для уменьшения размерности данных или фильтр пространственного сглаживания. Этот этап необязателен: многоспектральный снимок можно использовать напрямую, без предварительного выделения признаков.
  2. Обучение. Выбор пикселов для обучения классификатора распознаванию нужных тем (классов), а также для определения границ решения, которые сегментируют пространство признаков в соответствии со свойствами обучающих пикселов. Этот этап либо контролируется аналитиком (классификация с обучением), либо не контролируется им (классификация без обучения).
  3. Маркировка. Применение границ решения пространства признаков к снимку целиком для присвоения метки каждому пикселу. Если обучение было контролируемым, метки уже привязаны к областям пространства признаков; если неконтролируемым — аналитик присваивает метки областям на этом этапе. Итоговая карта содержит по одной метке на каждый пиксел.

В результате маркировки данные цифрового снимка преобразуются в описательные метки, распределяющие материалы поверхности или условия по категориям — то есть исходные данные приобретают информационную величину. Одновременно с классификацией происходит существенное сжатие данных: многоспектральное изображение (сотни каналов по 8 бит/пиксел/канал) уменьшается до размеров карты с не более чем несколькими десятками меток категорий, которую можно хранить в одном канале и менее чем 8 битами на пиксел.

Таблицу перехода от векторов пикселов к меткам называют кодировочной книгой: она используется для кодирования данных на передающей стороне и раскодирования на приёмной. Если каждый вектор признаков связан с уникальной меткой — это сжатие без потерь; если разные векторы получают одну метку — сжатие с потерями. Результативность такой системы сжатия оценивается по тому, насколько хорошо сохраняется исходное содержимое снимка — иными словами, классификатор выполняет аппроксимацию исходного снимка, а точность этой аппроксимации служит объективным численным критерием для сравнения разных классификаторов и оценки точности меток на подходящих испытательных полигонах.

Предварительная обработка изображения

Радиометрическая коррекция

Радиометрическая коррекция исправляет аппаратные искажения яркости, вызванные характеристиками съёмочного прибора, а не самой сценой. Для сканерных приборов такие дефекты видны как полосы, параллельные или перпендикулярные направлению сканирования. Основные причины:

  • неоднородность отклика отдельных детекторов и их элементов;
  • неисправность элементов детектора;
  • потеря данных при передаче, архивировании или извлечении из архива;
  • узкий динамический диапазон;
  • непостоянство параметров съёмки от снимка к снимку.

Коррекцию выполняют двумя способами: по известным параметрам и настройкам съёмочного прибора (корректировочным таблицам, полученным на наземных и лётных испытаниях) или статистически — по самому корректируемому изображению. Первый способ точнее, так как учитывает особенности конкретной системы; второй проще, но менее надёжен.

Отклики элементов сенсора нормализуются с помощью таблицы соответствия (LUT, Look-Up Table), опорным значением для которой обычно служит наименьшая интенсивность сигнала на снимке; та же таблица используется для обратного пересчёта нормализованных данных в исходные абсолютные единицы.

Пример статистического устранения полосчатости на изображении спектрорадиометра ASTER — удаление вертикальных полос:

устранение вертикальных полос на изображении ASTER статистическим методом: слева — исходное изображение с полосой, справа — результат коррекции
Рис. Л5.3 — устранение вертикальных полос на изображении ASTER статистическим методом: слева — исходное изображение с полосой, справа — результат коррекции

и удаление горизонтальных полос (выпадение строк):

устранение горизонтальных полос на изображении ASTER: слева заметны горизонтальные линии сканирования, справа они убраны
Рис. Л5.4 — устранение горизонтальных полос на изображении ASTER: слева заметны горизонтальные линии сканирования, справа они убраны

Для радиометрической коррекции изображений ASTER в видимом и ближнем ИК диапазоне (каналы VNIR) используется линейная формула:

\[ L = A_v \cdot \frac{V}{G_v} + D_v, \]

где \( L \) — скорректированное значение, \( V \) — исходное (корректируемое) значение, а \( A_v, G_v, D_v \) — коэффициенты коррекции, рассчитанные разработчиком датчика. Например, при \( A_v = 1{,}7965 \), \( V = 57 \), \( G_v = 2{,}472 \), \( D_v = -2{,}6339 \):

\[ L = 1{,}7965 \cdot \frac{57}{2{,}472} + (-2{,}6339) \approx 38{,}79. \]

Типичная ошибка

Результат такой коррекции, как в примере выше, — не целое число, и полученные значения уже не совпадают с исходными масштабированными цифровыми числами (Digital Numbers, DN) — это приближение к реальным физическим величинам (яркости или коэффициенту отражения), но ещё не точное их значение: атмосферная коррекция здесь не учтена.

Атмосферная коррекция

Излучение от земной поверхности проходит через атмосферу, прежде чем попасть на датчик спутника, поэтому данные ДЗЗ несут информацию не только о поверхности, но и о состоянии атмосферы в момент съёмки. Устранение влияния атмосферы требует знания количества водяного пара, распределения аэрозолей и видимости сцены; поскольку прямое измерение этих параметров редко доступно, применяют один из трёх подходов:

  • Стандартная абсолютная коррекция — способом эмпирической коррекции линий (по спектральным характеристикам объектов на выбранном участке строится спектральная кривая, из которой находится истинный коэффициент отражения для каждого пиксела) или способом вычитания абсолютно тёмных объектов (по гистограмме находят заведомо тёмные участки — тени в видимом диапазоне, глубокие чистые водоёмы в ближнем ИК — присваивают им коэффициент отражения 0 и вычитают полученную разницу яркости из всех пикселов снимка, предполагая, что она вызвана атмосферой). Недостаток первого способа — обычно неизвестны истинные значения коэффициентов отражения объектов; недостаток второго — риск переисправленных данных.
  • Стандартная относительная коррекция — калибровкой по плоской области спектральной кривой (нормализация значений в пределах участка снимка с ровным ходом спектральной кривой) или делением на среднее значение отражения по изображению. Первый способ неприменим при отсутствии на сцене плоских участков спектральной кривой; второй теряет точность при значительном покрытии растительностью.
  • Коррекция на основе заданных моделей — с помощью программ переноса излучения (например, MODTRAN), реализующих несколько моделей атмосферы и состава аэрозолей и рассчитывающих уникальное решение для каждого снимка; для термальных ИК-радиометров такая коррекция особенно важна, так как погрешность определения яркостной температуры без неё может достигать 10 К и более.

Снимок спутника Landsat без атмосферной коррекции и тот же снимок с атмосферной коррекцией:

снимок Landsat без атмосферной коррекции — заметна дымка, снижающая контраст и искажающая цвет
Рис. Л5.5 — снимок Landsat без атмосферной коррекции — заметна дымка, снижающая контраст и искажающая цвет
тот же снимок Landsat с атмосферной коррекцией — контраст и цветопередача ближе к реальным
Рис. Л5.6 — тот же снимок Landsat с атмосферной коррекцией — контраст и цветопередача ближе к реальным

Геометрическая коррекция и ортотрансформирование

Геометрическая коррекция исправляет искажения формы и масштаба снимка — как определяемые самой съёмочной системой, так и не зависящие от неё: вращение и форма Земли, геометрия фокальной плоскости сенсора, ориентация снимка относительно направления движения спутника, ошибки определения высоты орбиты и углов ориентации аппарата, рассинхронизация бортовой и наземной шкал времени. Для коррекции снимок трансформируют из собственной системы координат в систему координат наземной станции с помощью динамической модели съёмки; систематические искажения устраняются на самих станциях приёма, а искажения относительного позиционирования пикселов (полосчатость, выпадение строк) — статистическими методами, аналогичными радиометрической коррекции.

Геопривязка (математическая часть геометрической коррекции — приведение координат снимка к координатам местности) и передискретизация методами ближайшего соседа, билинейной интерполяции и бикубической свёртки подробно разбираются в Практике 4 вместе с выбором степени полинома и количества опорных точек. Отдельный вид геометрической коррекции — ортотрансформирование: компенсация смещений изображений объектов за счёт рельефа местности, для которой, помимо элементов внешнего и внутреннего ориентирования съёмочной системы, требуется цифровая модель рельефа (ЦМР).

Взаимодействие излучения с подстилающей поверхностью

Возможность различать объекты по данным ДЗЗ основана на том, что горные породы, почвы, вода, растительность и другие материалы по-разному отражают и поглощают электромагнитное излучение на разных длинах волн. При взаимодействии солнечного излучения с поверхностью происходит его поглощение и отражение; отражённое излучение может расходиться по всем направлениям (диффузное отражение) или отражаться преимущественно в одном направлении, как от зеркала (зеркальное отражение) — реальная картина обычно лежит между этими двумя крайними случаями.

Зеркальное отражение — гладкая поверхность, λ много больше шероховатости
Диффузное отражение — шероховатая поверхность, λ много меньше шероховатости

Рис. Л5.2 — два предельных случая отражения; реальные поверхности дают промежуточную картину

Излучение, проходящее через атмосферу, дополнительно ослабляется поглощением в характерных полосах: например, молекула водяного пара даёт линии поглощения на частотах 22,3; 183,4 и 323,8 ГГц, а молекула кислорода — на 119 ГГц и вблизи 60 ГГц. Именно поэтому съёмку ведут в окнах прозрачности атмосферы (см. Лекцию 1), а для количественного анализа отражения земной поверхности необходима атмосферная коррекция, описанная выше.

Улучшение качества изображения

«Пропущенные пикселы» и целые строки со сбойными значениями яркости (например, из-за отказа отдельных детекторов при съёмке или потерь при передаче) мешают тематической обработке снимка; их можно восстановить интерполяцией по соседним пикселам с определённой погрешностью.

Спектральные преобразования

Спектральные преобразования строятся на спектральной диаграмме (гистограмме) — графике зависимости количества пикселов изображения от значения их яркости. Контраст изображения — это разность между максимальным и минимальным значениями яркости; слабый контраст — самый распространённый дефект снимков, и для его повышения применяют несколько методов.

Минимаксное (линейное) растяжение присваивает всем значениям яркости новые значения так, чтобы охватить весь допустимый диапазон. Если исходные значения пикселов лежат в диапазоне \( [I_{\min}, I_{\max}] \), а привести их нужно к диапазону \( [D_{\min}, D_{\max}] \) (для 8-битного изображения обычно \( [0, 255] \)), то новое значение пиксела вычисляется по формуле:

\[ D = (I - I_{\min}) \cdot \frac{D_{\max} - D_{\min}}{I_{\max} - I_{\min}} + D_{\min}. \]

Нормализация гистограммы растягивает на весь диапазон значений яркости не всю гистограмму, а только её наиболее интенсивный (информативный) участок. Эквализация (выравнивание) гистограммы присваивает каждому пикселу новое значение яркости так, чтобы на каждый яркостный уровень приходилось примерно одинаковое количество пикселов — в отличие от линейного растяжения, это нелинейное преобразование, которое может заметно поменять относительный контраст между разными участками диапазона.

Ниже — интерактивная модель: один и тот же синтетический снимок с узким диапазоном яркостей обрабатывается тремя способами.

Интерактив

Способы повышения контраста

Слева — исходное изображение с узким диапазоном яркости и его гистограмма. Справа — результат выбранного преобразования.

Исходное изображение
Результат

Фильтрация

Пространственная фильтрация строится на понятии скользящего окна — квадратной матрицы весовых коэффициентов (обычно 3×3 или 5×5), которая проходит по изображению; пикселу в центре окна присваивается новое значение, вычисленное как функция от значений окружающих его пикселов. Формально это операция свёртки, разобранная в Лекции 4: сглаживающий фильтр (все веса равны и в сумме дают 1, например \( 1/9 \) в каждой из девяти ячеек) размывает шум, фильтры обнаружения линий и границ (веса с чередующимся знаком, например центральная строка или столбец с весом 2 и соседи с весом −1) подчёркивают резкие перепады яркости — опробовать оба типа фильтров и их отличие от нелинейной пороговой обработки можно на площадке операторов свёртки из Лекции 4.

Дешифрирование и классификация изображений

Дешифрирование — процесс обнаружения и распознавания объектов и явлений местности на снимке; оно бывает визуальным (ручным) и машинным (автоматическим). Машинное дешифрирование сводится к классификации: каждый пиксел представляется вектором в пространстве спектральных признаков, и по количественным связям спектральных яркостей разных объектов пикселы разделяются по классам. Классификация делится на классификацию с обучением, где аналитик заранее указывает эталонные участки («сигнатуры») известных классов, и классификацию без обучения, где границы классов ищутся автоматически без заранее заданных эталонов (см. Практику 6 и Практику 7 про формирование обучающей выборки и работу с сигнатурами в SCP).

При классификации с обучением, имея эталонные сигнатуры, новый пиксел относят к тому классу, до эталона которого он «ближе» в пространстве признаков — расстояние при этом можно определять разными способами:

  • Метод минимального расстояния. Между вектором пиксела и вектором каждого эталона вычисляется евклидово расстояние — корень из суммы квадратов покоординатных разностей. Пиксел относят к классу ближайшего эталона, если расстояние до него не превышает заданный порог.
  • Метод расстояния Махаланобиса. Похож на предыдущий, но вместо евклидова расстояния использует расстояние, учитывающее дисперсию (разброс) значений эталонной выборки по каждому направлению. Если евклидовы расстояния до двух эталонов равны, побеждает класс с большей дисперсией — то есть с более «размытой» и потому более вероятной областью значений.
  • Метод спектрального угла. Классы сравниваются не по расстоянию между векторами, а по углу между вектором-эталоном и вектором пиксела (см. также Практику 7, где этот угол используется для оценки разделимости сигнатур); пиксел относят к классу эталона, угол до которого меньше заданного порога.

Ниже — интерактивная модель, показывающая, почему эти три метода могут давать разные ответы для одного и того же пиксела.

Интерактив

Три способа сравнить пиксел с эталоном

Три готовых сценария показывают, когда какая метрика точнее. Кликните по полю, чтобы передвинуть тестовый пиксел (красный), и посмотрите, какой класс выберет каждый метод.

Проблема разделимости классов

Проблема разделимости возникает, когда объекты разных классов невозможно разделить в пространстве признаков линейно или другой сравнительно простой функцией. Если классы линейно разделимы, достаточно простых методов вроде описанных выше; если нет — нужны более сложные алгоритмы: методы с ядрами, такие как SVM (см. Практики 7–8), или нейронные сети.

Нейронные сети в обработке данных ДЗЗ

Искусственная нейронная сеть способна аппроксимировать любую непрерывную функцию при нелинейных функциях активации; её поведение определяется структурой сети и весами связей. В основе большинства нейросетевых архитектур для ДЗЗ лежат свёрточные слои: свёрточная нейронная сеть (СНС) состоит из свёрточных слоёв, субдискретизирующих (уменьшающих размерность) слоёв и слоёв «обычного» персептрона.

Свёрточный слой — это ровно та операция свёртки, что разобрана в Лекции 4: фильтр поэлементно умножается на значения пикселов под ним, произведения суммируются, и результат становится одним числом на выходной карте признаков; фильтр сдвигается по всему изображению, формируя карту отклика этого признака. Разница с классической фильтрацией в том, что веса фильтра в нейросети не задаются вручную, а подбираются в процессе обучения так, чтобы полученные признаки были полезны для конечной задачи (например, картирования засоления почв по обучающей выборке эталонных участков с известной степенью засоления).

Оценка точности классификации

Итоговую точность классификации оценивают матрицей ошибок (confusion matrix) — она подробно разобрана в Практике 8 вместе с ошибками первого и второго рода. На основе той же матрицы для класса «положительный» / «отрицательный» строится несколько сводных метрик.

Accuracy — доля правильных ответов алгоритма среди всех ответов:

\[ \text{Accuracy} = \frac{TP + TN}{TP + TN + FP + FN}. \]

Эта метрика интуитивно понятна, но вводит в заблуждение при несбалансированных классах: если, например, 95% территории — это один класс, классификатор, всегда выдающий этот класс, получит accuracy 95%, не различая объекты вовсе.

Precision (точность) и Recall (полнота) разбирают ошибки отдельно по направлениям:

\[ \text{Precision} = \frac{TP}{TP + FP}, \qquad \text{Recall} = \frac{TP}{TP + FN}. \]

Precision — это доля объектов, названных классификатором положительными, которые действительно положительны; Recall — доля реальных положительных объектов, которые классификатор нашёл. У этих метрик обратная зависимость: увеличивая порог уверенности, можно поднять Precision за счёт Recall и наоборот.

AUC-ROC и AUC-PR оценивают классификатор в целом, не привязываясь к конкретному порогу перевода вероятности в метку класса (естественный порог 0,5 не всегда оптимален, особенно при дисбалансе классов). ROC-кривая строится в координатах True Positive Rate и False Positive Rate при переборе всех возможных порогов, PR-кривая — в координатах Precision и Recall; AUC (area under curve) — площадь под соответствующей кривой, от 0,5 (случайное угадывание) до 1 (идеальный классификатор).

TPR FPR ROC случайный классификатор

Рис. Л5.3 — ROC-кривая: чем дальше от диагонали к верхнему левому углу, тем лучше классификатор при любом пороге

Коэффициент каппа (κ) оценивает согласие между классификацией и эталоном сверх того, что можно было бы получить случайным совпадением:

\[ \kappa = \frac{p_o - p_e}{1 - p_e}, \]

где \( p_o \) — наблюдаемая доля совпадений (аналог accuracy), а \( p_e \) — ожидаемая доля случайных совпадений при данном распределении классов в эталоне и в классификации. Чем выше κ, тем лучше алгоритм классификации согласуется с эталоном сверх случайного совпадения; в отличие от accuracy, каппа учитывает и то, насколько несбалансированы классы.

Типичная ошибка

Высокая accuracy сама по себе не доказывает качество классификации при несбалансированных классах (пример выше с 95% одного класса). Отчёт о точности стоит всегда сопровождать хотя бы одной метрикой, чувствительной к дисбалансу, — Recall для редкого класса, AUC-PR или коэффициентом каппа, а не только общей долей правильных ответов.

Контрольные вопросы