Commits

Casey Dunn committed f7e2b2b

moved convex minimization to a function

  • Participants
  • Parent commits 6390c84

Comments (0)

Files changed (2)

-
 #' regularization_by_thresholding
 #' 
 #' Regularizes a correlation matrix by thresholding
 	return( regularized )
 }
 
-regularization_by_convex_minimization <- function (cor_matrix) {
-
-	
-	return( luo )
-}
-
 
+#' regularization_by_convex_minimization
+#' 
+#' Regularizes a correlation matrix by convex minimization
+#' 
+#' @param cor_matrix A correlation matrix
+#' 
+#' @param n The number of independent observations that were used to construct 
+#' the 
+#' 
+#' @return The regularized correlation matrix
+#' 
+#' @references
+#' Luo, X. High Dimensional Low Rank and Sparse Covariance Matrix Estimation via 
+#' Convex Minimization. arXiv.org (2011). http://arxiv.org/abs/1111.1133
+#' 
+regularization_by_convex_minimization <- function ( cor_matrix, n ) {
 
-test_cov_matrix <- function() {
-	# Builds a simple character covariance matrix for simulating data
-	
-	G <- 100
+	p <- ncol( cor_matrix )
 	
-	trueCovariance = matrix(0,G,G)
-	trueCovariance[1:10,1:10] = 0.95
-	trueCovariance[11:80,11:80] = 0.3
-	trueCovariance[81:100,81:100] = 0.7
-	diag(trueCovariance) = 1
-
-	return( trueCovariance ) 
+	re.lorec <- lorec( cor_matrix, diag(1, 100), diag(1, 100), sqrt(p/n), 
+		sqrt(log(p)/n) )
+	re.lorec.eig <- eigen( re.lorec$L )
+	# threholding both
+	sel <- re.lorec.eig$values > sqrt( p/n )
+	if ( sum(sel) > 0 ) {
+		V <- re.lorec.eig$vectors
+		V <- V * ( abs(V)>sqrt(1/p) )
+		V <- V[ ,1:sum(sel) ]
+		Lhat <- V%*%diag( re.lorec.eig$values[1:sum(sel)] )%*%t(V)
+	} else {
+		Lhat <- matrix(0, p, p)
+	}
+	Shat <- re.lorec$S
+	Shat <- Shat*(abs(Shat)>sqrt(1/p))
+	regularized <- Lhat + Shat
+	return( regularized )
 }
-
-plot_matrix <- function(m, ... ) {
-	
-	nr <- nrow(m)
-	nc <- ncol(m)
-	image(1:nc, 1:nr, t(m[nr:1, ]), axes=F,xlab="", ylab="", ... )
-}
 
 set.seed(123456)
 
+
+# A couple functions to facilitate simulation
+
+test_cov_matrix <- function() {
+	# Builds a simple character covariance matrix for simulating data
+	
+	G <- 100
+	trueCovariance = matrix(0,G,G)
+	trueCovariance[1:10,1:10] = 0.95
+	trueCovariance[11:80,11:80] = 0.3
+	trueCovariance[81:100,81:100] = 0.7
+	diag(trueCovariance) = 1
+
+	return( trueCovariance ) 
+}
+
+plot_matrix <- function(m, ... ) {
+	# Basic plotting
+	nr <- nrow(m)
+	nc <- ncol(m)
+	image(1:nc, 1:nr, t(m[nr:1, ]), axes=F,xlab="", ylab="", ... )
+}
+
+
 # Build the true matrix
 trueCovariance = test_cov_matrix()
 
 bickel <- regularization_by_thresholding( contrastcor, n )
 
 # Regularization by convex Minimization
-# Luo, X. High Dimensional Low Rank and Sparse Covariance Matrix Estimation via 
-# Convex Minimization. arXiv.org (2011). http://arxiv.org/abs/1111.1133
-
-p <- ncol( contrastcor )
-
-re.lorec <- lorec( contrastcor, diag(1, 100), diag(1, 100), sqrt(p/n), 
-	sqrt(log(p)/n) )
-re.lorec.eig <- eigen( re.lorec$L )
-# threholding both
-sel <- re.lorec.eig$values > sqrt( p/n )
-if ( sum(sel) > 0 ) {
-	V <- re.lorec.eig$vectors
-	V <- V * ( abs(V)>sqrt(1/p) )
-	V <- V[ ,1:sum(sel) ]
-	Lhat <- V%*%diag( re.lorec.eig$values[1:sum(sel)] )%*%t(V)
-} else {
-	Lhat <- matrix(0, p, p)
-}
-Shat <- re.lorec$S
-Shat <- Shat*(abs(Shat)>sqrt(1/p))
-luo <- Lhat + Shat
+
+luo <- regularization_by_convex_minimization( contrastcor, n )
 
 
 # Set up plotting parameters