Indirect Robust Estimation of short - term interest rate processes(with Mathematica illustrations)
            by Nabeel Butt
                 Department of Applied Mathematics,UWO
                 28th November '2007

■  Generalized Method of Moment:
The generalized method of moments is a very general statistical method for obtaining estimates of parameters of statistical models.It is a generalization, developed by Lars Peter Hansen, of the method of moments.The term GMM is very popular among econometricians but is hardly used at all outside of economics, where the slightly more general term estimating equations is preferred.
  The idea of the generalized method of moments is to use moment conditions that can be found from the problem with little effort.We assume that the data are a stochastic process ((Subscript[Y, 1]),(Subscript[Y, 2]) ,..................).In mathematical language, we start out with a (vector valued) function f that depends both on the parameter and a single observation and that has mean zero for the true value of the parameter, θ=(Subscript[θ, o]), i.e.E[f((Subscript[Y, i]),(Subscript[θ, 0]))]=0.To turn this function into a parameter estimate, we minimize the associated chi - square statistic:
  (OverHat[θ])=arg (min_θ) (((f((Subscript[Y, i]),θ))^T)Af((Subscript[Y, i]),θ))
The construction of GMM estimator requires a set of orthogonality conditions.

where superscript T denotes transpose, and A is a positive definite "weighting" matrix.A may be known a priori or estimated from the sample.

http://www.mathstatica.com/book/

  The Mathematica Journal Volume 9,Issue
 " Skew Densities and Ensemble Inference for Financial Economics " by H. D. Vinod
 Usage of mathstatica package
■  Robust Generalized Method of Moment:  
■  Indirect Generalized Method of Moment:

Paper Review :
 "Indirect robust estimation of the short-term  interest rate process"
                             by  Veronika Czellar , G. Andrew Karolyi and Elvezio Ronchetti 
                              Journal of Empirical Finance 14 (2007) 546–563

Abstract of the paper :
  Indirect Robust Generalized Method of Moments (IRGMM) , a simulation - based estimation methodology is used to model the short - term interest rate processes.The primary advantage of IRGMM relative to classical estimators of continous - time short - rate diffusion processes is that it corrects errors due to discretization and the errors due to model mis - specification.The approach is applied to monthly US risk free rates and to various monthly Eurocurrency rates which provides extensive evidence of it' s predictive performances in a variety of settings. 
  
   Introduction :
 Understanding the dynamics of the short - term interest rate is important for building models of the term structure of interest rates, which, in turn, are used for pricing interest rate derivatives and instruments with embedded options such as callable bonds,options on bonds and mortgage - backed securities.The interest rate derivatives market is the largest derivatives market in the world. Market observers estimate that $60 trillion dollars by notional value of interest rate derivatives contract had been exchanged by May 2004.
Basic modeling issues remain unresolved. One issue stems from the difficulties associated with the statistical analysis of continuous-time processes. For example, in complex statistical models,like diffusion models described by stochastic differential equations (SDE) of the form:
d(Subscript[y, t])=η((Subscript[y, t]),θ)dt+σ((Subscript[y, t]),θ)d(Subscript[W, t])  -----------------------equation 1
where θ∈(R^p), η(·) and σ(·) are the drift and the volatility, respectively, and (Subscript[W, t]) a Wiener process,it is often difficult or impossible to carry out standard likelihood based estimation and inference.Such SDEs play an important role when modeling the short-term interest rate; see, for example,Vasicek (1977), Dothan (1978), Brennan and Schwartz (1977), Cox, Ingersoll and Ross (1981,1985), Chan, Karolyi, Longstaff and Sanders (CKLS, 1992), Brenner, Harjes and Kroner (1996),and Ahn and Gao (1999), among others.
Notable work done to estimate continous-time short-term process:
1)Efficient Method of Moments (Gallant and Tauchen, 1996; Andersen and Lund,1997)
2)Weighted least squares estimation (Chapman, Long and Pearson, 1999)
3)Simulated maximum likelihood estimation (Durham and Gallant, 2002; Durham, 2003)
4)Gibbs sampling-based Markov Chain Monte Carlo algorithms (Elerian et al., 2001; Eraker, 2001)
5)Retrospective sampling-based Monte Carlo methods (Beskos et al., 2006)

