Анализ данных и базоваямодель

Вступление

Эта статья основана на данных конкурса, который компания Driven Data опубликовала для решения проблем с источниками воды в Танзании.

Информация для конкурса была получена Министерством водных ресурсов Танзании с использованием платформы с открытым исходным кодом под названием Taarifa. Танзаниясамая большая страна в Восточной Африке с населением около 60 миллионов человек. Половина населения не имеет доступа к чистой воде, а 2/3 населения страдает от плохой санитарии. В бедных домах семьям часто приходится тратить несколько часов пешком, чтобы набрать воду из водяных насосов.

Для решения проблемы с пресной воды Танзании предоставляются миллиарды долларов иностранной помощи. Однако правительство Танзании не может разрешить эту проблему по сей день. Значительная часть водяных насосов полностью вышла из строя или практически не работает, а остальные требуют капитального ремонта. Министерство водных ресурсов Танзании согласилось с Taarifa, и они запустили конкурс в надежде получить подсказки от сообщества для выполнения стоящих перед ними задач.

Данные

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

Пункты водоснабжения разделены на исправные, нефункциональные и исправные, но нуждающиеся в ремонте. Цель конкурсапостроить модель, прогнозирующую функциональность точек водоснабжения.
Данные содержат 59400 строк и 40 столбцов. Целевая метка содержится в отдельном файле.
Показатель, используемый для этого соревнования,это classification rate, который вычисляет процент строк, в которых прогнозируемый класс совпадает с фактическим классом в тестовом наборе. Максимальное значение равно 1, а минимальное0. Цель состоит в том, чтобы максимизировать classification rate.

Анализ данных

Описания полей в таблице с данными:

• amount_tshобщий статический напор (количество воды, доступное для точки водоснабжения)

• date_recordedдата сбора данных

• funderкто спонсировал постройку колодца

• gps_heightвысота на которой находится колодец

• installerорганизация построившая колодец

• longitudeGPS координаты (долгота)

• latitudeGPS координаты (широта)

• wpt_nameназвание, если оно есть

• num_privateнет информации

• basinгеографический водный бассейн

• subvillageгеографическая локация

• regionгеографическая локация

• region_codeгеографическая локация (код)

• district_code географическая локация (код)

• lgaгеографическая локация

• ward географическая локация

• populationколичество населения около колодца

• public_meetingда/нет

• recorded_by кто собрал данные

• scheme_managementкто управляет колодцем

• scheme_nameкто управляет колодцем

• permitтребуется ли разрешение на доступ

• construction_yearгод постройки

• extraction_typeтип колодца

• extraction_type_groupгруппа типа колодца

• extraction_type_classкласс типа колодца

• managementкак управляется колодец

• management_groupгруппа управления колодца

• paymentстоимость воды

• payment_typeтип стоимости воды

• water_qualityкачество воды

• quality_groupгруппа качества воды

• quantityколичество воды

• quantity_groupгруппа количества воды

• sourceисточник воды

• source_typeтип источника воды

• source_classкласс источника воды

• waterpoint_typeтип колодца

• waterpoint_type_groupгруппа типа колодца

  Прежде всего, давайте посмотрим на целевую меткуклассы не имеют равномерного распределения:

  

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

  • уменьшать количество примеров у преобладающих классов(under-sampling)

  • дублировать строки с метками, которых мало(over-sampling)

  • генерировать синтетические выборкиэто случайная выборка атрибутов из экземпляров в классе меньшинства (SMOTE)

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

  Больше про тактики работы с несбалансированными датасетами можно почитать здесь.

  Посмотрим, как водяные насосы распределяются по территории страны.

  

  Некоторые данные содержат пустые значения.

  

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

  Следующая тепловая карта представляет наличие/отсутствие связи между переменными. Стоит обратить внимание на соотношение между permit, installer и funder.

  

  Давайте посмотрим на общую картину взаимоотношений на дендрограмме.

  

  В характеристиках водяных насосов есть та, которая показывает количество воды. Мы можем проверить, как количество воды связано с состоянием насосов (quantity_group).

  

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

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

  

  К сожалению, этот график не очень информативен, так как количество источников с хорошей водой преобладает. Попробуем сгруппировать только источники с менее качественной водой.

  

  Большинство колодцев с неизвестным качеством из quality_group не работают.

  Есть еще одна интересная характеристика водных точеких тип (waterpoint_type_group).

  

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

  

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

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

  

  Danidaсовместная организация Танзании и Дании по финансированию скважин, и хотя у них много работающих водозаборов, процент неисправных очень высок. Похожая ситуация с RWSSP (программа сельского водоснабжения и канализации), Dhv и некоторыми другими. Следует отметить, что большинство скважин, профинансированных Германией и частными лицами, находятся в рабочем состоянии. Напротив, большое количество скважин, которые финансируются государством, не функционируют. Большинство пунктов водоснабжения, установленных центральным правительством и районным советом, также не работают.

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

  

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

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

  

  Гипотеза полностью подтверждаетсяплата за воду помогает поддерживать источник в рабочем состоянии.

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

  

  Часть данных была заполнена значением 0 вместо реальных данных. Мы также можем видеть, что amount_tsh выше у исправных колодцев (label = 0). Также следует обратить внимание на выбросы в функции amount_tsh. В качестве особенности можно отметить перепад высот и тот факт, что значительная часть населения проживает на высоте 500 метров над средним уровнем моря.

Подготовка данных

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

• Фича installer содержит множество повторов с разными регистрами, орфографическими ошибками и сокращениями. Давайте сначала переведем все в нижний регистр. Затем по простым правилам уменьшаем количество ошибок и делаем группировку.

• После очистки мы заменяем все элементы, которые встречаются менее 71 раз (0,95 квантиля), на other.

• Повторяем по аналогии с фичей funder. Порог отсечки98.

• Данные содержат фичи с очень похожими категориями. Выберем только по одной из них. Поскольку в датасете не так много данных, мы оставляем функцию с наименьшим набором категорий. Удаляем следующие фичи: scheme_management, quantity_group, water_quality, payment_type, extraction_type, waterpoint_type_group, region_code.

• Заменим latitude и longitude значения у выбросов медианным значением для соответсвующего региона из region_code.

• Аналогично повторям для пропущенных значений для subvillage и scheme_name.

• Пропущенные значения в public_meeting и permit заменяем медианным значением.

• Для subvillage, public_meeting, scheme_name, permit, создаем дополнительные категориальные бинарные фичи, которы будут отмечать данные с пропущенными значениями. Так как мы заменяем их на медианные, то для модели оставим информацию, что пропуски были.

• Фичи scheme_management, quantity_group, water_quality, region_code, payment_type, extraction_type, waterpoint_type_group, date_recorded, и recorded_by можно удалить, так как там повторяются данные из других фичей, для модели они будут бесполезны.

Модель

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

Не будем подбирать оптимальные параметры, пусть это будет домашним заданием. Напишем функцию для инициализации и обучения модели.

def fit_model(train_pool, test_pool, **kwargs):    model = CatBoostClassifier(        max_ctr_complexity=5,        task_type='CPU',        iterations=10000,        eval_metric='AUC',        od_type='Iter',        od_wait=500,        **kwargs    )return model.fit(        train_pool,        eval_set=test_pool,        verbose=1000,        plot=False,        use_best_model=True)

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

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

def classification_rate(y, y_pred):    return np.sum(y==y_pred)/len(y)

Поскольку данных мало, разбивать полную выборку на обучающую и проверочную частине лучший вариант. В этом случае лучше использовать методику OOF (Out-of-Fold). Мы не будем использовать сторонние библиотеки; давайте попробуем написать простую функцию. Обратите внимание, что разбиение набора данных на фолды необходимо стратифицировать.

def get_oof(n_folds, x_train, y, x_test, cat_features, seeds):    ntrain = x_train.shape[0]    ntest = x_test.shape[0]              oof_train = np.zeros((len(seeds), ntrain, 3))    oof_test = np.zeros((ntest, 3))    oof_test_skf = np.empty((len(seeds), n_folds, ntest, 3))    test_pool = Pool(data=x_test, cat_features=cat_features)     models = {}    for iseed, seed in enumerate(seeds):        kf = StratifiedKFold(            n_splits=n_folds,            shuffle=True,            random_state=seed)                  for i, (train_index, test_index) in enumerate(kf.split(x_train, y)):            print(f'\nSeed {seed}, Fold {i}')            x_tr = x_train.iloc[train_index, :]            y_tr = y[train_index]            x_te = x_train.iloc[test_index, :]            y_te = y[test_index]            train_pool = Pool(data=x_tr, label=y_tr, cat_features=cat_features)            valid_pool = Pool(data=x_te, label=y_te, cat_features=cat_features)model = fit_model(                train_pool, valid_pool,                loss_function='MultiClass',                random_seed=seed            )            oof_train[iseed, test_index, :] = model.predict_proba(x_te)            oof_test_skf[iseed, i, :, :] = model.predict_proba(x_test)            models[(seed, i)] = modeloof_test[:, :] = oof_test_skf.mean(axis=1).mean(axis=0)    oof_train = oof_train.mean(axis=0)    return oof_train, oof_test, models

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

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

Посмотрев на важность фичей в модели, мы можем убедиться в отсутствии явной утечки информации (лика).

После усреднения прогнозов:

balanced accuracy: 0.6703822994494413classification rate: 0.8198316498316498

Попробуем загрузить результаты на сайт с конкурсом и посмотрим на результат.

Учитывая, что на момент написания этой статьи результат топ-5 был только на 0,005 лучше, мы можем сказать, что базовая модель получилась хорошей.

Проверим, что работы по анализу и очистке данных делались не зря построим модель исключительно на основе исходных данных. Единственное, что мы сделаем, это заполним пропущенные значения нулями.

balanced accuracy: 0.6549535670689709classification rate: 0.8108249158249158

Результат заметно хуже.

Итоги

В статье:

• познакомились с данными и поискали идеи для модели;

• очистили и подготовили данные для модели;

• приняли решение использовать CatBoost, так как основная масса фичей категориальные;

• написали функцию для OOF-предсказания;

• получили отличный результат для базовой модели.

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

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

Код из статьи можно посмотреть здесь.

В этой подборке, переводом которой мы решили поделиться к старту курса о машинном и глубоком обучении, по мнению автора, каждая библиотека заслуживает отдельной статьи. Всё начинается с самого начала: предлагается библиотека, которая сокращает шаблонный код импортирования; заканчивается статья пакетом удобной визуализации данных для исследовательского анализа. Автор также касается работы с картами Google, ускорения и упрощения работы с моделями ML и библиотеки, которая может повысить качество вашего проекта в области обработки естественного языка. Посвящённый подборке блокнот Jupyter вы найдёте в конце.

PyForest

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

Вот почему PyForest это одна из самых удобных библиотек, которые я знаю. С её помощью в ваш блокнот Jupyter можно импортировать более 40 популярнейших библиотек (Pandas, Matplotlib, Seaborn, Tensorflow, Sklearn, NLTK, XGBoost, Plotly, Keras, Numpy и другие) при помощи всего одной строки кода.

Выполните pip install pyforest. Для импорта библиотек в ваш блокнот введите команду from pyforest import *, и можно начинать. Чтобы узнать, какие библиотеки импортированы, выполните lazy_imports().

При этом с библиотеками удобно работать. Технически они импортируются только тогда, когда вы упоминаете их в коде. Если библиотека не упоминается, она не импортируется.

Emot

Эта библиотека может повысить качество вашего проекта по обработке естественного языка. Она преобразует эмотиконы в их описание. Представьте, например, что кто-то оставил в Твиттере сообщение I [здесь в оригинале эмодзи "красное сердце", новый редактор Хабра вырезает его] Python. Человек не написал слово люблю, вместо него вставив эмодзи. Если твит задействован в проекте, придётся удалить эмодзи, а значит, потерять часть информации.

