# Testing Yukawa-like potentials from $f(R)$-gravity in elliptical galaxies

N.R. Napolitano11affiliation: INAF – Osservatorio Astronomico di Capodimonte,
Salita Moiariello, 16, 80131 - Napoli, Italy , S. Capozziello22affiliation:
Dipartimento di Scienze Fisiche, Università di Napoli “Federico II”, Napoli,
Italy 33affiliation: Istituto Nazionale di Fisica Nucleare, Sez. di Napoli,
Italy , A.J. Romanowsky44affiliation: UCO/Lick Observatory, University of
California, Santa Cruz, CA 95064, USA , M. Capaccioli22affiliation:
Dipartimento di Scienze Fisiche, Università di Napoli “Federico II”, Napoli,
Italy , C. Tortora55affiliation: Universit$\rm\ddot{a}$t Z$\rm\ddot{u}$rich,
Institut f$\rm\ddot{u}$r Theoretische Physik, Winterthurerstrasse 190,
CH-8057, Z$\rm\ddot{u}$rich, Switzerland<EMAIL_ADDRESS>
###### Abstract

We present the first analysis of extended stellar kinematics of elliptical
galaxies where a Yukawa–like correction to the Newtonian gravitational
potential derived from $f(R)$–gravity is considered as an alternative to dark
matter. In this framework, we model long-slit data and planetary nebulae data
out to 7 $R_{\rm eff}$ of three galaxies with either decreasing or flat
dispersion profiles. We use the corrected Newtonian potential in a
dispersion–kurtosis Jeans analysis to account for the mass–anisotropy
degeneracy. We find that these modified potentials are able to fit nicely all
three elliptical galaxies and the anisotropy distribution is consistent with
that estimated if a dark halo is considered. The parameter which measures the
“strength” of the Yukawa–like correction is, on average, smaller than the one
found previously in spiral galaxies and correlates both with the scale length
of the Yukawa–like term and the orbital anisotropy.

###### Subject headings:

galaxies : kinematics and dynamics – galaxies : general – galaxies :
elliptical and lenticular – cosmology: theory

##  1\. Introduction

The “concordance” $\Lambda$CDM cosmological model, which includes some unseen
Cold Dark Matter (DM) and a cosmological constant ($\Lambda$) acting as a
repulsive form of Dark Energy (DE), has been remarkably successful in
explaining the formation and evolution of cosmological structures at different
scales (e.g., Springel et al. 2006).

However, at cosmological scales, the cosmological constant as “vacuum state”
of the gravitational field is about 120 orders of magnitude smaller than the
value predicted by any quantum gravity theory (Weinberg 1989) and comparable
to the matter density (coincidence problem), even if they evolved decoupled in
the history of the universe.

In addition, looking at the galaxy scales there are a few critical issues yet
to be solved, which are giving hard time to the whole $\Lambda$CDM framework.

Since the discovery of the flat rotation curves of spiral systems, galaxies
have been the most critical laboratory to investigate the gravitational
effects of the DM halos, to be compared against the expectation of the
cosmological simulations (Navarro et al. 1997, NFW hereafter, Burkert 1995,
Navarro et al. 2010). Here, the $\Lambda$CDM model is not able to fully
explain the shallow central density profile of spiral and dwarf galaxies
(Gilmore et al., 2007; Salucci et al., 2007; Kuzio de Naray et al., 2008).
Early-type galaxies (ETGs hereafter) have been proven only recently to be
consistent with $\Lambda$CDM predictions (and WMAP5 cosmological parameters,
e.g. Komatsu et al. 2009), from their centers (Napolitano et al. 2010) to
their peripheries (Napolitano et al. 2011, N+11 hereafter), although there are
also diverging results showing that ETGs in some cases have too high (Buote et
al. 2007) or too low (e.g. Mandelbaum et al. 2008) concentrations.

This very uncertain context has been a fertile soil for alternative approaches
to the so-called “missing mass”. The basic approach is that the Newtonian
Theory of Gravity, which has been tested only in the Solar System, might be
inaccurate on larger (galaxies and galaxy clusters) scales. The most popular
theory investigated so far, the Modified Newtonian dynamics (MOND) proposed by
Milgrom (1983), is based on phenomenological modifications of Newton dynamics
in order to explain the flat rotation curves of spiral galaxies, and passed a
number of observational tests (Ferreira & Starkman 2009), included ETG
kinematics (Milgrom & Sanders 2003; Tiret et al. 2007; Kroupa et al. 2010;
Cardone et al. 2010; Richtler et al. 2011). Only lately it has been derived in
a cosmological context (Bekenstein, 2004).

A new approach, motivated from cosmology and quantum field theories on a
curved space-time, has been proposed to study the gravitational interaction:
the Extended Theories of Gravity (Capozziello 2002; Capozziello & Faraoni
2011). In particular, the so called $f(R)-$gravity seem to have passed
different observational tests like spiral galaxies’ rotation curves, X-ray
emission of galaxy clusters and cosmic acceleration (see e.g. Capozziello et
al. 2007a, C+07 hereafter, Capozziello et al. 2009, C+09 hereafter,
Capozziello et al. 2008). This approach is based on a straightforward
generalization of Einstein theory where the gravitational action (the Hilbert-
Einstein action) is assumed to be linear in the Ricci curvature scalar $R$. In
the case of $f(R)$-gravity, one assumes a generic function $f$ of the Ricci
scalar $R$ (in particular analytic functions) and asks for a theory of gravity
having suitable behaviors at small and large scale lengths.

As shown in Capozziello et al. (2009), analytic $f(R)$-models give rise, in
general, to Yukawa–like corrections to the Newtonian potentials in the weak
field limit approximation (see also Lubini et al. 2011). The correction
introduces a new gravitational scale, besides the standard Schwarzschild one,
depending on the dynamical structure of the self-gravitating system.

Here we want to test these Yukawa-like gravitational potentials against a
sample of elliptical galaxies. This approach has been proposed earlier, in a
phenomenological scheme for anti–gravity, to model flat rotation curves of
spiral galaxies (Sanders, 1984), and recently, in $f(R)$ theories to model
disk galaxies combined with NFW haloes (see Cardone & Capozziello 2011). The
test we are proposing at galaxy scales is crucial: reproducing kinematics and
then dynamics of these very different classes of astrophysical systems in the
realm of the same paradigm is needed to test these new gravitational theory as
an alternative to DM which has not been definitely found out at fundamental
level.

The layout of the paper is the following. In §2, we sketch the main
ingredients of $f(R)$-gravity deriving, in the weak field limit, the Yukawa-
like corrected gravitational potential. §3 is devoted to the high-order Jeans
analysis suitable for ellipticals. The dispersion-kurtosis fitting and the
data sample are presented in §4. Discussion and conclusions are in §5.

##  2\. Post- Newtonian potentials from $f(R)$-gravity

We are interested in testing a class of modified potentials which naturally
arise in post-Newtonian approximation of $f(R)-$gravity for which no
particular choice of the Lagrangian has been provided.

The starting point is a general gravity action of the form

$\displaystyle\mathcal{A}\,=\,\int
d^{4}x\sqrt{-g}\biggl{[}f(R)+\mathcal{X}\mathcal{L}_{m}\biggr{]}\,,$ (1)

where $f(R)$ is an analytic function of Ricci scalar, $g$ is the determinant
of the metric $g_{\mu\nu}$, ${\displaystyle\mathcal{X}=\frac{16\pi G}{c^{4}}}$
is the gravitational coupling constant, and $\mathcal{L}_{m}$ describes the
standard fluid-matter Lagrangian. Such an action is the straightforward
generalization of the Hilbert-Einstein action obtained as soon as $f(R)=R$.

In Capozziello et al. (2009, and reference therein) it has been shown that if
one solves the field equations in the weak field limit under the general
assumption of an analytic Taylor expandable $f(R)$ functions of the form

$\displaystyle f(R)\simeq f_{0}+f_{1}R+f_{2}R^{2}+f_{3}R^{3}+...\,$ (2)

the following gravitational potential arises

$\displaystyle\Phi\,=\,-\left(\frac{GM}{f_{1}r}+\frac{L\delta_{1}(t)~{}e^{-\frac{r}{L}}}{6~{}r}\right),$
(3)

where $L\doteq\displaystyle-\frac{6f_{2}}{f_{1}}$, $f_{1}$ and $f_{2}$ are the
expansion coefficients obtained by Taylor expansion. We note that the $L$
parameter is related to the effective mass
$m=(-3/L^{2})^{-1/2}=(2f_{2}/f_{1})^{1/2}$ and can be interpreted also as an
effective length.

From Eq. 3, the standard Newton potential is recovered only in the particular
case $f(R)=R$. Furthermore, the parameters $f_{1}$ and $f_{2}$ and the
function $\delta_{1}$ represent the deviations with respect the standard
Newton potential. On the Solar system scale, it has been shown that
Yukawa–like deviation from the pure Newtonian potential are not in
contradiction with classical tests of General Relativity (see e.g. Capozziello
& Tsujikawa 2008), thanks to the so-called Chameleon mechanism (Khoury &
Weltman 2004). In particular, $f_{1}$ and $f_{2}$ parameters are expected to
allow the regular Newtonian potential, while at larger scales they can assume
non-trivial values (e.g. $f_{1}\neq 1,\,\delta_{1}(t)\neq 0,\,\xi\neq 1$, see
Capozziello et al. 2007b, Capozziello et al. 2009).

Eq. (3) can be recast as

$\Phi(r)=-\frac{GM}{(1+\delta)r}\left(1+\delta e^{-\frac{r}{L}}\right)\,,$ (4)

where the first term is the Newtonian–like part of the potential associated to
baryonic point–like mass $M/(1+\delta)$ (no DM) and the second term is a
modification of the gravity including a “scale length”, $L$ associated to the
above coefficient of the Taylor expansion. If $\delta=0$ the Newtonian
potential is recovered. Comparing Eqs. 3 and 4, we obtain that
$1+\delta=f_{1}$, and $\delta$ is related to $\delta_{1}(t)$ through

$\delta_{1}=-\frac{6GM}{L^{2}}\frac{\delta}{1+\delta}$ (5)

