Особенности дистанционной диагностики агрохимических свойств пахотных почв

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

Картографирование параметров почвенного плодородия проводится преимущественно на основе данных полевых работ и лабораторных исследований образцов почв. Работы ведутся в соответствии с утвержденными методиками (Методические указания…, 2003) и предполагают отбор и анализ большого количества почвенных образцов, что требует существенных временных и финансовых затрат. С другой стороны, пространственная частота пробоотбора при агрохимическом обследовании почв не доста- точна для отображения пространственного варьирования показателей плодородия пахотных почв, особенно элементов питания, содержание которых может сильно варьировать в пространстве (Wetterlind et al., 2010; Nyéki et al., 2022). Таким образом, большая трудоемкость полевых работ и достаточно сильная изменчивость параметров плодородия почв приводят к тому, что получить оперативно подобные данные для больших территорий традиционными методами практически невозможно и очень затратно.

Однако на протяжении последних десятилетий ведутся разработки в области создания новых методов картографирования параметров плодородия почв, основанных на использовании данных дистанционного зондирования (Савин, Симакова, 2012). Эти методы могут быть использованы для картографирования и оперативного мониторинга параметров плодородия пахотных почв (Ben-Dor, Banin, 1995; Viscarra Rossel et al., 2006; Савин и др., 2019; Prudnikova, Savin, 2021). Но, до сих пор они не нашли широкого применения в практике агрохимического обследования пахотных земель, что связано с целым рядом существующих ограничений и их недостаточной изученностью. Целью данной статьи был сравнительный анализ качества регрессионного моделирования параметров плодородия пахотных почв тестового участка в Тверской области с использованием данных дистанционного зондирования, различающихся уровнем пространственного обобщения информации.

ОБЪЕКТЫ И МЕТОДЫ

Тестовое поле расположено на землях Всероссийского научно-исследовательского института мелиорированных земель (филиал ФИЦ “Почвенный институт им. В.В. Докучаева”, н. п. Эммаусс, Тверская область) (рис. 1).

Размеры тестового участка составляют 90 на 120 м (1.08 га). Территория участка характеризуется выположенным рельефом с перепадом высот не более 1–2 м. В почвенном покрове господствуют осушенные закрытым дренажем дерново-подзолистые оглеенные почвы, развитые на супесях и легких суглинках, подстилаемых флювиогляциальными песками, часто с большим количеством валунных и галечниковых включений. Поверхность поля в момент отбора образцов и съемки была проборонована перед севом пшеницы, растительность отсутствовала.

Рис. 1. Территория исследований (справа красными границами выделено тестовое поле).

Fig. 1. Test field (within red boundary on the right).

На тестовом поле 15.05.2024 г. было отобрано 48 смешанных образцов поверхностного горизонта пахотных почв (0–5 см) (рис. 2).

В образцах определялись следующие показатели плодородия пахотных почв: рН в водной и солевой вытяжке, содержание органического вещества, массовая доля соединений фосфора, содержание обменного калия, обменных катионов по Шолленберге-ру (кальций, магний), массовая доля азота нитратов. Анализы образцов проводились в ИЛЦ ФИЦ “Почвенный институт им. В.В. Докучаева”.

В точках отбора образцов проводилось полевое спектромет-рирование поверхности почв. Для съемки спектральной отражательной способности (СОС) использовался полевой спектрора-диометр FieldSpec® HandHeld 2™ (“ASD Inc”, США), позволяю- щий получать спектры отражения в диапазоне длин волн от 300 до 1025 нм с интервалом 2 нм. Съемка проводилась в малооблачную погоду в дневное время. В каждой точке опробования измерения проводились в 5-кратной повторности. Для съемки спектров открытой поверхности почвы прибор располагался на высоте 90– 100 см над уровнем почвы. Полученные спектры осреднялись для каждой точки.

Рис. 2. Схема отбора образцов почв на тестовом поле.

Fig. 2. Soil sampling scheme.