Вот здесь и пригодится пакет emot, преобразующий эмодзи в слова. Для тех, кто не совсем понял, о чём речь, эмотиконы это способ выражения через символы. Например, :) означает улыбку, а :( выражает грусть. Как же работать с библиотекой?

Чтобы установить Emot, выполните команду pip install emot, а затем командой import emot импортируйте её в свой блокнот. Нужно решить, с чем вы хотите работать, то есть с эмотиконами или с эмодзи. В случае эмодзи код будет таким: emot.emoji(your_text). Посмотрим на emot в деле.

Выше видно предложение I [эмодзи "красное сердце"] Python, обёрнутое в метод Emot, чтобы разобраться со значениями. Код выводит словарь со значением, описанием и расположением символов. Как всегда, из словаря можно получить слайс и сосредоточиться на необходимой информации, например, если я напишу ans['mean'], вернётся только описание эмодзи.

Geemap

Говоря коротко, с её помощью можно интерактивно отображать данные Google Earth Engine. Наверное, вы знакомы с Google Earth Engine и всей его мощью, так почему не задействовать его в вашем проекте? За следующие несколько недель я хочу создать проект, раскрывающий всю функциональность пакета geemap, а ниже расскажу, как можно начать с ним работать.

Установите geemap командой pip install geemap из терминала, затем импортируйте в блокнот командой import geemap. Для демонстрации я создам интерактивную карту на основе folium:

import geemap.eefolium as geemapMap = geemap.Map(center=[40,-100], zoom=4)Map

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

Dabl

Позвольте мне рассказать об основах. Dabl создан, чтобы упростить работу с моделями ML для новичков. Чтобы установить её, выполните pip install dabl, импортируйте пакет командой import dabl и можно начинать. Выполните также строчку dabl.clean(data), чтобы получить информацию о признаках, например о том, есть ли какие-то бесполезные признаки. Она также показывает непрерывные, категориальные признаки и признаки с высокой кардинальностью.

Чтобы визуализировать конкретный признак, можно выполнить dabl.plot(data).

Наконец, одной строчкой кода вы можете создать несколько моделей вот так: dabl.AnyClassifier, или так: dabl.Simplefier(), как это делается в scikit-learn. Но на этом шаге придётся предпринять некоторые обычные шаги, такие как создание тренировочного и тестового набора данных, вызов, обучение модели и вывод её прогноза.

# Setting X and y variablesX, y = load_digits(return_X_y=True)# Splitting the dataset into train and test setsX_train, X_test, y_train, y_test = train_test_split(X, y, random_state=1)# Calling the modelsc = dabl.SimpleClassifier().fit(X_train, y_train)# Evaluating accuracy scoreprint(Accuracy score, sc.score(X_test, y_test))

Как видите, Dabl итеративно проходит через множество моделей, включая Dummy Classifier (фиктивный классификатор), GaussianNB (гауссовский наивный Байес), деревья решений различной глубины и логистическую регрессию. В конце библиотека показывает лучшую модель. Все модели отрабатывают примерно за 10 секунд. Круто, правда? Я решил протестировать последнюю модель при помощи scikit-learn, чтобы больше доверять результату:

Я получил точность 0,968 с обычным подходом к прогнозированию и 0,971 с помощью Dabl. Для меня это достаточно близко! Обратите внимание, что я не импортировал модель логистической регрессии из scikit-learn, поскольку это уже сделано через PyForest. Должен признаться, что предпочитаю LazyPredict, но Dabl стоит попробовать.

SweetViz

Это low-code библиотека, которая генерирует прекрасные визуализации, чтобы вывести ваш исследовательский анализ данных на новый уровень при помощи всего двух строк кода. Вывод библиотеки интерактивный файл HTML. Давайте посмотрим на неё в общем и целом. Установить её можно так: pip install sweetviz, а импортировать в блокнот строкой import sweetviz as sv. И вот пример кода:

my_report = sv.analyze(dataframe)my_report.show_html()

Вы видите это? Библиотека создаёт HTML-файл с исследовательским анализом данных на весь набор данных и разбивает его таким образом, что каждый признак вы можете проанализировать отдельно. Возможно также получить численные или категориальные ассоциации с другими признаками; малые, большие и часто встречающиеся значения. Также визуализация изменяется в зависимости от типа данных. При помощи SweetViz можно сделать так много, что я даже напишу о ней отдельный пост, а пока настоятельно рекомендую попробовать её.

Заключение

Все эти библиотеки заслуживают отдельной статьи и того, чтобы вы узнали о них, потому что они превращают сложные задачи в прямолинейно простые. Работая с этими библиотеками, вы сохраняете драгоценное время для действительно важных задач. Я рекомендую попробовать их, а также исследовать не упомянутую здесь функциональность. На Github вы найдёте блокнот Jupyter, который я написал, чтобы посмотреть на эти библиотеки в деле.

Этот материал не только даёт представление о полезных пакетах экосистемы Python, но и напоминает о широте и разнообразии проектов, в которых можно работать на этом языке. Python предельно лаконичен, он позволяет экономить время и в процессе написания кода, выражать идеи максимально быстро и эффективно, то есть беречь силы, чтобы придумывать новые подходы и решения задач, в том числе в области искусственного интеллекта, получить широкое и глубокое представление о котором вы можете на нашем курсе "Machine Learning и Deep Learning".

Узнайте, как прокачаться и в других специальностях или освоить их с нуля:

• Профессия Data Scientist

• Профессия Data Analyst

• Курс по Data Engineering

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

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

В результате родился проект Объясняем код. Посмотреть, что это такое можно на code-explained.com. Код проекта выложен на Гитхаб.

Чем я вдохновлялся

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

Например, если вы редактируете программный код, то результат должен отображаться в реальном времени, без повторной компиляции. Иначе, вы не увидите промежуточные состояния, и у вас не родятся связанные с ними идеи. Я не буду подробно пересказывать эту мысль, тем более, что на Хабре уже подробно о ней рассказали. Почитайте, если еще не успели.

Идея Брета Виктора о мгновенной обратной связи сильно меня вдохновила. Я подумал, что она могла бы изменить подход к изучению программирования, в частности алгоритмов и структур данных. Потому что сегодня этот процесс ужасно технозависим.

Как сегодня изучают алгоритмы

Преподаватель записывает на доске код, а потом по шагам объясняет, как он работает.

В видео на Youtube происходит нечто подобное ведущий берет код алгоритма и отрисовывает этапы его работы

На ИТ-ресурсах создают анимации.

А кто-то даже пытается объяснить работу алгоритма через танец.

Почему эти подходы казались мне неэффективными

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

Как я перешел от технозависимости к человечности

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

Вместе с командой мы разработали проект Объясняем код. Это такой интерактивный учебник, который с помощью визуализаций объясняет, как работают алгоритмы.

Выглядит это как-то так. Можно вводить свои данные. А затем с помощью плеера запускать и останавливать исполнение кода, перематывать его вперед и назад. При этом визуализация в реальном времени показывает, что делает алгоритм.

Как это технически реализовано

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

Чтобы сэкономить на копировании данных, я применил Immutable.js. Immutable.js это библиотека неизменяемых коллекций. Модификации таких коллекций всегда возвращают новую коллекцию, которую легко сохранить. Так я избегаю глубокого копирования. Не буду вдаваться в подробности, на Хабре про immutable.js уже писали. Для примера кусочек кода с итерированием по хеш-таблице:

while (true) { // Итерируемся по списку    this.addBP('check-not-found'); // Метод сохраняем состояние    if (this.newList.get(this.newListIdx) === null) {        // this.newList -- это немутабельный список        break;    }    this.addBP('check-found'); // Выполнена очередная строчка, сохраняем состояние    if (EQ(this.newList.get(this.newListIdx), this.number)) {        this.addBP('found-key');        return true;    }    this.fmtCollisionCount += 1; // Для динамических комментариев иногда нужно сохранять статистикуу    this.newListIdx = (this.newListIdx + 1) % this.newList.size; // Переходим к следующему индекксу    this.addBP('next-idx');}

В любой момент пользователь может отмотать время слайдером, и мы должны быть готовы анимировать переход. Поэтому анимации сделаны декларативно на React и CSS. На React описываются нужное состояние, а CSS transitions автоматически высчитывают и анимируют переход. Это позволяет избежать совсем сложного императивного кода с запутанным состоянием и ручными вычислениями. Получается тоже не просто, и в итоге анимациями заведует огромный класс на 1000 строк. Но я уверен, что императивно реализовать то же самое у меня вообще не получилось бы.

Много строк возникает из-за особенностей браузерных API и производительности в разных браузерах. Например, большой проблемой оказалось сделать так, чтобы браузеры не склеивали последовательные изменения друг с другом. Если добавить div с определённой начальной позицией, и потом сразу же поменять координаты на конечные, то браузер склеит эти два изменения в одно. Div сразу окажется в конечной позиции без анимации. Чтобы такое не происходило, приходится вставлять задержку в два фрейма анимации с помощью window.requestAnimationFrame().

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

Код проекта на гитхабе

Что дальше?

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

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

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

И нетрудно догадаться почему. Гистограммы весьма интуитивно наглядны: любой поймет их с первого взгляда. Более того, они объективно представляют реальность, не так ли? А вот и нет.

Гистограмма может ввести в заблуждение и привести к ошибочным выводам даже на простейшем наборе данных!

В этой статье мы на примерах рассмотрим 6 причин, почему, когда дело доходит до визуализации данных, гистограммы точно не является лучшим выбором:

1. Они слишком сильно зависят от количества интервалов.

2. Они слишком сильно зависят от максимума и минимума переменной.

3. Они не дают возможности заметить значимые значения переменной.

4. Они не позволяют отличить непрерывные переменные от дискретных.

5. Они делают сравнение распределений сложным.

6. Их построение затруднено, если в памяти находятся не все данные.

Ладно, я понял: гистограммы не идеальны. Но есть ли у меня выбор? Конечно есть!

В конце статьи я порекомендую другой график, называемый CDP, который минует эти недостатки.

Итак, что же не так с гистограммой?

1. Она слишком сильно зависит от количества интервалов.

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

Переменная представляет собой максимальную частоту сердечных сокращений (ударов в минуту), полученную у 303 людей во время некоторой физической активности (данные взяты из набора данных UCI по сердечным заболеваниям: источник).

Глядя на верхний левый график (который мы получим по умолчанию в Python и R), у нас сложится впечатление хорошего распределения с одним пиком (модой). Однако если бы мы рассмотрели бы другие варианты гистограммы, мы получили бы совершенно другую картину. Разные гистограммы одних и тех же данных могут привести к противоречивым выводам.

2. Она слишком сильно зависит от максимума и минимума переменной.

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

Например, давайте попробуем изменить максимум переменной, не меняя количество интервалов.

Отличается только одно значение, а весь график получается другим. Это нежелательное свойство, потому что нас интересует общее распределение: одно значение не должно так влиять на график!

3. Не дает возможности заметить значимые значения переменной.

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

Классическим примером является случай, когда отсутствующим значениям массово присваивается 0. В качестве примера давайте рассмотрим набор данных переменной, состоящий из 10 тысяч значений, 26% из которых нули.

График слева это то, что вы получаете по умолчанию в Python. Глядя на него, вы не заметите скопление нулей, и вы даже можете подумать, что эта переменная имеет плавную динамику.

График справа получен путем сужения интервалов и дает более четкое представление о реальности. Но дело в том, что как бы вы ни сужали интервалы, вы никогда не будете уверены, содержит ли первый интервал только 0 или какие-то другие значения.

4. Не позволяет отличить непрерывные переменные от дискретных.

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

Возьмем переменную Возраст (Age). Вы можете получить Возраст = 49 лет (когда возраст округлен) или Возраст = 49,828884325804246 лет (когда возраст рассчитывается как количество дней с момента рождения, деленное на 365,25). Первая дискретная переменная, вторая непрерывная.

Тот, что слева, непрерывен, а тот, что справа, дискретен. Однако на верхних графиках (по умолчанию в Python) вы не увидите никакой разницы между ними: они выглядят совершенно одинаково.

5. Сложно сравнивать распределения.

Часто бывает необходимо сравнить одну и ту же переменную в разных кластерах. Например, в отношении данных UCI о сердечных заболеваниях, приведенных выше, мы можем сравнить:

• все население (для справки)

• люди моложе 50 страдающие сердечными заболеваниями

• люди моложе 50 НЕ страдающие сердечными заболеваниями

• люди старше 60 лет страдающие сердечными заболеваниями

• люди старше 60 и НЕ страдающие сердечными заболеваниями.

Вот что мы получили бы в итоге:

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

6. Сложно построить, если в памяти находятся не все данные.

Если все ваши данные находятся в Excel, R или Python, построить гистограмму легко: в Excel вам просто нужно кликнуть по иконке гистограммы, в R выполнить команду hist(x), а в Python plt.hist(х).

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

| INTERVAL_LEFT | INTERVAL_RIGHT | COUNT |

|---------------|----------------|---------------|

| 75.0 | 87.0 | 31 |

| 87.0 | 99.0 | 52 |

| 99.0 | 111.0 | 76 |

| ... | ... | ... |

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

WITHSTATS AS (  SELECT     COUNT(*) AS N,    APPROX_QUANTILES(VARIABLE_NAME, 4) AS QUARTILES  FROM    TABLE_NAME),BIN_WIDTH AS (  SELECT    -- freedman-diaconis formula for calculating the bin width    (QUARTILES[OFFSET(4)]  QUARTILES[OFFSET(0)]) / ROUND((QUARTILES[OFFSET(4)]  QUARTILES[OFFSET(0)]) / (2 * (QUARTILES[OFFSET(3)]  QUARTILES[OFFSET(1)]) / POW(N, 1/3)) + .5) AS FD  FROM     STATS),HIST AS (  SELECT     FLOOR((TABLE_NAME.VARIABLE_NAME  STATS.QUARTILES[OFFSET(0)]) / BIN_WIDTH.FD) AS INTERVAL_ID,    COUNT(*) AS COUNT  FROM     TABLE_NAME,    STATS,    BIN_WIDTH  GROUP BY     1)SELECT   STATS.QUARTILES[OFFSET(0)] + BIN_WIDTH.FD * HIST.INTERVAL_ID AS INTERVAL_LEFT,  STATS.QUARTILES[OFFSET(0)] + BIN_WIDTH.FD * (HIST.INTERVAL_ID + 1) AS INTERVAL_RIGHT,  HIST.COUNTFROM   HIST,   STATS,   BIN_WIDTH

Немного громоздко, не правда ли?

Альтернатива: график кумулятивного распределения.

Узнав 6 причин, по которым гистограмма не является идеальным выбором, возникает естественный вопрос: Есть ли у меня альтернатива? Хорошие новости: существует лучшая альтернатива, которая называется График кумулятивного распределения (Cumulative Distribution Plot - CDP). Я знаю, что это название не такое запоминающееся, но гарантирую, оно того стоит.

График кумулятивного распределения это график квантилей переменной. Другими словами, каждая точка CDP показывает:

• по оси x: исходное значение переменной (как в гистограмме);

• по оси y: сколько наблюдений имеют такое же или меньшее значение.

Давайте посмотрим на пример с переменной максимальной частотой пульса.

Возьмем точку с координатами x = 140 и y = 90 (30%). По горизонтальной оси вы видите значение переменной: 140 ударов сердца в минуту. По вертикальной оси вы видите количество наблюдений, у которых частота сердцебиение равна или ниже 140 (в данном случае 90 человек, что означает 30% выборки). Следовательно, у 30% нашей выборки максимальная частота сердцебиения составляет 140 или менее ударов в минуту.

Какой смысл в графике, показывающем, сколько наблюдений равно или ниже заданного уровня? Почему не просто равно? Потому что в противном случае результат зависел бы от отдельных значений переменной. И это не сработает, потому что каждое значение имеет очень мало наблюдений (обычно только одно, если переменная непрерывна). Напротив, CDP полагаются на квантили, которые более стабильны, выразительны и легко читаются.

Вдобавок CDP намного полезнее. Если задуматься, вам часто приходится отвечать на такие вопросы, как у скольких из них от 140 до 160? Или у скольких из них больше 180?. Имея перед глазами CDP, вы можете дать немедленный ответ. С гистограммой это было бы невозможно.

CDP решает все проблемы, которые мы видели выше. Фактически, по сравнению с гистограммой:

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

2. Не страдает от выпадающих значений. Экстремальные значения не влияют на CDP, поскольку квантили не меняются.

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

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

5. Упрощает сравнение распределений. На одном графике легко сравнить два или более распределения, поскольку это просто кривые, а не области. Кроме того, ось y всегда находится в диапазоне от 0 до 100%, что делает сравнение еще более простым. Для сравнения, это пример, который мы видели выше:

6. Его легко построить, даже если у вас нет всех данных в памяти. Все, что вам нужно, это квантили, которые можно легко получить с помощью SQL:

SELECT   COUNT(*) AS N,  APPROX_QUANTILES(VARIABLE_NAME, 100) AS PERCENTILESFROM  TABLE_NAME

Как построить график кумулятивного распределения в Excel, R, Python

В Excel вам нужно построить два столбца. Первый с 101 числом, равномерно распределенными от 0 до 1. Второй столбец должен содержать процентили, которые могут быть получены по формуле: =PERCENTILE(DATA, FRAC), где DATA - это вектор, содержащий данные, а FRAC - это первый столбец: 0,00, 0,01, 0,02, 0,03,, 0,98, 0,99, 1. Затем вам просто нужно построить график по этим двум столбцам, разместив значения переменной на оси x.

В R это делается в одну строчку:

plot(ecdf(data))

В Python:

from statsmodels.distributions.empirical_distribution import ECDFimport matplotlib.pyplot as pltecdf = ECDF(data)plt.plot(ecdf.x, ecdf.y)

Спасибо за внимание! Надеюсь, эта статья оказалась для вас полезной.

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

Перевод материала подготовлен в рамках онлайн-курса "Machine Learning. Basic". Всех заинтересованных приглашаем на день открытых дверей курса, где можно будет узнать все подробности об обучении и пообщаться с преподавателем.

- Узнать подробнее о курсе "Machine Learning. Basic"

- Смотреть онлайн-встречу "День открытых дверей"

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

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

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

• "Теория вероятности и математическая статистика" Л. Н. Фадеева, А. В. Лебедев;

• "Теория вероятности и математическая статистика" Н. Ш. Кремер. Эти две книги на, мой взгляд, лучшие для первого знакомства с "тервером и статами". Вовсе не рекомендую вгрызаться в эти книги, продолжайте читать даже если что-то непонятно, просто для того что бы оценить "широту картинки" и привыкнуть к терминологии. Конечно, советовать продолжать читать, даже если что-то непонятно - это не комильфо. Но есть еще одна хорошая книга:

• "Теория вероятности и математическая статистика" В. П. Лисьев. В ней гораздо больше примеров и разъяснений. Так что если первые две книги не "побеждаются", то можете смело браться за Лисьева. Однако в первых двух книгах очень хорошее введение в теорию вероятности.

Ну и наконец книга "Статистика для всех" Сара Бослаф. Я бы порекомендовал для прочтения только эту книгу, но мне кажется, что если бы я ее читал, абсолютно не зная теорию вероятности и мат. статистику, то я бы очень много вообще не понял, ну или очень долго вникал в текст, формулы и картинки. Там конечно есть разделы с основами для конкретных новичков, но на мне такое не работает. Хотя, более чем вероятно, что этой книги для вас этого окажется более чем достаточно. В книге намеренно нет примеров с кодом, но это скорее плюс чем минус, потому что в книге есть примеры численных расчетов. Так что можно кодить и сравнивать свои результаты с результатами из книги.

Конечно в интернете море всякой литературы, но как пишет Сара в своей книге "Вода, вода, кругом вода, а мы не пьем", намекая на то, что книг море, но самоучки умирают от "жажды".

Вместо того что бы развозить математику и доказательства (кто-то скажет, что это кощунство) я предлагаю немного покодить и пойти по принципу "если что-то непонятно то моделируй и экспериментируй до тех пор пока не станет понятно". Само программирование потребует от нас хотя бы небольших познаний библиотек NumPy, matplotlib, seaborn и SciPy - открытых вкладок с документацией будет более чем достаточно. Итак, начнем.

Пример не для тех кто проголодался

Давайте представим... так постойте, прежде чем что-то представлять давайте сделаем все необходимые импорты и настройки:

import numpy as npfrom scipy import statsimport matplotlib.pyplot as plt#import seaborn as snsfrom pylab import rcParams#sns.set()rcParams['figure.figsize'] = 10, 6%config InlineBackend.figure_format = 'svg'np.random.seed(42)

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

norm_rv = stats.norm(loc=30, scale=5)samples = np.trunc(norm_rv.rvs(365))samples[:30]

array([32., 29., 33., 37., 28., 28., 37., 33., 27., 32., 27., 27., 31.,       20., 21., 27., 24., 31., 25., 22., 37., 28., 30., 22., 27., 30.,       24., 31., 26., 28.])

Конечно же теперь нам интересно выяснить среднее время доставки пиццы и его среднеквадратическое отклонение:

samples.mean(), samples.std()

(29.52054794520548, 4.77410133275075)

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

А еще, было бы интересно посмотреть на то, как распределены данные:

sns.histplot(x=samples, discrete=True);

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

norm_rv = stats.norm(loc=30, scale=5)beta_rv = stats.beta(a=5, b=5, loc=14, scale=32)gamma_rv = stats.gamma(a = 20, loc = 7, scale=1.2)tri_rv = stats.triang(c=0.5, loc=17, scale=26)x = np.linspace(10, 50, 300)sns.lineplot(x = x, y = norm_rv.pdf(x), color='r', label='norm')sns.lineplot(x = x, y = beta_rv.pdf(x), color='g', label='beta')sns.lineplot(x = x, y = gamma_rv.pdf(x), color='k', label='gamma')sns.lineplot(x = x, y = tri_rv.pdf(x), color='b', label='triang')sns.histplot(x=samples, discrete=True, stat='probability',             alpha=0.2);

Распределение точно нормальное?

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

• на перекрестке пришлось ждать две минуты пока светофор загорится зеленым;

• ударился ногой и из-за хромоты шел дольше обычного;

• доставщик оказался скейтбордистом и передвигался быстрее обычного;

• дорогу перебежала черная кошка и пришлось идти другим более долгим путем;

• развязались шнурки и пришлось тратить время на их завязывание;

• развязались шнурки и доставщик упал, поэтому пришлось тратить время на отряхивание грязи и завязывание шнурков.

Конечно, мы можем придумывать очень много таких событий, вплоть до самых невероятных (возможно одна из дорог в Ад вымощена не доставленными пицами, поэтому мы не будем выдумывать что-то обидное. Ок?). Тем не менее, для нас важно что бы эти события описывались такими переменными, значения которых равновероятны, т.е. распределены по равномерному закону. В качестве примера, можно придумать переменную , которая будет описывать время ожидания зеленого света светофора на перекрестке. Если это время заключено в промежутке от нуля до четырех минут, то сегодня это время может составлять:

unif_rv = stats.uniform(loc=0, scale=4)unif_rv.rvs()

0.8943833540778106

А завтра, после-завтра и после-после-завтра это время может быть равно:

unif_rv.rvs(size=3)

array([3.85289016, 0.0486179 , 3.87951531])

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

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

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

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

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

Но если значения каждой из переменных распределены равномерно, то как будет распределена их сумма - переменная ? Давайте сгенерируем 10 тысяч таких сумм и посмотрим на гистограмму:

Y_samples = [unif_rv.rvs(size=15).sum() for i in range(10000)]sns.histplot(x=Y_samples);

Не напоминает ли вам эта картинка тот самый колокол, того самого нормального распределения? Если напоминает, то вы только, что поняли смысл центральной предельной теоремы: распределение суммы случайных переменных стремится к нормальному распределению при увеличении количества слагаемых в этой сумме.

Конечно, пример с доставкой пиццы не совсем корректен для демонстрации предельной теоремы, потому что все события которые мы придумали носят условный характер:

• на перекрестке пришлось ждать две минуты пока светофор загорится зеленым. Но ждать придется только если мы подошли к светофору, который уже горит красным;

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

• доставщик оказался скейтбордистом и передвигался быстрее обычного. При условии что он не забыл скейтборд дома;

• дорогу перебежала черная кошка и пришлось идти другим более долгим путем. Существует ненулевая вероятность того, что доставщик не является фанатом бабы Нины;

• развязались шнурки и пришлось тратить время на их завязывание. Как часто развязываются шнурки?;

• развязались шнурки и доставщик упал, поэтому пришлось тратить время на отряхивание грязи и завязывание шнурков. Если дворник не халтурит, то и ни от какой грязи отряхиваться не надо.

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

А завтра это время может задаваться как:

Будет ли теперь распределена нормально? Учитывая что сумма нормально распределенных величин тоже имеет нормальное распределение, то можно дать утвердительный ответ. Именно поэтому, когда мы взглянули на распределение 365-и значений времени доставки, мы практически сразу решили, что перед нами нормальное распределение, даже несмотря на то что оно вовсе не похоже на идеальный колокол.

Z-значения

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

Почему наш сосед так уверен в долгой доставке? И вообще, оправдана ли его уверенность? Очевидно он, как и некоторые люди, думает, что минут означает, что доставка может длиться 27, 31, даже 35 минут, но никак не 23 или 38 минут. Однако, мы заказывали пицу 365 раз и знаем, что доставка может длиться и 20 и даже 45 минут. А фраза минут, означает лишь то, что какая-то значительная часть доставок бодет занимать от 25 до 35 минут. Зная параметры распределения, мы даже можем смоделировать несколько тысяч доставок и прикинуть величину этой части:

N = 5000t_data = norm_rv.rvs(N)t_data[(25 < t_data) & (t_data < 35)].size/N

0.7008

Где-то две трети значений укладываются в интервал от 25 до 35 минут. А сколько значений будет превосходить 40 минут?

t_data[t_data > 40].size/N

0.0248

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

0.02**3

8.000000000000001e-06

Конечно, компьютерное моделирование - это очень хорошо. Но в данном случае лучше воспользоваться так называемыми Z-значениями.

Вычислить Z-значение можно по следующей формуле:

Где - это время доставки, т.е. какое-то конкретное значение случайной переменной , которая имеет нормальное распределение, а и это мат. ожидание и среднеквадратическое отклонение, т.е. параметры распределения, в нашем случае они равны 30 и 5 минут соответственно. Давайте рассчитаем Z-значение для сорока минут:

Что мы сейчас сделали? В числителе мы вычислили на какую величину наше время доставки отличается от среднего времени доставки, а далее мы просто поделили это значение на стандартное отклонение времени доставки. Но как интерпретировать данный результат и зачем вообще использовать Z-значение? Что бы понять это придется немного "порисовать":

fig, ax = plt.subplots()x = np.linspace(norm_rv.ppf(0.001), norm_rv.ppf(0.999), 200)ax.vlines(40.5, 0, 0.1, color='k', lw=2)sns.lineplot(x=x, y=norm_rv.pdf(x), color='r', lw=3)sns.histplot(x=samples, stat='probability', discrete=True);

На данном графике мы изобразили гистограмму наших данных samples, но теперь высота прямоугольников показывает не количество вхождений каждого значения в выборку, а вероятность их появления в выборке. Красной линией мы нарисовали функцию распределения плотности вероятности значений времени доставки. Выше мы пытались экспериментальным путем определить какую долю составляют значения времени больше 40 минут. И данный график должен натолкнуть нас на мысль о том, что есть два способа сделать это. Первый - экспериментальный, т.е. мы моделируем, скажем, 5000 доставок, строим гистограмму и находим сумму высот прямоугольников, расположенных правее черной линии:

np.random.seed(42)N = 10000values = np.trunc(norm_rv.rvs(N))fig, ax = plt.subplots()v_le_41 = np.histogram(values, np.arange(9.5, 41.5))v_ge_40 = np.histogram(values, np.arange(40.5, 51.5))ax.bar(np.arange(10, 41), v_le_41[0]/N, color='0.8')ax.bar(np.arange(41, 51), v_ge_40[0]/N)p = np.sum(v_ge_40[0]/N)ax.set_title('P(Y>40min) = {:.3}'.format(p))ax.vlines(40.5, 0, 0.08, color='k', lw=1);

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

fig, ax = plt.subplots()x = np.linspace(norm_rv.ppf(0.0001), norm_rv.ppf(0.9999), 300)ax.plot(x, norm_rv.pdf(x))ax.fill_between(x[x>41], norm_rv.pdf(x[x>41]), np.zeros(len(x[x>41])))p = 1 - norm_rv.cdf(41)ax.set_title('P(Y>40min) = {:.3}'.format(p))ax.hlines(0, 10, 50, lw=1, color='k')ax.vlines(41, 0, 0.08, color='k', lw=1);

Значения вероятности в первом и втором случае получились довольно близкими. Но надо оговориться, что в первом случае мы использовали дискретную случайную величину, а во втором - непрерывную, поэтому код выглядит несколько "читерским", как буд-то бы все намеренно подгонялось так, что бы два результата были похожи. Конечно, возможно я что-то перемудрил, но вроде бы все верно. В первом случае мы считаем сколько раз доставка составляла более 40 минут, т.е. считали сколько доставок длились 41, 42, 43 и т.д. минут. Во втором случае мы просто считаем сколько непрерывных величин оказалось правее значения 41.0. В принципе, модуль SciPy позволяет создавать собственные распределения случайных переменных, в руководстве даже имеется пример, как создать дискретно-нормальное распределение. Но если задуматься, то наличие дискретных значений и использование непрерывных функций распределения, это своего рода - неизбежность, ведь подавляющее большинство измерений носит дискретный характер: рост, вес, доход и т.д. Так что далее мы не будем больше оговариваться по этому поводу, а лучше вернемся к Z-значениям.

Итак, допустим мы оказались в средиземье и каким-то образом выяснили что рост хобитов и гномов в сантиметрах распределен как и . Если рост Фродо равен 99 сантиметрам а рост Гимли 143 сантиметра, то как понять чей рост более типичен среди своих народов? Что бы выяснить это мы можем изобразить функцию распределения плотности вероятности для каждого народа с отмеченными значениями, а заодно определить долю тех, кто превышает эти значения:

fig, ax = plt.subplots(nrows=1, ncols=2, figsize = (12, 5))nrv_hobbit = stats.norm(91, 8)nrv_gnome = stats.norm(134, 6)for i, (func, h) in enumerate(zip((nrv_hobbit, nrv_gnome), (99, 143))):    x = np.linspace(func.ppf(0.0001), func.ppf(0.9999), 300)    ax[i].plot(x, func.pdf(x))    ax[i].fill_between(x[x>h], func.pdf(x[x>h]), np.zeros(len(x[x>h])))    p = 1 - func.cdf(h)    ax[i].set_title('P(H > {} см) = {:.3}'.format(h, p))    ax[i].hlines(0, func.ppf(0.0001), func.ppf(0.9999), lw=1, color='k')    ax[i].vlines(h, 0, func.pdf(h), color='k', lw=2)    ax[i].set_ylim(0, 0.07)fig.suptitle('Хобиты слева, гномы справа', fontsize = 20);

Эти графики, конечно не обладают свойством самоочевидности, но в принципе, можно сказать (и наверняка ошибиться), что рост Фродо несколько ближе к вершине распределения чем рост Гимли. А это значит, что вероятность встретить хобита с таким же ростом как у Фродо несколько больше вероятности встретить гнома с ростом как у Гимли. Именно это и понимается под словом "типичность".

Выполнить сравнение типичности гораздо проще и нагляднее если воспользоваться Z-значениями:

fig, ax = plt.subplots()N_rv = stats.norm()x = np.linspace(N_rv.ppf(0.0001), N_rv.ppf(0.9999), 300)ax.plot(x, N_rv.pdf(x))ax.hlines(0, -4, 4, lw=1, color='k')ax.vlines(1, 0, 0.4, color='r', lw=2, label='Фродо')ax.vlines(1.5, 0, 0.4, color='g', lw=2, label='Гимли')ax.legend();

Огромным преимуществом Z-значений является то, что они "стандартизированы", т.е. преобразованы так словно они взяты из стандартного нормального распределения , именно поэтому два Z-значения нарисованы на фоне единственной кривой. Однако, в общем случае, даже само рисование графиков вовсе не обязательно, потому что меньшие по модулю Z-значения обладают большей частотой появления. Одной записи:

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

Сравнение Z-значений нескольких величин из разных "нормальных" выборок с разными параметрами распределения возможно потому, что сами Z-значения измеряются в сигмах. Это становится более очевидным если еще раз взглянуть на знаменатель формулы:

Не важно, что вы пытаетесь сравнить: рост, вес, литры или доход; в чем бы не измерялись сравниваемые величины, после вычисления Z-значений они будут измеряться в сигмах. Чем меньше по модулю Z-значение, тем ближе оно к вершине распределения, а знак Z-значения укажет по какую сторону от вершины это значение находится.

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

Z-статистика

Теперь давайте снова вернемся к нашему соседу, который возмущен слишком долго доставкой пицы. Выше мы вычислили Z-значение для 40 минут:

Пока у нас нет опыта, достаточного для того что бы сразу понять много это или мало, лучше изображать эти значения графически:

fig, ax = plt.subplots()N_rv = stats.norm()x = np.linspace(N_rv.ppf(0.0001), N_rv.ppf(0.9999), 300)ax.plot(x, N_rv.pdf(x))ax.hlines(0, -4, 4, lw=1, color='k')ax.vlines(2, 0, 0.4, color='g', lw=2);

Это значение находится справа от вершины на расстоянии , что довольно много. Однако, сосед утверждает, что ждал пицу больше сорока минут три дня подряд. Возможно он и не вычислял среднее время ожидания своей пицы, но три значения - это уже статистика!

Характеристики генеральной совокупности называют параметрами, а характеристики выборки - статистиками. Замеряя время доставки на протяжении 365-и дней, мы сделали вывод о параметрах генеральной совокупности, т.е. времени всех возможных доставок пицы, решив что эти значения берутся из . А зная распределение мы можем поэкспериментировать. Например, наш сосед сделал всего три заказа, и по его ощущениям среднее время доставки, было больше сорока минут. А что если мы повторим этот эксперимент 5000 раз и посмотрим на то как распределено среднее трех заказов:

sns.histplot(np.trunc(norm_rv.rvs(size=(N, 3))).mean(axis=1),             bins=np.arange(19, 42));

Судя по графику, получить среднее значение трех заказов большее 40 минут - практически нереально. И вот тут возможны два предположения:

• либо сосед врет;

• либо что-то с генеральной совокупностью, т.е. по какой-то причине, доставка и правда длится дольше обычного.

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

1 - N_rv.cdf(2)

0.02275013194817921

Поэтому не удивительно, что вероятность получить среднее время трех доставок большее 40 минут исчезающе мала:

(1 - N_rv.cdf(2))**3

1.1774751750476762e-05

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

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

Давайте предположим что мы сделали три заказа и среднее значение оказалось равным 35 минутам, тогда Z-статистика будет вычисляться так:

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

который в нашем случае выглядит как [25; 35] минут. Как и ранее мы можем найти данную вероятность с помощью моделирования:

N = 10000means = np.trunc(norm_rv.rvs(size=(N, 3))).mean(axis=1)means[(means>=25)&(means<=35)].size/N

0.9241

N = 10000fig, ax = plt.subplots()means = np.trunc(norm_rv.rvs(size=(N, 3))).mean(axis=1)h = np.histogram(means, np.arange(19, 41))ax.bar(np.arange(20, 25), h[0][0:5]/N, color='0.5')ax.bar(np.arange(25, 36), h[0][5:16]/N)ax.bar(np.arange(36, 41), h[0][16:22]/N, color='0.5')p = np.sum(h[0][6:16]/N)ax.set_title('P(25min < Y < 35min) = {:.3}'.format(p))ax.vlines([24.5 ,35.5], 0, 0.08, color='k', lw=1);

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

x, n, mu, sigma = 35, 3, 30, 5z = abs((x - mu)/(sigma/n**0.5))N_rv = stats.norm()fig, ax = plt.subplots()x = np.linspace(N_rv.ppf(0.0001), N_rv.ppf(0.9999), 300)ax.plot(x, N_rv.pdf(x))ax.hlines(0, x.min(), x.max(), lw=1, color='k')ax.vlines([-z, z], 0, 0.4, color='g', lw=2)x_z = x[(x>-z) & (x<z)] # & (x<z)ax.fill_between(x_z, N_rv.pdf(x_z), np.zeros(len(x_z)), alpha=0.3)p = N_rv.cdf(z) - N_rv.cdf(-z)ax.set_title('P({:.3} < z < {:.3}) = {:.3}'.format(-z, z, p));

Обратите внимание, на то что Z-статистика зависит как от среднего выборки - , так и от величины выборки - . Если мы закажем пицу 5, 30 или 100 раз, то какова вероятность того, что среднее значение времени доставок окажется в интервале [29; 31]? Давайте взглянем:

fig, ax = plt.subplots(nrows=1, ncols=3, figsize = (12, 4))for i, n in enumerate([5, 30, 100]):    x, mu, sigma = 31, 30, 5    z = abs((x - mu)/(sigma/n**0.5))    N_rv = stats.norm()    x = np.linspace(N_rv.ppf(0.0001), N_rv.ppf(0.9999), 300)    ax[i].plot(x, N_rv.pdf(x))    ax[i].hlines(0, x.min(), x.max(), lw=1, color='k')    ax[i].vlines([-z, z], 0, 0.4, color='g', lw=2)    x_z = x[(x>-z) & (x<z)] # & (x<z)    ax[i].fill_between(x_z, N_rv.pdf(x_z), np.zeros(len(x_z)), alpha=0.3)    p = N_rv.cdf(z) - N_rv.cdf(-z)    ax[i].set_title('n = {}, z = {:.3}, p = {:.3}'.format(n, z, p));fig.suptitle(r'Z-статистики для $\bar{x} = 31$', fontsize = 20, y=1.02);

При 5 заказах среднее выборки попадает в интервал [29;31] скорее случайно, чем систематически. При 30 заказх около четверти средних значений так и не войдут в заданный интервал. И только при сотне заказов мы можем быть более-менее уверены в том что отклонение среднего выборки от среднего генеральной совокупности не будет больше 1 минуты.

С другой стороны, мы можем рассуждать и по другому: если среднее генеральной совокупности равно 30 минутам, то какова вероятность получить среднее выборки равное 31 минуте если мы сделаем 5, 30 или 100 заказов? Очевидно, что при n=5 среднее выборки, может отклоняться очень сильно, следовательно, вероятность получить минуте очень высока. Но при среднее выборки практически не отклоняется от среднего генеральной совокупности, поэтому получить случайным образом при практически невозможно. Что это значит? А это значит, что если мы сделали 100 заказов и получили среднее время доставки равное 31 минуте, то скорее всего мы ошибаемся насчет того, что среднее генеральной совокупности равно 30 минутам.

Но как быть с нашим соседом, ведь он сделал всего 3 заказа? Неужели он действительно прав? Судя по всему - да. Даже если среднее время ожидания составило 40 минут, то Z-статистика будет равна 3.81, а площадь под кривой будет практически равна 1:

x, n, mu, sigma = 41, 3, 30, 5z = abs((x - mu)/(sigma/3**0.5))N_rv = stats.norm()fig, ax = plt.subplots()x = np.linspace(N_rv.ppf(1e-5), N_rv.ppf(1-1e-5), 300)ax.plot(x, N_rv.pdf(x))ax.hlines(0, x.min(), x.max(), lw=1, color='k')ax.vlines([-z, z], 0, 0.4, color='g', lw=2)x_z = x[(x>-z) & (x<z)] # & (x<z)ax.fill_between(x_z, N_rv.pdf(x_z), np.zeros(len(x_z)), alpha=0.3)p = N_rv.cdf(z) - N_rv.cdf(-z)ax.set_title('P({:.3} < z < {:.3}) = {:.3}'.format(-z, z, p));

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

• либо наш сосед просто чрезвычайно "везучий" человек;

• либо доставка пицы и вправду теперь выполняется дольше обычного.

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

p-value

Мы смогли убедиться в том что Z-статистика позволяет оценить вероятность того, что среднее выборки размером , взятой из генеральной совокупности попадет в заданный интервал значений. Это удобно тем, что позволяет сделать вывод о случайности полученного . Чем меньше модуль значения Z-статистики, тем меньше достоверность среднего. Например выше мы видели, что вероятность попадания при в интервал [29;31] составляет всего около 0.35. В то время как вероятность непопадания в заданный интервал равна 10.35=0.65. Поэтому мы и сделали вывод о том, что значение при скорее обусловлено случайностью, чем какими-то объективными причинами.

Фактически, значение вероятности 0.65, полученное выше - это и есть то самое p-value которое как раз и показывает вероятность случайного выхода за границы интервала, задаваемого значением Z-статистики, что можно изобразить следующим образом:

x, n, mu, sigma = 31, 5, 30, 5z = abs((x - mu)/(sigma/n**0.5))N_rv = stats.norm()fig, ax = plt.subplots()x = np.linspace(N_rv.ppf(0.0001), N_rv.ppf(0.9999), 300)ax.plot(x, N_rv.pdf(x))ax.hlines(0, x.min(), x.max(), lw=1, color='k')ax.vlines([-z, z], 0, 0.4, color='g', lw=2)x_le_z, x_ge_z = x[x<-z], x[x>z]ax.fill_between(x_le_z, N_rv.pdf(x_le_z), np.zeros(len(x_le_z)), alpha=0.3, color='b')ax.fill_between(x_ge_z, N_rv.pdf(x_ge_z), np.zeros(len(x_ge_z)), alpha=0.3, color='b')p = N_rv.cdf(z) - N_rv.cdf(-z)ax.set_title('P({:.3} < z < {:.3}) = {:.3}'.format(-z, z, p));

Чем меньше p-value тем меньше вероятность того, что среднее выборки получено случайно. При этом p-value напрямую связано с двусторонними гипотезами, т.е. гипотезами о попадании величины в заданный интервал. Если мы получили какие-то результаты, но p-value оказалось довольно большим, то вряд ли эти результаты могут считаться значимыми. Причем, традиционно, уровень значимости, обозначаемый буквой равен 0.05, а это означает, что для подтверждения значимости результатов p-value должно быть меньше этого уровня. Однако, стоит обязательно отметить, что традиционный уровень значимости может быть непригоден в некоторых областях исследований. Например, в сфере образования наверняка можно обойтись , а вот в квантовой физике, запросто придется снизить этот уровень до 5 сигм, т.е. будет равна:

1 - (N_rv.cdf(5) - N_rv.cdf(-5))

5.733031438470704e-07

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

С другой стороны, что бы подтвердить небольшие отклонения от среднего генеральной совокупности, придется очень сильно увеличивать размер выборки. Так, например, если мы хотим заявить с уровнем значимости , что среднее время доставки пицы равно 31 минуте, а не 30 как считалось ранее, то придется сделать не менее 100 заказов.

На последок

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

Это кажется не очень правдоподобным, но давайте взглянем. Сгенерируем 1000 значений из равномерного, экспоненциального и Лапласова распределения, а затем, последовательно, для каждого распределения построим kde-графики распределений среднего значения выборок разного размера:

import matplotlib.animation as animationfig, axes = plt.subplots(nrows=2, ncols=3, figsize = (12, 7))unif_rv = stats.uniform(loc=10, scale=10)exp_rv = stats.expon(loc=10, scale=1.5)lapl_rv = stats.laplace(loc=15)np.random.seed(42)unif_data = unif_rv.rvs(size=1000)exp_data = exp_rv.rvs(size = 1000)lapl_data = lapl_rv.rvs(size=1000)title = ['Uniform', 'Exponential', 'Laplace']data = [unif_data, exp_data, lapl_data]y_max = [60, 310, 310]n = [3, 10, 30, 50]for i, ax in enumerate(axes[0]):    sns.histplot(data[i], bins=20, alpha=0.4, ax=ax)    ax.vlines(data[i].mean(), 0, y_max[i], color='r')    ax.set_xlim(10, 20)    ax.set_title(title[i])def animate(i):    for ax in axes[1]:        ax.clear()    for j in range(3):        rand_idx = np.random.randint(0, 1000, size=(1000, n[i]))        means = data[j][rand_idx].mean(axis=1)        sns.kdeplot(means, ax=axes[1][j])        axes[1][j].vlines(means.mean(), 0, 2, color='r')        axes[1][j].set_xlim(10, 20)        axes[1][j].set_ylim(0, 2.1)        axes[1][j].set_title('n = ' + str(n[i]))    fig.set_figwidth(15)    fig.set_figheight(8)    return axes[1][0], axes[1][1],axes[1][2]    dist_animation = animation.FuncAnimation(fig,                                       animate,                                       frames=np.arange(4),                                      interval = 200,                                      repeat = False)#  Сохраняем анимацию в виде gif файла:dist_animation.save('dist_means.gif',                 writer='imagemagick',                  fps=1)

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

В общем, спасибо за внимание. Жму F5 и жду ваших и комментариев.

Линейный график является неотъемлемой частью анализа данных. Он даёт нам представление о том, как величина изменяется при последовательных измерениях. В случае работы с временными рядами важность линейных графиков становится решающей. Тренд [направление], сезонность и корреляция вот некоторые характеристики, которые можно наблюдать на аккуратно сгенерированных линейных графиках. В этой статье мы будем создавать интерактивные линейные графики с помощью двух библиотек Python Pandas и Altair.

Мы уже затрагивали тему визуализаций при помощи библиотеки Altair на примере создания интерактивных карт, а сегодня, к старту курса о Data Science, решили поделиться простым руководством о том, как можно из множества графиков выделить самый важный; с этого руководства можно начать изучать Altair на практике.

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

Давайте начнём с генерации данных. Типичный пример использования линейных графиков это анализ цен на акции. Один из самых простых способов получения данных о ценах на акции предоставляет библиотека pandas-datareader. Сначала нам нужно импортировать её вместе с уже установленной в Google Colab Pandas.

import pandas as pdfrom pandas_datareader import data

Мы получим цены на акции трёх компаний за один год. Необходимо указать начальную дату, конечную дату и источник данных:

start = '2020-1-1'end = '2020-12-31'source = 'yahoo'

Необходимо знать ещё одну деталь название акции:

apple = data.DataReader("AAPL", start=start ,end=end, data_source=source).reset_index()[["Date", "Close"]]ibm = data.DataReader("IBM", start=start ,end=end, data_source=source).reset_index()[["Date", "Close"]]microsoft = data.DataReader("MSFT", start=start ,end=end, data_source=source).reset_index()[["Date", "Close"]]

Сейчас у нас есть цены на акции Apple, IBM и Microsoft в 2020 году. Лучше располагать их в одном фрейме данных. Перед объединением нам нужно добавить колонку, которая указывает, к какой акции относится та или иная цена. Следующий блок кода добавляет соответствующие столбцы, а затем объединяет фреймы данных с помощью функции concat:

apple["Stock"] = "apple"ibm["Stock"] = "ibm"microsoft["Stock"] = "msft"stocks["Month"] = stocks.Date.dt.monthstocks = pd.concat([apple, ibm, microsoft])

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

Altair

Altair это библиотека статистической визуализации на Python. Как мы увидим в примерах, её синтаксис чист и прост для понимания. Также с помощью Altair очень просто создавать интерактивные визуализации. Я кратко объясню структуру Altair, а затем сосредоточусь на создании интерактивных линейных графиков. Если вы новичок в Altair, вот туториал по Altair в виде серии из 4 частей:

Вот простой линейный график без какой-либо интерактивности:

alt.Chart(stocks).mark_line().encode(   x="Date",   y="Close",   color="Stock").properties(   height=300, width=500)

Основная структура начинается с объекта диаграммы верхнего уровня. Данные могут иметь формат фрейма данных Pandas, или можно написать строку с URL, указывающим на файл JSON или CSV. Затем указывается тип визуализации (mark_circle, mark_line и так далее).

Функция encode указывает Altair, что нужно построить в заданном фрейме данных. Таким образом, всё, что мы пишем в функции encode, должно быть связано с данными. Различать акции мы будем при помощи параметра color, который аналогичен параметру hue в Seaborn. Наконец, с помощью функции properties задаются определённые свойства графика.

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

selection = alt.selection_multi(fields=["Stock"], bind="legend")alt.Chart(stocks).mark_line().encode(   x="Date",   y="Close",   color="Stock",   opacity=alt.condition(selection, alt.value(1), alt.value(0.1))).properties(   height=300, width=500).add_selection(   selection)

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

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

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

hover = alt.selection(   type="single", on="mouseover", fields=["Stock"], nearest=True)

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

В коде ниже есть 3 компонента. Первый создаёт линейный график. Второй компонент это диаграмма рассеяния, построенная на линейном графике, эта диаграмма используется для определения ближайшей точки. Мы настроим непрозрачность так, чтобы диаграмма рассеяния не была видна. Третий компонент отвечает за выделение линии, содержащей захваченную точку на втором графике:

# line plotlineplot = alt.Chart(stocks).mark_line().encode(   x="Date:T",   y="Close:Q",   color="Stock:N",)# nearest pointpoint = lineplot.mark_circle().encode(   opacity=alt.value(0)).add_selection(hover)# highlightsingleline = lineplot.mark_line().encode(   size=alt.condition(~hover, alt.value(0.5), alt.value(3)))

Теперь, объединив второй и третий графики, можно создать интерактивный линейный график:

point + singleline

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

Заключение

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

Визуализация данных стоит усилий, потому что она помогает по-новому смотреть на данные в прямом, и переносном смыслах. Но визуализация это лишь один из аспектов работы с данными. Если хотите узнать об этом больше можете обратить внимание на наш курс о Data Science, в рамках которого вы узнаете как выгружать данные, очищать их, анализировать и использовать для решения своих задач.

Узнайте, как прокачаться и в других специальностях или освоить их с нуля:

• Профессия Data Scientist

• Профессия Data Analyst

• Курс по Data Engineering

Доброго времени суток, хабраледи и хабраджентельмены! В этой статье мы продолжим погружение в статистику вместе с Python. Если кто пропустил начало погружения, то вот ссылка на первую часть. Ну, а если нет, то я по-прежнему рекомендую держать под рукой открытую книгу Сары Бослаф "Статистика для всех". Так же рекомендую запустить блокнот, чтобы поэкспериментировать с кодом и графиками.

Как сказал Эндрю Ланг: "Статистика для политика все равно что уличный фонарь для пьяного забулдыги: скорее опора, чем освещение." Тоже самое можно сказать и про эту статью для новичков. Вряд ли вы почерпнете здесь много новых знаний, но надеюсь, эта статья поможет вам разобраться с тем, как использовать Python для облегчения самостоятельного изучения статистики.

Зачем выдумывать новые распределения?

Представим... так, прежде чем нам что-то представлять, давайте снова сделаем все необходимые импорты:

import numpy as npfrom scipy import statsimport matplotlib.pyplot as pltimport seaborn as snsfrom pylab import rcParamssns.set()rcParams['figure.figsize'] = 10, 6%config InlineBackend.figure_format = 'svg'np.random.seed(42)

А вот теперь давайте представим, что мы пивовары, придумавшие рецепт и технологию варки классного пива. Продажи идут хорошо, да и в целом все неплохо, но все-таки очень интересно, как потребители оценивают качество нашего пива. Мы решили опросить 1000 человек и оценить качество напитка по 100-бальной шкале в сравнении с другими сортами пива, которые им доводилось пробовать ранее. Данные опроса могли бы выглядеть так:

gen_pop = np.trunc(stats.norm.rvs(loc=80, scale=5, size=1000))gen_pop[gen_pop>100]=100print(f'mean = {gen_pop.mean():.3}')print(f'std = {gen_pop.std():.3}')

mean = 79.5std = 4.95

В данном случае, мы можем предположить, что оценки распределены нормально с мат. ожиданием и стандартным отклонением 80 и 5 баллов соответственно. Такой результат означает, что массовый потребитель, конечно же, не считает наше пиво лучшим в мире, но в среднем, видит его где-то в числе лучших.

Все могло бы быть ровно да гладко, особенно в нашем воображении. Но теперь давайте представим, что случилась внезапная неприятность - поставщик хмеля пострадал от неурожая и нам пришлось срочно покупать сырье у другого производителя. Очевидно, что это как-то скажется на качестве пива, но как? Единственный способ проверить это - сварить небольшую партию с новым хмелем и провести опрос у небольшой группы респондентов. Допустим, в оценке качества нового пива приняли участие всего 10 человек, а их оценки распределились следующим образом:

Теперь мы можем вычислить Z-статистику по уже известной нам формуле:

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

sample = np.array([89,99,93,84,79,61,82,81,87,82])sample.mean()

83.7

Вот значение Z-статистики:

z = 10**0.5*(sample.mean()-80)/5z

2.340085468524603

И еще обязательно вычисляем p-value:

1 - (stats.norm.cdf(z) - stats.norm.cdf(-z))

0.019279327322753836

Для нас, как для пивоваров, это крайне интересные цифры и вот почему: Z-значение отстоит от 0 более чем на 2 сигмы, т.е. отклонение среднего балла 10 респондентов находится очень далеко от среднего бала генеральной совокупности, а вероятность того, что это отклонение произошло случайно, всего около 0.02. Другими словами, если взять 10 человек из совокупности людей, которые оценивают наше "старое" пиво как , то вероятность того, что эти 10 человек могут случайно дать среднюю оценку "новому" пиву равную 83.7 баллам составляет всего около 2%. Скорее всего, если мы проведем полномасштабный опрос, то увидим, что потребители не просто не заметят изменения качества пива, а даже отметят его небольшое улучшение. Круто.

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

sample.std(ddof=1)

10.055954565441423

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

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

Для того, чтобы учитывать поправку в методе std() пакета NumPy есть параметр ddof, который по умолчанию равен 0, но если мы применяем std() к небольшой выборке, то необходимо явно указать, что ddof=1. Влияние данного параметра на расчеты может быть очень значимым. Например, если мы возьмем 10000 выборок по 10 элементов из и построим распределение стандартных отклонений, то увидим, что при ddof=0 пик распределения находится левее реального значения стандартного отклонения. При ddof=1 пик распределения так же расположен левее реального значения, но все-таки ближе к нему, чем при ddof=0:

fig, ax = plt.subplots(nrows=1, ncols=2, figsize = (12, 5))for i in [0, 1]:    deviations = np.std(stats.norm.rvs(80, 5, (10000, 10)), axis=1, ddof=i)    sns.histplot(x=deviations ,stat='probability', ax=ax[i])    ax[i].vlines(5, 0, 0.06, color='r', lw=2)    ax[i].set_title('ddof = ' + str(i), fontsize = 15)    ax[i].set_ylim(0, 0.07)    ax[i].set_xlim(0, 11)fig.suptitle(r'$s={\sqrt {{\frac {1}{10 - \mathrm{ddof}}}\sum _{i=1}^{10}\left(x_{i}-{\bar {x}}\right)^{2}}}$',             fontsize = 20, y=1.15);

Является ли этот факт важным и можно ли теперь быть уверенными в выводе, который мы сделали на основе Z-статистики? Первое, что приходит в голову - посмотреть, как может быть распределено стандартное отклонение большого количества выборок из генеральной совокупности. Давайте сделаем 5000 тысяч выборок по 10 человек, взятых из распределения и посмотрим на распределение стандартного отклонения баллов этих выборок:

deviations  = np.std(stats.norm.rvs(80, 5, (5000, 10)), axis=1, ddof=1)sns.histplot(x=deviations ,stat='probability');

Это похоже на колокол нормального распределения, но важнее то, что стандартное отклонение выборочных оценок не отклоняется до значения равного 10-и баллам. А это крайне настораживающе. Судите сами, если вероятность того, что отклонение среднего балла наших 10 респондентов от среднего генеральной совокупности составляет всего 2%, то вероятность того, что стандартное отклонение (уж простите за каламбур) может отклониться до значения в 10 баллов стремится к 0. А это означает, что для наших выводов об улучшении качества пива, есть вполне уместное замечание: может быть опрос 10-и респондентов свидетельствует не только об изменении среднего генеральной совокупности, может быть ее стандартное отклонение тоже изменилось.

В конце концов, добавление нового сорта хмеля могло привести к большей специфичности вкуса, следовательно, вызвать большую разнонаправленность в предпочтениях: кому-то этот вкус пива стал нравиться больше, кому-то меньше, но средняя оценка могла остаться прежней. Например, если генеральная совокупность оценок "нового" пива теперь распределена как , то Z-статистика и p-value для имеющихся 10-и оценок будет выглядеть не так обнадеживающе:

z = 10**0.5*(sample.mean()-80)/10p = 1 - (stats.norm.cdf(z) - stats.norm.cdf(-z))print(f'z = {z:.3}')print(f'p-value = {p:.4}')

z = 1.17p-value = 0.242

Оказывается, что для распределения резульататы оказались значимыми, но как только мы оценили стандартное отклонение генеральной совокупности тем же значением, что и у выборки, т.е. представили, что генеральная совокупность распределена как , то это очень сильно отразилось на значимости имеющихся оценок. Ведь вероятность получить такое же или даже еще более отклоненное среднее теперь равна не 2%, а почти 25%. Хотя, всего-то заменили на .

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

По сути мы просто придумали новую T-статистику, которая отличается от Z-статистики только тем, что в знаменателе стоит не генеральной совокупности, а выборки. А теперь сделаем 10000 выборок из , вычислим для каждой выборки Z- и T-статистику, а распределение значений изобразим в виде гистограмм:

fig, ax = plt.subplots(nrows=1, ncols=2, figsize = (12, 5))N = 10000samples = stats.norm.rvs(80, 5, (N, 10))statistics = [lambda x: 10**0.5*(np.mean(x, axis=1) - 80)/5,              lambda x: 10**0.5*(np.mean(x, axis=1) - 80)/np.std(x, axis=1, ddof=1)]title = 'ZT'bins = np.linspace(-6, 6, 80, endpoint=True)for i in range(2):    values = statistics[i](samples)    sns.histplot(x=values ,stat='probability', bins=bins, ax=ax[i])    p = values[(values > -2)&(values < 2)].size/N    ax[i].set_title('P(-2 < {} < 2) = {:.3}'.format(title[i], p))    ax[i].set_xlim(-6, 6)    ax[i].vlines([-2, 2], 0, 0.06, color='r');

Я все-таки не удержался и сделал гифку:

import matplotlib.animation as animationfig, axes = plt.subplots(nrows=1, ncols=2, figsize = (18, 8))def animate(i):    for ax in axes:        ax.clear()    N = 10000    samples = stats.norm.rvs(80, 5, (N, 10))    statistics = [lambda x: 10**0.5*(np.mean(x, axis=1) - 80)/5,                  lambda x: 10**0.5*(np.mean(x, axis=1) - 80)/np.std(x, axis=1, ddof=1)]    title = 'ZT'    bins = np.linspace(-6, 6, 80, endpoint=True)    for j in range(2):        values = statistics[j](samples)        sns.histplot(x=values ,stat='probability', bins=bins, ax=axes[j])        p = values[(values > -2)&(values < 2)].size/N        axes[j].set_title(r'$P(-2\sigma < {} < 2\sigma) = {:.3}$'.format(title[j], p))        axes[j].set_xlim(-6, 6)        axes[j].set_ylim(0, 0.07)        axes[j].vlines([-2, 2], 0, 0.06, color='r')    return axesdist_animation = animation.FuncAnimation(fig,                                       animate,                                       frames=np.arange(7),                                      interval = 200,                                      repeat = False)dist_animation.save('statistics_dist.gif',                 writer='imagemagick',                  fps=1)

Просто мне кажется, что так можно продемонстрировать, как едва заметные черты могут оказаться чрезвычайно важными для исследователя. На что я хочу обратить ваше внимание? Во-первых, и слева, и справа мы видим свиду два одинаковых колоколообразных распределения. Вполне уместно предположить, что это два стандартных нормальных распределения, верно? Но, вот что любопытно, в интервал должен содержать около 95.5% всех значений. Для Z-статистики это требование выполняется, а для придуманной нами T-статистики не очень, потому что только 92-93% ее значений укладываются в заданный интервал. Казалось бы, что отличие не столь велико, чтобы заподозрить какую-то закономерность, но мы можем провести больше экспериментов:

statistics = [lambda x: 10**0.5*(np.mean(x, axis=1) - 80)/5,              lambda x: 10**0.5*(np.mean(x, axis=1) - 80)/np.std(x, axis=1, ddof=1)]quantity = 50N=10000result = []for i in range(quantity):    samples = stats.norm.rvs(80, 5, (N, 10))    Z = statistics[0](samples)    p_z = Z[(Z > -2)&((Z < 2))].size/N    T = statistics[1](samples)    p_t = T[(T > -2)&((T < 2))].size/N    result.append([p_z, p_t])result = np.array(result)fig, ax = plt.subplots()line1, line2 = ax.plot(np.arange(quantity), result)ax.legend([line1, line2],           [r'$P(-2\sigma < {} < 2\sigma)$'.format(i) for i in 'ZT'])ax.hlines(result.mean(axis=0), 0, 50, color='0.6');

В каждом из 50 экспериментов мы видим одно и тоже. Можно было бы заподозрить, что есть ошибка в коде, но его не так много, и легко убедиться в том, что вычисления выполняются верно. Так в чем же дело? А дело в том, что мы получили совершенно новый тип распределения! Снова взгляните на гифку с гистограммами распределений значений Z- и T- статистик, присмотритесь к основанию колокола каждой из них. Вам не кажется, что у распределения T-статистик основание чуть шире? Это хорошо видно по выпирающим за красные линии, так называемым - хвостам. А то, что эти хвосты несколько больше, или как еще говорят - тяжелее, чем у нормального распределения, так же будет означать, что мы будем наблюдать несколько больше сильных отклонений от вершины распределения. Проще говоря, мы теперь можем учитывать дисперсию выборки при оценке параметров генеральной совокупности. Однако, мы так и не ответили на вопрос - хорошо ли, можно ли, да и вообще зачем оценивать генеральной совокупности значением стандартного отклонения выборки.

Задумаемся о дисперсии

Когда мы рассматривали распределение значений Z-статистик, мы обращали внимание только на распределение выборочного среднего и даже не задавались вопросом о том, что выборочное стандартное отклонение тоже может как-то распределяться. Давайте возьмем 10000 выборок из по 10 элементов в каждой и посмотрим, как будет выглядеть зависимость стандартного отклонения выборок от их среднего значения:

# если график строиться слишком долго,# то смените формат svg на png:#%config InlineBackend.figure_format = 'png'N = 10000samples = stats.norm.rvs(80, 5, (N, 10))means = samples.mean(axis=1)deviations = samples.std(ddof=1, axis=1)T = statistics[1](samples)P = (T > -2)&((T < 2))fig, ax = plt.subplots()ax.scatter(means[P], deviations[P], c='b', alpha=0.7,           label=r'$\left | T \right | < 2\sigma$')ax.scatter(means[~P], deviations[~P], c='r', alpha=0.7,           label=r'$\left | T \right | > 2\sigma$')mean_x = np.linspace(75, 85, 300)s = np.abs(10**0.5*(mean_x - 80)/2)ax.plot(mean_x, s, color='k',        label=r'$\frac{\sqrt{n}(\bar{x}-\mu)}{2}$')ax.legend(loc = 'upper right', fontsize = 15)ax.set_title('Зависимость выборочного стандартного отклонения\nот выборочного среднего',             fontsize=15)ax.set_xlabel(r'Среднее значение выборки ($\bar{x}$)',              fontsize=15)ax.set_ylabel(r'Стандартное отклонение выборки ($s$)',              fontsize=15);

Данный график интересен тем, что показывает зависимость среднего выборки от ее стандартного отклонения. На нем видно, что выборки, у которых и отклоняется очень сильно и одновременно практически не встречаются, т.е. точки стараются избегать углов графика. Это значит, что, если мы берем выборку из нормального распределения , то мы вряд ли увидим слишком большое значение при большом значении . А еще из данного графика следует, что экстремальные значения (красные точки) могут быть получены, только если выполняется следующее соотношение:

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

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

Как интерпретировать это расположение и цвет точек? Давайте вспомним пример, от которого мы отталкивались. Мы сварили пиво с новым сортом хмеля, нашли десять человек (наша выборка) и попросили их оценить качество нового пива по 100-больной шкале. На основании их оценок мы пытаемся понять, как могут оценить качество пива все, кто его купит в дальнейшем (генеральная совокупность). Допустим среднее полученных 10-и оценок равно 82-м баллам, а стандартное отклонение равно всего 2-м баллам. Могут ли такие оценки быть получены случайно, если предположить, что генеральная совокупность имеет мат.ожидание , а стандартное отклонение такое же как у выборки, т.е. ? Чтобы выяснить это нужно вычислить Z-статистику:

Чтобы оценить вероятность случайного получения такого значения из нужно вычислить p-value:

z = 10**0.5*(82-80)/2p = 1 - (stats.norm.cdf(z) - stats.norm.cdf(-z))print(f'p-value = {p:.2}')

p-value = 0.0016

Вероятность того что 10 человек случайно дадут средню оценку равную 82-м баллам очень мала и составляет менее 2%. Но данное утверждение является верным только если мы предполагаем что эти десять человек взяты из генеральной совокупности всех потребителей в которой их оценки распределены как . И если мы действительно предполагаем, что , то мы можем быть вполне уверены, что оценки обусловлены новым сортом хмеля в составе пива, а не случайностью.

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

С другой стороны мы могли бы пригласить крайне разношерстную публику в качестве 10 респондентов. Они могли бы оценить качество пива так же в 82 балла, но со стандартным отклонением равным, допустим, 9-и баллам. Может ли такая оценка быть получена случайно? Проверим:

z = 10**0.5*(82-80)/9p = 1 - (stats.norm.cdf(z) - stats.norm.cdf(-z))print(f'p-value = {p:.2}')

p-value = 0.48

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

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

import matplotlib.animation as animationfig, ax = plt.subplots(figsize = (15, 9))def animate(i):    ax.clear()    N = 10000        samples = stats.norm.rvs(80, 5, (N, i))    means = samples.mean(axis=1)    deviations = samples.std(ddof=1, axis=1)    T = i**0.5*(np.mean(samples, axis=1) - 80)/np.std(samples, axis=1, ddof=1)    P = (T > -2)&((T < 2))        prob = T[P].size/N    ax.set_title(r'зависимость $s$ от $\bar{x}$ при $n = $' + r'${}$'.format(i),                 fontsize = 20)    ax.scatter(means[P], deviations[P], c='b', alpha=0.7,               label=r'$\left | T \right | < 2\sigma$')    ax.scatter(means[~P], deviations[~P], c='r', alpha=0.7,               label=r'$\left | T \right | > 2\sigma$')    mean_x = np.linspace(75, 85, 300)    s = np.abs(i**0.5*(mean_x - 80)/2)    ax.plot(mean_x, s, color='k',            label=r'$\frac{\sqrt{n}(\bar{x}-\mu)}{2}$')    ax.legend(loc = 'upper right', fontsize = 15)    ax.set_xlim(70, 90)    ax.set_ylim(0, 10)    ax.set_xlabel(r'Среднее значение выборки ($\bar{x}$)',              fontsize='20')    ax.set_ylabel(r'Стандартное отклонение выборки ($s$)',              fontsize='20')    return axdist_animation = animation.FuncAnimation(fig,                                       animate,                                       frames=np.arange(5, 21),                                      interval = 200,                                      repeat = False)dist_animation.save('sigma_rel.gif',                 writer='imagemagick',                  fps=3)

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

Но постойте! Разве мы не вернулись к тому, с чего начали? Ведь мы снова так и не ответили на вопрос - почему при малых размерах выборки мы предполагаем, что стандартное отклонение генеральной совокупности должно быть таким же, как у выборки. В нашем конкретном примере, с которого все началось, оценки 10-и респондентов выглядели так:

Среднее этой выборки , а стандартное оклонение , причем мы приняли во внимание, что полномасштабный опрос, до введения в состав нового сорта хмеля, показал, что оценки распределены как . Но вместо Z-статистики, мы решили использовать придуманную T-статистику, единственное отличие которой от Z-статистики состоит в том, что мы просто заменили в формуле генеральной совокупности на выборки. По сути мы просто решили проверить гипотезу относительно генеральной совокупности, которая теперь по каким-то причинам распределена не как , а как - это вообще правильно? Разве не нужно было бы проверить это предположение для самых разных значений ?: Почему мы не предполагаем, что , , или любое другое значение ?

Чтобы ответить на этот вопрос, давайте снова вернемся к изначальному смыслу доказательства гипотез. Суть состоит в том, чтобы на основании выборки из генеральной совокупности сделать какие-то выводы о параметрах распределения этой генеральной совокупности. В нашем примере мы были уверены, что оценки его качества распределены как , однако, когда мы сварили пиво с новым сортом хмеля, то мы получили 10 оценок, стандартное отклонение которых составило 10 баллов. В самом начале мы видели, что получить выборку с таким отклонением из просто невозможно. А это значит, что мы просто вынуждены как-то по новому оценить стандартное отклонение генеральной совокупности, иначе все наши выводы не будут иметь никакого отношения к реальности.

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

N = 10000sigma = np.linspace(5, 20, 151)prob = []for i in sigma:    p = []    for j in range(10):        samples = stats.norm.rvs(80, i, (N, 10))        means = samples.mean(axis=1)        deviations = samples.std(ddof=1, axis=1)        p_m = means[(means >= 83) & (means <= 84)].size/N        p_d = deviations[(deviations >= 9.5) & (deviations <= 10.5)].size/N        p.append(p_m*p_d)    prob.append(sum(p)/len(p))prob = np.array(prob)fig, ax = plt.subplots()ax.plot(sigma, prob)ax.set_xlabel(r'Стандартное отклонение генеральной совокупности ($\sigma$)',              fontsize=20)ax.set_ylabel('Вероятность',              fontsize=20);

Как видите, максимум вероятности достигается при . Именно по этой причине, когда мы имеем дело с небольшими выборками, мы проверяем гипотезу о распределении генеральной совокупности не относительно известного нам значения ее стандартного отклонения , а относительно стандартного отклонения выборки . Если параметры генеральной совокупности как-то изменились, то мы заранее предполагаем, что ее стандартное отклонение теперь равно стандартному отклонению выборки, просто потому, что вероятность такого изменения максимальна. Вот и все.

А это точно мы изобрели T-статистику?

Давайте представим, что некий инженер в какой-нибудь компании по производству шариков для пинг-понга решил задачу по сверхплотной упаковке шаров. После этого в один грузовой контейнер стало влезать на 1% больше шариков, что в общем-то не густо. Но тем не менее, глава компании мог посчитать этот результат очень важным и запретить инженеру распространять свое решение, так как это решение теперь является коммерческой тайной. Вот так, важный для цивилизации теоретический результат мог бы остаться в тайне ради ничтожного увеличения прибыли какой-то компании по производству шариков для пинг-понга. Абсурд?

Несомненно - абсурд! Однако, что-то подобное чуть не произошло с Уильямом Госсетом, когда он работал на пивоварне "Гиннес" и открыл t-распределение. Оценить практическую значимость некоторых статистических открытий очень сложно, но не в тех случаях, когда на кону астрономические потенциальные прибыли и стратегическое развитие компании, или даже государства. Например, последовательный критерий отношений правдоподобия Вальда, открытый в 1943 году во время войны вообще засекретили, так как этот критерий позволял на 50% уменьшить среднее число наблюдений в выборке. Короче, статистика - это крайне полезная наука.

Итак, Уильям Госсет отвечал за контроль качества в пивоварне "Гиннес" и пытался решить проблему появления ошибок при использовании выборок ограниченного размера. Госсет проделал все то же самое, что мы проделали с вами выше (только без компьютера!) и сделал вывод о том, что если выборка берется из "нормальной" генеральной совокупности и что если использовать стандартное отклонение выборки для оценки стандартного отклонение этой самой генеральной совокупности, то распрделение выборочного среднего можно описать уже знакомой нам формулой:

Данная формула и является формулой так называемого t-распределения, оно же распределение Стьюдента. Учитывая, что история о том, как Госсет использовал псевдоним "Стьюдент", для обхода наложенным Гиннесом запрета на публикации работ, довольно широко известна, то, наверное, было бы правильнее называть это распределение - распределением Госсета. Но названия "распределение Стьюдента", или используемое немного реже "t-распределение" так сильно устоялись, что мы не будем нарушать эту общепринятую договоренность.

Еще один способ получить распределение Стьюдента заключается в использовании вот этой формулы:

которая позволяет прояснить понятие "степени свободы" выборки. Пусть случайные переменные берутся из стандартного нормального распределения, т.е. , тогда значение , т.е. объем выборки, как раз и определяет степень свободы распределения Стьюдента. А принадлежность какой-нибудь случайной величины к распределению Стьюдента с определенным числом степеней свободы обозначается, как:

Давайте попробуем взглянуть на то, как будет меняться распределение Стьюдента при увеличении числа степеней свободы:

import matplotlib.animation as animationfig, ax = plt.subplots(figsize = (15, 9))def animate(i):    ax.clear()    N = 15000        x = np.linspace(-5, 5, 100)    ax.plot(x, stats.norm.pdf(x, 0, 1), color='r')        samples = stats.norm.rvs(0, 1, (N, i))        t = samples[:, 0]/np.sqrt(np.mean(samples[:, 1:]**2, axis=1))    t = t[(t>-5)&(t<5)]    sns.histplot(x=t, bins=np.linspace(-5, 5, 100), stat='density', ax=ax)        ax.set_title(r'Форма распределения $t(n)$ при n = ' + str(i), fontsize = 20)    ax.set_xlim(-5, 5)    ax.set_ylim(0, 0.5)    return axdist_animation = animation.FuncAnimation(fig,                                       animate,                                       frames=np.arange(2, 21),                                      interval = 200,                                      repeat = False)dist_animation.save('t_rel_of_df.gif',                 writer='imagemagick',                  fps=3)

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

Но лучше воспользуемся готовой реализацией в SciPy:

import matplotlib.animation as animationfig, ax = plt.subplots(figsize = (15, 9))def animate(i):    ax.clear()    N = 15000        x = np.linspace(-5, 5, 100)    ax.plot(x, stats.norm.pdf(x, 0, 1), color='r')    ax.plot(x, stats.t.pdf(x, df=i))        ax.set_title(r'Форма распределения $t(n)$ при n = ' + str(i), fontsize = 20)    ax.set_xlim(-5, 5)    ax.set_ylim(0, 0.45)    return axdist_animation = animation.FuncAnimation(fig,                                       animate,                                       frames=np.arange(2, 21),                                      interval = 200,                                      repeat = False)dist_animation.save('t_pdf_rel_of_df.gif',                 writer='imagemagick',                  fps=3)

Здесь гораздо лучше видно, что при увеличении (параметр df в коде) распределение Стьюдента стремится к нормальному распределению. Что в общем-то и неудивительно, так как ранее мы уже видели, что при увеличении выборки, ее стандартное отклонение все меньше колеблется вокруг стандартного отклонения генеральной совокупности.

t-распределение на практике

Давайте сразу выполним t-тест из модуля SciPy для выборки из нашего примера с пивом:

sample = np.array([89,99,93,84,79,61,82,81,87,82])stats.ttest_1samp(sample, 80)

Ttest_1sampResult(statistic=1.163532240174695, pvalue=0.2745321678073461)

А теперь попробуем получить тот же самый результат самостоятельно:

T = 9**0.5*(sample.mean() -80)/sample.std()T

1.163532240174695

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

T = 10**0.5*(sample.mean() -80)/sample.std(ddof=1)T

1.1635322401746953

Хорошо, значение t-статистики мы вычислили, как вычислить значение p-value? Думаю, вы и до этого уже догадались - нужно проделать все то же самое, что мы делали при вычислении p-value для Z-статистики, но только с использованием t-распределения:

t = stats.t(df=9)fig, ax = plt.subplots()x = np.linspace(t.ppf(0.001), t.ppf(0.999), 300)ax.plot(x, t.pdf(x))ax.hlines(0, x.min(), x.max(), lw=1, color='k')ax.vlines([-T, T], 0, 0.4, color='g', lw=2)x_le_T, x_ge_T = x[x<-T], x[x>T]ax.fill_between(x_le_T, t.pdf(x_le_T), np.zeros(len(x_le_T)), alpha=0.3, color='b')ax.fill_between(x_ge_T, t.pdf(x_ge_T), np.zeros(len(x_ge_T)), alpha=0.3, color='b')p = 1 - (t.cdf(T) - t.cdf(-T))ax.set_title(r'$P(\left | T \right | \geqslant {:.3}) = {:.3}$'.format(T, p));

Мы видим, что p-value чуть больше 27%, т.е. вероятность того, что результат получен случайно - довольно велика, особенно в соответствии с самым распространенным уровнем значимости , который меньше p-value более, чем в 5 раз. Что ж, наши результаты не являются значимыми, но они пригодны для демонстрации вычисления доверительного интервала - границы значений математического ожидания генеральной совокупности, в пределах которых оно находится с заданной вероятностью , обычно равной 0.95:

Чтобы вычислить данный интервал в SciPy, достаточно воспользоваться методом interval с заданными параметрами loc (смещение) и scale (масштаб) для нашей выборки:

sample_loc = sample.mean()sample_scale = sample.std(ddof=1)/10**0.5ci = stats.t.interval(0.95, df=9, loc=sample_loc, scale=sample_scale)ci

(76.50640345566619, 90.89359654433382)

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

Напоследок

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

Конечно, хотелось бы в конце вставить какую-нибудь гифку, но хочу закончить фразой Герберта Спенсера: "Величайшая цель образования - не знание, а действие", так что запускайте свои анаконды и действуйте! Особенно это касается самоучек, вроде меня.

Жму F5 и жду ваших комментариев!

В этой статье, переводом которой мы решили поделиться специально к старту курса о Data Science, автор представляет новый пакет Python для генерации кластерограмм из решений кластеризации. Библиотека была разработана в рамках исследовательского проекта Urban Grammar и совместима со scikit-learn и библиотеками с поддержкой GPU, такими как cuML или cuDF в рамках RAPIDS.AI.

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

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

Маттиас Шонлау предложил другой подход кластерограмму. Кластерограмма это двухмерный график, отражающий потоки наблюдений между классами по мере добавления кластеров. Это говорит вам о том, как перетасовываются ваши данные и насколько хороши ваши сплиты. Тал Галили позже реализовал кластерограмму для k-средних в R. Я использовал реализацию Таля, перенёс ее на Python и создал clustergram пакет Python для создания кластерограмм.

clustergram в настоящее время поддерживает метод k-средних, использование scikit-learn (включая реализацию Mini-Batch) и RAPIDS.AI cuML (если у вас есть GPU с поддержкой CUDA), Gaussian Mixture Model (только scikit-learn) и иерархическую кластеризацию на основе scipy.hierarchy. В качестве альтернативы мы можем создать кластерограмму на основе меток и данных, полученных с помощью альтернативных пользовательских алгоритмов кластеризации. Пакет предоставляет API, подобный sklearn, и строит кластерные диаграммы с помощью matplotlib, что даёт ему широкий выбор вариантов оформления в соответствии со стилем вашей публикации.

Установка

Установить clustergram можно при помощи conda или pip:

conda install clustergram -c conda-forge

или

pip install clustergram

В любом случае вам нужно установить выбранный бэкенд (scikit-learn и scipy или cuML).

from clustergram import Clustergramimport urbangrammar_graphics as uggimport seaborn as snsimport matplotlib.pyplot as pltfrom sklearn.preprocessing import scalesns.set(style='whitegrid')

Давайте рассмотрим несколько примеров, чтобы понять, как выглядит кластерограмма и что с ней делать.

Набор данных о цветке ириса

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

iris = sns.load_dataset("iris")g = sns.pairplot(iris, hue="species", palette=ugg.COLORS[1:4])g.fig.suptitle("Iris flowers", y=1.01)

Похоже, что setosa относительно чётко определённая группа, тогда как разница между versicolor и virginica меньше, поскольку они частично перекрываются (или, в случае ширины чашелистика, полностью).

Итак, мы знаем, как выглядят данные. Теперь мы можем увидеть кластерограмму. Помните, мы знаем, что существует три кластера, и в идеале мы должны быть в состоянии распознать это по кластерограмме. Я говорю в идеале, потому что, даже если есть известные метки, это не значит, что наши данные или метод кластеризации способны различать эти классы.

Давайте начнём с кластеризации методом k-средних. Чтобы получить стабильный результат, мы можем запустить кластерную программу с 1000 инициализаций.

data = scale(iris.drop(columns=['species']))cgram = Clustergram(range(1, 10), n_init=1000)cgram.fit(data)ax = cgram.plot(    figsize=(10, 8),    line_style=dict(color=ugg.COLORS[1]),    cluster_style={"color": ugg.COLORS[2]},)ax.yaxis.grid(False)sns.despine(offset=10)ax.set_title('K-Means (scikit-learn)')

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

Мы ищем разделение, т. е. отвечаем на вопрос, принёс ли дополнительный кластер какое-либо значимое разделение? Шаг от одного кластера к двум большой хорошее и чёткое разделение. От двух до трёх свидетельство довольно хорошего раскола в верхней ветви. Но с 3 по 4 видимой разницы нет, потому что новый четвёртый кластер почти не отличается от существующей нижней ветви. Хотя сейчас она разделена на две части, это разделение не даёт нам много информации. Таким образом, можно сделать вывод, что идеальное количество кластеров для данных Iris три.

Мы также можем проверить некоторую дополнительную информацию, например оценку силуэта или оценку Калинского Харабазса.

fig, axs = plt.subplots(2, figsize=(10, 10), sharex=True)cgram.silhouette_score().plot(    xlabel="Number of clusters (k)",    ylabel="Silhouette score",    color=ugg.COLORS[1],    ax=axs[0],)cgram.calinski_harabasz_score().plot(    xlabel="Number of clusters (k)",    ylabel="Calinski-Harabasz score",    color=ugg.COLORS[1],    ax=axs[1],)sns.despine(offset=10)

По этим графикам можно предположить наличие 34 кластеров по аналогии с кластерограммой, но они не очень убедительны.

Набор данных о пингвинах со станции Палмера

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

penguins = sns.load_dataset("penguins")g = sns.pairplot(penguins, hue="species", palette=ugg.COLORS[3:])g.fig.suptitle("Palmer penguins", y=1.01)

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

data = scale(penguins.drop(columns=['species', 'island', 'sex']).dropna())cgram = Clustergram(range(1, 10), n_init=1000)cgram.fit(data)ax = cgram.plot(    figsize=(10, 8),    line_style=dict(color=ugg.COLORS[1]),    cluster_style={"color": ugg.COLORS[2]},)ax.yaxis.grid(False)sns.despine(offset=10)ax.set_title("K-Means (scikit-learn)")

Мы ищем разделения, и эта кластерограмма показывает достаточное их количество. На самом деле определить оптимальное количество кластеров довольно сложно. Однако, поскольку мы знаем, что происходит между различными вариантами, мы можем поиграть с этим. Если у нас есть причина быть консервативными, мы можем обойтись 4 кластерами (я знаю, это уже больше, чем первоначальный вид). Но и дальнейшее разделение также разумно, а это указывает на то, что даже более высокая детализация может дать полезную информацию, что могут существовать значимые группы.

Можно ли сказать, что их три? Поскольку мы знаем, что их должно быть три... Ну, не совсем. Разница между разделениями 23 и 34 незначительна. Однако здесь виновником является метод K ближайших соседей, а не кластерограмма. Он просто не может правильно кластеризовать эти данные из-за наложений и общей структуры. Давайте посмотрим, как работает смешанная Гауссова модель (Gaussian Mixture).

cgram = Clustergram(range(1, 10), n_init=100, method="gmm")cgram.fit(data)ax = cgram.plot(    figsize=(10, 8),    line_style=dict(color=ugg.COLORS[1]),    cluster_style={"color": ugg.COLORS[2]},)ax.yaxis.grid(False)sns.despine(offset=10)ax.set_title("Gaussian Mixture Model (scikit-learn)")

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

Подобная ситуация случается очень часто. Идеального случая не существует. В конечном счёте нам необходимо принять решение об оптимальном количестве кластеров. Clustergam даёт нам дополнительные сведения о том, что происходит между различными вариантами, как они расходятся. Можно сказать, что вариант с четырьмя кластерами в данных Iris не помогает. Также можно сказать, что пингвины Палмера могут быть сложными для кластеризации с помощью k-средних, что нет решающего правильного решения. Кластерограмма не даёт простого ответа, но она даёт нам лучшее понимание, и только от нас зависит, как мы её [кластерограмму] интерпретируем.

Установить clustergram можно с помощью conda install clustergram -c conda-forge или pip install clustergram. В любом случае вам всё равно придётся установить бэкенд кластеризации: либо scikit-learn, либо cuML. Документация доступна здесь, а исходный код здесь, он выпущен под лицензией MIT.

Если вы хотите поиграть с примерами из этой статьи, блокнот Jupyter находится на GitHub. Вы также можете запустить его в среде interactive binder в браузере. Более подробную информацию можно найти в блоге Тала Галили и оригинальных статьях Матиаса Шонлау.

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

Узнайте, как прокачаться и в других специальностях или освоить их с нуля:

• Профессия Data Analyst

• Курс по Data Engineering

Бесспорно, интерфейс Gravitee представляет достаточно наглядные и удобные средства визуализации работы шлюзов Gravitee. Но в любом случае, возникает потребность предоставить доступ к этим инструментам службе мониторинга, владельцам или потребителям API и при этом они могут находится вне закрытого контура, в котором расположен менеджер API. Да и иметь всю доступную информацию по различным API на одном экране всегда удобнее.
Видеть происходящее на шлюзах, при этом не вдаваясь в особенности пользовательского интерфейса Gravitee, а администраторам - не тратить время на создание пользователей и разделение ролей и привилегий внутри Gravitee.
На Хабре уже была пара статей посвященных APIM Gravitee, тут и тут. По этому, в своей заметке, буду подразумевать, что читатель уже знаком с процессом установки/настройки APIM Gravitee и Grafana, рассмотрю только процесс настройки их интеграции.

Почему нельзя пойти простым путём?

По умолчанию, хранилищем для аналитики Gravitee является ElasticSearch. Информация накапливается в четырёх различных индексах, с посуточной разбивкой:

• gravitee-request-YYYY.MM.DD - здесь хранится информация по каждому запросу (аналог access.log в nginx). Это наша основная цель;

• gravitee-log-YYYY.MM.DD - здесь уже хранится более подробная информация о запросе (при условии, что включена отладка, см. рисунок ниже). А именно полные заголовки запросов и ответов, а также полезная нагрузка. В зависимости от настроек, логироваться может как обмен между потребителем и шлюзом, так и/или шлюзом и поставщиком API;

  Экран включения/отключения расширенного логирования

• gravitee-monitor-YYYY.MM.DD - этот нас не интересует;

• gravitee-health-YYYY.MM.DD - этот нас не интересует.

И казалось бы, что может быть проще: подключай ElasticSearch в качестве источника данных в Grafana и визуализируй, но не всё так просто.
Во первых, в индексе хранятся только идентификаторы объектов, т.е. человеко-читаемых имён поставщиков и потребителей, вы там не увидите. Во вторых, получить полную информацию соединив данные из двух источников непосредственно в интерфейсе Grafana, крайне проблематично. Gravitee хранит информацию о настройках и статистику своей работы в разных местах. Настройки, в MongoDB или PostgreSQL, по сути статическая информация. Таким образом в одном месте у нас (в терминах Grafana) - таблица, в другом - временной ряд.

B как же быть?

Большим преимуществом СУБД PostgreSQL является богатый набор расширений для работы с внешними источниками данных, в том числе и с ElasticSearch (тут). Благодаря этому интеграция сводится к тому, что Grafana общается с единственным источником данных - СУБД PostgreSQL, которая в свою очередь получает данные из ElasticSearch и обогащает их информацией и делает читаемой для администратора или любого другого бенефициара.
Схематически это будет выглядеть следующим образом (рисунок ниже).

Ну что же, за дело!

Все ниже описанные действия актуальны для следующей конфигурации: CentOS 7, APIM Gravitee 3.6, СУБД PostgreSQL 11, ElasticSearch 7.+

Начнём с интеграции PostgreSQL и ElasticSearch. Сам процесс интеграции достаточно прост и делится на следующие шаги:

1. Устанавливаем расширение multicorn11 и если не установлен pip, то ставим и его:

   yum install multicorn11 python3-pip
   

2. Далее из pip-репозитория, устанавливаем библиотеку python3 для работы с ElasticSearch:

   pip3 install pg_es_fdw
   

3. Далее, переходим к настройке PostgreSQL. Подключаемся целевой БД и добавляем расширение multicorn и подключаем необходимую библиотеку:

   GRANT USAGE on FOREIGN DATA WRAPPER multicorn TO gatewaytest;GRANT USAGE ON FOREIGN SERVER multicorn_es TO gatewaytest;
   

    CREATE EXTENSION multicorn; CREATE SERVER multicorn_es FOREIGN DATA WRAPPER multicorn  OPTIONS (wrapper 'pg_es_fdw.ElasticsearchFDW');
   

4. Выдаём права, непривилегированному пользователю. В нашем случае это logreader:

   GRANT USAGE on FOREIGN DATA WRAPPER multicorn TO logreader;GRANT USAGE ON FOREIGN SERVER multicorn_es TO logreader;
   

5. Для удобства, создадим отдельную схему logging, владельцем которой будет наш пользователь logreader:

   CREATE SCHEMA logging AUTHORIZATION logreader;
   

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

   CREATE TABLE logging.requests (  id varchar(36),  "@timestamp" timestamp with time zone,  api varchar(36),  "api-response-time" int,  application varchar(36),  custom json,  endpoint text,  gateway varchar(36),  "local-address" varchar(16),  method int,  path text,  plan varchar(36),  "proxy-latency" int,  "remote-address" varchar(16),  "request-content-length" int,  "response-content-length" int,  "response-time" int,  sort text,  status int,  subscription varchar(36),  uri text,  query TEXT,  score NUMERIC) PARTITION BY RANGE("@timestamp");
   

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

7. Для подключения и отключения индексов, создадим небольшой shell-скрипт и будем запускать его раз в сутки через cron:

   #!/bin/shNEWPART=${1:-$(date +'%Y.%m.%d')}OLDPART=$(date --date='14 days ago' +'%Y.%m.%d')curl http://gateway.corp/testpsql gateway -U logreader -c "CREATE FOREIGN TABLE logging.\"requests_${NEWPART}\"PARTITION OF logging.requests   FOR VALUES FROM ('${NEWPART} 00:00:00') TO ('${NEWPART} 23:59:59')SERVER multicorn_esOPTIONS (host 'els-host',  port '9200',  index 'gravitee-request-${NEWPART}',  rowid_column 'id',  query_column 'query',  query_dsl 'false',    score_column 'score',  sort_column 'sort',  refresh 'false',  complete_returning 'false',  timeout '20',  username 'elastic-ro',  password 'Sup3rS3cr3tP@ssw0rd');"    psql gateway -U gatewaydev -c "drop foreign table logging.\"requests_${OLDPART}\""
   

   Немного пояснений:

   • NEWPART - текущая дата, для формирования имени партиции , при подключении нового индекса из ElasticSearch;

   • OLDPART - дата истекшего, неактуально индекса, здесь это 14 дней (определяется исходя из настроек ES Curator). Удалять партиции, ссылающиеся на несуществующие - обязательно. В противном случае запросы, к родительской таблице, будут прерываться с ошибками;

   • Вызов 'curl http://gateway.corp/test', необходим для того, что бы создавался индекс текущего дня, так как он создаётся в момент первого обращения к любому поставщику API. Если его не создать, то это будет приводить к ошибке, описанной выше. Такая проблема больше актуальна для тестовых стендов и стендов разработки;

   • Затем, создаём партицию на индекс текущего дня;

   • И на последнем шаге - удаляем неактуальный индекс.

   • Проверяем что всё работает

     TABLE logging.requests LIMIT 1;
     

     Если всё правильно, то должны получить похожий результат

   -[ RECORD 1 ]-----------+-------------------------------------id                      | 55efea8a-9c91-4a61-afea-8a9c917a6133@timestamp              | 2021-05-16 00:00:02.025+03api                     | 9db39338-1019-453c-b393-381019f53c72api-response-time       | 0application             | 1custom                  | {}endpoint                | gateway                 | 7804bc6c-2b72-497f-84bc-6c2b72897fa9local-address           | 10.15.79.29method                  | 3path                    | plan                    | proxy-latency           | 2remote-address          | 10.15.79.27request-content-length  | 0response-content-length | 49response-time           | 2sort                    | status                  | 401subscription            | uri                     | /testquery                   | score                   | 1.0
   

Рисуем графики

И вот мы подошли к тому, ради чего всё и делалось - визуализируем статистику Gravitee. Благодаря тому, что для доступа к аналитике используется единая точка входа, а именно СУБД PostgreSQL, это даёт дополнительные возможности. Например, выводить статическую информацию: количество поставщиков, количество потребителей и их статусы; количество и состояние подписок; параметры конфигурации для поставщика и многое другое, наряду с динамическими данными.
В том числе хотелось бы отметить, что у поставщиков и потребителей имеется раздел Metadata, которые можно заполнять кастомными данными и так же выводить в дашборды Grafana.

Вот тут:

А вот так это можно отобразить в Grafana:

SELECT  name "Наименование",  value "Значение"FROM  metadataWHERE  reference_id='${apis}'

Пример комплексного экрана

APIs (статика) - общее количество поставщиков и количество активных.

SELECT COUNT(*) AS "Всего" FROM apis;SELECT COUNT(*) AS "Активных" FROM apis WHERE lifecycle_state='STARTED';

Для Applications, запросы составляются по аналогично, только из таблицы applications

API Hits - количество вызовов по каждому поставщику. Тут уже немного по сложнее

SELECT  date_trunc('minute',"@timestamp") AS time,  apis.name,ee с Grafana  COUNT(*)FROM  logging.requests alJOIN  apis ON al.api = apis.idWHERE  query='@timestamp:[$__from TO $__to]'GROUP BY 1,2

Average response time by API - среднее время ответа, по каждому поставщику считается аналогичным способом.

SELECT  date_trunc('minute',"@timestamp") AS time,  apis.name,  AVG(al."api-response-time")FROM  logging.requests alJOIN  apis ON al.api = apis.idWHERE  query='@timestamp:[$__from TO $__to]'GROUP BY 1,2

Еще один интересный показатель Hits, by gateways, это равномерность распределения запросов по шлюзам. Считается так:

SELECT  date_trunc('minute',"@timestamp") as time,  al."local-address",  COUNT(*)FROM  logging.requests alWHERE  query='@timestamp:[$__from TO $__to]'GROUP BY 1,2

Заключение

Приведённое выше решение, по моему субъективному мнению, нисколько не уступает стандартным средствам визуализации APIM Gravitee, а ограничено лишь фантазией и потребностями.
Учитывая то, что Grafana, обычно является центральным объектом инфраструктуры мониторинга, то преимущества такого решения очевидны: более широкий охват, более высокая плотность информации и простая кастомизация визуальных представлений.

P.S.

В ближайшее время, планируется ещё статья по интеграции Gravitee с ActiveDirectory. Процесс достаточно прост, но как всегда, есть нюансы.

Конструктивная критика, пожелания и предложения приветствуются!

_{Соавтор:Кузнецова Юлия Андреевна - UX-писатель Экосистемы РСХБ}

Каким должен быть дизайнер в банке, чтобы и продукт хороший создавал, и коллеги не жаловались. Смотрим через призму софт-навыков вместе с UX-дизайнерами РСХБ.

#1 Коммуникабельность: не просто коллега человек

Дизайнер отличная профессия для интровертов. Рисуешь интерфейсы, слушаешь музыку, ни с кем не общаешься На самом деле, нет!

Работая в банке, ты участвуешь в большом количестве созвонов, слушаешь коллег, предлагаешь свои идеи. Споришь. Много споришь. Убеждаешь. Идешь на компромисс. Вместе с тобой над продуктом работают много других сотрудников продакт-менеджеры, аналитики, разрабы, редакторы, юристы ты должен уметь услышать каждого и найти общий язык.

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

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

#2 Презентация: готов объяснить свою работу

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

Крутой дизайнер умеет переключить фокус с творческого и ранимого художника на хладнокровного адвоката. Если кто-то говорит: Это плохо сделано, это не значит, что ты плохой. Это значит, что коллега видит в данном решении слабые стороны. Продвинутый дизайнер отключает эмоции и включает голову. Объясняет и обосновывает свою позицию. Задаёт вопросы, внимательно слушает, делает все, чтобы добраться до истины и хорошего продукта.

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

#3 Аналитика: думает не только о красоте

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

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

С одной стороны, ты должен на интуитивном уровне понимать психологию взаимодействия с интерфейсом. А с другой всегда остаешься скептиком и тестировать свои гипотезы. Порой пользователи воспринимают ту или иную функцию абсолютно непредсказуемо и не понимают твоей задумки. Одна только фраза Мне нравится. Думаю пользователям понравится тоже не поможет. Нужна аналитика. Много аналитики.

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

#4 Понимание бизнес-процесса: переводит с банковского языка на человеческий

Юный дизайнер: мастерски работает в Фигме, любит компоненты, боготворит филигранную верстку.

Продвинутый дизайнер: понимает, как устроен бизнес и финансы.

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

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

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

#5 Любознательность: постоянно учится новому

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

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

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

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

Автор статьи реализовал анализ опционов с помощью телеграмм-бота, что дает возможность быстрого получения информации об инвестиционном портфеле на смартфоне.

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

Запоминание состояния бота между вебхуками

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

$id_init = file_get_contents('php://input');$id=sbs($id_init, '"from":{"id":',',"is_bot":'); //в эту переменную записываем уникальный номер пользователяfunction sbs ($str,$m1,$m2){ //из строки str возвращает подстроку между двумя метками-словами m1 и m2$p1=strpos($str,$m1)+strlen($m1); //длина слова-метки слева$p2=strpos($str,$m2);return substr($str,$p1,$p2-$p1);}

Для каждого пользователя строится следующая структура данных:

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

2. файл данных по составу портфеля ( идентификационный номер ценной бумаги, тип, количество, цена, дата экспирации)

3. файлы для построения графика

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

Построение графика для анализа портфеля

График строится по состоянию на момент экспирации опционов, что дает возможность упрощенно представить график каждого инструмента в виде y=kx+b,

где y это размер прибыли/убытка на момент экспирации

х стоимость базового актива

После записи данных в портфель создаются файлы с расширением png с помощью библиотеки GD в несколько этапов:

1) определение Х-координат всех точек перелома на графике (очевидно, что сумма всех линий даст ломанную кривую)

2) определение масштаба изображения (средняя всех координат точек перелома по оси Х, максимальное отклонение от среднего, и максимальный размер по оси Y)

