fzero not working with nested functions and mvncdf

Hello everyone,

I have defined three functions which are (progressively) nested

function E = EBlackCox(V0, sigma, VB, r, d, L, TL)

a = (r - d - sigma.^2./2)./(sigma.^2);
D = @(x, lambda) (log(V0./x) + lambda.*sigma.^2.*TL)./(sigma.*sqrt(TL));
E = V0.*exp(-d.*TL).*(normcdf(D(L, a + 1)) - (VB./V0).^(2.*(a + 1)) .* ...
    normcdf(-D(VB.^2./L, -a - 1))) - L.*exp(-r.*TL).*(normcdf(D(L, a)) ...
    - (VB./V0).^(2.*a) .* normcdf(-D(VB.^2./L, -a)));

end

function c0 = CBlackCox(V0, sigma, VB, r, d, L, TL, K, T)

a = (r - d - sigma.^2./2)./(sigma.^2);
f = @(x) EBlackCox(x, sigma, VB, r, d, L, TL - T) - K;
Kstar = fzero(f, [K 10.*(K + L./r)]);

D = @(x, lambda, t) (log(V0./x) + lambda.*sigma.^2.*t)./(sigma.*sqrt(t));

rho = sqrt(T./TL);
mu = [0 0]; 
Rho_p = [1 rho; rho 1];
Rho_m = [1 -rho; -rho 1];

c0 = V0 .* exp(-d.*TL) .*(mvncdf([D(Kstar, a + 1, T) D(L, a + 1, TL)], mu, Rho_p) ...
    + mvncdf([-D(VB.^2./Kstar, a + 1, T) D(L, a + 1, TL)], mu, Rho_m) ...
    - (VB./V0).^(2.*(a + 1)) .* (mvncdf([D(Kstar, - a - 1, T) -D(VB.^2./L, - a - 1, TL)], mu, Rho_m) ...
    + mvncdf([-D(VB.^2./Kstar, - a - 1, T) -D(VB.^2./L, - a - 1, TL)], mu, Rho_p)))...
    - L .* exp(-r.*TL) .*(mvncdf([D(Kstar, a, T) D(L, a, TL)], mu, Rho_p) ...
    + mvncdf([-D(VB.^2./Kstar, a, T) D(L, a, TL)], mu, Rho_m) ...
    - (VB./V0).^(2.*a) .* (mvncdf([D(Kstar, - a, T) -D(VB.^2./L, - a, TL)], mu, Rho_m) ...
    + mvncdf([-D(VB.^2./Kstar, - a, T) -D(VB.^2./L, - a, TL)], mu, Rho_p)))...
    - K .* exp(-r.*T) .* (normcdf(D(Kstar, a, T)) - (VB./V0).^(2.*a) .* normcdf(-D(VB.^2./Kstar, -a, T)));

end

function s =  sigmaCBlackCox(V0, VB, r, d, L, TL, K, T, price)

option = @(sigma) CBlackCox(V0, sigma, VB, r, d, L, TL, K, T);
difference = @(sigma) (option(sigma) - price);
s = fzero(difference, [.01 .8]);

end

The first two functions work fine when inputting values of the parameter. If you want to try here is a testable set

r = 0.03;
d = 0.06;
sigma = 0.2;
V0 = 100;
L = 40;
VB = 0.6*L;
TL = 10;
T = 1;
a = (r - d - sigma.^2./2)./(sigma.^2);
K = 30;

If applied to the second function

c0 = CBlackCox(V0, sigma, VB, r, d, L, TL, K, T) 

it returns c0 = 4.3352.

The problem is when I use the third function which tries to find sigma given c0. An error occurs

FZERO cannot continue because user-supplied function_handle ==> @(s)(option(s)-c0)
failed with the error below.

 Function values at the interval endpoints must differ in sign.

Error in fzero (line 245)
    fa = localFirstFcnEval(FunFcn,FunFcnIn,a,varargin{:});

However, the function works fine as long as I do something like

option = @(sigma) CBlackCox(V0, sigma, VB, r, d, L, TL, K, T);
difference = @(sigma) (option(sigma) - c0);
difference(0.1)
difference(0.2)
difference(0.3)
...

which are difference(0.1) = -2.7101, difference(0.2) = 0, difference(0.3) = 3.3517 respectively

The purpose of the third function is indeed to return 0.2

Also, I cannot fplot the function difference as I get the following warning

