Как построить непрерывный график из дискретных данных?

Question

Есть массив значений значения абсцисс(начало и конец по времени или по расстоянию) и интеграл (сумма значений, например, объем переданных данных или количество человек) [(x1,x2,S1),(x3,x4,S2)....]

Как по нему построить гладкую кривую, такую чтоб площадь под кривой была примерно равна сумме площадей?

Ступенчатый график получил:

import heapq
l = [(1,5,10),(3,6,10),(4,5,10)]
heapq.heapify(l)
def agg(data):
    dec = []
    heapq.heapify(dec)
    y=0
    for x1,x2,s in data:
        while dec and dec[0][0] < x1:
            x2, dy = heapq.heappop(dec)
            yield x2,y
            y-=dy
            yield x2,y
        yield x1,y
        dy=s/(x2-x1)
        y+=dy
        yield x1,y
        dec.append((x2,dy))
    while dec:
        x2, dy = heapq.heappop(dec)
        y-=dy
        yield x2, y

for x,y in agg(l):
    print("({:.2f}, {:.2f})".format(x,y))

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

Другой способ придумал - разбить ступеньку на сектора и пересчитать площадь- построить новую ступеньку, а потом резать углы. Вот код с резанными углами.

import math
import heapq

def agg2(data, step = 1):
    step = step*2
    def monotonic():
        x=math.floor(data[0][0]/step)*step
        while True:
            yield x
            x+=step
    tail = []
    heapq.heapify(tail)
    for x in monotonic():
        s=0
        while tail and tail[0][0] < (x+step):
            x2,s2 = heapq.heappop(tail)
            s+=s2
        while data and data[0][0] < (x+step):
            x1,x2,s0 = data.pop()
            if (x2-(x+step)) > 0:
                s1=((x+step)-x1)/(x2-(x+step)) * s0
                heapq.heappush(tail, (x2,s0-s1))
            else:
                s1 = s0
            s+=s1

        yield x,x+step,s
        if not (data or tail):
            break


def diff(data):
    x = 0
    y = 0
    for x1,x2,s in data:
        dy = s/(x2-x1) - y
        x=x1
        y+=dy/2
        yield x, y
        x=(x2+x1)/2
        y+=dy/2
        yield x, y
        x=x2
        y=s/(x2-x1)/2
    yield x, y


data = [(1,5,10),(3,6,10),(4,5,10), (8,10,10)]
heapq.heapify(data)

for x,y in diff(agg2(data)):
    print("({:.2f}, {:.2f})".format(x,y))

Но неужели нет чего-то готового для таких интерполяций?

Лучше всего реализовать это через какую-то агрегацию на postgresql.

А пока остановлюсь на делении кубиков на монотонные отсчеты и сглажу колоколом на клиенте...

Answer 1

Замените каждый отрезок на функцию "колокола", и просуммируйте результат. Можно использовать разные функции "колокола" (производную от любого сигмоида). Например производная логистической функции:

#  integral( bell, -inf, +inf ) == 1.
bell = lambda x : math.exp(-x) * (1+ math.exp(-x) )**-2
bell_width = 1/bell(0.)
def bell_summ(x, data):
     y=0.
     for x1,x2,s in data:
          y+= s/(x2-x1) * bell( (x-(x1+x2)/2)/(x2-x1) * bell_width )
     return y

data = [(1,5,10),(3,6,10),(4,5,10)]
interpolated = lambda x : bell_summ(x, data)

for x in range(0., 10., 0.1):
     print("({:.2f}, {:.2f})".format( x , interpolated(x) )

Answer 2

На примере ваших данных сделал гладкое распределение средствами KDE - Kernel Density Estimation

Исходные данные: data=[(1.0,5.0,10.0),(3.0,6.0,10.0),(4.0,5.0,10.0)] Для пересекающихся интервалов я просуммировал значения. Получил вот такую нормализованную гистограмму:

Вот результат сглаживания:

Скрипт использует matplotlib.pyplot для генерации изображений и сохранения их в файле. На экран ничего не выводит, генерирует три файла:

bars.png - гистограмма,
pdf.png - сглаженная плотность вероятности (probability density function)
bars_pdf.png - гистограмма и плотность вероятности на одном графике

import numpy as np

# sample data
data=[(1.0,5.0,10.0),(3.0,6.0,10.0),(4.0,5.0,10.0)]

left=min(x1 for x1,_,_ in data)
right=max(x2 for _,x2,_ in data)

# to make the distributes smoother, resample the original distribution
step=0.1
points=np.arange(left,right+step/2,step)
slices = [(points[i],points[i+1]) for i in range(len(points)-1)]
points = points[:-1] # remove last point

def compute_weights(data, slices):
    for s1,s2 in slices:
        total = 0
        for x1,x2,d in data:
            if s2 <= x1: continue
            if x2 <= s1: continue
            total += d/(x2-x1)
        yield (s1, s2, total)

intervals = list(compute_weights(data, slices))
total_weight=sum(w for _,_,w in data)

# normalize the distribution
weights = np.array(list(w/total_weight for _,_,w in intervals))

# Draw the bars
import matplotlib.pyplot as plt

plt.bar(points, weights, step, align='edge')
# save to file rather than display on the screen
plt.savefig('bars.png')

# Approximate the distribution
from scipy import stats

kernel = stats.gaussian_kde(points, weights=weights)
X=np.arange(left-1,right+1,0.1)
Y=kernel(X)

# Erase the plot, draw smoothened distribution
plt.close()
plt.plot(X,Y,'r')
plt.savefig('pdf.png')

# Draw both bar and curve
fig, ax = plt.subplots()
ax.bar(points, weights, step, align='edge')
# red dashed
ax.plot(X,Y,'r--')
plt.savefig('bars_pdf.png')

БЛОГ НА HUSL

Как построить непрерывный график из дискретных данных?

Ответы (2 шт):