Применение функций к слоям данных

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

В качестве примера создадим слой, состоящий из максимальных значений пикселей исходного растра:

img <- rast("KANOPUS_20190829.tif")
# class       : SpatRaster 
# dimensions  : 1440, 2052, 4  (nrow, ncol, nlyr)
# names       : blue, green, red, nir

img_max <- max(img)
# class       : SpatRaster 
# dimensions  : 1440, 2052, 1  (nrow, ncol, nlyr)
# name        : max

(показаны только интересующие нас метаданные).

Выберем пиксели с одинаковым номером (ipix) и посмотрим, какие значения они содержат в разных слоях исходного растра и что получилось после применения функции max:

ipix <- 100000

img_val <- values(img)
img_val[ipix,]
# blue green   red   nir 
#    6    19    30    59
    
values(img_max)[ipix,]
# max 
#  59

Можно применять и другие подобные функции: mean, min, sum. Всякий раз при этом создается новый однослойный растр.

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

mean_img <- global(img, "mean")
#            mean
# blue   74.30830
# green  76.57495
# red    69.35716
# nir   102.35891

apply(mean_img, 2, mean)
# mean 
# 80.64983

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

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

sar <- rast("louisiana_SAR.tif")

# гистограмма
#   выборки значений пикселей
hist(sar, main = "Backscattering")
#   всех значений пикселей
hist(values(sar), main = "Backscattering")

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

В нашем случае все просто. На снимке изображена река Миссисипи в штате Луизиана

plot(sar, col = grey.colors(15))

Нижняя часть распределения, с пиком -22.5 дБ, соответствует водным объектам, а верхняя, с пиком -13.5 дБ, представляет сушу.

Разделяя гистограмму по минимуму между пиками (примерно: -20 дБ), получим маску водных объектов

water_mask <- sar < -20
plot(water_mask, col = rev(grey.colors(2)))

Конечно, это еще не готовая маска. На ней много ложных срабатываний. Но самый важный шаг в построении маски уже сделан.

Для выявления выбросов кроме hist можно использовать boxplot. А вот пример послойного центрирования значений растровых данных

img_centered <- scale(img, center = TRUE, scale = FALSE)

Здесь от значения каждого пикселя отнимается среднее значение пикселей соответствующего слоя.

Изучить взаимосвязь между спектральными каналами снимка позволяет матрица диаграмм рассеяния

pairs(img) 

Она считается заметно дольше всех рассмотренных выше функций.

Матрица диаграмм рассеяния показывает тесную попарную связь между каналами видимого света (коэффициент корреляции Пирсона превышает 0,9), и гораздо более слабую связь этих каналов с 4-м каналом — каналом ближнего ИК-излучения (NIR).

#R

Доступ к ячейкам данных

Для доступа к значениям ячеек (пикселей) растровых данных в пакете terra есть несколько функций.

values() используется, чтобы получить все значения растра или их часть. Например, значения одной строки или значения прямоугольного блока ячеек:

r <- rast("dem.tif")
v <- values(r)

v[1000:1005, ]
# [1] 125.9172 126.2032 126.4972 126.7908 127.0973 127.3987

values(r, row=100, nrow=1, col=35, ncol=6)
#       Layer_1
# [1,] 112.4701
# [2,] 112.5410
# [3,] 112.6132
# [4,] 112.6909
# [5,] 112.7424
# [6,] 112.7527

Обратите внимание, что значения растра представлены в виде вектора (одномерного массива). Для доступа к значениям ячеек используется стандартная индексация R.

Можно сначала узнать номер ячейки по номеру строки и столбца, при помощи cellFromRowCol, а затем получить значения выбранных ячеек:

cells <- cellFromRowCol(r, 100, 35:40)
cells
# [1] 288719 288720 288721 288722 288723 288724
r[cells]
#    Layer_1
# 1 112.4701
# 2 112.5410
# 3 112.6132
# 4 112.6909
# 5 112.7424
# 6 112.7527

Получим координаты выбранных ячеек

