Управление ходом вычислений: условные операторы

Скалярный условный оператор if в R записывается так:

if (условие) {
  выполняется, если условие истинно
} else {
  выполняются, если условие ложно
}

Конструкции else-if в R нет.

Внимание! Если в коде записано:

if (condition == TRUE) x <- TRUE
else x <- FALSE

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

##    Error: unexpected 'else' in " else"

Чтобы R интерпретировал else как часть предшествующего if’а, нужно указать при помощи фигурных скобок, что if еще не закончился:

if (condition == TRUE) { 
  x <- TRUE
} else {x <- FALSE}

Есть в R и скалярный switch, но используется он редко.

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

ifelse(условие, выполнить_если_ИСТИНА, выполнить_если_ЛОЖЬ)

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

x <- c(6:-4)
sqrt(x)

Значения аргументов необходимо предварительно отфильтровать, заменив отрицательные числа какими-нибудь допустимыми значениями, например, NA:

# фильтр для данных
x1 <- ifelse(x>=0,x,NA)
sqrt(x1)

Однако мы снова получим предупреждение в ситуации

ifelse(x>0,sqrt(x),NA)

потому что, несмотря на условие, sqrt() выполняет операцию над исходным "нефильтрованным" вектором.

#R

Управление ходом вычислений: циклы

Циклы в R используются довольно редко. "Винить" за это следует векторизацию. Тем не менее, иногда именно циклы являются самым простым путем решения задачи.

Синтаксис циклов в R похож на синтаксис циклов в C-подобных языках программирования. Циклы могут вкладываться друг друга в соответствии с обычными правилами вложения циклов.

Цикл со счетчиком записывается следующим образом:

for (i in start:stop) {
  <тело цикла>
}

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

for (i in 1:5) {
  print(i^2)
}

Функция print используется для вывода промежуточных данных из тела цикла.

Можно перебирать в цикле элементы вектора:

x <- c("Антон", "Арсений", "Дмитрий", "Сергей")
for(i in x) {
  print(paste("Имя", i, "состоит из", nchar(i), "символов."))
}

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

Тело цикла выполняется до тех пор, пока истинно условие:

while (условие) {
 <тело цикла>
}

Пример:

x <- 5
while (x < 10) {
 x <- x + 1
 print(x)
}

Для преждевременного выхода из цикла используется команда break. После нее управление передается во внешний цикл, команде, следующей за циклом из которого мы вышли. Команда next выполняет переход к следующей итерации цикла.

В следующем примере третий столбец матрицы A так и не будет заполнен:

A <- matrix(NA, nrow = 3, ncol = 3)
for (i in 1:3) {
  for (j in 1:3) {
    if (j > 2) break
      A[i,j] <- i + j
  }
}

#R

Спутник “СМОТР-Р”

Михаил Котов на 12-м Петербургском Международном Газовом Форуме сфотографировал несколько плакатов, посвященных спутнику “СМОТР-Р”.

“СМОТР-Р” — это радарный спутник, который планирует производить “Газпром СПКА”. Аппарат массой 650 кг, оснащенный радаром X-диапазона с активной фазированной антенной решеткой, должен быть изготовлен в 2025 году. Срок активного существования на орбите: 7 лет.

Параметры съемки (разрешение/ширина полосы захвата):

* детальный режим — 2 м / 11 км
* маршрутный режим — 5 м / 25 км
* обзорный режим — 15 м / 115 км

В перспективе планируется создание орбитальной спутниковой группировки из шести аппаратов “СМОТР-Р”. В таком случае период повторной съемки всего региона Северного морского пути составит 24 часа.

“Газпром СПКА” — дочерняя компания “Газпрома”, учрежденная в августе 2015 года для реализации проекта создания сборочного производства космических аппаратов (СПКА) в интересах “Газпрома” и других потенциальных заказчиков.

Ранее, производство спутников для системы дистанционного зондирования Земли "СМОТР" планировали начать на предприятии в подмосковном Щелкове уже в нынешнем году.

