I'm using the UnivariateSpline from SciPy, and it appears to me that the second formula it gives for the smoothing condition (for the s parameter) is incorrect. The first formula it gives to decide when to stop adding knots to the spline seems fine: sum((w * (y-spl(x)))**2, axis=0) <= s This is fine (despite the syntax not being proper Python code due to the i and using Numpy sum). After stating the above formula, they then state, "However, because of numerical issues, the actual condition is" and give this second formula: abs(sum((w * (y-spl(x)))**2, axis=0) - s) < 0.001 * s This does not seem correct at all to me. (If you care, I'm a mathematician.) The issue is as follows. Define the error as error = sum((w * (y-spl(x)))**2, axis=0) Their second formula is saying abs(error - s) < 0.001 * s This is the same (mathematically speaking) as saying that 0.999 * s < error < 1.001 * s This means that error is approximately s. Contrast this to their first formula, which is just saying error <= s. Am I missing something, or is their second formula not really the stopping condition for smoothing? It appears that the call to the absolute value function just shouldn't be there. Without it, my modified second formula would be equivalent to this (which makes a lot more sense to me): error < 1.001 * s I checked their claim on their stopping condition and found that the error is sometimes made to be outside the interval 0.999 * s < error < 1.001 * s: from scipy.interpolate import UnivariateSpline x = range(4) y = [0.01, 1.01, 2.01, 3.01] spl = UnivariateSpline(x, y, s=1) print(spl(x)) print(sum((yi-spl(xi))**2 for xi, yi in zip(x,y))) # prints this: #[0.01 1.01 2.01 3.01] #3.6412113011071178e-34 I tried looking up the source code and found the UnivariateSpline class, but it appears to me that the heart of their algorithm is the call to dfitpack.fpcurf0 which I haven't found the source code for. Continue reading...
