Nothing
deriv1LikeMult = function(theta,u1,u2,family,multFactor)
{
side = getSideIfParameterAtBound(theta,family)
result = grad(likeMultFactor,theta, method='simple',u1=u1,u2=u2,family=family,multFactor=multFactor, side=side)
return(result)
}
likeMult = function(params2,v1,v2,family2,params1,u1,u2,family1)
{
nPar1 = getNumbOfParameters(family1)
par1 = getParAsScalars(nPar1,params1)
nPar2 = getNumbOfParameters(family2)
par2 = getParAsScalars(nPar2,params2)
result = mean(log(BiCopPDF(u1,u2,family1,par1[1],par1[2]))*log(BiCopPDF(v1,v2,family2,par2[1],par2[2])))
return(result)
}
likeMultDeriv = function(params1,u1,u2,family1,params2,v1,v2,family2)
{
side = getSideIfParameterAtBound(params2,family2)
result = grad(likeMult,params2, method='simple',v1=v1,v2=v2,family2=family2,params1=params1,u1=u1,u2=u2,family1=family1, side=side)
return(result)
}
deriv2LikeMult = function(params1,u1,u2,family1,params2,v1,v2,family2)
{
side = getSideIfParameterAtBound(params1,family1)
result = jacobian(likeMultDeriv,params1, method='simple',u1=u1,u2=u2,family1=family1,params2=params2,v1=v1,v2=v2,family2=family2, side=side)
return(result)
}
## Asympt with ranks section
likeVar1 = function(u1,u2,family,params)
{
nPar = getNumbOfParameters(family)
par = getParAsScalars(nPar,params)
result = log(BiCopPDF(u1,u2,family,par[1],par[2]))
return(result)
}
likeVar2 = function(u2,u1,family,params)
{
nPar = getNumbOfParameters(family)
par = getParAsScalars(nPar,params)
result = log(BiCopPDF(u1,u2,family,par[1],par[2]))
return(result)
}
computeLikeWithCpitsDerivs = function(parVector, data, svcmDataFrame, cPitData, copulaInd)
{
d <- ncol(data)
nObs = nrow(data)
likeWithCpitsDerivs = array(0, dim = c(nObs, d))
family = svcmDataFrame$family[copulaInd]
if (copulaInd<d)
{
cPit1 = data[,svcmDataFrame$var1[copulaInd]]
cPit2 = data[,svcmDataFrame$var2[copulaInd]]
}
else
{
cPit1 = getCpit1(cPitData, svcmDataFrame, copulaInd)
cPit2 = getCpit2(cPitData, svcmDataFrame, copulaInd)
}
## with ranks
xxSide = vector(length=length(cPit1))
xxSide[] = NA
xxSide[cPit1>0.99] = -1
xxSide[cPit1<0.01] = 1
likeVar1Deriv = grad(likeVar1,cPit1, side = xxSide ,u2=cPit2,family=family,params=parVector, method='simple')
xxSide = vector(length=length(cPit2))
xxSide[] = NA
xxSide[cPit2>0.99] = -1
xxSide[cPit2<0.01] = 1
likeVar2Deriv = grad(likeVar2,cPit2, side = xxSide ,u1=cPit1,family=family,params=parVector, method='simple')
if (copulaInd < d) # first tree
{
likeWithCpitsDerivs[, svcmDataFrame$var1[copulaInd]] = likeVar1Deriv
likeWithCpitsDerivs[, svcmDataFrame$var2[copulaInd]] = likeVar2Deriv
}
else
{
condset = svcmDataFrame$condset[[copulaInd]]
cPit1Deriv = getCpit1Deriv(cPitData, svcmDataFrame, copulaInd, svcmDataFrame$var1[copulaInd])
likeWithCpitsDerivs[, svcmDataFrame$var1[copulaInd]] = likeVar1Deriv * cPit1Deriv
cPit2Deriv = getCpit2Deriv(cPitData, svcmDataFrame, copulaInd, svcmDataFrame$var2[copulaInd])
likeWithCpitsDerivs[, svcmDataFrame$var2[copulaInd]] = likeVar2Deriv * cPit2Deriv
for (condsetVariable in condset)
{
cPit1Deriv = getCpit1Deriv(cPitData, svcmDataFrame, copulaInd, condsetVariable)
cPit2Deriv = getCpit2Deriv(cPitData, svcmDataFrame, copulaInd, condsetVariable)
likeWithCpitsDerivs[, condsetVariable] = likeVar1Deriv * cPit1Deriv + likeVar2Deriv * cPit2Deriv
}
}
return(likeWithCpitsDerivs)
}
helpingfunctionBSspForCovWithRanks = function(parVector, data, svcmDataFrame, cPitData, copulaInd)
{
d <- ncol(data)
nObs = nrow(data)
family = svcmDataFrame$family[copulaInd]
w = array(0, dim = c(nObs, d))
likeWithCpitsDerivs = computeLikeWithCpitsDerivs(parVector, data, svcmDataFrame, cPitData, copulaInd)
orderingInds = apply(data,2,order, decreasing=TRUE)
for (iVar in 1:d)
{
xx = likeWithCpitsDerivs[orderingInds[,iVar], iVar]
w[orderingInds[,iVar], iVar] = cumsum(xx)/nObs
}
sumOfW = apply(w,1,sum)
return(sumOfW)
}
scoresForBSspWithRanks = function(params, data, svcmDataFrame, cPitData, copulaInd)
{
family = svcmDataFrame$family[copulaInd]
d <- ncol(data)
if (copulaInd<d)
{
cPit1 = data[,svcmDataFrame$var1[copulaInd]]
cPit2 = data[,svcmDataFrame$var2[copulaInd]]
}
else
{
cPit1 = getCpit1(cPitData, svcmDataFrame, copulaInd)
cPit2 = getCpit2(cPitData, svcmDataFrame, copulaInd)
}
result = likeVar1(cPit1, cPit2, family, params)
return(result)
}
likeMultWithRanks = function(params2,copulaInd2, params1,copulaInd1, data, svcmDataFrame, cPitData)
{
xx1 = helpingfunctionBSspForCovWithRanks(params1, data, svcmDataFrame, cPitData, copulaInd1)
xx2 = scoresForBSspWithRanks(params1, data, svcmDataFrame, cPitData, copulaInd1)
yy1 = helpingfunctionBSspForCovWithRanks(params2, data, svcmDataFrame, cPitData, copulaInd2)
yy2 = scoresForBSspWithRanks(params2, data, svcmDataFrame, cPitData, copulaInd2)
cov1 = mean((xx1-mean(xx1))*(yy1-mean(yy1)))
cov2 = mean((xx2-mean(xx2))*(yy1-mean(yy1)))
cov3 = mean((xx1-mean(xx1))*(yy2-mean(yy2)))
result = cov1 + cov2 + cov3
return(result)
}
likeMultDerivWithRanks = function(params1,copulaInd1, params2,copulaInd2, data, svcmDataFrame, cPitData)
{
xxFamily2 = svcmDataFrame$family[copulaInd2]
side = getSideIfParameterAtBound(params2,xxFamily2)
result = grad(likeMultWithRanks,params2, method='simple',
copulaInd2=copulaInd2,
params1=params1,copulaInd1=copulaInd1,
data=data, svcmDataFrame=svcmDataFrame, cPitData=cPitData, side=side)
return(result)
}
deriv2LikeMultWithRanks = function(params1,copulaInd1, params2,copulaInd2, data, svcmDataFrame, cPitData)
{
xxFamily1 = svcmDataFrame$family[copulaInd1]
side = getSideIfParameterAtBound(params1,xxFamily1)
result = jacobian(likeMultDerivWithRanks,params1, method='simple',
copulaInd1=copulaInd1,
params2=params2,copulaInd2=copulaInd2,
data=data, svcmDataFrame=svcmDataFrame, cPitData=cPitData, side=side)
return(result)
}
bSspForCovWithRanks = function(data, svcmDataFrame, cPitData, includeLastCopula = FALSE)
{
d <- ncol(data)
nObs = nrow(data)
if (includeLastCopula)
{
nCopulas = d*(d-1)/2
}
else
{
nCopulas = d*(d-1)/2-1
}
nParameters = sum(svcmDataFrame$nPar[1:nCopulas])
bSsp = matrix(0,nrow=nParameters,ncol=nParameters)
for (jCopula in 1:nCopulas)
{
if (svcmDataFrame$nPar[jCopula])
{
parVector1 = svcmDataFrame$par[[jCopula]]
for (lCopula in 1:jCopula)
{
if (svcmDataFrame$nPar[lCopula]) #&&
#((lCopula == jCopula) || !any(svcmDataFrame$copulaInd[lCopula]==c(svcmDataFrame$cPit1CopulaInd[[jCopula]],svcmDataFrame$cPit2CopulaInd[[jCopula]]))))
{
parVector2 = svcmDataFrame$par[[lCopula]]
bSsp[svcmDataFrame$parInd[[lCopula]], svcmDataFrame$parInd[[jCopula]]] =
deriv2LikeMultWithRanks(parVector1,jCopula,
parVector2,lCopula,
data, svcmDataFrame, cPitData)
}
}
}
}
xx = t(bSsp)
bSsp[lower.tri(bSsp)] = xx[lower.tri(xx)]
return(bSsp)
}
getOmegaWithLikesD = function(data, svcmDataFrame, cPitData, includeLastCopula = FALSE)
{
d <- ncol(data)
nObs = nrow(data)
if (includeLastCopula)
{
nCopulas = d*(d-1)/2
}
else
{
nCopulas = d*(d-1)/2-1
}
nParameters = sum(svcmDataFrame$nPar[1:nCopulas])
D = matrix(0,nrow=nParameters,ncol=nParameters)
for (jCopula in 1:nCopulas)
{
if (svcmDataFrame$nPar[jCopula])
{
family_1 = svcmDataFrame$family[jCopula]
par_1 = svcmDataFrame$par[[jCopula]]
if (jCopula<d)
{
cPit1_1 = data[,svcmDataFrame$var1[jCopula]]
cPit2_1 = data[,svcmDataFrame$var2[jCopula]]
}
else
{
cPit1_1 = getCpit1(cPitData, svcmDataFrame, jCopula)
cPit2_1 = getCpit2(cPitData, svcmDataFrame, jCopula)
}
for (lCopula in 1:jCopula)
{
if (svcmDataFrame$nPar[lCopula])
#((lCopula == jCopula) || !any(svcmDataFrame$copulaInd[lCopula]==c(svcmDataFrame$cPit1CopulaInd[[jCopula]],svcmDataFrame$cPit2CopulaInd[[jCopula]]))))
{
family_2 = svcmDataFrame$family[lCopula]
par_2 = svcmDataFrame$par[[lCopula]]
if (lCopula<d)
{
cPit1_2 = data[,svcmDataFrame$var1[lCopula]]
cPit2_2 = data[,svcmDataFrame$var2[lCopula]]
}
else
{
cPit1_2 = getCpit1(cPitData, svcmDataFrame, lCopula)
cPit2_2 = getCpit2(cPitData, svcmDataFrame, lCopula)
}
D[svcmDataFrame$parInd[[lCopula]],svcmDataFrame$parInd[[jCopula]]] = deriv2LikeMult(par_1,cPit1_1,cPit2_1,family_1,
par_2,cPit1_2,cPit2_2,family_2)
}
}
}
}
xx = t(D)
D[lower.tri(D)] = xx[lower.tri(xx)]
return(D)
}
computeEstEqWithCpitsDerivs = function(data, cPitData, svcmDataFrame, copulaInd, indGroup, theta)
{
d <- ncol(data)
nObs = nrow(data)
cPit1 = getCpit1(cPitData, svcmDataFrame, copulaInd)
cPit2 = getCpit2(cPitData, svcmDataFrame, copulaInd)
# Obtain the subsample
cPit1InGroup = cPit1[indGroup]
cPit2InGroup = cPit2[indGroup]
nObsPerGroup = length(cPit1InGroup)
mu1WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
mu2WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
sigma1WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
sigma2WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
rhoWithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
# Obtain the estimated parameters
mu1 = theta[1]
sigma2_1 = theta[2]
mu2 = theta[3]
sigma2_2 = theta[4]
rho = theta[5]
condset = svcmDataFrame$condset[[copulaInd]]
cPit1Deriv = getCpit1Deriv(cPitData, svcmDataFrame, copulaInd, svcmDataFrame$var1[copulaInd])
cPit1DerivInGroup = cPit1Deriv[indGroup]
mu1WithCpitsDerivs[, svcmDataFrame$var1[copulaInd]] = -cPit1DerivInGroup
sigma1WithCpitsDerivs[, svcmDataFrame$var1[copulaInd]] = -2*(cPit1InGroup-mu1)*cPit1DerivInGroup
rhoWithCpitsDerivs[, svcmDataFrame$var1[copulaInd]] = -cPit1DerivInGroup * (cPit2InGroup - mu2) / sqrt(sigma2_1*sigma2_2)
cPit2Deriv = getCpit2Deriv(cPitData, svcmDataFrame, copulaInd, svcmDataFrame$var2[copulaInd])
cPit2DerivInGroup = cPit2Deriv[indGroup]
mu2WithCpitsDerivs[, svcmDataFrame$var2[copulaInd]] = -cPit2DerivInGroup
sigma2WithCpitsDerivs[, svcmDataFrame$var2[copulaInd]] = -2*(cPit2InGroup-mu2)*cPit2DerivInGroup
rhoWithCpitsDerivs[, svcmDataFrame$var2[copulaInd]] = -cPit2DerivInGroup * (cPit1InGroup - mu1) / sqrt(sigma2_1*sigma2_2)
for (condsetVariable in condset)
{
cPit1Deriv = getCpit1Deriv(cPitData, svcmDataFrame, copulaInd, condsetVariable)
cPit2Deriv = getCpit2Deriv(cPitData, svcmDataFrame, copulaInd, condsetVariable)
cPit1DerivInGroup = cPit1Deriv[indGroup]
cPit2DerivInGroup = cPit2Deriv[indGroup]
mu1WithCpitsDerivs[, condsetVariable] = -cPit1DerivInGroup
sigma1WithCpitsDerivs[, condsetVariable] = -2*(cPit1InGroup-mu1)*cPit1DerivInGroup
mu2WithCpitsDerivs[, condsetVariable] = -cPit2DerivInGroup
sigma2WithCpitsDerivs[, condsetVariable] = -2*(cPit2InGroup-mu2)*cPit2DerivInGroup
rhoWithCpitsDerivs[, condsetVariable] = (-cPit1DerivInGroup * (cPit2InGroup - mu2) - cPit2DerivInGroup * (cPit1InGroup - mu1)) / sqrt(sigma2_1*sigma2_2)
}
return(list(mu1WithCpitsDerivs = mu1WithCpitsDerivs, sigma1WithCpitsDerivs = sigma1WithCpitsDerivs,
mu2WithCpitsDerivs = mu2WithCpitsDerivs, sigma2WithCpitsDerivs = sigma2WithCpitsDerivs,
rhoWithCpitsDerivs = rhoWithCpitsDerivs))
}
helpingfunctionBSspForEqualCorrWithRanks = function(data, cPitData, svcmDataFrame, copulaInd, indGroup, theta)
{
d <- ncol(data)
nObs = nrow(data)
dataInGroup = data[indGroup,]
nObsPerGroup = length(indGroup)
wMu1WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
wMu2WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
wSigma1WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
wSigma2WithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
wRhoWithCpitsDerivs = array(0, dim = c(nObsPerGroup, d))
result = array(0, dim = c(nObsPerGroup, 5))
estEqWithCpitsDerivs = computeEstEqWithCpitsDerivs(data, cPitData, svcmDataFrame, copulaInd, indGroup, theta)
orderingInds = apply(dataInGroup,2,order, decreasing=TRUE)
for (iVar in 1:d)
{
xx = estEqWithCpitsDerivs$mu1WithCpitsDerivs[orderingInds[,iVar], iVar]
wMu1WithCpitsDerivs[orderingInds[,iVar], iVar] = cumsum(xx)/nObs
xx = estEqWithCpitsDerivs$mu2WithCpitsDerivs[orderingInds[,iVar], iVar]
wMu2WithCpitsDerivs[orderingInds[,iVar], iVar] = cumsum(xx)/nObs
xx = estEqWithCpitsDerivs$sigma1WithCpitsDerivs[orderingInds[,iVar], iVar]
wSigma1WithCpitsDerivs[orderingInds[,iVar], iVar] = cumsum(xx)/nObs
xx = estEqWithCpitsDerivs$sigma2WithCpitsDerivs[orderingInds[,iVar], iVar]
wSigma2WithCpitsDerivs[orderingInds[,iVar], iVar] = cumsum(xx)/nObs
xx = estEqWithCpitsDerivs$rhoWithCpitsDerivs[orderingInds[,iVar], iVar]
wRhoWithCpitsDerivs[orderingInds[,iVar], iVar] = cumsum(xx)/nObs
}
result[,1] = apply(wMu1WithCpitsDerivs,1,sum)
result[,3] = apply(wMu2WithCpitsDerivs,1,sum)
result[,2] = apply(wSigma1WithCpitsDerivs,1,sum)
result[,4] = apply(wSigma2WithCpitsDerivs,1,sum)
result[,5] = apply(wRhoWithCpitsDerivs,1,sum)
return(result)
}
bsspEstEqEqualCorrSingleGroup = function(data, svcmDataFrame, cPitData, copulaInd, indGroup, theta, aInGroup, bInGroup)
{
nObs = nrow(data)
nObsPerGroup = length(indGroup)
lambdaInGroup = nObsPerGroup/nObs
xx = helpingfunctionBSspForEqualCorrWithRanks(data, cPitData, svcmDataFrame, copulaInd, indGroup, theta)
yy = cbind(aInGroup, bInGroup)
result = matrix(0, 5, 5)
for (iEstEq in 1:5)
{
for (jEstEq in iEstEq:5)
{
cov1 = mean((xx[,iEstEq]-mean(xx[,iEstEq]))*(xx[,jEstEq]-mean(xx[,jEstEq])))
cov2 = mean((yy[,iEstEq]-mean(yy[,iEstEq]))*(xx[,jEstEq]-mean(xx[,jEstEq])))
cov3 = mean((xx[,iEstEq]-mean(xx[,iEstEq]))*(yy[,jEstEq]-mean(yy[,jEstEq])))
result[jEstEq, iEstEq] = (cov1 + cov2 + cov3)/lambdaInGroup
}
}
xx = t(result)
result[upper.tri(result)] = xx[upper.tri(xx)]
return(result)
}
bsspEstEqEqualCorr = function(data, svcmDataFrame, indList, cPitData, copulaInd, theta, listOfMultipliers)
{
nObs = nrow(data)
# Obtain the subsamples
nGroups = length(indList)
bSsp = matrix(0,nrow = 4*nGroups+nGroups,ncol = 4*nGroups+nGroups)
for (iGroup in 1:nGroups)
{
aInGroup = listOfMultipliers$aInGroups[[iGroup]]
bInGroup = listOfMultipliers$bInGroups[[iGroup]]
# Obtain the estimated parameters
thetaInGroup = c(theta[(1 + 4*(iGroup-1)) : (4 + 4*(iGroup-1))], theta[(4*nGroups+iGroup)])
indGroup = indList[[iGroup]]
xx = bsspEstEqEqualCorrSingleGroup(data, svcmDataFrame, cPitData, copulaInd, indGroup, thetaInGroup, aInGroup, bInGroup)
thisGroupStartIndMuSigma = 4*(iGroup-1)
thisGroupMuSigmaIndexSet = (1 + thisGroupStartIndMuSigma):(4 + thisGroupStartIndMuSigma)
thisGroupRhoIndex = 4*nGroups + iGroup
thisGroupIndexSet = c(thisGroupMuSigmaIndexSet, thisGroupRhoIndex)
bSsp[thisGroupIndexSet, thisGroupIndexSet] = xx
}
xx = t(bSsp)
bSsp[upper.tri(bSsp)] = xx[upper.tri(xx)]
return(bSsp)
}
likeMultCrossTermWithRanks = function(params1,copulaInd1, data, svcmDataFrame, cPitData, yyEstEqWithCpitsDeriv, yyEstEq, indGroup)
{
nObs = nrow(data)
nObsPerGroup = length(indGroup)
lambdaInGroup = nObsPerGroup/nObs
xx1 = helpingfunctionBSspForCovWithRanks(params1, data, svcmDataFrame, cPitData, copulaInd1)
xx2 = scoresForBSspWithRanks(params1, data, svcmDataFrame, cPitData, copulaInd1)
xx1 = xx1[indGroup]
xx2 = xx2[indGroup]
cov1 = mean((xx1-mean(xx1))*(yyEstEqWithCpitsDeriv-mean(yyEstEqWithCpitsDeriv)))
cov2 = mean((xx2-mean(xx2))*(yyEstEqWithCpitsDeriv-mean(yyEstEqWithCpitsDeriv)))
cov3 = mean((xx1-mean(xx1))*(yyEstEq-mean(yyEstEq)))
result = (cov1 + cov2 + cov3)
return(result)
}
deriv1LikeMultWithRanks = function(params1,copulaInd1, data, svcmDataFrame, cPitData, yyEstEqWithCpitsDeriv, yyEstEq, indGroup)
{
xxFamily1 = svcmDataFrame$family[copulaInd1]
side = getSideIfParameterAtBound(params1,xxFamily1)
result = grad(likeMultCrossTermWithRanks,params1, method='simple',
copulaInd1=copulaInd1,
data=data, svcmDataFrame=svcmDataFrame, cPitData=cPitData,
yyEstEqWithCpitsDeriv = yyEstEqWithCpitsDeriv, yyEstEq = yyEstEq, indGroup = indGroup, side=side)
return(result)
}
bsspEstEqEqualCorrCrossTermsSinglegroup = function(data, svcmDataFrame, cPitData, copulaInd, indGroup, theta, aInGroup, bInGroup)
{
d <- ncol(data)
nObs = nrow(data)
nCopulas = d*(d-1)/2-1
nParameters = sum(svcmDataFrame$nPar[1:nCopulas])
bSsp = matrix(0,nrow=nParameters,ncol=5)
yyEstEqWithCpitsDeriv = helpingfunctionBSspForEqualCorrWithRanks(data, cPitData, svcmDataFrame, copulaInd, indGroup, theta)
yyEstEq = cbind(aInGroup, bInGroup)
for (jCopula in 1:nCopulas)
{
if (svcmDataFrame$nPar[jCopula])
{
parVector1 = svcmDataFrame$par[[jCopula]]
for (iTheta in 1:5)
{
bSsp[svcmDataFrame$parInd[[jCopula]],iTheta] =
deriv1LikeMultWithRanks(parVector1,jCopula, data, svcmDataFrame, cPitData,
yyEstEqWithCpitsDeriv[,iTheta], yyEstEq[,iTheta], indGroup)
}
}
}
return(bSsp)
}
bsspEstEqEqualCorrCrossTerms = function(data, svcmDataFrame, indList, cPitData, copulaInd, theta, listOfMultipliers)
{
d <- ncol(data)
nObs = nrow(data)
nCopulas = d*(d-1)/2-1
nParameters = sum(svcmDataFrame$nPar[1:nCopulas])
# Obtain the subsamples
nGroups = length(indList)
bSsp = matrix(0,nrow=4*nGroups+nGroups,ncol=nParameters)
for (iGroup in 1:nGroups)
{
aInGroup = listOfMultipliers$aInGroups[[iGroup]]
bInGroup = listOfMultipliers$bInGroups[[iGroup]]
# Obtain the estimated parameters
thetaInGroup = c(theta[(1 + 4*(iGroup-1)) : (4 + 4*(iGroup-1))], theta[(4*nGroups+iGroup)])
indGroup = indList[[iGroup]]
xx = bsspEstEqEqualCorrCrossTermsSinglegroup(data, svcmDataFrame, cPitData, copulaInd, indGroup, thetaInGroup, aInGroup, bInGroup)
thisGroupStartIndMuSigma = 4*(iGroup-1)
thisGroupMuSigmaIndexSet = (1 + thisGroupStartIndMuSigma):(4 + thisGroupStartIndMuSigma)
thisGroupRhoIndex = 4*nGroups + iGroup
thisGroupIndexSet = c(thisGroupMuSigmaIndexSet, thisGroupRhoIndex)
bSsp[thisGroupIndexSet, ] = t(xx)
}
return(bSsp)
}
getOmegaWithLikesE = function(data, svcmDataFrame, indList, cPitData, listOfMultipliers)
{
d <- ncol(data)
nObs = nrow(data)
nCopulas = d*(d-1)/2-1
nParameters = sum(svcmDataFrame$nPar[1:nCopulas])
nGroups = length(indList)
E = matrix(0,nrow=4*nGroups+nGroups,ncol=nParameters)
for (jCopula in 1:nCopulas)
{
if (svcmDataFrame$nPar[jCopula])
{
family_1 = svcmDataFrame$family[jCopula]
par_1 = svcmDataFrame$par[[jCopula]]
if (jCopula<d)
{
cPit1_1 = data[,svcmDataFrame$var1[jCopula]]
cPit2_1 = data[,svcmDataFrame$var2[jCopula]]
}
else
{
cPit1_1 = getCpit1(cPitData, svcmDataFrame, jCopula)
cPit2_1 = getCpit2(cPitData, svcmDataFrame, jCopula)
}
for (iGroup in 1:nGroups)
{
cPit1InGroup = cPit1_1[indList[[iGroup]]]
cPit2InGroup = cPit2_1[indList[[iGroup]]]
nObsPerGroup = length(cPit1InGroup)
lambdaInGroup = nObsPerGroup/nObs
E[1 + 4*(iGroup-1),svcmDataFrame$parInd[[jCopula]]] = deriv1LikeMult(par_1,cPit1InGroup,cPit2InGroup,family_1,listOfMultipliers$aInGroups[[iGroup]][,1])
E[2 + 4*(iGroup-1),svcmDataFrame$parInd[[jCopula]]] = deriv1LikeMult(par_1,cPit1InGroup,cPit2InGroup,family_1,listOfMultipliers$aInGroups[[iGroup]][,2])
E[3 + 4*(iGroup-1),svcmDataFrame$parInd[[jCopula]]] = deriv1LikeMult(par_1,cPit1InGroup,cPit2InGroup,family_1,listOfMultipliers$aInGroups[[iGroup]][,3])
E[4 + 4*(iGroup-1),svcmDataFrame$parInd[[jCopula]]] = deriv1LikeMult(par_1,cPit1InGroup,cPit2InGroup,family_1,listOfMultipliers$aInGroups[[iGroup]][,4])
E[(4*nGroups + iGroup),svcmDataFrame$parInd[[jCopula]]] = deriv1LikeMult(par_1,cPit1InGroup,cPit2InGroup,family_1,listOfMultipliers$bInGroups[[iGroup]])
}
}
}
return(E)
}
omegaMuSigmaRhoInSingleGroup = function(cPit1InGroup, cPit2InGroup, theta, nObs)
{
omega = matrix(0,nrow = 5,ncol = 5)
# Obtain the estimated parameters
mu1 = theta[1]
var1 = theta[2]
mu2 = theta[3]
var2 = theta[4]
rho = theta[5]
nObsPerGroup = length(cPit1InGroup)
lambdaInGroup = nObsPerGroup/nObs
bInGroups = rho - (cPit1InGroup - mu1) * (cPit2InGroup - mu2) / sqrt(var1*var2)
aInGroups = cbind(mu1-cPit1InGroup,
var1-(cPit1InGroup-mu1)^2,
mu2-cPit2InGroup,
var2-(cPit2InGroup-mu2)^2)
omega[1:4, 1:4] = 1/nObsPerGroup *(t(aInGroups) %*% aInGroups)/lambdaInGroup
omega[5, 1:4] = 1/nObsPerGroup *(t(aInGroups) %*% bInGroups)/lambdaInGroup
omega[5, 5] = mean(bInGroups^2)/lambdaInGroup
return(list(aInGroups = aInGroups, bInGroups = bInGroups, omega = omega))
}
omegaMuSigmaRho = function(data, indList, cPit1, cPit2, theta)
{
nObs = nrow(data)
# Obtain the subsamples
nGroups = length(indList)
omega = matrix(0,nrow = 4*nGroups+nGroups,ncol = 4*nGroups+nGroups)
aInGroups = vector("list",nGroups)
bInGroups = vector("list",nGroups)
for (iGroup in 1:nGroups)
{
# Obtain the subsample
cPit1InGroup = cPit1[indList[[iGroup]]]
cPit2InGroup = cPit2[indList[[iGroup]]]
# Obtain the estimated parameters
thetaInGroup = c(theta[(1 + 4*(iGroup-1)) : (4 + 4*(iGroup-1))], theta[(4*nGroups+iGroup)])
xx = omegaMuSigmaRhoInSingleGroup(cPit1InGroup, cPit2InGroup, thetaInGroup, nObs)
bInGroups[[iGroup]] = xx$bInGroups
aInGroups[[iGroup]] = xx$aInGroups
thisGroupStartIndMuSigma = 4*(iGroup-1)
thisGroupMuSigmaIndexSet = (1 + thisGroupStartIndMuSigma):(4 + thisGroupStartIndMuSigma)
thisGroupRhoIndex = 4*nGroups + iGroup
thisGroupIndexSet = c(thisGroupMuSigmaIndexSet, thisGroupRhoIndex)
omega[thisGroupIndexSet, thisGroupIndexSet] = xx$omega
}
listOfMultipliers = list(aInGroups=aInGroups,bInGroups=bInGroups)
xx = t(omega)
omega[upper.tri(omega)] = xx[upper.tri(xx)]
return(list(omega =omega, listOfMultipliers = listOfMultipliers))
}
omegaRvine = function(data, svcmDataFrame, indList, cPitData, theta, withRanks)
{
d <- ncol(data)
nObs = nrow(data)
nCopulas = d*(d-1)/2-1
nParameters = sum(svcmDataFrame$nPar[1:nCopulas])
nParametersFirstTree = sum(svcmDataFrame$nPar[1:(d-1)])
# Obtain the cPits to be tested
copulaInd = nrow(svcmDataFrame)
cPit1 = getCpit1(cPitData, svcmDataFrame, copulaInd)
cPit2 = getCpit2(cPitData, svcmDataFrame, copulaInd)
# Obtain the subsamples
nGroups = length(indList)
omega = matrix(0,nrow=nParameters+4*nGroups+nGroups,ncol=nParameters+4*nGroups+nGroups)
xx = omegaMuSigmaRho(data, indList, cPit1, cPit2, theta)
omegaCorrelations = xx$omega
listOfMultipliers = xx$listOfMultipliers
if(withRanks)
{
bSsp = bsspEstEqEqualCorr(data, svcmDataFrame, indList, cPitData, copulaInd, theta, listOfMultipliers)
omegaCorrelations = omegaCorrelations + bSsp
}
omega[(nParameters+1):(nParameters+5*nGroups), (nParameters+1):(nParameters+5*nGroups)] = omegaCorrelations
if (nParameters)
{
omegaWithLikesD = getOmegaWithLikesD(data, svcmDataFrame, cPitData)
omegasWithLikesE = getOmegaWithLikesE(data, svcmDataFrame, indList, cPitData, listOfMultipliers)
if(withRanks)
{
bSsp = bSspForCovWithRanks(data, svcmDataFrame, cPitData)
omegaWithLikesD = omegaWithLikesD + bSsp
bSsp = bsspEstEqEqualCorrCrossTerms(data, svcmDataFrame, indList, cPitData, copulaInd, theta, listOfMultipliers)
omegasWithLikesE = omegasWithLikesE + bSsp
}
omega[1:nParameters,1:nParameters] = omegaWithLikesD
omega[(nParameters+1):(nParameters+4*nGroups+nGroups),1:nParameters] = omegasWithLikesE
}
xx = t(omega)
omega[upper.tri(omega)] = xx[upper.tri(xx)]
return(omega)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.