Learning bases using the Lagrange dual (R)

LagrangeDualDict.r

# = = = = = include = = = = = #

library(MASS)

# = = = = = function = = = = = #

Obj_func <- function(Y,B,A,Lambda){
  MAT <- t(Y) %*% Y - Y %*% t(A) %*% solve(A%*%t(A)+Lambda) %*% t(Y%*%t(A)) - Lambda
  return <- sum(diag(MAT))
}

LagrangeDualDict <- function(Y,B,A,lambda=NULL){

	## norm constriant
	c <- 1

	## set lambda
	if(is.null(lambda)){
		lambda <- rep(1,ncol(B))
	}

	for(ite in 1:10){

		Lambda <- diag(lambda)

		print("Obj_value")
		print(Obj_func(Y, B, A, Lambda))
		
		StSpLinv <- solve(A %*% t(A) + Lambda)
		Box <- Y %*% t(A) %*% StSpLinv

		## Gradient
		Grad_func <- colSums( (Box)^2 ) - 1
		
		## Hessian
		H <- -2 * crossprod(Box) * StSpLinv

		lambda_new <- as.vector(lambda - solve(H) %*% Grad_func)
		if( sum( (lambda-lambda_new)^2 ) < 1e-4 ){
			break
		}
		lambda <- lambda_new
	}

	return( list(B=Box,lambda=lambda_new) )
}

## <<References>>
## [1] Efficient sparse coding algorithms
## 	Honglak Lee, Alexis Battle, Rajat Raina, and Andrew Y. Ng.
## 	http://ai.stanford.edu/~hllee/softwares/nips06-sparsecoding.htm