Cubic spline interpolation drops out halfway

Question

I am trying to make a cubic spline interpolation and for some reason, the interpolation drops off in the middle of it. It's very mysterious and I can't find any mention of similar occurrences anywhere online.

This is for my dissertation so I have excluded some labels etc. to keep it obscure intentionally, but all the relevant code is as follows. For context, this is an astronomy related plot.

from scipy.interpolate import CubicSpline
import numpy as np
import matplotlib.pyplot as plt

W = np.array([0.435,0.606,0.814,1.05,1.25,1.40,1.60])
sum_all = np.array([sum435,sum606,sum814,sum105,sum125,sum140,sum160])
sum_can = np.array([sumc435,sumc606,sumc814,sumc105,sumc125,sumc140,sumc160])

fall = CubicSpline(W,sum_all)
newallx=np.arange(0.435,1.6,0.001)
newally=fall(newallx)

fcan = CubicSpline(W,sum_can)
newcanx=np.arange(0.435,1.6,0.001)
newcany=fcan(newcanx)

#----plot

plt.plot(newallx,newally)
plt.plot(newcanx,newcany)
plt.plot(W,sum_all,marker='o',color='r',linestyle='')
plt.plot(W,sum_can,marker='o',color='b',linestyle='')
plt.yscale("log")
plt.ylabel("Flux S$_v$ [erg s$^-$$^1$ cm$^-$$^2$ Hz$^-$$^1$]")
plt.xlabel("Wavelength [n$\lambda$]")
plt.show()

The plot that I get from that comes out like this, with a clear gap in the interpolation:

And in case you are wondering, these are the values in the sum_all and sum_can arrays (I assume it doesn't matter, but just in case you want the numbers to plot it yourself):

sum_all:
[  3.87282732e+32   8.79993191e+32   1.74866333e+33   1.59946687e+33
   9.08556547e+33   6.70458731e+33   9.84832359e+33]
can_all:
[  2.98381061e+28   1.26194810e+28   3.30328780e+28   2.90254609e+29
   3.65117723e+29   3.46256846e+29   3.64483736e+29]

The gap happens between [0.606,1.26194810e+28] and [0.814,3.30328780e+28]. If I change the intervals from 0.001 to something higher, it's obvious that the plot doesn't actually break off but merely dips below 0 on the y-axis (but the plot is continuous). So why does it do that? Surely that's not a correct interpolation? Just looking with our eyes, that's clearly not a well-interpolated connection between those two points.

Any tips or comments would be extremely appreciated. Thank you so much in advance!

I do not understand the pgm language but I assume you use single cubic SPLINE which has 4 control points and t^3 as highest order polynomial part so that limits it up to 2 bumps and 1 inflex point. your dataset has more of each so either use higher order SPLINE (but that tends to oscillate too) or use piecewise interpolation. — Spektre
@Spektre scipy.interpolate.CubicSpline does perform a piecewise spline interpolatation. This means that the number of splines is one less than number of points to interpolate. You can therefore have 1 point of inflection per pair of points. — ImportanceOfBeingErnest
@ImportanceOfBeingErnest looks like you got it (+1) I did not spot the log scale .... — Spektre

ImportanceOfBeingErnest ImportanceOfBeingErnest · Accepted Answer · 2017-04-07T08:38:45

The reason for the breakdown can be better observed on a linear scale.

We see that the spline actually passes below 0, which is undefined on a log scale.

So I would suggest to first take the logarithm of the data, perform the spline interpolation on the logarithmically scaled data, and then scale back by the 10th power.

from scipy.interpolate import CubicSpline
import numpy as np
import matplotlib.pyplot as plt

W = np.array([0.435,0.606,0.814,1.05,1.25,1.40,1.60])
sum_all = np.array([  3.87282732e+32,   8.79993191e+32,   1.74866333e+33,   1.59946687e+33,
   9.08556547e+33,   6.70458731e+33,   9.84832359e+33])
sum_can = np.array([  2.98381061e+28,   1.26194810e+28,   3.30328780e+28,   2.90254609e+29,
   3.65117723e+29,   3.46256846e+29,   3.64483736e+29])

fall = CubicSpline(W,np.log10(sum_all))
newallx=np.arange(0.435,1.6,0.001)
newally=fall(newallx)

fcan = CubicSpline(W,np.log10(sum_can))
newcanx=np.arange(0.435,1.6,0.01)
newcany=fcan(newcanx)


plt.plot(newallx,10**newally)
plt.plot(newcanx,10**newcany)
plt.plot(W,sum_all,marker='o',color='r',linestyle='')
plt.plot(W,sum_can,marker='o',color='b',linestyle='')
plt.yscale("log")

plt.ylabel("Flux S$_v$ [erg s$^-$$^1$ cm$^-$$^2$ Hz$^-$$^1$]")
plt.xlabel("Wavelength [n$\lambda$]")
plt.show()

Cubic spline interpolation drops out halfway

1 Answers