CVA and vulnerable options in stochastic volatility models

TL;DR

This paper develops a general formula to evaluate the Credit Value Adjustment (CVA) for vulnerable options within stochastic volatility models, accounting for correlation with default risk, and validates it through numerical analysis.

Contribution

It introduces a unified representation formula for CVA in stochastic volatility models with default correlation, applicable to models like SABR, Hull-White, and Heston.

Findings

The formula explicitly captures the impact of correlation on CVA.

Numerical results show good agreement with Monte Carlo simulations.

The approach enhances accuracy in CVA evaluation for models with default risk.

Abstract

In this work we want to provide a general principle to evaluate the CVA (Credit Value Adjustment) for a vulnerable option, that is an option subject to some default event, concerning the solvability of the issuer. CVA is needed to evaluate correctly the contract and it is particularly important in presence of WWR (Wrong Way Risk), when a credit deterioration determines an increase of the claim's price. In particular, we are interested in evaluating the CVA in stochastic volatility models for the underlying's price (which often fit quite well the market's prices) when admitting correlation with the default event. By cunningly using Ito's calculus, we provide a general representation formula applicable to some popular models such as SABR, Hull \& White and Heston, which explicitly shows the correction in CVA due to the processes correlation. Later, we specialize this formula and construct its approximation for the three selected models. Lastly, we run a numerical study to test the formula's accuracy, comparing our results with Monte Carlo simulations.

Tables

Table 1. Table 1: Parameter sets for the Vasicek default intensity.

	${\lambda }_0$	$q$	$\mu$	$\sigma$
Set 1	0.09	0.3	0.4	0.1
Set 2	0.1	0.18	0.1	0.015

Table 2. Table 2: Parameter sets for the CIR default intensity.

	$\lambda$	$q$	$\mu$	$\sigma$
Set 1	0.03	0.02	0.161	0.08
Set 2	0.05	0.09	0.2	0.1
Set 3	0.01	0.8	0.02	0.2
Set 4	0.03	0.5	0.05	0.5

Equations

$$c^d(t,T)=\mathbf{E}[B(t,T)f(X_T)1_{\{\tau >T\}}+B(t,\tau )Z_{\tau }{1}_{\{t<\tau \le T\}}|{\mathcal{G}}_t],t\in [0,T],$$

$$\mathbf{E}\Big{[}\mathbf{1}_{\{\tau>t\}}Y|{\mathcal{G}}_{t}\Big{]}=P(\tau>t|{\mathcal{G}}_{t})\frac{\mathbf{E}\Big{[}\mathbf{1}_{\{\tau>t\}}Y|{\mathcal{F}}_{t}\Big{]}}{P(\tau>t|{\mathcal{F}}_{t})}.$$

$$\mathbf{E}[B(t,T)f(X_T)1_{\{\tau >T\}}|{\mathcal{G}}_t]$$

$$\mathbf{E}[B(t,\tau )Z_{\tau }{1}_{\{t<\tau \le T\}}|{\mathcal{G}}_t]$$

$$F_t=P(\tau \le t|{\mathcal{F}}_t),\forall t\ge 0,$$

$${\Gamma }_t:=-\ln (1-F_t)\Rightarrow F_t=1-{\mathrm{e}}^{-{\Gamma }_t}\forall t>0,{\Gamma }_0=0.$$

$$\mathbf{E}[B(t,T)f(X_T)1_{\{\tau >T\}}|{\mathcal{G}}_t]=1_{\{\tau >t\}}\mathbf{E}[B(t,T)f(X_T)1_{\{\tau >T\}}|{\mathcal{F}}_t]{\mathrm{e}}^{{\Gamma }_t}$$

$$=$$

$$\mathbf{E}[M_{\tau}\mathbf{1}_{\{t<\tau\leq T\}}|{\mathcal{F}}_{t}]=\mathbf{E}\Big{[}\int_{t}^{T}M_{u-}dF_{u}+\sum_{t<u\leq T}\Delta M_{u}\Delta H_{u}|{\mathcal{F}}_{t}\Big{]},$$

$$\mathbf{E}[M_{\tau }{1}_{\{t<\tau \le T\}}|{\mathcal{F}}_t]=\mathbf{E}[{\int }_t^T{M}_{u-}d{F}_u|{\mathcal{F}}_t].$$

$$\mathbf{E}[M_{\tau }{1}_{\{t<\tau \le T\}}|{\mathcal{F}}_t]=\mathbf{E}[{\int }_t^T{M}_{u}d{F}_u|{\mathcal{F}}_t].$$

$$\mathbf{E}[B(t,\tau )Z_{\tau }{1}_{\{t<\tau \le T\}}|{\mathcal{G}}_t]=1_{\{\tau >t\}}\mathbf{E}[B(t,\tau )Z_{\tau }{1}_{\{t<\tau \le T\}}|{\mathcal{F}}_t]{\mathrm{e}}^{{\Gamma }_t}$$

$$=$$

$$\!c^{d}(t,T)\!=\!\mathbf{1}_{\{\tau>t\}}\!\left[\mathbf{E}\Big{(}\mathrm{e}^{-\!\int_{t}^{T}(r_{s}\!+\lambda_{s})ds}f(X_{T})\!+\!\!\int_{t}^{T}\!\!\!\!Z_{s}\lambda_{s}\mathrm{e}^{-\!\int_{t}^{s}(r_{u}\!+\lambda_{u})du}ds|{\mathcal{F}}_{t}\Big{)}\right]\!,$$

$$\displaystyle\mathbf{E}[\int_{t}^{T}Z_{s}\lambda_{s}\mathrm{e}^{-\int_{t}^{s}(r_{u}+\lambda_{u})du}ds|{\mathcal{F}}_{t}]=R\mathbf{E}\Big{[}\int_{t}^{T}c(s,T)\lambda_{s}\mathrm{e}^{-\int_{t}^{s}(r_{u}+\lambda_{u})du}ds|{\mathcal{F}}_{t}\Big{]}$$

$$\displaystyle c^{d}(t,T)=\mathbf{1}_{\{\tau>t\}}\Big{[}R\mathbf{E}\Big{(}\mathrm{e}^{-\int_{t}^{T}r_{u}du}f(X_{T})|{\mathcal{F}}_{t}\Big{)}+(1-R)\mathbf{E}\Big{(}\mathrm{e}^{-\int_{t}^{T}(r_{u}+\lambda_{u})du}f(X_{T})|{\mathcal{F}}_{t}\Big{)}\Big{]},$$

$$CVA(t):=$$

$$=$$

$$\mathbf{E}\Big{[}\mathrm{e}^{-\int_{t}^{T}(r_{s}+\lambda_{s})ds}f(X_{T})|{\mathcal{F}}_{t}\Big{]}=\mathbf{E}\Big{[}\mathrm{e}^{-\int_{t}^{T}r_{s}ds}f(X_{T})|{\mathcal{F}}_{t}\Big{]}\mathbf{E}\Big{[}\mathrm{e}^{-\int_{t}^{T}\lambda_{s}ds}|{\mathcal{F}}_{t}\Big{]}.$$

$$d{X}_t$$

$$d{Y}_t$$

$$B_t^1$$

$$\text{SABR model: }\gamma \in (0,1),a(y)=y,b(t,y)=0,c(t,y)=cy,c>0$$

$$\text{Hull \& White model: }\gamma =1,a(y)=y,b(t,y)=b(t)y,c(t,y)=c(t)y,b,c:[0,T]\to {\mathbb{R}}^{+}\text{bounded}$$

$$\text{Heston model: }\gamma =1,a(y)=\sqrt{y},b(t,y)=k(\theta -y),c(t,y)=c\sqrt{y},k,\theta >0,c^2<2k\theta .$$

$$\mathbf{E}_{t}\left[(1-N^{t}_{T})U(T,X_{T},Y_{T})\right]\mathbf{1}_{\{\tau>t\}},\quad\textrm{where}\quad N^{t}_{s}:=\mathbf{E}_{s}\Big{(}\mathrm{e}^{-\int_{t}^{T}\lambda_{u}du}\Big{)}$$

$$dU(s,X_s,Y_s))$$