Two problems with estimation based upon dicretized SDE:
1)Estimated parametres are biased.
2)Model-misspecification can also lead to biased results.The theory of robust statistics can be used to mitigate this problem. Specifically, Dell'Aquila, Ronchetti and Trojani (2003) used the Robust Generalized Method of Moments (RGMM) technique, developed by Ronchetti and Trojani(2001), to estimate discrete-time interest rate processes.

Indirect estimation for the problem mentioned above:
Typically the estimation of Eq. (1) is performed by means of an auxiliary model which is a discretized version of the SDE. The resulting indirect estimation (Gouriéroux et al., 1993) is based on the following idea. Given a sample of observations  {(Subscript[y, t])}t=1,2,...........,n assumed to be generated from a probability model , θ∈(R^p) , define an auxiliary model (Subscript[(OverTilde[F]), μ])where the parameter μ∈(R^r) is easier to estimate than θ. For instance, (Subscript[F, θ]) could be the diffusion model (1) or a fine discretization of (1) and (Subscript[(OverTilde[F]), μ]) a crude discretization of (1).

Bias correction(Similar to control-variate technique for pricing financial options ?)-Indirect robust estimation technique:
 The main steps of indirect estimation are:
(i) an auxiliary estimator  (OverTilde[μ ]) is calculated with the original sample;
(ii) pseudo-observations are simulated from the model (Subscript[F, θ  ])generating S samples of pseudo-data
(Subscript[{((Subscript[y, t])^1)}, t=1,2,3,..............,n]) ,......................,(Subscript[{((Subscript[y, t])^S)}, t=1,2,3,..............,n])
(iii) an auxiliary estimator (Subscript[(OverTilde[μ]), S]) (θ) is calculated with the simulated pseudo-observations;
(iv) the indirect estimator (OverHat[θ]) of θ is obtained as that which minimizes the distance between the
auxiliary estimators (OverTilde[μ ])  and (Subscript[(OverTilde[μ]), S]) (θ):

(OverHat[θ])=arg (min_θ) ((((OverTilde[μ ]) -(Subscript[(OverTilde[μ]), S])(θ)))^T)Ω((OverTilde[μ ]) -(Subscript[(OverTilde[μ]), S])(θ))   (Different worlds have different notions of distances) 

Robust indirect estimation can reduce the bias due to discretization and model-misspecification:

In particular, the authors define the Indirect Robust Generalized Method of Moments(IRGMM), a simulation-based estimation of SDEs which reduces both the bias due to discretization and to contamination, and apply it to monthly US risk free and to various monthly Eurocurrency interest rate series. They also investigate three different diffusion models for shortterm interest rates proposed in the finance literature (CKLS, 1992; Ahn and Gao, 1999; Chapman and Pearson, 2000) to explore the sensitivity of our findings to different model specifications.We compare empirical results obtained with GMM, RGMM and IRGMM in both within- and out-of sample predictive tests. For most of the cases, we obtain superior predictive performances using the IRGMM estimation technique.

IRGMM analysis of the short-term interest rate process:
Robust Indirect inference
A necessary condition for the consistency of an indirect estimator defined in Eq. (2) is that the original data {(Subscript[y, t])} and the pseudo-data are generated from the same probability model, see Gouriéroux et al. (1993). In reality, the presence of large outliers and high kurtosis in the increments of interest rate data shows that a diffusion model with standard normal increments may be misspecified. As an illustration, let us consider the diffusion model (1) where the increments come from an ε-neighborhood of a standard normal distribution defined by:
(Subscript[W, t])-(Subscript[W, t-1]) =(Subscript[ϵ, t])~(Subscript[G, ϵ])=(1-ϵ)N(0,1)+ϵG
with 0≤ϵ<1 and G an unknown distribution.


