Skip to content

Commit

Permalink
Add Ridme example and changes examples from deadtime to tmin (#446)
Browse files Browse the repository at this point in the history
* Add basic RIDME example

* Change from deadtime to tmin

* Change RIDME stretch value

* Replaced all Acquisition Deadtime with Start time.

* Missed one example

* Bug fix with propagate and evaluate functions

The propagate and evaluate functions had a bug that stopped them working on functions.
  • Loading branch information
HKaras authored Aug 15, 2023
1 parent 92a3bd1 commit 4dcdcc2
Show file tree
Hide file tree
Showing 30 changed files with 310 additions and 103 deletions.
13 changes: 12 additions & 1 deletion deerlab/fit.py
Original file line number Diff line number Diff line change
Expand Up @@ -358,6 +358,8 @@ def fit(model_, y, *constants, par0=None, penalties=None, bootstrap=0, noiselvl=
Fitted parameter vector ordered according to the model parameter indices.
paramUncert : :ref:`UQResult`
Uncertainty quantification of the parameter vector ordered according to the model parameter indices.
paramlist : list
List of the fitted parameter names ordered according to the model parameter indices.
model : ndarray
Fitted model response.
regparam : scalar
Expand Down Expand Up @@ -505,6 +507,15 @@ def bootstrap_fcn(ysim):
FitResult_dict = {key: getattr(fitresults,key) for key in ['y','mask','param','paramUncert','model','cost','plot','residuals','stats','regparam','regparam_stats','__plot_inputs']}
_paramlist = model._parameter_list('vector')

param_idx = [[] for _ in _paramlist]
idxprev = 0
for islinear in [False,True]:
for n,param in enumerate(_paramlist):
if np.all(getattr(model,param).linear == islinear):
N = len(np.atleast_1d(getattr(model,param).idx))
param_idx[n] = np.arange(idxprev,idxprev + N)
idxprev += N

# Enforce normalization of the linear parameters (if needed) for the final output
FitResult_param_,FitResult_paramuq_ = FitResult_param.copy(),FitResult_paramuq.copy()
if normalization:
Expand All @@ -526,7 +537,7 @@ def _scale(x):
noiselvl = noiselvl[0]

# Generate FitResult object from all the dictionaries
fitresult = FitResult({**FitResult_param_,**FitResult_paramuq_, **FitResult_dict,'penweights':penweights,'noiselvl':noiselvl})
fitresult = FitResult({**FitResult_param_,**FitResult_paramuq_, **FitResult_dict,'penweights':penweights,'noiselvl':noiselvl,'paramlist':_paramlist, '_param_idx':param_idx})

fitresult._summary = _print_fitresults(fitresult,model)

Expand Down
109 changes: 75 additions & 34 deletions deerlab/fitresult.py
Original file line number Diff line number Diff line change
Expand Up @@ -87,45 +87,76 @@ def _extarct_params_from_model(self, model):
if callable(model):
try:
modelparam = model._parameter_list('vector')
function_type=False
except AttributeError:
function_type=True
modelparam = inspect.getfullargspec(model).args

if not hasattr(self,'param'):
raise ValueError('The fit object does not contain any fitted parameters.')

# Enforce model normalization
normfactor_keys = []
for key in model._parameter_list():
param = getattr(model,key)
if np.all(param.linear):
if param.normalization is not None:
normfactor_key = f'{key}_scale'
normfactor_keys.append(normfactor_key)
try:
model.addnonlinear(normfactor_key,lb=-np.inf,ub=np.inf,par0=1,description=f'Normalization factor of {key}')
getattr(model,normfactor_key).freeze(1)
except KeyError:
pass
# # Enforce model normalization
# normfactor_keys = []
# for key in modelparam:
# param = getattr(model,key)
# if np.all(param.linear):
# if param.normalization is not None:
# normfactor_key = f'{key}_scale'
# normfactor_keys.append(normfactor_key)
# try:
# model.addnonlinear(normfactor_key,lb=-np.inf,ub=np.inf,par0=1,description=f'Normalization factor of {key}')
# getattr(model,normfactor_key).freeze(1)
# except KeyError:
# pass


# Get some basic information on the parameter vector
modelparam = model._parameter_list(order='vector')
param_idx = [[] for _ in model._parameter_list('vector')]
idxprev = 0
for islinear in [False,True]:
for n,param in enumerate(model._parameter_list('vector')):
if np.all(getattr(model,param).linear == islinear):
N = len(np.atleast_1d(getattr(model,param).idx))
param_idx[n] = np.arange(idxprev,idxprev + N)
idxprev += N

params = {key : fitvalue if len(fitvalue)>1 else fitvalue[0] for key,fitvalue in zip(modelparam,[self.param[idx] for idx in param_idx])}
# # Get some basic information on the parameter vector
# modelparam = model._parameter_list(order='vector')
# param_idx = [[] for _ in model._parameter_list('vector')]
# idxprev = 0
# for islinear in [False,True]:
# for n,param in enumerate(model._parameter_list('vector')):
# if np.all(getattr(model,param).linear == islinear):
# N = len(np.atleast_1d(getattr(model,param).idx))
# param_idx[n] = np.arange(idxprev,idxprev + N)
# idxprev += N

fitparams = {key : fitvalue if len(fitvalue)>1 else fitvalue[0] for key, fitvalue in zip(self.paramlist,[self.param[idx] for idx in self._param_idx])}
# Check that all parameters are in the fit object
for param in modelparam:
if not param in params:
if not param in fitparams:
raise KeyError(f'The fit object does not contain the {param} parameter.')


params = {param : fitparams[param] for param in modelparam}
params_idx = [self._param_idx[self.paramlist.index(param)] for param in modelparam]

return modelparam, params, params_idx

def _extract_params_from_function(self,function):
"""
Extracts the fitted parameters from a callable function.
Assumes that all parameters are length 1.
"""
# Get the parameter names from the function definition
modelparam = inspect.getfullargspec(function).args

fitparam_idx = self._param_idx

# Get the parameter values from the fit object
fitparams = {param : self.param[i] for i,param in enumerate(self.paramlist)}
params = {param : fitparams[param] for param in modelparam}
params_idx = [fitparam_idx[self.paramlist.index(param)] for param in modelparam]
# fitparams = {key : fitvalue if len(fitvalue)>1 else fitvalue[0] for key, fitvalue in zip(modelparam,[self.param[idx] for idx in fitparam_idx])}




return modelparam, params, params_idx


return modelparam, params, param_idx

def evaluate(self, model, *constants):
# ----------------------------------------------------------------------------
Expand Down Expand Up @@ -154,9 +185,14 @@ def evaluate(self, model, *constants):
response : array_like
Model response at the fitted parameter values.
"""

modelparam, params, _ = self._extarct_params_from_model(model)
parameters = {param: params[param] for param in modelparam}
try:
modelparam = model._parameter_list('vector')
modelparam, fitparams, fitparam_idx = self._extarct_params_from_model(model)
except AttributeError:
modelparam, fitparams, fitparam_idx = self._extract_params_from_function(model)


parameters = {param: fitparams[param] for param in modelparam}

# Evaluate the input model
response = model(*constants,**parameters)
Expand Down Expand Up @@ -195,11 +231,16 @@ def propagate(self, model, *constants, lb=None, ub=None):
Uncertainty quantification of the model's response.
"""

modelparam,_, param_idx = self._extarct_params_from_model(model)
# Determine the indices of the subset of parameters the model depends on
subset = [param_idx[np.where(np.asarray(modelparam)==param)[0][0]] for param in modelparam]
try:
modelparam = model._parameter_list('vector')
modelparam, fitparams, fitparam_idx = self._extarct_params_from_model(model)

except AttributeError:
modelparam, fitparams, fitparam_idx = self._extract_params_from_function(model)


# Propagate the uncertainty from that subset to the model
modeluq = self.paramUncert.propagate(lambda param: model(*constants,*[param[s] for s in subset]),lb,ub)
modeluq = self.paramUncert.propagate(lambda param: model(*constants,*[param[s] for s in fitparam_idx]),lb,ub)
return modeluq


Expand Down
5 changes: 3 additions & 2 deletions examples/advanced/ex_comparing_uncertainties.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,15 +27,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
t = t - t[0] # Account for zerotime
t = t + tmin

# Distance vector
r = np.arange(2,6,0.05) # nm
Expand Down
6 changes: 3 additions & 3 deletions examples/advanced/ex_dipolarpathways_selection.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,16 +25,16 @@
# Experimental parameters
tau1 = 0.5 # First inter-pulse delay, μs
tau2 = 3.5 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime

t = t - t[0] # Account for zerotime
t = t + tmin
# Construct the distance vector
r = np.arange(2,5,0.05)

Expand Down
5 changes: 3 additions & 2 deletions examples/advanced/ex_extracting_gauss_constraints.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,15 +23,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
t = t - t[0] # Account for zerotime
t = t + tmin

# Construct the dipolar signal model
r = np.arange(2,6,0.02)
Expand Down
5 changes: 3 additions & 2 deletions examples/advanced/ex_forcefield_fit.py
Original file line number Diff line number Diff line change
Expand Up @@ -48,15 +48,16 @@ def forcefield_P(c0,c1,c2,c3):
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 5.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
t = t - t[0] # Account for zerotime
t = t + tmin

# Construct the energy and distance distribution models
forcefield_energymodel = dl.Model(forcefield_energy)
Expand Down
11 changes: 6 additions & 5 deletions examples/advanced/ex_global_twostates_parametric.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@
# Experimental parameters
tau1 = 0.4 # First inter-pulse delay, μs
tau2 = 4.5 # Second inter-pulse delay, μs
deadtime = 0.2 # Acquisition deadtime, μs
tmin = 0.2 # Start time, μs

Vmodels, ts, Vexps = [], [], []
for file in files:
Expand All @@ -39,10 +39,11 @@
t, Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp / np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime

Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp / np.max(Vexp) # Rescaling (aesthetic)
t = t - t[0] # Account for zerotime
t = t + tmin

# Put the datasets into lists
ts.append(t)
Vexps.append(Vexp)
Expand Down
5 changes: 3 additions & 2 deletions examples/advanced/ex_identifiability_analysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,15 +23,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
t = t - t[0] # Account for zerotime
t = t + tmin

# Truncate the signal
Vexp_truncated = Vexp[t<=2]
Expand Down
4 changes: 2 additions & 2 deletions examples/advanced/ex_long_threespin_analysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@
files = [f'../data/triradical_protein_deer_{dB}dB.DTA' for dB in [0,6,9]]

# Experiment information
t0 = 0.280 # Acquisition deadtime, μs
t0 = 0.280 # Start time, μs
tau1 = 0.40 # First interpulse delay, μs
tau2 = 9.00 # Second interpulse delay, μs

Expand All @@ -34,7 +34,7 @@
t,Vexp, descriptor = dl.deerload(file,full_output=True)
t = t[:-80]
Vexp = Vexp[:-80]
# Adjust the start time
# Adjust the Start time
t = t - t[0] + t0

# Pre-processing
Expand Down
9 changes: 5 additions & 4 deletions examples/advanced/ex_multigauss_fitting_4pdeer.py
Original file line number Diff line number Diff line change
Expand Up @@ -24,15 +24,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t, Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp / np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp / np.max(Vexp) # Rescaling (aesthetic)
t = t - t[0] # Account for zerotime
t = t + tmin

# Maximal number of Gaussians in the models
Nmax = 5
Expand Down
5 changes: 3 additions & 2 deletions examples/advanced/ex_profileanalysis.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,15 +19,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
t = t - t[0] # Account for zerotime
t = t + tmin

# Distance vector
r = np.arange(2,5,0.05) # nm
Expand Down
5 changes: 3 additions & 2 deletions examples/basic/ex_bootstrapping.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,15 +30,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime
t = t - t[0] # Account for zerotime
t = t + tmin

# Distance vector
r = np.linspace(2,5,100) # nm
Expand Down
6 changes: 3 additions & 3 deletions examples/basic/ex_fitting_4pdeer.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,16 +21,16 @@
# Experimental parameters
tau1 = 0.3 # First inter-pulse delay, μs
tau2 = 4.0 # Second inter-pulse delay, μs
deadtime = 0.1 # Acquisition deadtime, μs
tmin = 0.1 # Start time, μs

# Load the experimental data
t,Vexp = dl.deerload(path + file)

# Pre-processing
Vexp = dl.correctphase(Vexp) # Phase correction
Vexp = Vexp/np.max(Vexp) # Rescaling (aesthetic)
t = t + deadtime # Account for deadtime

t = t - t[0] # Account for zerotime
t = t + tmin
# Distance vector
r = np.arange(2.5,5,0.01) # nm

Expand Down
Loading

0 comments on commit 4dcdcc2

Please sign in to comment.