$$+a(Y_s){\mathrm{e}}^{-(1-\gamma )X_s}{\partial }_{x}U(s,X_s,Y_s)d{B}_s^1+c(s,Y_s){\partial }_{y}U(s,X_s,Y_s)d{B}_s^2$$

$$={\partial }_{x}U(s,X_s,Y_s)d{M}_s^X+{\partial }_{y}U(s,X_s,Y_s)d{M}_s^Y,$$

$$\begin{gathered} \mathcal{L}={\partial }_s+\frac{a^2(y){\mathrm{e}}^{-2(1-\gamma )x}}{2}({\partial }_{xx}^2-{\partial }_x)+\frac{c^2(s,y)}{2}{\partial }_{yy}^2 \\ +\eta a(y)c(s,y){\mathrm{e}}^{-(1-\gamma )x}{\partial }_{xy}^2+b(s,y){\partial }_y \end{gathered}$$

$$\{\begin{array}{ll}& \mathcal{L}U(s,x,y)=0\\ & U(T,x,y)=f(x).\end{array}$$

$$\boxed{\begin{aligned} CVA(t)=&\mathbf{E}_{t}\left[(1-N^{t}_{T})U(T,X_{T},Y_{T})\right]=(1\!-\!N^{t}_{t})U(t,x,y)\\ -&\mathbf{E}_{t}\Big{[}\!\int_{t}^{T}\!\!\!\!\partial_{x}U(s,X_{s},Y_{s})d\langle N^{t}\!,X\rangle_{s}\Big{]}\!-\mathbf{E}_{t}\Big{[}\!\int_{t}^{T}\!\!\!\partial_{y}U(s,X_{s},Y_{s})d\langle N^{t}\!,Y\rangle_{s}\Big{]},\end{aligned}}$$

$$\text{CVA}(t)={\mathbf{E}}_t\left[(1-N_T^t)u(T,X_T,v_T^2)\right]1_{\{\tau >t\}},$$

$$V_{t}^{2}=\frac{1}{T-t}\mathbf{E}_{t}\Big{(}\int_{t}^{T}\sigma^{2}ds\Big{)}=\frac{1}{T-t}\mathbf{E}_{t}\Big{(}\int_{t}^{T}a^{2}(Y_{s})ds\Big{)}.$$

Full text

Key words: Credit Value Adjustment, Vulnerable Options, Stochastic volatility model, Intensity approach

JEL Classification: E43, G12, G13.

Mathematics Subject Classification (2010): 91G60, 91G20, 60J75.

1 Introduction

Defaultable claims are derivatives that are subject to some default event, which concerns the solvability of the counterparty before the final settlement of the transaction. This is the so called ”Counterparty Credit Risk” (CCR), an immediate consequence of it being that the product’s price needs an adjustment to include in its quotation the possibility of default. This adjustment, which results in a price reduction, is commonly known as **Credit Value Adjustment **(CVA) and it was first introduced in a paper by Zhu and Pykhtin ([29]). The last financial crisis (2007-2008) has greatly increased the monitoring and pricing of CCR on OTC-markets products and many researchers and practitioners tried to develop a general framework for a better assessment of the CVA evaluation to compensate a derivatives holder for taking CCR. Indeed, along the years, other value adjustments have been additionally considered leading to the acronym (X)VA. An updated overview of the recent research directions under investigation is presented in [15].

Many works concentrate on CVA evaluation for swaps while others focus on European options. In this case , when the risk relates only to the issuer, these contracts are called vulnerable options and one can find a vast literature under this denomination (see e.g. [12], [8],[9], [4] and references therein, the last three focusing more on CVA).

Typically the default event is characterized by means of a random time, representing the time of default. At the time of default there might be a total loss for the investor or a partial recovery of the investment’s current value might be possible.

The difficulty in the evaluation of the CVA is twofold. First of all the default time might be not completely measurable with respect to the information generated by the market prices, since it reflects also other exogenous factors, secondly even under full knowledge of this random time, the derivative’s evaluation will call for the joint distribution of the random time and the price processes, usually very difficult to know.

In this framework we consider the CVA evaluation of a vulnerable European option under a variety of stochastic volatility models. The importance of this choice lies in the fact that stochastic volatility models better fit the market, reproducing the smiles and skews of the implied volatility. Pricing of vulnerable options under the Heston model was already discussed in some papers ( [21] and [28]) under the structural default framework (see [20]).

Alternatively, one might use the so called intensity approach (introduced in [24] to price defaultable bonds) in order to characterize the distribution of the default time conditionally to the information generated by the market prices. In a stochastic volatility setting a first result was obtained in [27] in discrete time, assuming a GARCH model for the underlying coupled with default intensity, but the literature is scarce in continuous time.

In the present paper, we describe the joint dynamics of the asset prices, of the default time and of the other stochastic factors as a Markovian system whose components may exhibit correlation. We remark that the presence of correlation among processes is crucial, indeed under independence between the intensity and the price/volatility processes, the evaluation formula reduces to two separate evaluations: the classical default-free derivative price multiplied by a factor (similar to a bond price paying $1$ at maturity) induced by the intensity process.

When in presence of correlation, computations do not come out as easily and we employ a keen integration by parts formula, inspired by a technique developed by E. Alos ([1]), to enucleate the contribution due to the correlations. The evaluation formula gets split into three terms: the first giving the zero correlations CVA value, i.e. the value for the independent case, the second and the third terms coming from the correlation of the intensity process respectively with the asset’s price and the stochastic volatility. We call this the “first order representation formula” since it involves the first order derivatives of the default free price. This expression points out the contribution of the correlation, but it does not identify it explicitly. Were it needed, one can enlighten the role played by each process by applying again Itô’s formula to the default-free price, getting a so called “second order representation formula” . The computations become more involved, but this second expression captures more accurately the behavior of the correlations when these are parametrized.

More in detail, we model the intensity either by a Vasicek or a CIR process each coupled with the SABR, Hull & White and Heston stochastic volatility models. By consequence, the representation formula is rewritten explicitly in terms of the correlation coefficients among the drivers of the SDE’s describing the model, in particular between the asset’s price and intensity ($\rho$) and between the intensity and the volatility ($\nu$). The first and second order formulas we obtain correspond to a first or second order expansion of the CVA with respect to the correlation parameters and this approach might be extended to any order at cost of a larger computational effort. These are approximated by “ad hoc” techniques for each case and the numerical accuracy and computational speed we achieve makes this technique a valid alternative to Monte Carlo methods.

To clarify the role of $\rho$, let us consider an investor $A$ buying from a bank $B$ a call option written on the asset of a bank $C$ competing with $B$ on a large market. If $B$ is going to default, the intensity $\lambda_{t}$ is going to have a large value and $\rho>0$ will reflect the propension of the asset $C$ to increase its value accordingly, due to gaining of market positions. Therefore $A$ will have in the portfolio a call contract which is deep in the money but with a counterparty proximal to default, describing a wrong way risk (WWR) situation.

Finally, we remark that our method can be straightforwardly extended to include a stochastic interest rate and its numerical accuracy and speed makes it interesting to apply it to more general XVA evaluation in a future research.

In the next two sections, we first introduce the theoretical framework for CVA evaluation and then we provide our general representation formulas for a Markovian setting. In section 4 we specialize the first order formula for the three stochstic volatility models associated with the two intensity models and we suggest an appropriate approximation technique for each of them. In section 5 we carry out the same discussion for the second order formula only for the SABR and Heston models coupled with a CIR intensity model. These are the most interesting cases and we decided to limit ourselves to those for the sake of exposition. The numerical discussion of our method is to be found in the last section.

2 CVA Evaluation of Defaultable European Claims

Let $[0,T]$ be a finite time interval and $(\Omega,\mathcal{F},P)$ a complete probability space endowed with a filtration $\{{\mathcal{F}}_{t}\}_{t\in[0,T]}$, augmented with the $P-$null sets and made right continuous. We also assume that all the processes have a cádlág version.

The market model will be described by the interest rate process $r_{t}$ and by a process $X_{t}$ representing an asset log-price whose dynamics will be specified later. The asset price may depend also on other stochastic factors. We shall assume

• •

that the filtration $\{{\mathcal{F}}_{t}\}_{t\in[0,T]}$ is rich enough to support all the aforementioned processes;

• •

to be in absence of arbitrage;

• •

that the given probability $P$ is a risk neutral measure, already selected by some criterion.

So we denote by $S_{t}=\mathrm{e}^{X_{t}}$ the asset price and by $B(t,s)=\mathrm{e}^{-\int_{t}^{s}r_{u}du}$ the forward discounting factor. In this market, a defaultable European contingent claim paying $f(X_{T})$ at maturity is traded, where $f$ is some function that will be specified later. We denote by $\tau$ (not necessarily a stopping time w.r.t. the filtration ${\mathcal{F}}_{t}$) the default time of the contingent claim and by $Z_{t}$ an ${\mathcal{F}}_{t}-$measurable bounded recovery process.

To properly evaluate this type of derivative we need to include the information generated by the default time. We denote by ${\mathcal{G}}_{t}$ the progressively enlarged filtration, that makes $\tau$ a ${\mathcal{G}}_{t}-$stopping time, that is ${\mathcal{G}}_{t}={\mathcal{F}}_{t}\lor\sigma(\{\tau\leq t\})$. Hence, denoting by $H_{t}=\mathbf{1}_{\{\tau\leq t\}}$ and by ${\mathcal{H}}_{t}$ its natural filtration, we choose ${\mathcal{G}}_{t}={\mathcal{F}}_{t}\lor{\mathcal{H}}_{t}$.

We make the fundamental assumption, known as the H-hypothesis (see e.g. [14] and [13] and the references therein), that

(H)        Every ${\mathcal{F}}_{t}-$martingale remains a ${\mathcal{G}}_{t}-$martingale.

Under this assumption $B(t,s)S_{s},\,s\geq t$ remains a ${\mathcal{G}}_{s}-$martingale under the unique extension of the risk neutral probability to the filtration ${\mathcal{G}}_{s}$. (To keep notation light, we do not indicate the probability we use for the expectations, assuming to be the corresponding one to the filtration in use).

In this setting, the price a defaultable claim, with final value $f(X_{T})$, default time $\tau$ and recovery process $\{Z_{t}\}_{t}$ is given by

[TABLE]

while the corresponding default free value is $c(t,T)=\mathbf{E}[B(t,T)f(X_{T})|{\mathcal{F}}_{t}].$

In many situations, investors do not know the default time and they may observe only whether it happened or not. The actual observable quantity is the asset price, therefore it is interesting to write the pricing formula (1) in terms of ${\mathcal{F}}_{t}$, rather than in terms of ${\mathcal{G}}_{t}$. For that we have the following key Lemma, see [7] or [5].

Lemma 2.1

For any integrable $\mathcal{G}-$measurable r.v. $Y$, the following equality holds

[TABLE]

Applying this lemma to the first and the second term of (1) and recalling that $1-H_{t}=\mathbf{1}_{\{\tau>t\}}$ is ${\mathcal{G}}_{t}-$measurable, we obtain

[TABLE]

which may be made more explicit by following the hazard process approach.

We denote the conditional distribution of the default time $\tau$ given ${\mathcal{F}}_{t}$ by

[TABLE]

whence, for $u\geq t$, $P(\tau\leq u|{\mathcal{F}}_{t})=\mathbf{E}(P(\tau\leq u|{\mathcal{F}}_{u})|{\mathcal{F}}_{t})=\mathbf{E}(F_{u}|{\mathcal{F}}_{t})$. We also assume that $F_{t}(\omega)<1$ for all $t>0$ to well define the so called hazard process

[TABLE]

With this notation, we rewrite (3) as

[TABLE]

Assuming that $B(t,\cdot)Z.$ is a bounded $\mathcal{F}-$ martingale (which is usually the case), we can treat the expectation in (4), applying an extension of Proposition 5.1.1 of [6], as developed in [3].

Proposition 2.1

Let $M$ be a bounded, ${\mathcal{F}}$-martingale. Then

1. (i)

for any $t\leq T$

[TABLE]

where by $\Delta X_{s}$ we denoted $X_{s}-X_{s-}$, for any process $X$. 2. (ii)

If moreover $M$ and $H$ do not have simultaneous jumps, then

[TABLE] 3. (iii)

*If either $M$ is *${\mathcal{F}}-$predictable or $M$ and $H$ do not have simultaneous jumps and $F$ is predictable, then

[TABLE]

Remark 2.1

The processes $H$ and $B(t,\cdot)Z.$ do not have simultaneous jumps, since $B(t,\cdot)Z.$ is cádlág and $\tau$ is not measurable with respect to $\{{\mathcal{F}}_{t}\}$, so it is not a ${\mathcal{F}}-$stopping time. Jump times of cádlág process are stopping times (see [26]), hence the two processes cannot jump simultaneously.

Keeping the previous remark in mind and assuming that the hazard process is differentiable with derivative $\lambda_{t}$, called the intensity process, $\Gamma_{t}=\int_{0}^{t}\lambda_{u}du,$ (4) becomes

[TABLE]

since $\Gamma$ is continuous and it is not going to charge any jump. Hence the pricing formula (1) becomes

[TABLE]

recovering Lando’s formulas (3.1) and (3.3) in [24].

This formula can be specialized even further assuming fractional recovery ([11]), $Z_{t}=Rc(t,T)$ for some $0\leq R<1$ and using the Optional Projection Theorem (see [25] Theorem 4.16)

[TABLE]

Putting the two pieces together, we finally obtain

[TABLE]

which was used also by Fard in [12] with $f(x)=(\mathrm{e}^{x}-K)^{+}$ and that can be interpreted as a convex combination of the default free price and the price with default.

As a consequence we have an expression also for the unilateral CVA, defined as difference between the default free price and the adjusted price

[TABLE]

It is immediate to notice that if $\lambda_{t}$ is independent from $(X_{t},r_{t})$, then

[TABLE]

Of course, the computability of the expectations will depend heavily on the modeling choices one makes for $\lambda$.

3 A representation formula for the CVA in stochastic volatility models

In this section we consider a family of stochastic volatility models and, for the sake of simplicity, from now on we assume zero fractional recovery ($R=0$) and risk free spot rate $r=0$. These are not restrictive assumptions, since the following discussion can be easily extended to the case of $0<R<1$ and $r$ a deterministic function of time. With the same techniques, increasing the dimensionality of the problem also a stochastic interest rate might be considered. Our market model is therefore given by the asset log-price $X_{t}$ and a stochastic factor, $Y_{t}$, determining the volatility process

[TABLE]

where $B_{t}$ and $Z_{t}$ are independent Brownian motions, $0\leq\eta^{2}<1$, $b,c:[0,T]\times\mathbb{R}^{+}\longrightarrow\mathbb{R}^{+}$ and $a:\mathbb{R}^{+}\longrightarrow\mathbb{R}^{+}$ are deterministic functions so that (14)-(15) admits a unique strong solution.

We remark that, because of the deterministic coefficients, the pair $(X_{t},Y_{t})$ as well as the stochastic factor $Y_{t}$ are Markovian processes, which implies that the claim’s price will be a deterministic function of these processes. We remark that, with the right choice of functions, several popular models are included in the above formulation. As a matter of fact

[TABLE]

From now on, we take the shorter notation $E_{t}$ to denote conditional expectation w.r.t. the filtration generated by the processes under consideration and let us consider the CVA problem given by formula (12), for $R=0$. In other words we have to evaluate the risk-neutral expectation

[TABLE]

for $t\leq s\leq T$ and by $U(s,X_{s},Y_{s})$ we denoted the default free price of the European claim whose payoff is $U(T,X_{T},Y_{T})=f(X_{T})$. The default free price $U$ is, by no arbitrage, a martingale and consequently applying Itô’s formula we have

[TABLE]

where $M^{X}$ and $M^{Y}$ are the respective martingale parts of $X$ and $Y$ and

[TABLE]

and $U$ verifies the PDE

[TABLE]

By integration by parts and making use of (18) we get the following basic representation formula of the difference between the ‘classical price’ and the adjusted price on the event $\{\tau>t\}$

[TABLE]

that is the uncorrelated term plus two extra terms that come from the correlation between the asset and the default and the volatility and the default.

Remark 3.1

We may view evaluation (16) from a different perspective:

[TABLE]

where $v^{2}_{t}:=V^{2}_{t}(T-t)$, $\sigma_{t}=a(Y_{t})$ is the stochastic volatility process and $V_{t}^{2}$ denotes the adapted averaged variance process (or zero-strike variance swap)

[TABLE]

The process $v^{2}$ is called the variance swap and it more easily estimated from the market data rather than non tradeable stochastic factor $Y$. These two proceses are strictly connected, and we assume

[TABLE]

and we may equivalently write the evaluations by means of $v^{2}$. Indeed, by (H1) $(X,v^{2})$ is still a Markovian pair and the function $u$ verifies an equivalent partial differential equation

[TABLE]

with

[TABLE]

arriving at

[TABLE]

Hypotheses (H1) is not restrictive, as a matter of fact it is verified by all the models we are interested in

[TABLE]

Remark 3.2

Since $\sigma_{t}=a(Y_{t})$, when $a$ is also an invertible and differentiable function, we have that $u(t,X_{t},v_{t}^{2})=U(t,X_{t},Y_{t})=U(t,X_{t},a^{-1}(\sigma_{t}))$, consequently

[TABLE]

whence we may give an approximation of the above expression (21), by freezing the integrands at the inital time $t$

[TABLE]

since

[TABLE]

When specializing to the three models we are interested in, assuming some specific dynamics for $\lambda$, the quantities appearing in (19) or (21) might be explicitly computed in terms of the correlation parameters between the intensity and the other two processes, which is what we are going to do in the next sections.

The representation formula (19) is quite crucial, because it allows to exploit an explicit expression of the martingale $N^{t}$, whenever possible. This happens, for instance when considering affine models for the dynamics of $\lambda$, in this case we may resort to the bond pricing theory and we may obtain that

[TABLE]

for some deterministic differentiable functions $\varphi$ and $\psi$ of time to maturity, which implies that

[TABLE]

where $M^{\lambda}$ denotes the martingale part of the process $\lambda$. Representation (19) hence becomes

[TABLE]

The quadratic covariations will express the correlation between the processes of the model. We remark, though, that correlations are still present between $\partial_{x}U(s,X_{s},Y_{s})$, $\partial_{y}U(s,X_{s},Y_{s})$ and the processes defined by

[TABLE]

so one might think of refining the above formula, by applying Itô’s formula and integration by parts once again.

By differentiating (18) with respect to $x$ and $y$, we may conclude that $\partial_{x}U$ and $\partial U_{y}$ verify the following system of PDE’s

[TABLE]

where

[TABLE]

Consequently, applying Itô calculus we have

[TABLE]

Similarly for the other term we obtain

[TABLE]

Plugging these expressions into (26), we obtain a second order representation formula for CVA

[TABLE]

In the next sections we will consider either the Vasicek or the CIR model for the intensity process and, after specializing the above representations in those cases for a call option, we are going to suggest appropriate approximation formulas, that are going to be discussed numerically in the last section.

We remark that, even though the Vasicek model has already been used in the literature for the intensity process (see [12]), it has the drawback that it can become negative with positive probability, making the process $\Gamma$ a bad candidate to represent a probability. This means that this model can be employed only if the estimated parameter values imply that the probability $\lambda$ may become negative is negligible.

4 The first order approximation

In this section we analyze the CVA first order representation formula (19) for a defaultable call option with strike price given by $K=\mathrm{e}^{\kappa},\,\kappa\in\mathbb{R}$, for our three models when either a Vasicek or a CIR dynamics is chosen for the intensity.

4.1 The Vasicek intensity model

We assume that the intensity process verifies

[TABLE]

with constants $q,\mu,\sigma>0$ and

[TABLE]

where ${W}_{s}=(Z_{s},B_{s},U_{s})$ are independent Brownian motions, $(\nu,\eta,\rho)$ are such that

[TABLE]

In this way we set up the correlations to be

[TABLE]

and to simplify notation we set

[TABLE]

Consequently, for $t\leq s\leq T$, we have

[TABLE]

and $N^{t}_{t}$ is the price of a zero coupon bond with spot rate $\lambda$.

4.1.1 SABR - V

Chosen $0<\gamma<1$, we have for $s>t$

[TABLE]

Therefore

[TABLE]

so that (19) becomes

[TABLE]

which might be easily approximated by freezing at the initial time all the terms containing $X$

[TABLE]

We remark that this choice forces the approximation to stay linear in $\rho$ and we will see later that this is well supported by the numerical simulations. Were it not the case, one might think of refining the calculations by building an approximation for the second order representation formula (28).

Therefore the problem is reduced to having an expression for the default free evaluation $U$. The power of the SABR model lies in the fact that an explicit formula for the call price can be derived by substituting in the Black’s formula the SABR implied volatility, $\sigma_{B}(t,x,y,\kappa)$ (see [16]) and in our case ($r=0$) we obtain

[TABLE]

The implied volatility may be efficiently approximated by truncating the formula (2.17a) in [16] up to the second order. In our setting this approximation is

[TABLE]

where

[TABLE]

We notice that actually $U(t,x,y)=U(t,x,\sigma_{B}(t,x,y))$, whence we may approximate $U$ by

[TABLE]

and we obtain

[TABLE]

We omit here the lengthy calculations of the two derivatives $\partial_{x}{\tilde{\sigma}}_{B}$ and $\partial_{y}{\tilde{\sigma}}_{B}$, which we performed with the help of Matlab. Inserting the above expressions in (33), we obtain

[TABLE]

Arrived at this stage, we might obtain a very handy approximation by simply freezing also all the other processes at the initial time. As well, we might think of assuming $\nu=0$ (often done by practitioners) so that $\mathbf{E}_{t}(N^{t}_{s}Y_{s})=\mathbf{E}_{t}(N^{t}_{s})\mathbf{E}_{t}(Y_{s})=yN_{t}^{t}$. Both choices are rather crude, so we decided to push computations a little further to get a more accurate evaluation of this expectation. Recalling the explicit expressions of $Y_{s}$ and $\lambda_{s}$

[TABLE]

the definition of $\psi(T-s)$ and $\varphi(T-s)$ and that $B^{2}_{u}=B_{s}$ and $B^{3}_{u}=\nu B_{u}+\alpha Z_{u}+\beta U_{u}$ and gathering the alike terms together, after some calculations we arrive at the following expression

[TABLE]

where we used that the independence of the Brownian motions $B,Z,U$ implies that the exponent is Gaussian with density ${\mathcal{N}}(0,f_{2}^{2}(t,s))$, with

[TABLE]

4.1.2 Hull & White - V

Here we have

[TABLE]

and we remark that we assume $\eta<0$, since the martingale nature of the price under a risk neutral probability is guaranteed if and only if the correlation coefficient is non positive, see [18]. Since the martingale $N^{t}_{s}$ is defined by (30), for this model we have

[TABLE]

with (19) that might be approximated by

[TABLE]

As before, the problem is reduced to evaluating $E_{t}(N^{t}_{s}Y_{s})$, easily done since both processes are lognormally distributed

[TABLE]

Instead, it does not exist an explicit pricing formula for the call price in a Hull & White stochastic volatility model, hence to approximate those $U(t,X_{t},Y_{t})$ and its derivatives we employ the power series expansion introduced in [2], leading to the following approximation

[TABLE]

where

[TABLE]

4.1.3 Heston - V

Here the market model is given by

[TABLE]

and from (30), for this model we have

[TABLE]

so that (19) might be approximated by

[TABLE]

By using Fourier transform techniques, we can evaluate $U$ and its derivatives. Indeed being the Heston model affine, from [17] we have

[TABLE]

where for $j=1,2$, setting $\delta_{1}=1,\delta_{2}=0$, we have

[TABLE]

and, passing the derivative under the integral sign, we may conclude that

[TABLE]

It remains to compute $\mathbf{E}_{t}(N^{t}_{s}\sqrt{Y_{s}})$. Indeed it does not exists an explicit expression of $Y$, so we need a manageable approximation of this process to be able to compute the previous expectation. Since our choice of parameters ($c^{2}<2k\theta$) guarantees the strict positivity of $Y$, we suggest to approximate it by a geometric Brownian motion

[TABLE]

where the deterministic functions $\gamma_{1}$ and $\gamma_{2}$ are chosen appropriately. To determine such functions we perform a moment matching between the two random variables, because we can compute explicitly the first two moments of $Y_{s}$ given $Y_{t}=y$

[TABLE]

while on the other hand we have

[TABLE]

Matching the expectations, we obtain the system

[TABLE]

explicitly solved by

[TABLE]

Consequently, employing this approximation we have

[TABLE]

and when evaluating $\mathbf{E}_{t}(N^{t}_{s}\sqrt{Y_{s}})$, we may approximate it by the same formula as the first of (39), choosing

[TABLE]

4.2 The CIR intensity model

In this section we consider a CIR dynamics for the intensity process, the advantage being that, under the Feller condition, it guarantees the positivity of the process preserving the mean reversion property. Hence we choose

[TABLE]

given that $\sigma^{2}<2q\mu$, with $\sigma,\mu,q>0$.

We know that, by Fourier inversion, the martingale $\displaystyle N^{t}_{s}=\mathrm{e}^{-\int_{t}^{s}\lambda_{u}du}\mathrm{e}^{\varphi(T-s)\lambda_{s}+\psi(T-s)}$ has an explicit expression (see[23]) with

[TABLE]

with $p^{2}=q^{2}+2\sigma^{2}$ and it has dynamics $dN^{t}_{s}=-\sigma\varphi(T-s)N^{t}_{s}\sqrt{\lambda_{s}}dB^{3}_{s}.$

4.2.1 SABR - C I

Combining the above with the SABR model (31), we obtain that

[TABLE]

so that (19) may be approximated by

[TABLE]

It remains to compute $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}}Y_{s})$. Due to the empirical experience that the dependence between default and stochastic volatility is rather weak, practioners often consider them as independent, hence a rough, but handy way, to approximate this factor might be

