ADMM example

scaled dual ADMM algorithm

• Given x, z and , to some initial value.

• Repeat:

  • $x:=\arg \underset{x}{max}(f(x)+\frac{\rho }{2}\parallel Ax+Bz-c+\mu {\parallel }_2^2)$
  • $x:=\arg \underset{x}{max}(g(z)+\frac{\rho }{2}\parallel Ax+Bz-c+\mu {\parallel }_2^2)$
  • $\mu :=\mu +(Ax+Bz-c)$
  • Stopping criterion : quit $\parallel r{\parallel }_2<ϵ$ and $\parallel s{\parallel }_2<ϵ$
    
    
    

Stopping criterion

We can define the primal and dual residuals in ADMM at step k+1.
* Primal residuals : $r^{k+1}=A{x}^{k+1}+B{z}^{k+1}-c$
* Dual residuals : $s^{k+1}=\rho A^{T}B(z^{k+1}-z^k)$

Therefore stopping criterion satisfies that $\parallel r{\parallel }_2$ and $\parallel s{\parallel }_2$ are smaller than any $ϵ$

Lasso example

$$\begin{array}{rlrl}& \underset{\beta }{\text{minimize}}& & \sum _{i=1}^n(y_i-{\beta }_0-x_i^t\beta {)}^2+\lambda \sum _{j=1}^p|{\beta }_j|\end{array}$$

$\iff$

$$\begin{array}{rlrl}& \underset{\beta }{\text{minimize}}& & f(\beta )+f(z)\\ & \text{subject to}& & I\beta -IZ=0\end{array}$$

$$\begin{array}{r}r=I\beta -IZ\\ L_{\rho }(\beta ,z,v)& =f(\beta )+g(z)+v^{t}r+\frac{\rho }{2}||r|{|}_2^2\\ & =f(\beta )+g(z)+\frac{\rho }{2}||r+\frac{1}{\rho }v|{|}_2^2-\frac{\rho }{2}||v|{|}_2^2\\ & =f(\beta )+g(z)+\frac{\rho }{2}||r+\mu |{|}_2^2-constan{t}_v,\mu =\frac{1}{\rho }v\end{array}$$

$$\begin{array}{rl}{\beta }^{k+1}& :=\underset{\beta }{\mathrm{argmin}}(f(\beta )+\frac{\rho }{2}||I\beta -I{Z}^k+{\mu }^k|{|}_2^2)\\ & =\underset{\beta }{\mathrm{argmin}}(y-X\beta {)}^t(y-X\beta )+\frac{\rho }{2}||I\beta -I{Z}^k+{\mu }^k|{|}_2^2)\\ & \Rightarrow -2X^{t}y+2X^{t}X\beta +\rho \beta -\rho Z^k+\rho {\mu }^k=0\\ & \iff (2X^{t}X+\rho I)\beta =2X^{t}y+\rho (Z^k-{\mu }^k)\\ & \therefore {\beta }^{k+1}=(2X^{t}X+\rho I{)}^{-1}(2X^{t}y+\rho (Z^k-{\mu }^k))\end{array}$$

$$\begin{array}{rl}{\beta }^{k+1}& :=\underset{\beta }{\mathrm{argmin}}(f(\beta )+\frac{\rho }{2}||I\beta -I{Z}^k+{\mu }^k|{|}_2^2)\\ & =\underset{\beta }{\mathrm{argmin}}(y-X\beta {)}^t(y-X\beta )+\frac{\rho }{2}||I\beta -I{Z}^k+{\mu }^k|{|}_2^2)\\ & \Rightarrow -2X^{t}y+2X^{t}X\beta +\rho \beta -\rho Z^k+\rho {\mu }^k=0\\ & \iff (2X^{t}X+\rho I)\beta =2X^{t}y+\rho (Z^k-{\mu }^k)\\ & \therefore {\beta }^{k+1}=(2X^{t}X+\rho I{)}^{-1}(2X^{t}y+\rho (Z^k-{\mu }^k))\end{array}$$

The prox operatior for $g(z)=\lambda ||z|{|}_1$

$$\begin{array}{rl}pro{x}_{\lambda ,g}(z)& =\underset{v}{\mathrm{argmin}}(\lambda ||z|{|}_1+\frac{1}{2}||z-v|{|}_2^2)\\ & =\underset{v}{\mathrm{argmin}}(||v|{|}_1+\frac{1}{2\cdot \lambda }||z-v|{|}_2^2)\\ & \therefore \underset{v_i}{\mathrm{argmin}}(\frac{1}{2}(v_i-z_i{)}^2+\lambda |v_i|)\end{array}$$

