Skip to content

Commit

Permalink
gauss: change width to std (#280)
Browse files Browse the repository at this point in the history 
* gauss: change width to std

* Fix typo

Co-authored-by: Luis Fabregas <48292540+luisfabib@users.noreply.github.com>
  • Loading branch information
@stestoll @luisfabib
stestoll and luisfabib committed Mar 10, 2022
1 parent 4954618 commit 4855675
Show file tree
Hide file tree
Showing 18 changed files with 257 additions and 258 deletions.
49 changes: 24 additions & 25 deletions deerlab/dd_models.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -154,14 +154,14 @@ def _nonparametric():
where `\left<r\right>` is the mean distance and `\sigma` the standard deviation.
"""
def _gauss(r,mean,width):
return _multigaussfun(r,mean,width)
def _gauss(r,mean,std):
return _multigaussfun(r,mean,std)
# Create model
dd_gauss = Model(_gauss,constants='r')
dd_gauss.description = 'Gaussian distribution model'
# Parameters
dd_gauss.mean.set(description='Mean', lb=1.0, ub=20, par0=3.5, units='nm')
dd_gauss.width.set(description='Standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss.std.set(description='Standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
# Add documentation
dd_gauss.__doc__ = _dd_docstring(dd_gauss,notes) + docstr_example('dd_gauss')

Expand All @@ -177,16 +177,16 @@ def _gauss(r,mean,width):
where `\left<r\right>_i` are the mean distances, `\sigma_i` the standard deviations, and `a_i` are the amplitudes of the Gaussians.
"""
def _gauss2(r,mean1,width1,mean2,width2):
return _multigaussfun(r,[mean1,mean2],[width1,width2])
def _gauss2(r,mean1,std1,mean2,std2):
return _multigaussfun(r,[mean1,mean2],[std1,std2])
# Create model
dd_gauss2 = Model(_gauss2,constants='r')
dd_gauss2.description = 'Sum of two Gaussian distributions model'
# Parameters
dd_gauss2.mean1.set(description='1st Gaussian mean', lb=1.0, ub=20, par0=2.5, units='nm')
dd_gauss2.mean2.set(description='2nd Gaussian mean', lb=1.0, ub=20, par0=4.5, units='nm')
dd_gauss2.width1.set(description='1st Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss2.width2.set(description='2nd Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss2.std1.set(description='1st Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss2.std2.set(description='2nd Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss2.addlinear('amp1',description='1st Gaussian amplitude', lb=0, par0=1, units='')
dd_gauss2.addlinear('amp2',description='2nd Gaussian amplitude', lb=0, par0=1, units='')
# Add documentation
Expand All @@ -204,18 +204,18 @@ def _gauss2(r,mean1,width1,mean2,width2):
where `\left<r\right>_i` are the mean distances, `\sigma_i` the standard deviations, and `a_i` are the amplitudes of the Gaussians.
"""
def _gauss3(r,mean1,width1,mean2,width2,mean3,width3):
return _multigaussfun(r,[mean1,mean2,mean3],[width1,width2,width3])
def _gauss3(r,mean1,std1,mean2,std2,mean3,std3):
return _multigaussfun(r,[mean1,mean2,mean3],[std1,std2,std3])
# Create model
dd_gauss3 = Model(_gauss3,constants='r')
dd_gauss3.description = 'Sum of three Gaussian distributions model'
# Parameters
dd_gauss3.mean1.set(description='1st Gaussian mean', lb=1.0, ub=20, par0=2.5, units='nm')
dd_gauss3.mean2.set(description='2nd Gaussian mean', lb=1.0, ub=20, par0=3.5, units='nm')
dd_gauss3.mean3.set(description='3rd Gaussian mean', lb=1.0, ub=20, par0=5.0, units='nm')
dd_gauss3.width1.set(description='1st Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss3.width2.set(description='2nd Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss3.width3.set(description='3rd Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss3.std1.set(description='1st Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss3.std2.set(description='2nd Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss3.std3.set(description='3rd Gaussian standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gauss3.addlinear('amp1',description='1st Gaussian amplitude', lb=0, par0=1, units='')
dd_gauss3.addlinear('amp2',description='2nd Gaussian amplitude', lb=0, par0=1, units='')
dd_gauss3.addlinear('amp3',description='3rd Gaussian amplitude', lb=0, par0=1, units='')
Expand All @@ -239,15 +239,15 @@ def _gauss3(r,mean1,width1,mean2,width2,mean3,width3):
where `\left<r\right>` is the mean distance,`\sigma` is the standard deviation, and `\beta` is the kurtosis of the distribution.
"""
def _gengauss(r,mean,width,kurt):
P = kurt/(2*width*spc.gamma(1/kurt))*np.exp(-abs(r-mean)/width**kurt)
def _gengauss(r,mean,std,kurt):
P = kurt/(2*std*spc.gamma(1/kurt))*np.exp(-abs(r-mean)/std**kurt)
return _normalize(r,P)
# Create model
dd_gengauss = Model(_gengauss,constants='r')
dd_gengauss.description = 'Generalized Gaussian distribution model'
# Parameters
dd_gengauss.mean.set(description='Mean', lb=1.0, ub=20, par0=3.5, units='nm')
dd_gengauss.width.set(description='Standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gengauss.std.set(description='Standard deviation', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_gengauss.kurt.set(description='Kurtosis', lb=0.25, ub=15, par0=5.0, units='')
# Add documentation
dd_gengauss.__doc__ = _dd_docstring(dd_gengauss,notes) + docstr_example('dd_gengauss')
Expand All @@ -265,20 +265,20 @@ def _gengauss(r,mean,width,kurt):
<img src="../_images/model_scheme_dd_skewgauss.png", style="width: 50%">
<br><br><br>
:math:`P(r) = \frac{1}{\sqrt{2\pi}}\exp\left(-\frac{(r-\left<r\right>)^2}{\sqrt(2)\sigma^2}\right)\frac{1}{2}\left(1 + \mathrm{erf}\left(\frac{(r-\left<r\right>)}{\sqrt{2}\sigma}\right) \right)`
:math:`P(r) = \frac{1}{\sqrt{2\pi}}\exp\left(-\frac{(r-\left<r\right>)^2}{2\sigma^2}\right)\left(1 + \mathrm{erf}\left(\frac{(r-\left<r\right>)}{\sqrt{2}\sigma}\right) \right)`
where `\left<r\right>` is the center distance,`\sigma` is the spread, and `\alpha` is the skewness of the distribution.
"""
def _skewgauss(r,center,width,skew):
x = (r-center)/width/np.sqrt(2)
P = 1/np.sqrt(2*np.pi)*np.exp(-x**2)*(1 + spc.erf(skew*x))
"""
def _skewgauss(r,center,std,skew):
x = (r-center)/std
P = 1/np.sqrt(2*np.pi)*np.exp(-x**2/2)*(1 + spc.erf(skew*x/np.sqrt(2)))
return _normalize(r,P)
# Create model
dd_skewgauss = Model(_skewgauss,constants='r')
dd_skewgauss.description = 'Skew Gaussian distribution model'
# Parameters
dd_skewgauss.center.set(description='Center', lb=1.0, ub=20, par0=3.5, units='nm')
dd_skewgauss.width.set(description='Spread', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_skewgauss.std.set(description='Spread', lb=0.05, ub=2.5, par0=0.2, units='nm')
dd_skewgauss.skew.set(description='Skewness', lb=-25, ub=25, par0=5, units='')
# Add documentation
dd_skewgauss.__doc__ = _dd_docstring(dd_skewgauss,notes) + docstr_example('dd_skewgauss')
Expand Down Expand Up @@ -1018,12 +1018,11 @@ def _wormchain(r,contour,persistence):
Radial Distribution Function of Semiflexible Polymers
Phys. Rev. Lett. 77(12), 2581-2584, 1996
"""
def _wormgauss(r,contour,persistence,width):
sigma = width
def _wormgauss(r,contour,persistence,std):
P = _wlc(r,contour,persistence)
# Compute Gaussian convolution window
idx = np.argmax(P)
gauss = np.exp(-((r - r[idx])/sigma)**2)
gauss = np.exp(-((r - r[idx])/std)**2)
# Convolution with size retention
P = np.convolve(gauss,P,mode='full')
# Adjust new convoluted axis
Expand All @@ -1040,6 +1039,6 @@ def _wormgauss(r,contour,persistence,width):
# Parameters
dd_wormgauss.contour.set(description='Contour length', lb=1.5, ub=10, par0=3.7, units='nm')
dd_wormgauss.persistence.set(description='Persistence length', lb=2, ub=100, par0=10, units='nm')
dd_wormgauss.width.set(description='Gaussian standard deviation', lb=0.01, ub=5, par0=0.2, units='nm')
dd_wormgauss.std.set(description='Gaussian standard deviation', lb=0.01, ub=5, par0=0.2, units='nm')
# Add documentation
dd_wormgauss.__doc__ = _dd_docstring(dd_wormgauss,notes) + docstr_example('dd_wormgauss')
4 changes: 2 additions & 2 deletions docsrc/source/basics.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,14 +31,14 @@ Parametric distance distributions

>>>print(dl.dd_gauss)
Description: Gaussian distribution model
Signature: (r, mean, width)
Signature: (r, mean, std)
Constants: [r]
Parameter Table:
======= ======= ======= ======== ======== ======= ====================
Name Lower Upper Type Frozen Units Description
======= ======= ======= ======== ======== ======= ====================
mean 1 20 nonlin No nm Mean
width 0.05 2.5 nonlin No nm Standard deviation
std 0.05 2.5 nonlin No nm Standard deviation
======= ======= ======= ======== ======== ======= ====================

Expand Down
6 changes: 3 additions & 3 deletions docsrc/source/fitting_guide.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -177,8 +177,8 @@ In this example, let us construct a penalty that imposes the smoothness of a Gau
# Differential operator of second order
L = dl.regoperator(x,2)
# Define the penalty function
def smoothness_fcn(center,width):
gaussian = gauss(center,width)
def smoothness_fcn(center,std):
gaussian = gauss(center,std)
penalty_vector = L@gaussian
return penalty_vector
# Construct the penalty
Expand Down Expand Up @@ -219,7 +219,7 @@ A summary of the fit can be accessed by printing the ``FitResult`` object as ret
Parameter Value 95%-Confidence interval Units Description
=========== ======= ========================= ======= ================
center 5.002 (4.995,5.010) None None
width 0.204 (0.196,0.212) None None
std 0.204 (0.196,0.212) None None
scale 0.999 (0.999,0.999) None Scaling factor
=========== ======= ========================= ======= ================
Expand Down

0 comments on commit 4855675

Please sign in to comment.