【中金固收·固收＋】久期测算的探索：细节处理与Python实践

债基的久期一直不是容易计算的问题，否则市场上也不会有如此多的讨论和尝试。实际上这不是一个学术问题——模型不会设计得很复杂，债券投资者也往往不太接受这种尝试——而是存在不少细节问题需要处理。比如，基金净值本就是锯齿状的“准连续数据”，再如债券的波动性要比股票低得多，导致噪声的比例要远比股票仓位模型大（例如我们之前介绍过的固收+基金风格分解模型，在这里就无法照搬）。其中一些问题经常被忽略，比如有的算法不稳定，于是用N日均线或者只看全市场中位数来掩盖——这也是为何，中位数的走势比均值更常见。好在这些问题也都并非无解，在此我们介绍一种稳定性与时效性兼顾的方法。

图表: 久期算法框架

资料来源：万得资讯、中金公司研究部

在解释算法前，首先我们要明确可能遇见的问题：

1、如前所述，债券波动要比股市小得多，这也意味着所有可能作为分母的数据，我们都要更加小心——有的算法可能会出现市场波动越大（往往在债券牛市）基金久期越短，进而出现了一种“市场择时能力很强”的假象，就是因为波动越大时，不容易出现“分母太小”的情况；

图表: 指数波动小带来的数字陷阱，呈现“逆势交易”假象

资料来源：万得资讯、中金公司研究部

2、债基净值虽然每天发布，但最多只能精确到小数点后四位（亦有精确到后三位的），再加上债券本身波动就很小，导致这个数据实际介于连续和离散数据之间的形态；

图表: 基金净值曲线举例

资料来源：万得资讯、中金公司研究部

3、为了区分不同产品的风格差异，我们要用到多条债券指数及其久期数据，但这些指数之间有着很强的多重共线性——在此前固收+基金风格分解中，我们用到了分步过滤法来解决这个问题；

4、基于回归的方法都会出现“拖尾问题”，即当日的结果是用过去很长的交易日数据拟合的，这样存在时效性问题。

综上考虑，我们选择了基于分步过滤和波动比较的方法缓解上述问题。简单来说，步骤包括：

1、获取基金净值和债券指数的数据：为了缓解“锯齿净值”的问题，我们对这些数据的涨跌幅都做了滚动窗口加总（即过去n日加总）。在考虑滚动窗口时，我们观察了数据的稳定性，并最终选择8日。为方便调用，我们在程序实现时，先把这些数据（连同后面的计算公式）都封装进fundDuration类中，程序逻辑如下：

图表: 初始化数据

class fundDuration (object):

    def __init__(self, codes, start, end, n=30, rol=8):
        # 基金代码，起始日、结束日，回归窗口区间、滚动长度
        self.codes = codes
        self.start, self.end = start, end
        self.rol = rol
        self.n = n

        self.dfNavRt = self.getNavRt(n) #净值的原始值
        self.dfIndex = self.getBondIndex(n) #指数的原始值

        self.dfNavRtRoll = self.dfNavRt.rolling(rol).sum().dropna() # 获取基金净值波幅表，去空值，rol日滚动窗口化
        self.dfIndexRoll = self.dfIndex.rolling(rol).sum().dropna() # 获取债券指数变动波幅表，去空值，rol日滚动窗口化

        self.dfDur = self.getBondDur(n) # 债券指数久期表

    def _bondDict(self):
        bondDict = {"CBA05821.CS": u"1到3年利率债",
        "CBA05831.CS": u"3到5年利率债",
         "CBA05841.CS": u"5到7年利率债",
         "CBA05851.CS": u"7到10年利率债",
        "CBA02711.CS": u"1年以内信用债",
        "CBA02721.CS": u"1到3年信用债",
         "CBA02731.CS": u"3到5年信用债",
         "CBA02741.CS": u"5到7年信用债",
         "CBA02751.CS": u"7到10年信用债",
         "CBA05801.CS":u"利率债综合",
         "CBA02701.CS": u"信用债综合",
         "CBA01901.CS":u"高等级信用",
         "CBA03801.CS" :u"高收益信用"
        }

        return bondDict

资料来源：万得资讯、中金公司研究部

这里首先涉及到三个获取数据的子方法，都是直接调用万得数据即可，没有进一步的处理。不过，由于有的基金净值披露不连续可能导致计算无法进行，我们在这里加入了一个检查步骤，即数据不全者直接剔除（见下图）。在随后的dfNavRtRoll表，我们用到了pandas自带的滚动窗口加总的方法。