[TABLE]

being $Y$ a martingale. It remains to compute $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}})$ and we remark that applying the integration by parts formula, we may write

[TABLE]

whence, by considering that the martingale part gives null contribution, we have

[TABLE]

and we decide to approximate this expectation by freezing the $\frac{1}{\sqrt{\lambda}_{s}}$ factor at the initial value and solving the resulting ordinary differential equation we obtain

[TABLE]

4.2.2 Hull and White - C I

In this case we have

[TABLE]

then the first order approximation becomes

[TABLE]

As before, we approximate $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}}\sqrt{Y_{s}})$ by $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}})\mathbf{E}_{t}(\sqrt{Y_{s}})$, with $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}})$ approximated as in (49) and

[TABLE]

4.2.3 Heston - C I

In this case we have

[TABLE]

which implies the approximation formula

[TABLE]

As before, we approximate $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}}\sqrt{Y_{s}})$ by $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}})\mathbf{E}_{t}(\sqrt{Y_{s}})$, with $\mathbf{E}_{t}(N^{t}_{s}\sqrt{\lambda_{s}})$ approximated as in (49) and the factor \mathbf{E}_{t}\Big{[}\sqrt{Y_{s}}\Big{]} can be computed by lognormal approximation exactly as before by (44).

5 The second order approximation

