Статистическая регуляризация некорректных обратных задач им. Турчина (часть 1)

    Привет, Хабр! Сегодня мы хотим рассказать, чем занимается лаборатория методов ядерно-физических экспериментов, входящая в JetBrains Research.

    Где JetBrains и где ядерная физика, спросите вы. Мы сошлись на почве любви к Kotlin, хотя в данном посте о нем речи не пойдет. Наша группа ориентируется на развитие методик анализа данных, моделирования и написание софта для ученых, и поэтому ориентирована на сотрудничество и обмен знаниями с IT-компаниями.

    В этой статье мы хотим поговорить о популяризуемом нами методе статистической регуляризации, предложенном В.Ф.Турчиным в 70-х годах XX века, и его реализации в виде кода на Python и Julia.

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


    Возникновение проблемы: зачем вообще кого-то регуляризировать?


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

    $f(y) = \int \limits_a^b dx K(x,y)\varphi(x) $

    Фактически это уравнение описывает следующее: наш измерительный прибор представлен здесь своей аппаратной функцией $K(x,y)$, которая действует на исследуемый спектр или иной входной сигнал $\varphi$, в результате чего исследователь наблюдает выходной сигнал $f(y)$. Цель исследователя — восстановить сигнал $\varphi$ по известным $f(y)$ и $K(x,y)$. Также можно сформулировать это выражение в матричной форме, заменив функции векторами и матрицами:

    $f_m = K_{mn}\varphi_n$

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

    $\varphi(x) = 2*N(2, 0.16) + N(4, 0.04)$

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

    $K_{mn} = \theta(x_m-y_n)$

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

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

    Восстанавливать сигнал мы будем методом наименьших квадратов:

    $\varphi^{МНК} = (K^TK)^{-1}K^Tf$

    И в результате получим:

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

    Но давайте сначала разберемся, из-за чего сия неудача постигла нас? Очевидно дело в ошибках измерений, но на что же они влияют? Дело в том, что еще Жак Адамар (тот самый, который добавил черточку в формулу Коши — Адамара) разделил все задачи на корректно поставленные и некорректно.

    Вспомнив классиков: «Бессмыслица — искать решение, если оно и так есть. Речь идет о том, как поступать с задачей, которая решения не имеет. Это глубоко принципиальный вопрос…» — мы не будем говорить о корректных задачах и сразу возьмемся за некорректные. Благо мы уже встретили такую: написанное выше уравнение Фредгольма является некорректной обратной задачей — даже при бесконечно малых флуктуациях во входных данных (а уж наши ошибки измерения далеко не бесконечно малые) решение уравнения, полученное точным аналитическим образом, может сколь угодно отличаться от истинного.

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

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

    Теоретическое описание статистической регуляризации


    Стратегия


    Сформулируем нашу задачу в терминах математической статистики: по известной реализации $f$ (которую мы измеряем в эксперименте) нам нужно оценить значение параметра $\varphi$. Функционал $\hat{S}$ вычисляющий $\varphi$ на основе $f$ мы будем называть стратегией. Для того чтобы определить, какие стратегии более оптимальны, введем квадратичную функцию потерь. Реальная функция потерь может быть любой, почему же мы выбираем именно квадратичную? Потому что любая функция потерь вблизи своего минимума может быть аппроксимирована квадратичной функцией:

    $L(\varphi,\hat{S}[f]) = ||\hat{\varphi}-\hat{S}[f])||_{L_2},$

    где $\hat{\varphi}$ — наилучшее решение. Тогда потери для выбранной нами стратегии задаются функцией риска:

    $R_{\hat{S}[f]}(\varphi) \equiv E[L(\varphi,\hat{S}[f])] = \int L(\varphi,\hat{S}[f])P(f|\varphi)df.$

    Здесь $P(f|\varphi)$ определяет плотность вероятности нашего ансамбля, по которому производится усреднение потерь. Этот ансамбль образован гипотетическим многократным повторением измерений $f$ при заданном $\varphi$. Таким образом, $P(f|\varphi)$ — это та самая известная нам плотность вероятности $f$, полученная в эксперименте.

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

    $r_{\hat{S}}(\varphi) \equiv E_{\varphi}E_{f}[L(\varphi,\hat{S}[f])|\varphi]$

    В этом случае оптимальная стратегия хорошо известна:

    $\hat{S}[f] = E[\varphi|f] = \int \varphi P(\varphi|f)d\varphi,$

    где апостериорная плотность $P(\varphi|f)$ определяется по теореме Байеса:

    $P(\varphi|f)= \frac{P(\varphi)P(f|\varphi)}{\int d\varphi P(\varphi)P(f|\varphi)}$

    Такой подход позволит определить дисперсию (корреляционную функцию) полученного решения:

    $D(x_1,x_2) = E[\varphi(x_1) - \hat{S}[f](x_1)][\varphi(x_2) - \hat{S}[f](x_2)]$

    Итак, мы получили оптимальное решение нашей задачи, введя априорную плотность $P(\varphi)$. Можем ли мы сказать что-либо о том мире функций $\varphi(x)$, который задается априорной плотностью?

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

    Именно в этом и заключается метод статистической регуляризации — регуляризация решения за счет введения дополнительной априорной информации о $\varphi(x)$. Если исследователь уже обладает какой-либо априорной информацией (априорной плотностью $P(\vec{\varphi})$), он может просто вычислить интеграл и получить ответ.

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

    Априорная информация


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

    Иначе говоря, он назначает более гладким $\varphi(x)$ более высокую априорную плотность вероятности. Так давайте попробуем ввести априорную вероятность, основанную на гладкости. Для этого вспомним, что введение априорной информации — это некоторое насилие над миром, принуждающее законы природы выглядеть удобным для нас образом.

    Это насилие следует свести к минимуму, и, вводя априорную плотность вероятности, необходимо, чтобы информация Шеннона относительно $\varphi(x)$, содержащаяся в $P(\vec{\varphi})$, была минимальной. Формализуя вышесказанное, выведем вид априорной плотности, основанной на гладкости функции. Для этого мы будем искать условный экстремум информации:

    $I[P(\vec{\varphi})] = \int \ln{P(\vec{\varphi})} P(\vec{\varphi}) d\vec{\varphi} \to min $

    При следующих условиях:

    1. Условие на гладкость $\varphi(x)$. Пусть $\Omega$ — некоторая матрица, характеризующая гладкость функции. Тогда потребуем, чтобы достигалось определенное значение функционала гладкости:

      $\int(\vec{\varphi},\Omega\vec{\varphi}) P(\vec{\varphi}) d\vec{\varphi} = \omega$

      Внимательный читатель должен задать вопрос об определении значения параметра
      $\omega$. Ответ на него будет дан далее по тексту.
    2. Нормированность вероятности на единицу: $ \int P(\vec{\varphi}) d\vec{\varphi} = 1$
      При этих условиях доставлять минимум функционалу будет следующая функция:

      $P_{\alpha}(\vec{\varphi}) = \frac{\alpha^{Rg(\Omega)/2}\det\Omega^{1/2}}{(2\pi)^{N/2}} \exp(-\frac{1}{2} (\vec{\varphi},\alpha\Omega\vec{\varphi})) $

      Параметр $\alpha$ связан с $\omega$, но поскольку у нас нет информации о конкретных значениях функционала гладкости, выяснять, как именно он связан, бессмысленно. Что же тогда делать с $\alpha$, спросите вы. Здесь перед вами раскрываются три пути:

      1. Подбирать значение параметра $\alpha$ вручную и тем самым фактически перейти к регуляризации Тихонова
      2. Усреднить (проинтегрировать) по всем возможным $\alpha$, предполагая все возможные $\alpha$ равновероятными
      3. Выбрать наиболее вероятное $\alpha$ по его апостериорной плотности вероятности $P(\alpha|\vec{f})$. Этот подход верен, поскольку дает хорошую аппроксимацию интеграла, если в экспериментальных данных содержится достаточно информации об $\alpha$.

    Первый случай нам мало интересен. Во втором случае мы должны посчитать вот такой вот страшненький интеграл:

    $\left\langle \varphi_i \right\rangle = \frac{\int d\varphi\, \varphi_i P(f|\varphi) \int\limits d\alpha\,P(\alpha) \alpha^{\frac{Rg(\Omega)}{2}} \exp(-\frac{\alpha}{2} (\vec{\varphi},\Omega\vec{\varphi}))}{\int d\varphi P(f|\varphi) \int\limits d\alpha\,P(\alpha) \alpha^{\frac{Rg(\Omega)}{2}} \exp(-\frac{\alpha}{2} (\vec{\varphi},\Omega\vec{\varphi}))} $

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

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

    Дискретизация


    Мы говорили о функциях, но любой реальный прибор не может измерить не то что континуальное, но даже счетное множество точек. Мы всегда проводим измерение в конечном наборе точек, поэтому вынуждены проводить процедуры дискретизации и перехода от интегрального уравнения к матричному. В методе статистической регуляризации мы поступаем следующим образом: будем раскладывать $\varphi(x)$ по некоторой системе функций $\{T_n\}$:

    $\varphi(x) = \sum \limits_n \varphi_n T_n(x).$

    Таким образом, коэффициенты этого разложения образуют некоторый вектор $\vec{\varphi}$, который является вектором в функциональном пространстве.

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

    Тогда элементы матрицы $K$ вычисляются как:

    $K_{mn} = (\hat{K}T_n(x))(y_m),$

    где $y_m$ — точки, в которых производились измерения. Элементы матрицы $\Omega$ будем вычислять по формуле:

    $\Omega_{ij} = \int\limits_a^b \left(\frac{d^pT_i(x)}{dx}\right)\left(\frac{d^pT_j(x)}{dx}\right)dx,$

    где $a$ и $b$ — границы интервала, на котором определена функция $\varphi(x)$.

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

    $D[\varphi(x)] = D[\sum \limits_n \varphi_n T_n(x)] = \sum\limits_{i,j} \varphi_i\varphi_j cov(T_i(x), T_j(x)).$

    При этом нужно учитывать, что в некоторых случаях представление функции при помощи вектора конечной размерности приводит к частичной потере или изменению информации. Фактически мы можем считать алгебраизацию разновидностью регуляризации, однако слабой и недостаточной для превращения некорректной задачи в корректную. Но, так или иначе, мы теперь перешли от поиска $\varphi(x)$ к поиску вектора $\vec{\varphi}$ и в следующем разделе найдем-таки его.

    Случай гауссовых шумов


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

    $\vec{\varphi} = (K^T\Sigma^{-1}K + \alpha^*\Omega)^{-1}K^T\Sigma^{-1^{T}}\vec{f}$

    $\Sigma_{\vec{\varphi}} = (K^T\Sigma^{-1}K+\alpha^*\Omega)^{-1}, $

    где $\Sigma$ — ковариационная матрица многомерного распределения Гаусса, $\alpha^*$ — наиболее вероятное значение параметра $\alpha$, которое определяется из условия максимума апостериорной плотности вероятности:

    $P(\alpha|\vec{f}) = C'\alpha^{\frac{Rg(\Omega)}{2}}\sqrt{|(K^T\Sigma^{-1}K+\alpha\Omega)^{-1}|}\exp(\frac{1}{2} \vec{f}^T\Sigma^{-1}K^{T}(K^T\Sigma^{-1}K+\alpha\Omega)^{-1}K^T\Sigma^{-1^{T}}\vec{f}) $


    А если у меня не гауссовы ошибки?


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

    $\left\langle \varphi_i \right\rangle = \frac{\int d\varphi\, \varphi_i P(f|\varphi) \int\limits d\alpha\,P(\alpha) \alpha^{\frac{Rg(\Omega)}{2}} \exp(-\frac{\alpha}{2} (\vec{\varphi},\Omega\vec{\varphi}))}{\int d\varphi P(f|\varphi) \int\limits d\alpha\,P(\alpha) \alpha^{\frac{Rg(\Omega)}{2}} \exp(-\frac{\alpha}{2} (\vec{\varphi},\Omega\vec{\varphi}))} $


    Основная проблема в том, что этот страшный интеграл, во-первых многомерный, а во-вторых в бесконечных пределах. Причем он сильно многомерный, вектор $\vec{\varphi}$ спокойно может иметь размерность $m = 30-50$, а сеточные методы вычисления интегралов имеют сложность типа $O(n^m)$, поэтому малоприменимы в данном случае. При взятии многомерных интегралов хорошо работает интегрирование методами Монте-Карло.

    Причем поскольку у нас пределы бесконечные, мы должны использовать метод существенной выборки (important sampling), однако тогда нам потребуется подбирать функцию для сэмплирования. Что бы сделать все более автоматизированным следует использовать Markov Chain Monte Carlo (MCMC), который может самостоятельно адаптировать семплирующую функцию под вычисляемый интеграл. Про применение MCMC мы поговорим в следующей статье.

    Практическая часть


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

    Python


    Установка

    Устанавливаем через pip:

    pip install statreg

    или качаем исходный код.

    Примеры

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

    Импортируем нужные научные пакеты.

    import numpy as np
    import matplotlib.pyplot as plt
    from scipy.stats import norm
    from scipy.integrate import quad
    %matplotlib inline

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

    a = 0
    b = 5
    # Восстанавливаемый сигнал
    phi = lambda x: 4*norm.pdf(x-2, scale=0.4) + 2*norm.pdf(x-4, scale = 0.5)
    x = np.linspace(a, b,100)
    plt.plot(x, phi(x));


    Определим ядро и операцию свертки функций (Примечание: np.convolution определенно для массивов):

    kernel = lambda x,y : np.heaviside(x-y, 1) # Ядро уравнения
    convolution =  np.vectorize(lambda y: quad(lambda x: kernel(x,y)*phi(x), a,b)[0])

    Генерируем измеренные данные и зашумляем их с помощью нормального распределения:

    y = np.linspace(a, b, 50)
    ftrue = convolution(y)
    sig = 0.05*ftrue +0.01 # Ошибка измерения
    f = norm.rvs(loc = ftrue, scale=sig)
    plt.errorbar(y, f, yerr=sig);



    Решаем интегральное уравнение

    Импортируем решатель и вспомогательный класс для дискретизации:

    from statreg.model import GaussErrorUnfolder
    from statreg.basis import CubicSplines

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

    basis = CubicSplines(y, boundary='dirichlet')
    model = GaussErrorUnfolder(basis, basis.omega(2))

    Решаем уравнение:

    phi_reconstruct = model.solve(kernel, f, sig, y)

    Строим график решения:

    plt.plot(x,phi(x))
    phir = phi_reconstruct(x)
    phiEr = phi_reconstruct.error(x)
    plt.plot(x, phir, 'g')
    plt.fill_between(x, phir-phiEr, phir + phiEr, color='g', alpha=0.3);



    Решаем матричное уравнение

    Импортируем решатель и вспомогательный класс для дискретизации:

    from statreg.model import GaussErrorMatrixUnfolder
    from statreg.basis import CubicSplines

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

    cubicSplines = CubicSplines(y, boundary='dirichlet')
    omega = cubicSplines.omega(2)
    Kmn = cubicSplines.discretizeKernel(kernel,y)

    Решаем матричное уравнение:

    model = GaussErrorMatrixUnfolder(omega)
    result = model.solve(Kmn, f, sig)

    Строим график:

    phir = lambda x: sum([p*bf(x) for p, bf in zip(result.phi,cubicSplines.basisFun)])
    plt.plot(x,phir(x))
    plt.plot(x,phi(x));



    Julia


    Как мы упоминали, для дальнейшего развития методики требуется продвинутое Монте-Карло интегрирование. Мы могли бы использовать какой-либо модуль на Python (например, мы работали с PyMC3), однако мы, помимо прочего, участвуем в совместном проекте с институтом Макса Планка в Мюнхене.

    Этот проект называется Bayesian Analysis Toolkit. Его цель — создание фреймворка с инструментами для байесовских методов анализа, в первую очередь включая инструменты для MCMC. Сейчас команда работает над второй версией фреймворка, которая пишется на Julia (первая написана на плохом C++). Одна из задач нашей группы — продемонстрировать возможности этого фреймворка на примере статистической регуляризации, поэтому мы написали реализацию на Julia.

    using PyCall
    include("../src/gauss_error.jl")
    include("../src/kernels.jl")
    
    a = 0.
    b = 6.
    
    function phi(x::Float64)
        mu1 = 1.
        mu2 = 4.
        n1 = 4.
        n2 = 2.
        sig1 = 0.3
        sig2 = 0.5
    
        norm(n, mu, sig, x) = n / sqrt(2 * pi*sig^2) * exp(-(x - mu)^2 / (2 * sig^2))
        return norm(n1, mu1, sig1, x) + norm(n2, mu2, sig2, x)
    end
    x = collect(range(a, stop=b, length=300))
    
    import PyPlot.plot
    
    myplot = plot(x, phi.(x))
    savefig("function.png", dpi=1000)


    На этот раз используем другое ядро, будем брать не интегрирующую ступеньку, а свертку с гауссом, которая фактическим наводит некий «блюр» на наши данные:

    function kernel(x::Float64, y::Float64)
        return getOpticsKernels("gaussian")(x, y)
    end
    
    convolution = y -> quadgk(x -> kernel(x,y) * phi(x), a, b, maxevals=10^7)[1]
    y = collect(range(a, stop = b, length=50))
    ftrue = convolution.(y)
    sig = 0.05*abs.(ftrue) +[0.01 for i = 1:Base.length(ftrue)]
    using Compat, Random, Distributions
    noise = []
    for sigma in sig
        n = rand(Normal(0., sigma), 1)[1]
        push!(noise, n)
    end
    f = ftrue + noise
    plot(y, f)


    Аналогично возьмем базис сплайнов с закрепленными концами:

    basis = CubicSplineBasis(y, "dirichlet")
    Kmn = discretize_kernel(basis, kernel, y)
    model = GaussErrorMatrixUnfolder([omega(basis, 2)], "EmpiricalBayes", nothing, [1e-5], [1.], [0.5])
    result = solve(model, Kmn, f, sig)
    phivec = PhiVec(result, basis)
    
    x = collect(range(a, stop=b, length=5000))
    plot(x, phi.(x))
    
    phi_reconstructed = phivec.phi_function.(x)
    phi_reconstructed_errors = phivec.error_function.(x)
    
    plot(x, phi_reconstructed)
    fill_between(x, phi_reconstructed - phi_reconstructed_errors, phi_reconstructed + phi_reconstructed_errors, alpha=0.3)



    Пример на реальных данных

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

    Вот так выглядит исходный интегральный спектр:

    А так — результат восстановления:

    Анализ с помощью фита имеет три основных недостатка:

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

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

    Выводы и анонс следующей части


    Мы рассказали вам о полезной методике, которую можно адаптировать под многие задачи анализа данных (в том числе машинного обучения), причем получить честную «подгонку» ответа — наиболее рациональное решение уравнение в условиях неопределенности, вызванной ошибками измерения. Как приятный бонус мы получаем значения для ошибки решения. Желающие участвовать в развитии или применять метод статистической регуляризации, могут внести свой вклад в виде кода на Python, Julia или на чем-нибудь еще.

    В следующей части мы поговорим о:

    • Использовании MCMC
    • Разложении Холецкого
    • В качестве практического примера рассмотрим применение регуляризации для обработки сигнала с модели орбитального детектора протонов и электронов

    Ссылки



    Автор статьи: Михаил Зелёный, исследователь лаборатории методов ядерно-физических экспериментов в JetBrains Research.

    Комментарии 8

      0
      Какой хороший язык АЛГОЛ-60, до сих пор код исключительно понятен)
      И документацию тогда умели писать не хуже, чем сейчас (а во многом лучше).
      А исходник — это Приложение к чему?
      0
      Это тот Турчин, который изобрел Рефал?
        0
        Да, это он.
          0
          я фанат его теории эволюции
            0
            Да, книжка (Феномен науки) рекомендуется к прочтению всем, кто видит себя в должности архитектора сложных систем. Вот запись семинара на эту тему: youtu.be/fpBGnYMm3PM (качество к сожаление посредственное).
              0
              Я бы еще для архитекторов рекомендовал бы нашего эволюциониста Северцова почитать.

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

      Самое читаемое