Статистическое программирование на R: Часть 3. Повторное использование кода и объектное программирование

Обзор основных возможностей R

R — это мощная бесплатно распространяемая статистическая среда, которая включает в себя язык программирования, интерактивную оболочку и обширные графические возможности. Эта статья продолжает две предыдущие публикации Дэвида (написанные совместно с Брэдом Хантингом), рассматривая объектно-ориентированное программирование в R, а также некоторые общие концепции программирования в R.

Девид Мертц (David Mertz), Разработчик, Gnosis Software, Inc.

Девид Мертц считает, что искусственные языки вполне естественны, а естественные кажутся немного искусственными. Вы можете связаться с ним по mertz@gnosis.cx; вы можете исследовать все стороны его жизни на его личной web-странице. Посмотрите его книгу, Text Processing in Python. Пожелания и предложения по поводу прошлых и будущих статей только приветствуются.



13.08.2008

Первые две статьи этой серии рассматривали R в применении к «реальному миру». Мы исследовали различные возможности статистического анализа и графического отображения, используя большие наборы данных о температуре, которые собрал один из авторов этих статей. Как было упомянуто в предыдущих статьях, мы фактически обследовали лишь малую часть всей глубины и богатства статистических библиотек R.

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

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

Назад к основам

Перед рассмотрением объектной модели R давайте изучим и уясним объекты и функции R. Главное, что нужно запомнить об объектах R — это то, что «все — вектор». Даже объекты, которые на первый взгляд отличаются от векторов — матрицы, массивы, структуры данных и остальное — на самом деле векторы с дополнительными (измененными) атрибутами, которые «указывают» R обрабатывать их специальным образом.

Массив (обозначаемый dim) — это один из наиболее важных атрибутов, которые имеют (некоторые) вектора R. Функции matrix(), array() и dim() - просто удобные функции для установки размерности вектора. Система ООП (объектно-ориентированного программирования) R также уходит корнями к атрибутам объектных (или других) классов.

Давайте рассмотрим объявление массива на примере кода в Листинге 1.

Листинг 1. Создание векторов и установка размерности
> v = 1:1000
> typeof(v)
[1] "integer"
> attributes(v)
NULL
> dim(v) = c(10,10,10)  # Установка размерности
> attributes(v)
$dim
[1] 10 10 10
> v2 = matrix(1:1000, nrow=100, ncol=10)
> typeof(v2)
[1] "integer"
> attributes(v2)
$dim
[1] 100  10
> attr(v2,'dim') = c(10,10,10)  # Изменение размерности
> attributes(v2)
$dim
[1] 10 10 10

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

Единственное, что может смутить в положении «все — вектор» — это построчные и постолбцовые операции, которые могут быть не так интуитивно понятны. Например, достаточно легко создать двумерный массив (матрицу) и поработать с одним столбцом или строкой:

Листинг 2. Построчные операции над матрицей
> m = matrix(1:12, nrow=3, ncol=4)
> m
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12
> sum(m)  # суммирование всех элементов матрицы m
[1] 78
> sum(m[1,])  # суммирование по первой строке
[1] 22

Но если бы вы захотели создать вектор сумм по строкам, вам, возможно, захотелось бы сделать что-либо вроде следующего:

Листинг 3. Неправильный способ выполнения построчных операций
> sum(m[c(1,2,3),])  # это НЕ суммирование строк
[1] 78

Здесь можно было бы использовать цикл, но это идет вразрез с функциональными и векторно-ориентированными операциями R. Правильное решение - использовать функцию apply():

Листинг 4. Функция apply() для построчных операций
> apply(m, 1, sum) # по строке
[1] 22 26 30
> apply(m, 2, sum) # по столбцу
[1]  6 15 24 33
> apply(m, c(1,2), sum) # по столбцу и строке (суммирование каждой отдельной ячейки)
     [,1] [,2] [,3] [,4]
[1,]    1    4    7   10
[2,]    2    5    8   11
[3,]    3    6    9   12
# Бесполезно суммировать каждую отдельную ячейку, но можно сделать так:
> a = array(1:24,c(3,4,2))
> apply(a, c(1,2), sum)  # суммирование вглубь в 3-х мерном массиве
     [,1] [,2] [,3] [,4]
[1,]   14   20   26   32
[2,]   16   22   28   34
[3,]   18   24   30   36
> apply(a, c(3), sum)    # суммирование по каждому направлению вглубь
[1]  78 222

Бесконечная последовательность