When choosing a CIR model for the intensity process, as shown in [9], even in the simpler Black and Scholes model the effect of the correlation parameter between asset’s price and intensity in the CVA evaluation is more marked and it requires an approximation that accounts for terms of order $\rho^{2}$. This is justified also by the Monte Carlo simulations run in this context. Instead, the effect due to correlation between intensity and volatility always seems to be quite irrelevant.

To capture this behavior, we are going to employ the second order representation formula (28), which leads to a second order expansion in $\rho$. Indeed the first order approximation applied in the previous section, due to the freezing of the terms containing $X$, will work efficiently only when the dependence upon $\rho$ is roughly linear. The computations needed to exploit representation formula (28) are rather lengthy, hence, for the sake of exposition we are going to specify them only for the more commonly used SABR and Heston models coupled with CIR and under the assumption that $\nu=0$, as we remarked that the behaviour in this parameter is usually very well captured by the first order approximation.

For $\nu=0$, $\langle Y,\lambda\rangle_{s}=0$ and consequently ${\partial_{s}}H^{y}_{s}=N^{t}_{s}\frac{d\langle Y,\lambda\rangle_{s}}{ds}=0$, so (28) becomes

[TABLE]

where $\displaystyle H^{x}_{s}=N^{t}_{s}\frac{d\langle X,\lambda\rangle_{s}}{ds}.$

