信用风险建模 in Python 系列 2 - 独立模型上

本文含 5120 字，44 图表截屏

建议阅读 39 分钟

0

引言

本文是「信用风险建模 in Python」系列的第二篇，其实在之前的 Cufflinks 那篇已经埋下了信用风险的伏笔，

1. 信用组合可视化
2. 信用风险 101
3. 玩具模型 - 伯努利模型

违约指的就是债务人无能力或者不愿意偿还债务。如果让你来分析哪些因素会造成违约，你会怎么写？最明显的两类因素应该是：

• 债务人类型：个人、家庭、公司和主权实体的行为会不同。
• 财务状态：盈利能力、收入支出、国家或地区税收。

此外债务水平、未来前景、商业模式、政治环境、行业竞争力、地理位置和整体规模都是可考量的因素；一些宏观经济因子比如利率和汇率水平，商品价格可会影响违约；最后可能还会考量公司的管理能力、国家领导人的政治领导力等。

能影响违约的因素太多了，没有一个模型能捕捉到现实世界所有复杂性，模型只能对复杂世界做降维处理。我们首先分析的是最简单但也能挖出价值的信用风险模型 - 伯努利（Bernoulli）模型。

该系列是理论和代码相结合，首先引入所需的 Python 包。

1

简介

对于 N 个债务人，定义第 n 个的违约事件为

，相对应的违约指示函数（default indicator）为

信用组合的损失 L 可写成每个债务人的损失 Ln 之和

其中 EAD 是exposure-at-default 的缩写，代表违约时暴露；LGD 是 loss-give-default的缩写，代表违约时损失率，而 LGD 等于 1 减去恢复率（recovery rate）。

本帖先不深入讨论 EAD 和 LGD，至少目前不是重点。当 EAD 和 LGD 之间没有相关性时，为了之后推导简化，我们可以

很明显

是个伯努利随机变量。有随机性就能求出其期望、方差、协方差和相关系数等统计指标。

我们发现协方差为零，因此上述模型没有考虑违约相关（default dependence），因此在本系列后面的文章中，我们要改进模型使其考虑违约相关。

对于损失变量 L，我们可以计算出它的各项统计指标，首先看其期望

因此一个组合的预期损失（expected loss）等于违约概率加权（default probability-weighted）的损失暴露（loss exposure）的总和。

风险价值 (value-at-risk, VaR) 是组合在持有期间内在给定置信区间内（或给定概率水平下）由于市场变动所导致的最大损失值。比如在一年内有 1% 的可能性损失一百万。VaR 衡量的是极端损失，而期望损失（expected shortfall, ES）更进一步，它衡量的是当损失超过 VaR 阈值时的平均损失，因此 ES 考虑的是更加极端的损失。

通俗来讲，VaR 用来衡量坏事，而 ES 是说如果坏事发生那么到底有多坏。下面来看严谨的数学定义：

上贴一直强调一点，要计算这些统计指标，计算出损失分布即可。通常有两种方法来计算损失分布，即损失变量 L：

1. 数值法：蒙特卡洛模拟，更通用，用来处理实际问题
2. 解析法：简化假设硬推，更快速，更直观的理解模型

后面两节分别从理论和代码的角度来阐明。

2

数值解

1.1

理论

回顾损失变量里面成分，只有违约指标是随机的，我们假设它服从伯努利分布。假设有 M 个模拟路径，N 个借贷人，那么对 n =1, 2, …, N 和 m= 1, 2, …, M, 我们需要模拟出 NM 个违约指标。

模拟方法如上式和下图所示：

剩下的操作就简单了，对于第 m 个模拟情境，计算出组合损失

将上面过程重复 M 遍得到 L(1), L(2), …, L(M)，再根据均值和方差的定义来计算它们（用 hat 表示它们是估计量而不是数学定义）

只要 M 够大，上面这些计算出来的值可以近似当成真实值。

要计算 VaR 和 ES，我们首先需要对 L(1), L(2), …, L(M) 按升序排序，即排序后 L(1)最小而 L(M) 最大。要计算 VaRq，当 q = 1%，我们计算索引 Mq，在返回对应该索引的损失值 L(Mq)。如果 Mq 不是整数，假设 M = 250, q = 99%，那么 Mq = 247.5 不是整数，有两种方法返回合理的损失值：