Warning: Error updating FunctionLine.

 The following error was reported evaluating the function in FunctionLine update:
 Operands to the logical AND (&&) and OR (||) operators must be convertible to logical
 scalar values. Use the ANY or ALL functions to reduce operands to logical scalar
 values.

and the figure appears blank.

I suspect this errors occurs because of the nesting toghether with using the function mvncdf (multivariate density of a normal vector) because I have a similar code (again three series of nested functions) which does not require the use of mvncdf and the fzero works propertly.

Can someone suggest a way to modify my code so that fzero works?

Thank you

Answer 1

[Torsten]

r = 0.03;
d = 0.06;
sigma = 0.2;
V0 = 100;
L = 40;
VB = 0.6*L;
TL = 10;
T = 1;
a = (r - d - sigma.^2./2)./(sigma.^2);
K = 30;
price = 4.3352;
[s,fval] = sigmaCBlackCox(V0, VB, r, d, L, TL, K, T, price)

s = 0.2360

fval = 0

function [s,fval] =  sigmaCBlackCox(V0, VB, r, d, L, TL, K, T, price)
option = @(sigma) CBlackCox(V0, sigma, VB, r, d, L, TL, K, T);
difference = @(sigma) arrayfun(@(x)(option(x) - price),sigma);
[s,fval] = fzero(difference, [.1 .8]);
end

function c0 = CBlackCox(V0, sigma, VB, r, d, L, TL, K, T)
a = (r - d - sigma.^2./2)./(sigma.^2);
f = @(x) EBlackCox(x, sigma, VB, r, d, L, TL - T) - K;
Kstar = fzero(f, [K 10.*(K + L./r)]);
D = @(x, lambda, t) (log(V0./x) + lambda.*sigma.^2.*t)./(sigma.*sqrt(t));
rho = sqrt(T./TL);
mu = [0 0] ;
Rho_p = [1 rho; rho 1];
Rho_m = [1 -rho; -rho 1];
c0 = V0 .* exp(-d.*TL) .*(mvncdf([D(Kstar, a + 1, T) D(L, a + 1, TL)], mu, Rho_p) ...
    + mvncdf([-D(VB.^2./Kstar, a + 1, T) D(L, a + 1, TL)], mu, Rho_m) ...
    - (VB./V0).^(2.*(a + 1)) .* (mvncdf([D(Kstar, - a - 1, T) -D(VB.^2./L, - a - 1, TL)], mu, Rho_m) ...
    + mvncdf([-D(VB.^2./Kstar, - a - 1, T) -D(VB.^2./L, - a - 1, TL)], mu, Rho_p)))...
    - L .* exp(-r.*TL) .*(mvncdf([D(Kstar, a, T) D(L, a, TL)], mu, Rho_p) ...
    + mvncdf([-D(VB.^2./Kstar, a, T) D(L, a, TL)], mu, Rho_m) ...
    - (VB./V0).^(2.*a) .* (mvncdf([D(Kstar, - a, T) -D(VB.^2./L, - a, TL)], mu, Rho_m) ...
    + mvncdf([-D(VB.^2./Kstar, - a, T) -D(VB.^2./L, - a, TL)], mu, Rho_p)))...
    - K .* exp(-r.*T) .* (normcdf(D(Kstar, a, T)) - (VB./V0).^(2.*a) .* normcdf(-D(VB.^2./Kstar, -a, T)));
end

function E = EBlackCox(V0, sigma, VB, r, d, L, TL)
a = (r - d - sigma.^2./2)./(sigma.^2);
D = @(x, lambda) (log(V0./x) + lambda.*sigma.^2.*TL)./(sigma.*sqrt(TL));
E = V0.*exp(-d.*TL).*(normcdf(D(L, a + 1)) - (VB./V0).^(2.*(a + 1)) .* ...
    normcdf(-D(VB.^2./L, -a - 1))) - L.*exp(-r.*TL).*(normcdf(D(L, a)) ...
    - (VB./V0).^(2.*a) .* normcdf(-D(VB.^2./L, -a)));
end

[Federico Maglione]

Thank you so much @Torsten! That works perfectly. Would you mind explaining me why an arrayfun is needed?

[Torsten]

Your function CBlackCox doesn't accept array inputs (inputs with more than one element at a time) for sigma. With arrayfun, it is called with only one value for sigma at a time.