where $6GM/L^{2}$ and $\delta_{1}$ can be assumed quasi-constant. From Eq. 5,
it turns out that $L\propto\sqrt{-\delta/(1+\delta)}$. Due to the arbitrarity
of $\delta_{1}(t)$, the actual value of the $\delta$ parameter can assume any
values, however, in order to have a Yukawa potential with a non imaginary
exponent (i.e. $L$ must be real) it is required that $\xi<0$ or $-1<\delta<0$.
As comparison, Sanders (1984) adopted the same potential as in Eq. 4 under the
assumption of anti–gravity generated by massive particles (of mass $m_{0}$)
carrying the additional gravitational force. In this case a typical scale
length would naturally arise ($L=h/m_{0}c$ being a Compton length) and a
$-1<\delta<0$ would provide a repulsive term to the Newtonian–like term
producing flat rotation curves at $r\gg L$ as observed in spiral galaxies. In
particular, for a small sample of spiral systems Sanders (1984) found
$-0.95\mathrel{\hbox to0.0pt{\lower 3.5pt\hbox{\hskip 0.5pt$\sim$}\hss}\raise
0.5pt\hbox{$<$}}\delta\mathrel{\hbox to0.0pt{\lower 3.5pt\hbox{\hskip
0.5pt$\sim$}\hss}\raise 0.5pt\hbox{$<$}}-0.92$.

Figure 1.— Circular velocity produced by the modified potential in Eq. 4 for
the two galaxies N4494 (top) and N4374 (bottom). In both cases the $M/L_{*}$
has been fixed to some fiducial value (as expected from stellar population
models and Kroupa 2001 IMF): $M/L_{*}=4.3\Upsilon_{\odot,B}$ for NGC 4494 and
$M/L_{*}=5.5\Upsilon_{\odot,V}$ for NGC 4374. The potential parameters adopted
are: $L=250^{\prime\prime}$ and $\delta$=0, -0.65, -0.8, -0.9 (lighter to
darker solid lines) and $L=180^{\prime\prime}$ and $\delta$=-0.8 (dashed
lines). The dotted line is a case with positive coefficient of the Yukawa–like
term and $L=5000^{\prime\prime}$ which illustrates that positive $\delta$
cannot produce flat circular velocity curves. Finally some reference
Navarro–Frenk–White (NFW) models are shown as dot-dashed lines.

Here we want to test the modified potential as in Eq. 4 in elliptical galaxies
and check whether it is able to provide a reasonable match to their kinematics
and how the model parameters compare with the results obtained from spiral
systems. We construct equilibrium models based on the solution of the radial
Jeans equation (see §3) to interpret the kinematics of planetary nebulae (PNe,
see Napolitano et al. 2002, 2005; Romanowsky et al. 2003; Coccato et al 2009)
which are the only stellar–like tracers for galaxy dynamics available in ETGs
out to $\sim 5-10$ effective radii ($R_{\rm eff}$).

We will use the inner long slit data and the extended PN kinematics for three
galaxies which have published dynamical analyses within DM halo framework: NGC
3379 (Douglas et al. 2007; De Lorenzi et al. 2009, DL+09 hereafter), NGC 4494
(Napolitano et al. 2009, N+09), NGC 4374 (N+11). The decreasing velocity
dispersion profiles of the first two galaxies have been modeled with an
intermediate mass halo, $\log M_{\rm vir}\sim 12-12.2M_{\odot}$, with
concentration $c_{\rm vir}=6-8$ and a fair amount of radial anisotropy in the
outer regions. For NGC 4374, having a rather flat dispersion profile, a more
massive (adiabatically contracted) halo with $\log M_{\rm vir}\sim
13.4M_{\odot}$ and $c_{\rm vir}\sim 7$ was required with a negligible amount
of anisotropy in the outer regions. These models turned out to be in fair
agreement with the expectation of WMAP5 $c_{\rm vir}-M_{\rm vir}$ relation
(N+11) and with a Kroupa (2001) IMF, making this sample particularly suitable
for a comparison with an alternative theory of gravity with no-DM as we want
to propose here.

Before we go on with detailed stellar dynamics, we show in Fig. 1 the circular
velocity of the modified potential as a function of the potential parameters
$L$ and $\delta$ for NGC 4494 and NGC 4374. As for the spiral galaxies,
negative values of the $\delta$ parameter make the circular velocity more and
more flat also reproducing the typical dip (e.g. NGC 4374) of the circular
velocity found for the DM models (dot–dashed curves) of the most massive
systems. On the contrary, positive $\delta$ values cannot produce flat
circular velocity curves (see Fig. 1).

Figure 2.— Dispersion in km s-1 (top) and kurtosis (bottom) fit of the galaxy
sample for the different $f(R)$ parameter sets: the anisotropic solution
(solid lines) is compared with the isotropic case (dashed line – for NGC 4374
and NGC 4494 this is almost indistinguishable from the anisotropic case). From
the left, NGC 4494, NGC 3379 and NGC 4374 are shown with DM models as gray
lines from N+09, DL+09 (no kurtosis is provided), and N+11 respectively.

###  2.1. A consistency check with galaxy scaling relations

To conclude the inspection of the modified potential as in Eq. 4 here we want
to show that, beside flat rotation curves, this also naturally accounts for
fundamental scaling relations of galaxies: the Tully-Fisher (TF) relation for
spirals and Faber-Jackson (FJ) relation for ETGs.

Both relations connect the total mass $M$ of galaxies with some characteristic
velocity defining the kinetic energy of the systems (i.e. the maximum rotation
velocity, $v_{\rm max}$, for spirals and the central velocity dispersion,
$\sigma_{0}$, for ETGs). In either cases the kinematical quantities involved
are proportional to the circular velocity of the systems through some
“structure” constant, thus the arguments below apply to galaxies in general.

Although the point–like version of the potential implies that the circular
velocity $v_{c}$ scales with mass as $M\sim v_{\rm c,max}^{2}$ (as pointed out
by Sanders 1984), if one derives the circular velocity for an extended galaxy
this can be generalized as

$v_{c}^{2}(r)=(GM_{\rm tot}/r_{*})\times f(r/r_{*};\delta,L/r_{*})$ (6)

where $r_{*}$ is a characteristic radius (e.g. the disk length for spirals or
the effective radius encircling half of the galaxy light for ETGs),
$f(r/r_{*};~{}\delta,~{}L/r_{*})$ is a generic function which includes the
radial dependence of the enclosed mass and the above Yukawa-like term. This
function is defined such as, for $\delta=0$, it gives $v_{c}^{2}(r)=GM(r)/r$
as the usual Newtonian expression. It is easy to show that if galaxies are
homologous the maximum of $v_{c}$ is reached at the same $r/r_{*}$, for a
given $\delta$ and $L/r_{*}$ and this maximum can be written as

$v_{\rm c,max}^{2}={\rm K}M_{\rm tot}/r_{*}$ (7)

where the constant $K$ depends on the set of parameters
{$\delta,~{}L,~{}r_{*}$} adopted. In Eq. 7, though, $M_{\rm tot}$ and $r_{*}$
are linked by the size–mass relation which is generally written as
$r_{*}\propto M_{\rm tot}^{\alpha}$, from which Eq. 7 can be written as

$v_{\rm c,max}^{2}\propto M_{\rm tot}^{1-\alpha}.$ (8)

The size-mass relation of spiral galaxies can be found in Persic et al. (1996,
see also ) to be $r_{*}\propto M_{\rm tot}^{0.4}$, while it is $r_{*}\propto
M_{\rm tot}^{0.6}$ for ETGs (e.g. Shen et al. 2003, Napolitano et al. 2005).
This would give a TF slope of $3.33$ and FJ slope of $5$ which are both in the
range of the observed relations (see e.g. McGaugh 2005 and Nigoche-Netro et
al. 2010 respectively) with the remaining discrepancy being mainly due to the
conversion factor to the observed quantities and non homologies.

We finally remark that the TF relation has been found not to be conflicting
with $f(R)$ potentials in Capozziello et al. (2006), although the potentials
from $f(R)\propto R^{n}$ adopted there are just a series expansion of the
Yukawa–like potential coming out from a more general polynomial f(R) as in Eq.
3.

##  3\. High-order Jeans analysis

From the model point of view, the problem of fitting a modified potential as
in Eq. 4 (which is formally self-consistent since the source of the potential
is the only mass of the dynamical tracers, i.e. stars) implies the same kind
of degeneracies between the anisotropy parameter,
$\beta=1-\sigma_{\theta}^{2}/\sigma_{r}^{2}$ (where $\sigma_{\theta}$ and
$\sigma_{r}$ are the azimuthal and radial dispersion components in spherical
coordinates), and the non–Newtonian part of the potential (characterized by
two parameters like typical dark haloes) in a similar way of the classical
mass-anisotropy degeneracy. We have shown (N+09, N+11) that these degeneracies
can be alleviated via higher-order Jeans equations including in the dynamical
models both the dispersion111For the slow–rotating models we use the velocity,
$v_{\rm rms}=\sqrt{v^{2}+\sigma^{2}}$ as a measure of the velocity dispersion.
($\sigma_{\rm p}$) and the kurtosis ($\kappa$) profiles of the tracers.

