Постановка задачи
Не буду вдаваться в подробности о том, откуда берутся миллионы временных серий и почему они умудряются изменяться еженедельно. Просто возникла задача еженедельно сделать прогноз на 2-8 недель по паре миллионов временных серий. Причем не просто прогноз, а с кроссвалидацией и выбором наиболее оптимальной модели (ARIMA, нейронная сеть, и т.п.).
Имеется свыше терабайта исходных данных и достаточно сложные алгоритмы трансформации и чистки данных. Чтобы не гонять большие массивы данных по сети решено было реализовать прототип на одном сервере.
Первый блин комом
На 32-х ядерный виртуальный сервер с 256 ГБ оперативной памяти был установлен MS SQL 2016 Standart и R service. В принципе, вместо R вполне можно было использовать Python или даже CLR. Но аналога auto.arima в Python тогда найдено не было (может появился уже?). А так как сервер должен был использоваться не только для этой задачи, рисковать стабильностью его работы вынося код в CLR не хотелось.
Трансформация и очистка данных заняла на нем порядка 6 часов. После чего, средствами sp_execute_external_script вызывался код на R, выполняющий уже восстановление пропусков в временных сериях фильтром Калмана, обучающий несколько моделей и сравнивающий результаты их прогнозов кроссвалидацией. Возвращался результат прогноза выигравшей модели.
При этом возник целый ряд проблем с производительностью.
Если вызывать sp_execute_external_script, передавая ему по одной временной серии, то только на вызов любой функции сервиса внешних языков уходило около 30 мс. Вроде бы немного, но для пары миллионов вызовов это получалось уже свыше 16 часов. При этом редакция MS SQL Standard поддерживает только 2 ядра при интеграции с R.
Можно вызывать sp_execute_external_script, передавая ему множество временных серий. Но для этого надо объединять их через PIVOT в матрицу, что MS SQL делал очень и очень медленно.
Попытка реализовать параллельное обучение моделей в R так и не увенчалась успехом. Какими бы средствами я не пользовался, периодически один из рабов просто отваливался. А на миллионах вызовов - отваливался всегда. Даже просто включение более одного ядра в функции auto.arima опять таки, раз в ~100 тыс. вызовов, приводило к молчаливому возврату не обученной модели.
Как результат, на двух ядрах процесс прогнозирования пары миллионов временных серий не успевал уложиться в неделю, при том что прогнозировать надо еженедельно.
Сложности выбора
Первой возникшей мыслью был переход на MS SQL Enterprise. Было произведено даже исследование, показавшее, что на MS SQL Developer результат можно получить, примерно, за трое суток. Вот только столь существенное увеличение бюджета (Enterprise стоит намного дороже Standard) согласовать не удалось.
Решение проблемы
Так как уже были готовы алгоритмы трансформации и очистки данных на T-SQL и код обучения моделей и кроссвалидации на R, то логичным вариантом решения оказался PostgreSQL, напрямую, в своей среде, поддерживающий выполнение функций на R. Поэтому на ту же виртуальную машину был установлен Oracle Linux, PostgreSQL и R.
И действительно, временные затраты на вызов R функций у PostgreSQL оказались ниже 1 мс. При этом проблемы параллельных вычислений средствами R сохранились. Поэтому решено было передавать функции на R по одной временной серии при каждом вызове.
Просто так выполнять запрос с функцией на R на множестве ядер PostgreSQL не захотел. Поэтому, средствами dblink() запрос вручную параллелится на заданное количество ядер.
Подробности реализации
Была написана довольно простая процедура
CREATE OR REPLACE PROCEDURE perform_query_async_sp (
slave_query text,
cursor_query text,
slaves smallint=16,
batch_size int=64 ) AS $proc$
<<proc>>
DECLARE
dbname text=current_database();
batches_cur refcursor;
slave_list text[];
slave_connection text;
conn_status text;
is_not_dispatched boolean;
in_list text;
BEGIN
FOR i IN 1..slaves LOOP
slave_connection='Slave'||i::text;
conn_status=dblink_connect(slave_connection,'dbname='||dbname);
IF conn_status<>'OK' THEN
INSERT INTO PBD_ExecutionLog (LogLevel, SourceName, LogMessage, LogCharData)
VALUES (32, 'perform_query_async_sp', format('Slave %s connection error',slave_connection), slave_connection);
COMMIT AND CHAIN;
RAISE EXCEPTION 'Slave % connection error', slave_connection;
END IF;
slave_list=array_append(slave_list,slave_connection);
END LOOP;
OPEN batches_cur FOR EXECUTE cursor_query USING batch_size;
<<cursor_loop>>
LOOP
FETCH batches_cur INTO in_list;
EXIT cursor_loop WHEN NOT FOUND;
is_not_dispatched=TRUE;
WHILE is_not_dispatched LOOP
<<slaves_loop>>
FOREACH slave_connection IN ARRAY slave_list LOOP
IF dblink_is_busy(slave_connection)=0 THEN
conn_status='';
WHILE conn_status IS NOT NULL LOOP
SELECT res FROM dblink_get_result(slave_connection) AS R (res text) INTO conn_status;
END LOOP;
IF dblink_send_query(slave_connection,format(slave_query,in_list))=0 THEN
INSERT INTO PBD_ExecutionLog (LogLevel, SourceName, LogMessage, LogCharData)
VALUES (32, `{OBJECT_NAME`}, format('Slave %s send query error',slave_connection), slave_connection);
COMMIT AND CHAIN;
RAISE EXCEPTION 'Slave % send query error', slave_connection;
END IF;
is_not_dispatched=FALSE;
EXIT slaves_loop;
END IF;
END LOOP;
IF is_not_dispatched THEN
PERFORM pg_sleep(0.01);
END IF;
END LOOP;
END LOOP;
CLOSE batches_cur;
FOREACH slave_connection IN ARRAY slave_list LOOP
WHILE dblink_is_busy(slave_connection)=1 LOOP
PERFORM pg_sleep(0.01);
END LOOP;
conn_status='';
WHILE conn_status IS NOT NULL LOOP
SELECT res FROM dblink_get_result(slave_connection) AS R (res text) INTO conn_status;
END LOOP;
SELECT res FROM dblink_get_result(slave_connection) AS R (res text) INTO conn_status;
conn_status=dblink_disconnect(slave_connection);
END LOOP;
END; $proc$ LANGUAGE plpgsql;В параметре slave_query процедуре передается произвольный SQL запрос, обязанный содержать один и только один placeholder %s. Обычно, в виде конструкции IN (%s). В параметре cursor_query передается текст запроса для курсора, который обязан возвращать одну строку (varchar), которая будет подставляться вместо %s при передаче запроса из slave_query очередному рабу. Текст курсора обязан содержать placeholder $1, в который подставляется значение параметра batch_size. Таким образом проще управлять количеством элементов списка, возвращаемых курсором в строке.
Пример строки cursor_query:
SELECT array_to_string(array_agg(Q.VectorId),',') AS in_list
FROM (
SELECT V.VectorId, ROW_NUMBER() OVER (ORDER BY V.VectorId) AS rn
FROM tmp_global_statistics_data_time_series_arrays V
ORDER BY V.VectorId) Q
GROUP BY Q.rn/$1 ORDER BY Q.rn/$1 DESC;Следует отметить, что права на вызов процедуры perform_query_async_sp остаются только у административной роли. Соответственно, вызывать эту процедуру могут только процедуры созданные с правами администратора и с опцией SECURITY DEFINER. Все же через dblink() без аутентификации исполнение запросов происходит под административными правами и это надо учитывать.
Хочу обратить внимание, что если таким путем множеством рабов заполнять одну и ту же таблицу, то желательно, чтобы она была нежурналируемая (UNLOGGED). Если таблица создается в том же соединении перед вызовом процедуры, то обязательно нужно между созданием таблицы и вызовом процедуры выполнить COMMIT.
Кроме того, благодаря поддержке массивов в PostgreSQL удалось существенно упростить трансформацию и очистку данных, а так же в разы уменьшить объемы промежуточных таблиц.
В результате сейчас свыше двух с половиной миллионов прогнозов вычисляются, примерно, за 16 часов с ограничением на использование 24 ядер. Так как за воскресенье, когда нагрузка от пользователей минимальна и 8 оставшихся ядер им достаточно, один сервер вполне управляется с этой задачей, то прототип как то сам собой превратился в продуктивное решение.
Комментарии (23)
dgoncharov
26.08.2023 14:53Интересно, но как-то не очень понятно изложено. Я правильно понял, что у вас есть 2*10^6 временны'х рядов (какой длины - не понял) и вы для каждого ряда строите несколько моделей, выбираете из них лучшую и потом делаете прогноз на 8 недель? Ну и какой метод (ARIMA или что-то еще) показывает себя лучше? Правильно я понял, что ряды независимы, так что с распараллеливанием никаких принципиальных проблем нет? Технические детали реализации, честно говоря, не очень интересны, т.к. воспроизводить такие вычисления вряд ли кто-то тут станет.
свыше двух с половиной миллионов прогнозов вычисляются
Эти два миллиона надо еще умножить на число методов, вероятно?
ptr128 Автор
26.08.2023 14:53+1Так цель статьи была не рассказать о самом прогнозировании по миллионам временных серий, а показать, как особенности PostgreSQL позволяют решать эту задачу существенно эффективней и дешевле, чем средства MS SQL или даже Spark (терабайты по сети гонять медленней, чем в оперативке).
Да, для каждого из 2.5 миллионов временных рядов выполняется прогнозирование на 2-8 недель (параметризируемо). Собственно говоря именно из-за обилия параметров, влияющих на количество временных серий и фильтрацию фактов, и не получилось просто дообучать модели свежими данными.
Ну и какой метод (ARIMA или что-то еще) показывает себя лучше?
Зависит от временного ряда. На длинных временных рядах в лидерах ARIMA и нейронная сеть, примерно поровну. На средних - явно лидирует ARIMA. На коротких - ARIMA начинает иногда проигрывать уже ETS (экспоненциальному сглаживанию), HW (модель Хольта-Винтера) или даже простейшим методам (медиане, среднему между 0.25. и 0.75 перцентилем, наивному прогнозу между 0.25. и 0.75 перцентилем).
Эти два миллиона надо еще умножить на число методов, вероятно?
Почти. Еще плюс два-четыре. Сначала фильтром Калмана заполняются пропуски, что уже требует обучения ARIMA. Затем делается прогноз без учета факта за последние N недель (где N - требуемое количество недель прогноза). Выполняется вычисление прогноза по каждой модели и сверяется с фактом за эти N недель. После чего уже выигравшая модель обучается по полной временной серии. Если результат прогноза вышел за пределы 0.05 и 0.95 перцентиля, то прогноз бракуется и обучается всей временной серией уже следующая по качеству кроссвалидации модель.В редких случаях, возможен аналогичный переход уже к третьей по очереди модели. Глубже не видел.
В планах увеличить глубину кроссвалидации для оценки сезонности. Тогда количество расчетов существенно возрастет.
dgoncharov
26.08.2023 14:53Спасибо за ответ, очень интересно. Почитал также другие ветки, стало еще интереснее. "Модели обучаются не всеми фактам, а только признанными нормальным" - т.е. вы сначала как-то отделяете сигнал от шума и исключаете outliers, так? Согласен с другим комментатором, что распределение вероятностей было бы более информативным, чем мат.ожидание, особенно если вам надо создавать резерв ресурсов (вагонов) с заданной наперед вероятностью, что ресурсов хватит. Я, конечно, не знаю достоверно, но предполагаю, что тут запросто может быть распределение с "черными лебедями".
ptr128 Автор
26.08.2023 14:53как-то отделяете сигнал от шума и исключаете outliers, так?
Не совсем так. У меня же на входе результат трекинга вагонов не только с временем и их местоположением, но и с множеством других аналитик. Например, если вагон выехал из пункта А в пункт Б и посередине в пункте С была операция бросания, отцепа или перевода в неисправный парк, то я и без математических методов могу забраковать все события с этим вагоном до отправления с пункта С после бросания или отцепа. Хотя, конечно, это не отменяет выявления нештатной ситуации с вагоном просто по его длительной задержке на промежуточной станции.
распределение вероятностей было бы более информативным, чем мат.ожидание, особенно если вам надо создавать резерв ресурсов (вагонов) с заданной наперед вероятностью
Распределение вероятностей, естественно, тоже есть, так как на его основе считается мат. ожидание - насколько выгодно заключать конкретный контракт. Резерв порожняка - тоже деньги. Бизнес волнует рентабельность. Для ряда грузоотправителей по условиям договора и ПСДЦ может оказаться нерентабельным вообще резерв порожняка иметь. Выгодней оформить накладные на подачу только после того, как что-то случилось.
Ad_fesha
26.08.2023 14:53+1Добрый день
У Питона тоже есть автоарима, но скорость расчетов существенно выше (по крайней мере была в 2020-2021, как сейчас - уже не скажу) с разницей примерно +- минута против 1-2 сек у R на 1 ts. Кстати, реализованы алгоритмы по разному, потому как Арима питона, хоть и дольше по времени но скор показывала ощутимо выше чем арима на R
Про prophet кстати упоминания не увидел. Его частенько ставят в пару к ариме - для выбора лучшей модельки
Еще вопрос про масштабируемость и 100500 моделек, подскажите, вариант с "запихнуть все в df и скормить бустингам" наверняка же рассматривали. Подскажите - почему решили, в итоге, в пользу описанного в статье варианта?ptr128 Автор
26.08.2023 14:53У Питона
У Python в PostgrSQL, к сожалению, до сих пор есть неустранимый недостаток - он не может быть безопасным (не иметь прав доступа ни к чему на сервере), как R. Поэтому использовать хранимые процедуры на Python для продуктивного сервера PostgreSQL пока слишком рискованно.
Про prophet кстати упоминания не увидел. Его частенько ставят в пару к ариме - для выбора лучшей модельки
Так auto.arima, и есть обертка над arima, выбирающая лучшую модель с учетом заданных ограничений.
запихнуть все
Как только начинается "запихнуть", то это сразу переливание терабайтов по сети, о чем я писал. Вся прелесть решения как раз в том, что если данные и копируются, то только в оперативке одного сервера.
mikko_kukkanen
Как ни обучай модели, они все равно прогнозируют по принципу "завтра будет так как сегодня или вчера". Это все равно полезно, но можно ли, например, предсказать "черный вторник"?
tessob
Нельзя предсказывать хвосты распределений. Нет такой математики. Хвосты обычно оцениваются через призму рисковых моделей.
mikko_kukkanen
Не очень понял Ваш коммент. Хвост распределения - это ближайшее прошлое или ближайшее будущее (раз его можно или нельзя предсказывать)? А математика есть любая. Вопрос - в надежности предсказания.
tessob
Хвост распределения — это не будущее или прошлое, а редкие события, например вероятность цунами, землетрясения или извержения вулкана. Вот тут «завтра будет как вчера» вообще не работает, т.к. на кануне события тоже будет «как вчера». И математики, которая может сказать, что завтра вероятность такого события больше чем сегодня и чем послезавтра банально нет и не предвидится.
mikko_kukkanen
Я понял, это еще называют "длинные хвосты". А насчет отсутствия на эту математики я бы не был так категоричен. Нет проблем написать модель с "накоплением энергиии".
ptr128 Автор
Сезонность тоже неплохо прогнозируется. То бишь, не только "как сегодня или вчера", но и с учетом, как было в этом периоде в прошлые годы.
Риски же прогнозируются через их математическое ожидание.
mikko_kukkanen
Сезонность "спрогнозировать" просто - зимой холоднее, чем летом, днем теплее, чем ночью, люди больше работают днем, чем ночью (с поправкой на временной пояс), поэтому днем больше транзакций и т.п. И скорее всего, так и будет повторяться, с точность до форс-мажора.
А вот насчет мат ожидания риска можно поподробней? Хотя бы в осях. Как известно, мат ожидания может просто не быть, если соответствующий интеграл расходится.
Хотя по Вашим ответам на комменты вижу, что скорее это - не Ваша тема. У Вас, скорее, обеспечение тех поддержки (в широком смысле).
ptr128 Автор
Так сезонность пропускной способности участков железных дорог Вы не спрогнозируете. Где то выше загрузка летом, где-то - зимой. Ремонты на Кубани начинают уже в марте, а на БАМ - с июня.
Мат. ожидание рассчитывается при кластеризации событий. Модели обучаются не всеми фактам, а только признанными "нормальным" на основании ряда аналитик. "Ненормальные" - разбиваются на кластеры. На основании количества "нормальных" и "ненормальных" событий уже несложно посчитать мат. ожидание того, в какой кластер попадет событие.
Подробней описать не могу, так как получится уже совсем другая тема и объем намного больший этой статьи. Только методика прогнозирования, без уточнения алгоритмов, листов на 100.
mikko_kukkanen
Наверное, все-таки не мат ожидание, а распределение вероятностей? В смысле с какой вероятностью будет попадание в какой кластер.
ptr128 Автор
Так было бы, если бы кластеризация была в разрезе тех же аналитик, по которым формируются разные временные серии. Но это далеко не так.
mikko_kukkanen
В рамках теории нечетких множеств сумма функций принадлежности вроде может быть больше 1 (а также меньше). Но боюсь, что Вы опять скажате, что это не то. Тогда интересно было бы получить ссылку на пример.
ptr128 Автор
В моем случае риски планируется через подачу излишка порожняка на станции и количество вагонов в предварительных заявках операторам. То есть, мне надо посчитать сколько порожняка нужно иметь доступным для подачи и когда. Отсюда и математическое ожидание.
mikko_kukkanen
Конечно, мат ожидание. Только мне кажется, это регрессия, а не кластеризация.
ptr128 Автор
Была бы регрессией, если бы выявлялась только средствами мат. статистики. А так как выявляется, преимущественно, на основе сопутствующих аналитик - то кластеризация.
mikko_kukkanen
Вот уже не хочу спорить о терминах. Как говорится, хоть горшком назови только в печь не ставь. А задача конкретная, респект.