xy <- xyFromCell(r, cells)
#              x        y
# [1,] -77.84894 39.49693
# [2,] -77.84890 39.49693
# [3,] -77.84887 39.49693
# [4,] -77.84884 39.49693
# [5,] -77.84881 39.49693
# [6,] -77.84878 39.49693

Получить значения ячеек можно по их номерам (cells) или координатам (xy), с помощью функции extract:

extract(r, cells)
extract(r, xy)
#    Layer_1
# 1 112.4701
# 2 112.5410
# 3 112.6132
# 4 112.6909
# 5 112.7424
# 6 112.7527

Координаты и значения всех ячеек растра в виде таблицы:

dt <- as.data.frame(r, xy=TRUE)

Функцию values можно использовать для замены значений в объекте SpatRaster или для присвоения новых значений ячейкам. Если вы заменяете значение в объекте SpatRaster, созданном на основе файла, то связь с этим файлом теряется, поскольку теперь значения объекта отличаются от файла

r
# class       : SpatRaster 
# dimensions  : 1620, 2916, 1  (nrow, ncol, nlyr)
# resolution  : 3.08642e-05, 3.08642e-05  (x, y)
# extent      : -77.85, -77.76, 39.45, 39.5  (xmin, xmax, ymin, ymax)
# coord. ref. : lon/lat WGS 84 (EPSG:4326) 
# source      : dem.tif 
# name        :  Layer_1
# min value   :  87.12177 
# max value   : 158.5774

values(r)[1000:1005, ] <- 0
r
# class       : SpatRaster 
# dimensions  : 1620, 2916, 1  (nrow, ncol, nlyr)
# resolution  : 3.08642e-05, 3.08642e-05  (x, y)
# extent      : -77.85, -77.76, 39.45, 39.5  (xmin, xmax, ymin, ymax)
# coord. ref. : lon/lat WGS 84 (EPSG:4326) 
# source(s)   : memory
# name        :  Layer_1 
# min value   :   0.0000 
# max value   : 158.5774

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

#R

Из полученного растра исчезнут единицы. Если задать порог численности равным 9, то исчезнут все скопления, кроме троек (и, естественно, нулей).

#R

*Гонсалес Р., Вудс Р. Цифровая обработка изображений. 3-е издание, исправленное и дополненное. — М.: Техносфера, 2012. — 1104 с.

Отображение растровых данных

Основными функциями для отображения растровых данных в terra являются функция plot, которая отображает слои данных, и plotRGB, которая создает и отображает RGB-композитные изображения.

Функция plot отображает растровые данные по слоям. Основные ее возможности видны на примере знакомой нам функции plot_raster:

plot_raster <- function(r) {
  plot(r, axes = F, legend = F)
  plot(as.polygons(r, na.rm = F, dissolve = F, trunc = F), add = T)
  text(r, digits = 2)
}

Первый вызов plot отображает растр без осей координат (axes = F) и легенды (legend = F).

Второй plot добавляет поверх первого изображения (add = T) набор полигонов, созданный из ячеек исходного растра. Опция add = T очень удобна для добавления изображений векторных объектов, например границ, поверх растровых данных. Можно накладывать друг на друга и растровые данные, задавая их прозрачность параметром alpha, который изменяется от 0 (прозрачный) до 1 (непрозрачный).

Функция text отображает… текст поверх исходного растра r. В нашем примере она отображает значения соответствующих ячеек r.

При отображении пространственных данных можно использовать знакомые по базовой графике функции points и lines. В terra они служат для отображения точечных и линейных векторных объектов, а функция polys отображает полигоны.

plot_raster() рассчитана на работу с однослойными данными. Нам же, чаще всего, предстоит работать с данными многослойными. В следующем примере plot отображает 9 каналов снимка, полученного прибором OLI спутника Landsat 8:

lc <- rast("LC08_044034_20170614.tif")
# class       : SpatRaster 
# dimensions  : 1245, 1497, 9  (nrow, ncol, nlyr)
# resolution  : 30, 30  (x, y)
# extent      : 594090, 639000, 4190190, 4227540  (xmin, xmax, ymin, ymax)
# coord. ref. : WGS 84 / UTM zone 10N (EPSG:32610)