3) создание ассоциативного массива точек, в котором координата X это ключ, координата Y величина, для всей цифровой полуплоскости:

$typ опцион колл, пут или фьючерс

$q количество ( отрицательное продажа)

$cena цена приобретения ценной бумаги

$strike страйк для опционов

$x0 начальная координата по оси Х

$sx масштаб по оси Х

function pparr($typ, $q, $cena, $strike,$x0,$sx){ //функция выдает одномерный массив - координаты x=>y точек по //типу цб, направлению (покупка продажа), цене приобретения и страйку (для опционов)  if ($q<0) { $q=-$q; $drct='-'; }  else $drct='+'; $a=array(); $b=array(); $delta=$sx; //расстояние между точками равно масштабу $scalx for ($i=0;$i<740;$i++){ //кол во точек 740 определено заранее $xkk=$x0+$delta*$i; //значение по оси X if ($typ=='fut') { if ($drct=='+') $a[$xkk]=($xkk-$cena)*$q; else $a[$xkk]=(-$xkk+$cena)*$q; } if ($typ=='call'){ if ($drct=='+') { if ($xkk<=$strike) $a[$xkk]=-$cena*$q; else $a[$xkk]=$q*($xkk-$strike-$cena);} else { if ($xkk<=$strike) $a[$xkk]=$q*$cena; else $a[$xkk]=(-$xkk+$strike+$cena)*$q;}  } if ($typ=='put'){ if ($drct=='+') { if ($xkk<=$strike) $a[$xkk]=(-$xkk+$strike-$cena)*$q; else $a[$xkk]=-$q*$cena;} else { if ($xkk<=$strike) $a[$xkk]=($xkk-$strike+$cena)*$q; else $a[$xkk]=$cena*$q;}  } $b[(string)$xkk]=(string)$a[$xkk]; }return $b; };

4) создание файлов изображений, при этом одновременно строится график для каждой строки из портфеля ( зеленый цвет) и результирующий для портфеля ( красный цвет). Кроме того, know-how заключается в том, что одновременно строится еще четыре изображения для увеличения/уменьшения изображения по оси Х и по оси Y. За счет этого достигается эффект работы он-лайн с клавишами X+,X-,Y+,Y- под графиком. Таким образом, для каждого пользователя в каждый момент времени существует пять файлов изображения.

