DCC composite likelihood

Hi everyone,

I want to estimate a standard DCC model by maximizing a composite 
likelihood function by optimizing over the function "comp" as shown in 
the code below. I want to get the optimal parameters alpha and beta for 
a set of 29 assets, i.e. 406 pairwise correlations. Each asset has 1019 
observations.

So what I am trying is basically the double for-loop in the function. 
Calculating the pairwise rho for each point in time and for each pair. 
CL here denotes the composite likelihood. I only want to have one alpha 
and one beta. The functional form is shown in eq. 2.8 on page 9/44 in 
the following paper: http://fic.wharton.upenn.edu/fic/papers/11/11-20.pdf

I guess it must be the for-loop since my optimization does not come to 
an end.

Thanks in advance.


----------------------------------------------------------------------------------------------------------------

library("stats")
library("rugarch")
library("moments")
library("tseries")
library("fGarch")
library("nlme")
library("fOptions")
library("NMOF")

res <- matrix(rnorm(1019*29, 1019, 29)   ## artificial data

comb <- combn(1:29, 2, FUN = NULL, simplify = TRUE)
comb.sq <- matrix(,1019,406)
comb.cr <- matrix(,1019,406)
for(j in 1:406) { comb.sq[,j] <- res[,comb[1,j]]^2 + res[,comb[2,j]]^2
                   comb.cr[,j] <- res[,comb[1,j]]*res[,comb[2,j]] 
}         ## only data manipulation

comp <- function(pars) {

alpha <- pars[1]
beta <- pars[2]
q11 <- matrix(,1019,406)
q12 <- matrix(,1019,406)
q22 <- matrix(,1019,406)
rho <- matrix(,1019,406)
rho12b <- matrix(,1,406)
like <- matrix(, 1019, 406)
CL <- NULL
comb <- combn(1:29, 2, FUN = NULL, simplify = TRUE)

q11[1,] <- 1         ## starting values for each qii, qij
q22[1,] <- 1
q12[1,] <- apply(comb.cr, 2, mean)
rho12b[1,] <- apply(comb.cr, 2, mean)
rho[1,] <- q12[1,]/sqrt(q11[1,]*q22[1,])

for(j in 1:406)  {     ## calculating time-varying pairwise correlation 
for each pair j and each point in time i
for(i in 2:1019) {


q11[i,j] <-  1 + alpha*(res[i-1,comb[1,j]]^2-1)+beta*(q11[i-1,j]-1)

q12[i,j] <- rho12b[1,j] + alpha*(comb.cr[i-1,j]-rho12b[1,j]) + 
beta*(q12[i-1,j]-rho12b[1,j])

q22[i,j] <- 1 + alpha*(res[i-1,comb[2,j]]^2-1) + beta*(q22[i-1,j]-1)

rho[i,j] <- q12[i,j]/sqrt(q11[i,j]*q22[i,j])

CL <- c(CL, 0.5*(log(1-rho[i,j]^2) + (comb.sq[i,j] - 
2*rho[i,j]*comb.cr[i,j])/(1-rho[i,j]^2)))

                   }
                   }
CL <- sum(CL)
return(CL)
                    }



const1 <- matrix(c(1,0,-1,0,1,-1), 3, 2)   # constraints for nonlinear 
optimization
const2 <- c(0.0001, 0.0001, -0.999)
constrOptim(c(0.1, 0.85), comp, NULL, ui=const1, ci=const2)


---------------------------------------------------------------------------------------------------

R version 2.14.1 (2011-12-22)
Platform: i386-pc-mingw32/i386 (32-bit)

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United 
States.1252    LC_MONETARY=English_United States.1252
[4] LC_NUMERIC=C                           LC_TIME=English_United 
States.1252

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods base

other attached packages:
  [1] NMOF_0.25-4           fOptions_2160.81 nlme_3.1-104          
fGarch_2110.80        fBasics_2160.85
  [6] MASS_7.3-16           timeSeries_2160.95 timeDate_2160.95      
tseries_0.10-29       zoo_1.7-7
[11] quadprog_1.5-4        moments_0.13 rugarch_1.0-11        
Rsolnp_1.12           truncnorm_1.0-5
[16] chron_2.3-42          numDeriv_2012.3-1 RcppArmadillo_0.3.4.0 
Rcpp_0.9.10           R.utils_1.16.2
[21] R.oo_1.9.9            R.methodsS3_1.4.2

loaded via a namespace (and not attached):
[1] grid_2.14.1      lattice_0.20-10  stabledist_0.6-4 tools_2.14.1