Блог

На r языке программирования следующие

 require(stats)
require(splines)
knots = quantile(women$height, seq(0.1,0.9,length.out = 5))
bs(women$height, knots=knots, degree=3)
 

ВОЗВРАТ

 1   2   3   4   5   6   7   8
0.0000000000    0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.0000000000    0.00000000
0.6284418529    0.323939099 0.024295432 0.000000000 0.000000000 0.000000000 0.0000000000    0.00000000
0.2155814707    0.599894720 0.182883868 0.001639942 0.000000000 0.000000000 0.0000000000    0.00000000
0.0349854227    0.495626822 0.438289602 0.031098154 0.000000000 0.000000000 0.0000000000    0.00000000
0.0001619695    0.245586330 0.620809038 0.133442663 0.000000000 0.000000000 0.0000000000    0.00000000
0.0000000000    0.072886297 0.584548105 0.338678328 0.003887269 0.000000000 0.0000000000    0.00000000
0.0000000000    0.009110787 0.384718173 0.561892614 0.044278426 0.000000000 0.0000000000    0.00000000
0.0000000000    0.000000000 0.166666667 0.666666667 0.166666667 0.000000000 0.0000000000    0.00000000
0.0000000000    0.000000000 0.044278426 0.561892614 0.384718173 0.009110787 0.0000000000    0.00000000
0.0000000000    0.000000000 0.003887269 0.338678328 0.584548105 0.072886297 0.0000000000    0.00000000
0.0000000000    0.000000000 0.000000000 0.133442663 0.620809038 0.245586330 0.0001619695    0.00000000
0.0000000000    0.000000000 0.000000000 0.031098154 0.438289602 0.495626822 0.0349854227    0.00000000
0.0000000000    0.000000000 0.000000000 0.001639942 0.182883868 0.599894720 0.2155814707    0.00000000
0.0000000000    0.000000000 0.000000000 0.000000000 0.024295432 0.323939099 0.6284418529    0.02332362
0.0000000000    0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.0000000000    1.00000000
 

Существует ли эквивалент Python? Я пробовал BSpline от scipy , но для этого требуется, чтобы коэффициенты уже были известны и передавались.

Как я могу просто сгенерировать базовую матрицу B-сплайна, передав массив, узлы и степень?

Чтобы воспроизвести входной Python, вы можете сделать:

 import numpy as np

women_height = np.array([58,59,60,61,62,63,64,65,66,67,68,69,70,71,72])
knots = array([59.4, 62.2, 65. , 67.8, 70.6])
 

Превращение комментария в ответ- BSpline.design_matrix это построение того, что вам нужно, в разреженном формате csr. Он будет доступен с scipy 1.8, когда он будет выпущен. До тех пор вы можете либо захватить главную ветвь scipy, либо воспользоваться обходным путем, предложенным в документах (https://scipy.github.io/devdocs/reference/generated/scipy.interpolate.BSpline.design_matrix.html#scipy.интерполировать.BSpline.design_matrix) :

 t = ... 
c = np.eye(len(t) - k - 1)
design_matrix_gh = BSpline(t, c, k)(x)
 

ПРАВКА: документация R, https://www.rdocumentation.org/packages/splines/versions/3.6.2/topics/bs, утверждает, что knots аргумент является внутренним узлом. scipy’s BSpline не накладывает узлы автоматически, поэтому вам нужно сделать это самостоятельно. Использование данных OP:

 In [22]: women_height = np.array([58,59,60,61,62,63,64,65,66,67,68,69,70,71,72])
    ...: knots = np.array([59.4, 62.2, 65. , 67.8, 70.6])

In [23]: t = np.r_[(58,)*4, knots, (72,)*4]   # <<<<<< here

In [24]: m = BSpline.design_matrix(women_height, t, k=3)

In [25]: with np.printoptions(linewidth=120, precision=5):
    ...:     print(m.toarray())
    ...:
[[1.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00]
 [2.33236e-02 6.28442e-01 3.23939e-01 2.42954e-02 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 2.15581e-01 5.99895e-01 1.82884e-01 1.63994e-03 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 3.49854e-02 4.95627e-01 4.38290e-01 3.10982e-02 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 1.61970e-04 2.45586e-01 6.20809e-01 1.33443e-01 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 0.00000e 00 7.28863e-02 5.84548e-01 3.38678e-01 3.88727e-03 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 0.00000e 00 9.11079e-03 3.84718e-01 5.61893e-01 4.42784e-02 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 1.66667e-01 6.66667e-01 1.66667e-01 0.00000e 00 0.00000e 00 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 4.42784e-02 5.61893e-01 3.84718e-01 9.11079e-03 0.00000e 00 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 3.88727e-03 3.38678e-01 5.84548e-01 7.28863e-02 0.00000e 00 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 1.33443e-01 6.20809e-01 2.45586e-01 1.61970e-04 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 3.10982e-02 4.38290e-01 4.95627e-01 3.49854e-02 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 1.63994e-03 1.82884e-01 5.99895e-01 2.15581e-01 0.00000e 00]
 [0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 2.42954e-02 3.23939e-01 6.28442e-01 2.33236e-02]
 [0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 0.00000e 00 1.00000e 00]]

 

что похоже на вывод R из операции по модулю первого столбца. Из документов R не сразу понятно, как именно он заполняет вектор узлов, но если вы хотите получить тот же результат, вы можете просто отрезать первый столбец ( m.toarray()[1:, :] или что-то в этом роде).