Анализ рисков опционного портфеля

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

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

Гарантийное обеспечение для фьючерсов считается по формуле, которая предложена самой Московской биржей:

ГО=БГО+(Цена-Расчетная_Цена)*БП;

где Расчетная_Цена определяется по результатам клиринга и практически равна цене базового актива на момент клиринга, определяется в пунктах

БП стоимость пункта в рублях, зависит, в частности, от курса доллара,

БГО базовое гарантийное обеспечение, определяется биржей, может произвольно увеличиваться, например, перед праздниками

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

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

Вот как выглядит покупка колла в нашем боте:

Гарантийное обеспечение=23.

Что нам показывает график: если стоимость базового актива ( в данном случае стоимость фьючерса на курс рубля к доллару) [вечером] 03.06.2021 будет 75000 и ниже, то наш убыток составит 23 . При повышении этой стоимости до 75023 мы выйдем в безубыток, при дальнейшем росте получим прибыль.

Что мы имеем с точки зрения риска: не при каких обстоятельствах наш убыток не превысит сумму 23. Следовательно, наш опцион совершенно не похож на фьючерс, и в расчете ГО мы можем записать просто сумму 23.

Покупка пута примерно та же картина.

Продажа пута.

Гарантийное обеспечение= 5436.

Можно убедиться, что при стоимости базового актива выше 72750 мы имеем прибыль 44. При снижении стоимости БА до 72706 мы выходим в ноль. При дальнейшем падении стоимости БА наш убыток НИЧЕМ НЕ ОГРАНИЧИВАЕТСЯ.

