Я пытаюсь создать какое-то упрощенное приложение Oracle Crystal Ball для своих геологических исследований, которое будет использовать значения P90 (90% достоверности) и P10 (10% достоверности) в качестве входных данных и обратное распределение различных вероятностных сценариев. Похоже на распределение Монте-Карло. Я новичок в Python, только недавно начал, кстати:)
Эта тема будет разделена на четыре ключевые части:
- Общее описание объема работ.
- Псевдокодирование (хотя никогда раньше не пытался).
- Фактический код Python.
- Причина, по которой я здесь, или проблемы с логикой / кодом.
ЧАСТЬ 1. Общее описание объема работ.
Для простоты предположим, что у нас есть только три категории, каждая с параметрами P90 и P10 без каких-либо шагов междуих:
- cat_1: [1, 2]
- cat_2: [2, 4]
- cat_3: [3, 6]
Используя декартово произведение, мы получаем следующие 8 списков с возможными сценариями:
- [1, 2, 3], [1, 2, 6], [1, 4, 3], [1, 4, 6], [2, 2, 3], [2, 2, 6], [2, 4, 3], [2, 4, 6]
Умножение параметров в каждом списке приводит к следующим продуктам:
- [6, 12, 12,24, 12, 24, 24, 48]
Измерение частоты каждого продукта приводит к:
{6: 1, 12: 3, 24: 3, 48: 1} или с учетом процентов:
{6: 12,5%, 12: 37,5%, 24: 37,5%, 48: 12: 5%,}, что означает, что вероятность появления 12 или 24 выше, чем 6 или 48.
Вот результат, который я хотел бы получить: знаявероятность того, что продукты смогут получить средние, срединные и модовые значения.
Трудная часть для моего оборудования - это огромное количество возможных сценариев в реальном случае. Всего имеется шесть категорий с небольшими шагами между значениями P90 и P10. С учетом метрической системы диапазон значений P90 и P10 может быть следующим:
- площадь квадрата: 0,01 - 100,00 км2, шаг 0,01;
- толщина слоя: 0,10 - 100,00 м,шаг 0,1;
- пористость: 0,01-1,00 пу, шаг 0,01;
- насыщение: 0,01-1,00 пу, шаг 0,01;
- давление: 1 - 2000 атм, шаг 1атм;
- поверхность: 0,01 - 1,00 пу, шаг 0,01.
Обычно в реальных примерах используются гораздо более узкие диапазоны, скажем, 0,1 - 2,0 км2 дляплощадь квадратного сечения, 1 - 10 м для толщины, 8 - 15 для пористости и т. д. Тем не менее, даже в этом случае это звучит как количество возможных сценариев "google" с учетом упомянутых шагов. В результате я получаю следующее уведомление, которое является ключевой проблемой:
Процесс завершен с кодом выхода 137 (прерван сигналом 9: SIGKILL).
Это происходит, когда общий объем вычислений превышает ~ 10 мм и ~ 1 минуту (проверено экспериментально, следовательно, цифры приблизительны).
ЧАСТЬ 2. Псевдокодирование.
Хорошая практика говорит, что не следует абстрагироваться при псевдокодировании, однако у меня нет опыта в этой области, поэтому я буду стараться изо всех сил.
User inputs minimum possible values (P90) for total 6 categories
User inputs maximum possible values (P10) for total 6 categories
Total 6 list are created (square area, layer thickness, porosity etc.), 1 per each category that contain a range of possible values and indicated step (P90_category1, P10_category1, step1)
Use a Cartesian product to create a list_of_tuples with possible scenarios
Convert list_of_tuples to the list_of_lists
Create empty_list
for each element in the list_of_lists:
calculate its product
append to the empty_list
Round values in the empty_list
Create a dictionary that counts similar values in the empty_list
Calculate a probability of each value according to its repetition frequency in the dictionary
Вот и все. Некоторые базовые статистические данные и графики также применяются, однако это не является ключевым моментом здесь.
ЧАСТЬ 3. Фактический код Python.
При первых значениях P90 (90%достоверность):
P90_area = float(input('P90 area: '))
P90_thickness = float(input('P90 thickness: '))
P90_porosity = float(input('P90 porosity: '))
P90_saturation = float(input('P90 saturation: '))
P90_pressure = float(input('P90 pressure: '))
P90_surface = float(input('P90 surface: '))
Затем значения P10 (достоверность 10%):
P10_area = float(input('P10 area: '))
P10_thickness = float(input('P10 thickness: '))
P10_porosity = float(input('P10 porosity: '))
P10_saturation = float(input('P10 saturation: '))
P10_pressure = float(input('P10 pressure: '))
P10_surface = float(input('P10 surface: '))
Создание диапазона значений от P90 до P10 с определенным шагом
area_values = np.arange(P90_area, P10_area + 0.01, 0.01)
thickness_values = np.arange(P90_thickness, P10_thickness + 0.1, 0.1)
porosity_values = np.arange(P90_porosity, P10_porosity + 0.01, 0.01)
saturation_range = np.arange(P90_saturation, P10_saturation + 0.01, 0.01)
pressure_range = np.arange(P90_pressure, P10_pressure + 1, 1)
surface_range = np.arange(P90_surface, P10_surface + 0.01, 0.01)
Объединить все списки в декартово произведение (т.е. [(area1, толщина1, пористость1), (area1, толщина1, пористость2) и т. Д.]):
list_of_tuples = list(itertools.product(area_values, thickness_values, porosity_values, saturation_range, pressure_range, surface_range)
Преобразовать список кортежей в список списков:
list_of_lists = [list(elem) for elem in list_of_tuples]
Создать список с умноженными значениями и отсортировать их (np.prod возвращает продукт для каждого списка):
multiplied_values = []
for i in list_of_lists:
i = np.prod(np.array(i))
multiplied_values.append(i)
multiplied_values = sorted(multiplied_values)
Значения округлости:
rounded_values = [float(Decimal('%.2f' % elem)) for elem in multiplied_values]
Создать словарь, который подсчитывает все похожие / уникальные объекты:
counts = Counter(rounded_values)
Вычислить вероятность путем деления значения на общее количество элементов в списке:
probability_mass = {k: v/total for k, v in counts.items()}
Этоработает, здесь идет простая статистика и диаграмма для конкретного случая:
- Всего вычислений: 4899510
- P90 составляет: 5,60
- P10 составляет: 43,41
- P50 (значение с максимальной вероятностью): 15,24
- Среднее значение: 23.80
Рис. Диаграмма распределения вероятностей
Первый вопрос имеет решающее значение, поскольку он блокирует вычисление большого стека данных:
ЧАСТЬ 4. Ключевые проблемы.
Q1. Ключевая проблема:
В результате я получаю следующее уведомление, которое является ключевой проблемой:
Процесс завершен с кодом выхода 137 (прерывается сигналом 9: SIGKILL).
По схожим темам, скорее всего, мой скрипт был убит ОС из-за чрезмерной загрузки процессора. Я проверил загрузку ЦП с помощью команды 'top' во время выполнения кода, и ЦП был загружен до 100%, когда он мог обрабатывать входные параметры, и в некоторых моментах до 110% при прерывании.
Характеристики: ноутбук Asus G531GU |i7-9750H CPU 2.60GHz |GeForce GTX 1660 TI, 6 Гб |16 Гб DDR4 |Ubuntu 18 |IDE сообщества PyCharm.
Вопрос: В любом случае, как я могу избавиться от этого прерывания и позволить скрипту работать столько, сколько необходимо? Я умею ждать столько, сколько нужно, чтобы получить правильное распределение для большого стека данных. Увеличение шага для каждого параметра является жесткой опцией, я бы не стал этого делать.
Q2 . Диаграмма распределения вероятностей не похожа на классическое нормальное распределение, тогда как разница между максимально вероятными и средними значениями значительна. Как вы думаете, могут ли быть какие-либо проблемы с логикой кода?
PS Я понимаю, что этот скрипт выглядит довольно ухабистым, надеюсь, ваши глаза не будут кровоточить)