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