С точки зрения рисков это фьючерс, купленный по цене 72706. Подставляем это число в формулу ГО для фьючерса и получаем ГО для опциона! Этот ГО достаточно велик (5436), но может превратится в прибыль в течение нескольких дней.

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

С продажей колла будет аналогичная ситуация.

А если одновременно продать пут и кол?

Гарантийное обеспечение не изменилось!

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

с наибольшим ГО , оно и будет мерилом риска.

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

Интерфейс телеграмм-бота

На экране отображаются следующие группы данных:

• подробная инструкция по работе с ботом

• текущие котировки ближайших фьючерсов по трем базовым активам

• таблица Портфель, отображающая состав портфеля, сделанная с помощью интерфейса телеграмм ( можно редактировать)

• значение гарантийного обеспечения

• таблица, дублирующая состав портфеля, построенная как изображение формата png, которую можно копировать и сохранять

• нижняя клавиатура, которую можно скрывать, с кнопками: Добавить позицию в портфель, Анализ портфеля с помощью графика, Обновить котировки

Таблица Портфель создана с помощью InlineKeyboard.

При нажатии на клавиши в этой таблице происходят следующие действия:

• клавиша в столбце К-во редактирует количество ценных бумаг, соответствующее строке, в которой нажата клавиша

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

• клавишей Х можно удалить строку.

