Economica: Mathematica с MATLAB-акцентом

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

Так что решил перенести пару функций из MATLAB в Mathematica.

linspace[x,y,n] – дает набор равномерных цифр от x до y размером n. x и y могут быть векторами или матрицами (для чисел отсюда, для векторов и матриц дописал сам):

[wlcode]linspace[x_?NumericQ, y_?NumericQ, n_Integer] := Array[# &, n, {x, y}];
linspace[x_List, y_List, n_Integer] :=
Module[{temp},
If[Dimensions@x != Dimensions@y, Abort[]];
temp = MapThread[linspace[#1, #2, n] &, {x, y}, 1];
temp[Transpose]
];[/wlcode]
Для получения нулевой матрицы можно использовать zeros, для единичной матрицы – eye:
[wlcode]zeros[x : Sequence[Integer__]] := ConstantArray[0, {x}];
eye[x_Integer] := IdentityMatrix[x];[/wlcode]Для того, чтобы избавиться от списков с одним елементом нужно использовать squeeze (отсюда):

[wlcode]squeeze[A_] := ArrayReshape[A, Dimensions[A]~DeleteCases~1][/wlcode]

Все они доступны в Economica.

Economica: Mathematica с MATLAB-акцентом

Economica: Данные Quandl, сезонное сглаживание и тест Грейнджера

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

Немногое изменилось и сегодня, но в моем запасе теперь несколько советов Леонида Шифрина  и Дэвида Вагнера. 
Новость же в том, что сегодня добавлены функции, которые могут быть действительно полезными (в ссылках инструкции по их использованию в pdf):
1. QuandlData – получение данных с сайта Quandl. Если вы зарегистрируетесь и получите собственный токен, то можете использовать с помощью функции QuandlToken
2. Сезонное сглаживание на основе официального алгоритма SeasonalAdjustmentX13 от Бюро переписей США (мой ответ на mma.se на эту тему).
3. Тест Грейнджера для двух рядов GrangerCausalityTest.
Для установки пакета следует распаковать архив в вашей $UserBaseDirectory. В справке должны появиться страницы, соответствующие функциям пакета.

Все это находиться на очень раннем этапе, и я буду рад любым откликам. 
Economica: Данные Quandl, сезонное сглаживание и тест Грейнджера

MMA: Данные по генерации и потреблению электроэнергии

Аналитики зачастую используют статистику по генерации и потреблению электричества для построения оценок интенсивности активно экономической деятельности. Нередко данные по электричеству могут улучшить оценки индекса промышленного производства или ВВП.

Статистика по генерации и потреблению доступна на сайте системного оператора ЕЭС здесь.

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

Например:

Более того, решил, что начну систематизировать накопившиеся подобные кусочки в пакет Economica вот здесь. Он пока содержит одну эту функцию, то есть находится на уровне задумки.

MMA: Данные по генерации и потреблению электроэнергии

Mathematica: ARIMA-X13

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

Основная запись опубликована на замечательном сайте Романа Осипова, но комментарии к ней можно оставить и здесь.
Немного позже мне хотелось бы собраться и с помощью этого инструмента повторить результаты впечатлившей меня работы Бессонова и Петроневич о проблемах, к которым может привести использование сезонного сглаживания при его использовании в режиме “не приходя в сознание”.
Mathematica: ARIMA-X13

График для индексов с особой легендой

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

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

Попробуем воспроизвести подобный график в Mathematica.  Здесь довольно легко добиться близкого вида, используя BarLegend, но, после несколько попыток поместить легенду в необходимой близости от графика, мне показалось сделать собственную легенду (см. пример в приложении).

Например, следующий график показывает динамику (лог-)выпуска на человека в долларах  нормированного на значение в 2000 г. для каждой из стран (данные МВФ):

Номинальный ВВП на душу в текущих долларах к уровню 2000 г., МВФ, собств. расчет.
График для индексов с особой легендой

Heapsters gonna heapst

На прошлой неделе нашел немого времени просмотреть материалы UDACITY CS215 Crunching Social Networks и должен сказать, что курс, если и неплохой, то не оправдывает своего названия. Вместо того, чтобы действительно поделиться подходами к анализу структуры отношений между вершинами графа (будь они люди или, скажем, публикации) Михаил Литтман рассказывает об алгоритмах поиска и сортировки, об их сложности и скорости, что на самом деле не так интересно, особенно если набор данных не астрономических масштабов.

Ниже хочу поделиться примером решения игрушечной задачки, которую мне бы хотелось уметь решать лучше:
Возьмем данные о ссылках американских политических блогов друг на друга[1] со страницы Марка Ньюмана здесь и попробуем найти самых влиятельных блоггеров.

Решение
Попробуем напрямую импортировать данные, надеясь, что Mathematica сама все поймет:
lesmis = Import[“polblogs.gml”];

Программа покажет ошибку, которая говорит о том, что граф сложный. Обычно это значит, что одно из ребер повторяется два раза*. Предлагаю избавиться от дублей, импортировав граф в виде текста и удалив все задвоения. Если открыть сам файл “.gml” в блокноте можно увидеть, что устроен он несложно. Итак:
lesmis = Import[“polblogs.gml”, “Text”];

s1 = Reverse[
   StringCases[lesmis,
    RegularExpression[“(?s)(node|edge) [.*? ]”]]];
s2 = StringJoin[Riffle[Union[s1], “n”]];

Export[FileNameJoin[{NotebookDirectory[], “better.gml”}],
 StringJoin[
  “Creator “Lada Adamic on Tue Aug 15 2006″ngraph [n  directed 1n
“,
  s2, “n]”], “Text”]

Попытаемся импортировать граф снова:
lesmis2 = Import[“better.gml”];
…и все впорядке.

Теперь возможно увидеть на кого ссылались много (out-degree) и какие блоги сами часто ссылались на других (out-degree):

(нажать, чтобы увеличить)

Из графика видим, что безусловным лидером был blogsforbush.com в обеих номинациях, но довольно много ссылок на себя получил и atrios.blogspot.com, на который частенько ссылается Кругман.

У меня есть ощущение, что было бы интересно (и не сложно) собрать данные и сделать граф ссылок между дневниками наших местных экономистов, но проблема в том, что мне не известно так уж много хороших дневников с которых можно было бы начать. Есть идеи? Это не должны быть обязательно очень популярные блоги, взгяните на график – большинство из дневников имеют в сумме менее 50 ссылок.

Скачать расчет и данные.

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

Источник
[1] Lada A. Adamic and Natalie Glance, “The political blogosphere and the 2004 US Election”, in Proceedings of the WWW-2005 Workshop on the Weblogging Ecosystem (2005).

Heapsters gonna heapst

Python 101


 В личных целях нужно научиться пользоваться COM из Python, к чему должна привести книжка справа —>

Однако, как оказалось, нужно начинать с самого начала – мой навык письма на Python на столько плох, что однажды пришлось переписать готовый код для раскрашивания карт только потому, что не смог понять как работает ‘import’.
В части “Problem solving” попавшегося под руку/замечательного “Byte of Python” предлагается использовать  GNUwin32 для того, чтобы архивировать содержание каталога в zip с указанием даты и других подробностей. 
Решение этой задачки дается в BoP, но так же легко возможно проделать это и в Mathematica, используя   Run[], которая походит на os.system():

dirA = "d:\a\" <> DateString[Date[], {"Year", "Month", "Day"}] <> 
".zip"
dirB = "d:\b"
cmd = "zip -qr " <> dirA <> " " <> dirB
Run[cmd]

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

Python 101

Неслучайные связи

На прошлой неделе написал функцию для тестирования наличия причинность по Грейнджеру для пары переменных. Код и ссылка на пример ниже. 
На всякий случай проверил результаты с Eviews на выборке, которую предлагает Jamie Monogan из Washington University здесь: результат одинаковый:

Пример

Для наглядности приведу игрушечный пример использования теста Грейнджера.
Из опыта нам известно, что с ростом безработицы для многих выпускников аспирантура становится привлекательнее (скажем, из-за снижения альтернативных потерь: более низкой оплаты труда/вероятности трудоустройства на рынке).
Проверим это предположение. Данные по безработице из Вышки и данные по статистике образования. Я перехожу к приростам, так данные должны быть стационарны.
Действительно, вероятность того, что безработица не (Грейджер-)приводит к изменению в приростах в количестве поступающих в аспиратуру чуть более пяти процентов и не стоило бы ее отклонять. Так же логично, что обратное не верно, и мы можем уверенно отклонить гипотезу о том, что изменения в приростах поступающих могут влиять на изменения в приростах уровня безработицы. 

Было бы интересно:
1. Сделать пример с поиском потенциально интересных причинно-следственных связей в массиве переменных.
2. Сделать тест для >2 переменных.
Неслучайные связи

Рынок услышан

Петр обратил внимание на популярный пост smart-lab с предложением слушать котировки рынка.  Если постараться не задумываться о полезности (или последствиях) подобного времяпрепровождения, то остается только найти простой способ идею воплотить в жизнь.

К счастью, для этого не нужно “гуглить спектральный анализ, преобразование Фурье, вейвлеты” [1].

1. Используем данные, которые предлает Финам. Для примера я выбрал исторические котировки индекса ММВБ.

2. В Mathematica используем ListPlay[x] для того, чтобы получить звук из котировок.

3. (по желанию) Этот звук возможно экспортировать, напримерв в wav:
Export[“micex.wav”,ListPlay[x]]

Скачать .nb.

1.Rorschach, Надоело смотреть на графики, решил их послушать, smart-lab.ru, 2012

Рынок услышан

Замкнутый круг

Процедура под названием “усреднение” позволяет кредитным организациям держать не все обязательные резервы в Банке России постоянно, а лишь их часть (например, сейчас 40% от размера обязательствых резервов).

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

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

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

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

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

Угол – (номер дня усреднения/количество дней усреднения в данном периоде)*2 Pi,
расстояние от центра – (Ставка RUONIA в данный день)/(Средняя ставка  RUONIA за данный период усреднения).
Всего периодов усреднения – 26.

Черной линией показана оценка (устойчивой к выбросам) линейной регрессии. Казалось бы сезонности нет. Для того чтобы проверить себя, отложим наши данные в привычной декартовой плоскости:

Хочу отметить, что периоды усреднения не совпадают с календарными месяцами, а проходят с 10 по 10 число месяцев.
Синей линией показана скользящая средняя по процентной ставке.

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

Мне наглядное доказательство наличия сезонности кажется убедительным.

Какие формальные или наглядные тесты Вы используете для того, чтобы проверить наличие сезонности?

Приложение

  1. Расчеты в .nb
  2. Данные в .xlsx
Замкнутый круг