Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
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