$$\begin{array}{rl}{Z}^{k+1}& :=\underset{Z}{\mathrm{argmin}}(g(Z)+\frac{\rho }{2}||I{\beta }^{k+1}-IZ+{\mu }^k|{|}_2^2)\\ & =\underset{Z}{\mathrm{argmin}}(g(Z)+\frac{\rho }{2}||{\beta }^{k+1}+{\mu }^k-Z)|{|}_2^2\\ & =\underset{Z}{\mathrm{argmin}}(g(Z)+\frac{1}{2\cdot \frac{1}{\rho }}||{\beta }^{k+1}+{\mu }^k-Z)|{|}_2^2\\ & \therefore pro{x}_{\frac{1}{\rho },g}({\beta }^{k+1}+{\mu }^k)\end{array}$$

• Given $\beta$, $z$, $\mu$, $\rho$ to some initial value

• Repeat:

  • $\therefore {\beta }^{k+1}=(2X^{t}X+\rho I{)}^{-1}(2X^{t}y+\rho (Z^k-{\mu }^k))$
  • $\therefore Z^{k+1}=pro{x}_{\frac{1}{\rho },g}({\beta }^{k+1}+{\mu }^k)$
  • ${\mu }^{k+1}:={\mu }^k+({\beta }^{k+1}-Z^{k+1})$
  • Stopping criterion : quit $\parallel r{\parallel }_2<ϵ$ and $\parallel s{\parallel }_2<ϵ$

• prime resdual : $r^{k+1}={\beta }^{k+1}-z^{k+1}$
• dual resdual : $s^{k+1} = -(z^{k+1} - z^k) $

R code

nll <- function(X, Y, beta) {
                A <- Y - X %*% beta
                loglike <- crossprod(A)
                return(loglike)
}

# Proximal operator
prox.l1 <- function(u, lambda) {
                uhat <- abs(u) - lambda
                prox <- sign(u) * pmax(rep(0, length(u)), uhat)
                return(prox)
}

l2norm <- function(x) sqrt(sum(x^2))


ADMM <- function(X,Y,rho=5,lambda=.1,iter=100, eps = 0.0001){
                
                n <- nrow(X)
                p <- ncol(X)
                
                beta <- matrix(0, nrow=iter, ncol=p) 
                beta[1,] <- rep(0, p)
                
                obj <- rep(0, iter)
                obj[1] <- nll(X, Y, beta[1,]) + lambda * sum(abs(beta[1,]))
                
                z <- matrix(0, nrow=iter, ncol=p)
                v <- rep(0, p)     
                
                invmat <- solve(2*crossprod(X) + diag(rho, p))
                
                s <- 0    
                r <- 0    
                t <- 0
                
                for (t in 2:iter){
                                
                                beta[t,] <- invmat %*% (2*crossprod(X, Y) + rho * (z[t-1,]-v))
                                z[t,] <- prox.l1(beta[t,] + v, lambda/rho)
                                v <- v + beta[t,] - z[t,]
                                obj[t] <- nll(X, Y, beta[t,]) + lambda * sum(abs(beta[t,]))
                                
                                r <- beta[t,] - z[t,]
                                s <- -rho * (z[t,] - z[t-1,])
                                
                                r.norm <- l2norm(r)
                                s.norm <- l2norm(s)
                                
                                if (r.norm < eps & s.norm < eps) {
                                                break
                                }
                }
                beta <- beta[-c(t+1:iter),]
                obj <- obj[-c(t+1:iter)]
                result <- list("beta.hat" = beta[nrow(beta),], "beta"=beta, "objective"=obj, "iter"=t)
                return(result)
}

x <- cbind(1, matrix(rnorm(1000*4), ncol = 4))
beta <- c(1.4, -2, -3, 4, 5)

eps <- rnorm(1000*1)　
y <- x%*%beta + eps 

ADMM(X = x, Y = y)

$beta.hat
[1]  1.379323 -1.990279 -2.998234  4.018089  5.044086

$beta
         [,1]      [,2]      [,3]     [,4]     [,5]
[1,] 0.000000  0.000000  0.000000 0.000000 0.000000
[2,] 1.376414 -1.985059 -2.991330 4.008202 5.031837
[3,] 1.379269 -1.990212 -2.998172 4.018014 5.044009
[4,] 1.379323 -1.990278 -2.998234 4.018089 5.044086
[5,] 1.379323 -1.990279 -2.998234 4.018089 5.044086

$objective
[1] 58940.123  1022.952  1022.613  1022.613  1022.613

$iter
[1] 5

Citation

BibTeX citation:

@online{don2021,
  author = {Don, Don and Don, Don},
  title = {ADMM Example},
  date = {2021-05-30},
  url = {https://dondonkim.netlify.app/posts/2021-06-12-admm/admm.html},
  langid = {en}
}

For attribution, please cite this work as:

Don, Don, and Don Don. 2021. “ADMM Example.” May 30, 2021. https://dondonkim.netlify.app/posts/2021-06-12-admm/admm.html.