In the following, we will use the assumption of spherical symmetry since
galaxies in the sample are all E0–E1 for which, if one exclude the singular
chance that they are all flattened systems seen face–on (see discussion in
Sect. 8.1. of Douglas et al. 2007), the spherical approximation is good at 10%
(Kronawitter et al. 2000)222 The effect of non–spherical models is outside the
scope of this paper, but details for NGC 3379 and NGC 4494 can be found in
DL+09 and N+09.. Under spherical assumption, no-rotation and $\beta=\rm const$
(corresponding to the family of distribution functions
$f(E,L)=f_{0}L^{-2\beta}$, see Łokas 2002 and references therein333 Here,
there is the caveat that the solution of Jeans Equations does not ensure that
the final distribution function is non negative and thus fully physical (see
e.g., An & Evans 2006)., the 2-nd and 4-th moment radial equations can be
compactly written as:

$s(r)=r^{-2\beta}\int_{r}^{\infty}x^{2\beta}H(x)dx$ (9)

where $s(r)=\\{\rho\sigma_{r}^{2};\rho\overline{v_{r}^{4}}\\}$, $\beta$ is the
anisotropy parameter, and

$H(r)=\left\\{\rho\frac{d\Phi}{dr};3\rho\frac{d\Phi}{dr}\overline{v_{r}^{2}}\right\\}$

respectively for the dispersion and kurtosis equations, being the latter
$\kappa(r)={\overline{v_{r}^{4}}}/\sigma_{r}^{4}$. In the same equations,
$\Phi(r)$ is the spherical extended source version of the point–like potential
as in Eq. 4444This is obtained assuming the onion shell approximation:
$\Phi(r)=\int_{0}^{r}\int_{0}^{2\pi}\int_{0}^{\pi}\phi(r)~{}r^{2}\sin\theta~{}d\theta
d\varphi dr$, see also Eq. 18 of C+09. and $\rho(r)$ is the 3D density of the
tracer obtained by multiplying the deprojection of the stellar surface
brightness profile, $j_{\star}(r)$, by some constant stellar mass-to-light
ratio, $M/L_{\star}$.

This $M/L_{\star}=const$ might be a strong assumption to check further in a
separate paper as it neglects the presence of stellar population gradients
(see e.g. Tortora et al. 2010). However, colour (and $M/L$) gradients are
generally stronger within $R_{\rm eff}$ (see e.g. Tortora et al 2010) and
might mainly drive the best fit in the central regions, while they are
possibly shallower outside (e.g., Tamura & Ohta 2003) where the $f(R)$
parameters should be better constrained. In the following, $j_{\star}(r)$ is
derived by photometry presented in previous dynamical studies (i.e. DL+09,
N+09, N+11 for NGC 3379, NGC 4494 and NGC 4374 respectively).

Eqs. 9 are the ones interested by the potential modification and include four
free parameters to be best-fitted: the $f(R)$ parameters {$\delta,L$}, the
“dynamically inferred” stellar mass-to-light ratio $M/L_{\star}$ and the
constant anisotropy $\beta$ (see also §4). The solutions of Eqs. 9 on a
regular grid in the parameter space are then projected to match the observed
line-of-sight kinematical profile via ordinary Abel integrals (see N+09 for
details).

As mentioned earlier, Eqs. 9 are written under the assumption of a constant
$\beta$ with radius, which provides a average global anisotropy distribution
over all the galaxy. As seen in previous analyses (e.g. N+09, DL+09 and N+11),
it is likely that this might not be a fair assumption, as $\beta$ turns out to
be constant somewhere in the outer regions, but strongly varying in the
central radii. In this preliminary test we will skip this implementation of
the models since we expect this to possibly improve the fit to the data in the
central part only, where we do not expect the overall dynamics of being
strongly ruled by the $f(R)$ potential, whose parameters are the main focus of
this work. Furthermore, we have shown previously (see e.g. N+09 and N+11) that
the assumption of the constant or radial varying anisotropy did not strongly
affect the determination of the other important parameter, the dynamically
based stellar $M/L$. In the following we will take the $\beta=$const as fair
estimate of the average galaxy anisotropy.

##  4\. Dispersion-kurtosis fitting

In Fig. 2 we show the dispersion and kurtosis profiles of the three galaxies
with the $f(R)$ models superimposed (solid lines). The fitting procedure is
based on the simultaneous $\chi^{2}$ minimization of the dispersion and
kurtosis profiles over a regular grid in the parameter space. The best–fit
parameters are summarised in Table 1 together with some info of the galaxy
sample.

Overall the agreement of the model curves with data is remarkably good and it
is comparable with models obtained with DM modeling (gray lines in Fig. 2).

In all cases, the $f(R)$ models allow to accommodate a constant orbital
anisotropy $\beta$ which is very close to the estimates from the DM models
(see e.g. Table 1555For NGC 4374 only to be nicely fitted at all radii we
needed to include some radial anisotropy in the very central regions,
following the $\beta(r)$ distribution adopted in N+11 (see Eq. 5, whith best
fit $r_{a}=22.5$).). This is mainly guaranteed by the fit to the $\kappa(r)$
which does not respond much differently to the modified potential with respect
the DM models. Thus, an important result of the analysis is that the orbital
anisotropy is fairly stable to the change of the galaxy potential. In
particular, the use of the kurtosis profiles has allowed us to solve the
degeneracy of the models and favor the anisotropic solutions for NGC 3379 and
NGC 4494 (NGC 4374 being almost isotropic everywhere). Although the isotropy
solutions provide also a good fit for the dispersion profile only (see e.g.
the dashed lines in Fig. 2), they not correctly match the observed $\kappa$.
This produces a significantly worse total $\chi^{2}$/dof (NGC 3379: 45/26;
NGC4374: 35/40; NGC 4494: 27/44) with respect to the best–fit in the Table 1,
although still close to $\chi^{2}$/dof$~{}\sim 1$ mainly because of the large
error bars.

Table 1Model parameters for the $f(R)$ potential.

Galaxy | Mag (band) | $R_{\rm eff}$ | $M/L_{\star}$ | $L$ | $\delta$ | $\beta$ |  $\chi^{2}$/dof
---|---|---|---|---|---|---|---
NGC3379 | -19.8(B) | 2.2 | 6 (7) | 6 | -0.75 | 0.5($<$0.8) | 14/25
NGC4374 | -21.3(V) | 3.4 | 6 (6) | 24 | -0.88 | 0.01(0.01) | 14/39
NGC4494 | -20.5(B) | 6.1 | 3 (4) | 20 | -0.79 | 0.5(0.5) | 18/43

NOTES – Galaxy ID, total magnitude, effective radius and model parameters for
the unified solution. DM–based estimates for $M/L_{\star}$ and $\beta$ (NGC
3379: DL+09; NGC 4374: N+11; NGC 4494: N+09) are shown in parentheses for
comparison. $M/L_{\star}$ are in solar units, $R_{\rm eff}$ and $L$ in kpc.
Typical errors on $M/L_{\star}$ are of the order of 0.2$M/L_{\odot}$ and on
$\beta$ of 0.2 (see also Fig. 3). The small $\chi^{2}$ values are mainly due
to the large data error bars.

Finally, the best fit $M/L_{\star}$ in Table 1 are very similar to the values
found for DM models (reported between brackets) in all cases, generally
consistent with a Kroupa (2001) IMF.

Looking at the $f(R)$ parameters, in Fig. 3 we show the marginalized
confidence contours of the main two potential parameters for the three
galaxies. As also reported in Table 1, the $\delta$ parameter has a mean value
$\delta=-0.81\pm 0.07$ which is inconsistent with the one previously found for
spiral galaxies (e.g. Sanders 1984, also shown in Fig. 3). On the contrary,
$\delta$ seems nicely correlated with the other potential parameter, $L$, as
expected from Eq. 5. In the same figure the correlation is supported by the
tentative fit into the $\delta-L$ plane (whether or not the spiral galaxy
sample is included in the fit), although the sample is too small to drive any
firm conclusion.

Interestingly, there seems to be a possible increasing trend of $\delta$ with
the orbital anisotropy: this is also shown in Fig. 3 where we have added the
fiducial value obtained for the spiral sample (having assumed a reference
$\beta=-1$ for fiducial tangential anisotropy for late-type systems, see e.g.
Battaglia et al. 2005). This evidence leaves room for an interpretation of
$\delta$ and the physics of the galaxy collapse (e.g. the spherical infall
model, Gunn & Gott 1972; Gunn 1977).

In fact, as discussed in §2, $\delta$ is linked to $\delta_{1}$, which is an
arbitrary function that comes out because the field equations in the post-
Newtonian approximation, depending only on the radial coordinate. From a
physical point of view, such a function could be related to second order
effects related to anisotropies and non-homogeneities which could trigger the
formation and the evolution of stellar systems. To take into account such a
situation, one should perform the post-Newtonian limit of the theory not only
in the simple hypothesis of homogeneous spherical symmetry (Schwarzschild
solution) but also considering more realistic situation as
Lema$\hat{i}$tre–Tolman–Bondi solutions (see e.g. Herrera et al. 2010).

Figure 3.— Top: 1- and 2-$\sigma$ confidence levels in the $\delta-L$ space
marginalized over $M/L_{\star}$ and $\beta$ (see also Table 1). Spiral galaxy
results from Sanders (1984) are shown as empty triangle with error bars. Solid
(dashed) curve shows the tentative best-fit to the data including (excluding)
the spiral galaxies and assuming a $L\propto\sqrt{\delta/(1+\delta)}$
correlation as expected from Eq. 5. Bottom: the anisotropy and the $\delta$
parameters turn out to be correlated for the elliptical sample (full squares).
This correlation seems to include also the spiral sample cumulatively shown as
the empty triangle (here we have assumed $\beta=-1.0\pm 0.5$ as a fiducial
value for spiral galaxies to draw a semi-quantitative trend across galaxy
types).

##  5\. Discussion and Conclusions

There is a growing attention to alternative model to the $\Lambda$CDM paradigm
as the latter is still suffering some discrepancies at the galaxy scales and,
most importantly, is based on the assumption of the existence of two
ingredients (DM and DE) whose nature is still unknown.

Different attempts have tried to circumvent the problem by introducing a
modified dynamics, e.g. with the MOND theory (see Sanders & McGaugh 2002,
Swaters et al. 2010, Cardone et al. 2010), but this seems still needing some
DM at least to cluster scales which might be still consistent with the
primordial nucleosynthesis (e.g. via high energy neutrinos, Angus et al. 2010)
and does not provide an explanation for the DE.

Lately $f(R)$–gravity models have made their step out as a natural explanation
for the two dark ingredients of the Universe assuming that they are related to
the fact that gravitational interaction could present further degrees of
freedom whose dynamical effects emerge at large scales (Capozziello & Faraoni
2011). In this paper we have checked the Yukawa-like modification to the
Newtonian potential obtained as post-Newtonian approximation of $f(R)-$gravity
for which no particular choice of the Lagrangian has been provided, with the
only assumption that $f(R)$ is analytic function.

We have used a combination of long-slit spectroscopy and planetary nebulae
kinematics out to $\sim$7 $R_{\rm eff}$ in three systems (NGC 3379, NGC 4374,
NGC 4494) for which $\Lambda$CDM models turned out to be fairly consistent
with WMAP5 measurements (see N+11 for a discussion).

Due to the small galaxy sample, the spirit of this analysis has been to check
whether $1)$ the modified potential introduced by the $f(R)$-gravity allowed a
fit to the galaxy kinematics comparable to the DM models; $2)$ the three
galaxies returned a parameter $\delta$ which is comparable with spiral
galaxies (Sanders 1984).