图表: 数据获取与剔除

     def _qstart(self, n):
        # 由于频繁用到日期回跳，这里做一个私有函数
        n += self.rol
        if not w.isconnected(): w.start() 
        return w.tdaysoffset(-n, self.start).Data[0][0].strftime("%Y%m%d")


    def getNavRt(self, n=30):
        # 基金净值，需要从start_date向前读 n + rol天，所以先调整start日
        qstart = self._qstart(n)
        _, df = w.wsd(','.join(self.codes), "NAV_adj_return1", qstart, self.end, usedf=True)
        df.index = [pd.to_datetime(x).strftime("%Y%m%d") for x in df.index]

        srsCount = df.count()
        codes = list(srsCount[srsCount == srsCount.max()].index)
        self.codes = codes
        return df.loc[:, codes]




    
    def getBondIndex(self, n=30):
        # 债券指数的表现，需要从start_date向前读 n + rol天
        qstart = self._qstart(n)
        indexCodes = self._bondDict().keys() 
        _, dfIndex = w.wsd(",".join(indexCodes), "pct_chg", qstart, self.end, usedf=True)
        dfIndex.index = [pd.to_datetime(x).strftime("%Y%m%d") for x in dfIndex.index]
        return dfIndex 

    def getBondDur(self, n=30):
        # 指数久期数据
        qstart = self._qstart(n)
        indexCodes = self._bondDict().keys()

        _, dfDur = w.wsd(",".join(indexCodes), "duration", qstart, self.end, usedf=True)
        dfDur.index = [pd.to_datetime(x).strftime("%Y%m%d") for x in dfDur.index]

        return dfDur

资料来源：万得资讯、中金公司研究部

2、初步回归选出“最优指数”：为避免多重共线性问题，我们在这里要用分布过滤法首先找到对该基金拟合效果最好的指数曲线（即“最优指数”），虽然再用其他曲线去解释“最优指数”无法刻画的部分。因此我们首先去用每个指数逐个与基金净值涨跌做一阶线性回归，得到它们的拟合优度，继而选择拟合优度最大的那一个指数。这一步与我们在固收+基金风格分解时的做法类似，实现方法如下，可以返回最优指数的代码及其对应的alpha、beta。

注意在此前还有一个_getSlice函数，用来辅助寻找对应日期所需要的净值和指数切片数据，在最后组合拼接的时候会用到。这里，为了避免因为大额申赎、信用违约等等原因造成的净值跳动，我们将波动过大的交易日做了剔除处理（阈值是当日净值波幅超过了当日波幅最大指数的3倍，见_getSlice）。

图表: 最优指数选择




    
     def _getSlice(self, code, date, n):
        i = self.dfNavRtRoll.index.get_loc(date)
        srsNav = self.dfNavRtRoll[code].iloc[i-n:i+1]
        dfIndex = self.dfIndexRoll.iloc[i-n:i+1]

        _t = srsNav.apply(pd.np.abs) <= dfIndex.applymap(pd.np.abs).max(axis=1) * 3
        ids = _t[_t].index

        if any(srsNav.isnull()):
            print u"非完整数据基金产品，跳过"
            return None, None
        else:
            return srsNav.loc[ids], dfIndex.loc[ids]

    def findLevel1(self, dfIndex, srsNav):
        # 用回归的方式寻找出对srsNav回归效果最好的指数，返回指数代码，一阶线性回归下的alpha与beta
        L = LinearRegression(fit_intercept=True)
        indexCodes = dfIndex.columns

        score,alpha,beta = 0, None, None
        for idCode in indexCodes:

            srsIndex = dfIndex[idCode]
            L.fit(srsIndex.values.reshape(-1,1), srsNav.values.reshape(-1, 1))
            if L.score(srsIndex.values.reshape(-1,1), srsNav.values.reshape(-1, 1)) >= score:
                score = L.score(srsIndex.values.reshape(-1,1), srsNav.values.reshape(-1, 1)) 
                level1,alpha,beta = idCode,L.intercept_[0


    
], L.coef_[0][0] 
        return level1,alpha,beta

资料来源：万得资讯、中金公司研究部

3、接下来要用其他指数去解释上一步留下的残差：但在此之前，我们仍然不能把剩余所有的指数都放进来，因为这样仍然会形成较大的多重共线性，造成典型的“精确错误”。于是我们先用聚类（这里用了最简单的KMeans法）把指数按照行为模式分为两组。然后找到两组里，各自与“最优指数”相关系数最低的一个指数，从而得到一共两个在最优指数之外的备选指数。下面这个函数返回两个备选指数的代码。

图表: 备选指数选择

      def _getGroup(self, dfIndex, level1):
        indexCodes = list(dfIndex.columns)
        indexCodes.remove(level1)

        km = KMeans(n_clusters=2)
        t = dfIndex / dfIndex.std()
        km.fit(t.loc[:,indexCodes]. values.transpose())
        srsRet = pd.Series(km.labels_, index=indexCodes)
        lstGroup0, lstGroup1 = list(srsRet[srsRet==0].index), list(srsRet[srsRet==1].index)

        srsCorr = dfIndex.corr()[level1]
        select0 = pd.to_numeric(srsCorr[lstGroup0]).argmin()
        select1 = pd.to_numeric(srsCorr[lstGroup1]).argmin()

        return select0, select1

资料来源：万得资讯、中金公司研究部