Предполагается, что в систему СМОТР войдут также спутники высокодетальной оптической съемки "СМОТР-В". Помимо сенсоров для мультиспектральной и панхроматической съемки, "СМОТР-В" оснащен газоанализатором, предназначенным для выявления источников выбросов метана и других парниковых газов.

📸Фото

#россия

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

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

Базовая графика в R

Рассмотрим основы работы с базовой двумерной графикой R. Базовая она, потому что реализована в пакете graphics, который поставляется вместе с R и подключается автоматически. В будущем мы познакомимся с другими графическими пакетами R, в первую очередь, со знаменитым ggplot2. Но сейчас — про основы, ибо с функцией plot нам предстоит встречаться очень часто.

Построим график синусоиды, координаты которой хранятся в векторах x и y:

x <- seq(-pi,pi,.1)
y <- sin(x)

График можно составить из отдельных элементов:

plot.new() # пустое окно
plot.window(xlim = c(-pi,pi), ylim = c(-1,1))
points(x,y)
axis(1)
axis(2)
box()
title(xlab = "x")
title(ylab = "y")

Здесь:

* plot.new — создает пустое графическое окно;
* plot.window — задает пределы изменения x- и y-координат на графике, а также другие графические параметры;
* points — строит график в виде точек;
* axis — добавляет на график оси координат: 1 – снизу, 2 – слева, 3 – сверху, 4 – справа;
* box — строит рамку вокруг графика;
* title — выводит пояснения к графикам. Параметры xlab, ylab задают обозначения осей координат.

Тот же график можно построить при помощи всего лишь одной функции:

plot(x,y)

Раньше мы имели дело с низкоуровневыми функциями, которые создают график из отдельных элементов, как из кубиков. plot — функция высокого уровня. Она создает нужный график целиком.

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

plot(x, y, main="2d-Graph", xlab="Axis X", ylab="Axis Y")

По умолчанию, точки данных изображаются в виде незакрашенных кружков черного цвета. Символы для отображения точек (маркеры) настраиваются параметром pch (plot character). Цвет задается параметром col. Цвета по названиям цветов на английском языке ("red", "green" и т. п.), по номерам цветов, а также при помощи шестнадцатиричных строк #RRGGBB, как HTML.

Пусть точки данных отображаются в виде красных "звездочек":

plot(x, y, pch="*", col="red")

Маркеры можно задавать по номерам:

plot.new()
plot.window(xlim=c(0,2.1), ylim=c(.8,1.1), yaxs="i")
n <- 20
x <- seq(0.1,2,length=n)
points(x,rep(1,length(x)),pch=1:n,cex=2)
text(x,rep(.9,length(x)),1:n)

Тип отображения кривых определяется параметром type: "p" — точки (по умолчанию), "l" — линия, "b" — точки и линия (both).

x <- seq(-pi,pi,.3)
y <- sin(x)
plot(x, y, col="#009900", type="b")

Для отображения кривых существуют специальные функции: points, lines и т. п. Они добавляют кривые в уже существующий график:

lines(x+.3, y, col="#0000FF")

Наиболее употребительные графические параметры:

* main — заголовок графика
* xlab, ylab — подписи к осям x и y
* xlim, ylim — диапазон изменения данных по осям x и y
* pch — маркер, которым отображается точка данных
* cex — размер маркера
* col — цвет
* lty — тип линии (solid, dashed, dotted...)
* lwd — толщина линии

Для глобальной настройки параметров графиков используется функция par. Вызов par без аргументов возвращает текущие значения графических параметров. Указав графический параметр в качестве аргумента par, получим его значение:

par("col")

По умолчанию, для рисования используется черный цвет.

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

x <- seq(-pi,pi,.3)
y <- sin(x)
par(col="red")
plot(x,y,type="b")
lines(x+.3, y)
points(x-.3, y, col="blue")

