diff --git a/main.tex b/main.tex index 37ba9cb..99acdbe 100644 --- a/main.tex +++ b/main.tex @@ -219,8 +219,8 @@ и \(x0\) --- не скаляр, возвращаемое значение равно \verb|(x, converged, zero_der)|, где: \begin{itemize} - \item converged --- ndarray из значений bool. Указывает, какие элементы сошлись успешно. - \item zero\_der --- ndarray из значений bool. Указывает, какие элементы имеют нулевую производную. + \item converged --- \verb|ndarray| из значений bool. Указывает, какие элементы сошлись успешно. + \item zero\_der --- \verb|ndarray| из значений bool. Указывает, какие элементы имеют нулевую производную. \end{itemize} \item \(disp\) --- bool, необязательный @@ -904,7 +904,8 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения и ограничений поставленной перед исследователем задачи. \subsubsection{Реализации метода в библиотеках numpy, scipy} В библиотеке scipy существует множество реализаций данного метода, -которые отличаются используемым видом сплайна. В модуле +которые отличаются используемым видом сплайна а также поддерживаемой +размерностью функции. В модуле \textbf{scipy.interpolate} описаны следующие объекты \cite{links:scipy_doc}: \begin{enumerate} @@ -912,8 +913,9 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения \item класс \textbf{PchipInterpolator}; \item класс \textbf{CubicHermiteSpline}; \item класс \textbf{Akima1DInterpolator}; - \item класс \textbf{RectBivariateSpline}, используется для многомерной - \item интерполяции; + \item классы \textbf{RectBivariateSpline}, + \textbf{LSQBivariateSpline}, используются для многомерной + интерполяции; \item функция \textbf{interp1d}, которая может интерполировать таблично заданную функцию разными методами, в т.ч. сплайн-интерполяцией. Признана устаревшей и поэтому не будет рассмотрена; @@ -964,13 +966,11 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения Одномерный массив, содержащий значения независимой переменной. Значения должны быть действительными, конечными числами и находиться в строго возрастающем порядке. - \item \(y\) --- array\_like Массив, содержащий значения зависимой переменной. Он может иметь произвольное количество измерений, но длина должна совпадать с длиной \(x\). Значения должны быть конечными. - \item \(axis\) int, необязательный Ось, вдоль которой предполагается изменение y. Это означает, @@ -989,19 +989,16 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения сегменты на конце кривой представляют собой один и тот же полином. Это хороший вариант по умолчанию, когда нет информации о граничных условиях. - \item \verb|"periodic"|: Предполагается, что интерполируемые функции являются периодическими с периодом \(x[-1] \ \mbox {---}\ x[0]\). Первое и последнее значение y должны быть идентичными: \(y[0] == y[-1]\). Это граничное условие приведет к тому, что \(y'[0] == y'[-1]\) и \(y''[0] == y''[-1]\). - \item \verb|"clamped"|: Первая производная на концах кривых равна нулю. Предполагая, что y одномерный, \(bc\_type=((1, 0.0), (1, 0.0))\) --- то же самое условие. - \item \verb|"natural"|: Вторая производная на концах кривой равна нулю. Предполагая, что y одномерный, \(bc\_type=((2, 0.0), (2, 0.0))\) --- то же самое @@ -1017,7 +1014,6 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения \begin{itemize} \item \(order\): порядок производной --- 1 или 2. - \item \(deriv\_value\): array\_like Содержит значения производной, форма должна быть такой @@ -1028,7 +1024,6 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения \(axis=1\), то \(deriv\_values\) должно быть двухмерным и иметь форму \((n0, n2)\). \end{itemize} - \item \(extrapolate\) --- {bool, \verb|"periodic"|, \verb|None|}, необязательный @@ -1041,29 +1036,385 @@ LU-разложение \cite[с. 259]{book:levitin}. Для получения значение \verb|True| в противном случае. \end{enumerate} -Конструктор экземпляра класса \textbf{PchipInterpolator} имеет следующие -параметры: +Класс \textbf{PchipInterpolator} представляет из себя кубический +эрмитов сплайн, поэтому является дочерним классом +\textbf{CubicHermiteSpline}. Конструктор экземпляра класса имеет +следующие параметры: \begin{enumerate} - \item + \item \(x\) --- ndarray, вида (xlen,1), то есть одномерный + + Одномерный массив монотонно возрастающих действительных + значений. Не может включать повторяющиеся значения (иначе + сплайн будет слишком конкретизированный) + \item \(y\) --- ndarray, вида (...,xlen,...) + + N-мерный массив действительных значений. Длина \(y\) по оси + интерполяции должна быть равна длине \(x\). Используйте + параметр \(axis\), чтобы выбрать ось интерполяции. + \item \(axis\) int, необязательный + + См. описание одноименного метода родительского класса. + \item \(extrapolate\) --- bool, необязательный + + Определяет, следует ли экстраполировать точки, + которые выходят за пределы границ, установленные первым и + последним интервалами или возвращать \verb|NaN|. + По умолчанию \verb|None|, значение см. в пункте + \ref{list:CubicHermiteSpline-extrapolate} описания + конструктора родительского класса. \end{enumerate} Конструктор экземпляра класса \textbf{CubicHermiteSpline} имеет следующие параметры: \begin{enumerate} - \item + \item \(x\) --- array\_like, вида (n,), то есть одномерный + + Одномерный массив, содержащий значения независимой переменной. Значения должны быть действительными, конечными и находиться в строго возрастающем порядке. + \item \(y\) --- array\_like + + Массив, содержащий значения зависимой переменной. + Он может иметь произвольное количество измерений, но длина по + оси \(axis\) должна совпадать с длиной \(x\). + Значения должны быть конечными. + \item \(dydx\) --- array\_like + + Массив, содержащий производные зависимой переменной. Он может + иметь произвольное количество измерений, но длина по оси + \(axis\) должна совпадать с длиной \(x\). Значения должны быть + конечными. + \item \(axis\) --- int, необязательный + + Ось, вдоль которой предполагается изменение y. Это означает, что для \(x[i]\) соответствующие значения равны + \verb|np.take(y, i, axis=axis)|. + По умолчанию --- 0. + \item \(extrapolate\) --- {bool, \verb|"periodic"|, + \verb|None|}, необязательный + + Если bool, определяет, следует ли экстраполировать точки, + которые выходят за пределы границ, установленные первым и + последним интервалами или возвращать \verb|NaN|. Если + \verb|"periodic"|, используется периодическая экстраполяция. + Если \verb|None| (по умолчанию), используйте метод \textbf{extrapolate} + этого класса. + \label{list:CubicHermiteSpline-extrapolate} \end{enumerate} Конструктор экземпляра класса \textbf{Akima1DInterpolator} имеет следующие параметры: \begin{enumerate} - \item + \item \(x\) --- ndarray, вида (n, 1) + + Одномерный массив монотонно возрастающих действительных + значений. + \item \(y\) --- ndarray, вида (..., n, ...) + + N-мерный массив реальных значений. Длина \(y\) по оси + интерполяции должна быть равна длине \(x\). Используйте + параметр \(axis\), чтобы выбрать ось интерполяции. + \item \(axis\) --- int, необязательный + + Ось в массиве \(y\), соответствующая значениям координаты + \(x\). По умолчанию --- 0. \end{enumerate} - \subsection{Сглаживание. Метод наименьших квадратов} +С помощью метода наименьших квадратов (далее, МНК) осуществляется +интерполяция функции с помощью сглаживания. \subsubsection{Описание метода} -\subsubsection{Реализации метода в библиотеках numpy, scipy} +Задача метода --- минимизировать отклонение от исходных данных: +\begin{equation} + S = \sum_{i=1}^{n}(\widetilde{y}(x_i,a)-y_i)^2 \rightarrow \min + \label{formula:mnk} +\end{equation} +где \(a\) --- неизвестные параметры; \(n\) --- количество известных +пар значений (\(x,y\)); \(x_i,y_i\) --- значения \(i\) пары. +В качестве \(\widetilde{y}(x)\) необходимо выбрать функцию вида, +который лучше всего подходит под исходные данные. Этот процесс +называют выбором зависимости. + +Примеры зависимостей (\(a,b,c,d\) --- неизвестные параметры): +\begin{alignat*}{2} + & ax + b & \ & \text{--- линейная}; \\ + & ax^2 + bx + c & \ & \text{--- полиномиальная}; \\ + & a \ln(x) + b & \ & \text{--- логарифмическая}; \\ + & be^{ax} & \ & \text{--- экспоненциальная}; +\end{alignat*} + +В результате подстановки выбранной зависимости и известных данных +в (\ref{formula:mnk}) задача минимизации сводится к нахождению +экстремумов функции \(S(a)\), и выбору подходящего из них. +\subsubsection{Реализации метода в библиотеках numpy, scipy} +В модуле \textbf{scipy.optimize} библиотеки scipy существует +несколько реализаций данного метода, каждая из которых имеет свою +специализацию. +\begin{enumerate} + \item \textbf{curve\_fit} --- функция, предназначенная для поиска + значения параметра \(a\) произвольной зависимости \(\widetilde + {y}(x_i,a)\), при котором выполняется условие + (\ref{formula:mnk}). + \item \textbf{least\_squares} --- функция, которая отличается от + \textbf{curve\_fit} возможностью задать собственный критерий + близости, область допустимых решений; позволяет задать больше + параметров точности. + \item \textbf{leastsq} --- функция, которая является устаревшей + версией \textbf{least\_squares}, поэтому ее подробного + рассмотрения не будет. + \item \textbf{lsq\_linear} --- функция, предназначенная для + минимизации выражения \(||ax-b||^2\), где \(a,b\) --- матрицы + параметров; функция \textbf{nnls} отличается от нее + ограничением на положительность искомых значений. +\end{enumerate} + +Функция \textbf{curve\_fit} имеет следующий набор параметров: +\begin{enumerate} + \item \(f\) --- callable + + Модельная функция \(f(x,\dots)\). Она должен принимать независимую переменную в качестве первого аргумента, а параметры --- в качестве оставшихся. + \item \(xdata\) --- array\_like + + Значения независимой переменной. Обычно это должна быть последовательность длины \(M\) или массив \((k,M)\)-образной формы для функций с \(k\) предикторами, и каждый элемент должен быть конвертируемым в \(float\), если это объект, подобный массиву. + \item \(ydata\) --- array\_like + + Значения зависимой переменной, массив длиной M. + При найденных параметрах \(f\) совпадает с значением + \(f(xdata, \dots)\). + \item \(p0\) --- array\_like, необязательный + + Начальное значение параметров \(f\) (длины N). Если + \verb|None| (значение по умолчанию), то все начальные + значения будут равны \verb|1| (если количество параметров + функции можно определить с помощью интроспекции, в противном + случае генерируется исключение \verb|ValueError|). + \item \(sigma\) --- \verb|None| / sequence длины M / массив MxM, необязательный + + Задает неопределенность в \(ydata\). Если мы определим + остатки как \(r = ydata - f(xdata, *popt)\), то интерпретация + значения \(sigma\) будет зависеть от ее количества измерений: + \begin{itemize} + \item одномерная \(sigma\) должна содержать значения + стандартных отклонений ошибок в \(ydata\). Используется + при вычислении параметра + \verb|chisq = sum((r/sigma)**2)|; + \item двухмерная \(sigma\) должна содержать ковариационную + матрицу ошибок в \(ydata\). Используется + при вычислении параметра + \verb|chisq = r.T @ inv(sigma) @ r|. + \end{itemize} + \item \(absolute\_sigma\) --- bool, необязательный + + Если значение \verb|True|, \(sigma\) используется в абсолютном смысле, и предполагаемая ковариация параметра \(pcov\) отражает эти абсолютные значения. + + Если значение \verb|False| (по умолчанию), имеют значение + только относительные величины значений \(sigma\). + Ковариационная матрица возвращаемого параметра \(pcov\) + основана на масштабировании \(sigma\) постоянным + коэффициентом. Эта константа устанавливается требованием, чтобы + приведенный выше \(chisq\) для параметров + \(popt\) при использовании масштабированной \(sigma\) был равен + единице. Другими словами, \(sigma\) масштабируется так, чтобы + соответствовать выборочной дисперсии остатков после подгонки. + Математически, + \verb|pcov(absolute_sigma=False) = pcov(absolute_sigma=True)| + \verb| * chisq(popt)/(M-N)|. + \item \(check\_finite\) --- bool, необязательный + + Если \verb|True|, то входные массивы проверяются на содержание + значений \verb|NaN| или \verb|Inf|, и если они содержат, то + генерируется исключение \verb|ValueError|. + + Значение \verb|False| может привести к выдаче бессмысленных + результатов, если входные массивы содержат \verb|NaN|. + + По умолчанию установлено значение \verb|True|, если + \(nan\_policy\) не указано явно, и значение \verb|False| в + противном случае. + \item \(bounds\) --- tuple из двух array\_like / \verb|Bounds|, + необязательный + + Нижняя и верхняя границы параметров \(f\). По умолчанию нет ограничений. Есть два способа указать границы: + \begin{itemize} + \item через экземпляр класса \verb|Bounds|; + + \item \verb|tuple| из двух array\_like, где каждый элемент + кортежа должен быть либо массивом с длиной, равной + количеству параметров, либо скаляром (в этом случае + граница считается одинаковой для всех параметров). + Используйте \verb|np.inf| с соответствующим знаком, + чтобы отключить ограничения для всех или некоторых + параметров. + \end{itemize} + \item \(method\) --- \{\verb|"lm"|, \verb|"trf"|, \verb|"dogbox"|\}, необязательный + + Метод, используемый для вычисления параметров. См. описание + аналогичного параметра функции \textbf{least\_squares} для + получения более подробной информации. По умолчанию используется + \verb|"lm"| для задач без ограничений и \verb|"trf"|, если + указано значение \(bounds\). + \item \(jac\) --- callable / string / \verb|None|, необязательный + + Функция с сигнатурой \(jac(x, ...)\), которая вычисляет матрицу + Якоби модельной функции относительно параметров как плотную + array\_like структуру. Он будет масштабироваться в + соответствии с предоставленным параметром \(sigma\). Если + \verb|False| (по умолчанию), якобиан будет оцениваться + численно. Строковые ключевые слова для методов \verb|"trf"| и + \verb|"dogbox"| можно использовать для выбора + конечно-разностной схемы, см. описание + \textbf{least\_squares}. + \item \(full\_output\) --- boolean, необязательный + + Если значение = \verb|True|, эта функция возвращает дополнительную информацию: \(infodict\), \(mesg\) и \(ier\). По умолчанию + \verb|False|. + \item \(nan\_policy\) --- {\verb|"raise"|, \verb|"omit"|, \verb|None|}, необязательный + + Определяет, как действовать, если входные данные содержат nan. Доступны следующие параметры (по умолчанию — \verb|None|): + \begin{itemize} + \item \verb|"raise"|: выдает ошибку + + \item \verb|"omit"|: выполняет вычисления, игнорируя + значения \verb|NaN|. + + \item \verb|None|: никакая специальная обработка \verb|NaN| + не выполняется (кроме того, что делается с помощью + \(check\_finite\)); поведение при наличии + \verb|NaN| зависит от реализации и может измениться. + + \end{itemize} + Обратите внимание: если значение указано явно (не \verb|None|), + для \(check\_finite\) будет установлено значение \verb|False|. + \item \(**kwargs\) + + Именованные параметры, которые передаются в \textbf{leastsq} + для \(method\)='lm' или \textbf{least\_squares} в противном + случае. +\end{enumerate} +Данная функция возвращает результат в виде +\((popt,pcov,infodict,mesg,ier)\), где +\begin{enumerate} + \item \(popt\) --- вычисленные параметры + \(f\), при которых выполняется условие + \(f(xdata, *popt) - ydata \rightarrow \min\). + \item \(pcov\) --- 2-D array + + Вычисленная приблизительная ковариация \(popt\). + \item \(infodict\) --- dict + + По ключу \(nfev\) + хранится значение количества вызовов функции. Методы \("trf"\) + и \("dogbox"\) не учитывают вызовы функций для численной + аппроксимации якобиана, в отличие от метода \("lm"\). + По ключам \(fvec\), \(fjac\), \(ipvt\), \(qtf\) можно получить + дополнительную информацию о работе алгоритма. + \item \(mesg\) --- str + + Cтроковое сообщение с информацией о решении. + \item \(ier\) --- int + + Целочисленный флаг. Если он равен 1, 2, 3 или 4, решение + найдено. В противном случае решение не было найдено. \(mesg\) + предоставляет дополнительную информацию о решении. +\end{enumerate} +Разработчики отмечают, что входные данные \(xdata\), \(ydata\) и +выходные данные \(f\) должны иметь формат \(float64\), иначе данная +функция может вернуть неверные результаты \cite{links:scipy_doc}. + +Функция \textbf{least\_squares} предназначена для решения следующей +задачи минимизации \(F(x)\), где +\verb|F(x) = 0.5 * sum(rho(f_i(x)**2), i = 0, ..., m - 1)| +при ограничениях \(lb \leq x \leq ub\). То есть, \(f\) должна принимать +на вход значение параметров \(a\) и возвращать значение \(f(x,a) - y\). + +Данная функция имеет следующий набор параметров: +\begin{enumerate} + \item \(fun\) --- callable + Функция, которая вычисляет вектор остатков, имеет сигнатуру + \(fun(x, *args, **kwargs)\), т.е. минимизация продолжается + относительно ее первого аргумента. Аргумент \(x\), + передаваемый этой функции, представляет собой N-мерный массив размерности \((N,K)\) (никогда не скаляр, даже для \(N=1\)). \(fun\) должна вернуть одномерный массив типа array \((m,)\) или скаляр. Если \(x \in \textbf{C}\) или функция \(fun\) возвращает комплексные числа, ее необходимо обернуть действительнозначной функцией от действительнозначных аргументов. + \item \(x_0\) --- array\_like размерности (N,K) / float + + Начальное значение \(x\). Если параметр имеет тип float, он + будет рассматриваться как одномерный массив с одним элементом. + \item \(jac\) --- \{\verb|"2-point"|, \verb|"3-point"|, \verb|"cs"|, callable\}, необязательный + + Метод вычисления матрицы Якоби (матрица размером \(m \times n\), где элемент \((i, j)\) является частной производной f[i] по x[j]). Ключевые слова выбирают схему конечных разностей для числовой оценки. Схема \verb|"3-point"| более точная, но требует в два раза больше операций, чем \verb|"2-point"| (по умолчанию). Схема \verb|"cs"| использует комплексные шаги и, хотя потенциально является наиболее точной, она применима только тогда, когда \(fun\) правильно обрабатывает комплекснозначные входные данные и может быть аналитически продолжена на комплексную плоскость. \(method\) \verb|"lm"| всегда использует схему \verb|"2-point"|. Если параметр имеет тип \verb|callable|, он используется как \verb|jac(x, *args, **kwargs)| и должен возвращать хорошее приближение (или точное значение) для якобиана в виде array\_like (к выводу \(jac\) будет применен \verb|np.atleast_2d|), разреженной матрицы (предпочтительна \verb|csr_matrix| для производительности) или \verb|scipy.sparse.linalg.LinearOperator|. + \item \(bounds\) --- tuple из 2 array\_like / \verb|Bounds|, необязательный + + Есть два способа указать границы: + \begin{itemize} + \item через экземпляр класса \verb|Bounds|; + + \item \verb|tuple| из двух array\_like, где каждый элемент + кортежа должен быть либо массивом с размерностью + \(x0\), либо скаляром (в этом случае + граница считается одинаковой для всех параметров). + Используйте \verb|np.inf| с соответствующим знаком, + чтобы отключить ограничения для всех или некоторых + параметров. + \end{itemize} + \item \(method\) --- \{\verb|"lm"|, \verb|"trf"|, \verb|"dogbox"|\}, необязательный + + Допустимое отклонение на завершение по норме градиента. + По умолчанию --- 1e-8. Точное состояние зависит от используемого + значения \(method\): + \begin{itemize} + \item Для \verb|"trf"|: \(norm(g\_scaled, ord=np.inf) < gtol\), где \(g\_scaled\) --- значение градиента, масштабированное для учета наличия границ. + \item Для \verb|"dogbox"|: \(norm(g_free, ord=np.inf) < gtol\), где \(g_free\) --- градиент по отношению к переменным, которые не находятся в оптимальном состоянии на границе. + \item Для \verb|"lm"|: максимальное абсолютное значение косинуса углов между столбцами якобиана и вектором невязки меньше \(gtol\), или вектор невязки равен нулю. + \end{itemize} + + Если \verb|None| и \(method\) не \verb|"lm"|, завершение по этому условию отключено. Если \(method\) --- \verb|"lm"|, этот допуск должен быть выше, чем машинный эпсилон. + \item \(x\_scale\) --- array\_like / \verb|"jac"| + + Характеристический масштаб каждой переменной. + Установка \(x\_scale\) эквивалентна переформулировке проблемы + в масштабированных переменных \(xs = x/x\_scale\). + Альтернативная точка зрения состоит в том, что размер + доверительной области по j-му измерению пропорционален + \(x\_scale[j]\). Улучшения сходимости можно достичь, установив + \(x\_scale\) таким образом, чтобы шаг заданного размера по + любой из масштабируемых переменных имел аналогичный эффект на + функцию стоимости. Если установлено значение \verb|"jac"|, + масштаб итеративно обновляется с использованием обратных норм + столбцов матрицы Якоби. + \item \(loss\) --- str / callable, необязательный + + Определяет функцию потерь. Допускаются следующие значения ключевых слов: + \begin{itemize} + \item \verb|"linear"| (по умолчанию) : \verb|rho(z) = z|. Дает стандартную задачу наименьших квадратов. + \item \verb|"soft_l1"| : \verb|rho(z) = 2 * ((1 + z)**0.5 - 1)|. Гладкая аппроксимация потерь l1 (абсолютное значение). Обычно хороший выбор для надежного метода наименьших квадратов. + \item \verb|"huber"| : \verb|rho(z) = z if z <= 1 else 2*z**0.5 - 1|. Работает аналогично \verb|"soft_l1"|. + \item \verb|"cauchy"| : \verb|rho(z) = ln(1 + z)|. Сильно ослабляет влияние выбросов, но может вызвать трудности в процессе оптимизации. + \item \verb|"arctan"| : \verb|rho(z) = arctan(z)|. Ограничивает максимальный убыток по одному остатку, имеет свойства, аналогичные \verb|"cauchy"|. + \end{itemize} + + Если параметр имеет тип \verb|callable|, он должен принимать + одномерный \verb|ndarray| \verb|z=f**2| и возвращать + array\_like с формой \((3, m)\), где строка 0 содержит значения функции, строка 1 содержит первые производные, а строка 2 содержит вторые производные. \(method\) \verb|lm| поддерживает только \verb|"linear"| \(loss\). + + \item \(f\_scale\) --- float, необязательный + + Значение границы между ошибками данных и выбросами в данных, + по умолчанию - 1,0. Значение функции потерь оценивается + следующим образом: \verb|rho_(f**2) = C**2 * rho(f**2 / C**2)|, + где C --- это \(f\_scale\), а \(rho\) определяется параметром + \(loss\). Этот параметр не имеет никакого эффекта при + \(loss\)=\verb|"linear"|, но для других значений потерь он имеет решающее + значение. + \item \(max\_nfev\) --- \verb|None| / int, необязательный + + Максимальное количество выполнений функции перед завершением работы. Если \verb|None| (по умолчанию), значение выбирается автоматически: + \begin{itemize} + \item Для \verb|"trf"| и \verb|"dogbox"|: \(100 \cdot n\). + + \item Для \verb|"lm"|: \(100 \cdot n\), если \(jac\) является вызываемым, и \(100 \cdot n \cdot (n + 1)\) в противном случае (поскольку \verb|"lm"| учитывает вызовы функций в оценке Якобиана). + \end{itemize} + \item \(diff\_step\) --- \verb|None| / array\_like, необязательный + + + \item +\end{enumerate} \section{Численное интегрирование} \subsection{Метод прямоугольников} \subsubsection{Описание метода}