如果未做特别说明，文中的程序都是 python3 代码。

QuantLib 金融计算——收益率曲线之构建曲线（3）

载入 QuantLib 和其他包：

import QuantLib as ql
import seaborn as sb
import numpy as np
import pandas as pdprint(ql.__version__)

1.15

概述

本文展示利用 quantlib-python 根据样本券的交易数据估算出即期利率的期限结构的完整流程，并指出当前实现所存在的问题。

示例所用的样本券交易数据来自专门进行期限结构分析的 R 包——termstrc。具体来说是数据集 govbonds 中的 GERMANY 部分，包含 2008-01-30 这一天德国市场上 52 只固息债的成交数据。

注意：为了适配 QuantLib，实际计算中删除了两只债券的数据，以保证所有样本券的到期时间均不相同。样本券数据在附录中列出。

估算期限结构的步骤

QuantLib 中估算期限结构的核心流程有两步：

1. 配置 *Helper 对象，描述样本券信息，包括付息时间表（schedule）、价格（默认用净价）、票息等；
2. 配置期限结构模型，可以额外提供样本券权重、优化方法、参数正则化条件等参数辅助计算。

读取样本券数据

govBond = pd.read_csv('GERMANY_INFO.csv',parse_dates=['MATURITYDATE', 'ISSUEDATE'])numberOfBonds = govBond.shape[0]PRICE = [ql.QuoteHandle(ql.SimpleQuote(p)) for p in govBond['PRICE']]
MATURITYDATE = [ql.Date(m.day, m.month, m.year) for m in govBond['MATURITYDATE']]
ISSUEDATE = [ql.Date(i.day, i.month, i.year) for i in govBond['ISSUEDATE']]
COUPONRATE = govBond['COUPONRATE'].values

一些基本配置

# 查看 govbonds 数据集可知样本券均为每年付息一次
frequency = ql.Annual
# termstrc 的日期计算并不如 QuantLib 精细，
# 为了和 termstrc 的算法保持一致，示例使用如下天数计算规则
dc = ql.Actual365Fixed(ql.Actual365Fixed.Standard)
paymentConv = ql.Unadjusted
terminationDateConvention = ql.Unadjusted
convention = ql.Unadjusted
redemption = 100.0
faceAmount = 100.0
# 其实我不知道样本券所在的交易所，
# 所以不确定是不是该用这个日历 :)
calendar = ql.Germany(ql.Germany.Eurex)# 估值日期 2008-01-30
today = calendar.adjust(ql.Date(30, 1, 2008))
ql.Settings.instance().evaluationDate = today# 为了和 termstrc 的算法保持一致，示例采用 T+0 结算
bondSettlementDays = 0
bondSettlementDate = calendar.advance(today,ql.Period(bondSettlementDays, ql.Days))

配置 *Helper 对象

instruments = []for j in range(numberOfBonds):# 配置付息时间表schedule = ql.Schedule(ISSUEDATE[j],MATURITYDATE[j],ql.Period(frequency),calendar,convention,terminationDateConvention,ql.DateGeneration.Backward,False)# 配置 Helper 对象# 因为样本券均为固息债，所以采用 FixedRateBondHelper 类# 对于其他金融工具，需要使用对应的 Helper 类helper = ql.FixedRateBondHelper(PRICE[j],bondSettlementDays,faceAmount,schedule,[COUPONRATE[j]],dc,paymentConv,redemption)instruments.append(helper)

配置期限结构

tolerance = 1.0e-6
max = 5000# 即期利率的 Svensson 模型
sf = ql.SvenssonFitting()
# 即期利率的 Nelson Siegel 模型
nsf = ql.NelsonSiegelFitting()
# 用指数样条函数拟合贴现因子
esf = ql.ExponentialSplinesFitting()
# 用简单多项式函数拟合贴现因子
spf = ql.SimplePolynomialFitting(8)
# 用三次 B-样条函数拟合贴现因子
knots = [-20.0, -10.0, 0.0, 0.25, 0.5, 1, 3, 5, 10, 20, 30, 40, 50]
cbsf = ql.CubicBSplinesFitting(knots)tsSvensson = ql.FittedBondDiscountCurve(bondSettlementDate, instruments, dc,sf,tolerance, max)tsNelsonSiegel = ql.FittedBondDiscountCurve(bondSettlementDate, instruments, dc,nsf,tolerance, max)tsExponentialSplines = ql.FittedBondDiscountCurve(bondSettlementDate, instruments, dc,esf,tolerance, max)tsSimplePolynomial = ql.FittedBondDiscountCurve(bondSettlementDate, instruments, dc,spf,tolerance, max)tsCubicBSplines = ql.FittedBondDiscountCurve(bondSettlementDate, instruments, dc,cbsf,tolerance, max)

估算期限结构

sv = []
ns = []
es = []
sp = []
cbs = []
matList = []matDate = bondSettlementDatewhile matDate <= bondSettlementDate + ql.Period(31, ql.Years):matDate = matDate + ql.Period(1, ql.Days)matList.append(dc.yearFraction(bondSettlementDate, matDate))sv.append(tsSvensson.zeroRate(matDate, dc, ql.Continuous, frequency).rate() * 100)ns.append(tsNelsonSiegel.zeroRate(matDate, dc, ql.Continuous, frequency).rate() * 100)es.append(tsExponentialSplines.zeroRate(matDate, dc, ql.Continuous, frequency).rate() * 100)sp.append(tsSimplePolynomial.zeroRate(matDate, dc, ql.Continuous, frequency).rate() * 100)cbs.append(tsCubicBSplines.zeroRate(matDate, dc, ql.Continuous, frequency).rate() * 100)