Функция par устанавливает единый для всех графиков цвет вывода кривых и осей координат (красный). Локальное значение параметра (синий цвет последнего графика) перекрывает действие глобального.

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

# узнать текущие настройки
par()
# сделать их копию
old_par <- par()
# создать новые настройки
par(...)
# построить необходимые графики
<plot>
# восстановить прежние настройки       
par(old_par)

Функция dev.off() закрывает текущее графическое окно, а graphics.off() закрывает все графические окна в текущем сеансе.

#R

Глобальные карты высоты леса

В (Potapov et al., 2020) описан метод создания глобальной карты высоты лесного полога с пространственным разрешением 30 м, основанной на данных космического лидара GEDI и разновременных данных Landsat.

Для создания карты использованы данные GEDI (апрель–октябрь 2019 года) и данные Landsat за 2019 год. Высота леса моделировалась ансамблем деревьев регрессии с движущимся окном. Калибровка модели осуществлялась локально, с применением метрики GEDI RH95 (относительная высота на уровне 95%). В качестве признаков взяты разновременные данные Landsat, призванные отразить фенологические изменения на поверхности. Поскольку лидар GEDI, работающий на борту МКС, позволяет проводить измерения лишь в полосе широт от 51.6° с.ш. до 51.6° ю.ш., то для создания глобальной карты пришлось экстраполировать построенную модель в бореальные регионы (за пределы диапазона данных GEDI).

🌍 Карту можно получить на GEE: Global Forest Canopy Height from GEDI & Landsat.

Карта подготовлена лабораторией Global Land Analysis and Discovery (GLAD) департамента Географических наук университета штата Мэриленд. GLAD известна своими картами Global Forest Change и данными Harmonized Landsat Sentinel-2.

Следующая глобальная карта высоты лесного полога построена по данным 2020 года с пространственным разрешением 10 м (Lang et al., 2023). Для ее создания использованы данные GEDI, спутниковые снимки Sentinel-2 и ансамбль моделей на основе сверточных нейросетей (CNN). Итоговая модель позволяет получить в любой точке Земли высоту лесного полога и погрешность оценки высоты полога.

Применявшиеся ранее подходы для создания карт высоты полога на основе данных GEDI и Landsat-8, в частности (Potapov et al., 2020) реализовывали попиксельное отображение пространства признаков в высоту полога, без учета свойств окрестности пикселя и текстуры изображения. Применение сверточных нейросетей позволило учесть эти локальные особенности, что является существенной новизной предложенного подхода.

Проблема насыщения при расчете высоты полога, которая приводит к ухудшению точности оценок с увеличением высоты полога, в (Lang et al., 2023) решалась чисто механически — за счет огромного объема обучающих данных. Точность оценки высоты по-прежнему падает с увеличением высоты деревьев, хотя и несколько лучше, чем у более ранних подходов.

Важно, что предложенный (Lang et al., 2023) подход к расчету высоты полога может переноситься на другой временной период. Напомним, что лидар GEDI будет работать на орбите по крайней мере до конца 2024 года.

Исходный код и обученные модели доступны на Github: github.com/langnico/global-canopy-height-model.

🌍 Карта на GEE: ETH Global Sentinel-2 10m Canopy Height (2020).

Код примера

Представляется, что данная карта высот лесного полога является лучшей на сегодняшний день среди глобальных карт с точки зрения пространственного разрешения и точности. В то же время, она не использует лидарные данные для оценки высоты за пределами полосы измерений GEDI, что открывает возможности для ее совершенствования, например, с привлечением данных лидара ICESat-2.


P. Potapov, X. Li, A. Hernandez-Serna, A. Tyukavina, M.C. Hansen, A. Kommareddy, A. Pickens, S. Turubanova, H. Tang, C. E. Silva, J. Armston, R. Dubayah, J. B. Blair, M. Hofton (2020). https://doi.org/10.1016/j.rse.2020.112165

