Как выровнять два массива разной длины в Python (используя NaN, где нет соответствующего элемента)

Question

У меня есть массив y, который содержит значения, наблюдаемые в данный день месяца. День месяца в массиве х.

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

Если я хочу построить оригинал y и интерполированный y (т.е. yd), мне нужно выровнять их по одной оси. Это ось, учитывающая весь день месяца, xd.

Существует ли эффективный способ быстрого создания нового массива y, который содержит точно исходный элемент y в нужном месте на основе новой оси x со всеми другими элементами, заполненными нулями или NaN (предпочтительно)?

Так, например, мой первый y доступен только во второй день, поэтому в новом массиве y мне нужен первый элемент, чтобы показать 0 / NaN. Тогда второй элемент покажет оригинал y = 11, третий покажет NaN и т. Д.

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

import numpy as np
import scipy.interpolate as sp

x = [2, 5, 7, 11, 13, 16, 19, 23, 25, 30]
y = [11, 10, 12, 14, 16, 19, 17, 14, 18, 17]

xd = np.linspace(0, max(x), int(max(x))+1) # create the new x axis
ipo = sp.splrep(x, y, k=3) # cubic spline
yd = sp.splev(xd, ipo) # interpolated y values

newY = np.zeros((1, len(yd)), dtype=float) # preallocate for the filled y values

for i in x: 
    if(i in xd): 
        idx, = np.where(xd == i) # find where the original x value is in the new x axis
        idx2, = np.where(np.array(x) == i)
        newY[0, int(idx)] = y[int(idx2)] # replace the y value of the new vector with the y value from original set

РЕДАКТИРОВАТЬ:

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

Я знаю, что в этом случае сплайн всегда будет проходить через все точки, которые я даю в качестве входных данных, поэтому различия будут равны нулю, но нижеприведенная функция построения графика должна работать с любым видом сравнения (т. Е. С любым типом интерполированных значений против реальный вклад). Используемая ниже функция построения графика:

def plotInterpolatedVsReal(xaxis, yaxis1, yaxis2, xlab='Dates', mainTitle='', width=25, zero2nan=True):
    if(zero2nan):
        yaxis1[yaxis1 == 0] = np.nan
        yaxis2[yaxis2 == 0] = np.nan

    fix, (ax1, ax2, ax3) = plt.subplots(3, sharex=True, figsize=(10, 10))
    ax1.plot(xaxis, yaxis1, label='Interpolated')
    ax1.plot(xaxis, yaxis2, 'ro', label='Input')
    ax1.set_ylabel('Prices')
    ax1.legend(loc=0)
    ax2.bar(xaxis, yaxis1 - yaxis2, width=width)
    ax2.axhline(y=0, linewidth=1, color='k')
    ax2.set_ylabel('Errors [diff]')
    ax3.bar(xaxis, 100*(yaxis1/yaxis2 - 1), width=width)
    ax3.axhline(y=0, linewidth=1, color='k')
    ax3.set_ylabel('Errors [%]')
    ax3.set_xlabel(xlab);
    plt.suptitle(mainTitle)

РЕДАКТИРОВАТЬ 2:

Добавление показателей эффективности для данного предложения. Мой цикл (метод A) быстрее, потому что он зацикливается только по вектору x, в то время как другие 2 метода зацикливаются по xd, который может быть значительно больше. В моем случае здесь x имеет 23 элемента, а xd имеет 3655 элементов.

def A():
    for i in x: 
        if(i in xd): 
            idx, = np.where(xd == i) # find where the original x value is in the new x axis
            idx2, = np.where(np.array(x) == i)
            newY[int(idx)] = y[int(idx2)] # replace the y value of the new vector with the y value from original set 

def B():
    for i, date in enumerate(xd):
        if date in x:
            new_y[i] = date

def C(): 
    known_values = dict(zip(x, y))

    for i,u in enumerate(xd):
        if u in known_values:
            newY[i] = known_values[u]

% времени A () 219 мкс ± 8,8 мкс на цикл (среднее ± стандартное отклонение из 7 циклов, 1000 циклов в каждом)

% время B () 8,87 мс ± 95,3 мкс на цикл (среднее ± стандартное отклонение из 7 циклов, по 100 циклов в каждом)

% времени этого C () 408 мкс ± 11,3 мкс на цикл (среднее ± стандартное отклонение из 7 циклов, по 1000 циклов в каждом)

Я также пытался передать свою функцию A () Numba для компиляции JIT:

A_nb = numba.jit(A)

получение:

% timeit A_nb () 226 мкс ± 610 нс на цикл (среднее ± стандартное отклонение 7 циклов, 1000 циклов каждый)

Answer 1

Извините, если я полностью неправильно понял ваш код, но разве np.linspace(0, max(x), int(max(x))+1) не является просто окольным способом написания np.array(range(1+max(x)))? Выглядит так, как будто вы просто берете 1+max(x) выборок с линейным интервалом в диапазоне от 0 до max(x) включительно, что равносильно тому, что вы берете целые числа от 0 до max (x).

И в таком случае, нужно ли это делать?

if(i in xd): 
    idx, = np.where(xd == i) # find where the original x value is in the new x axis

Если xd действительно является просто списком целых чисел от 0 до max (x) включительно, то все элементы в x будут в xd по определению, а idx всегда должно быть равно i.

(Предполагается, что x, конечно, содержит только неотрицательные целочисленные значения.)

xd = np.array(range(1+max(x)))
newY = np.zeros(len(xd))

for i,j in zip(x, y):
    newY[i] = j

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

known_values = dict(zip(x, y))

xd = [... your new axis ...]
newY = np.zeros(len(xd))

for i,x in enumerate(xd):
    if x in known_values:
        newY[i] = known_values[x]

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

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

k = 0
for i,j in zip(x,y):
    while k < len(xd) and xd[k] < i:
        k += 1
    if xd[k] == i:
        newY[k] = j

Answer 2

Я понимаю, что смысл всего этого - отображать значения y на одном графике, почему бы не сделать это напрямую?Оси могут легко обрабатывать разные оси x на одном и том же графике, например:

import numpy as np
import scipy.interpolate as sp
import matplotlib.pyplot as plt

x = [2, 5, 7, 11, 13, 16, 19, 23, 25, 30]
y = [11, 10, 12, 14, 16, 19, 17, 14, 18, 17]

xd = np.linspace(0, max(x), int(max(x)) + 1)  # create the new x axis
ipo = sp.splrep(x, y, k=3)  # cubic spline
yd = sp.splev(xd, ipo)  # interpolated y values

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x, y, label='Original')
ax.plot(xd, yd, label='Interpolated')
plt.legend()
plt.grid()

plt.show()

Как вы и хотели, все данные 'y' выровнены по своей оси X без предварительной обработки,Единственная интерполяция, выполняемая здесь, - это та, которую Matplotlib выполняет для отображения.

Поскольку вам действительно нужно заполнить свой массив с помощью Nan, вот рабочий способ сделать это:

new_y = np.NAN * np.zeros(yd.shape)
for i, date in enumerate(xd):
    if date in x:
        new_y[i] = date

вероятно, может быть уменьшен в некотором причудливом одном вкладыше