Иногда, из чисто практических соображений, полезно использовать такую конструкцию, как бесконечная числовая последовательность. Например, мой соавтор по предыдущей статье делал некий анализ с помощью метода интегрирования Монте-Карло, и для его задачи было полезно использовать бесконечно длинную последовательность случайных чисел. Необходимо понять, что бесконечная последовательность — это не просто возможность генерации нового числа в случае необходимости; это также необходимость иметь возможность обращаться к любому предшествующему элементу и получать его прежнее значение.

Очевидно, нет ни компьютерного языка, ни компьютера, способного хранить бесконечную последовательность — все, что они могут хранить, - это «ленивые» (lazy) и неограниченные (unbounded) последовательности. Новые элементы добавляются к уже реализованному списку тогда и только тогда, когда это необходимо. Например, в Python можно выполнить такое создание спископодобных объектов методом .__getitem__(), который расширяет внутренние списки по необходимости. В Haskell «ленивые» последовательности встроены глубоко внутрь языка — в результате все последовательности являются «ленивыми». В моем руководстве по Haskell (Ресурсы) я использовал пример создания списка всех простых чисел:

Листинг 5. Список простых чисел в Haskell, полученных с помощью решета Эратосфена
primes :: [Int]
primes = sieve [2 .. ]
sieve (x:xs) = x : sieve [y | y <- xs, (y `rem` x)/=0]

В том, что касается бесконечности, R ближе к Python, чем к Haskell. Вам необходимо явно создавать новые элементы в случае необходимости. Подождем до раздела ООП, где мы увидим, как работает скрытый механизм векторного индексирования; здесь еще не так все запутано.

Листинг 6. Объявление вектора и способа его динамического расширения
inf_vector = rnorm(10, 0, 1)   # произвольный старт w/ ,
                               # инициализация 10 элементов
assure <- function(index) {
  extend_by = max(index-length(inf_vector), 0)
  extras = rnorm(extend_by, 0, 1)
  v <- c(inf_vector, extras)
  assign("inf_vector", v, env=.GlobalEnv)
  return(index)
}
getRand <- function(index) {
  assure(index)
  return(inf_vector[index])
}

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

Листинг 7. Использование функции-контейнера в качестве посредника для бесконечного виртуального вектора
> getRand(3)                # единичный индекс
[1] 0.5557101
> getRand(1:5)              # диапазон
[1] -0.05472011 -0.30419695  0.55571013  0.91667175 -0.40644081
> getRand(sqrt(c(4,16)))    # расчет индекса коллекции
[1] -0.3041970  0.9166717
> getRand(100)              # принудительное расширение вектора
[1] 0.6577079

Если необходимо, перед использованием элементов можно создать достаточно большой вектор с помощью функции assure():

Листинг 8. Расширение вектора (при необходимости) перед работой с ним
> assure(2000)
[1] 2000
> inf_vector[1500]
[1] 1.267652

Объектно-ориентированное программирование в R

R полностью отвечает требованиям объектно-ориентированного программирования, но чтобы понять это, нужно вспомнить о том, что такое ООП. Пользователи таких языков, как Java™ и C++ и даже Python, Ruby или Smalltalk, могут иметь несколько ограниченное представление об объектно-ориентированном программировании. Не ошибочное, но ограниченное одной моделью.

Принципы ООП в R в большей степени основываются на обобщенных функциях, чем на иерархиях классов. Эта концепция будет близка читателям, использующим CLOS Lisp или тем, кто читал мои дискуссии по множественной диспетчеризации в Python (Ресурсы). К сожалению, подход R — это единичная диспетчеризация; в этом отношении он эквивалентен «традиционным» языкам, упомянутым ранее (C++, Java и другие).

Необходимо заметить, хотя это и не будет обсуждаться подробно в данной статье, что наиболее новая версия R сопровождается пакетом methods, который определяет и использует так называемые «формальные методы». Применение этих формальных методов во многом накладывает некую дисциплину (и ограничения), знакомые по традиционным языкам ООП. В любом случае формальное ООП в R строится поверх «неформального ООП», о котором будет рассказано в этой статье. Насколько я знаю, пакет methods не имеет окончательного статуса, но некая слегка модифицированная версия этого пакета, по-видимому, будет сохранена в будущих версиях R. В разделе Ресурсы можно найти более подробную информацию по этому вопросу.

Необходимо помнить, что суть концепции ООП на самом деле не в наследовании, а в более общем принципе - решениях по диспетчеризации . Например, вызов obj.method() в традиционных языках ООП будет использовать порядок разрешения методов (method resolution order - MRO) объекта для поиска «первого» класса-предка obj, который имеет метод .method().

