SciPy - это мощнейшая библиотека для Python, позволяющая с легкостью решать сложные научные и инженерные задачи. Хотите узнать, как с ее помощью расширить возможности Python для анализа данных, моделирования и визуализации? Тогда эта статья для вас!
Что такое SciPy и зачем она нужна
SciPy - это библиотека Python с открытым исходным кодом, предназначенная для выполнения научных и инженерных расчетов . Она построена на базе NumPy и позволяет управлять данными, а также визуализировать их с помощью разных высокоуровневых команд.
Целевая аудитория SciPy - пользователи продуктов MATLAB и Scilab. Для визуализации результатов расчетов часто применяется библиотека Matplotlib, являющаяся аналогом средств вывода графики MATLAB.
Преимущества SciPy:
- Открытый исходный код
- Интеграция с экосистемой Python
- Высокая производительность за счет использования С, C++ и Fortran в критичных к скорости вычислениях
Основной структурой данных в SciPy является многомерный массив, реализованный модулем NumPy (более старые версии SciPy использовали модуль Numeric).
Установка и импорт SciPy
Установить SciPy можно с помощью менеджеров пакетов pip или conda:
pip install scipy conda install scipyТакже есть "научные дистрибутивы" Python, куда уже входит SciPy вместе с другими полезными для анализа данных библиотеками. Примеры таких дистрибутивов:
- Anaconda
- Canopy
Чтобы импортировать SciPy, достаточно написать:
import scipyНо чаще импортируют отдельные пакеты в зависимости от решаемой задачи:
from scipy import linalg # для линейной алгебры from scipy.optimize import minimize # для оптимизации import scipy.stats as st # статистика 
Исследуем возможности SciPy на примерах
Возможности библиотеки SciPy распределены по нескольким модулям, или пакетам, которые объединены общим назначением.
Базовые функции
SciPy предоставляет функции help(), info() и source() для получения информации о других функциях и методах:
help(scipy.integrate.quad) # информация о функции quad info(scipy.signal) # инфо о модуле обработки сигналов source(scipy.spatial.Delaunay) # исходный код функции Делоне 
Специальные математические функции
В SciPy есть набор специальных функций, используемых в математической физике: эллиптические, гамма, бета и так далее. Например, для вычисления экспоненты используется scipy.special.exp10:
import scipy.special as sc sc.exp10(3) # 1000.0 Другой пример - вычисление синуса в радианах:
x = 2.0 sc.sin(x) # 0.909297426826 Работа с интегралами и дифференциальными уравнениями
Для вычисления определенного интеграла в SciPy используется функция quad():
from scipy import integrate integrate.quad(lambda x: x**2, 0, 4) # (21.333333333333332, 2.3684757858670003e-13) Результат - кортеж со значением интеграла и оценкой погрешности.
А для решения обыкновенных дифференциальных уравнений (ОДУ) используется scipy.integrate.solve_ivp. Например, решим уравнение Лоренца:
from scipy.integrate import solve_ivp import numpy as np def lorenz(t, z): # Уравнения Лоренца как система ОДУ x, y, z = z return [-8/3*x + y*z, -10*y + 28*x - y*z, -z + x*y] # Начальные условия z0 = [8, 27, -7] sol = solve_ivp(lorenz, [0, 4], z0) print(sol.y) # решение - массив значений x, y, z Оптимизация и поиск экстремумов
Для оптимизации в SciPy используется модуль scipy.optimize. Он содержит различные алгоритмы:
- метод сопряженных градиентов
- метод Нелдера-Мида
- дифференциальная эволюция
Например, найдем минимум функции Розенброка:
from scipy.optimize import minimize def rosenbrock(x): return (1 - x[0])**2 + 105*(x[1] - x[0]**2)**2 x0 = [-3, -1] res = minimize(rosenbrock, x0, method='Nelder-Mead') print(res.x) # [1, 1] - координаты минимума Интерполяция данных
Интерполяция в SciPy может быть как одномерной, так и многомерной. Для одномерной интерполяции используется функция interp1d:
import numpy as np from scipy.interpolate import interp1d x = np.linspace(0, 10, 6) y = np.sin(x) interp_func = interp1d(x, y) interp_func(3.5) # интерполированное значение для x=3.5 А для многомерной интерполяции применяется, например, interp2d:
from scipy.interpolate import interp2d import numpy as np x = [1, 2, 3, 4] y = [10, 20, 25, 50] z = [[1, 2, 1, 2], [3, 4, 3, 4], [9, 4, 3, 8], [3, 7, 8, 5]] interp_func = interp2d(x, y, z) interp_func(2.5, 20) # интерполированное значение для x=2.5, y=20 Таким образом, SciPy предоставляет гибкие средства для интерполяции данных произвольной размерности.
Быстрое преобразование Фурье
Анализ Фурье часто используется для обработки сигналов и изображений. В SciPy для этого есть функции fft и ifft:
from scipy import fft import numpy as np signal = np.array([1, 2, 1, -1, 4, 10]) fft_signal = fft(signal) ifft_signal = ifft(fft_signal) print(ifft_signal) # восстановленный исходный массив Цифровая фильтрация и обработка сигналов
Пакет scipy.signal предоставляет различные возможности для обработки и анализа сигналов.
Например, функция correlate вычисляет кросс-корреляцию двух сигналов:
from scipy import signal import numpy as np sig1 = np.array([1, 2, 3, 4, 5]) sig2 = np.array([3, 4, 5, 6]) signal.correlate(sig1, sig2) # array([ 9, 14, 19, 24, 14]) А convolve - свертку сигналов:
signal.convolve(sig1, sig2[::-1]) # array([ 3, 10, 22, 28, 22, 8]) Для фильтрации шумов из сигнала используются различные цифровые фильтры: bessel, butter, cheby1 и другие:
from scipy import signal import numpy as np import matplotlib.pyplot as plt t = np.linspace(0, 1, 1000) sig = np.sin(2*np.pi*10*t) + 0.2*np.random.rand(len(t)) f_cut = 15 sos = signal.cheby1(10, 5, f_cut, 'hp', output='sos', fs=1000) filtered = signal.sosfilt(sos, sig) plt.plot(t[:30], sig[:30]) plt.plot(t[:30], filtered[:30]) plt.legend(['Исходный сигнал', 'Отфильтрованный']) Здесь Chebyshev фильтром высоких частот убран шум с частотой ниже 15 Гц.
SciPy для математического моделирования
Благодаря наличию различных специальных функций, дифференциальных уравнений и оптимизаторов, с помощью SciPy можно создавать математические модели для имитации физических, химических и других процессов.
Использование констант и специальных функций
Например, константы из подмодуля scipy.constants позволяют задавать параметры среды:
from scipy import constants liquid_n2_temp = constants.nitrogen_boiling_temperature # 77.355 K room_temp = constants.zero_Celsius + 25.0 # 298.15 K А специальные математические функции могут описывать физические законы. К примеру, закон Аррениуса для скорости химической реакции:
import scipy.special as sc A = 2.5e9 # предэкспоненциальный множитель Ea = 40e3 # энергия активации R = 8.31 # газовая постоянная T = 298.15 # температура k = A * np.exp(-Ea / (R*T)) # константа скорости реакции Решение дифференциальных уравнений
Рассмотрим модель химического реактора как систему обыкновенных дифференциальных уравнений (ОДУ). SciPy позволяет численно решить такую систему:
from scipy.integrate import odeint import numpy as np # скорости реакций k1 = 0.13 k2 = 0.19 # система ОДУ def reactor_model(C, t): Ca, Cb = C dCa_dt = -k1*Ca - k2*Ca*Cb dCb_dt = k1*Ca - k2*Ca*Cb return [dCa_dt, dCb_dt] # начальные условия Ca0, Cb0 = 1, 0 t = np.linspace(0, 10) C0 = [Ca0, Cb0] C = odeint(reactor_model, C0, t) Здесь концентрации реагентов Ca и Cb получены решением ОДУ для разных моментов можно визуализировать полученное решение, например с помощью Matplotlib:
import matplotlib.pyplot as plt plt.plot(t, C[:,0], label='C_a') plt.plot(t, C[:,1], label='C_b') plt.legend() plt.xlabel('Время') plt.ylabel('Концентрация') plt.show() Таким образом, используя SciPy, можно решать дифференциальные уравнения, описывающие кинетику химических или физических процессов.
Оптимизация параметров моделей
Часто требуется подобрать параметры модели, при которых ее поведение наилучшим образом соответствует экспериментальным данным. Это можно сделать с помощью оптимизационных алгоритмов из SciPy.
Например, для нашей модели реактора можно варьировать константы скоростей реакций k1 и k2 и минимизировать невязку между расчетными и измеренными концентрациями:
from scipy.optimize import minimize def optimize(params): k1, k2 = params # решение модели C_model = odeint(reactor_model, C0, t, args=(k1, k2)) # сравнение с экспериментальными данными error = np.sqrt(np.sum((C_exp - C_model)**2)) return error res = minimize(optimize, x0=[0.1, 0.1]) opt_k1, opt_k2 = res.x # оптимальные k1 и k2 Благодаря этому модель становится более точной.
Верификация моделей
Полученную модель всегда нужно верифицировать, сравнивая ее поведение с реальными данными. SciPy может помочь вычислить различные метрики:
- Цреднеквадратичная ошибка
- Коэффициент детерминации \((R^2)\)
- Среднее абсолютное отклонение
Например:
from scipy import stats MSE = np.square(np.subtract(C_model, C_exp)).mean() R2 = stats.pearsonr(C_model[:,0], C_exp[:,0])[0] ** 2 Эти метрики позволяют численно оценить, насколько хорошо модель соответствует данным.
Визуализация данных с SciPy и Matplotlib
Хотя в SciPy есть базовые средства визуализации, обычно для построения графиков используют специализированные библиотеки вроде Matplotlib, Plotly или Bokeh.
Интеграция SciPy и Matplotlib
Matplotlib это фактически стандарт де-факто для визуализации данных с Python. SciPy и Matplotlib прекрасно работают вместе.
Например, результаты моделирования можно визуализировать следующим образом:
import matplotlib.pyplot as plt plt.plot(t, C[:,0], 'r') plt.plot(t, C_exp[:,0], 'bo') plt.legend(['Модель', 'Эксперимент']) plt.xlabel('t, c') plt.ylabel('Концентрация, моль/л') plt.show() 2D и 3D графики
С помощью Matplotlib можно строить различные 2D графики: линейные, столбчатые, круговые и т.д. А для 3D визуализации данных используется модуль mpl_toolkits.mplot3d.
Например, посторим поверхность по рассчитанным значениям Ца и Cb:
from mpl_toolkits import mplot3d import numpy as np import matplotlib.pyplot as plt fig = plt.figure() ax = plt.axes(projection='3d') ax.plot_surface(t, t, C[:,0].reshape(10,10), rstride=1, cstride=1) ax.set_xlabel('t') ax.set_ylabel('t') ax.set_zlabel('Ca') plt.show() Интерактивные графики
Для создания интерактивных дашбордов используется библиотека Plotly, которая умеет работать с данными из SciPy.
Пример простого интерактивного графика концентраций:
import plotly.graph_objects as go fig = go.Figure(data=[ go.Scatter(y=C[:,0], name='Ca'), go.Scatter(y=C[:,1], name='Cb') ]) fig.update_layout( title='Концентрации', xaxis_title='t, c', yaxis_title='C, моль/л' ) fig.show() Расширение SciPy своими алгоритмами
Хотя SciPy содержит очень много готовых алгоритмов, иногда возникает необходимость реализовать свои собственные, отвечающие специфике решаемой задачи.
Создание пользовательских функций
Благодаря открытой архитектуре SciPy, можно создавать собственные функции и модули на основе ее API.
Например, реализуем алгоритм Винера для фильтрации зашумленных сигналов:
from scipy.signal import convolve import numpy as np def wiener(signal, noise): signal_fft = fft(signal) noise_fft = fft(noise) H = signal_fft / (np.abs(noise_fft)**2 + np.abs(signal_fft)**2) filtered = ifft(H * signal_fft) return np.real(filtered) filtered_sig = wiener(signal, noise) Использование Numba для оптимизации
Чтобы ускорить работу алгоритмов на Python можно задействовать JIT-компилятор Numba.
Декоратор @numba.njit позволяет это сделать очень просто:
from numba import njit @njit def wiener_fast(signal, noise): # тот же код, что и выше filtered_sig = wiener_fast(signal, noise) За счет этого выполнение функции ускорится в разы.
Интеграция с моделями машинного обучения
SciPy часто используют совместно с библиотеками машинного обучения, такими как sklearn, PyTorch, TensorFlow.