As an illustrative example, Fig. 2.1 presents simulated data of sample size n=300 generated from a contaminated geometric Brownian motion with drift:
d(Subscript[y, t])=β(Subscript[y, t])dt+σ(Subscript[y, t])d(Subscript[W, t])
where the increments follow (Subscript[G, ϵ])=(1-ϵ)N(0,1)+ϵt3 with ε=0.05, β=0.1, σ=0.5 and t3 a t-distribution with 3 degrees of freedom. Fig. 2.1 shows that a contaminated diffusion model
can generate large outliers in the data. The kurtosis of (Subscript[y, t]) -(Subscript[y, t-1])for the plotted sample is approximately 10.50, which is much higher than 3, the kurtosis under the normal model.

nor[mu_,sigma_]:=Random[NormalDistribution[mu,sigma]];
student[v_]:=Random[StudentTDistribution[v]];
deltaY[y_,beta_,sigma_,epsilon_]:=beta*y +
	sigma*y*((1-epsilon)*nor[0,1]+epsilon*student[3]);
shortY[y_,beta_,sigma_,epsilon_]:=y+deltaY[y,beta,sigma,epsilon];
SeedRandom[10]
tt=NestList[shortY[#,0.1,0.5,0.05]&,0.01,300];
ListPlot[tt,Joined->True,PlotLabel->"Simulated data from GBM diffusion with 5% contamination",AxesLabel->{"Observation","data"}]

Our goal is to construct an indirect estimator of the diffusion model (1) which is robust to misspecifications of the underlying stochastic structure of the model described in Eq. (3). Genton
and Ronchetti (2003) have shown that indirect estimators are robust if the auxiliary estimator is a robust estimator of the parameter of the auxiliary model. In our case, the auxiliary model (OverTilde[(Subscript[F, μ])]) is a
discretization of the diffusion model (1):
(Subscript[y, t])=(Subscript[y, t-1])+η((Subscript[y, t-1]),μ)+σ((Subscript[y, t-1]),μ)(Subscript[ϵ, t  ])-----------------------------------equation5
Where (Subscript[ϵ, t]) are i.i.d standard normal variables.Equation 5 is referred to as by authors as crude discretization.
In order to use Genton and Ronchetti's result to construct a robust indirect estimator of the parameter θ of the diffusion model, we need a robust estimator of the parameter μ. The robust version of the GMM estimator (see Ronchetti and
Trojani, 2001) is a GMM estimator where the classical H-dimensional orthogonality condition h is replaced by a truncated orthogonality condition:

(Subscript[(h^(A,τ)), c])((Subscript[y, t]),(Subscript[y, t-1]),μ)=A[h((Subscript[y, t]),(Subscript[y, t-1]),μ)-τ](Subscript[w, c])((Subscript[y, t]),(Subscript[y, t-1]),μ)--------------------------------equation6

Where c>(Sqrt[H ]) is a tuning constant and (Subscript[w, c])((Subscript[y, t]),(Subscript[y, t-1]),μ)=min(1,(c/(||A(h((Subscript[y, t]),(Subscript[y, t-1]),μ)-τ)||)))
The non-singular matrix A∈(R^H)×(R^H) and the vector τ∈(R^H)  are determined by the implicit equations:

(Subscript[E, (Subscript[μ, o])])((Subscript[h, c])^(A,τ))=0--------------------equation7
(1/n)(∑_(t=1)^n)((((Subscript[h, c])^(A,τ))(((Subscript[h, c])^(A,τ))))^T)=I---------------------------------------------------------------equation8

If τ=0 equation 7 can be eliminated.

Denote by (OverTilde[μ]) the robust GMM(RGMM) estimator associated with the truncated orthogonality function(6).By construction,the RGMM estimator is defined by:
(Overscript[μ, `])=arg (min_μ) g({(Subscript[y, t])},μ),
g({(Subscript[y, t])},μ)=(1/n)(∑_(t=1)^n)(((Subscript[(h^(A,τ)), c])((Subscript[y, t]),(Subscript[y, t-1]),μ))^T)(1/n)(∑_(t=1)^n)(Subscript[(h^(A,τ)), c])((Subscript[y, t]),(Subscript[y, t-1]),μ)

The authors define the indirect robust GMM(IRGMM) estimator as an indirect estimator computed with an auxiliary RGMM estimator.The authors then present an algorithm to compute the IRGMM estimator.

The IRGMM technique:

 We define the IRGMM estimator (OverHat[θ])  according to the following steps:
 
 1) Compute the covariance matrix A,vector τ, and the RGMM estimator (OverTilde[μ]) defined in equations(7)-equations(10) using the real data
 2)Simulate S pseudo-data sets using discretized version of the model
 3)For each s ,calculate the function g({((Subscript[y, t])^s)},μ) defined in equation (10) using A and τ obtained in 1) and the auxilliary estimator for the pseudo-data:
(Subscript[(OverTilde[μ]), S])(θ)=arg (min_μ)(∑_s)g({((Subscript[y, t])^s)},μ)
4)The indirect robust estimator is one which minimizes the distance between the auxilliary estimators (OverTilde[μ]) and (Subscript[(OverTilde[μ]), S])(θ):
(OverHat[θ])=arg (min_θ) ((((OverTilde[μ ]) -(Subscript[(OverTilde[μ]), S])(θ)))^T)Ω((OverTilde[μ ]) -(Subscript[(OverTilde[μ]), S])(θ)) 

where Ω is a positive definite symmetric matrix e.g the inverse of the assymptotic covariance matrix of the auxilliary estimator.



Computational complexity of the algorithm.Is the usage of GRID MATHEMATICA a suitable solution?

2.3 The author's choice of the parametres S and δ to achieve convergence in estimates:
To calibrate their procedure authors generate data of sample size n=100 from the exact discretization of Eq. (4) with β=0.1, σ=0.5 and standard normal increments. The computed indirect
estimators with auxiliary maximum likelihood estimators of the crude discretization of Eq. (4).In Figs. 2.3–2.4 they plot the indirect parameter estimates of β and σ obtained with S=10, 100,
1000 pseudo-data sets and with increasing values of 1 / δ. The plots show that there is no significant gain in stability of the indirect estimator computed with S=1000 compared to the one
computed with S=100 if taken into account the fact that computation time is increasing linearly with S.Fig. 2.4 shows that the indirect estimator settles down for 1/ δ≥10, so a choice of δ≤1/10 is
necessary. In the next section, author's estimate monthly rates and in this case it is natural to choose δ=1/22 for the simulation of daily rates (given there are about 22 trading days per month).


3.Model estimation:
The data considered by the authors is the one month risk free rates from the Monthly CRSP Fama Risk Free Rates File covering the period from June 30, 1964 to December 31, 2004, which is an extension of the
dataset used by Chan et al. (1992).
Interest rate model:d(Subscript[y, t])=((Subscript[β, 0])+(Subscript[β, 1])(Subscript[y, t])+(Subscript[β, 2])((Subscript[y, t])^2)+((Subscript[β, 3])/(Subscript[y, t])))dt+σ((Subscript[y, t])^γ)d(Subscript[W, t])

The estimate our dataset using the three models and three different techniques: classical GMM,RGMM and IRGMM, as defined in Section 2. For the estimation of the unrestricted model by GMM, we consider Chapman and Pearson's (2000) orthogonality conditions:

E[h((Subscript[y, t]),(Subscript[y, t-1]),μ)]=0

h((Subscript[y, t]),(Subscript[y, t-1]),μ)=(({{(Subscript[v, t])}, {(Subscript[v, t])(Subscript[y, t-1])}, {(Subscript[v, t])((Subscript[y, t-1])^2)}, {(Subscript[v, t])((Subscript[y, t-1])^(-1))}, {((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^(2γ))}, {(((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^(2γ)))(Subscript[y, t-1])}}))



where  (Subscript[v, t])=(Subscript[y, t])-(Subscript[y, t-1])-(Subscript[β, 0])-(Subscript[β, 1])(Subscript[y, t-1])-(Subscript[β, 2])((Subscript[y, t-1])^2)-((Subscript[β, 3])/(Subscript[y, t-1]))

For the CKLS model, we use the orthogonality conditions proposed by Chan et al. (1992):


h((Subscript[y, t]),(Subscript[y, t-1]),μ)=(({{(Subscript[v, t])}, {(Subscript[v, t])(Subscript[y, t-1])}, {((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^(2γ))}, {(((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^(2γ)))(Subscript[y, t-1])}}))
where  (Subscript[v, t])=(Subscript[y, t])-(Subscript[y, t-1])-(Subscript[β, 0])-(Subscript[β, 1])(Subscript[y, t-1])

Finally for the GMM estimation of the Ahn and Gao's model, we
use the orthogonality conditions proposed by Ahn and Gao (1999):
h((Subscript[y, t]),(Subscript[y, t-1]),μ)=(({{(Subscript[v, t])}, {(Subscript[v, t])(Subscript[y, t-1])}, {(Subscript[v, t])((Subscript[y, t-1])^2)}, {(((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^3))}, {(((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^3))(Subscript[y, t-1])}, {(((Subscript[v, t])^2)-(σ^2)((Subscript[y, t-1])^3))((Subscript[y, t-1])^3)}}))
where  (Subscript[v, t])=(Subscript[y, t])-(Subscript[y, t-1])-(Subscript[β, 0])-(Subscript[β, 1])(Subscript[y, t-1])-(Subscript[β, 2])((Subscript[y, t-1])^2)


Estimation and goodness of fit-tests:
Table 3.3 presents the parameter estimates obtained with the three models and each model is estimated with the three different techniques (t-statistics of the parameter estimates are in
parentheses). The IRGMM estimation is performed using S=100 pseudo-data sets generated from the fine discretization (11) with δ=1/ 22.4 The weighting matrix in Eq. (13) is replaced by
the inverse of the asymptotic matrix of the RGMM estimator, i.e. Ω= E(([((∂ ((Subscript[h, c])^(A,τ)))/(∂ μ))])^T)E[((∂ ((Subscript[h, c])^(A,τ)))/(∂ μ))].We use the Broyden–Fletcher–Goldfarb–Shanno (BFGS) variant of the Davidon–Fletcher–Powell minimization algorithm for which we compute the gradient using Ridders method of polynomial extrapolation. This minimization algorithm consists of finding a local minimum in the
neighborhood of the initial point, for which we choose the parameter estimates obtained with RGMM. Since each evaluation of the function in Eq. (13) needs an optimization over S=100
simulated datasets, an IRGMM estimation is computationally intensive. The authors had written there entire program in C.It is interesting to see if we could have used the mathematica enviroment for programmin purposes.
In order to compare numerically the different estimation techniques we define the following goodness-of-fit measures(in-sample and out-sample goodness of fit) :

RAMSE=((((1/n)(∑_(t=1)^n)(1/S)(∑_(s=1)^S)(((((Subscript[y, t])^((s)))-(Subscript[y, t])))^2)))^0.5)
AMAD=((1/n)(∑_(t=1(Underscript[median, s]))^n)(|((Subscript[y, t])^((s)))-(Underscript[median(((Subscript[y, t])^((l)))), l])|))
AMB=(1/n)(∑_(s=1)^n)|(Underscript[median, s])(((Subscript[y, t])^((s))))-(Subscript[y, t])|

The smallest values of goodness-of-fit measures for each model are italicized. The smallest value of bRAMSE for all models and estimation techniques is marked by a star.




Can be extended to exchange rates data.Out-of-sample forecasts of eurocurrency data with UK against Switzerland.
For forcasting purposes IRGMM outperforms RGMM and GMM.
 
Conclusion:
In this paper the authors have proposed an application of the indirect robust GMM (IRGMM) simulation-based estimation technique for models of short-term interest rates. We presented 
an empirical comparison of three estimation techniques (GMM, RGMM and IRGMM) for three different interest-rate diffusion models that have received attention in the finance literature. They evaluated the within-sample performance of these techniques for US interest rates between 1964 and 2004. They then compared the out-of-sample predictive power of the three techniques on monthly Eurocurrency rates for the US dollar, UK Sterling and Swiss Franc. They founnd that the improvement in goodness-of-fit provided by IRGMM estimation tends to occur most consistently for those models with more complex specifications, such as those with non-linear drift.A possible future research direction would be to go beyond single factor models and to investigate the performance of IRGMM in more complex diffusion models.

Term structure of interest rates

    Quantitative financial modelling using Mathematica:
    
    "You can write your code,write down a mathematical model,test them,plot the results,and produce a report for risk managment/trading/regulatory bodies all within one enviroment.................." William.T.Shaw (Phd Oxford)

Vasicek's model
This is one of the most widely used models for the pricing of bonds.Vasicek (1977) uses the Ornstein-Uhlenbeck process for the spot interest rater:dr=α(γ-r)dt+σdB

Simulating the Vasicek model
The Vasicek model is aone-factormodel:All rates ultimately depend on the shortest-term interest rate,which we call thespot interest rate(or simply thespot rate) and denote byr.To simulate this rate we discretize the basic Vasicek equation by considering changes in the interest rate over a short period Δt:Δr=α(γ-r)Δt+σZ(Sqrt[Δt])Note that-as opposed to stock price models,which are multiplicative-term structure models areadditive.This means that ifrt is the spot rate at timet,then the spot rate at timet + Δtis given by
(Subscript[r, t+Δ])=(Subscript[r, t])+Δr=(Subscript[r, t])+α(γ-r)Δt+σZ(Sqrt[Δt])

Clear[nor,deltaR,shortTermR, alpha, sigma, r]
nor[mu_,sigma_]:=Random[NormalDistribution[mu,sigma]];
deltaR[alpha_,gamma_,r_,sigma_,deltaT_]:=alpha*(gamma-r)*deltaT +
	sigma*nor[0,1]*Sqrt[deltaT]
shortTermR[alpha_,gamma_,r_,sigma_,deltaT_]:=r+
	deltaR[alpha,gamma,r,sigma,deltaT]

Set α=3%,γ=4%,σ=0.12,(Subscript[δ, t])=0.0001.This means that the time step is about 1 hour and that(Subscript[r, today])=3%.We now simulate 300 steps.

SeedRandom[2]
tt=NestList[shortTermR[0.03,0.04,#,0.12,0.0001]&,0.03,300];
ListPlot[tt,Joined->True,PlotLabel->"Vasicek Simulation"]

The Cox-Ingersoll-Ross term structure model
CIR consider an interest rate process of the type:Δr=α(γ-r)Δt+σ(Sqrt[r])dBwhere α>0,γ=long-run mean interest rate,r =current interest rate.

Simulating the CIR spot interest rate process
We can repeat the simulation we performed above for Vasicek, simulating the CIR short - term interest rate process :

Clear[ deltaR, shortTermR, alpha, gamma, sigma ];
deltaR[alpha_,gamma_,r_,sigma_,deltaT_]:=
	alpha*(gamma-r)*deltaT +sigma*nor[0,1]*Sqrt[deltaT]*Sqrt[r];
shortTermR[alpha_,gamma_,r_,sigma_,deltaT_]:=
	r+deltaR[alpha,gamma,r,sigma,deltaT];
SeedRandom[3];
lst=NestList[shortTermR[1,0.04,#,0.02,0.08333]&,0.04,300];
ListPlot[lst,Joined->True,PlotLabel->"CIR Simulation"]

In addition to theVasicek and the Cox-Ingersoll-Ross models of the term structure,there are a variety of other models:
1. Merton 1973 dr=b dt+σ dB
2. Vasicek 1977 dr=α (γ-r) dt+σ dB
3. Dothan 1978 dr=b r dt+σ r dB
4. Cox,Ingersoll,Ross (CIR) 1985 dr=α (γ-r) dt+σ(Sqrt[r])dB
5. Ho,Lee 1986 dr=θ(t) dt+σ dB
6. Hull,White (extended Vasicek) 1990 dr=(θ(t)-βr) dt+σ dB
7. Hull,White (extended CIR) 1990 dr=(θ(t)-βr) dt+σ(Sqrt[r])dB
8. Black,Karasinski 1991 dlog r=(θ(t)-βlog r) dt+σdB

Mathematica program illustrating estimation by generalized method of moments for Vasicek and CIR model of interest rates

Nested Optimization 
Global Optimization
Local Optimization

h[ylist_,t_,c_,β0_,β1_,σ_]:={{v[ylist,t,c,β0,β1,σ]},{(ylist[[t-1]])*v[ylist,t,c,β0,β1,σ]},{(((v[ylist,t,c,β0,β1,σ])^2)-(σ^2))},{(((v[ylist,t,c,β0,β1,σ])^2)-(σ^2))*ylist[[t-1]]}};
v[ylist_,t_,c_,β0_,β1_,σ_]:=ylist[[t]]-ylist[[t-1]]-β0-β1*ylist[[t-1]];
wc[A_,ylist_,t_,c_,β0_,β1_,σ_]:=Min[1,(c/(Norm[A.h[ylist,t,c,β0,β1,σ]]))];
hca[A_,ylist_,t_,c_,β0_,β1_,σ_]:=(A.h[ylist,t,c,β0,β1,σ])*wc[A,ylist,t,c,β0,β1,σ];
c=5.85;
A={{a11,a12,a13,a14},{a21,a22,a23,a24},{a31,a32,a33,a34},{a41,a42,a43,a44}};
f[y_]:=NMinimize[{Flatten[(((∑_(s=2)^(Length[y]))Transpose[hca[A,y,s,c,β0,β1,σ]]).((∑_(u=2)^(Length[y]))hca[A,y,u,c,β0,β1,σ]))][[1]],β0>0,β1<0,σ>0}
,{β0,β1,σ,a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34,a41,a42,a43,a44},MaxIterations->10000]
Clear[β0,β1,σ]
nor[mu_,sigma_]:=Random[NormalDistribution[mu,sigma]];
deltay[y_,β0_,β1_,σ_,deltaT_]:=(β0+β1*y)*deltaT +
	σ*nor[0,1]*Sqrt[deltaT];
shortTermy[y_,β0_,β1_,σ_,deltaT_]:=y+
	deltay[y,β0,β1,σ,deltaT];
SeedRandom[2]
(tt=Table[{NestList[shortTermy[#,0.0012,-0.03,0.01,0.01]&,0.04,10000][[1;;10000;;100]]},{1}];)
f@@@tt

ListPlot[tt[[1]],Joined->True,AxesLabel->{"time","interest rate"},PlotLabel->"Vasicek Simulation"]

Estimating CIR interest rate process by Generalized method of Moments:

h[ylist_,t_,c_,β0_,β1_,σ_]:={{v[ylist,t,c,β0,β1,σ]},{(ylist[[t-1]])*v[ylist,t,c,β0,β1,σ]},{(((v[ylist,t,c,β0,β1,σ])^2)-(σ^2)*ylist[[t-1]])},{(((v[ylist,t,c,β0,β1,σ])^2)-(σ^2)*ylist[[t-1]])*ylist[[t-1]]}};
v[ylist_,t_,c_,β0_,β1_,σ_]:=ylist[[t]]-ylist[[t-1]]-β0-β1*ylist[[t-1]];
wc[A_,ylist_,t_,c_,β0_,β1_,σ_]:=Min[1,(c/(Norm[A.h[ylist,t,c,β0,β1,σ]]))];
hca[A_,ylist_,t_,c_,β0_,β1_,σ_]:=(A.h[ylist,t,c,β0,β1,σ])*wc[A,ylist,t,c,β0,β1,σ];
c=5.85;
A={{a11,a12,a13,a14},{a21,a22,a23,a24},{a31,a32,a33,a34},{a41,a42,a43,a44}};
f[y_]:=NMinimize[{(Flatten[(((∑_(s=2)^(Length[y]))Transpose[hca[A,y,s,c,β0,β1,σ]]).((∑_(u=2)^(Length[y]))hca[A,y,u,c,β0,β1,σ]))][[1]]),σ≥0,β0>0,β1<0}
,{β0,β1,σ,a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34,a41,a42,a43,a44},MaxIterations->10000]
Clear[β0,β1,σ]
nor[mu_,sigma_]:=Random[NormalDistribution[mu,sigma]];
deltay[y_,β0_,β1_,σ_,deltaT_]:=(β0+β1*y)*deltaT +
	σ*nor[0,1]*(Sqrt[y])Sqrt[deltaT];
shortTermy[y_,β0_,β1_,σ_,deltaT_]:=y+
	deltay[y,β0,β1,σ,deltaT];
SeedRandom[2]
(tt=Table[{NestList[shortTermy[#,0.0012,-0.03,0.01,0.01]&,0.05,10000][[1;;10000;;100]]},{1}];)
ListPlot[tt[[1]],Joined->True,AxesLabel->"CIR Simulation"]
f@@@tt

GMM when applied to a general diffusion process in mathematica of the form:
d(Subscript[y, t])=((Subscript[β, 0])+(Subscript[β, 1])(Subscript[y, t])+(Subscript[β, 2])((Subscript[y, t])^2)+((Subscript[β, 3])/(Subscript[y, t])))dt+σ((Subscript[y, t])^γ)d(Subscript[W, t])

h[ylist_,t_,c_,β0_,β1_,β2_,β3_,σ_,γ_]:={{v[ylist,t,c,β0,β1,β2,β3,σ,γ]},{(ylist[[t-1]])*v[ylist,t,c,β0,β1,β2,β3,σ,γ]},{((ylist[[t-1]])^2)*v[ylist,t,c,β0,β1,β2,β3,σ,γ]},{(((v[ylist,t,c,β0,β1,β2,β3,σ,γ])^2)-(σ^2)((ylist[[t-1]])^3))},{(((v[ylist,t,c,β0,β1,β2,β3,σ,γ])^2)-(σ^2)((ylist[[t-1]])^3))*ylist[[t-1]]},{(((v[ylist,t,c,β0,β1,β2,β3,σ,γ])^2)-(σ^2)((ylist[[t-1]])^3))*(ylist[[t-1]]^3)}};
v[ylist_,t_,c_,β0_,β1_,β2_,β3_,σ_,γ_]:=ylist[[t]]-ylist[[t-1]]-β0-β1*ylist[[t-1]]-β2*(ylist[[t-1]])^2;
wc[A_,ylist_,t_,c_,β0_,β1_,β2_,β3_,σ_,γ_]:=Min[1,(c/(Norm[A.h[ylist,t,c,β0,β1,β2,β3,σ,γ]]))];
hca[A_,ylist_,t_,c_,β0_,β1_,β2_,β3_,σ_,γ_]:=(A.h[ylist,t,c,β0,β1,β2,β3,σ,γ])*wc[A,ylist,t,c,β0,β1,β2,β3,σ,γ];
c=5.85;
A={{a11,a12,a13,a14,a15,a16},{a21,a22,a23,a24,a25,a26},{a31,a32,a33,a34,a35,a36},{a41,a42,a43,a44,a45,a46},{a51,a52,a53,a54,a55,a56},{a61,a62,a63,a64,a65,a66}};
f[y_]:=NMinimize[{Flatten[(((∑_(s=2)^(Length[y]))Transpose[hca[A,y,s,c,β0,β1,β2,β3,σ,γ]]).((∑_(u=2)^(Length[y]))hca[A,y,u,c,β0,β1,β2,β3,σ,γ]))][[1]],β0>0,β1<0,σ≥0}
,{β0,β1,β2,β3,σ,γ,a11,a12,a13,a14,a15,a16,a21,a22,a23,a24,a25,a26,a31,a32,a33,a34,a35,a36,a41,a42,a43,a44,a45,a46,a51,a52,a53,a54,a55,a56,a61,a62,a63,a64,a65,a66},MaxIterations->10000];
Clear[β0,β1,β2,β3,σ,γ]
nor[mu_,sigma_]:=Random[NormalDistribution[mu,sigma]];
deltay[y_,β0_,β1_,β2_,β3_,σ_,γ_,deltaT_]:=(β0+β1*y+β2*(y^2)+((β3)/y))*deltaT +
	σ*(y^γ)*nor[0,1]*Sqrt[deltaT];
shortTermy[y_,β0_,β1_,β2_,β3_,σ_,γ_,deltaT_]:=y+
	deltay[y,β0,β1,β2,β3,σ,γ,deltaT];
SeedRandom[2]
(tt=Table[{NestList[shortTermy[#,0.0012,-0.03,-0.08,0.035,0.03,0.5,0.01]&,0.04,10000][[1;;10000;;100]]},{1}];)
f@@@tt