有了这些准备，我们就可以得到基金净值波动的回归式了。由于准备充分，最终的回归处理反而显得更简洁，如下。这里会返回一个Series变量，其中包括基金净值波动的alpha、最优指数代码及其对应的beta、两个备选指数的代码及其对应的beta值。

图表: 最终回归

     def regForSingleDaySingleCode(self, code, date, n=55):
        
        srsNav, dfIndex = self._getSlice(code, date, n)
        if srsNav is None: 
            return


    
 None                
        
        # 寻找哪个才是该基金的主要指数
        level1,alpha, beta = self.findLevel1(dfIndex, srsNav) 
        
        # 将各指数用KMeans分2组，并找到各组与Level相关系数最低的那一个指数 
        select0, select1 = self._getGroup(dfIndex, level1)
        
        # 用剩余的两种风格去回归掉残差
        
        dfX = dfIndex.loc[:, [select0, select1]]
        for sel in (select0, select1):
            dfX[sel] -= dfIndex[level1]
            
        y = srsNav - beta*dfIndex[level1] - alpha 
        
        LineFinal = LinearRegression(fit_intercept=True)
        LineFinal.fit(dfX.values, y)
        # 整理最后结果，输出alpha，beta1，beta2，beta3
        srsRet = pd.Series(index=["alpha",level1, select0, select1])
        
        srsRet["alpha"] = alpha + LineFinal.intercept_
        srsRet[level1] = beta - sum(LineFinal.coef_)
        srsRet.loc[[select0, select1]] = LineFinal.coef_
        
        return srsRet

资料来源：万得资讯、中金公司研究部

虽然为了避免多重共线，前后有分步，但最后的结果仍然是简单容易理解的，即：

净值涨跌 = alpha + beta * 最优指数波动 + beta0 * 备选指数0波动 + beta1 * 备选指数1波动

4、但回归完成并不是终点，我们要考虑长期模式和当下久期的区别——虽然债券本身比股票波动小，但实操中投资者会做一些久期的选择乃至波段交易。而在上面的回归中，我们需要一个比较长的时间窗口（我们选用的55个交易日，为斐波那契数列中的一个值），与真实情况相比存在比较明显的时滞。这里我们需要用“有效波动比”来矫正，以反映近期信息。这个方法不难理解，即用前一步产生的回归结果，把“净值涨跌-alpha”看做最优指数和备选指数的组合（称作拟合组合）。然后用“净值涨跌 – alpha” / 拟合组合涨跌得到“波动比”，再以这个比值乘以拟合组合的久期（指数有直接披露的值），从而得到最终的久期。

但什么是“有效波动”？我们的考虑是，债券指数虽然没有“锯齿型”的问题，但波动也很小，此时把其波幅作为“分母”是不太合适的——我们需要波动达到一定幅度、信用足够有用时，再计算波动比。算法上，从T日向前回溯，仅当累计波动达到阈值（我们设为0.4%）时，我们才认定“有效波动”，从而计算久期。这里的算法很简单，这里输出最终的久期，如下。

图表: 有效波动的认定和最终输出久期

     def calcDur(self, code, srsRet, date):
        # srsRet是regForSingleDaySingleCode的返回值，date为日期数据
        lstIndexes = list(srsRet.index[1:])
        durS = (self.dfDur.loc[date, lstIndexes] * srsRet[lstIndexes]).sum()
        
        absFenmu = 0.0
        for num in range(13, 60):
            # 从13开始寻找有效波动值
            i = self.dfNavRt.index.get_loc(date)
            if i >= num:
                numDaysBefore = self.dfNavRt.index[i-num]
                fenmu = (self.dfIndex.loc[numDaysBefore:date, lstIndexes] * srsRet[lstIndexes]).sum().sum() 



    
                if abs(fenmu) > 0.4:
                    ratio = (self.dfNavRt.loc[numDaysBefore:date, code].sum() - srsRet.alpha*num/float(self.rol))/ fenmu 
                    
                    break
                elif abs(fenmu) > absFenmu:
                    ratio = (self.dfNavRt.loc[numDaysBefore:date, code].sum() - srsRet.alpha*num/float(self.rol))/ fenmu 
                    absFenmu = abs(fenmu)
        else:
            print u"疑似没有有效波动，ratio强制设为期间波动最大时的水平"
            
        durRet = durS * ratio
        return durRet

资料来源：万得资讯、中金公司研究部

久期测算值怎样用？简单来说，一方面用来观察市场，即以公募纯债债基为考察对象，观察市场的选择、情绪与分歧，当然这里也要涉及主观的判断。不过，我们并不建议用样本的均值或者中位数——这样做主要的作用是掩盖算法的不稳定性——即便我们通过算法设计解决了多数问题，但我们也要承认稳定与波动同在，否则测算久期也失去了意义。我们更建议观察“扩散指标”，比如样本中有多少个基金的久期位于自己过去1年的85%以上分位数或15%以内分位数。

图表: 样本内测算久期均值和中位数示例