We have found that the modified potentials allow to nicely model the three
galaxies with a distribution of the $\delta$ parameters which turned out to be
inconsistent with the results found in spiral systems. We have shown some
hints that $\delta$ might be correlated with the galaxy anisotropy, $\beta$,
and the scale parameter, $L$, with elliptical and spiral galaxies following
the same pattern.

This evidence can have interesting implication on the ability of the theory to
make predictions on the internal structure of the gravitating systems after
their spherical collapse (e.g. Gunn 1977) which has to be confirmed on a
larger galaxy sample which we expect to do in a near future.

Despite of some simplifications on the model adopted (e.g. constant $M/L$ and
anisotropy across the galaxy) and the degeneracies between the model
parameters, the results are very encouraging. The fit to the data is very good
in all cases and both the stellar $M/L$ (with Kroupa IMF generally favored)
and orbital anisotropy turn out to be similar to the one estimated if a dark
halo is considered.

Getting a modified gravity to work self-consistently for all gravitating
systems in general, and all galaxy families in particular, is a very non-
trivial challenge that has foiled other theories (e.g. MOND).

We thank the anonymous referee for constructive comments which allowed us to
significantly improve the paper. AJR was supported by National Science
Foundation Grants AST-0808099 and AST-0909237. CT was supported by the Swiss
National Science Foundation.

## References

  * An & Evans (2006) An, J. H., & Evans, N. W. 2006, AJ, 131, 782
  * Angus et al. (2010) Angus, G. W., Famaey, B., & Diaferio, A. 2010, MNRAS, 402, 395
  * Battaglia et al. (2005) Battaglia, G., Helmi, A., Morrison, H., et al. 2005, MNRAS, 364, 433
  * Bekenstein (2004) Bekenstein, J. D. 2004, Phys. Rev. D, 70, 083509
  * Buote et al. (2007) Buote, D. A., et al., 2007, ApJ, 664, 123
  * Burkert (1995) Burkert, A. 1995, ApJL, 447, L25
  * Capozziello (2002) Capozziello, S. 2002, International Journal of Modern Physics D, 11, 483
  * Capozziello et al. (2006) Capozziello, S., Cardone, V. F., & Troisi, A. 2006, JCAP, 8, 1
  * Capozziello et al. (2007a) Capozziello, S., Cardone, V. F., Troisi, A., 2007, MNRAS, 375, 1423 (C+07)
  * Capozziello et al. (2007b) Capozziello, S., Stabile, A., & Troisi, A. 2007, Phys. Rev. D, 76, 104019
  * Capozziello et al. (2008) Capozziello, S., et al., 2008a, Phys. Rev. D 78, 063504
  * Capozziello & Tsujikawa (2008) Capozziello, S., & Tsujikawa, S. 2008, Phys. Rev. D, 77, 107501
  * Capozziello et al. (2009) Capozziello, S., Stabile, A., & Troisi, A. 2009, Modern Physics Letters A, 24, 659
  * Capozziello et al. (2009) Capozziello, S., de Filippis, E., & Salzano, V. 2009, MNRAS, 394, 947 (C+09)
  * Capozziello & Faraoni (2011) Capozziello, S. & Faraoni V. 2011, Beyond Einstein Gravity: A Survey Of Gravitational Theories For Cosmology And Astrophysics, Springer, New York .
  * Cardone et al. (2010) Cardone, V. F., Angus, G., Diaferio, A., Tortora, C., & Molinaro, R. 2011, MNRAS, in press, arXiv:1011.5741
  * Cardone & Capozziello (2011) Cardone, V. F., & Capozziello, S. 2011, MNRAS, in press, arXiv:1102.0916
  * Coccato et al (2009) Coccato, L., et al. 2009, MNRAS, 394, 1249
  * De Lorenzi et al. (2009) De Lorenzi, F., et al. 2009, MNRAS, 395, 76 (DL+09)
  * Douglas et al. (2007) Douglas, N. G., et al. 2007, ApJ, 664, 257
  * Ferreira & Starkman (2009) Ferreira, P. G., & Starkman, G. D. 2009, Science, 326, 812
  * Gilmore et al. (2007) Gilmore, G., et al., 2007, ApJ, 663, 948
  * Gunn & Gott (1972) Gunn, J. E., & Gott, J. R., III 1972, ApJ, 176, 1
  * Gunn (1977) Gunn, J. E. 1977, ApJ, 218, 592
  * Herrera et al. (2010) Herrera L., et al., 2010, PRD 82, 024021
  * Khoury & Weltman (2004) Khoury, J., & Weltman, A. 2004, Phys. Rev. D, 69, 044026
  * Komatsu et al. (2009) Komatsu, E., et al. 2009, ApJS, 180, 330
  * Kronawitter et al. (2000) Kronawitter, A., Saglia, R. P., Gerhard, O., & Bender, R. 2000, A&AS, 144, 53
  * Kroupa (2001) Kroupa, P., 2001, MNRAS, 322, 231
  * Kroupa et al. (2010) Kroupa, P., Famaey, B., de Boer, K. S., et al. 2010, A&A, 523, 32
  * Kuzio de Naray et al. (2008) Kuzio de Naray, R., McGaugh, S. S., & de Blok, W. J. G., 2008, ApJ, 676, 920
  * Łokas (2002) Łokas, E. L. 2002, MNRAS, 333, 697
  * Lubini et al. (2011) Lubini, M., Tortora, C., Näf, J., Jetzer, Ph., Capozziello, S. 2011, accepted for publication in EPJC, preprint: arXiv:1104.2851
  * Mandelbaum et al. (2008) Mandelbaum, R., Seljak, U., & Hirata, C. M. 2008, JCAP, 8, 6
  * McGaugh (2005) McGaugh, S. S. 2005, ApJ, 632, 859
  * Milgrom (1983) Milgrom, M. 1983, ApJ, 270, 365
  * Milgrom & Sanders (2003) Milgrom, M., & Sanders, R. H. 2003, ApJ, 599, L25
  * Moore et al. (1999) Moore, B., Quinn, T., Governato, F., Stadel, J., Lake, G., 1999, MNRAS, 310, 1147
  * Navarro et al. (1997) Navarro, J.F., Frenk, C.S., White, S.D., 1997, ApJ, 490, 493 (NFW)
  * Navarro et al. (2010) Navarro, J. F., et al. 2010, MNRAS, 402, 21
  * Napolitano et al. (2002) Napolitano, N.R., Arnaboldi, M., Capaccioli, M., 2002, A&A, 383, 791
  * Napolitano et al. (2005) Napolitano, N. R., et al. 2005, MNRAS, 357, 691
  * Napolitano et al. (2009) Napolitano, N. R., et al. 2009, MNRAS, 393, 329 (N+09)
  * Napolitano et al. (2010) Napolitano, N. R., Romanowsky, A. J., & Tortora, C. 2010, MNRAS, 405, 2351
  * Napolitano et al. (2011) Napolitano, N. R., et al. 2011, MNRAS, 411, 2035 (N+11)
  * Nigoche-Netro et al. (2010) Nigoche-Netro, A., Aguerri, J. A. L., Lagos, P., et al. 2010, A&A, 516, 96
  * Persic et al. (1996) Persic, M., Salucci, P., & Stel, F. 1996, MNRAS, 281, 27
  * Richtler et al. (2011) Richtler, T., Famaey, B., Gentile, G., & Schuberth, Y. 2011, A&A, 531, 100
  * Romanowsky et al. (2003) Romanowsky, A. J., et al., 2003, Science, 301, 1696
  * Sanders (1984) Sanders, R. H. 1984, A&A, 136, L21
  * Sanders (1986) Sanders, R. H. 1986, MNRAS, 223, 539
  * Sanders & McGaugh (2002) Sanders, R. H., & McGaugh,S. S. 2002, ARA&A, 40, 263
  * Salucci et al. (2007) Salucci, P., et al., 2007, MNRAS, 378, 41
  * Shen et al. (2003) Shen, S., Mo, H. J., White, S. D. M., et al. 2003, MNRAS, 343, 978
  * Springel et al. (2006) Springel, V., Frenk, C. S., & White, S. D. M. 2006, Nature, 440, 1137
  * Swaters et al. (2010) Swaters, R. A., Sanders, R. H., & McGaugh, S. S. 2010, ApJ, 718, 380
  * Tamura & Ohta (2003) Tamura, N., & Ohta, K. 2003, AJ, 126, 596
  * Thomas et al. (2009) Thomas, J., Saglia, R. P., Bender, R., et al. 2009, ApJ, 691, 770
  * Tiret et al. (2007) Tiret, O., Combes, F., Angus, G. W., Famaey, B., & Zhao, H. S. 2007, A&A, 476, L1
  * Tortora et al. (2009) Tortora, C., Napolitano, N. R., Romanowsky, A. J., Capaccioli, M., & Covone, G. 2009, MNRAS, 396, 1132
  * Tortora et al (2010) Tortora, C. et al., 2010, MNRAS, 407, 144
  * Weinberg (1989) Weinberg S., 1989, Rev. Mod. Phys. 61, 1


In this paper, we focus on finding the optimal hedging strategy of a credit index option using reinforcement learning.
We take a practical approach, where the focus is on realism i.e. discrete time, transaction costs; even testing our policy on real market data.
We apply a state of the art algorithm, the Trust Region Volatility Optimization (TRVO) <cit.> algorithm and show that the derived hedging strategy outperforms the practitioner's Black & Scholes delta hedge.

§ INTRODUCTION

Hedging consists in investing to reduce the risk of adverse price movements of financial instruments, and it is one of the main concerns in finance. In this paper we focus on the concept of option hedging, where an option is a contract which offers the buyer the opportunity to buy or sell the underlying asset at a predefined strike price in the future. In particular, the options considered here are credit index options i.e., the underlying is a Credit Default Swap (CDS) index. Option hedging is based on a mathematical theory started with Black & Scholes (B&S) <cit.>. This theory is motivated by a strong set of assumptions which tend to be unrealistic <cit.>. In particular, hedging is assumed to be done costlessly and continuously.

Several approaches have been proposed to extend the B&S model to account for transaction costs, starting with <cit.> and more recently <cit.> which uses stochastic optimal control. The main difference with respect to these approaches is that is data-driven and model-free i.e., it thus does not require any assumptions on the dynamics of the assets.