5.1 SABR - C II

In the case of SABR model we have

[TABLE]

so that we have

[TABLE]

Recalling (48) and

[TABLE]

and applying the multidimensional Itô’s formula with $f(w,y,z)=wyz$ respectively to the processes $W_{s}=N^{t}_{s}\sqrt{\lambda_{s}}$, $Y_{s}$, $Z_{s}=\mathrm{e}^{-(1-\gamma)X_{s}}$, we have that

[TABLE]

whence

[TABLE]

So freezing at the initial time the processes $\partial_{x}U(s,X_{s},Y_{s}),\partial^{2}_{xx}U(s,X_{s},Y_{s}),\partial^{2}_{xy}U(s,X_{s},Y_{s})$, $\frac{1}{\sqrt{\lambda}_{s}}$ and $\mathrm{e}^{-(1-\gamma)X_{s}}$, keeping in mind that the martingale parts give no contributions, we arrive at (after some lengthy calculations) the following approximation formula

[TABLE]

where, setting

[TABLE]

we have

[TABLE]

it only remains to evaluate

[TABLE]

and the integrand $\mathbf{E}(N^{t}_{s}\sqrt{\lambda_{s}})$, to be done as in (49).

5.1.1 Heston - C II

In this case we have $\gamma=1$ and $H^{x}_{s}=\rho\sigma N^{t}_{s}\sqrt{\lambda_{s}}\sqrt{Y_{s}}$, so that we have