Что такое «первый» - это более тонкий вопрос, чем кажется на первый взгляд (в разделе Ресурсы представлено хорошее обсуждение развития MRO в Python). R принимает те же решения, но выворачивает идею наследования наизнанку — вместо набора классов, которые могут объявлять и аннулировать различные методы внутри себя, R порождает семейство обобщенных (generic) функций, которые имеют метки, указывающие тип объекта, которым они хотят оперировать.

Обобщенные функции

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

Листинг 9. Создание обобщенной функции и помеченных методов
#------------- Создание обобщенного  метода
> whoami <- function(x, ...) UseMethod("whoami")
> whoami.foo <- function(x) print("I am a foo")
> whoami.bar <- function(x) print("I am a bar")
> whoami.default <- function(x) print("I don't know who I am")

Ключевая идея состоит в том, что каждый объект в R может принадлежать нулю, одному или большему числу классов. MRO любого заданного объекта (относительно конкретного метода) — это просто вектор именованных классов (если они есть) в атрибуте его class. Например:

Листинг 10. Назначение объектам меток членства в классе
> a = 1:10
> b = 2:20
> whoami(a)                 # Нет соответствующего класса
[1] "I don't know who I am"
> attr(a,'class') <- 'foo'
> attr(b,'class') <- c('baz','bam','bar')
> whoami(a)
[1] "I am a foo"
> whoami(b)                 # поиск MRO для описываемого вектора
[1] "I am a bar"
> attr(a,'class') <- 'bar'  # изменение класса 'a'
> whoami(a)
[1] "I am a bar"

Как и в традиционных наследующих языках, объект не обязан использовать один и тот же класс для всех вызываемых методов. Традиционно, если Child наследуется от Mom и Dad, объект типа Child может использовать .meth1() от Mom и .meth2() от Dad. Все это можно сделать и в R, но Mom и Dad ничего не значат — это просто имена:

Листинг 11. Разрешение методов
> meth1 <- function(x) UseMethod("meth1")
> meth1.Mom <- function(x) print("Mom's meth1")
> meth1.Dad <- function(x) print("Dad's meth1")
> meth2 <- function(x) UseMethod("meth2")
> meth2.Dad <- function(x) print("Dad's meth2")
> attr(a,'class') <- c('Mom','Dad')
> meth1(a)  # несмотря на существование meth1.Dad, Mom используется первым для a
[1] "Mom's meth1"
> meth2(a)
[1] "Dad's meth2"

Включение предков

Необходимость явного указания MRO объекта вместо неявного разрешения через синтаксис наследования может показаться ограничивающей. В действительности можно легко реализовать основанный на наследовании синтаксис MRO, используя минимальное количество функций-контейнеров. MRO, используемый в Листинге 11, вероятно, не лучший из всех возможных (см. очерк Симионато в Ресурсах), но он демонстрирует идею:

Листинг 12. Реализация основанного на наследовании MRO с использованием минимального количества функций-контейнеров
char0 = character(0)
makeMRO <- function(classes=char0, parents=char0) {
    # Создание MRO из опционального явного списка
    # и опционального списка предков
    mro <- c(classes)
    for (name in parents) {
        mro <- c(mro, name)
        ancestors <- attr(get(name),'class')
        mro <- c(mro, ancestors[ancestors != name])
    }
    return(mro)
}
NewInstance <- function(value=0, classes=char0, parents=char0) {
    # Создание нового объекта, основываясь на первоначальном значении,
    # явных классах и предках (все опционально)
    obj <- value
    attr(obj,'class') <- makeMRO(classes, parents)
    return(obj)
}
MaternalGrandma <- NewInstance()
PaternalGrandma <- NewInstance()
Mom <- NewInstance(classes='Mom', parents='MaternalGrandma')
Dad <- NewInstance(0, classes=c('Dad','Uncle'), 'PaternalGrandma')
Me <- NewInstance(value='Hello World', 'Me', c('Mom','Dad'))

В действии код Листинга 12 выглядит следующим образом:

Листинг 13. Объект с основанным на наследовании MRO
> print(Me)
[1] "Hello World"
attr(,"class")
[1] "Me"              "Mom"             "MaternalGrandma" "Dad"
[5] "Uncle"           "PaternalGrandma"

Если следовать традиционному подходу отношений класс/наследование, надо включать имя создаваемого класса (подобно Mom в аргументе classes). Фактически, учитывая, что каждый класс — это нормальный объект, рассмотренная выше система ближе к прототипному ООП, чем к основанному на классах.

Опять же, вся система обладает достаточной гибкостью для реализации всех вариантов. Можно при желании совершенно свободно отделить объекты классов от объектов экземпляров — можно различать классы по соглашениям именования (например, Mom в отличие от mom), присваивая другие атрибуты (например, type может быть class или instance; функция обработки должна проверять тип) или другими способами.

