## The function in this file is the R function for testing the corresponding c function
GroupStandardize2=function(x,y)
{
##declaration
##scale each collum of x
n=dim(x)[1]
p=dim(x)[2]
scale=rep(0,p)
tempx=x
for(i in 1: p)
{
scale[i]=sqrt(x[,i]%*%x[,i]/n)
if(scale[i]!=0)
{
tempx[,i]=x[,i]/scale[i]
}
}
## Return list
list(tempx,y,scale)
}
##calculate the value of Xj'Y/n
CrossProduct2=function(x,y,begin,end)
{
n=length(y)
t(x[,begin:end])%*%y/n
}
CrossProductL2Norm2=function(x,y,begin,end)
{
out=CrossProduct2(x,y,begin,end)
sqrt(t(out)%*%out)
}
##find the max ||Xj'Y/n|| for max lamda
MaxProduct2=function(x,y,groupInfo)
{
##declaration
p=length(groupInfo)
n=length(y)
maxVal=0
begin=0
end=0
for(j in 1:p)
{
begin=end+1;
end=end+groupInfo[j];
val= CrossProductL2Norm2(x,y,begin,end);
if(maxVal<val)
{
maxVal=val;
}
}
maxVal
}
GCDReg2=function(x, y, groupInfo, penalty,gamma, lamda, delta, maxIter)
{
##declaration
n=length(y)
m=dim(x)[2]
p=length(groupInfo)
L=length(lamda)
begin=0
end=0
lstart=1 ##since lamda[1] will give all beta 0
##reslut to be returned
beta=matrix(ncol=m,nrow=L,0)
loss=rep(0,L)
iter=rep(0,L)
##temp
r=rep(0,n)
betaPre=rep(0,m)
betaShift=rep(0,m)
z=rep(0,m)
tempb=rep(0,m)
##initial
r=y
for(j in 1:p) ##initial z
{
begin=end+1
end=end+groupInfo[j]
tempz=CrossProduct2(x,y,begin,end) ## division by n included
for(i in begin:end)
{
z[i]=tempz[i-begin+1]
}
}
loss[1]=t(y)%*%y ##initial loss[1]
##iteration for each lamda
for(l in lstart:L)
{
if(l>=2)
{
betaPre=beta[l-1,] ##assign previous beta to betaPre
}
##iteration for all covariates
while(iter[l]<maxIter)
{
iter[l]=iter[l]+1
begin=end=0
##iteration for each covariate group
for(j in 1:p)
{
begin=end+1
end=end+groupInfo[j]
##(1)calculate z j group
tempz=CrossProduct2(x,r,begin,end)
for(i in begin:end)
{
z[i]=tempz[i-begin+1]+betaPre[i]
}
##(2)update beta j group
if (penalty=="MCP")
{
tempb=McPGroup2(z,begin,end,lamda[l],gamma)
beta[l,begin:end]=tempb
}
##(3)update r
betaShift[begin:end]=beta[l,begin:end]-betaPre[begin:end]
r=r-x[,begin:end]%*%betaShift[begin:end]
}
##update betaPre for next iteration
betaPre=beta[l,]
## Check for convergence
if(t(betaShift)%*%betaShift<delta)
{
break;
}
}
##compute square of loss
loss[l]=t(r)%*%r
}
list(beta,loss,iter)
}
McPGroup2=function(z, begin,end,lamda,gamma)
{
##error checking
#if(end<begin||end<0||begin<0)
size=end-begin+1
tempz=rep(0,size)
multiplier=sqrt(size)
znorm=sqrt(t(z[begin:end])%*%z[begin:end])
lamda2=multiplier*lamda
if(znorm<=lamda2) ## set all to zero
{
tempz=rep(0,size)
}
else if(lamda2<znorm && lamda2*gamma>=znorm)
{
s=gamma/(gamma-1)*(1-lamda2/znorm)
tempz=z[begin:end]*s
}
else
{
tempz=z[begin:end]
}
tempz
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.