Прочие клавиши в таблице заблокированы и служат только для передачи информации.

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

Заключение

В телеграмм-боте реализовано запоминание портфеля ценных бумаг для каждого пользователя. Под ценными бумагами понимаются опционы и фьючерсы, базовым активом для которых являются: курс рубля к доллару (Si), стоимость нефти брент (BR), а также индекс РТС (RI). Это самые высоколиквидные деривативы на московской бирже.

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

С помощью телеграмм-бота можно анализировать опционный портфель на графике прибылей/убытков (P/L график).

Протестировать телеграмм-бота можно по ссылке t.me/@test09062020bot. Или попробовать найти в телеграмме по названию опционный портфель.

Если ваша жизнь DBA, сопровождающего PostgreSQL, наполнена вопросами "а почему так медленно?" и "как сделать, чтобы запрос не тормозил?", наш сервис анализа и визуализации планов запросовexplain.tensor.ru сделает ее немного легче за счет привлечения коллег и обновленных подсказок.

"Ландон из зе кепитал оф Грейт Британ"

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

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

Подсказки к плану

Про базовый набор подсказок и способов ими воспользоваться я уже рассказывал в статье "Рецепты для хворающих SQL-запросов" - теперь мы сделали их еще больше и удобнее!

• цветовая и текстовая кодировка разных видов

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

