Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
This post shows how to use constrOptim.nl() R function to solve non-linear optimization problem with or without equality or inequality constraints. Nelson-Siegel yield curve model is used as an target example. Its calculation time is faster than nloptr() function.

# Nelson-Siegel model using constrOptim.nl() R function

In this post, the non-linear least squares problem is solved by using constrOptim.nl() function from alabama R package. Like the previous posts, Nelson-Siegel model is also used for a testing purpose.

For your information, this blog has dealt with Nelson-Siegel model as follows.

### Nelson-Siegel model as a non-linear optimization problem

Nelson-Siegel model is a non-linear least square problem with 4 parameters with some inequality constraints.
\begin{align} y(\tau) &= \beta_1 + \beta_2 \left( \frac{1-e^{- \tau \lambda }}{\tau \lambda }\right) \\ &+ \beta_3 \left(\frac{1-e^{- \tau \lambda }}{\tau \lambda }-e^{- \tau \lambda }\right) \end{align} \begin{align} \beta_1 > 0, \beta_1 + \beta_2 >0, \lambda > 0 \end{align}
Here, $$\tau$$ is a maturity and $$y(\tau)$$ is a continuously compounded spot rate with $$\tau$$ maturity. $$\beta_1, \beta_2, \beta_3$$ are coefficient parameters and $$\lambda$$ is the time decay parameter.

### constrOptim.nl() function

After installing alabama R package, constrOptim.nl() function can be used. Its arguments are so straightforward to understand that I just point out two arguments for constraints. The following two arguments denote inequality and equality constraint functions respectively. In particular, unlike nloptr() function, the inequality constraint function (hin) returns a vector of greater than zero constraints. You can easily find their functional form in the R code which will be shown below

 # hin : a vector function specifying inequality constraints such that hin[j] > 0 for all j # heq : a vector function specifying equality constraints such that heq[j] = 0 for all j

### R code

The following R code is similar to the previous post in that it estimates parameters of Nelson-Siegel model with yield curves at two dates.

 123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123 #========================================================## Quantitative ALM, Financial Econometrics & Derivatives # ML/DL using R, Python, Tensorflow by Sang-Heon Lee ## https://kiandlee.blogspot.com#——————————————————–## Nelson-Siegel using constrOptim.nl function#========================================================# graphics.off(); rm(list = ls()) library(alabama)    #———————————————–# constrOptim.nl objective function#———————————————–    objfunc_co <– function(para, y, m) {        beta <– para[1:3]; la <– para[4]        C  <– cbind(            rep(1,length(m)),                (1–exp(–la*m))/(la*m),                         (1–exp(–la*m))/(la*m)–exp(–la*m))        return(sqrt(mean((y – C%*%beta)^2)))    }    #———————————————–# constrOptim.nl constraint function ( > : larger )#———————————————–    # beta0 + beta1 > 0    constfunc_co <– function(para, y, m) {        beta  <– para[1:3]; la <– para[4]                h <– rep(NA, 4)        h[1] <– beta[1]         # b1 > 0        h[2] <– beta[1]+beta[2] # b1+b2 > 0        h[3] <– la–0.001        # lambda > 0.001        h[4] <– 0.1–la          # lambda < 0.1                           return(h)    } #=======================================================# 1. Read data#=======================================================        # b1, b2, b3, lambda for comparisons    ns_reg_para_rmse1 <– c(        4.26219396, –4.08609206, –4.90893865,          0.02722607,  0.04883786)    ns_reg_para_rmse2 <– c(        4.97628654, –4.75365297, –6.40263059,          0.05046789,  0.04157326)        str.zero <– “        mat     rate1      rate2        3       0.0781     0.0591        6       0.1192     0.0931        9       0.1579     0.1270        12      0.1893     0.1654        24      0.2669     0.3919        36      0.3831     0.8192        48      0.5489     1.3242        60      0.7371     1.7623        72      0.9523     2.1495        84      1.1936     2.4994        96      1.4275     2.7740        108     1.6424     2.9798        120     1.8326     3.1662        144     2.1715     3.4829        180     2.5489     3.7827        240     2.8093     3.9696″        df <– read.table(text = str.zero, header=TRUE)    m  <– df$mat y1 <– df$rate1; y2 <– df$rate2 #==============================================================# 2. constrOptim.nl # : Nonlinear least squares with constraints#============================================================== #——————————————— # NS estimation with 1st data #——————————————— y <– y1 x_init <– c(y[16], y[1]–y[16], 2*y[6]–y[1]–y[16], 0.0609) m_coptim <– constrOptim.nl(par = x_init, y = y, m = m, fn = objfunc_co, hin = constfunc_co, control.optim = list(maxit=50000, trace=0, reltol = 1e–16, abstol = 1e–16) , control.outer = list(eps = 1e–16, itmax = 50000, method = “Nelder-Mead”, NMinit = TRUE)) ns_coptim_out1 <– c(m_coptim$par, m_coptim$value) #——————————————— # NS estimation with 2nd data #——————————————— y <– y2 x_init <– c(y[16], y[1]–y[16], 2*y[6]–y[1]–y[16], 0.0609) m_coptim <– constrOptim.nl(par = x_init, y = y, m = m, fn = objfunc_co, hin = constfunc_co, control.optim = list(maxit=50000, trace=0, reltol = 1e–16, abstol = 1e–16) , control.outer = list(eps = 1e–16, itmax = 50000, method = “Nelder-Mead”, NMinit = TRUE)) ns_coptim_out2 <– c(m_coptim$par, m_coptim\$value)    #=======================================================# 3. Results and Comparisons#=======================================================        ns_reg_para_rmse1    ns_coptim_out1        ns_reg_para_rmse2    ns_coptim_out2    Colored by Color Scripter cs

The following estimation results show that non-linear optimization results from constrOptim.nl() function are nearly similar to the answers.

Using these estimated parameters at two selected dates, two fitted term structure of interest rates can be generated as follows.

### Concluding Remarks

This post shows how to use constrOptim.nl() R function to perform non-linear optimization problem with or without equality or inequality constraints. The reason why I introduce this optimization function is that its computation time tends to be a little faster than nloptr() function. $$\blacksquare$$