[TABLE]

We recall that the function $\partial_{x}U$ verifies the same PDE as $U$, which implies that $\partial_{x}U(s.X_{s},Y_{s})$ is again a martingale verifying

[TABLE]

We have

[TABLE]

substituting in the above we obtain

[TABLE]

As before, by freezing $U$ and its derivatives, $\frac{1}{\sqrt{\lambda_{u}}}$ $\frac{1}{\sqrt{Y_{u}}}$ at the initial points and by exploiting the independence between $Y$ and $\lambda$, this formula may be approximated by

[TABLE]

where

[TABLE]

Again we have to approximate the expectations inside the integrals:

[TABLE]

6 Numerical results

To show the numerical efficiency of our method, in this section we present some numerical implementations for the SABR and Heston stochastic volatility models coupled with either a Vasicek or a CIR intensity model. The Hull-White model shares essentially the same behavior, but we are not going to present it, since in this case the efficiency of our correlation method is less evident due to the additional error introduced by the approximation of the zero-th term of the expansion. Numerical experiments were implemented in MatLab (version 9.2.0 - R2017a) on an Intel Core i7 2.40GHZ with 8GB RAM. The CVA was computed for different values of $\rho$ (the underlying-intensity correlation) spanning the interval $(-1,1)$ and we compare the results of our approximation formulas with a full Monte Carlo estimation of (12) with $R=0$, taken as the benchmark CVA value. To reduce the variance of these latter estimates, we used the default-free option value as a control variate: the high correlation between default free and defaultable option price (up to $99\%$ in our experiments) reduces the simulation error by one order of magnitude.

