# Variational Multiscale Nonparametric Regression: Algorithms and Implementation
## 1 Introduction
### Computational Challenges and Scope of the Paper  
However, the practical and theoretical superiority of regularised multiscale methods over purely dictionary thresholding methods comes at a cost: instead of simple thresholding, a complex optimisation problem has to be solved, with an increase in computation time and the need of improved optimisation methods. In practice, the problem (3) is non-smooth (e.g., if \(R\) is chosen as the TV seminorm) and is furthermore high-dimensional due to a huge number of constraints. For instance, in the case of images of \(256\times 256\) pixels, the number of constraints is 1,751,915 for the dictionary of cubes of edge length \(\leq 30\) and 3,211,264 for that of shearlets.  
Besides our theoretical review in Section 2, the goal of this paper is to discuss different particularly successful algorithmic approaches for the solution of (3) in Section 3. Among the first attempts is the usage of the alternating direction method of multipliers (ADMM) in [31; 32; 36], including a problem-specific convergence analysis. The computational disadvantage of this approach is the necessity to compute projections to the constraint set, which is most generally done by [PERSON]'s algorithm. Even though this can be efficiently implemented using GPUs (cf. [37; 38]), the projection step provides a computational bottleneck of the overall algorithm. Several other convex optimisation algorithms such as Douglas-Rachford-Splitting (see [39]) or the Proximal Alternating Predictor Corrector algorithm (cf. [40]) also require the computation of such projections, and hence their computational performance is similar to the ADMM-based version. In view of this, we furthermore discuss two new approaches based on the Chambolle-Pock algorithm [41] on the one hand and a semismooth [PERSON] method (cf. [42; 43]) on the other hand. The former allows avoiding computing projections to the constraint sets, but requires only the computation of the corresponding resolvent operator, which reduces to a high-dimensional soft-thresholding problem. The latter allows solving a regularised version of the optimality conditions by [PERSON]'s method and then apply a path-continuation scheme to decrease the amount of regularisation. Both algorithms hence avoid [PERSON]'s projection step, but still have favourable convergence properties.  
In Section 4, we finally discuss the practical advantages and disadvantages of the previously described algorithms along different numerical examples. For all examples checked in this paper, we found that the Chambolle-Pock algorithm is superior compared to the others. Finally, we provide a Matlab implementation of the Chambolle-Pock algorithm for this problem and code to run all examples in [[https://github.com/housenli/MIND](https://github.com/housenli/MIND)]([https://github.com/housenli/MIND](https://github.com/housenli/MIND)).  
Figure 1: Comparison of MIND with TV regularisation and the soft-thresholding (ST) estimator on the “building” image from Darmstadt Noise Dataset [34] with different dictionaries. Details of dictionaries can be found in Section 4.2. The thresholds for MIND and ST are chosen the same as the 50% quantile of \(\|\tilde{\sigma}\epsilon\|_{\text{MS}}\) in (5) (see Section 2.2). The noise level \(\sigma\) is estimated by a second order difference-based estimator [35] (see Section 4.3). (**a**) Truth; (**b**) Noisy image; (**c**) MIND, wavelets; (**d**) MIND, curvelets; (**e**) MIND, shearlets; (**f**) ST, wavelets; (**g**) ST, curvelets; (**h**) ST, shearlets.