calibrator.java

金融资产定价,随机过程,MONTE CARLO 模拟 JAVA 程序和文档资料

 *
 *                      READ SWAPTION PRICES
 *
 ******************************************************************************/       
       
   /** Reads the swaptions into the vector <code>swaptions</code>. 
    *  Missing swaptions are ignored during calibration.
    *  The file containing the swaptions must adhere to the following
    *  <a href=#swaptionFile>file format</a>.
    */ 
    private void readSwaptions()
    {
        // the current swaption and its fields
        int p,q; double kappa, price;
        
        while(!swaptionReader.eof()){
            
            p=swaptionReader.readInt();
            q=swaptionReader.readInt();
            kappa=swaptionReader.readDouble();
            price=swaptionReader.readDouble();
            
            if(!swaptionReader.eof())
            swaptions.add(new Swpn(p,q,kappa,price));
        }
        swaptionReader.close();
        
    } // end readSwaptions
    
                           
    
       
   
/*******************************************************************************
 *
 *                             CONSTRUCTOR
 *
 ******************************************************************************/    

    
    /** <p>This constructor does not allocate the file readers to read the
     *  caplet and swaption prices. It is used for unit testing of other
     *  methods.</p>
     *
     * @param n dimension (number of Libors including <code>L_0</code>)
     * @param l initial Libors <code>l[j]=L_j(0)</code>.
     * @param c the <code>c[j]=c_j</code> in the definition of the volatilities
     * <code>sigma_j(t)</code>.
     * @param Tc tenor structure, <code>Tc[j]=T_j</code>.
     */
    public Calibrator(int n, double[] l, double[] c, double[] Tc) 
    {
        
        this.n=n;
        this.l=l;
        this.c=c;
        this.Tc=Tc;
        
        // allocate delta array and compute deltas
        delta=new double[n];
        for(int j=0;j<n;j++)delta[j]=Tc[j+1]-Tc[j];
        
        // initial X-Libors X_j(0)=delta_jL_j(0)
        x=new double[n];
        for(int j=0;j<n;j++)x[j]=delta[j]*l[j];

        factorLoading=null;
        
        // allocate and initialize the zero coupon bonds
        Bi=new double[n+1];
        for(int i=0; i<n+1;i++)Bi[i]=B0(i);
        
        // allocate and initialize the annuities and swap rates,
        // slightly wasteful (unused cells) but less akward
        B_pq=new double[n][];
        S_pq=new double[n][];
        for(int p=0; p<n;p++){ 
            
            B_pq[p]=new double[n-p+1];
            S_pq[p]=new double[n-p+1];
            for(int q=p+1; q<n+1;q++){ B_pq[p][q-p]=B_pq(p,q);
                                       S_pq[p][q-p]=swapRate(p,q); }
        }
        
        // array of Cholesky roots of the covariation matrices
        tCholeskyRoots=new ColtMatrix[n-1];
        
        // array of volatility factors c_j from sigma_j(t)=c_j*g(1-t/T_j)
        c=new double[n];
        
        // caplet array, cell 0 not used but makes indexing more natural
        capletSigma=new double[n];
        
    } // end constructor
    
    
    
    /** <p>This constructor allocates everything necessary to carry out the
     *  calibration procedure.</p>
     *
     * <p>The file <a name="capletFile"><code>capletFile<code></a> 
     * must contains the prices of the caplets
     * <code>cplt([T_i,T_{i+1}],kappa)</code> in precisely the following 
     * form: <code>i kappa price</code>, where this pattern is repeated
     * separated only by spaces. Calibration will try to match all caplet
     * prices. Missing caplets will be interpolated from nearest neighbors.</p>
     *
     * <p>The file <a name="swaptionFile"><code>swaptionFile<code></a>
     * must contains the prices of the swaptions
     * <code>swpn([T_p,T_q],kappa)</code> in precisely the following 
     * form: <code>p q kappa price</code>, where this pattern is repeated
     * separated only by spaces. Calibration will try to match the existing
     * swaption prices by minimizing the squared error.</p>
     *
     * @param n dimension (number of Libors including <code>L_0</code>)
     * @param l initial Libors <code>l[j]=L_j(0)</code>.
     * @param c the <code>c[j]=c_j</code> in the definition of the volatilities
     * <code>sigma_j(t)</code>.
     * @param Tc tenor structure, <code>Tc[j]=T_j</code>.
     * @param capletFile file containing the caplet prices.
     * @param swaptionFile file containing the swaption prices.
     */
    public Calibrator(int n, double[] l, double[] c, double[] Tc, 
                      String capletFile, String swaptionFile) 
    {
        this(n,l,c,Tc);
        capletReader=new EasyReader(capletFile);
        swaptionReader=new EasyReader(swaptionFile);
    }
        

             
/*******************************************************************************
 *
 *        ZERO COUPON BONDS AND ANNUITY NUMERAIRE AT TIME ZERO
 *
 ******************************************************************************/

        
 
    /** The zero coupon bond <code>B_i(0)=B(0,T_i)</code>.
     *  
     * @param i bond matures at time <code>T_i<code>
     */
     public double B0(int i)
     { 
         double f=1;
         // accumulate 1 from time t=0 to time t=T_i                    
         for(int j=0;j<i;j++)f*=(1+x[j]); 
         return 1/f;                                 // B_i(0)=1/f  
     }   

     

/*******************************************************************************
 *
 *                        FORWARD SWAP RATES
 *
 ******************************************************************************/                         
             
   /** <p>The forward swap rate <code>S_pq(t)=k(t,[T_p,T_q])</code> at time 
    *  <code>t=0</code>.</p>
    *
    * @param p swap begins at <code>T_p</code>.
    * @param q swap ends at <code>T_q</code> (settled).
    */ 
     public double swapRate(int p, int q)
     { 
        double f=1+x[q-1], S=delta[q-1];
        for(int k=q-2;k>=p;k--){ S+=delta[k]*f; f*=(1+x[k]); }
 
        return (f-1)/S;
     } //end Swap_Rate

     
