SERVING THE QUANTITATIVE FINANCE COMMUNITY

Cuchulainn
Posts: 53306
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: System of nonlinear equations

It could happen to a bishop! Seriously, people can get stuck in a rut and then always have several solutions.

https://en.wikipedia.org/wiki/How_to_Solve_It

I am surprised DE does not converge in your case. DE is very robust. So, it might be bad input or a bug in Matlab , or other?
http://www.datasimfinancial.com

“Sir Walter Scott created rank & caste in the South and also reverence for and pride and pleasure in them"

Mark Twain

Alan
Posts: 9171
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: System of nonlinear equations

@Daniel, if you still want the paraphrase: the basic model (Hull example), models the equity as a call option -- with the asset value $V_t$ playing the role of the underlying and a call strike equal to the debt value L. In that model, the company goes bankrupt if and only if the asset value is less than the debt at expiration T. But paths where the intermediate asset values (prior to expiration) fall below the strike (=debt) are fine (do not result in bankruptcy) as long as $V_T > L$.

The second model tightens up the bankruptcy criterion so that the company is bankrupt if the asset value i) falls below the barrier H at *any* time prior to T  or ii) is below L at expiration. The natural choice is H = L, but you could argue that H < L could make sense, too. In any event, If H is much smaller than L (H << L), the likelihood of the particle reaching the barrier is small, so it's as if the barrier does not exist, and so you should recover the basic model results.

Cuchulainn
Posts: 53306
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: System of nonlinear equations

Gencho wrote:
P.S. I am very surprised, to say at least, that I could not guess myself the following:
1) I could not understand that scaling is so important.
2) I could not even think that the root solver will work worse than constrained optimization. I thought it is the most robust approach.
So again, you really helped me!

Once a given approach becomes accepted then people just tend to use and if necessary tweak it until it works (or not). In the current case the usual suspect is a variation on Newton Raphson.
Of course the "dropping penny" is the realisation $f(x) = 0 \Leftrightarrow min f(x)^2.$

There was a discussion on this very topic with many (ingenious) solutions. AFAIR minimisaton was the most robust.

viewtopic.php?f=34&t=97812

//
Another robust method to resolve NR is homotopy/continuation/Davidenko (from topology) is to embed f(x) = 0 in a lager space to get H(x,t) = where 0 <= t <= 1 and solve the corresponding ODE.
http://www.lehigh.edu/~wes1/apci/24feb00-s.pdf

Here's code in 1d for Kepler's equation; In 2d you would us the Jacobian. AFAIR Matlab has a ODE solver ode45??
#include <iostream>#include <vector>#include <cmath>#include <functional>#include <boost/numeric/odeint.hpp>#include <boost/numeric/odeint/stepper/bulirsch_stoer.hpp>// The type of container used to hold the state vector typedef std::vector<double> state_type;typedef std::function<double (double)> FunctionType;namespace Bode = boost::numeric::odeint;using namespace boost::numeric::odeint;double p = 2.0;double c = 10.0;const double eps = 13.5;const double M = 0.8;double Kepler(double y){ double d = y - eps*std::sin(y) - M; return d;}// NR function stuffdouble f(double x){ return Kepler(x);}double fd(double x){  // Kepler return 1.0 - eps*std::cos(x);}//double x0 = 138.3377; // far away//double x0 = 12.3; // far awaydouble x0 = 2.3000001; // try -20.0 as wellclass WriteOutput{public: std::vector<double> tValues; std::vector<double> yValues;public:    WriteOutput()  : tValues(std::vector<double>()), yValues(std::vector<double>()) { }  void operator ()( const state_type &x , const double t ) { tValues.push_back(t); yValues.push_back(x[0]);  std::cout << "Time and value (FO): " << t << " " << x[0] << std::endl; } };// The rhs of x' = f(x) defined as a class class Davidenko{public:    Davidenko() { }    void operator() (const state_type& x, state_type& dHdt, const double t)    {  dHdt[0] = -f(x0)/fd(x[0]);    }};void write_out( const state_type& x , const double t ){ if (t == 1.0) { std::cout << "Time and value: " << t << " " << x[0] << std::endl; }}int main(){    namespace Bode = boost::numeric::odeint;    // Initial condition    state_type x(1); x[0] = x0; double L = -1.0; double U = 2.0; double dt = 1.0e-3; // Integration_class, V2 based on Runge Kutta54 Cash Karp stepper (5th order) Davidenko ode;    std::size_t steps = Bode::integrate(ode, x, L, U, dt, write_out); std::cout << "Number of steps " << steps << std::endl; std::cout << x[0] << ", " << f(x[0])<< std::endl;  WriteOutput wo;// Choice B: Fast adaptive typedef runge_kutta_cash_karp54< state_type > error_stepper_type; typedef Bode::controlled_runge_kutta< error_stepper_type > controlled_stepper_type; Bode::runge_kutta_fehlberg78< state_type, double> myStepper;// steps = Bode::integrate_const(myStepper, ode, x, L, U, dt, boost::ref(wo));// double abs_err = 1.0e-3 , rel_err = 1.0e-3, a_x = 1.0 , a_dxdt = 1.0; double dtODE = 1.0e-3;  //   controlled_stepper_type controlled_stepper(default_error_checker< double >( abs_err , rel_err , a_x , a_dxdt ) );// steps = Bode::integrate_adaptive( controlled_stepper,ode, x , L, U, dtODE, boost::ref(wo) );  // std::cout << "Steps " << steps << ", " << f(x[0]) << std::endl;  return 0;}
Last edited by Cuchulainn on June 20th, 2017, 12:25 pm
http://www.datasimfinancial.com

“Sir Walter Scott created rank & caste in the South and also reverence for and pride and pleasure in them"

Mark Twain

Cuchulainn
Posts: 53306
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: System of nonlinear equations

Alan wrote:

@Daniel, if you still want the paraphrase: the basic model (Hull example), models the equity as a call option -- with the asset value $V_t$ playing the role of the underlying and a call strike equal to the debt value L. In that model, the company goes bankrupt if and only if the asset value is less than the debt at expiration T. But paths where the intermediate asset values (prior to expiration) fall below the strike (=debt) are fine (do not result in bankruptcy) as long as $V_T > L$.

The second model tightens up the bankruptcy criterion so that the company is bankrupt if the asset value i) falls below the barrier H at *any* time prior to T  or ii) is below L at expiration. The natural choice is H = L, but you could argue that H < L could make sense, too. In any event, If H is much smaller than L (H << L), the likelihood of the particle reaching the barrier is small, so it's as if the barrier does not exist, and so you should recover the basic model results.

Clear. Thanks.
There is a PDE to describe first exit time. For Merton this is OK. For double barriers it should also be feasible? The company's equity would be the derivative.
http://www.datasimfinancial.com

“Sir Walter Scott created rank & caste in the South and also reverence for and pride and pleasure in them"

Mark Twain

Alan
Posts: 9171
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: System of nonlinear equations

Well, you'd have to have some plausible boundary condition at the putative upper barrier U. Maybe you could say the company gets acquired there, boosting the equity value to some  value $A$, so $E(t,V) = A$ if $V = U$. I'd want $A \ge U - L$, where $L$ is the debt; otherwise the problem sounds wonky. (So, this should keep the equity value (= call price) monotone increasing in the asset value $V$ (= underlying price) at a fixed time slice, which might be needed for a solution to the nonlinear system to exist).

Gencho
Topic Author
Posts: 11
Joined: June 14th, 2017, 7:12 am

### Re: System of nonlinear equations

Be careful using this forum! I just wrote a huge message and I was asked to login again. of course all the text I wrote got destroyed.
Dear all, thanks again for help! The DPs are in reasonable bounds now, as I said previously.
Nevertheless, from the Question 2 pdf you may see that the DPs are still not stable everywhere. You may see windows of time where our DPs are aligned with the benchmark and when they are very different. Therefore I have a question:
What is the best way to find the parameter/issue that causes divergence of my results compared to the benchmark? ( I believe that divergence is because of the drift parameter)
Subquestions:
1) Is MLE a good estimate for \mu_V? ( Acknowledging that I find \sigma_V,V_0 from the system of nonlinear equations)
2) How should I estimate V_0 if I use MLE for finding \mu_V, \sigma_V?
Initially I thought of solving 3 equations together, hence finding 3 roots: \mu_V,V_0,\sigma_V. But it is seems wrong(details in pdf).
Additionally, I have been reading about non-linearity of the barrier option itself, and of its delta greek, on which I heavily rely on. I have checked as well sensitivity of DPs to the change of maturity time( more detailed info in the file), so I deal with very non linear problem
Thank you!