plot(lc)

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

plotRGB(x, r=1, g=2, b=3)

x — объект типа SpatRast. Следующие три аргумента — номера слоев данных (спектральных каналов), использующихся в качестве красного (r), зеленого (g) и синего (b) цветов итогового изображения. Аргументов у plotRGB намного больше, мы рассматриваем лишь самые основные.

При заданных по умолчанию номерах спектральных каналов, для снимка Landsat 8 мы получим изображение в комбинации 1-2-3, где в роли красного выступает Coastal/Aerosol (канал 1), в роли зеленого — синий (канал 2), а синим будет зеленый (канал 3). Не слишком удачный выбор. Комбинация “естественные цвета” 4-3-2 (красный-зеленый-синий) задается так

plotRGB(lc, 4, 3, 2)

В результате получим темный прямоугольник. Дело в том, что по умолчанию значения пикселей в каждом канале находятся в интервале от 0 до 255, тогда как на снимке эти значения лежат в интервале от 0 до 1. Исправить ситуацию может параметр scale — целое число, которое задает максимальное значение пикселя в каждом из трех каналов

plotRGB(lc, 4, 3, 2, scale = 1)

Второй вариант решения проблемы — параметр stretch. Он задает способ растяжения изображений для повышения контрастности: "lin" (линейный) или "hist" (гистограмма).

plotRGB(lc, 4, 3, 2, stretch = "lin")

Комбинацию 5-4-3 мы видели здесь. Для изучения сельскохозяйственных земель, водно-болотных угодий, а также при анализе пустынь может пригодиться комбинация 7-5-3 (SWIR, NIR, Green):

plotRGB(lc, 7, 5, 3, stretch = "lin")

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

#R

Сигнальная информация или автоматизация сбора аннотаций статей

В древности, каждые две недели в наш домашний почтовый ящик приходила сигнальная информация ВИНИТИ — ксерокопии оглавлений журналов по выбранной тематике. Сигнальная информация, увы, давно уже не приходит, но потребность в ней осталась. Попробуем решить эту проблему при помощи нескольких строк кода в R и библиотеки rcrossref.

Ключевой фрагмент кода:

start_date <- "2024-01-01"
end_date   <- "2024-01-30"

jname <- "Remote Sensing"

dt_list <- cr_works(filter = c(from_pub_date = start_date,
                               until_pub_date = end_date,
                               container_title = jname,
                               type = "journal-article"), 
                    limit = 1000)

cat("Найдено статей:", dt_list$meta$total_results)
dt.0 <- dt_list$data

cols <- c("title","url","abstract")
dt <- dt.0[, cols]

Самое главное делает функция cr_works. Она отбирает из базы данных CrossRef статьи (type = "journal-article"), опубликованные в журнале Remote Sensing (jname) в период от start_date до end_date.

На выходе cr_works возвращает список с метаданными поиска (dt_list$meta) и таблицей найденных статей (dt_list$data). В этой таблице несколько десятков колонок.

Дальше мы отбираем из таблицы dt.0 нужные колонки ("title","url","abstract") и работаем с ними. Оставшаяся часть кода посвящена отображению результатов поиска в виде таблицы с гиперссылками. Полный код и пример построенной таблицы прилагаются ⬇️.

В CrossRef попадают аннотации статей далеко не ото всех журналов. Почему так происходит — неясно. Однако название и DOI есть для каждой публикации. Вообще, CrossRef — неплохой инструмент для наукометрических исследований. Но это уже другая история.

Задержка с помещением сведений о публикации в базу CrossRef для в журнала Remote Sensing составляет около двух суток.

#R #справка

⭐️ СТРАНЫ / КОМПАНИИ / СПУТНИКИ

