Сгладить кривую в питоне без ошибок на границах?

Самый простой способ сделать это — детрендировать сигнал перед фильтрацией. Краевые эффекты, которые вы видите, в основном связаны с тем, что сигнал не является стационарным (т.е. он имеет наклон).

Прежде всего, давайте продемонстрируем проблему:

import numpy as np
import matplotlib.pyplot as plt
from scipy.ndimage import gaussian_filter1d

x, y = np.loadtxt('spectrum.dat').T

# Smooth with a guassian filter
smooth = gaussian_filter1d(y, 10)

fig, ax = plt.subplots()
ax.loglog(x, y, color='black')
ax.loglog(x, smooth, color='red')
plt.show()

Ой! Краевые эффекты особенно плохи в конце (правый размер) ваших данных, потому что там наклон самый крутой. Если бы у вас был более крутой склон в начале, вы бы также увидели более сильный краевой эффект.

Хорошей новостью является то, что есть несколько способов исправить это. Ответ @ChristianK. показывает, как использовать сглаживающие сплайны для эффективного предварительного формирования фильтра нижних частот. Я продемонстрирую использование нескольких других методов обработки сигналов для достижения той же цели. Что является «лучшим», все зависит от ваших потребностей. Сглаживание сплайнов выполняется прямолинейно. Использование более «причудливых» методов обработки сигналов дает вам контроль над тем, какая именно частота отфильтровывается.

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

В качестве быстрого примера:

import numpy as np
import matplotlib.pyplot as plt
from scipy.ndimage import gaussian_filter1d

x, y = np.loadtxt('spectrum.dat').T

# Let's detrend by fitting a second-order polynomial in log space
# (Note that your data looks like a parabola in log-log space.)
logx, logy = np.log(x), np.log(y)
model = np.polyfit(logx, logy, 2)
trend = np.polyval(model, logx)

# Smooth with a guassian filter
smooth = gaussian_filter1d(logy - trend, 10)

# Add the trend back in and convert back to linear space
smooth = np.exp(smooth + trend)

fig, ax = plt.subplots()
ax.loglog(x, y, color='black')
ax.loglog(x, smooth, color='red')
plt.show()

Обратите внимание, что у нас все еще есть некоторые краевые эффекты. Это связано с тем, что фильтр Гаусса, который я использовал, вызывает фазовый сдвиг. Если бы мы действительно хотели пофантазировать, мы могли бы убрать тренд, а затем использовать фильтр с нулевой фазой, чтобы еще больше минимизировать краевые эффекты.

import numpy as np
import matplotlib.pyplot as plt
import scipy.signal as signal

def main():
    x, y = np.loadtxt('spectrum.dat').T

    logx, logy = np.log(x), np.log(y)
    smooth_log = detrend_zero_phase(logx, logy)
    smooth = np.exp(smooth_log)

    fig, ax = plt.subplots()
    ax.loglog(x, y, 'k-')
    ax.loglog(x, smooth, 'r-')
    plt.show()

def zero_phase(y):
    # Low-pass filter...
    b, a = signal.butter(3, 0.05)

    # Filtfilt applies the filter twice to avoid phase shifts.
    return signal.filtfilt(b, a, y)

def detrend_zero_phase(x, y):
    # Fit a second order polynomial (Can't just use scipy.signal.detrend here,
    # because we need to know what the trend is to add it back in.)
    model = np.polyfit(x, y, 2)
    trend = np.polyval(model, x)

    # Apply a zero-phase filter to the detrended curve.
    smooth = zero_phase(y - trend)

    # Add the trend back in
    return smooth + trend

main()

Если все ваши данные медленно изменяются в пространстве журнала, я бы сделал следующее:

1. сильно понижать логарифмические данные в линейной шкале
2. рассчитать сглаживающий сплайн
3. преобразовать обратно в линейный масштаб

e.g.:

import numpy as np
from scipy.interpolate import interp1d, splrep, splev
import pylab

x = np.log10(x)
y = np.log10(y)

ip = interp1d(x,y)
xi = np.linspace(x.min(),x.max(),10)
yi = ip(xi)

tcl = splrep(xi,yi,s=1)
xs = np.linspace(x.min(), x.max(), 100)
ys = splev(xs, tcl)

xi = np.power(10,xi)
yi = np.power(10,yi)
xs = np.power(10,xs)
ys = np.power(10,ys)

f = pylab.figure()
pl = f.add_subplot(111)
pl.loglog(aset.x,aset.y,alpha=0.4)
pl.loglog(xi,yi,'go--',linewidth=1, label='linear')
pl.loglog(xs,ys,'r-',linewidth=1, label='spline')
pl.legend(loc=0)
f.show()

Это дает: