project logo

Начало работы с R

R — это свободное программное обеспечение для статистических вычислений и графики.

Этот документ описывает как:

  • использовать R для простых вычислений;
  • загружать данные из шейп-файла и наносить их на карту;
  • производить преобразования координат;
  • добавлять точки на карту.

Запуск R

Чтобы запустить программу нужно:

  • Выбрать R Statistics из раздела меню Spatial Tools — появится окон терминала с запущенным R. Или
  • ввести R в окне консоли. R запустится в этом окне.

Не бойтесь командной строки — она является очень мощным средством. Используйте кнопки вверх и вниз, чтобы вызывать уже введенные команды и исправлять ошибки. Используйте CTRL-C, если застряли и нужно вернуться обратно в консоль.

Выход из R

Почти всё в R является функциями, включая функцию для выхода из программы. Введите q() и нажмите “Enter”. Если вы наберёте просто q, то увидите исходный код для функции q.

R спросит, хотите ли вы сохранить рабочую среду в специальный файл (.RData). Если вы снова запустите R из папки, в которой уже есть этот файл, то он автоматически восстановит оттуда все данные.

Начало работы с R

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

Вы можете начать с простой арифметики:

> 3*2
[1] 6

> 1 + 2 * 3 / 4
[1] 2.5

> sqrt(2)
[1] 1.414214

> pi * exp(-1)
[1] 1.155727

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

Основная системная “подсказка” в консоли выглядит так >, но есть еще “подсказка” для продолжения — +, она появляется, когда R ожидает от вас ещё что-либо для завершения выражения. С ней можно столкнуться, например, если вы забыли скобку или кавычку.

> sqrt(
+ 2
+ )
[1] 1.414214

Структура данных

Вы могли задаться вопросом, что же в результате показывает загадочная единица в квадратных скобках? Единица показывает, что результат — одно число. R может хранить данные в одномерных векторах, двумерных матрицах и многомерных массивах. Они создаются в процессе различных операций. Вот простой пример:

> seq(1,5,len=10)
[1] 1.000000 1.444444 1.888889 2.333333 2.777778 3.222222 3.666667 4.111111
[9] 4.555556 5.000000

Вы можете видеть, что [9] говорит нам, что 4.555 является девятым элементом вектора.

Если вы сделаете матрицу, что получите подписанные ряды и колонки:

> m=matrix(1:12,3,4)
> m
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12

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

> m[2,4]
[1] 11

> m[2,]
[1]  2  5  8 11

> m[,3:4]
     [,1] [,2]
[1,]    7   10
[2,]    8   11
[3,]    9   12

Фреймы данных являются структурами, очень похожими на структуры данных СУБД, таких, как PostgreSQL или MySQL. Каждый ряд может рассматриваться как запись, а колонка — как поле в базе данных. Так же, как и в базе данных, значение в поле должно быть одинакового типа для каждой записи.

Во многом эти структуры работают как матрицы, но вы можете также задавать и обращаться к колонкам, используя знак $:

> d = data.frame(x=1:10,y=1:10,z=runif(10)) # z это 10 случайных чисел
> d
        x  y          z
    1   1  1 0.44128080
    2   2  2 0.09394331
    3   3  3 0.51097462
    4   4  4 0.82683828
    5   5  5 0.21826740
    6   6  6 0.65600533
    7   7  7 0.59798278
    8   8  8 0.19003625
    9   9  9 0.24004866
    10 10 10 0.35972749

> d$z
 [1] 0.44128080 0.09394331 0.51097462 0.82683828 0.21826740 0.65600533
 [7] 0.59798278 0.19003625 0.24004866 0.35972749

> d$big = d$z > 0.6  # d$big это бинарное булево значение true/false
> d[1:5,]
  x y          z   big
1 1 1 0.44128080 FALSE
2 2 2 0.09394331 FALSE
3 3 3 0.51097462 FALSE
4 4 4 0.82683828  TRUE
5 5 5 0.21826740 FALSE

