Кратко — практический пошаговый протокол с пояснениями и ключевыми формулами.
1) Подготовка данных (геометрия и радиометрия)
- СфМ/фотограмметрия: получить точно орторектифицированные ортофото и цифровую модель рельефа (DSM/DTM) для каждой съёмки. Обязательно точная ко-регистрация растров между сезонами.
- Конвертация в поверхностную отражательную способность (surface reflectance): либо через калибровочную панель и эмпирическую линию, либо атмосферную коррекцию (6S/ATCOR/Sen2Cor) — чтобы сравнение NDVI было физически обосновано.
2) Учёт топографии (топографическая коррекция)
- Вычислить локальный угол падения луча (угол инциденции) $i$ по DEM, склону $\alpha$, аспекту $\varphi$ и положению Солнца (зенитный угол ${\theta }_s$, азимут ${\varphi }_s$):
$$\cos i=\cos {\theta }_s\cos \alpha +\sin {\theta }_s\sin \alpha \cos ({\varphi }_s-\varphi ).$$ - Наиболее распространённые коррекции:
- Косинусная (простая): скорректированная отражательная способность ${\rho }_c=\rho \frac{\cos {\theta }_s}{\cos i}$ (работает плохо в тенях/крутых склонах).
- C-correction (рекомендую): для каждой полосы оценивается параметр $C$ по регрессии $\rho$ на $\cos i$, затем
$${\rho }_c=\rho \frac{\cos {\theta }_s+C}{\cos i+C}.$$ - Minnaert: ${\rho }_c=\rho {\left(\frac{\cos {\theta }_s}{\cos i}\right)}^k$, где $k$ — параметр, оцениваемый по данным.
Выбор: C-correction чаще даёт стабильные результаты для мультиспектральных снимков.
3) Обработка теней и угла обзора
- Сформировать маску теней по порогу NIR или по $\cos i$ (низкие значения указывают на недоосвещённые пиксели). В тенях NDVI корректируется с осторожностью или исключается.
- Для БПЛА: учёт вариации угла обзора (BRDF) — если есть данные с разных углов, выполнить BRDF-нормализацию (Ross-Li или эмпирические модели), иначе планировать съёмки с одинаковой высотой и солнцегеометрией.
4) Сравнение по фенологическим стадиям
- Если возможно, сравнивать съёмки, снятые в одинаковых фазах вегетации (посев, вегетация, уборка). Если даты различаются, применять временную нормализацию/моделирование сезонности (см. п.6).
5) Многосезонный анализ NDVI и выделение сигналов
- Посчитать NDVI для каждой съёмки после топокоррекции: $\mathrm{NDVI}=\frac{NIR-RED}{NIR+RED}$.
- Первичный анализ: разности NDVI между аналогичными стадиями (год1 vs год2) — уменьшает влияние сезонности.
6) Разделение сезонности, рельефа и реального изменения (методы)
- Временная декомпозиция: модель сезонности + тренд + остаток. Примеры:
- Гармоническая регрессия (период $T$, например год):
$$NDVI(t)=a_0+\sum _{n=1}^N\left[a_n\cos \left(\frac{2\pi nt}{T}\right)+b_n\sin \left(\frac{2\pi nt}{T}\right)\right]+ϵ(t).$$ - STL / LOESS / TIMESAT / BFAST для разделения тренда и сезонной компоненты.
- Регрессионный способ (учесть влияние рельефа статистически): строите модель
$$NDVI=f_{season}(t)+{\beta }_{1}slope+{\beta }_{2}aspect+{\beta }_{3}illumination+\epsilon ,$$ где $f_{season}(t)$ — сезонная компонента (напр. гармоническая), а остаток $\epsilon$ содержит реальный сигнал изменения растительности. Можно использовать GLM/GAM для нелинейных связей.
- Простая проверка корреляции: вычислить пространственную корреляцию изменения NDVI с картами наклона/аспекта/инсоляции; сильная корреляция — сигнал топографии.
- Флаги ненадёжности: помечать участки с крутым уклоном/постоянными тенями/низкой освещённостью как сомнительные.
7) Валидация и пороги значимости
- Использовать поля/контрольные точки или PIF (pseudo-invariant features) для контроля удачности радиометрической нормализации.
- Статистически тестировать изменения (t-тест/Bootstrap) на уровне пикселя или объектов, учитывать пространственную корреляцию.
8) Практические рекомендации и порядок работ
1. Орто/DEM и ко-регистрация.
2. Привести все ортофото к surface reflectance.
3. Топографическая коррекция (C-correction) + маска теней.
4. Вычислить NDVI.
5. Нормализация между датами (PIF/histogram matching), либо моделирование сезонности.
6. Применить временную декомпозицию / регрессию с переменными рельефа.
7. Отобрать значимые остатки/персистентные изменения и верифицировать полем.
Коротко: коррекция отражательной способности и топографии (формулы выше) + сезонная модель (гармоники/Timesat/BFAST) + регрессия NDVI на параметры рельефа дает разделение вкладов рельефа, сезонности и реального изменения растительности.