Credit index options market makers have the target of making profit without keeping open risk positions. The most straightforward strategy is to buy and sell the same amount of the same option in order to have a return from the difference between the two prices (the difference between the buy and sell price of a security is called bid-ask spread). Given the low liquidity of these options, most of the time this is not possible, so the market maker's portfolio results in a combination of different options, and she needs to hedge at least the risk given by the underlying instrument i.e. the delta risk.
It is possible to hedge this risk following blindly the B&S delta hedge, perhaps with an automated software connected directly to the market, but, specifically in cases with high transaction such as in the case of the CDS index which we are analyzing, this can quickly become quite expensive; the resulting costs may obfuscate the returns of the market maker. With other asset classes, one could concentrate on optimizing the execution thus reducing transaction costs, but in the CDS index markets this is not the case as execution costs and impact are known. So, the only way to reduce transaction costs is to minimize the transaction amount, ideally without increasing the risk related to an open delta exposure. Similar types of behavior can be found in other OTC instruments such as interest rate swaptions, so while focusing on a specific instrument, our approach remains general.
Finally, this approach is certainly interesting for XVA traders and more generically bank resource managers, who typically have to deal with portfolios with hybrid and convex risks, and experience high rebalancing costs.

The contribution of this paper is a robust instrument capable of giving the trader a hedging signal, or even capable of autonomously trading if fitted with a market access, which is more accurate than the B&S delta which is currently widely used, as it is optimized in discrete time and with transaction costs. Such an instrument can be created through the use of Reinforcement Learning (RL), specifically by applying TRVO <cit.>, an algorithm capable of optimizing together the hedging (i.e. risk reduction) and Profit and Loss (p&l) objectives. By controlling the risk adverseness parameter, we are capable of creating a frontier, thus the job of the trader can be reduced to simply deciding on which point of the frontier to place himself.

To our knowledge, this is the first time the problem of hedging credit index options is analyzed from a RL perspective and the first time this approach is tested using real data for the underlying instrument.

Related Works

The issue of delta hedging using RL has been analyzed by various authors. Among the most recent approaches we mention <cit.>.
These papers can be subdivided into two categories, one addresses the problem from a practitioner's perspective and is focused on the details of the hedging strategies chosen by the agent; the other builds on the formal mathematical structure of option pricing and uses machine learning techniques to overcome the problems posed by realistic features such as transaction costs.
The distinction is faint as a hedging strategy implies a price, and vice versa.

The first category includes <cit.> and is also pertinent for this chapter.
The most comparable, regarding the financial environment, are <cit.>, which use the same MDP formulation considered in this dissertation.
The main difference consists in the use of an approximate variance formulation in the RL objective, compared to the full variance used in this chapter.
Furthermore, <cit.>, uses a one-step SARSA update, a value based approach, instead of a policy search method, while <cit.> considers both DQN <cit.> and PPO <cit.>.

<cit.> also consider an environment very close to ours, but with a transaction costs size that is larger than what we considered.
Regarding the RL algorithm, they use value-function methods and, in particular, risk-averse deep Q-learning.
It is an advanced approach taken from the risk-averse reinforcement learning literature <cit.>.
They consider two Q functions, one for the first moment and another for the second moment.

The paper then focuses on the agent's efficiency as a function of the rebalancing frequency.
Differently to this chapter where we also analyze what happens when changing the risk-aversion parameter, in <cit.>, only a single value of risk aversion is tested.

The second category includes <cit.>.
In <cit.>, the problem of option pricing in discrete time has been addressed from a machine learning perspective, neglecting hedging costs.
In <cit.>, the option pricing problem is undertaken by considering a class of convex risk measures and embedding them in a deep neural network environment. Initially, the dependence of the option price and hedge on the risk aversion parameter is studied in the absence of transaction costs. Then, a study of the option price dependence on transaction costs is discussed and the functional dependence of the price on the cost parameter is reconstructed.

What distinguishes our approach is the algorithm we considered: the risk-averse policy search algorithm TRVO. One of the advantages of TRVO compared to value based algorithms like the ones used by <cit.> is the fact that being policy search, TRVO is natively compatible with continuous states and actions and thus does not suffer from the problems of using a function approximator. Furthermore, being risk-averse, it is not necessary to apply any transformation to the reward differently from what is done for example in <cit.> and it is able to create a policy specific on the risk aversion of the user. Moreover, an advantage of model free RL algorithms, is that the policy learned is independent from the model used to generate the data. Thus TRVO can be used as is in an option hedging framework, and only requires the standard hyperparameter tuning typical of RL algorithms.

Paper Outline
In Section <ref> we present in detail the financial framework, describing CDS indexes and options before explaining how the hedging problem can be described using Markov Decision Processes (MDPs). In Section <ref> we will describe which reinforcement learning algorithm was used and for which reasons. Finally in Section <ref> we will see and evaluate the experimental performance.

§ FINANCIAL ENVIRONMENT

In this section we introduce the relevant financial instruments, the Markit iTraxx Europe Senior Financial index and the credit index options built on it. We show how to price them in a standard Black & Scholes environment and how options can be managed via standard delta hedging.

We will consider a financial environment where interest rates are set to zero for simplicity. Since both the options and the underlying are derivatives, trading them does not attract significant cash needs, and we can assume a substantial decoupling between rates and credit without loss of generality.

§.§ Markit iTraxx Europe Senior Financial index

The Markit iTraxx Europe Senior Financial index is a basket of credit default swaps on 30 european financial institutions (banks and insurances), equally weighted, with standardized maturities, coupons and payment dates.

Every 6 months, on the 20/09 and 20/03, or the Business Day immediately thereafter if it is not a Business Day, a new Series of the index is originated (or "rolled"). The new Series will be called "on-the-run", until a new one is generated.  Different maturities are traded for this CDS index (3, 5 and 10 years) with the maturity date that is the 20/12 or 20/06, respectively. For our purposes we consider the CDS index with 5Y maturity because it is the most liquid and there are many more options compared to the other maturities.
The index composition may be different from one Series to the other either in the number of constituents [For index versions originated before March 2015 the number of constituents was 25.] or in the CDS reference entities considered.
At the present time the Markit iTraxx Europe Senior Financial index on-the-run is the Series 35, started on 22 March 2021 and with maturity date 20 June 2026.

Each CDS index has a premium leg and a protection leg. The premium leg has standardized coupon dates: 20/03, 20/06, 20/09 and 20/12 (or the Business Day immediately thereafter if it is not a Business Day).
The coupon equals $\text{N}\,1\%\,\tau(t_{i-1},t_i) $, where N is the notional, expressed in Euro, $\tau(t_{i-1},t_i)$ is the year fraction, equal to the number of days between the present $t_i$ and the previous $t_{i-1}$ coupon date, divided by 360, while 1% is the standardized coupon; the coupon is paid on $t_i$.[The only caveat is about the last coupon date, which corresponds to the index maturity equal to the 20/06 or 20/12 even in case that day is a holiday, with a year fraction including an extra day.]
The protection leg pays, in case of a default of one of the $j$-th Series constituents occurred before the Series' maturity, an amount equal to $\text{LGD}_j\, \text{N} \, 1/n_j$, where $\text{LGD}_j$ is the loss given default[$\text{LGD}_j$ is equal to $1-R_j$, where $R_j$ is the recovery rate determined at the end of the ISDA CDS auction triggered by the credit event.] and $1/n_j$ is the constituent weight with $n_j$ the number of constituents at the default time (at the first default $n_1=30$).
Upon default of a constituent and settlement of the relative protection leg a new version of the Series is spin-offed including the surviving constituents, the notional N is rescaled accordingly.

Since the premium leg has a standardized $1\%$ coupon, the two legs are unbalanced by an amount that is exchanged at inception as a premium, this is referred to as upfront.

Even though the upfront amount is precisely the price of the derivative, the market does not quote it directly. Rather, following the standard single name CDS convention, what is traded is the running coupon of a par (i.e. upfront equal to zero) CDS. The relation between the traded spread $S$ and the upfront, assuming the latter to be received by the protection buyer from the protection seller, is:

\begin{equation}
\label{eq:upfront}
\text{Upf}(t,S_t) = \left(1\%-S_t\right)A_S(t)+1\% \tau(t_{acc},t),
\end{equation}

where t is the evaluation date, $\tau(t_{acc},t)$ is the year fraction, $t_{acc}$ the coupon date immediately before $t$, and $A_S(t)$ the annuity at time $t$.[In the computation of the accrual term the year fraction is modified adding an extra day.]
The latter quantity is defined as:

\begin{equation}
\label{eq:annuity}
A_S(t) = \sum_{t^+< \left\{t_i\right\} \le t_n} \tau(\max(t_{i-1}, t),t_{i}) \frac{P_S(t,t_{i-1})+P_S(t,t_i)}{2},
\end{equation}

where $t^+ = t + 1$ day, $\{t_i\}$ is the strip of index coupon dates, $t_n$ is the index maturity, $P_S(t,\theta)$ the survival probability between the present time $t$, and any future time $\theta$, given the current credit spread $S=S_t$ ($A_S(t)$ does not depend on $S_t$ directly but through $P_S(t)$).[Notice that if $t$ is the day before a coupon date, this coupon is excluded from the strip.]
The survival probability can be approximated as in <cit.>:

\begin{equation}
P_S(t,\theta) = e^{-S_t\tau(t,\theta) \text{LGD}^{-1}},
\end{equation}

with LGD usually set to 60$\%$ by convention.
Making trading decisions based on the credit spread is convenient as the upfront amount has jumps at the coupon dates due to $\tau(t_{acc},t)$, while the credit spread maintains a smoother behavior.

In the following we will consider the traded spread $S$ as a sort of underlying, having it's own specific dynamics, which we will simulate with Geometric Brownian Motion (GBM). The dynamics of the index will be inherited by the dynamics of $S$.
Thus, let $S_t$ be the underlying at time $t$, then it can be described as:

\begin{equation}\label{eq:GBM}
\mathrm{d}S_t = \mu S_t \mathrm{d}t + \sigma S_t \mathrm{d}W_t
\end{equation}

where $W_t$ is Brownian motion, $\mu$ the drift (which we assume to be 0 throughout the paper without loss of generality) and $\sigma$ the volatility.

For an initial value $S(0)$, the SDE has the analytic solution:

