FindSparse: Implements the Fast Laplace Algorithm

In bcs: Bayesian Compressive Sensing Using Laplace Priors

Description

Implements the fast Laplace algorithm given a measurement matrix and samples from a signal or function. Note that this function is a convenient wrapper for the function FastLaplace.

Usage

FindSparse(PHI, y, eta = 1e-08, roundit = FALSE, verbose = FALSE)

Arguments

PHI	typically equals the product of a measurment matrix and basis representation matrix, such as the wavelet basis. The solution vector w (see Details below) is assumed to be sparse in the chosen basis.
y	CS measurements, samples from the signal or function.
eta	tolerance level in determining convergence of marginal likelihood.
roundit	whether or not to round the marginal likelihood, in order to avoid machine precision error when comparing across platforms.
verbose	print to screen which basis are added, re-estimated, or deleted.

Details

This code implements the fast Laplace algorithm. The fast Laplace algorithm is a method used to solve the compressive sensing problem, or in general, a highly underdetermined system of equations. It does this by taking the system of equations

y = Φ w + n

and converting it into a minimization problem where we minimize the error with a constraint on w (the vector we are solving for) that enforces sparsity. The fast Laplace method uses a Bayesian framework, and in particular, uses a Laplace prior to enforce sparsity on w. See [1] for more information.

Value

The sparse signal w as found by the fast Laplace algorithm.

References

[1] S. D. Babacan, R. Molina and A. K. Katsaggelos, "Bayesian Compressive Sensing Using Laplace Priors," in IEEE Transactions on Image Processing, vol. 19, no. 1, pp. 53-63, Jan. 2010.

Examples

# size of the basis function expansion
N <- 64

# generate sparse coefficient vector
w <- rep(0,N)
w[sample(1:N,10)] <- runif(10,-1,1)

# create wavelet basis trasform matrix
wavelet.basis <- WaveletBasis(N)
# generate actual signal
signal <- wavelet.basis%*%w

# now we try and recover 'w' and 'signal' from samples
num_samps <- 25

# create random measurement matrix
measure.mat <- matrix(runif(num_samps*N),num_samps,N)
measure.mat <- measure.mat/matrix(rep(sqrt(apply(measure.mat^2,2,sum)),
                                      num_samps),num_samps,N,byrow=TRUE);

PHI <- measure.mat%*%wavelet.basis
# actual samples we see
y <- measure.mat%*%signal

# use fast Laplace algorithm
w_est <- FindSparse(PHI, y)

# compare plots of the sparse vector and the estimated sparse vector
plot(w,type='h')
lines(w_est,type='h',col='red')

# estimate signal
signal_est <- wavelet.basis%*%w_est

# Root mean squared error of estimate
error <- sqrt(mean((signal - signal_est)^2))

Example output

bcs documentation built on May 29, 2017, 11:58 a.m.