Further, sigma = 0.01 does not work in CBlackCox ; I had to choose sigma = 0.1 as lower bound in the call to "fzero".

[Federico Maglione]

Thank you!

Answer 2

[Jan]

I tried to reproduce your problem, but get a different error message:

r     = 0.03;
d     = 0.06;
sigma = 0.2;
V0    = 100;
L     = 40;
VB    = 0.6*L;
TL    = 10;
T     = 1;
a     = (r - d - sigma.^2./2) ./ (sigma.^2);
K     = 30;

option     = @(sigma) CBlackCox(V0, sigma, VB, r, d, L, TL, K, T);
difference = @(sigma) (option(sigma) - price);
s          = fzero(difference, [0.01, 0.8]);

Error using fzero>localFirstFcnEval
FZERO cannot continue because user-supplied function_handle ==> @(sigma)(option(sigma)-price) failed with the error below.

Function values at interval endpoints must be finite and real.

Error in fzero (line 245)
fa = localFirstFcnEval(FunFcn,FunFcnIn,a,varargin{:});

function c0 = CBlackCox(V0, sigma, VB, r, d, L, TL, K, T)

a = (r - d - sigma.^2./2) ./ (sigma.^2);
f = @(x) EBlackCox(x, sigma, VB, r, d, L, TL - T) - K;
Kstar = fzero(f, [K, 10.*(K + L./r)]);

D = @(x, lambda, t) (log(V0./x) + lambda.*sigma.^2.*t) ./ (sigma.*sqrt(t));

rho = sqrt(T./TL);
mu = [0 0]; 
Rho_p = [1 rho; rho 1];
Rho_m = [1 -rho; -rho 1];

c0 = V0 .* exp(-d.*TL) .*(mvncdf([D(Kstar, a + 1, T), D(L, a + 1, TL)], mu, Rho_p) ...
    + mvncdf([-D(VB.^2./Kstar, a + 1, T), D(L, a + 1, TL)], mu, Rho_m) ...
    - (VB./V0).^(2.*(a + 1)) .* (mvncdf([D(Kstar, - a - 1, T), -D(VB.^2./L, - a - 1, TL)], mu, Rho_m) ...
    + mvncdf([-D(VB.^2./Kstar, - a - 1, T) -D(VB.^2./L, - a - 1, TL)], mu, Rho_p)))...
    - L .* exp(-r.*TL) .*(mvncdf([D(Kstar, a, T), D(L, a, TL)], mu, Rho_p) ...
    + mvncdf([-D(VB.^2./Kstar, a, T), D(L, a, TL)], mu, Rho_m) ...
    - (VB./V0).^(2.*a) .* (mvncdf([D(Kstar, - a, T), -D(VB.^2./L, - a, TL)], mu, Rho_m) ...
    + mvncdf([-D(VB.^2./Kstar, - a, T), -D(VB.^2./L, - a, TL)], mu, Rho_p)))...
    - K .* exp(-r.*T) .* (normcdf(D(Kstar, a, T)) - (VB./V0).^(2.*a) .* normcdf(-D(VB.^2./Kstar, -a, T)));
end

function E = EBlackCox(V0, sigma, VB, r, d, L, TL)
a = (r - d - sigma.^2./2)./(sigma.^2);
D = @(x, lambda) (log(V0./x) + lambda.*sigma.^2.*TL)./(sigma.*sqrt(TL));
E = V0.*exp(-d.*TL).*(normcdf(D(L, a + 1)) - (VB./V0).^(2.*(a + 1)) .* ...
    normcdf(-D(VB.^2./L, -a - 1))) - L.*exp(-r.*TL).*(normcdf(D(L, a)) ...
    - (VB./V0).^(2.*a) .* normcdf(-D(VB.^2./L, -a)));
end

[Federico Maglione]

I don't know (we are using the same version of matlab) and I still get this error

Error using fzero>localFirstFcnEval

FZERO cannot continue because user-supplied function_handle ==>

@(sigma)(option(sigma)-price) failed with the error below.

Function values at the interval endpoints must differ in sign.

Error in fzero (line 245)

fa = localFirstFcnEval(FunFcn,FunFcnIn,a,varargin{:});

[Federico Maglione]

Did you replace price in difference with c0 = 4.3352?

[Jan]

"Did you replace price in difference with c0 = 4.3352?" - I've used copy&paste with code from your question and did not replace anything.