Страны: #австралия #германия #индия #иран #испания #канада #китай #португалия #россия #США #япония и т. п.
Но:
#корея обозначает Северную и Южную Кореи
#РБ — Республика Беларусь
#UK — Великобритания

Компании: #planet #maxar

Спутники: #landsat #sentinel1 #sentinel2

⭐️ ДЗЗ

Методы и приборы
#альтиметр
#гиперспектр — гиперспектральная оптическая съемка
#лидар
#оптика — мультиспектральная оптическая съемка
#радиометр — микроволновой радиометр
#dnb — ночная съёмка (day / night band)
#SIF — солнечно-индуцированная флуоресценция хлорофилла
#ro — радиозатменный метод
#SAR — радарная съемка
#InSAR — радарная интерферометрия
#LST — съемка в тепловом инфракрасном диапазоне
#GNSSR — ГНСС-рефлектометрия
#sigint — радиоэлектронная разведка

Виды орбит: #ГСО — геостационарная, #VLEO — сверхнизкая

#основы — обучающие материалы по ДЗЗ
#обучение курсы, обучающие сервисы и т. п.
#история — в основном, история ДЗЗ
#индексы — спектральные индексы
#комбинация — комбинации каналов

Данные
#данные — коллекции данных ДЗЗ, наземных данных, карты и т.п.
#датасет — набор данных для машинного обучения
Дополнительные хештеги, описывающие данные:
#LULC — Land Use & Land Cover
#осадки
#SST — Sea Surface Temperature
#nrt — (near real time) изображения, получаемые в режиме, близком к реальном времени
#debris — космический мусор
#границы — административные границы
#DEM — цифровая модель рельефа (ЦМР)
#keyhole — рассекреченные снимки разведспутников

Литература, справочная информация
#справка — спектральные каналы, орбиты спутников, поиск данных и т.п.
#обзор
#книга — текст книги прикреплён к сообщению.
Дополнительные хештеги:
#наблюдение — ресурсы для наблюдения спутников и орбиты спутников
#космодромы

#конференции — анонс конференций/семинаров/школ, посвященных ДЗЗ и анализ их материалов.
#конкурсы — анонс конкурсов/чемпионатов/олимпиад.
#МВК — материалы заседаний Межведомственной комиссии (МВК) по использованию результатов космической деятельности.

#снимки — поучительные (хоть в чем-то интересные) снимки, первые снимки

Программные инструменты / Языки
#нейронки #софт #GEE #R #tool #python
#ИИ #FM — Foundation Model (Remote Sensing Foundation Model)

⭐️ ОТРАСЛИ / ТЕМАТИЧЕСКИЕ ЗАДАЧИ

#археология #атмосфера #вода #война #засуха #климат #лед #лес #нефть #океан #оползни #наводнение #пожары #почва #растительность #севморпуть #сельхоз #снег
#AGB — надземная биомасса
#ЧС — мониторинг стихийных бедствий и катастроф
#GHG — парниковые газы
Отдельные газы: #CO2 #NO2
#энергетика — космическая энергетика
#SSA — Space Situational Awareness

И еще — к ⬆️

Пакет R cropCalendars — для моделирования календарей сельскохозяйственных культур в соответствии с подходами Waha et al. (2012) и Minoli et al. (2019).

#R

Гармонизированные наземные данные Land Use/Cover Area Frame Survey (LUKAS)

📖 d’Andrimont, R. et al. (2020). Harmonised LUCAS in-situ land cover and use database for field surveys from 2006 to 2018 in the European Union. Scientific Data, 7(1). https://doi.org/10.1038/s41597-020-00675-z

В Европейском союзе (ЕС) с 2006 года проводятся трехгодичные обследования почвенно-растительного покрова и землепользования в рамках программы Land Use/Cover Area Frame Survey (LUCAS). В ходе пяти обследований LUCAS было собрано 1351293 наблюдений в 651780 уникальных точках для 106 переменных, а также 5,4 млн фотографий. В работе все эти данные были объединены в единую гармонизированную базу данных, которая теперь является наиболее полным натурным набором данных о почвенно-растительном покрове и землепользовании в ЕС.

