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

Часто я сталкиваюсь с ситуацией, когда оценка может привести к вмешательству в дорогостоящее оборудование. Баланс между надежностью оборудования и капитальными решениями — дело тонкое. Здесь может пригодиться байесовский анализ. Оставшаяся часть этого поста посвящена «реальному» бизнес-кейсу и может служить введением в использование байесовского анализа в контексте IIoT. Будет некоторый код, однако я буду максимально абстрагироваться от математики байесовской оценки.

Какова цель этого анализа

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

Почему это важно?

Функция потока паровой турбины — параметр, по которому можно оценить изменения в паровом тракте турбины, т.е. изменения/повреждения лопаток. В этом случае функция потока строится с течением времени. Устойчивые изменения функции потока обычно указывают на необходимость вмешательства и замены лезвий. Хотя предметом данного анализа является турбина, этот метод применим к любому типу оборудования и измерений: отклонение от кривой насоса или параметр Отношение скоростей/безразмерность для аортальных клапанов. Подводя итог, мы заинтересованы в понимании и количественной оценке изменений технических параметров и, следовательно, работоспособности оборудования с течением времени.

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

Подождите, но почему байесовский

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

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

Но почему бы просто не использовать обычный t-тест нулевой гипотезы?

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

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

Приступаем к делу

Байесовский анализ здесь выполняется с использованием pymc3, в будущем будут обновления с использованием TensorFlow Probability (TFP) и Pyro.

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

plt.style.use('fivethirtyeight')
flowdf.hist('Flow_function', by='group', figsize=(12, 8))

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

Немного практического опыта байесовского оценивания

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

В классическом машинном обучении у нас есть y=f(X), где f(X) может быть нейронной сетью, линейной регрессией, SVM и т. д. Мы начнем с определения f() иначе известный как модель, затем мы используем алгоритм оптимизации вместе с целевой функцией, чтобы найти параметры модели, алгоритмом оптимизации может быть градиентный спуск. в основном то, что я только что описал, это то, что происходит за кулисами, когда вы вводите следующие две строки в scikit-learn:

model= LinearRegression()
model.fit(X,y)

В байесовском анализе мы начинаем с определения модели для наших параметров, мы называем ее априорной. Обычно это какое-то распределение, например. нормальное распределение или распределение Пуассона. Затем мы определяем модель для наших данных по параметрам, это может быть линейная регрессия или нейронная сеть. Затем мы «подгоняем» данные, запуская тысячи образцов для параметров из предыдущего распределения. Каждая выборка будет исходить из предыдущего распределения, а также будет соответствовать данным, которые мы моделируем. Эта выборка обычно называется MCMC (цепь Маркова Монте-Карло). Эта выборка параметров по всему набору данных дает нам окончательный результат , апостериорный, который представляет собой модель с обновленными параметрами. Априорное и апостериорное распределения представлены в виде распределений, подобных нормальному распределению.

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

В двух словах, теорема гласит, что если у кого-то есть линейная модель: y=ax+b, мы начинаем с начального предположения о параметрах a и b (априор), эти параметры могут быть обновлены при представлении данных, конечным результатом является апостериор.

Сложность обычно возникает из-за кажущегося простым термина p(данные)(доказательства). Чтобы вычислить p(data), необходимо вычислить интеграл от p(data|q).p(q) по всей модели параметры. поэтому в простом y=ax+b интеграция будет выглядеть примерно так.

Любой, кто выполнял интеграцию двух или трех параметров, может подтвердить, насколько это сложно для некоторых функций и во многих случаях не поддается расчету. В нашем примере использования у нас есть довольно сложная функция/модель с 9 параметрами (объяснено позже). Интегрировать более 9 параметров непросто. Следовательно, мы используем MCMC в качестве численного метода для аппроксимации p(data). Помимо MCMC, существуют и другие способы создания апостериорной , но я не буду описывать их в этой статье.