Simulation of the Heston paths were realized through an Euler scheme with full-truncation for the volatility component [22]. For the SABR model, we implemented an exact simulation of the volatility component and a log-Euler scheme for the asset’s price. We considered a uniform time grid $0=t_{0}<t_{1}<\cdots<t_{N}=T$, with $N=10^{3}$ and $M=10^{6}$ sample paths.

The implementation of our approximation formulas is quite straightforward: since the integrands resulted generally well-behaved in all the considered cases, the time (iterated) integrals were evaluated by using a trapezoidal rule with step $\Delta t=10^{-2}$.

In order to obtain consistent results, we fixed the parameters of the two models to fit approximately the same (given) set of default-free call prices at time $t=0$: $\kappa=1.15$, $\theta=0.04$, $c=0.39$, $\eta=-0.34$ and $y=0.034$ for the Heston model, $\gamma=.7367$, $c=0.7356$, $\eta=-0.3$ and $y=.5887$ for the SABR. Without loss of generality we took $r=0$. To keep formulas simple, we always took $\nu=0$ (the volatility-intensity correlation) in the models, since our main interest lied on the interaction between default and asset’s price (WWR). Moreover, extensive numerical simulations confirmed that the contribution due to the terms multiplied by $\nu$ was usually quite modest.

As already mentioned at the end of section 3 the Vasicek model, even though used in the literature (see e.g. [12]), has the drawback of assuming negative values with positive probability. To contain this effect, we arbitrarily chose the Vasicek parameters so that the asymptotic value of the intensity has a small probability ($<0.1\%$) to become negative (Table (1)). In this context, since Monte Carlo simulations seemed to exhibit always a linear behavior in $\rho$, we implemented only the first order approximation, which performed satisfactorily in both SABR and Heston model (Figs. 1, 2, 3 and 4).

