Source code for regularization_methods

# -*- coding: utf-8 -*-
"""
Created on Wed Dec 04 09:00:13 2013

Regularization of RichLucy is carried out by damping the update in each iteration.

The function RichLucy.RichLucy therefore uses

    Shat_new = regFun(Shat,delta)*(RL_update)

where RL_update is conv(M/(conv(Shat,b)),b(-lambda)).



Example: RichLucy with Tikhonov regularization

    Shat = RichLucy.RichLuy(b,delta,M,regFun=lambda Sh,delta: Tikhonov(Sh,reglambda,delta))

with reglambda the chosen regularization parameter (e.g. reglambda=1e-5)


Implemented regularization methods

* maxEntropy
* Tikhonov
* TotalVar

@author: Sascha Eichstaedt
"""

import numpy as np

[docs]def maxEntropy(Shat,deltaS,reglambda): """ maximum entropy regularization Criterion introduced for RichLucy in J Boutet de Monvel, S Le Calvez, and M Ulfendahl Image Restoration for Confocal Microscopy: Improving the Limits of Deconvolution, with Application to the Visualization of the Mammalian Hearing Organ Biophysical Journal Volume 80 May 2001 2455–2470 """ inds = np.nonzero(Shat>1e-18)[0] Shat[inds] -= reglambda*Shat[inds]*np.log(Shat[inds]) return Shat
[docs]def Tikhonov(Shat,reglambda,deltaS): """ Regularization using Tikhonov regularization Therefore, each estimated spectrum in SH is weighted with the regularizer: 1+lambda*Laplace(SH(k-1,:)) i.e., the smoothness of the previous iteration. Criterion introduced for RichLucy in G van Kempen, Image Restoration in Flourescence Microscopy, PhD Thesis, Technische Universiteit Delft, Holland, 1999 """ DSH = -np.gradient(np.gradient(Shat,deltaS),deltaS) # Shat = Shat + (1-1.0/(1+2*reglambda*DSH)) Shat *= 1.0/(1+2*reglambda*DSH) Shat[Shat<0]=0.0 return Shat
[docs]def TotalVar(Shat,reglambda,deltaS): """ Determination of optimal number of iterations by using total variation regularization Therefore, each estimated spectrum in SH is weighted with the regularizer: 1-lambda*divergence(gradient(SH(k-1,:)/abs(gradient(SH(k-1,:))) Criterion introduced for RichLucy in N Dey, L Blanc-Feraud, C Zimmer, P Roux, Z Kam, J C Olivo-Marin and J Zerubia 3D Microscopy Deconvolution using Richardson-Lucy Algorithm with Total Variation Regularization Rapport de Recherce, No 8272, July 2004, INRIA """ DHS = np.gradient(Shat,deltaS) inds = np.nonzero(np.abs(DHS)>0)[0] Shat[inds] *= 1.0/(1.0 - reglambda*np.gradient( DHS[inds]/np.abs(DHS[inds]), deltaS) ) Shat[Shat<0]=0.0 return Shat