Некоторые быстрые методы работы с формулой Блэка-Шоулса

При торговле опционами весьма неплохо знать и понимать теорию Блэка и Шоулса. Можно, конечно, смотреть профили позиций и прочее на многочисленных специальных сервисах типа www.option.ru, но, как известно, хочешь сделать хорошо--сделай все сам. В применении к опционам это вполне правильная вещь--не стоит доверять сторонним сервисам. Не потому, что они плохи (они обычно вполне корректно все рассчитывают), а потому, что опционы надо чувствовать.

Краткая и лаконичная суть теории Блэка и Шоулса изложена здесь: http://anatoly-utkin.livejournal.com/2835.html . Ничего сложного в ней нет, это просто теория эффективного рынка в применении к опционам, не более. В настоящей заметке я хотел бы привести некоторые быстрые расчетные методы для работы с формулой Блэка-Шоулса, позволяющие быстро находить цены опционов и IV.


Итак, формула Блэка-Шоулса имеет вид: C=KN(d1)-SN(d2)  ( Wikipedia ). Первое, что тут есть из нетривиального--это функция N(x)--функция нормального распределения. В трейдерской тусовке модно аппроксимировать N(x) полиномом, однако мне это режет глаз, поскольку при этом не выполнено экспоненциальное стремление N(x) к единице на плюс бесконечности и к нулю на минус бесконечности. Поэтому такая метода мне абсолютно не нравится.



Вообще, N(x)--это интеграл от простой и имеющейся во всех приличных языках программирования экспоненциальной функции exp(-x^2). Однако, влобовую брать интегралы через интегральные суммы--дело накладное в плане вычислительных мощностей. Для нормального распределения существует множество асимптотических формул. Одну из них--так называемый T алгоритм--я и использую. Математическое описание T алгоритма: http://algolist.manual.ru/maths/matstat/NormalDF/NormalDF1.php
Код VBA:

Function FERT(x As Single) As Single
Dim i%, n%
Dim S!, T!, s1!, x2!, epsilon!
epsilon = 0.00001
s1 = 0
T = Abs(x)
s1 = T
x2 = x ^ 2
n = 3
For i = 0 To 100
   T = x ^ 2 * T / n
   S = s1
   s1 = S + T
   n = n + 2
   If Abs(S — s1) < epsilon Then GoTo 10
   Next i
10   If x >= 0 Then FERT = 0.5 + s1 * NORM(x)
If x < 0 Then FERT = 0.5 — s1 * NORM(x)
'FERT = FERT — 0.5
If x < -10 Then FERT = 0
If x > 10 Then FERT = 1
End Function

Здесь NORM(x)=1/sqrt(2pi)*exp(-x^2/2)--плотность единичного нормального распределения, ее код:

Function NORM(x As Single) As Single
NORM = Exp(-x ^ 2 / 2) / (2 * 3.1415926) ^ (0.5)
End Function


Вторая нетривиальная и часто нужная вещь в формуле БШ--это как найти волатильность при известной цене опциона. Иными словами, нужно обратить формулу Блэка-Шоулса и найти в ней величину z=sigma*sqrt(T) по известным C, K, S. Для этого я использую метод Ньютона ( Wikipedia ). Введем вспомогательную функцию g(z)=KN(d1)-SN(d2)-C. Тогда искомая величина z будет решением уравнения g(z)=0. g(z)--это хорошая, аналитическая, достаточно плавная функция своего аргумента. Поэтому метод Ньютона для нее будет работать хорошо. Суть метода Ньютона неплохо поясняет картинка ниже:
















а его итерационная формула имеет простой вид: z_next=z-f(z)/(df/dz(z)).Особую прелесть именно методу Ньютона в применении к формуле Блэка-Шоулса придает тот факт, что производная от C по волатильности (так называемая вега) имеет простой вид S*NORM(d1)*sqrt(T). Код VBA:

Function IV_C(S As Single, K As Single, C As Single) As Single
Dim x!, x_next!, f!, df!, epsilon!
Dim i%, i_max%
epsilon = 0.00001
x = 0.5
i_max = 15
For i = 1 To i_max
   f = C — S * FERT((Log(S / K) + x ^ 2 / 2) / x) + K * FERT((Log(S / K) — x ^ 2 / 2) / x)
   df = -S * NORM((Log(S / K) + x ^ 2 / 2) / x)
   x_next = x — f / df
   If (Abs(x — x_next)) < epsilon Then GoTo 20
   x = x_next
   Next i
20 IV_C = x
End Function

Ну и напоследок общий финальный код VBA, объединяющий все эти вещи и позволяющий вертеть формулу Блэка-Шоулса во все стороны:

'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NORM(x As Single) As Single
NORM = Exp(-x ^ 2 / 2) / (2 * 3.1415926) ^ (0.5)
End Function
'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function FERT(x As Single) As Single
Dim i%, n%
Dim S!, T!, s1!, x2!, epsilon!
epsilon = 0.00001
s1 = 0
T = Abs(x)
s1 = T
x2 = x ^ 2
n = 3
For i = 0 To 100
   T = x ^ 2 * T / n
   S = s1
   s1 = S + T
   n = n + 2
   If Abs(S — s1) < epsilon Then GoTo 10
   Next i
10   If x >= 0 Then FERT = 0.5 + s1 * NORM(x)
If x < 0 Then FERT = 0.5 — s1 * NORM(x)
'FERT = FERT — 0.5
If x < -10 Then FERT = 0
If x > 10 Then FERT = 1
End Function

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function BS_C(K As Single, S As Single, T As Single, V As Single) As Single
Dim x!
x = V * (T / 365) ^ (0.5)
BS_C = S * FERT((Log(S / K) + (x ^ 2) / 2) / x) — K * FERT((Log(S / K) — (x ^ 2) / 2) / x)
End Function
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function IV_C(S As Single, K As Single, C As Single) As Single
Dim x!, x_next!, f!, df!, epsilon!
Dim i%, i_max%
epsilon = 0.00001
x = 0.5
i_max = 10
For i = 1 To i_max
   f = C — S * FERT((Log(S / K) + x ^ 2 / 2) / x) + K * FERT((Log(S / K) — x ^ 2 / 2) / x)
   df = -S * NORM((Log(S / K) + x ^ 2 / 2) / x)
   x_next = x — f / df
   If (Abs(x — x_next)) < epsilon Then GoTo 20
   x = x_next
   Next i
20 IV_C = x
End Function
'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub ggg()
MsgBox (BS_C(120000, 110000, 10, 0.76))
MsgBox (IV_C(110000, 120000, 2091.91) / (10 / 365) ^ (0.5))
End Sub
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''

Функция BS_C  рассчитывает цену опциона по известным страйку K, ЦБА S, волатильности V и времени до экспирации T, волатильность измерена в 1/корень из года, время до экспирации в днях.

Функция IV_C  рассчитывает величину «волатильность*корень из времени до экспирации» по известным страйку K, ЦБА S и цене опциона C.

Все это написано для коллов, для путов обобщение тривиально.