When assuming a CIR default intensity, we ran our numerics with three different sets of parameters given in Table (2), we chose the first two arbitrarily, while the third was taken as in [8], [9] where it was considered consistent with observed patterns of CDS spreads and implied volatilities.

Simulated paths were generated with the same full-truncation scheme used for the volatility component of the Heston model.

Figures 5, 6 refer to the first and second order formulas (47) and (54) for the SABR model, while figures 7 and 8 correspond to (51) and (57) for the Heston model. We notice that for Sets 1 and 2 both the first and the second order approximations behave quite well for both models and maturities, uniformly in $\rho$, instead when choosing Set 3 the CVA pattern shows a curvature and the second order approximation seems to be more appropriate. Finally we remark that for small values of $|\rho|$, approximations are very satisfactory for all the models, while they always tend to worsen for growing maturities at large values of $|\rho|$ (only for the negative ones for the Heston model). This effect might be due to the freezing we used to produce the approximation formulas.

Finally, we analyzed the impact of the intensity model parameters on the CVA. Taking the CIR model as an example, Figures 9 and 10 show the impact of the parameters $\sigma$ and $\mu$. In particular, the linear part in (47) and (51) represents the contribution to the CVA due to correlation and we notice that increasing $\sigma$ and $\mu$ determines and increase of the WWR effect ($\rho>0$) in all the considered set of parameters.

Set 4 deserves for a separate discussion. Indeed, even though it strongly violates Feller’s condition for the intensity, we decided to employ it in the simulations, since it was another set consistent with the observed patterns of CDS spreads and implied volatilities (see [9]). This is the most troublesome set and even if it shows still acceptable results for short maturities in both SABR and Heston models, so capturing the Right and Wrong Way risks, it becomes unsatisfactory for larger maturities. One might argue that CIR model is not appropriate for this set of parameters as a noticeable probability to assume negative values is a remarkable modeling flaw. Perhaps, in this context, it could be worth to make alternative choices to describe the default intensity and/or to sharpen the approximations used in our representation formulas. Here we employed the most classical intensity models found in the literature and we postpone to future work more detailed discussion and analysis on the impact of the default model for CVA evaluation.