• несколько подсказок у одного узла

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

• на самом видном месте

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

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

Масштабируемая диаграмма

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

Пользуйтесь! Возникнут идеи или замечания - прошу в комментарии.

В школе все мы решали задачки вида едет из пункта А в пункт Б. Речь преимущественно шла о скорости и времени как быстро доберётся транспортное средство? Реальность, однако, подбрасывает задачки значительно интереснее: Существует масштабная ритейл-сеть по продаже товаров, которой необходимо, чтобы огромное количество номенклатурных позиций доезжало в каждый из 17000 магазинов, расположенных на половине площади самой большой страны в мире, вовремя и в нужном количестве. Для решения такой задачи в X5 Group существует ряд реализованных решений, и одним из самых важных является процесс автозаказа товаров.

Техническую поддержку этого направления в X5 Group обеспечивает команда 2-SAP Логистики.

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

Автозаказ это комплекс процессов управления запасами и заказами

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

Планирование заказа.

Формирование заказа.

Отправка заказа.

Экономическое обоснование поддерживаемого уровня.

Контроль за состоянием запасов.

Управление ассортиментом.

Для этого используются разные механизмы прогнозирования в зависимости от формата магазина и гибких настроек в системе.

Прогноз продаж, построенный SAP это прогноз спроса на период 42 дня, который строится на исторических данных по продажам. Период анализа продаж и модель прогноза определяется из настройки Профиля прогноза.