Одновременно поверхность поля была отснята с использованием БПЛА DJI Matrice 200 с установленной на гиростабилизи-руемой подвеске камерой. Съемку осуществляли камерой Zenmuse X4S в видимом спектральном диапазоне (RGB). Основные характеристики камеры: матрица CMOS, 1"; число эффективных пикселей – 20 Мп; объектив F/2.8-11, 8.8 мм, FOV 84°. Для предполетной настройки коптера и камеры использовали приложение DJI GO 4, для автоматизированного построения маршрута полета и осуществления аэрофотосъемки – приложение DroneDeploy.

Обработку данных аэрофотосъемки проводили с использованием специального программного обеспечения Agisoft Metashape Professional

.

Полученное цветное изображение с пространственным разрешением 5 см на местности было разложено в ГИС ILWIS v.3.3 на три канала: R, G и B. После этого вокруг каждой точки отбора образцов (их положение было определено по характеру изображения, полученного с БПЛА) была построена буферная зона радиусом около 3 м. Для каждой зоны были рассчитаны средние и преобладающие значения тона изображения в каналах R, G и B.

Для анализа также использовались спутниковые данные, полученные 18 мая 2024 г. Атмосферно скорректированные спутниковые данные Sentinel-2 (продукты Level-2A) с пространственным разрешением 10 м на местности были получены с сайта: В исследование были включены 4-й, 3-й и 2-й каналы, приблизительно соответствующие каналам R, G и B используемого БПЛА.

Предобработка данных о спектральной отражательной способности почв заключалась в усреднении, удалении зашумленных участков и сглаживании. Для каждой точки на основе 5 повторностей рассчитывалась средняя спектральная кривая. При этом из-за зашумленности для анализа использовался только участок 350– 900 нм, который захватывает видимый и ближний инфракрасный участки спектра. Сглаживание спектральных кривых осуществлялось фильтром Савицкого-Голея в программе R с помощью функции savitzkyGolay пакета prospectr .

Далее для всех собранных дистанционных данных были рассчитаны спектральные индексы и спектральные отношения (табл. 1): 1) для спектрально генерализованных данных полевого спектрорадиометра; 2) для данных БПЛА на основе среднего и преобладающего значений отражения; 3) для спутниковых данных Sentinel-2.

Таблица 1 . Использованные спектральные индексы

Table 1 . Spectral indices used for the analysis

Индекс	Уравнение для данных спек-трорадиометра и БПЛА	Уравнение с учетом каналов Sentinel-2	Ссылка на уравнение по данным Sentinel-2
BI	\ (B 2 + G 2 +R 2 )	j (Band2 2 + Band32 + B and 4 2 )	Mathieu et al., 1998
	\ 3	J             3
CI	(R-G) (R + G)	(Band4 - Band3~) (Band4 + Band3~)	Mathieu et al., 1998
HI	(2*R-G-B)	(2 * Band4 - Band3 - Band2)	Mathieu et al., 1998
	(G-B)	(Band3 - Band2)
