2018-04-16 开胃学习数学系列 - Antithetic

For this homework, we only apply the variance reduction in the common market factor z, you should not change the random number e that were drew with in the drawDefaultIndicator function, i.e., only modify the simCDO code, re-use but do not modify the CDO class. Unless explicitly mentioned, keep the simulation path the same as the base case above.

• 对于这个作业只对 common market factor z 进行variance reduction
• 不改变在drawDefaultIndicator 函数中的随机数e，e是一个随机 size 125
• z 是 z循环随机
• 即只修改simCDO 的代码，重用但不修改 CDO类。
• 除非明确提到，否则保持模拟路径与上述基本情况相同。

在z上使用Antithetic Variate 来降低Variance
我不确定这个是正确的答案，等下期更新会确认！

原函数

## price the tranches using simulation
def simCDO(cdo, rho, disc, paths) :
    
    # 1000条随机path
    zs = np.random.normal(size=[paths])
    pv = np.zeros(np.shape(cdo.a))
    pv2 = np.zeros(np.shape(cdo.d))
    for z in zs:
        thisPV, _ = cdo.drawPV(z, rho, discf)
        pv += thisPV
        pv2 += thisPV*thisPV
        
    v = pv/paths
    var = pv2/paths - v**2
    return pv/paths, np.sqrt(var/paths)

使用Antithetic Variate 函数1

## price the tranches using simulation
def simCDO02(cdo, rho, disc, paths) :   
    
    # 1000条随机path
    ys = np.random.normal(size=[paths]) 
    
    pv = np.zeros(np.shape(cdo.a))
    pv2 = np.zeros(np.shape(cdo.d))
    for z in ys:
        thisPV, _ = cdo.drawPV(z, rho, discf)
        pv += thisPV
        pv2 += thisPV*thisPV              
    for z in -ys:
        thisPV, _ = cdo.drawPV(z, rho, discf)
        pv += thisPV
        pv2 += thisPV*thisPV
        
        
    v = pv/(2*paths)
    var = pv2/(2*paths) - v**2
    return pv/(2*paths), np.sqrt(var/(paths*2))

使用Antithetic Variate 函数2 (实际是一样的)

## price the tranches using simulation
def simCDO01(cdo, rho, disc, paths) :   
    
    # 1000条随机path
    ys = np.random.normal(size=[paths]) 
    zs = concatenate((ys, -ys))
    pv = np.zeros(np.shape(cdo.a))
    pv2 = np.zeros(np.shape(cdo.d))
    for z in zs:
        thisPV, _ = cdo.drawPV(z, rho, discf)
        pv += thisPV
        pv2 += thisPV*thisPV
        
    v = pv/(2*paths)
    var = pv2/(2*paths) - v**2
    return pv/(2*paths), np.sqrt(var/(paths*2))

做表对比一下

pv_0, err_0 = simCDO(cdo, rho, discf, npath)
df = pd.DataFrame(np.array([cdo.a, cdo.d, pv_0, err_0]), 
                  index=['Attach', 'Detach', 'PV', 'MC err'])
​
fmt.displayDF(df, fmt='4g')

MC Err 有一定的下降。
我不确定这个是正确的答案，等下期更新会确认！

CDO 类

from scipy.stats import norm

class CDO(object) :
    def __init__(self, w, defProbs, recovery, a, d) :
        self.w = w/np.sum(w)
        self.p = defProbs
        self.rec = recovery
        self.rho = rho
        self.a = a
        self.d = d

    
    # 返回的是累计分布和DefProbs 的差值是否为正的 Bool Array，前者大就是False
    def drawDefaultIndicator(self, z, rho) :
        '''return a list of default indicators given common factor z, using one factor Gaussian Copula
        '''
        e = np.random.normal(size=np.shape(self.p))
        x = z*np.sqrt(self.rho) + np.sqrt(1-self.rho)*e
        return np.less(norm.cdf(x), self.p)
    

    # 累计出来的违约概率计算出 portfolio 的 损失总值
    def portfolioLoss(self, defIndicator) :
        '''compute portfolio loss given default indicators'''
        return np.sum(defIndicator*self.w*(1-self.rec))

    
    def tranchePV(self, portfLoss, discf) :
        '''compute tranche PV from portfolio loss
        Args:
            portfLoss: the total portfolio loss
            discf: discount factor
        Returns:
            tranche PVs'''
        
        sz = self.d - self.a
        return discf/sz*np.minimum(np.maximum(portfLoss - self.a, 0), sz)

    def drawPV(self, z, rho, discf) :
        ''' compute PV and portfolio Loss conditioned on a common factor z'''
        di = self.drawDefaultIndicator(z, rho)
        pfLoss = self.portfolioLoss(di)
        return self.tranchePV(pfLoss, discf), pfLoss
    
    
cdo = CDO(w, defProbs, recovery, attachements, detachements)