Lang, N., Jetz, W., Schindler, K. _et al._ A high-resolution canopy height model of the Earth. _Nat Ecol Evol_ (2023). https://doi.org/10.1038/s41559-023-02206-6

#лидар #лес #данные #GEE

Функции

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

Синтаксис функции имеет вид:

имя <- function(список_аргументов) {
  тело функции
}

Например:

sqr <- function(x) {
  x^2
}

или

sqr <- function(x) x^2

Вызов функции по имени:

sqr(2)

Обратите внимание, что присваивание имени функции выглядит так же, как присваивание значения переменной. И действительно, переменная sqr относится к классу функций ("function"):

class(sqr)

Имя функции без круглых скобок возвращает код функции.

Можно создавать анонимные функции:

(function(x) x^2)(3)

Список аргументов функции может быть пуст, тогда после имени функции указываются лишь круглые скобки ().

Можно задавать значения аргументов по умолчанию в виде x = a, y = b и т. д. Обратите внимание: для присваивания значений аргументов используется только знак равенства, стрелка <- тут недопустима.

Список аргументов переменной длины формируется добавлением многоточия ... в конце списка. Так открывается широкий простор для создания функций-"оберток".

Пусть нам нужна функция, строящая графики линиями красного цвета. Создадим функцию-обертку вокруг plot, то есть функцию, основу которой составляет plot, но с некоторыми "доработками" — в данном случае, с фиксированным значением параметра цвета (col):

plot.red <- function(x, y, ...) {
  plot(x, y, col = "red", ...)
}

Вот как это работает:

x <- seq(-pi,pi,.1)
y <- sin(x)
plot.red(x,y)

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

Функции в R могут возвращать несколько объектов, относящихся к разным типам данных. Такие объекты нужно предварительно “завернуть” в список:

sqr <- function(x) {
  val <- x^2
  txt <- "Squaring"
  ret <- list(value=val,text=txt)
  return(ret)
}
sqr(2)

Созданные внутри функции переменные являются локальными, то есть их значения доступны (видимы) только внутри данной функции или вложенных в нее функций. Соответственно и функция "видит" переменные, заданные вне ее пределов:

a <- 5

foo <- function() {
  b <- 10  # локальная переменная foo
  a <- 3
  print(a) # foo видит a
  print(b) # foo видит b
}

foo()
# a в основной программе не изменилась
print(a)
# b в основной программе не определена
print(b)

Хотя a задана вне функции и не передается в нее через входные параметры, функция все равно видит эту переменную 1️⃣. Изменения внешних переменных внутри функции не сказываются на значении этих переменных по выходе из функции.

Области видимости вкладываются друг в друга 2️⃣.

Глобальные переменные внутри функции можно задать при помощи специального оператора присваивания <<-:

foo <- function() {
   b <<- 10
}
foo()
# теперь b видима из основной программы
b

Локальность и глобальность считается по отношению к конкретному окружению. Так, в следующем примере функция g() является глобальной по отношению к h():

g <- function(x) {
  h <- function(x) {x-1}
  print(environment(h))
  print(h(x))
  return(x+1)
}
g(2)

#R

Функциональное программирование

По сути, R является языком функционального программирования. Практически все в нем делается при помощи функций. Например, элементы векторов выбираются при помощи функции-квадратной скобки:

x <- c(1,2,3)
x[1] 

Ничего необычного? Тогда запишем функцию-квадратную скобку в “более функциональном" виде:

`[`(x,1)

За следующей записью:

`+`(1,`*`(2,3))

кроется обычное

1+2*3

а операцию присваивания можно выполнить так:

assign("a",1)

В качестве примера создадим оператор +=, которого нет в базовом R:

`%+=%` <- function(o1,o2) eval.parent(substitute(o1 <- o1 + o2))

x <- 1
x %+=% 2

Внутри новый оператор остается функцией.

Возможности функций в R огромны, мы лишь намекнули на них. Подробности можно найти во второй части книги Хэдли Уикэма: Wickam H. Advanced R. — CRC Press, 2015.

#R

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