> d$name = letters[1:10] # создать новое символьное поле
> d[1:5,]
  x y          z   big name
  1 1 1 0.44128080 FALSE    a
  2 2 2 0.09394331 FALSE    b
  3 3 3 0.51097462 FALSE    c
  4 4 4 0.82683828  TRUE    d
  5 5 5 0.21826740 FALSE    e

Загрузка геоданных

Существует множество модулей для управления пространственными данными и статистического анализа. Некоторые из них рассмотрены здесь, а другие можно скачать из CRAN.

Давайте загрузим в R два шейп-файла — границы стран и населённые пункты из набора данных Natural Earth. Мы используем два пакета для работы с геоданными — sp и maptools:

> library(sp)
> library(maptools)

> countries = readShapeSpatial("/usr/local/share/data/natural_earth/10m_admin_0_countries.shp")
> places = readShapeSpatial("/usr/local/share/data/natural_earth/10m_populated_places_simple.shp")
> plot(countries)

Мы увидим простую карту Земли:

../../_images/r_plot112.png

Когда OGR-совместимый набор данных загружается в R таким образом, мы получаем объект, который ведёт себя во многом как фрейм данных. Мы можем использовать поле ADMIN для выбора данных, например, только по Великобритании:

> uk = countries[countries$ADMIN=="United Kingdom",]
> plot(uk); axis(1); axis(2)
../../_images/r_plot212.png

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

> proj4string(uk)
[1] NA

NA обозначает, что данные отсутствуют. Нам нужно назначить систему координат объекту до того, как трансформировать его в другую систему координат с помощью функции spTransform из пакета rgdal. Мы будем трансформировать данные в систему координат EPSG:27700, которую использует Ordnance Survey Великобритании:

> proj4string(uk)=CRS("+init=epsg:4326")
> library(rgdal)
> ukos = spTransform(uk,CRS("+init=epsg:27700"))
> proj4string(ukos)
[1] " +init=epsg:27700 +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs
+towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894"

> plot(ukos);axis(1);axis(2)

Результат — карта, показывающая трансформированные данные. Теперь нам нужно добавить несколько точек из набор данных о населённых пунктах. Мы снова отберём нужные нам точки и переведём их в нужную систему координат:

> ukpop = places[places$ADM0NAME=="United Kingdom",]
> proj4string(ukpop)=CRS("+init=epsg:4326")
> ukpop = spTransform(ukpop,CRS("+init=epsg:27700"))

Мы добавим эти точки к нашей карте, масштабировав их по размеру в зависимости от квадратного корня населения (это сделает символ пропорциональным населению), мы также сделаем их красного цвета и используем в качестве символа “залитый” кружок:

> points(ukpop,cex=sqrt(ukpop$POP_MAX/1000000),col="red",pch=19)
> title("UK Population centre sizes")

и наша конечная карта будет выглядеть так:

../../_images/r_plot312.png

Виньетки

В прошлом документация по пакетам R обычно представляла собой скупые описания для каждой функции. Теперь авторам предлагается также писать так называемые “виньетки” — “дружественное к пользователю” введение в использование пакета. Если ввести команду vignette() без аргументов, то можно получить список виньеток, доступных на данный момент в системе. Попробуйте ввести vignette("sp") и вы получите техническое введение в структуры пространственных данных R, или vignette("spdep") чтобы почитать про статанализ пространственных автокорреляций. vignette("gstat") познакомит с использованием пакета gstat для пространственной интерполяции с использованием кригинга.

Дополнительная информация

Общую информацию об R можно найти в официальном Введении в R или в любой другой документации с домашней страницы Проекта R Project.

Подробную информацию о пространственных возможностях R можно найти в R Spatial Task View

Также может оказаться полезным страница R-Spatial на sourceforge, где можно найти полезные ссылки и информацию о листе рассылки R-sig-Geo.