\begin{equation} \label{eq:spread_definition}
S_t=S_0\exp\left(\left(\mu-\frac{\sigma^2}{2}\right)t+\sigma W_t\right)
\end{equation}

where: $W_{t+u}-W_t \sim N(0,u) = N(0,1) \sqrt{u}$.

§.§ Options on the CDS index

In this section, we consider options on the CDS index.
A receiver option gives the buyer the possibility of selling protection on the index at the expiry date at a spread equal to the strike. Conversely, a payer option gives the buyer the choice of buying protection at the expiry date at a spread equal to the strike.
Upon exercise in case of a payer (receiver) option, the option seller (buyer) physically delivers the underlying.
In terms of the strike K and the traded spread $S(T)$ at expiry, the payoff at expiry is:

\begin{align}
\label{eq:opt_payoff_c1}
{\max \,((S_T\, A_S(T) - K\,A_K(T)), \,0)}\\
\label{eq:opt_payoff_c2}
{\max \,((K\,A_K(T)-S_T\, A_S(T)), \,0)}
\end{align}

respectively for a payer (Pay) and a receiver (Rec) option and where $ A_K(T)$ is the same expression as $ A_S(T)$ that considers $S_t=K$ in $P_S(t)$. In this paper, for simplicity we consider the payoff

\begin{align}
\label{eq:opt_payoff}
{\max \,((S_T-K)\, A_S(T), \,0)},\\
{\max \,((K-S_T)\, A_S(T), \,0)},
\end{align}

which allows a treatment à la Black & Scholes on $S_t$[We focus on this simplification since the extension to the payoff of Equation (<ref>) and (<ref>), which is trivial from a numerical/RL perspective, complicates the analytical treatment in a way beyond our interest.], since the payoff of Equation (<ref>) can be seen as a call on the underlying $S_t$.

Considering an option traded at time $t$ with expiry $T$ and strike $K$, where for ease of notation we may write $S_t$ instead of $S(t)$:

\begin{align}
\label{eq:opt_pricing payer}
& \text{Pay}(t,S_t) =\left[\Phi(d_t)S_t(T) - \Phi(e_t)K\right] A_S(T), \\
&\text{Rec}(t,S_t)=\left[\Phi(-e_t)K- \Phi(-d_t)S_t(T)\right] A_S(T), \label{eq:opt_pricing receiver}\\
& d_t = \frac{1}{\sigma\sqrt{{\tau}(t,T)}}\left[\log\left(\frac{S_t(T)}{K}\right)+\left(\frac{\sigma^2}{2}\right){\tau}(t,T)\right], \nonumber\\
&e_t =d_t-\sigma\sqrt{{\tau}(t,T)}, \nonumber
\end{align}

where $S_T(t)$ is the forward value of $S_t$, $\sigma$ is the volatility and $\bar{\tau}(t,T)$, the number of days between $t$ and $T$ divided by 365[ACT/365 convention.][$T$ for the annuity is the settlement date, $T$ for $d_1$, $d_2$ is the expiry date.].

Of course, as is common for options, Equation (<ref>) and Equation (<ref>) are a way of mapping the option price into a volatility surface, which is convenient since the latter is a smoother function of the expiries and the strikes than the price is.

When trading a payer option, the buyer pays the option premium upfront to the seller, which delivers the underlying in case of exercise. There are no extra payments.

In case a name in the index defaults between $t$ and $T$, the option doesn't knock out i.e. there is not an automatic close out and, if exercised, it delivers both the protection leg of the defaulted name and the spin-offed index. In this sense the underlying of the option remains unchanged even if a default occurs, so that the relation between the underlying and the option is default-neutral.
Hence, one can neglect, as we will do, jump-to-default effects in modeling the underlying and option dynamics.

Finally, since the buyer of a payer (receiver) index option receives (pays) protection substantially from trading time $t$ and not from expiry $T$, the option price needs to be adjusted consequently in order to consider any losses due to the default before $T$. This is done by a proper adjustment to the forward spread.

Assuming zero interest rates for simplicity, the adjusted forward $S_T(t)$ is

\begin{align}
\label{eq:opt_fwd_adj_definition_1}
S_t(T) & = S_t +  \text{LGD}(1-P_S(t,T))\frac{1}{A_S(T)}.
\end{align}

In the limit  $S(t)\rightarrow 0$, assuming the option is traded at time $t$

\begin{align}
\label{eq:opt_fwd_adj_definition_2}
S_T(t) & \sim  S_t\left(1+\frac{\tau(t,T)}{\tau(T,t_n)}\right).
\end{align}

§.§ Trading the Index

The evolution of the iTraxx Europe Financial Senior 5y mid spread (in basis points, left axis) and bid/ask (in basis points, right axis) from mid 2020 to mid 2021.

The iTraxx Europe Senior Financial index is traded on the so-called "over-the-counter" (OTC) market. One of the most important differences with regulated exchanges is the trade execution; in regulated markets anonymous orders are sorted and matched through an order book managed directly by the exchange. In OTC markets, CDS indexes like the Senior Financial are traded against dealers through a Multilateral Trade Facility (MTF). Dealers contribute continuously a bid and an ask price for a given notional, ranging from 10 to up to 400 mln Eur, with most of the contributions ranging between 50 and 200 mln Eur. Trading times are not strictly regulated but generally go between 9am CET and 18pm CET.
Each contributor has a typical bid/ask span, which depends on the market conditions and on the level of the index spread (a wide spread is typically correlated with a larger bid-ask).
The market spread published by the dealers can be applied, but acceptance from the dealer is not always ensured: some dealers ensure that the spread level will be always confirmed, some retain the right to review it. Moreover, dealers can also decide to publish bid and ask spreads which cannot be executed at all, and are often off market. This is another difference from the market makers' quotes in the regulated markets, which are binding. Another consequence is that it is difficult to define an order book for OTC traded objects, to distinguish if a quote is applicable or not and what would have been the maximum size executable.

Defining trading costs
We approach the trading costs problem from a statistical perspective. The starting point is a dataset containing the most recent bid and ask spread quoted by all the dealers (about 20 in the dataset) every 30 minutes, during the most liquid trading hours (9:30am CET and 17:30pm CET).

In order to use this data it is necessary to discard quotes that are most likely typos, not executable or technological problems. Thus for each time-step and for both the bid and ask we consider the mean and standard deviation of the quotes of all the dealers and discard from the set the spreads which differ from the mean by more than two standard deviations.

Considering the processed data, we define as applicable bid the average of the remaining bid spreads, and as applicable ask the average of the remaining ask spreads. Finally we obtain the mid spread as the average of the applicable bid and applicable ask, and the bid/ask spread as the difference between them.

An alternative approach we considered to calculate trading costs was to use the median of the unfiltered bid and ask quotes. The resulting bid/ask spreads did not differ significantly from the first method, which we ultimately considered robust enough for our purpose.

In Figure <ref> we show the mid and bid/ask spreads from the cleaned dataset, on a one-year time horizon, considering intra-day data with 30-minutes time-steps. In the rest of the paper, we identify the mid with the spread $S$ introduced in Equation <ref>. The unit of measurement of bid-ask spread and mid spread is the Basis Point (bp) i.e. $1\text{bp} = \frac{1}{10000}$.

The way in which applicable bid and ask spreads are built ensures that notionals up to hundreds of millions of the index can be traded at that spread, so that we can discard execution-related issues, such as slippage and assume that the trading/rebalancing costs can be computed from the bid/ask shown in Figure <ref>, as a linear functional of the traded notional N:

\begin{equation} \label{eq:costdef}
c(\text{N})= \text{N} \left|\text{Upf}_t\left(S_t\pm\frac{\text{ba}}{2}\right)-\text{Upf}_t(S_t)\right|,
\end{equation}

where $ba$ is the bid/ask, Upf$(S,t)$ the index upfront at the execution time $t$, as per Equation <ref>. In Equation <ref> the $+$ ($-$) sign should be considered when buying (selling) protection.

§.§ Embedding in a Markov Decision Process

The hedging problem is a sequential decision problem, where the trader needs to decide at each time-step how much of the underlying instrument to trade based on information coming from the market. This sequential decision problem can be described through a Markov Decision Process (MDP).
The p&l in one timestep of a trader long a payer option and holding $h(t)$ in the hedging portfolio is:

\begin{align}
p\&l & = \text{Pay}_{t+1}-\text{Pay}_t   \label{eq:delta_pl}\\  &- h_t\cdot(\text{Upf}_{t+1} - \text{Upf}_t) - c(a_t - a_{t-1}) \nonumber
\end{align}

We can define the delta hedge as:

\begin{equation*}\label{eq:CDS_delta}
N_h(t) = \left(\frac{\partial P_S(t,T)}{\partial S}\right)\left(\frac{\partial \text{ Upf}(t,S_t)}{\partial S}\right)^{-1}.
\end{equation*}

The B&S model assures that $p\&l \rightarrow 0$ when $p\rightarrow 0$, $h(t) = N_h(t)$ and there are no transaction costs ($c(N)=0$).

From now on, we consider p>0, in particular we take as a reference point 17 rebalances per day and c(N) as defined in Equation (<ref>).

We can now transition to a reinforcement learning scenario, which will be rigorously defined in the next section. We shall define this hedging environment as a sequential decision problem, specifically as a Markov Decision Problem (MDP):

* the action $a_t=h_t \in [0,1]$

* the state $s_t = (S_t,\text{Pay}_t,N_h(t),a_{t-1})$

* the reward is Equation (<ref>)

The above formulation is similar to what is used in <cit.>, and is called accounting P&L formulation in <cit.>.

§ REINFORCEMENT LEARNING

In this section, we give a brief introduction of reinforcement learning, focusing on the algorithm which we based our analysis on.

A discrete-time Markov Decision Process (MDP) is defined as a tuple $\langle\Sspace,\Aspace, \mathcal{P}, \mathcal{R}, \gamma, \mu\rangle$, where $\Sspace$ is the state space, $\Aspace$ the (continuous) action space, $\mathcal{P}(\cdot|s,a)$ is a Markovian transition model that assigns to each state-action pair $(s,a)$ the probability of reaching the next state $s'$, $\mathcal{R}(s,a)$ is a bounded reward function, $\gamma\in[0,1)$ is the discount factor, and $\mu$ is the distribution of the initial state. The policy of an agent is characterized by $\pi(\cdot|s)$, which assigns to each state $s$ an action with a certain probability.

