-:math:`P(r) = \frac{\beta}{2\sigma\Gamma(1/\beta)}\exp\left(-\left(\frac{(r-\left)}{\sigma}\right)^\beta \right)`
+:math:`P(r) = \frac{\beta}{2\sigma\Gamma(1/\beta)}\exp\left(-\left(\frac{|r-\left|}{\sigma}\right)^\beta \right)`
-where `\left` is the mean distance,`\sigma` is the standard deviation, and `\beta` is the kurtosis of the distribution.
-"""
-def _gengauss(r,mean,std,kurt):
- P = kurt/(2*std*spc.gamma(1/kurt))*np.exp(-abs(r-mean)/std**kurt)
+where `\left` is the mean distance,`\sigma` is the standard deviation, and `\beta` determines the shape of the distribution.
+"""
+def _gengauss(r,mean,std,beta):
+ P = beta/(2*std*spc.gamma(1/beta))*np.exp(-(abs(r-mean)/std)**beta)
return _normalize(r,P)
# Create model
dd_gengauss = Model(_gengauss,constants='r')
@@ -248,7 +248,7 @@ def _gengauss(r,mean,std,kurt):
# Parameters
dd_gengauss.mean.set(description='Mean', lb=1.0, ub=20, par0=3.5, unit='nm')
dd_gengauss.std.set(description='Standard deviation', lb=0.05, ub=2.5, par0=0.2, unit='nm')
-dd_gengauss.kurt.set(description='Kurtosis', lb=0.25, ub=15, par0=5.0, unit='')
+dd_gengauss.beta.set(description='beta parameter', lb=0.25, ub=15, par0=5.0, unit='')
# Add documentation
dd_gengauss.__doc__ = _dd_docstring(dd_gengauss,notes) + docstr_example('dd_gengauss')
diff --git a/deerlab/deerload.py b/deerlab/deerload.py
index 88e49409..e6c294b8 100644
--- a/deerlab/deerload.py
+++ b/deerlab/deerload.py
@@ -28,17 +28,17 @@ def deerload(fullbasename, plot=False, full_output=False, *args,**kwargs):
Returns
-------
- t : ndarray
+ t : ndarray, or list of ndarray
Time axis in microseconds. Its structure depends on the dimensionality of the experimental datasets:
- * 1D-datasets: the time axis of ``N`` points is returned as a one-dimensional ndarray of shape ``(N,)``
- * 2D-datasets: the ``M`` time-axes of ``N`` points are returned as a two-dimensional ndarray of shape ``(N,M)``. The i-th axis can be accessed via ``t[:,i]``.
+ * 1D datasets: the time axis is returned as a one-dimensional ndarray of shape ``(N,)``
+ * 2D datasets: the axes are returned as elements of a list, with the first axis typically being the time axis of shape ``(N,)``
V : ndarray
Experimental signal(s). Its structure depends on the dimensionality of the experimental datasets:
- * 1D-datasets: the signal of ``N`` points is returned as a one-dimensional ndarray of shape ``(N,)``
- * 2D-datasets: the ``M`` signals of ``N`` points are returned as a two-dimensional ndarray of shape ``(N,M)``. The i-th signal can be accessed via ``V[:,i]``.
+ * 1D datasets: the signal of ``N`` points is returned as a one-dimensional ndarray of shape ``(N,)``
+ * 2D datasets: the ``M`` signals of ``N`` points are returned as a two-dimensional ndarray of shape ``(N,M)``. The i-th signal can be accessed via ``V[:,i]``.
pars : dict
Parameter file entries, returned if ``full_output`` is ``True``.
@@ -217,7 +217,7 @@ def deerload(fullbasename, plot=False, full_output=False, *args,**kwargs):
abscissa = np.squeeze(abscissa)
abscissas = []
# Convert to list of abscissas
- for absc in abscissa.T:
+ for absc in abscissa.T:
# Do not include abcissas full of NaNs
if not all(np.isnan(absc)):
# ns -> µs converesion
diff --git a/deerlab/dipolarmodel.py b/deerlab/dipolarmodel.py
index e5cc3e52..fec3691c 100644
--- a/deerlab/dipolarmodel.py
+++ b/deerlab/dipolarmodel.py
@@ -728,7 +728,7 @@ def ex_3pdeer(tau, pathways=[1,2], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
@@ -792,7 +792,7 @@ def ex_4pdeer(tau1, tau2, pathways=[1,2,3,4], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
@@ -860,7 +860,7 @@ def ex_rev5pdeer(tau1, tau2, tau3, pathways=[1,2,3,4,5], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
@@ -935,7 +935,7 @@ def ex_fwd5pdeer(tau1, tau2, tau3, pathways=[1,2,3,4,5], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
@@ -979,7 +979,7 @@ def ex_sifter(tau1, tau2, pathways=[1,2,3], pulselength=0.016):
r"""
Generate a 4-pulse SIFTER dipolar experiment model.
- The figure below shows the dipolar pathways in 4-pulse SIFTER. The observer (black) and pump (grey) pulses and their interpulse delays are shown on the top.
+ The figure below shows the dipolar pathways in 4-pulse SIFTER. The pulses and the interpulse delays are shown on the top.
The middle table summarizes all detectable modulated dipolar pathways `p` along their dipolar phase accumulation factors `\mathbf{s}_p`,
harmonics `\delta_p` and refocusing times `t_{\mathrm{ref},p}`. The most commonly encountered pathways are highlighted in color.
The bottom panel shows a decomposition of the dipolar signal into the individual intramolecular contributions (shown as colored lines).
@@ -1005,7 +1005,7 @@ def ex_sifter(tau1, tau2, pathways=[1,2,3], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
@@ -1044,7 +1044,7 @@ def ex_ridme(tau1, tau2, pathways=[1,2,3,4], pulselength=0.016):
r"""
Generate a 5-pulse RIDME dipolar experiment model.
- The figure below shows the dipolar pathways in 5-pulse RIDME. The observer (black) and pump (grey) pulses and their interpulse delays are shown on the top.
+ The figure below shows the dipolar pathways in 5-pulse RIDME. The pulses and their interpulse delays are shown on the top.
The middle table summarizes all detectable modulated dipolar pathways `p` along their dipolar phase accumulation factors `\mathbf{s}_p`,
harmonics `\delta_p` and refocusing times `t_{\mathrm{ref},p}`. The most commonly encountered pathways are highlighted in color.
The bottom panel shows a decomposition of the dipolar signal into the individual intramolecular contributions (shown as colored lines).
@@ -1069,7 +1069,7 @@ def ex_ridme(tau1, tau2, pathways=[1,2,3,4], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
@@ -1110,7 +1110,7 @@ def ex_dqc(tau1, tau2, tau3, pathways=[1,2,3], pulselength=0.016):
r"""
Generate a 6-pulse DQC dipolar experiment model.
- The figure below shows the dipolar pathways in 6-pulse DQC. The observer (black) and pump (grey) pulses and their interpulse delays are shown on the top.
+ The figure below shows the dipolar pathways in 6-pulse DQC. The pulses and their interpulse delays are shown on the top.
The middle table summarizes all detectable modulated dipolar pathways `p` along their dipolar phase accumulation factors `\mathbf{s}_p`,
harmonics `\delta_p` and refocusing times `t_{\mathrm{ref},p}`. The most commonly encountered pathways are highlighted in color.
The bottom panel shows a decomposition of the dipolar signal into the individual intramolecular contributions (shown as colored lines).
@@ -1139,7 +1139,7 @@ def ex_dqc(tau1, tau2, tau3, pathways=[1,2,3], pulselength=0.016):
pulselength : float scalar, optional
Length of the longest microwave pulse in the sequence in microseconds. Used to determine the uncertainty in the
- boundaries of the pathay refocusing times.
+ boundaries of the pathway refocusing times.
Returns
-------
diff --git a/docsrc/source/dipolar_guide_modeling.rst b/docsrc/source/dipolar_guide_modeling.rst
index cb791a1f..f8a421e0 100644
--- a/docsrc/source/dipolar_guide_modeling.rst
+++ b/docsrc/source/dipolar_guide_modeling.rst
@@ -3,86 +3,103 @@
Modeling
=========================
-DeerLab provides a very flexible framework to model dipolar signals originating from any dipolar EPR spectroscopy experiments. Choosing a model that properly describes your sample and experiment is of paramount importance. The DeerLab function ``dipolarmodel`` already defines the core model structure based on dipolar pathways, with the following components to be chosen:
+DeerLab provides a very flexible framework for simulating/modeling signals from a wide range of dipolar EPR spectroscopy experiments. For this, the central task is to use the DeerLab function ``dipolarmodel`` to define the model for the signal: ::
-* **Distance range**: Also called the interspin distance axis, is the range of distances where the distribution is defined.
+ Vmodel = dl.dipolarmodel(t,r,Pmodel,Bmodel,experiment)
-* **Distribution model**: Describes the intra-molecular distance distribution in either a parametric (e.g. a Gaussian distribution) or a non-parametric way.
+The following components need to be provided:
-* **Background model**: Describes the dipolar background signal arising from the inter-molecular contributions.
+* **Time range** ``t``: This is the range of dipolar evolution times for the desired signal.
-* **Number of pathways**: Sets the number of dipolar pathways contributing to the dipolar signal.
+* **Distance range** ``r``: This is the range of intra-molecular spin-spin distances where the distribution is defined.
-For each of these four components, a choice needs to be made:
+* **Distribution model** ``Pmodel``: Describes the intra-molecular distance distribution in either a parametric form (e.g. a Gaussian distribution) or a non-parametric form (as a histogram).
-Choosing a distance range
-*************************
+* **Background model** ``Bmodel``: Describes the dipolar background signal arising from the inter-molecular spin pairs.
+
+* **Type of experiment** ``experiment``: This sets the number of dipolar pathways contributing to the dipolar signal.
-The distance range :math:`[r_\mathrm{min},r_\mathrm{max}]` is an important choice, as any distance distribution is truncated to this range, i.e. :math:`P(r)=0` for :math:`rr_\mathrm{max}`. The lower limit of the distance range is determined by the bandwidth of the pulses, and also by the time increment. Typically, 1.5 nm is a reasonable choice. The upper limit depends on the distances in your sample. The number of points in ``r`` is usually set to a certain resolution (typically 0.01-0.05nm). Such a distance-axis is usually defined as ``r`` is most easily defined using the ``linspace`` function from NumPy: ::
+For each of these components, a choice needs to be made:
- r = np.linspace(1.5,6.5,100) # define distance range from 1.5nm to 6.5nm with a resolution of 0.05nm
+The time range
+*************************
-Choosing a distribution model
-******************************
+The time range depends the expected range of distances, but is typically beteen about 0.5 and several microseconds. If you are planning to use this model to fit experimental data, the time vector should be the same as the one used in the experiment. Otherwise, to construct a time vector, use the ``linspace`` or ``arange`` functions from NumPy: ::
-A non-parametric distribution is specified by setting the choice of ``Pmodel`` keyword in ``dipolarmodel`` to ``None``. In a non-parametric distribution, each element :math:`P_i` of the distribution is a linear parameter. Non-parametric distributions are obtained via methods such as Tikhonov regularization. If there are reasons to believe that the distance distribution has a specific shape (e.g. Gaussian, Rice, random-coil, etc.), or if there is very little information in the data, use a parametric distance distribution model from the :ref:`list of available models`.
+ # time from 0.2 µs to 3.0 µs with a resolution of 0.01 µs, 321 points
+ t = np.linspace(0.2,3.0,321) # start, stop, number of points
+ t = np.arange(0.2,3.0,0.01) # start, stop, step size
+
+Note that DeerLab interprets ``t`` such that its zero point is right at the end of the preceding pulse (e.g. right after the second observer pulse in 4-pulse DEER) rather than at the refocusing point of the signal (after `\tau_1` after the second observer pulse).
-Choose a background model
+The distance range
*************************
-Typically, a background model of a homogenous 3D distribution of spins is appropriate. The associated parametric model function is :ref:`bg_hom3d`. In some cases, depending on the properties of your sample, other background models might be needed, such as backgrounds arising from distributions of spins in fractal dimensions or when accounting for volume-exclusion effects. In such cases, use the associated parametric background models from the :ref:`list of available models`. If there is no inter-molecular background in your sample, or it is negligible, set the background model to ``None``.
+The distance range is a vector that goes from a minimum distance :math:`r_\mathrm{min}` to a maximum distance :math:`r_\mathrm{max}`. Any distance distribution is truncated to this range, i.e. it is assumed to be zero outside that range. The lower limit of the distance range is determined by the bandwidth of the pulses, and also by the time increment. Typically, 1.5 nm is a reasonable choice. The upper limit depends on the distances in your sample. The number of points in ``r`` is usually set to a certain resolution (typically 0.01-0.05 nm). To construct a distance vector, use the ``linspace`` or ``arange`` functions from NumPy: ::
+
+ # define distance range from 1.5 nm to 6.5 nm with a resolution of 0.05 nm
+ r = np.linspace(1.5,6.5,101)
+ r = np.arange(1.5,6.5,0.05)
+
+The distance distribution model
+**********************************
+A key aspect of the model is the type of distance distribution. A non-parametric distribution is specified by setting the choice of ``Pmodel`` keyword in ``dipolarmodel`` to ``None``. In a non-parametric distribution, `P` is a vector of the same length as ``r``, and each element :math:`P_i` is a parameter. Basically, the distribution is represented as a histogram. When analyzing data, such distributions are used with methods such as Tikhonov regularization. If there are reasons to believe that the distance distribution has a specific shape (e.g. Gaussian, Rice, random-coil, etc.), or if there is very little information in the data, use a parametric distance distribution model from the :ref:`list of available models`.
+
+The background model
+*************************
+
+Typically, a background model of a homogenous 3D distribution of spins is appropriate. The associated parametric model function is :ref:`bg_hom3d`. In some cases, depending on the properties of your sample, other background models might be needed, such as backgrounds arising from distributions of spins in fractal dimensions or when accounting for volume exclusion effects. In such cases, use the associated parametric background models from the :ref:`list of available models`. If there is no inter-molecular background in your sample, or if it is negligible, set the background model to ``None``.
-Choosing the number of dipolar pathways
-***************************************
+The number of dipolar pathways
+***************************************
-This decision should be based on the experiment you used to acquire the data and the type of pulses you used. In the case of 4-pulse DEER data, when analyzing a standard 4-pulse DEER signal without 2+1 component at the end a single pathway suffices. If the 2+1 component (appearing at the right edge of the time trace) is present, then it should be fitted as well, including its counterpart appearing at negative times, making a total of three dipolar pathways. Experiments such as 5-pulse DEER typically require at least two dipolar pathways to be properly modelled.
+What distibguishes different dipolar EPR experiments within DeerLab is the number of dipolar pathways and their refocusing times. How many to include depends on the experiment you used to acquire the data and the type of pulses you used. For example, for analyzing a standard 4-pulse DEER signal without 2+1 component at the end, a single pathway suffices. If the 2+1 component (appearing at the right edge of the time trace) is present, then it should be fitted as well, including its counterpart appearing at negative times, making a total of three dipolar pathways. Experiments such as 5-pulse DEER typically require at least two dipolar pathways to be correctly modelled.
-Using experimental pulse delays
+Experimental pulse delays
********************************
The dipolar pathways of a newly constructed dipolar model are initialized at arbitrary refocusing times and fully unconstrained. The refocusing times can be strongly constrained by knowing the experimental pulse sequence delays used to acquire the data. If the experiment used to acquire the data is known, as well as its pulse delays, then it is strongly recommended do so.
DeerLab provides a selection of experimental information generators for some of the most widely employed experimental methods (see the of :ref:`list of available experiments `). These are functions that take the pulse sequence delays, and return an ``ExperimentInfo`` object. This can be passed to the ``dipolarmodel`` function via the ``experiment`` keyword argument, to incorporate the experiment information on the model and constrain some of its parameters. These experimental information generators can also take information on the duration of the longest microwave pulses to more accurately constraint the parameters when using long pulses such as in frequency-swept or shaped microwave pulses.
-When using experimental time delays and the ``experiment`` argument, the model assumes that the experimental time axis ``t`` has its zero time at the beginning of the interpulse delay (see the illustrations of the individual experiment models for details). However, experimentally it is common not to record the first few hundred nanoseconds of signal. This results in a so-called start time (often also referred to as a deadtime), which many commercial spectrometers (such as Bruker) do not account for when storing the time-vector. It is very important to account for that time shift in the model via ::
+When using an experimental time axis ``t`` and the ``experiment`` argument, the model assumes that ``t`` is zero at the beginning of the interpulse delay (see the illustrations of the individual experiment models for details). For example, for 4-pulse DEER, this zero point would be immediately after the second pulse. However, experimentally, it is common that the time axis is defined and stored differently. For 4-pulse DEER, a commercial spectrometer usually defines ``t`` such that it is zero ``tau1`` after the second pulse. It is important to correct for this offset (also called start time or deadtime): ::
- t0 = 0.4 # Experimental start time of 400ns, in μs
- t = t - t0 # Shift the time axis to account for the start time
+ t0 = 0.4 # Experimental time offset of 400 ns, in μs
+ # To convert from experimental time axis to DeerLab time axis
+ t = t + t0 # Shift the experimental time axis to match DeerLab's definition
+
+ # To convert from DeerLab time axis to experimental time axis
+ t = t - t0
+
+Without this shift, an incorrectly defined time vector will result in wrong simulated or fitted signals, and in wrong plots.
-If the time vector ``t`` does not have any time shift, this step can be skipped. Otherwise, an incorrectly defined time vector will results in wrong results.
-
Constructing the dipolar model
*******************************
Once all the decisions above have been made, the dipolar model can be constructed using the ``dipolarmodel`` function. The models that have an associated parametric function, e.g. ``bg_hom3d``, must be passed directly as inputs to ``dipolarmodel``. In Python, functions can be passed as inputs to other functions. See the :ref:`details ` on ``dipolarmodel`` for more information.
-Example: Single-pathway 4-pulse DEER model
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-For example, a 4pDEER signal with non-parametric distance distribution and homogenous 3D background can be constructed using ::
+For example, a model for 4-pulse DEER with a non-parametric distance distribution and a homogenous 3D background can be constructed using ::
expinfo = dl.ex_4pdeer(tau1=0.5, tau2=5.5, pathways=[1])
Vmodel = dl.dipolarmodel(t, r, Pmodel=None, Bmodel=dl.bg_hom3d, experiment=expinfo)
-By default, the function ``dipolarmodel`` assumes a non-parametric distance distribution, a homogenous 3D background and a single pathway. Thus the above is equivalent to ::
+By default, ``dipolarmodel`` assumes a non-parametric distance distribution, a homogenous 3D background and a single pathway. Thus the above is equivalent to ::
expinfo = dl.ex_4pdeer(tau1=0.5, tau2=5.5, pathways=[1])
Vmodel = dl.dipolarmodel(t, r, experiment=expinfo)
-
-Example: Two-pathway 5-pulse DEER model
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-For example, a 5pDEER signal with non-parametric distance distribution and homogenous 3D background can be constructed using ::
+To construct a model for 5-pulse DEER with non-parametric distance distribution and homogenous 3D background, use ::
expinfo = dl.ex_rev5pdeer(tau1=0.5, tau2=5.5, tau3=0.2, pathways=[1,2])
Vmodel = dl.dipolarmodel(t, r, Pmodel=None, Bmodel=dl.bg_hom3d, experiment=expinfo)
-Manipulating the model
-***********************
+Using and modifying the model
+**********************************
-A full summary of the constructed model(s) can be inspected by printing the model object ::
+To obtain a summary of the constructed model, print it: ::
>>> print(Vmodel)
Description: Dipolar signal model
@@ -99,5 +116,4 @@ A full summary of the constructed model(s) can be inspected by printing the mode
========= ======= ======= ======= ======== ======== ====== ======================================
-From this point on, the model can be modified, manipulated and expanded freely as any other DeerLab model. Check out the :ref:`modeling guide ` for more details and instructions on model manipulation.
-
+Once defined, the model can be modified, manipulated and expanded freely as any other DeerLab model. Refer to the :ref:`modeling guide ` for more details and instructions on model manipulation.
diff --git a/docsrc/source/fitting_guide.rst b/docsrc/source/fitting_guide.rst
index 483a5348..69a718b2 100644
--- a/docsrc/source/fitting_guide.rst
+++ b/docsrc/source/fitting_guide.rst
@@ -150,9 +150,9 @@ The first step is to construct the penalty function. First, it must be a callabl
Third, a selection functional for the optimal selection of the penalty weight must be chosen from a list of built-in functionals:
- Informational complexity criterion (``icc``)
-- Akaike complexity criterion (``aic``)
-- Bayesian complexity criterion (``aic``)
-- Corrected Akaike complexity criterion (``aicc``)
+- Akaike information criterion (``aic``)
+- Bayesian information criterion (``bic``)
+- Corrected Akaike information criterion (``aicc``)
Last, the penalty can be constructed using the ``Penalty()`` constructor passing the penalty function and selection functional ::
diff --git a/docsrc/source/images/ex_dqc_pathways.png b/docsrc/source/images/ex_dqc_pathways.png
index fad90f52..abff9384 100644
Binary files a/docsrc/source/images/ex_dqc_pathways.png and b/docsrc/source/images/ex_dqc_pathways.png differ
diff --git a/docsrc/source/images/ex_dqc_pathways.svg b/docsrc/source/images/ex_dqc_pathways.svg
index 0e7111e5..c413d7c0 100644
--- a/docsrc/source/images/ex_dqc_pathways.svg
+++ b/docsrc/source/images/ex_dqc_pathways.svg
@@ -1,15 +1,15 @@