Таким образом, с точки зрения непрофессионала, 1) мы начинаем наш анализ с функции/модели для анализа, точно так же, как мы делаем с классическим машинным обучением, скажем, LinearRegression() , 2) мы устанавливаем начальное предположение о параметрах модели (априор), 3 ) мы наблюдаем данные, а затем меняем диапазон возможных значений этих параметров с помощью алгоритма MCMC, чтобы получить апостериорную

Итак, есть 3 основных шага:

  1. Определите модель данных, которые будут генерировать апостериорную
  2. Определите предварительное распределение, из которого будут выбираться параметры модели a и b.
  3. Создайте новый достоверный диапазон параметров с помощью MCMC и, следовательно, получите апостериорную

Шаг 1: Определение модели данных, f(x)

Цель анализа — оценить, как изменилась функция потока между группами 1, 2, 3 и 4. Чтобы сделать эту оценку, мы можем сравнить статистику распределения, такую ​​как среднее значение, стандартное отклонение между каждой группой. Для выполнения этого сравнения я выбрал распределение Стьюдента-t в качестве описательной модели для распределений в каждой группе. Я мог бы выбрать распределение Гаусса, однако распределение Стьюдента-t более устойчиво к выбросам, присутствующим в данных. Важно отметить, что распределение Стьюдента -t здесь используется для описания распределения данных, это не t-распределение для t-теста . Иными словами,распределение Стьюдента-t — это f(x): модель, используемая для нахождения описательная статистика каждой группы.

Распределение Стьюдента-t имеет три параметра: среднее μ, точность (обратная дисперсия) λ и параметр степеней свободы ν:

Здесь предполагается, что переменная ν одинакова для всех 4 групп. Среднее значение и параметры дисперсии будут разными для каждой группы. Это приводит к модели с 9 параметрами (параметры среднего значения и стандартного отклонения для каждой из 4 групп и один ν для всех групп). Большие значения ν делают распределение похожим на нормальное, а малые значения приводят к тяжелым хвостам, т. е. выбросам.

Шаг 2: Определение исходных значений параметров модели, также известных как Prior

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

а. Средние параметры модели

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

import pandas as pd
import numpy as np
import pymc3 as pm
μ_m = flowdf['Flow_function'].values.mean()
μ_s = flowdf['Flow_function'].values.std() * 2
with pm.Model() as model:
group1_mean = pm.Normal("group1_mean", mu=μ_m, sd=μ_s)
group2_mean = pm.Normal("group2_mean", mu=μ_m, sd=μ_s)
group3_mean = pm.Normal("group3_mean", mu=μ_m, sd=μ_s)
group4_mean = pm.Normal("group4_mean", mu=μ_m, sd=μ_s)

б. Параметры стандартного отклонения модели

Предыдущие значения параметра стандартного отклонения каждой группы будут выбраны из равномерного распределения. Равномерное распределение будет иметь нижний предел (1/100) * (стандартное отклонение набора данных) и верхний предел (100 * стандартное отклонение набора данных):

σ_low = (μ_s/2)/100
σ_high = (μ_s/2)*100
with model:
group1_std = pm.Uniform("group1_std", lower=σ_low, upper=σ_high)
group2_std = pm.Uniform("group2_std", lower=σ_low, upper=σ_high)
group3_std = pm.Uniform("group3_std", lower=σ_low, upper=σ_high)
group4_std = pm.Uniform("group4_std", lower=σ_low, upper=σ_high)

в. Параметр степеней свободы модели (ν)

В оригинальной статье Крушке априорные значения для ν выбираются из экспоненциального распределения. В этом анализе также используется экспоненциальное распределение со средним значением 30; малое значение ν означает, что данные имеют тяжелые хвосты, большее значение ν означает, что данные нормально распределены. Использование экспоненциально распределенной выборки ν означает, что модель распределения Стьюдента, f(x), будет адаптироваться к диапазону от нормальных данных (без выбросов) до данных с тяжелыми хвостами (с выбросами). Для ν>30 модель примет форму нормального распределения, для ν‹30 модель примет форму нормального распределения с длинными хвостами.