1. 向上取整 [Mq]。这时 [247.5] = 248，返回 L(248)
2. 找相邻索引，再做线性插值。247.5 的相邻索引为 247 和 248，那么在 L(247) 和 L(248) 线性插值。

第一种方法更保守（因为返回一个更大的损失值）；第二种方法更精确（因为返回一个内插值），不同系统商实现的方式不同，但第二种更常见。

上面公式初看很吓人，但是把 Mq = 247.5 和 [Mq] = 248 带进去就好理解了。

计算出 VaRq 后，再计算 ESq 也非常简单。

计算 ES 不需要像计算 VaR 那样分情况，因为不管什么情况，L[Mq] 总是第一个大于 VaR 的损失值。那么总共大于 VaR 的损失有 M – [Mq] 个，将它们求个平均就可以了。下面一图胜千言。

总结：在若干组情境下模拟出来的组合损失值可以当成损失分布，由损失分布套公式就能很容易的计算出 EL, UL, VaR 和 ES 了。

1.2

代码

样本组合

我们还是用之前提到的样本组合（sample portfolio），它包含 100 个不同的借贷人，有如下三个假设：

1. 组合的总规模为 1000，意味着平均每个借贷人的敞口（exposure）为 10。
2. 实际敞口是根据韦伯分布（Weibull）模拟得出，范围从小于 1 到 50。
3. 借贷人的无条件违约概率（unconditional default probability）根据卡方分布（chi-square）模拟得出，均值设为 1％。

我模拟好违约率和敞口存成两个 numpy 格式文件 expFile 和 dpFile，加载存储成变量 p 和 c，此外

• N 为借贷人数，等于 100
• M 为模拟次数，设为 1000000
• q 为百分数的列表，想看 VaR 和 ES 在不同置信度下的表现
• q_style, money_fmt 和 number_fmt 只是为了打印出来的 DataFrame 好看些

c = np.load(expFile)  p = np.load(dpFile)N = len(c)M = 1000000q = [0.95, 0.97, 0.99, 0.995, 0.999, 0.9997, 0.9999]q_style = [str(i*100) +'%' for i in q]money_fmt = '${0:,.2f}'number_fmt = '{0:,.2f}'

编写三个函数，分别计算损失分布（binomial_LD），计算风险指标（risk_measure）和整体模拟（binomial_simulation）。代码很简单，按照 1.1 的公式和逻辑就能轻易实现。

运行来生成 EL, UL, VaR 和 ES。

EL, UL, VaR, ES = binomial_simulation( N, M, p, c, q )

打印出在不同置信度下的 VaR 和 ES 值，都是递增的。

打印出该组合的 EL 和 UL 值，它们对于不同的置信度的都是一样的。我们可以看出极端损失（VaR, ES）要比 UL 大，因此损失波动率并不是一个可能捕捉注组合风险的好指标。

再看看损失分布，这里只模拟了 100 条，主要只想看看大概的样子。我发现用 1000000 条画 cufflinks 的图 Notebook 会卡死。

LD = binomial_LD( N, M=100, p, c, q=0.99 )df_LD = pd.DataFrame( LD[::-1] )df_LD.iplot( kind='bar',             histnorm='percent',             title='Loss Exposure',             xTitle='USD',             yTitle='Relative Frequency',             color='rgb(220,38,36)',             theme='ggplot' )

最后再可视化在不同置信度下的 VaR 和 ES，可以很清楚的看出 ES 总比 VaR 大。

df1_IB.iplot( kind='scatter',              mode='lines+markers',              size=10,               title='Tail Loss Distribution',               xTitle='Quantile',               yTitle='USD',              color=color,               theme='ggplot')

2

解析法

2.1

理论

推导解析解时需要做进一步模型假设，即假设所有借贷人的违约概率和损失暴露都有相同的 p 和 c。这种假设只在借贷人很多的“大型风险分散”组合才合理，但我们的目标并不是复现在“真实组合”上用蒙特卡洛计算出来的值，而是去深入了解模型背后的关键假设。

其中

是整个组合违约的个数（不一定是整数）。

由于每个借贷人的损失暴露都为 c，因此损失分布可以在一组 {0c, 1c, 2c, …, Nc} 离散点上构建，组合最小损失为 0，最大损失为 Nc。对于离散的损失值，对应的概率质量函数（probability mass function, PMF）为

