Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

snlls: allow extrapenalty to take linear parameters as input #175

Merged
merged 2 commits into from
May 10, 2021
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
2 changes: 1 addition & 1 deletion deerlab/fitmodel.py
Original file line number Diff line number Diff line change
Expand Up @@ -568,7 +568,7 @@ def separable_nonlinear_lsq_analysis(Vexp):
prescales = [1 for V in Vexp]
Vexp_ = [Vexp[i]/prescales[i] for i in range(nSignals)]

def scale_constraint(nonlinpar):
def scale_constraint(nonlinpar,linpar):
# --------------------------------------------------------
penalty = np.zeros(nSignals)
for i in range(nSignals):
Expand Down
16 changes: 9 additions & 7 deletions deerlab/snlls.py
Original file line number Diff line number Diff line change
Expand Up @@ -87,9 +87,10 @@ def snlls(y, Amodel, par0, lb=None, ub=None, lbl=None, ubl=None, nnlsSolver='cvx
The regularization parameter can be manually specified by passing a scalar value
instead of a string. The default ``'aic'``.

custom_penalty: callable
Custom penalty function to impose upon the solution. Must return a vector to be
added to the residual vector.
extrapenalty: callable
Custom penalty function to impose upon the solution. Must take two inputs, a vector of non-linear parameters
and a vector of linear parameters, and return a vector to be added to the residual vector (``pen = fcn(pnonlin,plin)``).
The square of the penalty is computed internally.

alphareopt : float scalar, optional
Relative parameter change threshold for reoptimizing the regularization parameter
Expand Down Expand Up @@ -383,7 +384,7 @@ def ResidualsFcn(p):

# Compute residual from custom penalty
if callable(extrapenalty):
penres = extrapenalty(p)
penres = extrapenalty(p,linfit)
penres = np.atleast_1d(penres)
res = np.concatenate((res,penres))

Expand Down Expand Up @@ -434,10 +435,11 @@ def ResidualsFcn(p):
Jnonlin = Jacobian(ResidualsFcn,nonlinfit,lb,ub)

# Jacobian (linear part)
Jlin = np.zeros((len(res),len(linfit)))
Jlin[:len(y),:] = scales_vec[:,np.newaxis]*Amodel(nonlinfit)
Jlin = scales_vec[:,np.newaxis]*Amodel(nonlinfit)
if callable(extrapenalty):
Jlin = np.concatenate((Jlin, Jacobian(lambda plin: extrapenalty(nonlinfit,plin),linfit,lbl,ubl)))
if includeRegularization:
Jlin[len(res)-Nlin:,:] = reg_penalty(regtype, alpha, L, linfit, huberparam, Nnonlin)[1]
Jlin = np.concatenate((Jlin, reg_penalty(regtype, alpha, L, linfit, huberparam, Nnonlin)[1]))

# Full Jacobian
J = np.concatenate((Jnonlin,Jlin),axis=1)
Expand Down