Skip to content
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
45 changes: 21 additions & 24 deletions splinecalib/calib_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,31 +8,28 @@
from scipy.stats import binom
from loss_fun_c import pen_ll_fun, pen_ll_fun_grad

def _natural_cubic_spline_basis_expansion(xpts, knots):
def _natural_cubic_spline_basis_expansion(z, knots):
"""Does the natural cubis spline bases for a set of points and knots"""
num_knots = len(knots)
num_pts = len(xpts)
outmat = np.zeros((num_pts,num_knots))
outmat[:, 0] = np.ones(num_pts)
outmat[:, 1] = xpts

def make_func_H(k):
def make_func_d(k):
def func_d(x):
denom = knots[-1] - knots[k-1]
numer = (np.maximum(x-knots[k-1], np.zeros(len(x))) ** 3 -
np.maximum(x-knots[-1], np.zeros(len(x))) ** 3)
return numer/denom
return func_d

def func_H(x):
d_fun_k = make_func_d(k)
d_fun_Km1 = make_func_d(num_knots-1)
return d_fun_k(x) - d_fun_Km1(x)
return func_H
for i in range(1, num_knots-1):
curr_H_fun = make_func_H(i)
outmat[:, i+1] = curr_H_fun(xpts)
nr_knots = knots.shape[0]
nr_data = z.shape[0]

# we explicitly make z a column vector to facilitate broadcasting, namely in the
# subtraction z - knots[:nr_knots - 1] -> array shape [nr_data, nr_knots - 1]
z = z[:, np.newaxis]

outmat = np.zeros((nr_data, nr_knots), dtype=np.float64)
outmat[:, 0] = np.ones(nr_data, dtype=np.float64)
outmat[:, 1] = z.squeeze()

final_cubic_term = non_negative((z - knots[-1])) ** 3
numerator = (
non_negative((z - knots[: nr_knots - 1])) ** 3 - final_cubic_term
) # array shape (nr_data, nr_knots-1)
denominator = knots[-1] - knots[: nr_knots - 1] # array shape (nr_knots-1, )
d = numerator / denominator # array shape (nr_data, nr_knots-1)
phi = d[:, :-1] - d[:, [-1]] # n.b. the [] around -1 req'd to retain the dimension

outmat[:, 2:] = phi
return outmat


Expand Down