SI	(R-B) (R + B)	(Band4 - Band2) (Band4+Band2)	Mathieu et al., 1998
RI	R2 (B*G 3 >	Band42 (Band2 * Band3 3 )	Mathieu et al., 1998
GSI	(R-B)	(Band4 - Band2)	Xiao et al., 2006
	(R + B + G)	(Band4+Band3+Band2)
CRUST	(R-B) (R + B)	(Band4 - Band2~) (Band4 + Band2~)	Karnieli, 1997

Примечание. B – спектральная отражательная способность в синем участке спектра; G – спектральная отражательная способность в зеленом участке спектра; R – спектральная отражательная способность в красном участке спектра.

Note. B – the spectral reflectance in the blue region of the spectrum; G – the spectral reflectance in the green region of the spectrum; R – the spectral reflectance in the red region of the spectrum.

Моделирование зависимостей анализируемых свойств осуществлялось методом множественной линейной регрессии с помощью функции lm пакета stats программы R и проводилось отдельно для данных, полученных с помощью полевого спектро-радиометра, БПЛА и Sentinel-2 по следующей схеме.

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

Отобранные таким образом информативные параметры затем использовались для моделирования параметров плодородия пахотных почв.

Для проверки полученных моделей проводилась процедура кросс-валидации. Выборка случайным образом делилась на 10 блоков, из которых один использовался для проверки, а остальные для обучения. Для устойчивости значений параметров оценки качества подобное деление повторялось 100 раз. Для кроссвалидации использовался пакет caret программы R .

Основными параметрами для оценки качества моделей по результатам кросс-валидации являлись следующие: коэффициент детерминации (R2cv), среднеквадратическая ошибка предсказания (RMSEcv) и отношение межквартильного размаха (IQR) к среднеквадратической ошибке предсказания (RPIQ).

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

В отобранных образцах pH солевой вытяжки варьирует от 4.20 до 5.06, а рН водной вытяжки – от 4.70 до 5.75 (табл. 2). Для данных свойств отмечается наименьший коэффициент вариации: 13.69 и 9.87% соответственно. Сильнее всего варьирует содержание обменных катионов магния – от 0.46 до 13.98 мг-экв/100 г с коэффициентом вариации 81.48%, а также содержание подвижного фосфора – от 24.10 до 760.20 мг/кг с коэффициентом вариации 79.75%.

Таблица 2. Параметры плодородия пахотных почв тестового поля

Table 2. Arable soil fertility parameters

Параметры плодородия	Минимум	Среднее	Максимум	Стандартное отклонение	Коэффициент вариации
pH водной вытяжки	4.70	5.75	7.10	0.57	9.87
pH солевой вытяжки	4.20	5.06	6.70	0.69	13.69
Подвижный фосфор, мг/кг	24.10	146.35	760.20	116.71	79.75
Обменный калий, мг/кг	87.58	228.82	376.94	64.36	28.13
Массовая доля азота нитратов, млн-1	7.80	11.72	15.50	1.97	16.78
Органическое вещество, %	1.08	2.31	3.34	0.52	22.31
Обменные катионы кальция, мг-экв/100 г	3.13	11.04	41.57	5.85	52.95
Обменные катионы магния, мг-экв/100 г	0.46	3.61	13.98	2.94	81.48

Содержание обменных катионов кальция также колеблется достаточно сильно: от 3.13 до 41.57 мг-экв/100 г с коэффициентом вариации 52.95%. Содержание обменного калия варьирует в пределах от 87.58 до 376.94 мг/кг с коэффициентом вариации 28.13%, массовая доля азота нитратов – от 7.80 до 11.72 млн-1 с коэффициентом вариации 16.78%. Содержание органического вещества изменяется в пределах от 1.08 до 3.34% с коэффициентом вариации 22.31%.

В целом, вне зависимости от используемых дистанционных данных, модели для тестового поля получились не очень высокого качества (табл. 3). R2cv варьировал от 0.24 до 0.40, RPIQ – от 0.86 до 1.62.

Что касается содержания органического вещества, статистически достоверная модель была получена только при использовании данных полевого спектрорадиометра. R2cv модели равен 0.28, RPIQ – 1.14. Такая же ситуация отмечается для рН солевой вытяжки и обменных катионов магния. Для рН солевой вытяжки R2cv модели, полученной по данным полевого спектрорадиометра, составил 0.36, RPIQ – 1.62. Для обменных катионов магния R2cv модели, полученной по данным полевого спектрорадиометра, составил 0.36, RPIQ – 0.86.

В случае содержания обменного калия статистически достоверные модели были получены только для данных БПЛА. При этом при использовании преобладающих значений отражения предсказательная способность модели оказалась чуть выше, чем при использовании средних значений. R2cv модели составил 0.36, RPIQ – 1.52.

Для содержания подвижного фосфора статистически достоверные модели удалось получить по данным полевого спектрора-диометра и БПЛА. Лучшая модель по величине R2cv и RPIQ была получена при использовании преобладающих значений отражения, рассчитанных по данным БПЛА (R2cv = 0.40, RPIQ = 1.33).

Что касается рН водной вытяжки, лучшая модель по величине R2cv была получена по данным полевого спектрорадиометра (R2cv = 0.37), а по величине RPIQ – по спутниковым данным Sentinel-2 (RPIQ = 1.32).

Таблица 3. Моделирования параметров плодородия пахотных почв по дистанционным данным для тестового поля

Table 3. Models for arable soil parameters detection using remote sensing data

Свойства	Информативные параметры	R2cv	RMSEcv	RPIQ
Данные, полученные с помощью полевого спектрорадиометра (генерализованные)
Органическое вещество	R, G, B, B/R, SI, BI	0.28	0.51	1.14
Обменный калий	не удалось получить статистически достоверной модели
Подвижный фосфор	G, CI, BI, GSI	0.38	101.23	1.19
рН солевой вытяжки	B/G, B/R, CI	0.36	0.65	1.62
рН водной вытяжки	B/G, B/R, HI	0.37	0.58	1.21
Обменные катионы кальция	не удалось получить статистически достоверной модели
Обменные катионы магния	CI	0.36	0.21	0.86
Массовая доля азота нитратов	не удалось получить статистически достоверной модели

Продолжение таблицы 3

Table 3 continued

Свойства	Информативные параметры	R2cv	RMSEcv	RPIQ
Данные, полученные с помощью БПЛА (среднее значение отражения)
Органическое вещество	не удалось получить статистически достоверной модели
Обменный калий	SI	0.30	60.3	1.42
Подвижный фосфор	BI, RI	0.39	95.99	1.26
рН солевой вытяжки	не удалось получить статистически достоверной модели
рН водной вытяжки	SI	0.33	0.53	1.32
Обменные катионы кальция	не удалось получить статистически достоверной модели
Обменные катионы магния	не удалось получить статистически достоверной модели
Массовая доля азота нитратов	не удалось получить статистически достоверной модели

Продолжение таблицы 3

Table 3 continued

Свойства	Информативные параметры	R2cv	RMSEcv	RPIQ
Данные, полученные с помощью БПЛА (преобладающее значение отражения)
Органическое вещество	не удалось получить статистически достоверной модели
Обменный калий	B, G, B/G, B/R, CI	0.36	56.26	1.52
Подвижный фосфор	G, RI	0.40	90.98	1.33
рН солевой вытяжки	не удалось получить статистически достоверной модели
рН водной вытяжки	R	0.26	0.54	1.30
Обменные катионы кальция	не удалось получить статистически достоверной модели
Обменные катионы магния	не удалось получить статистически достоверной модели
Массовая доля азота нитратов	B/G, SI	0.31	1.83	1.48

Продолжение таблицы 3

Table 3 continued

Свойства	Информативные параметры	R2cv	RMSEcv	RPIQ
Данные, полученные с помощью Sentinel-2
Органическое вещество	не удалось получить статистически достоверной модели
Обменный калий	не удалось получить статистически достоверной модели
Подвижный фосфор	не удалось получить статистически достоверной модели
рН солевой вытяжки	не удалось получить статистически достоверной модели
рН водной вытяжки	BI	0.34	0.53	1.32
Обменные катионы кальция	не удалось получить статистически достоверной модели
Обменные катионы магния	не удалось получить статистически достоверной модели
Массовая доля азота нитратов	B/R, SI	0.27	1.94	1.39

Примечание. R, G, B – данные в соответствующих каналах БПЛА (соответствуют 4, 3 и 2 каналу Sentinel-2). В случае данных полевого спектрорадиометра – спектрально генерализованные данные в соответствующих каналах БПЛА. B/G, B/R – спектральные отношения отражения в соответствующих каналах. Формулы для расчета остальных параметров представлены в таблице 1.

Note. R, G, B – data in the corresponding UAV bands (corresponding to Sentinel-2 bands 4, 3, and 2). For field spectroradiometer data – spectrally generalized data in the corresponding UAV bands. B/G, B/R – spectral reflectance ratios in the corresponding bands. Formulas for calculating the remaining parameters are presented in Table 1.

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

Для массовой доли азота нитратов лучшая модель была получена по данным БПЛА на основе преобладающего значения отражения. R2cv модели составил 0.31, RPIQ – 1.48. По данным БПЛА на основе среднего значения отражения и по данным полевого спектрометрирования для данного свойства не удалось получить статистически достоверной модели.

На рисунке 3 представлены карты отдельных свойств почвы, построенные на основе полученных регрессионных зависимостей (БПЛА и спутниковые данные) и на основе пространственной интерполяции данных измерений в точках методом скользящего среднего (данные спектрометрирования).

Судя по данным рисунка 3, данные полевого спектрометри-рования позволили построить больше статистически значимых моделей с отдельными свойствами почв. Но на точность полученных на основе этих моделей карт, кроме качества самих моделей, сильное влияние оказывает используемый метод пространственной интерполяции точечных данных. Кроме того, практическое использование данного подхода сильно ограничивается его большой трудоемкостью. Для построения хороших моделей необходимы сотни точек с полевыми обследованиями (Докучаев и др., 2017).

Точность моделей, полученных по данным БПЛА, оказалась достаточно близкой к точности моделей, полученных по данным полевого спектрорадиометра. Но количество свойств, для которых удалось получить статистически достоверные модели сократилось. Поскольку съемка БПЛА выполнялась одновременно с полевым спектрометрированием, это позволило избежать влияния таких факторов, как состояние открытой поверхности пахотных почв (влажность, шероховатость, наличие корки (Persson, 2005; Прудникова и др., 2021)), а также разницы в условиях освещения (Xu et al. 2019). Так (Heil et al., 2022) отмечают, что данные факторы обуславливают разницу в качестве моделей, получаемых на основе данных БПЛА и спектральных данных, полученных в лабораторных условиях.

Дистанционные данные	Параметр плодородия
	Органическое вещество	Обменный калий	Подвижный фосфор	pH солевой вытяжки	pH водной вытяжки	Обменные катионы магния	Массовая доля азота нитратов
Полевое спектрометр иро-ванне	«	модель отсутствует	1	Я	1	Я	модель отсутствует
Данные с БПЛА	модель отсутствует		*	модель отсутствует	•	модель отсутствует
Спутниковые данные Sentmel-2	модель отсутствует	модель отсутствует	модель отсутствует	модель отсутствует	-л	модель отсутствует	Й

Рис. 3. Карты параметров плодородия почв тестового поля, построенные на основе дистанционных данных.

Fig. 3. Maps of soil fertility parameters of the test field, modelled on the basis of remote sensing data.

В то же время избыточная детальность данных, получаемых с БПЛА, вероятно, негативно сказалась на качестве получаемых результатов. Влияние на изображение полос вспашки, а также растительных остатков, крупных камней на поверхности почвы приводит к падению точности получаемых карт. Именно поэтому в случае использования данных БПЛА чуть лучше оказались модели, полученные на основе преобладающих значений тона изображения в отдельных каналах, а не среднего тона.

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

Наименьшее количество моделей было построено с использованием спутниковых данных, что связано с их генерализованно-стью и сильным пространственным варьированием свойств почв на тестовом поле (Лозбенев и др., 2024). Конечно же, не может не вносить свой негативный вклад и влияние атмосферных условий съемки и качества предобработки спутниковых данных.

Разница в состоянии открытой поверхности пахотных почв во время полевых работ и во время получения спутниковых данных также может повлиять на результаты моделирования, как это было показано ранее (Prudnikova, Savin, 2021). Однако поскольку спутниковые данные были получены достаточно близко к дате полевых работ, вероятнее всего, различия были не очень существенными.

Отмечаемое невысокое качество регрессионных моделей может быть связано с ограниченным набором каналов, используемых при анализе, охватывающих только видимый диапазон спектра. Использование на БПЛА камеры, работающей также в ближнем инфракрасном или тепловом участке спектра, может позволить увеличить информативность получаемых данных (Грубина и др., 2020).

ВЫВОДЫ

В ходе проведенных исследований было установлено, что наиболее точные модели для предсказания агрохимических свойств почв могут быть получены на основе полевого спектро- метрирования, менее качественные модели получены с использованием данных БПЛА и наименее качественные модели – на основе спутниковых данных. Основной причиной этого результата, по-видимому, является пространственное варьирование агрохимических показателей почв и генерализованность изображения их поверхности на данных БПЛА и спутниковых снимках.