with model:
ν = pm.Exponential("ν_minus_one", 1 / 29.0) + 1
# The way we arrive to the exponential function above is:
# p(ν|λ)=(1/λ)*exp[-(ν-1)/λ], ν≥1 and λ=29 
pm.kdeplot(np.random.exponential(30, size=10000), fill_kwargs={"alpha": 0.5})

Шаг 3: Получите апостериорную модель, model.fit()

До него мы устанавливали вероятностную модель. Этот шаг позволит подобрать модель с помощью MCMC и получить заднюю часть. Оценки интересов представляют собой сравнения между группами, потому что мы хотим быть уверены в влиянии эксплуатационных событий на функцию потока паровой турбины. Сравнение между группами будет происходить в виде разницы в средних значениях, разницы в стандартном отклонении и величине эффекта. Эти меры будут настроены как детерминированные узлы. Присвоение Deterministic сообщает PyMC, что мы хотим, чтобы выборочные значения возвращались как часть вывода.

Еще одним интересным показателем является величина эффекта.

Размер эффекта поможет нам понять значение изменения.

with model:
diff_of_means_pre = pm.Deterministic("difference of means pre", group2_mean - group1_mean)
diff_of_stds_pre = pm.Deterministic("difference of stds pre", group2_std - group1_std)
effect_size_pre = pm.Deterministic(
"effect size pre", diff_of_means_pre / np.sqrt((group2_std ** 2 + group1_std ** 2) / 2))
diff_of_means_post = pm.Deterministic("difference of means post", group3_mean - group2_mean)
diff_of_stds_post = pm.Deterministic("difference of stds post", group3_std - group2_std)
effect_size_post = pm.Deterministic(
"effect size post", diff_of_means_post / np.sqrt((group3_std ** 2 + group2_std ** 2) / 2))
diff_of_means_new = pm.Deterministic("difference of means new", group4_mean - group3_mean)
diff_of_stds_new = pm.Deterministic("difference of stds new", group4_std - group3_std)
effect_size_new = pm.Deterministic(
"effect size new", diff_of_means_new / np.sqrt((group4_std ** 2 + group3_std ** 2) / 2))

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

with model:
trace = pm.sample(draws=2000,tune=1000)

Чтобы лучше понять, как работает сэмплирование MCMC и пошаговый алгоритм NUTS (No-U-Turn Sampler), используемый под капотом PyMC, я рекомендую прочитать этот учебник.

Мы наносим апостериорные результаты для каждой группы. Графики суммируют апостериорные распределения для каждого параметра модели. Эти графики представляют нам интервал 95% и средние значения для каждого параметра.

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

pm.plot_posterior(
trace,
var_names=["difference of means pre", "difference of stds pre", "effect size pre",
"difference of means post", "difference of stds post", "effect size post",
"difference of means new", "difference of stds new", "effect size new"],
ref_val=0, color="#87ceeb")

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

При сравнении группы 2 с группой 3 (послеоперационное изменение) и группы 3 с группой 4 (событие «Новые операции»), представленной ниже, наблюдается очевидное изменение среднего значения, однако дисперсия также увеличилась. Величина эффекта, хотя и остается положительной, остается стабильной и значительно сокращается по сравнению с первоначальным изменением между группами 1 и 2. Это указывает на более медленное, но продолжающееся ухудшение характеристик паровой турбины.

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

Закрывать

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

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

Ссылки:

  1. Байесовская оценка заменяет t-критерий, Джон К. Крушке
  2. Байесовская оценка заменяет t-критерий, PyMC Notebook
  3. Байесовский анализ данных, 3-е издание, А. Гельман, Дж. Б. Карлин, Х. С. Стерн, Д. Б. Дансон, А. Вехтари, Д. Б. Рубин