对于 k = 0, 1, 2, …, N。

现在随机变量是

，可能的取值是 0, 1, 2, …, N。首先它的期望和方差：

由上面结果可知

服从二项分布（binomial distribution），它的 PMF 和累积分布函数（cumulative distribution function, CDF）为

由于

，因此可以得到

2.2

代码

编写一个函数，计算二项分布的 PMF, CDF, VaR, ES 以及组合违约总个数 DN，代码也不难，

• PMF 和 CDF 直接用 scipy.stats 里面的函数
• DN 跟 quantile 有关，把 CDF 当做自变量，个数当做变量，线性插出就行
• VaR 就是 c 乘上 DN ，因为我们把损失整数离散化了
• ES 稍微麻烦点，但对于每个 quantile qi，将 qi 到 1 分 1000 个点，然后求出均值

根据之前的假设，p 和 c 对于所有借贷人都是一样的，再根据真实样本组合的平均违约率为 1%，总敞口为 1000 美元，有 100 个借贷人，那么平均每个敞口为 10 美元，因此将 p 和 c 设为 0.01 和 10。

(p, c) = (0.01, 10)PMF, CDF, VaR, ES, D_N = binomial_analytical( N, p, c, q, True )

打印出在不同置信度下的违约总个数 DN，VaR 和 ES 值，都是递增的。但这值明显比蒙特卡洛模拟出来的值要小。但这么比也不公平，因为结果基于两个不同的组合，一个是每个借贷人一样的 p 和 c（同质组合，风险当然会小）；一个是每个借贷人不同的 p 和 c（异质组合，风险当然会大）。

打印出该组合的 EL 和 UL 值，记注它们，和下面用蒙特卡洛在同质组合上模拟的结果对比。

为了苹果比苹果，我们用蒙特卡洛在相同的同质组合上模拟，再和之前的解析解结果对比。

M = 1000000EL, UL, VaR, ES = binomial_simulation( N, M, p*np.ones((N)), c*np.ones((N)), q )

打印出在不同置信度下的 VaR 和 ES 值，都是递增的，和解析解比较接近，和异质组合下的蒙特卡洛模拟结果相差很远。

打印出 EL 和 UL 值，解析解几乎一样。

把上面结果全部整理到一张表里，我们可以得出很清楚的看出，异质组合的风险要比同质组合的风险大很多（高 VaR 和 ES 值）。

那么我们费那么多功夫研究同质组合有什么用呢？请看下章（前方高能，数学不好的请忽略）。

3

硬核讨论

同质组合的特点就是所有借贷人的违约概率一样，记为 p，那么我们来看看当 N - 借贷人的个数 - 趋近无穷大时模型的表现。具体而言，我们想看看 DN /N，组合的违约个数比，是否收敛于 p。收敛的定义分为两种：

• 概率收敛（convergence in probability）
• 几乎收敛（convergence almost surely）

证明概率收敛需要切比雪夫不等式（Chebyshev’s inequality），推导如下

证明几乎收敛需要马可尔夫不等式（Markov’s inequality）和二项分布四阶中心矩（4th central moment），推导如下

现在我们证出了两个不等式，

因此，当借贷人个数非常多时，组合的违约个数比 DN /N 可以认为是 p，那么你可以用一个参数 p 来模拟该组合的损失分布，这个时候模型的表现会非常稳健。

4

总结

诚然，虽然二项分布的模型过于简单，但深入研究它我们至少知道为什么不好，原因有二：

1. 违约独立性的假设，这和实际情况不符。
2. 二项式分布的极限收敛于正态分布，而正态分布通常不能捕捉到尾部损失。

从保守派风险管理者看来，其他所有条件都一样，我们希望将更多的概率分配给极端事件。而这在信用风险尤其重要，因为我们可以完全将注意力集中在极端事件。二项模型在现实中的组合表现通常不会太好，但是，只要组合里人数够多，模型完全由单参数来描述。

因此，明确了模型缺点之后，我们就要改进它，这将激发该系列之后考虑的许多模型，比如混合型模型，阈值型模型，它们都以各自的方法考虑了违约相关而因此可以捕捉到尾部损失。虽然模型有很各种，但本质都是为了克服独立模型的缺点：

违约独立

带着这个主思路继续学习，不会乱。

Stay Tuned！