В периоде анализа продажи могут быть Плановые регулярные продажи и Внеплановые промо продажи.

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

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

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

Около 04:30 эти данные по потребности из JDA поступают в SAP, где обрабатываются фоновым заданием с интервалом запуска каждые 15 минут, в результате чего создаются Автозаказы.

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

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

Через автозаказ пополняется до 80% основного ассортимента магазинов, а ручное пополнение выполняется для товаров in-out.

Схема реализации товара:

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

Так выглядит децентрализованная цепочка:

- в 00:00 стартует задание по АЗ, в период его работы по нашему товару рассчиталась потребность, сформировался автозаказ, и он бы отправлен из SAP в программу магазина GK(Пятёрочка), там у сотрудников магазина есть время отредактировать рассчитанное кол-во (в большую или меньшую сторону) до наступления времени автосогласования (например, до 11 утра, время зависит от тайминга поставщика, т.е. до какого момента он должен получить заказ, чтобы собрать и вовремя привезти заказ в магазин);

- после внесения корректировок эти изменения отправляются в SAP, и дальше происходит процесс автосогласования заказ согласовывается и автоматически отправляется поставщику по почте или его провайдеру (в зависимости от того, по какой схеме поставщик работает). Например, если поставщик работает по EDI (через электронный документооборот), то после доставки заказа в его систему, нам провайдер передает информацию о размещении заказа в своем ПО и в системе ER2 заказ переходит в статус Размещен;

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

Автозаказ (Ночной профиль) фоновое задание, имеющее внутри себя несколько шагов, выполняющих определённые функции. Результатом задания является создание заказов на закупку товаров. Запускается параллельно по магазинам.

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

• Продажи за период в прошлом

• Прогноз на период в будущее (на основании продаж строится прогноз SAP)

• Остатки на объекте получателе (остаток в SAP на момент расчёта)

• Открытые заказы (заказы на поставку и возвраты)

Продажи приходят из магазинов каждый вечер до расчёта фонового задания Автозаказа.

Факт получения продаж запускает расчет Автозаказа. Если до 03:00 продажи не получены SAP ERP, то происходит безусловный запуск задания Автозаказа.

Автозаказ формируется строго по графику заказа/поставки согласованному с поставщиком, так категории FRESH и ULTRA FRESH, как правило, заказываются и поставляются в магазины ежедневно.

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

Для того, чтобы заказать оптимальное кол-во товара, учитывается прогноз, который строится в SAP в момент расчета Автозаказа.

Если по этому товару уже есть запас или размещенный заказ, который еще не закрыт, то система также учтёт это при определении размера финального заказа.

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

Запускается задание ежедневно в 00:00. Первыми выполняются расчеты по регионам Сибирь и Урал, далее по регионам Волга, ЮГ, Москва, Центральные регионы и Северо-запад. Последовательность выполнения крайне важна, т.к. наши магазины расположены в разных часовых поясах, и пока в Москве все сладко спят, в Сибири уже в разгаре рабочий день.

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

Для этого мы реализовали различные системы мониторинга, инструменты для анализа, а в случае необходимости и возможность отправки в магазины резервных шаблонов заказа.<o:p>

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

Также предусмотрен инструмент для выгрузки резервных заказов, рассчитанных накануне, в случае сбоев в системе JDA по прогнозированию и пополнению магазинов больших форматов Перекрёсток и Карусель.