汇总结果

# 以 termstrc 的估算结果作为比较基准
beta0 = 5.017052
beta1 = -1.117214
beta2 = -3.173622
tau = 2.443936
termstrc = [beta0 + \beta1 * (1 - np.exp(-m / tau)) / (m / tau) + \beta2 * ((1 - np.exp(-m / tau)) / (m / tau) - np.exp(-m / tau)) for m in matList]df = pd.DataFrame(dict(maturity=matList * 6,rate=sv + ns + es + sp + cbs + termstrc,type=np.repeat(['Svensson', 'NelsonSiegel', 'ExponentialSplines','SimplePolynomial', 'CubicBSplines', 'termstrc'], len(matList))))print(tsSvensson.fitResults().solution())
print(tsNelsonSiegel.fitResults().solution())
print(tsExponentialSplines.fitResults().solution())
print(tsSimplePolynomial.fitResults().solution())
print(tsCubicBSplines.fitResults().solution())print(tsSvensson.fitResults().minimumCostValue())
print(tsNelsonSiegel.fitResults().minimumCostValue())
print(tsExponentialSplines.fitResults().minimumCostValue())
print(tsSimplePolynomial.fitResults().minimumCostValue())
print(tsCubicBSplines.fitResults().minimumCostValue())sb.relplot(x='maturity', y='rate', kind='line', hue='type',size='type', sizes=[2, 2, 2, 2, 2, 4],data=df, height=5, aspect=1.6)

[ -134.509; 134.548; 136.375; -0.0234826; 0.000905134; 0.422332 ]
[ -9.7222; 9.75752; 38.0528; 3.14624e-05 ]
[ -31238.6; -95904.4; 46747.5; 53652; 21350.9; 2438.62; -77566.2; 36236.2; 0.000302241 ]
[ -0.0373576; 0.00165862; -0.000203106; 9.83472e-06; -1.47693e-07 ]
[ 1.13572; 0.977843; 0.948002; 0.90216; 0.798218; 0.616569; 0.355499; 0.293158 ]0.00032255014803192505
0.0020683101877706977
0.0019235041138322383
0.0008108667989534623
0.00037198366450668654

图 1：QuantLib 的结果

图 2：termstrc 的结果

注意：尽管以 termstrc 的结果作为基准，并不意味着基准就是正确答案。

NelsonSiegel、SimplePolynomial 和 ExponentialSplines 的结果与基准相去甚远。Svensson 和 CubicBSplines 的结果在短端与基准非常接近，但在长端依然有明显差距，Svensson 和 CubicBSplines 的结果要略低于基准。

考虑到基准似乎在长端高估了真实利率水平，Svensson 和 CubicBSplines 的结果可能要好于基准。另外，CubicBSplines 甚至顾及到了 0 附近的两个“异常值”。

当前实现存在的问题与对策

粗看结果似乎还可以接受，但实际上经不起推敲。

• 首先，根据 Nelson Siegel 模型和 Svensson 模型的经济意义（参考文献 1），估计结果绝对值的数量级应该和利率处于同一水平，通常是 10 以内的某个数。尤其是模型中的 \(\beta_0\)，大致应该等于利率的平均值。
• 其次，三次 B-样条的结果在远端出现了 S 形的扭曲，可能是 knot 选择不当的结果，最终结果对 knot 的选择其实非常敏感（参考文献 3）。另外，QuantLib 采用了 \(d(t) = \sum_{i=0}^{n} c_i \times N_{i,3}(t)\) 的格式（\(N_{i,3}(t)\) 是基本样条函数，\(d(t)\) 是贴现因子），因为 \(N_{i,3}(t)\) 在最外测的两个 knot 上取值非常小（MATLAB 的演示），这使得使用者必须提供而外的 knot 完全覆盖当前的期限范围才能有合理的估计，相当反人类的设计。
• 第三，指数样条等方法的参数过于极端。

所有问题的根源是相同的，因为估算期限结构本质上是一个优化问题。以 Nelson Siegel 模型和 Svensson 模型为例，参数估计本身是一个相当有挑战性的非凸优化问题（参考文献 2），可能需要借助一些特殊的技术手段（参考文献 2），而不是依赖于某个优化算法。但是 quantlib-python 在封装期限结构接口的时候只保留了样本券权重一个自由度，优化算法、正则化条件等选项均被忽略，特别是优化算法，统一使用比较原始的单纯形算法。

若要改良当前的结果，一种方法是编写 C++ 程序使用其他优化算法，并配置正则化条件；另一种方法是自定义 swig 的接口文件，修改 quantlib-python 期限结构类的接口，使其能使用其他优化算法，并接受正则化条件。

如果无法实现以上两种方法，在当前有限的条件下推荐使用三次 B-样条函数估算期限结构，但要注意 knot 的选择。

参考文献

1. 《收益率曲线的建模和预测——基于 DNS 方法创新》，东北财经大学出版社
2. R. Ferstl and J. Hayden, "Zero-Coupon Yield Curve Estimation with the Package termstrc" Journal of Statistical Software, August 2010, Volume 36, Issue 1.
3. James, J. and N. Webber, "Interest Rate Modelling" John Wiley, 2000.

附录

样本券数据。