Снова о бесконечном векторе

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

Операторы в R – это просто сокращенный способ вызова функций; вы можете свободно дифференцировать поведение операторов на основе классов, так же, как и для вызовов других функций. Попутно исправим еще несколько недостатков первой системы:

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

Сделаем все это:

Листинг 14. Объявление индексируемого бесконечного случайного вектора
"[.infinite_random" <- function(v, index) {
    name <- attr(v, 'name')
    rfunc <- attr(v, 'rfunc')
    extend_by = max(index-length(v), 0)
    extras = rfunc(extend_by)
    new <- c(v, extras)
    makeInfiniteRandomVector(name, v=new, rfunc)
    return(new[index])
}
unitnorm <- function(n) return(rnorm(n,0,1))
empty <- vector('numeric', 0)
makeInfiniteRandomVector <- function(name, v=empty, rfunc=unitnorm) {
    # Создание бесконечного вектора
    # можно расширить существующий вектор, настраиваемая функция rand
    attr(v,'class') <- 'infinite_random'
    attr(v,'name') <- name
    attr(v,'rfunc') <- rfunc
    assign(name, v, env=.GlobalEnv)
}
makeInfiniteRandomVector('v')
# makeInfiniteRandomVector('inf_poisson', rfunc=my_poisson)
# Использование: v[1]; v[10]; v[9:12]; и так далее.

Индексирование уже определено в R как обобщенная функция, так что для его настройки не нужно вызывать метод UseMethod(); можно просто определить столько новых специализаций, сколько нужно. Аналогично, встроенная функция print() тоже является обобщенной. Ею можно воспользоваться так:

Листинг 15. Печать бесконечного вектора
print.infinite_random <- function(v) {
    a_few = 5
    len = length(v)
    end_range = (len-a_few)+1
    cat('* Infinite Random Vector *\n')
    cat('[1] ', v[1:a_few], '...\n')
    cat('[')
    cat(end_range)
    cat('] ', v[end_range:len], '...\n')
}

In action, the code produces the following:

Листинг 16. Пример печати бесконечного вектора
> v[1000]
[1] -1.341881
> print(v)
* Бесконечный случайный вектор *
[1]  -0.6247392 1.308057 1.654919 1.691754 -2.251065 ...
[996]  1.027440 0.8376 -0.7066545 -0.7778386 -1.341881 ...

Заключение

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

Вся хитрость такого понимания ООП в том, чтобы мыслить в терминах "какой код вызван?" и "как сделан выбор?". Не привязывайтесь к синтаксису, который применяется в конкретных языках - C++, Objective C, Java, Ruby или Python; сконцентрируйтесь на концепции диспетчеризации.

Ресурсы

Научиться

Получить продукты и технологии

Комментарии

developerWorks: Войти

Обязательные поля отмечены звездочкой (*).


Нужен IBM ID?
Забыли Ваш IBM ID?


Забыли Ваш пароль?
Изменить пароль

Нажимая Отправить, Вы принимаете Условия использования developerWorks.

 


Профиль создается, когда вы первый раз заходите в developerWorks. Информация в вашем профиле (имя, страна / регион, название компании) отображается для всех пользователей и будет сопровождать любой опубликованный вами контент пока вы специально не укажите скрыть название вашей компании. Вы можете обновить ваш IBM аккаунт в любое время.

Вся введенная информация защищена.

Выберите имя, которое будет отображаться на экране



При первом входе в developerWorks для Вас будет создан профиль и Вам нужно будет выбрать Отображаемое имя. Оно будет выводиться рядом с контентом, опубликованным Вами в developerWorks.

Отображаемое имя должно иметь длину от 3 символов до 31 символа. Ваше Имя в системе должно быть уникальным. В качестве имени по соображениям приватности нельзя использовать контактный e-mail.

Обязательные поля отмечены звездочкой (*).

(Отображаемое имя должно иметь длину от 3 символов до 31 символа.)

Нажимая Отправить, Вы принимаете Условия использования developerWorks.

 


Вся введенная информация защищена.


  • Bluemix

    Узнайте больше информации о платформе IBM Bluemix, создавайте приложения, используя готовые решения!

  • Библиотека документов

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

static.content.url=http://www.ibm.com/developerworks/js/artrating/
SITE_ID=40
Zone=Linux
ArticleID=330589
ArticleTitle=Статистическое программирование на R: Часть 3. Повторное использование кода и объектное программирование
publish-date=08132008