In TFOCS I’m minimizing the L1 norm of the derivative of x. So the argument looks like ||Dx|| where D is the derivative matrix and the norm is L1. So I decided to make a new nonsmooth function prox_l1_d.m based upon the prox_l1.m script already in TFOCS.
So take a look at prox_l1.m: I change the input argument q (and all it’s subsequent occurrences) to D (not strictly necessary…but it keeps my head clear), comment out the nargin stuff, and change line 41 from,
v = norm( qq(:).*x(:), 1 );
to
v = norm( qq*x(:), 1 );
But I have no idea what to do about lines 44-48 as I don’t understand the idea of a proximal function. So I have two questions:
Is this a good way to implement L1 minimization of the derivative of x in TFOCS? How do I change lines 44-48 to make it work?
As you know, TFOCS is not as “convenient” as CVX in the fact that you are required to implement the core computations for each function that is used (unless, of course, your functions happen to be among the ones we have already coded up). For a “smooth” function, the code must compute the function value and, when requested, the gradient, at any point.
For a prox functionh(z), the requirements are different. It must be able to compute the function value, of course, but also, when requested, it must be able to solve the proximity minimization
$$\Phi_h(x,t) = \mathop{\textrm{argmin}}_z h(z) + \tfrac{1}{2}t^{-1}|z-x|_2^2$$The prox_l1 function, for instance, solves
$$\Phi_h(x,t) = \mathop{\textrm{argmin}}_z |z|_1 + \tfrac{1}{2}t^{-1}|z-x|_2^2$$
which happens to be relatively simple, because it is decoupled: each z_i can be solved separately. In your case, you need to solve
$$\Phi_h(x,t) = \mathop{\textrm{argmin}}_z |Dz|_1 + \tfrac{1}{2}t^{-1}|z-x|_2^2$$which is more difficult, because the z values are coupled together. I will be honest, I do not know if an analytic solution exists, and I do not have time at the moment to figure this out. Perhaps Stephen Becker will chime in with his thoughts on the matter. This does look very similar to a total variation (TV) norm, so perhaps there are some tricks to borrow.
If a prox function cannot be constructed (actually, it always can, the question is whether it is simple enough to be wortwhile), then the alternative is to use the tfocs_SCD approach to solving the problem.
Everything I’ve offered here is in the documentation.
In general computing the prox of ||Dx||_1 will not be easy. If D has periodic boundary conditions, then you may be able to do it somewhat efficiently using the FFT, but in general it is not easy. You can use certain stencil patterns in D (e.g., not just the usual (1,-2,1) ) that can then decompose D in to the sum of several operators, each of which is easy to deal with.
But the simplest solution is not to do any of this, and just use D as a linear operator and use the usual prox_l1.m. TFOCS was designed exactly for this situation! Now, you have to switch to the TFOCS_SCD formulation, not TFOCS.m. The operator “D” you want is already implemented in linop_TV.m. Pass this in as the linear operator part.
Terrific! I’m having trouble figuring out syntax, can you offer some guidance? The full problem is:
min_x ||Dx||_1 subject to ||Ax-b||_2 <= epsilon
This is VERY similar to solver_sBPDN.m (and even CLOSER to solver_sBPDN_W.m) but my attempts still don’t seem to work. I’m currently using:
x = tfocs_SCD( [], {linop_TV({n,1}),0;A,-b}, {prox_l1(1),prox_l2( epsilon )}, 0 );
This return NaN. Am I translating the problem correctly?
Furthermore, I’d ACTUALLY like the constraint,
||Ax-b||_2 <= epsilon
to look like this instead:
||a.*fft(x)-b||_2 <= epsilon.
Where a is a vector and it is dot multiplied by the fft of x. If there isn’t much penalty in computation time I’d be willing to write it as a diagonal matrix multiplication to get rid of the matlab specific idea of a ‘dot multiply’.
You can see I’m struggling with syntax here. Your guidance will be very much appreciated! Thank you!
I think the prox_l1(1) should be proj_linf(1) since you need to dualize it (for example, the constraint is proj_l2 but the code uses prox_l2, as it should, since this is the dual function). So that may help some issues.
For the second issue, it shouldn’t be hard either. First, create the dot-multiply linear action you want. You could create your own linop file, or simpler, call linop_matrix and then give it spdiags(a,0,n,n) (I forget the exact syntax, but something like this), and that will be almost as efficient as the dot-multiply.
Then create a FFT linop using the premade linop_fft.