Message Boards Message Boards

0
|
4923 Views
|
1 Reply
|
1 Total Likes
View groups...
Share
Share this post:
GROUPS:

Optimization problem : FindMinimum function very slow

Posted 9 years ago

Dear all,

I am currently trying to minimize a function of 5 variables with FindMinimum. The latter function should minimize the square between a set of data and an other function that i have defined. The problem is that the minimization takes too much time. In the case provided bellow, it takes 4min22sec while I only have 2 observations. When I add more observations I do not have results after 30 minutes.

Does someone knows what could imply this run time ? What should I correct to my code to increase the performance?

Please find de code bellow (divided into two parts)

part one : I define some variables and a function that is called in the FindMinimum :

discountcurve[T_] := (0.98)^(T/2);
forwardcurve[
   T_] := (1/(12/12))*((discountcurve[(T - (12/12))]/
       discountcurve[T]) - 1);
instforwardrate[
   T_] := -(Log[discountcurve[T + 0.0001]] - Log[discountcurve[T]])/
   0.0001;

(******************************)
(*  Swaption pricing formula  *)
(******************************)
\
Clear[swaptionpricefunction]
swaptionpricefunction[kappa_?NumericQ, sigma_?NumericQ, eta_?NumericQ,
    alpha_?NumericQ, beta_?NumericQ, tenor_, swaptionmaturity_, 
   swapmaturity_, strike_, psi_] := Block[{},

   (* Expectation and Variance of the OIS (risk-
   free) short rate process under the T-forward measure *)

   oisprocessexpectation = (-((sigma^2)/(kappa^2))*(1 - 
         Exp[-kappa*swaptionmaturity])) + (((sigma^2)/(2*
           kappa^2))*(1 - Exp[-2*kappa*swaptionmaturity]));
   oisprocessvariance = ((sigma^2)/(2*kappa))*(1 - 
       Exp[-2*kappa*swaptionmaturity]);

   (* Density function of the OIS process *)

   oisprocessdensity[y_] := (1/((2*Pi*oisprocessvariance)^0.5))*
     Exp[-(((y + oisprocessexpectation)^2)/(2*oisprocessvariance))];

   (* Bond price formula in the Hull-White model *)

   bondb[t_, T_] := (1/kappa)*(Exp[-kappa*(T - t)] - 1);
   bonda[t_, T_] := -instforwardrate[t]*bondb[t, T] + 
     Log[discountcurve[T]/
       discountcurve[t]] + ((sigma^2)/4*
        kappa)*((bondb[t, T])^2)*(Exp[-2*kappa*t] - 1);
   bondprice[t_, T_, y_] := Exp[bonda[t, T] + bondb[t, T]*y]; 

   (* Daycount convetion definition *)
   daycountfloat = tenor/12;
   daycountfix = tenor/12;

   (* Variables definition  *)

   epsilon[T_] := 
    forwardcurve[T]*1.2- (1 + alpha)*
      forwardcurve[T] - beta;

   cvalue[y_] := 
    Sum[daycountfloat*bondprice[swaptionmaturity, k, y]*beta, {k, 
      swaptionmaturity + daycountfloat, 
      swaptionmaturity + swapmaturity, daycountfloat}];
   dvalue[y_] := 
    strike*Sum[
       daycountfix*bondprice[swaptionmaturity, j, y], {j, 
        swaptionmaturity + daycountfix, 
        swaptionmaturity + swapmaturity, daycountfix}] - 1 - 
     alpha + (1 + alpha - 
        daycountfloat*epsilon[swaptionmaturity + swapmaturity])*
      bondprice[swaptionmaturity, swaptionmaturity + swapmaturity, 
       y] - Sum[(alpha - alpha + daycountfloat*epsilon[k])*
       bondprice[swaptionmaturity, k, y], {k, 
       swaptionmaturity + daycountfloat, 
       swaptionmaturity + swapmaturity - daycountfloat, 
       daycountfloat}];
   avalue[y_] := cvalue[y]*psi;
   bvalue[y_] := dvalue[y]*psi;

   (* Standard deviation of the logarithm of the independent \
martingale X *)
   v = eta*(swaptionmaturity^0.5);

   (* Definition of the integrand *)

   blackformula[F_, K_, V_, w_] := 
    w*F*CDF[NormalDistribution[0, 1], w*(Log[F/K] + (V^2)/2)/V] - 
     w*K*CDF[NormalDistribution[0, 1], w*(Log[F/K] - (V^2)/2)/V];
   integrandvalue[y_] := 
    If[avalue[y] > 0 && bvalue[y] > 0, 
     blackformula[avalue[y], bvalue[y], v, 1],
     If[avalue[y] < 0 && bvalue[y] < 0, 
      blackformula[-avalue[y], -bvalue[y], v, -1],
      If[avalue[y] >= 0 && bvalue[y] <= 0, avalue[y] - bvalue[y] ,
       If[avalue[y] <= 0 && bvalue[y] >= 0, 0 , Null
        ]
       ]
      ]
     ];


   (* Integral discretization - Method used : Trapèze *)

   inf = -3*(oisprocessvariance)^0.5;
   sup = 3*(oisprocessvariance)^0.5;
   step = (sup - inf)/200;
   integralnodepoints = 
    Table[integrandvalue[y]*oisprocessdensity[y], {y, inf, sup, step}];
   integralvalue = 
    Sum[step*((integralnodepoints[[j + 1]] + integralnodepoints[[j]])/
        2), {j, 1, Length[integralnodepoints] - 1}];

   (* Determination of the Swaption price *)

   swaptionprice = discountcurve[swaptionmaturity]*integralvalue

   ];

Part two : I define the observations and the function that I have to minimize :

(* Parematers linked to the swaption contract *)
tenor = 6;
strike = 0;
psi = -1;

(*swaptionpricefunction[1.24,1.02,1.07,0.49,1.09,tenor,\
swaptionmaturity,swapmaturity,strike,psi]*)


test = {{5, 5, 2}, {7, 5, 3}};
minimizingfunction = 
  FindMinimum[{Sum[((test[[i, 3]] - 
         swaptionpricefunction[kappa, sigma, eta, alpha, beta, tenor, 
          test[[i, 2]], test[[i, 3]], strike, psi])^2), {i, 1, 
      Length[test]}]}, {kappa, 1}, {sigma, 1}, {eta, 1}, {alpha, 
    1}, {beta, 1}]

Sorry if the code is not "clean", I am just starting with mathematica ...

Thank you very much !

POSTED BY: redone norie
Posted 9 years ago

without having gone through the code in detail, I have a couple observations that might be useful:

  1. defining f[x]:=..., g[x]:=..., h[x]:=... and then doing something like h[g[f[x]]] might be slow(er) than using a definition for h[x] that encapsulates the definitions of f[x] and g[x] directly (without needing to define those functions too).

  2. Using nested Table[] and probably Sum[] structures might slow the whole computation down. Perhaps it would be more helpful if you could figure out what a Table[...,Table[..., {...}],{...}] looks like and construct it beforehand (perhaps using Outer[] or Table[]). The idea is that passing a big complicated table as a parameter might be better than reconstructing the table in every iteration.

  3. Optimization in mathematica is sometimes a pain (in my personal experience). Most of the times, though, using N[] or perhaps something like N[10^-x Round[...10^x]] in places where numerical quantities are expected does the trick. Also playing with the precision options of the routine might help too.

POSTED BY: nik tsak
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract