Автономное картирование геологии карьерных уступов с использованием гиперспектральных данных: сравнительный анализ трёх алгоритмов классификации
Авторы: S. Schneider, R.J. Murphy, A. Melkumyan, E. Nettleton
Источник: Proceedings of APCOM 2011 (35th International Symposium on Application of Computers and Operations Research in the Minerals Industry)
Организация: Australian Centre for Field Robotics (ACFR), University of Sydney; Rio Tinto Centre for Mine Automation
Аннотация
Идентификация геологии и минералогии непосредственно на уступе карьера является необходимым условием для автоматизации горных работ — от разведки до сортировки руды и учёта запасов. Гиперспектральные данные, получаемые с наземных платформ, идеально подходят для этих задач. Однако изменчивость солнечного освещения вносит в сигнал вариации, не связанные с геологией, что усложняет анализ. Алгоритмы должны работать в вероятностном фреймворке (для последующего слияния данных) и быть устойчивыми к переменному освещению.
В работе сравнены три алгоритма классификации гиперспектральных данных: Spectral Angle Mapper (SAM), гауссовский процесс с квадратичной экспоненциальной ковариационной функцией (GP‑SE) и гауссовский процесс с ковариационной функцией, зависящей от угла наблюдения (GP‑OAD). Два эксперимента — на лабораторной спектральной библиотеке и на реальных изображениях вертикального уступа — показали: лабораторная оценка может вводить в заблуждение (GP‑SE хорош в контролируемых условиях, но не работает в поле). GP‑OAD сочетает устойчивость к изменениям освещённости (благодаря использованию спектрального угла), вероятностный выход и учёт внутриклассовой изменчивости, что делает его предпочтительным для автономных горных систем.
1. Введение
В открытых карьерах многих стран мира (Австралия, США, Канада, Чили, ЮАР) добываются золото, железная руда, уголь, медь. Карьеры часто расположены в труднодоступных или опасных для человека районах. Одна из наиболее рискованных операций — картирование геологических границ на вертикальных уступах перед началом выемки. Нестабильность уступов создаёт прямую угрозу для геологов. Следовательно, существует острая потребность в автоматизации таких задач.
Дистанционное зондирование позволяет картировать геологию с безопасного расстояния. Гиперспектральные данные предпочтительнее мультиспектральных, поскольку они обеспечивают идентификацию минералов (Kurz et al., 2008) и оценку их содержаний (Cudahy & Ramanaidou, 1997). Гиперспектральные сенсоры измеряют отражённый свет во многих узких непрерывных каналах, что позволяет разрешать диагностические полосы поглощения железа, глин, карбонатов.
Съёмка с наземных платформ необходима для достижения требуемого пространственного разрешения и ориентации на вертикальных уступах. Однако открытые карьеры создают серьёзные вызовы:
сложная геометрия уступов приводит к глубоким теням;
изменение положения солнца в течение дня вызывает сильную вариабельность падающего освещения.
Для использования гиперспектральных данных в автономных горных системах алгоритмы должны удовлетворять четырём требованиям:
Вероятностный выход — каждый пиксель карты должен иметь ассоциированную вероятность для последующего слияния с данными других сенсоров.
Устойчивость к изменчивости освещения (включая тени и затенение).
Независимость от конкретного уступа — классификация должна быть возможна с использованием независимо полученной спектральной библиотеки (без априорного знания об уступе).
Вычислительная эффективность для работы в реальном времени.
Ранее (Schneider et al., 2010) было показано, что метод GP‑OAD (Gaussian Process с Observation Angle Dependent ковариационной функцией) превосходит классический SAM для данных, полученных с невизуализирующим спектрометром в контролируемых условиях. В настоящей статье эта работа расширена: проведено прямое сравнение GP‑OAD, GP‑SE и SAM на лабораторных данных и на реальных гиперспектральных изображениях уступа, снятых в разное время суток.
2. Материалы и методы
2.1. Район исследований и геология
Все образцы отобраны в действующем железорудном карьере провинции Хамерсли (Западная Австралия). Район характеризуется обширными полосатыми железными формациями (BIF), которые в некоторых местах подверглись рудной минерализации за счёт выветривания и выщелачивания грунтовыми водами. При этом кремнезём (основной компонент BIF) выщелачивается, а железо концентрируется в форме гётита и мартита (разновидность гематита).
Другие основные породы:
сланец (включая тонкие прослои вулканического сланца и обширные отложения West Angeles Shale);
марганцовистый сланец, содержащий каолинит и/или галлуазит с прослоями пиролюзита;
водореактивная глина (каолинит) — небольшие изолированные пятна.
2.2. Гиперспектральные данные
2.2.1. Спектральная библиотека
Независимая библиотека отражения (0,35–2,5 мкм) создана с использованием полевого спектрорадиометра FieldSpec3 Pro (ASD Inc., Боулдер, Колорадо, США) по репрезентативным образцам из провинции Хамерсли. Образцы двух типов:
керны алмазного бурения (ширина 7,5 см) — отобраны как наиболее чистые спектры для каждого типа породы;
целые образцы с уступов — для учёта изменчивости (запылённость, выветрелость).
Искусственное освещение обеспечивало оптимальное отношение сигнал/шум и стандартные условия. Калибровка по панели Spectralon (99% отражения). Спектры преобразованы к полосам пропускания визуализирующих сенсоров (Specim) с помощью свёртки Гаусса. Библиотека содержит 605 спектров, неравномерно распределённых по семи типам пород (Таблица 1).
Таблица 1. Состав спектральной библиотеки
Тип породы
Акроним
Количество спектров
Полосчатая железная формация
BIF (purple)
87
Кремнистый BIF
CHT (gray)
40
Мартит-гётит (тонкозернистый, твёрдый)
MAR (yellow)
206
Водореактивная глина
WRC (white)
95
Гётит-лимонит
GOL (blue)
61
Сланец
SHL (red)
69
Марганцовистый сланец
SHN (green)
47
2.2.2. Гиперспектральные изображения уступа
Изображения получены с двух гиперспектральных сканеров Specim (Финляндия), установленных на вращающейся платформе:
VNIR (0,4–0,97 мкм), спектральное разрешение FWHM = 2,22 нм;
SWIR (0,97–2,5 мкм), FWHM = 6,35 нм.
Калибровочная панель Spectralon (99%) размещена перед уступом. Изображения калиброваны на отражение по стандартному протоколу (Murphy et al., 2008). Диапазоны с атмосферными поглощениями (1,45 и 1,85 мкм) удалены.
Съёмка выполнена в разное время суток:
11:30 — «без теней» (Рис. 1a в оригинале);
13:30 — «затенённое» (Рис. 1b).
Уступ имеет ориентацию и геометрию, вызывающую значительные изменения угла и интенсивности падающего света. Большая часть уступа покрыта рыхлым материалом (осыпи), поэтому геологические границы не очень чёткие. Выделены две крупные зоны: руда (мартит) и пустая порода (сланец), граница определена традиционной геологической съёмкой.
2.3. Алгоритмы классификации
2.3.1. Spectral Angle Mapper (SAM)
SAM вычисляет угол между целевым и эталонным спектрами, рассматривая каждый спектр как D-мерный единичный вектор в ℝᴰ:
Обычно полагают μ(x)=0 (центрирование данных). Совместное распределение обучающих и тестовых значений имеет вид:
Здесь K — ковариационная матрица, вычисленная с помощью ядра, σn2 — дисперсия шума, I — единичная матрица.
Условное распределение для тестовых точек:
Обучение GP заключается в нахождении гиперпараметров ядра путём максимизации логарифма маргинального правдоподобия:
2.3.3. Квадратичная экспоненциальная ковариационная функция (GP‑SE)
Наиболее распространённое ядро:
где σ0 — масштабный параметр (амплитуда), li — характеристический масштаб по каждому измерению. Из-за высокой размерности гиперспектральных данных (десятки-сотни каналов) оптимизация всех li вычислительно сложна, поэтому в данной работе принято li=l для всех i.
Важно: ядро использует евклидово расстояние между векторами. Это делает его чувствительным к изменениям общей яркости (смещению спектра как вектора).
2.3.4. Ковариационная функция, зависящая от угла наблюдения (GP OAD)
Предложена Melkumyan и Nettleton (2009). Ядро зависит только от угла между спектрами α(x,x'):
Здесь σ0, l и ϕ — гиперпараметры, подлежащие обучению. Угол α вычисляется как в SAM (формула 1). Поскольку ядро использует только угол, оно инвариантно к яркости (нормам векторов). При этом GP OAD, в отличие от SAM, моделирует внутриклассовую изменчивость и даёт вероятностный выход.
2.4. Процедура классификации
2.4.1. SAM для изображений
Каждый пиксель изображения сравнивается со всеми 605 спектрами библиотеки. Класс назначается по минимальному углу (без порога). Такой подход позволяет оценить «чистую» разделяющую способность без влияния порога.
2.4.2. GP методы: обучение и предсказание
Реализована схема «один против всех» (One vs All). Для каждого из 7 типов пород строится бинарный классификатор: спектры данного типа — положительный класс (+1), все остальные — отрицательный (−1). Гиперпараметры для каждого бинарного GP находятся максимизацией логарифма правдоподобия (итерации градиентного спуска с 10 случайными инициализациями для избежания локальных минимумов).
Для тестового пикселя x* каждый из 7 GP выдаёт вещественное предсказание E[f*] и дисперсию. Для получения дискретной метки класса вычисляется вероятность принадлежности к классу c через сигмоидальное преобразование или, как в данной работе, назначается класс с максимальным значением предсказания. Для классификации изображений метка присваивается по максимальной вероятности среди 7 классов.
2.4.3. Перекрёстная проверка (эксперимент 1)
10-кратная перекрёстная проверка с сохранением стратификации. Для каждого разбиения вычисляются точность (accuracy) и F мера:
Для SAM перекрёстная проверка невозможна (нет этапа обучения). Вместо этого каждый спектр по очереди используется как эталон для классификации всех остальных; накопление результатов в единой матрице ошибок даёт средние показатели для сравнения.
2.5. Оценка влияния освещения (эксперимент 2)
Классификация изображений «без теней» и «затенённого» выполняется с использованием одной и той же спектральной библиотеки. Сравниваются: 1) визуальное соответствие геологии; 2) доля пикселей, изменивших класс между двумя изображениями. Идеальный алгоритм должен давать одинаковую карту, независимо от времени съёмки.
3. Результаты
3.1. Эксперимент 1 (лабораторные данные)
Таблица 2. Средняя точность и F‑мера по 7 классам (в процентах).
Алгоритм
Точность (accuracy)
F мера
SAM (порог 0,1)
80%
40%
GP SE
95%
81%
GP OAD
97%
88%
Различие между GP‑SE и GP‑OAD статистически незначимо (p = 0,32 по парному t‑тесту). Оба метода существенно превосходят SAM (p < 0,01). Низкая F‑мера SAM объясняется большим числом ложных срабатываний: порог 0,1 слишком жёсткий для одних классов и слишком мягкий для других.
3.2. Эксперимент 2 (изображения уступа)
SAM (Рис. 3 в оригинале):
На изображении 11:30 — зоны руды и сланца идентифицированы верно, но присутствует шум (пиксели GOL в центре, BIF слева, что не подтверждено полевой съёмкой).
На изображении 13:30 — шум усиливается, особенно в тенях. 29% пикселей изменили класс между двумя изображениями (Таблица 3).
GP‑SE (Рис. 4):
Качественно неверно: вся рудная зона и сланец отнесены к MAR (мартит). Лишь небольшое пятно SHL (сланец) выделено внутри сланцевой зоны. Это полное несоответствие реальной геологии.
Формально только 6% пикселей изменили класс, но это следствие того, что классификация была систематически ошибочной на обоих изображениях.
GP‑OAD (Рис. 5):
Карты на обоих изображениях качественно соответствуют геологии: рудная зона (MAR) и сланцевая зона (SHL, SHN) разделены корректно. Шум значительно ниже, чем у SAM.
16% пикселей изменили класс между изображениями. Важно: изменения происходят в основном внутри сланцевой зоны (переход SHL ↔ SHN) и не затрагивают границу руда/сланец. Рудная зона остаётся стабильной.
Таблица 3. Доля пикселей, изменивших класс между изображениями 11:30 и 13:30
Алгоритм
Абсолютное число пикселей
Относительная доля
GP OAD
51449
0,1587 (15,9%)
GP SE
20637
0,0637 (6,4%)
SAM
93793
0,2892 (28,9%)
4. Обсуждение
4.1. Почему лабораторные тесты вводят в заблуждение?
В эксперименте 1 все спектры получены при искусственном стандартизированном освещении с фиксированной геометрией. В этих условиях яркость варьирует главным образом из-за внутренних свойств образцов, а не из-за внешних факторов. Поэтому GP SE, использующий евклидово расстояние, работает почти так же хорошо, как GP OAD.
В полевых условиях (эксперимент 2) яркость пикселя может измениться в 2–3 раза только из-за тени, при этом форма спектра (диагностические полосы поглощения) может оставаться неизменной. Евклидово расстояние между затенённым полевым спектром и ярким лабораторным эталоном становится большим, даже если это одна и та же порода. GP SE не способен это компенсировать и выдаёт ошибочную классификацию.
4.2. Преимущества угловой метрики
SAM и GP OAD используют спектральный угол, который инвариантен к масштабированию вектора (яркости). Именно поэтому они справляются с задачей, а GP SE — нет. Однако SAM имеет принципиальные ограничения:
дискретный выход (нет оценки неопределённости);
использует только один эталонный спектр на класс, игнорируя внутриклассовую изменчивость;
требует задания порога, который невозможно выбрать универсально для всех типов пород и условий съёмки.
GP OAD снимает эти ограничения, оставаясь в вероятностном фреймворке и используя всю обучающую библиотеку.
4.3. Почему GP OAD лучше SAM, несмотря на общую угловую метрику?
Моделирование изменчивости: GP OAD обучается на множестве спектров каждого класса, оценивая гиперпараметры σ0,l,ϕ, которые адаптируются к разбросу углов внутри класса и между классами.
Вероятностный выход: Для каждого пикселя выдаётся не только наиболее вероятная метка, но и мера уверенности (дисперсия). Это критически важно для последующего слияния с данными LiDAR, радара, фотограмметрии в автономных системах.
Отсутствие порога: Классификация происходит автоматически на основе максимальной вероятности, без ручного подбора углового порога.
4.4. Остаточная чувствительность к теням GP OAD (16% изменений)
Большинство изменений пикселей между изображениями в GP OAD происходит внутри сланцевой зоны между двумя различными типами сланцев (SHL и SHN). Это указывает на то, что спектры SHL и SHN очень похожи (близкие углы), и изменение освещения может смещать классификацию с одного на другой. Рудная зона (MAR) при этом остаётся стабильной, что подтверждает робастность метода для основного производственного различия «руда vs. пустая порода».
5. Дальнейшие исследования (продолжение работы)
Авторы указывают следующие направления для улучшения классификации затенённых областей и повышения общей надёжности:
Адаптивное моделирование тени
Разработка методов, которые явно детектируют затенённые пиксели (по интенсивности или по углу съёмки) и применяют к ним отдельную модель освещения. Например, можно ввести дополнительный гиперпараметр в ковариационную функцию, учитывающий угол падения солнечного света.
Слияние с данными LiDAR
Лазерное сканирование даёт высокоточную 3D-геометрию уступа. Зная ориентацию каждого пикселя относительно солнца и сенсора, можно компенсировать вариации освещения аналитически. GP OAD предоставляет вероятностный выход, что делает его естественным кандидатом для байесовского слияния с LiDAR.
Расширение спектральной библиотеки
Библиотека должна включать спектры, измеренные при разных углах освещения и разных уровнях яркости (не только в стандартных лабораторных условиях). Это позволит GP OAD обучаться на более репрезентативных данных.
Переход от пиксельной к объектной классификации
Вместо классификации каждого пикселя независимо — использование пространственного контекста (сегментация, CRF). Ожидается, что это снизит шум и улучшит стабильность границ между геологическими зонами.
Масштабирование на большие площади
Вычислительная сложность GP составляет O(N3) для обучения и O(N2) для предсказания, что при N=605 приемлемо, но при увеличении библиотеки до тысяч спектров требует аппроксимаций (индуктивные точки, разреженные GP).
Интеграция в реальное время
Реализация GP OAD на бортовом компьютере для автоматических буровых станков или роботизированных экскаваторов. Требуется оптимизация вычислений и, возможно, использование GPU.
6. Заключение
В работе проведено систематическое сравнение трёх алгоритмов классификации гиперспектральных данных для картирования геологии вертикальных карьерных уступов. Показано, что:
Оценка только на лабораторных данных недостаточна: GP SE демонстрирует высокую точность в контролируемых условиях, но полностью неработоспособен в полевых условиях из-за чувствительности к изменениям яркости.
Метрики на основе спектрального угла (SAM, GP OAD) критически необходимы для обеспечения инвариантности к освещению.
GP OAD сочетает преимущества угловой метрики, вероятностного выхода и учёта внутриклассовой изменчивости, что делает его превосходящим как SAM (дискретный, невероятностный), так и GP SE (неустойчив к освещению).
GP OAD успешно классифицирует зоны руды и пустой породы на изображениях, полученных в разное время суток, с приемлемой стабильностью (84% пикселей не меняют класс).
Рекомендуется использовать GP OAD в качестве базового алгоритма для систем автономного картирования в горной промышленности, а также для других задач дистанционного зондирования, где освещение вариабельно (лесное хозяйство, экологический мониторинг, сельское хозяйство).