We consider infinite-horizon problems in which future rewards are exponentially discounted with $\gamma$.  Following a trajectory $\tau \coloneqq (s_0, a_0, s_1, a_1, s_2, a_2, ...$), let the returns be defined as the discounted cumulative reward:
$G = \sum_{t=0}^\infty \gamma^t \mathcal{R}(s_t,a_t).$
For each state $s$ and action $a$, the action-value function is defined as:

\begin{equation}
Q_\pi(s,a) \coloneqq \EV_{\substack{s_{t+1}\sim \mathcal{P}(\cdot|s_{t},a_{t})\\a_{t+1}\sim\pi(\cdot|s_{t+1})}}\left[\sum_{t=0}^\infty \gamma^t \mathcal{R}(s_t,a_t)|s_0 = s, a_0 = a\right], %\equiv \EV_{s_t, a_t \sim\tau}\big[\sum_{t=0}^\infty \gamma^t \mathcal{R}(s_t,a_t)|s_0 = s, a_0 = a\big],
\label{eq:Q_fun}
\end{equation}

The typical RL objective is to maximize the action value function, given the initial state distribution.
This objective can be maximized in two main ways, the first is by learning the action-value function for each state and action, in general using the Bellman Equation. Once the action-value function is known, the policy is: $\pi(a|s) = \mathrm{argmax}_{a} Q(s,a)$, these algorithms are called value-based <cit.>, and are used in <cit.>. This approach becomes cumbersome in a hedging environment where both states and actions are (almost) continuous. There are approaches which use function approximation to interpolate the action value function, but it is necessary to discretize the state-action space loosing precision.

The other family is instead policy search methods <cit.>, which optimize the objective by searching directly in the policy space. They can easily handle continuous actions, learn stochastic policies in partially observable, non-Markovian environments and are robust when working in datasets with large amounts of noise <cit.>.

For all these reasons, we focused on policy search algorithms <cit.>.

Risk Averse Reinforcement Learning

The typical objective is to maximize the expected cumulative rewards, which in our context means maximizing the expected cumulative p&l. But maximizing this quantity is not the correct objective for this type of problem, in fact in an ideal B&S model, this quantity is as close as possible to zero, which translates to optimizing a risk averse objective.

Given the great experimental results achieved in <cit.>, we decided to use the Trust Region Volatility Optimization (TRVO) algorithm defined in in <cit.>.

The risk averse objective is:
$ \eta = J+\lambda \nu^2$,  where:

\begin{equation}
\nu^2_\pi \coloneqq (1-\gamma)\EV\left[\sum_{t=0}^\infty \gamma^t \left(\mathcal{R}(s_t,a_t)-J_\pi\right)^2\right].
\label{eq:reward_vola}
\end{equation}

One of the interesting things of this risk metric, to which we will refer to as reward-volatility is the it bounds the return-variance <cit.>. We would like to bring the reader's attention to the meaning of this reward-volatility term: it is minimizing the variations between one step and the next, in contrast to the return-variance which is minimizing the variance at the end of each path.

In this paper we aim at training agents with different risk aversions, in order to find target balances between risk (volatility) and reward. In a static environment, this can be achieved by training each agent with a specific value of $\lambda$ and algorithm will find a minimum with a specific, $\lambda$-dependent risk-reward ratio.
Instead, in an evolving environment with variable bid ask spread, a given specific value for $\lambda$ may induce different risk-reward targets, due to the fact that the terms in the risk averse objective $\nu$ will change in value even if the market conditions remain equal (i.e. the same action will induce different transaction costs).
This is the case in our problem as we can see from Figure <ref>, where the bid/ask varies significantly.
Intuitively, we can see that the dependence of the $J$ on $ba$ will be at most linear, and typically sublinear. Thus, distortions could come from a different scaling for the variance term $\nu^2$, but as will be apparent from the experiments, also in this case the distortion is sublinear, so there is no need to implement modifications/rescaling to take the issue into account.

The hedging strategy chosen by agents trained at different value of  the risk aversion parameter $\lambda$ is compared with the delta hedging strategy in a zero-cost environment. In the vertical axis the hedging notional as a percentage of the option notional.

§ EXPERIMENTS

In this section we present the experimental results. Once described the data generation and training parameters, we can will see the results obtained on a GBM simulated market, a heston simulated market and finally also real market data.

§.§ Data generation and agent training

We trained our agents on generated data, with episodes of 40 working days, with 17 observations per day, beginning at 9.30 and ending at 17.30. We simulated only the traded spread $S$, by using the GBM described in Equation (<ref>) with $\sigma$, the annualized volatility, equal to 60% and neglecting the drift term. We did not consider the possibility of a default of one of the components as no default has been observed in recent times for the instrument in consideration.
In each simulation, the underlying spread starts from an initial value of 100 bps; we define the stochastic evolution on the actual time span between the time-steps: 30 minutes during the day, 16 hours between the last step of one trading day and the first step of the next trading day in case of two contiguous trading days, a span of $16 + 24 n$ hours in the case of trading days separated by $n$ holidays or weekend days.

We trained our agents to hedge a position short a payer (but any other position would have been equivalent) option with 2 months maturity, thus maturing at the end of each episode. The strike $K$ was 100bps, equal to the initial value of the underlying at the begging of the episode.  We assumed an option notional of 100 mln Eur, which implies an hedging portfolio containing an underlying notional between 0 and 100 mln Eur. Given the market structure, our results are valid even assuming an option notional 10 times larger . We also assumed continuous underlying trading, which is reasonable given the option size and the fact that in the market small clips (down to 100KEur or less) can be traded. Assuming a risk neutral volatility equal to 60% the option has initial value of 530KEur.

The hedging strategy chosen by agents trained at different value of  the risk aversion parameter $\lambda$ is compared with the delta hedging strategy in an environment including hedging costs. In the vertical axis the hedging notional as a percentage of the option notional. In the upper (lower) plot the bid/ask equals 0.5 (2) basis points.

We built a training set of 40,000 episodes, and trained our agents varying two parameters: the risk aversion parameter $\lambda$ and the bid/ask parameter $ba$. We considered $\lambda $ following <cit.>, in order to span an efficient frontier in the risk/reward space; we considered $10^{-6}\lesssim\lambda\lesssim10^{-3}$, for better interpretability, in the rest of the paper we will rescale $\lambda$ by $10^5$, so to have bounds between 0.1 and 100. The choice of $ba$ as an extra parameter comes from the observation that the bid/ask of the instrument considered here shows a highly dynamic pattern (see Figure <ref>).

We considered ba ranging from 0.5 to 2 basis points (bp) as per Figure <ref>. We also considered the case with low values of $ba$, even $ba=0$ in order to further test our algorithms and to check that the standard delta-hedging strategy is smoothly recovered in the limit $ba\rightarrow 0$.

Each dot represents the performance of an agent on a GBM-simulated market in terms of p&l (w.r.t. delta hedging) and p&l volatility, depending on $\lambda $ (annotated next to each dot) and the $ba$ parameter.

§.§ Testing on a GBM-simulated market

We tested our agents on a dataset of 2,000 episodes with the underlying spreads generated by the GBM, with the same parameters of the training dataset. We performed different tests varying the $ba$ spread in order to monitor agents' performances comparing to the delta-hedging strategy.

In the $ba=0$ case, the trained agent perfectly replicates the delta-hedge, this can be seen in Figure <ref> where, for a specific testing scenario, the delta hedging strategy (in red) is compared with the action chosen by the agents trained with different values of the risk aversion parameter (in green, blue and purple). Given the absence of trading costs, all the agents replicate the same strategy, which is the optimal one[Indeed, neglecting hedging costs, the Black&Scholes paradigm is violated only by the assumed time discretization.], minimizing risks. Under the $ba=0$ assumption, the strategy has in average zero cumulated p&l.

Introducing hedging costs $ba>0$ the average cumulated p%l of the delta hedging strategy is shifted to negative values, depending linearly on $ba$, specifically, considering a $ba$ of 1 bp the cumulated p&l is on average -136kEur. The presence of hedging costs during training induces a smoother strategy for the agent, in terms of underlying allocation changes. Since each action becomes more expensive as $ba$ increases, the agent cuts costs through the reduction of portfolio rebalances. The downside of this approach is an increase the variability of the rewards, since the option is not continuously hedged. The desired balance between cost reduction and low reward volatility can be achieved changing the lambda aversion parameter of the model. This relationship is plotted in Figure <ref>, where different degrees of smoothness in the variation of the hedging portfolio can be seen to be dependent on $\lambda$.
The smoothness degree depends also on the size of the hedging cost: defining a certain risk aversion, a higher $ba$ implies a higher smoothness, as it is apparent by comparing the upper and lower plot.

The distribution of the p&l of the $\lambda=4$ agent relative to the p&l of delta hedge assuming a the $ba$ parameter equal to 1.5 bps.

The performance of the agents w.r.t. the delta hedging strategy in terms of cumulated p&l for different values of $\lambda$ and the $ba$ parameter is summarized in Figure <ref>.
In the figure, each dot represents the performance of an agent having the $\lambda$ indicated by the nearby annotated number and acting in an environment with $ba$ depending on the color (red for 0.5 basis point, orange for 1 basis point etc.).
The position on the vertical axis indicates the average p&l performance of the agent w.r.t. to the delta hedging strategy in an environment having the same $ba$.
The average is taken with respect to the terminal p&l measured on the 2,000 testing scenarios.
The position on the horizontal axis, instead, indicates the square root of the variance of the terminal p&l (the p&l volatility) on the same testing sample.
The colored dots laying on the horizontal axis indicate the variance performance of the delta hedging strategy in terms of p&l volatility at different values of the $ba$ parameter.
As se can see, all the agents perform better than the corresponding delta-hedging strategies in terms of p&l, while a certain number of agents (those lying left of the corresponding colored vertical line) perform better than the delta hedging strategy also in terms of p&l volatility. In this sense, all the frontiers dominate the corresponding delta hedge, and it is striking to notice that the level of dominance depends on the $ba$ parameter: at low costs the dominance is mild (as it was also experienced in <cit.>, where the very low hedging costs of listed equity products have been considered), at high costs the delta hedging is barely reasonable a strategy. As an example, one can consider the $\lambda =4$ point of the blue frontier (which assumes very large costs and beats delta hedging both in terms of p&l and p&l volatility) and observe from Figure <ref> how smooth its action is.
Another thing to notice is the $\lambda$ parametrization of the different frontiers, there is a shift of $\lambda$ to the left at the increase of the $ba$ parameter. This $\lambda$-scaling in $ba$, which is very mild, is in agreement with the considerations made at the end of Section <ref>.
The benefit of adopting our approach instead of the delta hedging strategy is apparent also from Figure <ref>, where we show the distribution of the p&l of the $\lambda =4$ agent relative to the p&l of delta hedging in the realistic case of $ba = 1.5$ bp: the agent essentially performs always better.

