# SABR Model for Normal Volatilities

The SABR model takes as input the parameters α, β, ρ and ν. Small time expansion yields an approximation formula for the implied normal volatility.
The general structure of the approximation formula is given as

kill(all);

sigma(S,K,T):=nu*(S-K)/chi(zeta)*(1+I*T);

Here S is the underlying forward, K is the option strike and T is the time to maturity.

## Denominator Computation

The denominator χ(ζ) is given by

chi(zeta):=log((sqrt(1-2*rho*%zeta+%zeta^2)-rho+%zeta)/(1-rho));

ζ is calculated based on the Lamperti transform

nu/alpha*integrate(1/C(s), s, K, S);

For the SABR model the local volatility function is

C(S):=S^beta;

This gives the representation of ζ (we give it explicitly to simplify Maxima evaluations) as

%zeta:nu / alpha*(S^(1-beta)-K^(1-beta)) /(1-beta);

As a result we get a representation of χ(ζ) in terms of input variables

chi(zeta);

We will use this representation later on to derive limiting cases of the formula.

## Computation of the Correction Term

The correction term I is given as

I:(2*%gamma_2(S_av)-%gamma_1(S_av)^2)/24*alpha^2*'C(S_av)^2+rho*nu*alpha*%gamma_1(S_av)/4*'C(S_av)+(2-3*rho^2)/24*nu^2;

The auxilliary parameters γ1 and γ2 are given by

define(%gamma_1(s),diff('C(s),s)/'C(s));

and

define(%gamma_2(s),diff('C(s),s,2)/'C(s));

The correction terms become

ev(%gamma_1(S), nouns, S=S_av);

and

ev(%gamma_2(S), nouns, S=S_av);

Consequently, the correction term I becomes

I2:ev(I, nouns);

The average Sav is calculated based on arithmetic average. That is

(S+K)/2;

This gives

ev(I2, nouns, S_av=(S+K)/2);

## ATM Volatility

ATM volatility requires special treatment because

chi2:chi(zeta);

ev(chi2,S=K);

We apply l'Hospital's rule manually to

(S-K)/'chi(zeta);

N:diff(S-K,S);

D:diff(chi(zeta),S);

D:ratsimp(ev(D,K=S));

This gives the ATM volatility formula

sigma_A(S,T):=nu*N/D*(1+'I*T);

ev(sigma_A(S,T),nouns);

## Normal SABR

Normal SABR formula is derived for β=0. Then we get

sigma_N(S,K,T):=nu*(S-K)/chi(zeta)*(1+I*T);

with χ(ζ) becomming

ev(chi(zeta),eval,beta=0);

and I becomming

ev(I,nouns,eval,beta=0);

## Log-Normal SABR

Log-normal SABR formula is derived for β→1. Again, we need to carefully take limit for χ(ζ)

chi_L:logcontract(tlimit(chi(zeta),beta,1,minus));

The correction term I simplifies to

I_L:ev(I,nouns,eval,beta=1);

## Boundary Case CEV

For ν=0 we get the boundary case of a CEV model approximation. Again, we need to apply l'Hospital's rule to

nu*(S-K)/'chi(zeta);

CEV1:factor(tlimit(nu*(S-K)/chi(zeta),nu,0,plus));

The correction term I becomes

CEV2:ev(I,nouns,nu=0);

This gives the representation of the CEV model normal volatility

sigma_C(S,K,T):=CEV1*(1+CEV2*T);

sigma_C(S,K,T);

## Boundary Case Normal Model

The normal model boundary case can be derived by setting β=0 in the CEV boundary case. That is

ratsimp(ev(sigma_C(S,K,T),eval,beta=0));

## Plotting

We plot example normal smile of the SABR model with ~100bp ATM volatility

sigma_0:0.01;
beta:0.5;
rho:0.1;
nu:0.05;
T:5;
S_0:0.05;
alpha:sigma_0/C(S);

SABR:ev(sigma(S,K,T), nouns,eval);

SABR:ev(SABR,S_av=(S+K)/2,S=S_0,nouns,eval);

CEV:ev(sigma_C(S,K,T),S_av=(S+K)/2,S=S_0,nouns,eval,eval);

We need to re-calibrate α for lognormal SABR...

alpha:sigma_0/S;

LSABR:ev(nu*(S-K)/chi_L*(1+I_L*T),S_av=(S+K)/2,S=S_0,nouns,eval);

plot2d([SABR,CEV,LSABR],[K,-0.01,0.09], [legend, "SABR", "CEV","LSABR"] ,[y,0,0.02])$