Для доступа к данным создан специальный R-пакет lukas.

🔗Страница данных Harmonised LUCAS in-situ land cover and use database for field surveys from 2006 to 2018 in the European Union, со ссылкой на скачивание
🔗 Данные на figshare

📸 Примеры данных LUCAS по одной точке, которую посетили пять раз в период 2006–2018 гг.

#данные #R

Применение функций к слоям данных: app и lapp

app применяет ко всем слоям растрового объекта SpatRaster заданную функцию fun, подобно тому как это делает apply для матрицы или таблицы.

r <- rast(ncols=3, nrows=3)
values(r) <- 1:ncell(r)
x <- c(r, sqrt(r), r+50)
s <- app(x, fun=sum)
s
# Для функций вроде
# "sum", "mean" и "max" можно сделать так:
# sum(x)

Можно посмотреть, как суммируются значения соответствующих пикселей всех слоёв x:

values(x)
values(s)

lapp (от layer-apply) позволяет объединить несколько слоёв растровых объектов. Его можно использовать как альтернативу растровой алгебре:

s <- rast(system.file("ex/logo.tif", package="terra")) + 1  
ss <- s[[2:1]]

fvi <- function(x, y) {(x - y ) / (x + y)} 
x <- lapp(ss, fun=fvi)
# Другой вариант:
y <- fvi(s[[2]], s[[1]])

Количество аргументов функции fun должно соответствовать количеству слоёв растра.

f2 <- function(x, y, z){ (z - y + 1) / (x + y + 1) } 
p1 <- lapp(s, fun=f2)
p2 <- lapp(s[[1:2]], f2, z=200)

#R

Применение функций к слоям данных: tapp и xapp

tapp применяет функцию fun к подмножествам слоёв SpatRaster, аналогично tapply и aggregate.

Слои объединяются на основе индекса. Например, для растра с шестью слоями можно использовать index=c(1,1,1,2,2,2). Он означает, что три первых слоя принадлежат одному подмножеству, а три следующих — другому. Если задана функция fun=sum, то на выходе получим SpatRaster с двумя слоями: суммой первых трёх слоев входного растра и суммой последних трёх слоев его же:

r <- rast(ncols=10, nrows=10)
values(r) <- 1:ncell(r)
s <- c(r, r, r, r, r, r)
s <- s * 1:6
b1 <- tapp(s, index=c(1,1,1,2,2,2), fun=sum)
b1
b2 <- tapp(s, c(1,2,3,1,2,3), fun=sum)
b2

Индексы могут повторяться, так что index=c(1,2) для растра с шестью слоями снова вернёт растр, состоящий их двух слоёв: один основан на нечетных слоях (1,3,5), другой — на чётных (2,4,6).

А вот пример поинтересней. В нём используется новый атрибут растровых данных — время, которое задаётся функцией time. С помощью tapp мы применяем заданные функции к слоям, сгруппированным по месяцам:

r <- rast(system.file("ex/elev.tif", package="terra"))
rst <- rep(r, 6) * 1:6
time(rst) <- as.POSIXct(c("2023-01-16","2023-01-17","2023-01-18","2023-02-17","2023-02-18","2023-02-19"),tz="UTC")

# Максимум по месяцам
x <- tapp(rst, "months", max) 
# 95-й процентиль по месяцам
y <- tapp(rst, "months", fun = function(i) quantile(i, probs = 0.95, na.rm = T))

Применяет функцию к ячейкам двух растров, в том числе многослойных.

Вычислим попиксельную корреляцию между двумя многослойными растрами:

r <- rast(ncols=3, nrows=3, nlyr=5)
set.seed(1)
r <- init(r, runif)
s <- init(r, runif)
x <- xapp(r, s, fun=cor)

#R

SpatRasterDataset и SpatRasterCollection

До сих пор мы работали с объектами класса SpatRaster.

SpatRaster представляет собой прямоугольную часть мира, разделённую на прямоугольные ячейки одинаковой площади (в единицах заданной системы координат). Для каждой ячейки может быть несколько значений (“слоёв”).

