13

Say I have some_data of shape (1, n). I have new incoming_data of shape (1, n±x), where x is some positive integer much smaller than n. I would like to squeeze or stretch incoming_data such that it is of the same length as n. How might this be done, using the SciPy stack?

Here's an example of what I'm trying to accomplish.

# Stretch arr2 to arr1's shape while "filling in" interpolated value
arr1 = np.array([1, 5, 2, 3, 7, 2, 1])
arr2 = np.array([1, 5, 2, 3, 7, 1])
result
> np.array([1, 5, 2, 3, 6.x, 2.x 1])  # of shape (arr1.shape)

As another example:

# Squeeze arr2 to arr1's shape while placing interpolated value.
arr1 = np.array([1, 5, 2, 3, 7, 2, 1])
arr2 = np.array([1, 5, 2, 3, 4, 7, 2, 1])
result
> np.array([1, 5, 2, 3.x, 7.x, 2.x, 1])  # of shape (arr1.shape)
ericmjl
  • 13,541
  • 12
  • 51
  • 80
  • FWIW, I have tried using `scipy.interpolate`, but I think I'm missing something to get the interpolation done right, because I keep getting an error because the arrays are not of the same length. – ericmjl Jun 27 '16 at 23:19
  • 1
    Use `interpolate` with `x=np.arange(arr2.size)` and `arr2` as basis points, and interpolate to new `x` values given by `np.linspace(0,arr2.size-1,arr1.size)`. These elements should be your interpolated values. – Andras Deak -- Слава Україні Jun 27 '16 at 23:21
  • Ooh! I'm feeling real close to the answer :D. I just tried doing what you wrote, but still hitting a block. What goes into the `scipy.interp1d(....)` function call? – ericmjl Jun 27 '16 at 23:27
  • Look at its manual. Try either `help(scipy.interpolate.interp1d)`, or [the same online](http://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.interp1d.html). – Andras Deak -- Слава Україні Jun 27 '16 at 23:28
  • Ok, got it - so with your solution, I end up having to keep track of which is the longer array... unless I'm doing something wrong? (I'm pretty sure I am doing something wrong...) I have a vector which is of a reference length (in the above example, it is `arr1`), and I'd like to interpolate everything to be of the reference length. – ericmjl Jun 28 '16 at 00:05
  • I'm a noob, please ignore my previous comments. I will post below my answer, and acknowledge your help, @AndrasDeak! – ericmjl Jun 28 '16 at 00:09
  • OK, I was just about to explain in an answer:) Go ahead with yours. – Andras Deak -- Слава Україні Jun 28 '16 at 00:11
  • Actually, I'd love to see your answer instead :). Please go ahead! – ericmjl Jun 28 '16 at 00:13
  • OK, I answered, thanks:) I made a picture to make it worth your while;) – Andras Deak -- Слава Україні Jun 28 '16 at 00:45

2 Answers2

23

You can implement this simple compression or stretching of your data using scipy.interpolate.interp1d. I'm not saying it necessarily makes sense (it makes a huge difference what kind of interpolation you're using, and you'll generally only get a reasonable result if you can correctly guess the behaviour of the underlying function), but you can do it.

The idea is to interpolate your original array over its indices as x values, then perform interpolation with a sparser x mesh, while keeping its end points the same. So essentially you have to do a continuum approximation to your discrete data, and resample that at the necessary points:

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

arr_ref = np.array([1, 5, 2, 3, 7, 1])  # shape (6,), reference
arr1 = np.array([1, 5, 2, 3, 7, 2, 1])  # shape (7,), to "compress"
arr2 = np.array([1, 5, 2, 7, 1])        # shape (5,), to "stretch"
arr1_interp = interp.interp1d(np.arange(arr1.size),arr1)
arr1_compress = arr1_interp(np.linspace(0,arr1.size-1,arr_ref.size))
arr2_interp = interp.interp1d(np.arange(arr2.size),arr2)
arr2_stretch = arr2_interp(np.linspace(0,arr2.size-1,arr_ref.size))

# plot the examples, assuming same x_min, x_max for all data
xmin,xmax = 0,1
fig,(ax1,ax2) = plt.subplots(ncols=2)
ax1.plot(np.linspace(xmin,xmax,arr1.size),arr1,'bo-',
         np.linspace(xmin,xmax,arr1_compress.size),arr1_compress,'rs')
ax2.plot(np.linspace(xmin,xmax,arr2.size),arr2,'bo-',
         np.linspace(xmin,xmax,arr2_stretch.size),arr2_stretch,'rs') 
ax1.set_title('"compress"')
ax2.set_title('"stretch"')

The resulting plot:

result

In the plots, blue circles are the original data points, and red squares are the interpolated ones (these overlap at the boundaries). As you can see, what I called compressing and stretching is actually upsampling and downsampling of an underlying (linear, by default) function. This is why I said you must be very careful with interpolation: you can get very wrong results if your expectations don't match your data.

  • 1
    It works perfectly, Andras! And FWIW, I ended up going with the [Pchipinterpolator](http://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.PchipInterpolator.html#scipy.interpolate.PchipInterpolator), simply because it allows for extrapolation, which specifically for my application, is what I need. – ericmjl Jun 28 '16 at 01:12
  • 2
    @ericmjl I'm glad I could help:) Thanks for the feedback. Be especially careful when extrapolating from interpolating functions, that can lead to uncontrolled errors. (Also, for completeness: as of scipy version 0.17 `interp1d` supports extrapolation for linear interpolation, with the keyword `full_value='extrapolate'`.) – Andras Deak -- Слава Україні Jun 28 '16 at 01:16
  • Also worth mentioning the `resampy` package. If your arrays are representations of audio, that package performs good interpolation faster than scipy. – Jake Stevens-Haas Feb 06 '20 at 21:14
  • 1
    @JakeStevens-Haas I'm not familiar with that package at all, but if it can do what the question is asking for you should consider adding another answer with it. – Andras Deak -- Слава Україні Feb 06 '20 at 23:26
7

There's another package that works very well for upsampling and downsampling: resampy. It has a simpler command than scipy.interpolate.interp1d but only uses a single interpolation function. As @Andras Deak said, you have to be careful in choosing interpolation functions.

MWE:

import numpy as np
import resampy
from matplotlib import pyplot as plt

x_mesh = np.linspace(0,1,10)
short_arr = np.sin(x_mesh*2*np.pi)
plt.plot(short_arr)

coarse_plot

interp_arr = resampy.resample(short_arr, 20, 100)
plt.plot(interp_arr)

fine_plot
Two words of caution:

  1. resampy uses a "band-limited sinc interpolation". Check the documentation for more info. It works best if your array originally came from data with local frequency components, e.g. sound, images, and other time-series data. It's used in some of the tensorflow examples on audio, which is what I use. I'm not sure whether your example array was small for demonstration purposes, but if that truly is the size of your array, interpolating may be bad whatever method you use, linear, spline, or otherwise.

  2. Your examples demonstrated more than interpolation. It seems you found a portion of the arrays that matched (e.g. [1,5,2,3]) then interpolated the rest. Depending on whether you want to match the beginning of the array or an arbitrary number of patches, you may be asking for a two methods: one to identify the correct portions of an array to interpolate, and one to interpolate those portions. If that's the case, look at numpy.isin for a basic method or levenshtein distance for more generally matching a set of substrings.

Jake Stevens-Haas
  • 1,186
  • 2
  • 14
  • 26