/*******************************************************************************
 *
 *                         ANNUITY NUMERAIRE
 *
 ******************************************************************************/                    
 
     
   /** <p>The annuity <code>B_pq(t)=sum_{k=p}^{q-1}delta_kB_{k+1}(t)</code>
    *  at time <code>t=0</code>.
    *  </p>
    *
    * @param p annuity begins at <code>T_p</code>.
    * @param q annuity ends at <code>T_q</code> (settled).
    */ 
     public double B_pq(int p, int q)
     { 
         double S=0, F=B0(q);
         for(int k=q-1;k>=p;k--){ S+=delta[k]*F; F*=(1+x[k]); }
         return S;
     } //end B_pq

     
/*******************************************************************************
 *
 *               LIBOR VOLATILITIES, CAPLET IMPLIED VOLS AND PRICES
 *
 ******************************************************************************/            

      
     /** <p>The quantity <code>Sigma=sigma*sqrt(T_i)</code>, where 
      *  <code>sigma</code> is the annualized volatility 
      *  of the caplet <code>cplt([T_i,T_{i+1}],kappa)</code>.</P>
      *  
      *  <p>In other words <code>Sigma</code> is the squareroot of the quadratic 
      *  variation <code>&lt;log(L_i)&gt;_0^{T_i}</code>. 
      *  Remember that this is the quantity that must be entered into the
      *  Black caplet formula.</p>
      *
      * @param i caplet caps Libor on <code>[T_i,T_{i+1}]</code>
      */
     public double capletSigma(int i)
     {
         double
         sigmasqr=factorLoading.integral_sgi_sgj_rhoij(i,i,0,Tc[i]);
         
         return Math.sqrt(sigmasqr);
     }
     
     /** <p>The {@link #capletSigma} computed from price 
      *  <code>c</code>, strike <code>kappa</code> and time
      *  <code>T_i</code> by inversion the Black caplet formula.</p>
      *
      * @param i caplet caps Libor on <code>[T_i,T_{i+1}]</code>
      * @param c caplet cash price at time <code>t=0</code>
      * @param kappa strike rate
      */
     public double capletImpliedSigma(int i, double kappa, double c)
     throws NoSolutionException
     {
         double L_i0=l[i],                                // L_i(0)
                // black-scholes function part in black caplet price
                bsf=c/(delta[i]*B0(i+1));     
         return FinMath.BisectionSolveBSF(L_i0,kappa,bsf);
     }
     
     
     /** <p>Black caplet cash price at time zero.</p>
      *
      * @param i caplet caps Libor on <code>[T_i,T_{i+1}]</code>.
      * @param kappa strike rate.
      */
     public double capletPrice(int i, double kappa)
     {
         double  L_i0=l[i],              // L_i(t)
                 cSigma=capletSigma(i),  // <log(L_i)>_0^{T_i}
                 Nplus,
                 Nminus;                     
    
         Nplus=FinMath.N(FinMath.d_plus(L_i0,kappa,cSigma));
         Nminus=FinMath.N(FinMath.d_minus(L_i0,kappa,cSigma));
     
         return delta[i]*B0(i+1)*(L_i0*Nplus-kappa*Nminus);
   
     } // end capletPrice
     
/*******************************************************************************
 *
 *                   SWAP-RATES AND SWAPTION PRICES
 *
 ******************************************************************************/             
     
     
     /** <p>The weight <code>w^{p,q}_j(t)</code> in the representation of the 
      *  swap rate <code>S_pq(t)</code> as a convex combination of Libors at 
      *  time <code>t=0</code>. See <i>LiborProcess.ps</i>.
      */
     public double wpq(int p, int q, int j)
     { 
         return delta[j]*Bi[j]/B_pq[p][q-p];
     }
         
   
     /** <p>The vector <code>x_pq(0)</code> used in the computation 
      *  of the approximation of the swap rate volatility. See 
      *  <i>LiborProcess.ps</i>. Note that the vector is allocated in size
      *  <code>n-p</code> and set to zero beyond coordiante <code>q-1</code>
      *  to fit with the <code>n-p by n-p</code> matrix 
      *  <code>choleskyRoots[p]</code>.
      *
      * @param p swap begins at <code>T_p</code>
      * @param q swap ends at <code>T_q</code>
      */
     public ColtVector xpq(int p, int q)
     { 
         ColtVector x=new ColtVector(n-p);
         for(int j=p;j<q;j++){
             
             double x_j=wpq(p,q,j)*l[j];
             x.setQuick(j-p,x_j);
         }
         return x;
     } // end xpq
    
    
    /** <p>Deterministic approximation to the aggregate volatility (square root 
     *  of the quadratic variation <code>&lt;log(S_pq)&gt;_0^{T_p}</code>) to 
     *  expiration of the swap rate <code>S_pq</code>.</p>
     *
     * <p>This quantity is used in the analytic approximation to the swaption
     * price. See <i>LiborProcess.ps</i></p>
     *
     * @param p swap starts at <code>T_p</code>
     * @param q swap stops at <code>T_q</code>
     */
    public double Sigma(int p, int q)
    {
        ColtMatrix Q=tCholeskyRoots[p];
        double nrm=xpq(p,q).timesEquals(Q).L2Norm();
        
        return nrm/S_pq[p][q-p];
    }
    
 
          
   /** <p>Analytic approximation to the swaption price. Based on a log-normal
    *  distribution of the swap rate.