在Python中获取“bs”（样条曲线）等价项

Question

在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等效的版本？我试过从scipy开始BSpline，但它要求系数已经知道并传入

在传递了数组、节点和度之后，如何生成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])

Answer 1

将评论转化为答案，BSpline.design_matrix是以csr稀疏格式构建您所追求的内容。它将在发布时从scipy 1.8中提供。在此之前，您可以获取scipy的主分支，或者使用docs（https://scipy.github.io/devdocs/reference/generated/scipy.interpolate.BSpline.design_matrix.html#scipy.interpolate.BSpline.design_matrix）建议的解决方法：

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

EDIT:R文档https://www.rdocumentation.org/packages/splines/versions/3.6.2/topics/bs指出knots参数是内部节点。scipy的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:, :]或类似的内容）