SpatRaster может указывать на один или несколько файлов на диске, в которых хранятся значения ячеек, и/или хранить эти значения в памяти. Эти объекты могут быть созданы с помощью метода rast.

Из объектов SpatRaster можно создавать новые объекты того же класса, а кроме них — SpatRasterDataset и SpatRasterCollection.

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

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

SpatRasterDataset создаётся функцией (методом) sds:

r <- rast(system.file("ex/logo.tif", package="terra"))   
x <- sds(r, r/2)
names(r) <- c("first", "second")
r

# Узнаем длину SpatRasterDataset
length(x)

# Извлечём 2-й SpatRaster
x[2]

SpatRasterCollection создаётся функцией sprc:

x <- rast(xmin=-110, xmax=-50, ymin=40, ymax=70, ncols=60, nrows=30)
y <- rast(xmin=-80, xmax=-20, ymax=60, ymin=30)
res(y) <- res(x)
values(x) <- 1:ncell(x)
values(y) <- 1:ncell(y)

z <- sprc(x, y)
z
z[1]

#R

Введение в геопространственный анализ и визуализацию в R

Журнал MDPI Data, наряду с данными, публикует статьи о методах обработки и анализа данных. В нём нам попалась забавная учебная статья

📖 Maesen, P., & Salingros, E. (2024). Introduction to Reproducible Geospatial Analysis and Figures in R: A Tutorial Article. Data, 9(4), 58. https://doi.org/10.3390/data9040058

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

В статье кратко представлены основные понятия растровых и векторных данных, системы координат, также описан базовый рабочий процесс для проведения воспроизводимых исследований в R. Приведены примеры создания различных типов карт (scatter, bubble, choropleth, hexbin и faceted) на основе открытых экологических данных. На этих примерах демонстрируются основные операции с геопространственными векторными данными (чтение, преобразование систем координат, создание геометрий, буферных зон вокруг существующих геометрий, определение пересечений между геометриями и т. п.).

Данные и скрипты есть в Supplementary Materials (скачать https://www.mdpi.com/article/10.3390/data9040058/s1).

📸 Карта средних значений органического углерода почвы (soil organic carbon, SOC) в Германии, измеренных в период 2011–2018 гг.

#R

Пакет mlhrsm для картографирования влажности почвы с высоким пространственным разрешением

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

В работе

📖 Peng, Y., Yang, Z., Zhang, Z., & Huang, J. (2024). A Machine Learning-Based High-Resolution Soil Moisture Mapping and Spatial–Temporal Analysis: The mlhrsm Package. Agronomy, 14(3), 421. https://doi.org/10.3390/agronomy14030421

предлагается очередная модель машинного обучения для картографирования влажности почвы. Она основана на алгоритме квантильного случайного леса (quantile random forest) и использует данные наземных датчиков влажности, параметры поверхности земли (растительность, рельеф и почву), а также оценки влажности почвы на поверхности и в прикорневой зоне, полученные по спутниковым данным. В работе используются данные спутников SMAP, Sentinel-1, Landsat, а также данные приборов MODIS. Область исследования: CONUS (contiguous USA), где существуют открытые данные наземных датчиков влажности почвы.

Модель позволяет создавать карты влажности почвы высокого разрешения (от 30 до 500 м, от ежедневных до ежемесячных) и строить оценки неопределенности на участках по территории CONUS на уровнях 0–5 см и 0–1 м.

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

Привлекает в работе то, что весь расчёт оформлен в виде пакета mlhrsm на языке R с открытым исходным кодом. По сути, статья — это руководство пользователя mlhrsm, где показан расчёт влажности почвы на примере одного поля. Её можно использовать для анализа сильных и слабых сторон подобных моделей, а также как основу для создания собственных моделей.

#почва #сельхоз #R

1️⃣ Скрипт для обработки ежемесячных данных.
2️⃣ Список космических и суборбитальных запусков в первом полугодии 2024 года.

#R