Eeach dot represents the performance on a Heston-simulated market of an agent in terms of p&l (w.r.t. delta hedging) and p&l volatility, depending on $\lambda $ (the number close to the dot) and the $ba$ parameter.

§.§ Testing on a Heston-simulated market

In order to make a further step towards realism, we challenge the assumption of the GBM constant volatility, as we know it does not hold in the financial markets. We thus generated a new testing set of 2,000 episodes with spreads derived from the Heston model, which introduces a dynamic for the volatility:

\begin{align}
dS_t & = \sqrt{\nu_t}\, S_t\, dW_t^S \\
d\nu_t & = \kappa\, (\theta - \nu_t)\,dt + \xi \sqrt{\nu_t}dW^\nu_t
\end{align}

with $\nu_0 = 60\%^2$, so to recover the initial volatility used in training, $\kappa = 2$, $\theta = \nu_0$, $\xi = 0.9$, and no correlation between the stochastic terms $dW_t^S$ and $dW^\nu_t$. With this configuration $\nu_t$ oscillates significantly reaching values as high as $\sim 120\%$ and as low as $\sim 0\%$.
When pricing the option we maintained the B&S formulation with $\sigma=60\%$.

Even if the agents were trained on a dataset generated with GBM, they are able to achieve very good performance over the heston dataset (see Figure <ref>). The reason could be that the hedging of an option is a task that implies a deep knowledge of the relationship between the underlying price and the option premium, but the way in which the underlying evolves is probably a secondary aspect.

The delta hedging strategy (red line) is compared to the strategy selected by agents trained with $ba=1$, at various risk aversions (other colored lines), on a real episode (the black line shows the underlying spread $S$ observed between July and September 2020.

The delta hedging strategy (red line) is compared to the strategy selected by agents trained with $ba=1$, at various risk aversions (other colored lines), on a real episode (the black line shows the underlying spread $S$ observed between September and November 2020.

§.§ Testing on real market data

In order to move a further step towards a realistic setup, we consider now real market data for iTraxx Europe Senior Financial index. We use the dataset constructed in Section <ref>, thus considering real market prices and real transaction parameters $ba$ as seen in Figure <ref>, and simulate the option price with $\sigma=60\%$. The available data is sufficient for 5 episodes of 40 days, which we used as a test set for agents trained with different values of $\lambda$ with $ba=1$ bp.

In Figures <ref> and <ref> we show the action of the the various agents compared with the delta hedge. We also show the market data dynamics (in black), on the right vertical axis. We can see as in the previous Figures, how lower values of $lambda$ generate smoother hedging policies.

Table <ref> summarizes the performance of the various agents (all the figures in kEur): in all the scenarios, all the considered agents
overperform the delta hedging strategy in terms of p&l.  Considering risk, given the low number of scenarios at hand, the cumulated p&l volatility previously considered is a very noisy estimator, thus, we considered the volatility of the p&l along each scenario, a measure similar to the reward volatility defined in Equation (<ref>) and used in  Section <ref> to define the objective of the TRVO algorithm, as described in <cit.>. Using this measure no agent outperforms delta hedging, but the volatility increase is anyway very small when compared with the cost reduction obtained by adopting our agents.

1c|     10c|Scenario

1c|     2c|1     2c|2     2c|3     2c|4     2c|5
$\lambda$     p&l     vol      p&l     vol      p&l     vol      p&l     vol      p&l     vol

1    -177   4.1   183   8.7   -257   3.8   -174   2.6   -282   2.4

1.25    -165   3.6   221   8.9   -224   3.6   -185   2.4   -272   2.4

1.5    -174   3.5   236   8.9   -212   3.6   -200   2.1   -270   2.4

2    -181   3.2   219   7.1   -188   3.5   -223   1.8   -285   2.3

3    -188   3.5   154   5.3   -144   3.6   -215   1.8   -294   2.3

4    -212   3.6   95   3.9   -129   3.7   -234   1.7   -308   2.4

5    -224   3.5   62   3.3   -134   3.4   -242   1.7   -317   2.4

6    -229   3.4   40   3.0   -138   3.3   -247   1.7   -323   2.4

7.5    -231   3.0   18   2.7   -145   3.0   -252   1.6   -329   2.4

10    -237   2.8   -15   2.2   -154   2.8   -259   1.6   -338   2.4

15    -264   2.4   -53   1.9   -191   2.7   -269   1.6   -349   2.3

25    -293   2.2   -80   1.7   -224   2.6   -278   1.6   -356   2.3

$\delta$     -376   2.1   -149   1.6   -310   2.5   -299   1.6   -372   2.3

The performance of the agents trained with $ba=1$ at various risk aversions on the real market data of Figure <ref>. On each scenario the p&l and path volatility of the delta hedging strategy (last line) is compared with those of the agents (the other lines). The agents always outperform delta hedging in terms of p&l, at the cost of an increase of the path volatility which is null or negligible for high values of $\lambda$ and it's very moderate in the other cases.

§ CONCLUSION

In this paper we tackled the credit index option hedging problem with the use of reinforcement learning. As we are in a dealer market scenario, there is no market impact, thus the only way to reduce costs when trading is by optimizing the trading policy. We showed that through the use of a state of the art RL algorithm TRVO, it is possible to learn a strategy which beats the practitioner's delta hedge in terms of risk, reward and generates lower transaction costs. This result was obtained not only on data generated through a GBM, but also when generating the underlying with a heston process and using real market data.
Interesting future works would be to consider a portfolio of options, or to complicate the financial environment by considering hybrid options.

, Sabbioni
, Vittori
, Papini

Bisi, L.
, Sabbioni, L.
, Vittori, E.
, Papini, M.
Restelli, M.

Risk-Averse Trust Region Optimization for
Reward-Volatility Reduction Risk-averse trust region optimization for
reward-volatility reduction.
arXiv preprint arXiv:1912.03193.

Black, F.
Scholes, M.

The pricing of options and corporate liabilities The
pricing of options and corporate liabilities.
Journal of political economy813637–654.

, Gonon
, Teichmann

Buehler, H.
, Gonon, L.
, Teichmann, J.
Wood, B.

Deep hedging Deep hedging.
Quantitative Finance1–21.

, Chen
, Hull

Cao, J.
, Chen, J.
, Hull, JC.
Poulos, Z.

Deep Hedging of Derivatives Using Reinforcement Learning
(December 20, 2019) Deep hedging of derivatives using reinforcement
learning (december 20, 2019).

, Neumann
, Peters

Deisenroth, MP.
, Neumann, G.
, Peters, J.
A survey on policy search for robotics A survey on
policy search for robotics.
Foundations and Trends® in

Du, J.
, Jin, M.
, Kolm, PN.
, Ritter, G.
, Wang, Y.
Zhang, B.

Deep Reinforcement Learning for Option Replication and
Hedging Deep reinforcement learning for option replication and
The Journal of Financial Data Science2444–57.

Guéant, O.
Pu, J.

Option pricing and hedging with execution costs and
market impact Option pricing and hedging with execution costs and market
Mathematical Finance273803–831.

Halperin, I.

QLBS: Q-learner in the Black-Scholes (-Merton) worlds
Qlbs: Q-learner in the black-scholes (-merton) worlds.
arXiv preprint arXiv:1712.04609.

Halperin, I.

The QLBS Q-Learner goes NuQLear: fitted Q iteration,
inverse RL, and option portfolios The qlbs q-learner goes nuqlear: fitted q
iteration, inverse rl, and option portfolios.
Quantitative Finance1–11.

Jarrow, RA.
Turnbull, SM.

Pricing derivatives on financial securities subject to
credit risk Pricing derivatives on financial securities subject to credit
The journal of finance50153–85.

Kolm, PN.
Ritter, G.

Dynamic replication and hedging: A reinforcement
learning approach Dynamic replication and hedging: A reinforcement learning
The Journal of Financial Data Science11159–171.

Leland, HE.

Option pricing and replication with transactions costs
Option pricing and replication with transactions costs.
The journal of finance4051283–1301.

Mnih, V.
, Kavukcuoglu, K.
, Silver, D.
, Graves, A.
, Antonoglou, I.
, Wierstra, D.
Riedmiller, M.

Playing atari with deep reinforcement learning Playing
atari with deep reinforcement learning.
arXiv preprint arXiv:1312.5602.

Moody, J.
Saffell, M.

Learning to trade via direct reinforcement Learning to
trade via direct reinforcement.
IEEE transactions on neural Networks124875–889.

Peters, J.
Schaal, S.

Reinforcement learning of motor skills with policy
gradients Reinforcement learning of motor skills with policy
Neural networks214682–697.

, Wolski
, Dhariwal
, Radford

Schulman, J.
, Wolski, F.
, Dhariwal, P.
, Radford, A.
Klimov, O.

Proximal Policy Optimization Algorithms Proximal
policy optimization algorithms.

Sutton, RS.
Barto, AG.

Introduction to Reinforcement Learning Introduction to
reinforcement learning (1st ).
Cambridge, MA, USAMIT Press.

, Di Castro

Tamar, A.
, Di Castro, D.
Mannor, S.

Learning the variance of the reward-to-go Learning the
variance of the reward-to-go.
The Journal of Machine Learning

, Trapletti

Vittori, E.
, Trapletti, M.
Restelli, M.

Option Hedging with Risk Averse Reinforcement Learning
Option hedging with risk averse reinforcement learning.
arXiv preprint arXiv:2010.12245.

Yalincak, OH.

Criticism of the Black-Scholes Model: But Why is It
Still Used?: (The Answer is Simpler than the Formula) Criticism of the
black-scholes model: But why is it still used?: (the answer is simpler than
the formula).
The Answer is Simpler than the Formula (July 22,