Constrain General Eqm Parameters

robertdkirkby · June 9, 2025, 10:01am

You can now constrain the general eqm parameters (the ones named in GEPriceParamNames). This works for both infinite and finite horizon models (and is now explained in an Appendix in the Intro to OLG Models, following is pretty much copy-paste from there).

Say you are finding the wage, w, in general equilibrium. You know this has to be positive, and if you try to solve the model with a negative wage then it gives errors. You can set constraints on the general equilibrium parameters to deal with these situations. There are three possible constraints, all of which are applied by specifying the name of the parameter you want to constrain. In the following, we will use w as the name of the parameter (in GEPriceParamNames) that we wish to constain.

You can constrain a parameter to be positive using
heteroagentoptions.constrainpositive={‘w’}
You can constrain a parameter to be between 0 and 1 (e.g., because it is a probability) using
heteroagentoptions.constrain0to1={‘w’}
You can constrain a parameter to be between A and B using
heteroagentoptions.constrainAtoB={‘w’}
and then you specify the values of A and B, e.g. that w is between 3 and 5, with
heteroagentoptions.constrainAtoBlimits.w=[3,5]
[Use name of parameter as field, ‘.w’, and then use a vector of two values for the range to apply to that parameter]

robertdkirkby · June 9, 2025, 10:06am

Note that this is the exact same setup that you use to constrain parameters when doing calibration or GMM estimation, except of course that in those cases you put these into caliboptions and estimoptions, respectively, rather than heteroagentoptions.

aledinola · June 26, 2025, 1:47pm

Small question: why the case 0-1 is specified separately from A-B?

robertdkirkby · June 27, 2025, 9:25am

Essentially because A-B works by converting to 0-1 and then using that. So I had to implement 0-1 as part of implementing A-B, and then seemed no reason not to just give user direct access to 0-1 as even though it is a subcase of A-B it is a commonly used one, e.g. for probabilities.

Rest of this just gives unnecessary detail that explains why implementing A-B works via 0-1…

Internally, ‘constrain positive’ replaces \theta (which is \in [0,Inf)) with \hat{\theta}\equiv log(\theta) (which is \in \mathbb{R}), and then performs the unconstrained optimization using \hat{\theta} and in each iteration inside the optimization just evaluates \theta=exp(\theta) (which is \in (0,Inf)) before passing it to the model.

Similarly, ‘constrain 0 to 1’ replaces \theta (which is \in [0,1]) with \hat{\theta}\equiv log(\frac{\theta}{1-\theta}) (which is \in \mathbb{R}) and then performs the unconstrained optimization using \hat{\theta} and in each iteration inside the optimization just evaluates \theta=1/(1+exp(-\hat{\theta})) (which is \in (0,1)) before passing it to the model.

Finally, ‘constrain A to B’ first replaces \theta (which is \in [A,B]) with \tilde{\theta}\equiv \frac{\theta - A}{B-A} (which is \in [0,1]), and then does the same as the constrain 0 to 1 (so \hat{\theta}\equiv log(\frac{\tilde{\theta}}{1-\tilde{\theta}}) (which is \in \mathbb{R}) and does the unconstrained optimization using \hat{\theta} and in each iteration inside the optimization just first evaluates \tilde{\theta}=1/(1+exp(-\hat{\theta})) (which is \in (0,1)) and then \theta=A+(B-A)\hat{\theta} (which is \in (A,B)) before passing it to the model.

Note, one minor weakness of all this, is that the constrained parameters can’t quite hit the end points. So if, e.g., you constrain 0-to-1 and the true solution is 0 (is 1) you will get a solution like 0.00001 (like 0.99999). You can see this in each paragraph above as the sets are closed (square-brackets [,]) to begin, and open (round-brackets (,)) at the end.

aledinola · June 27, 2025, 10:14am

Thanks for the answer!

"in each iteration inside the optimization just evaluates \hat{\theta} = \exp(\theta) "

What about possible overflows? I used the exp log transformation in the past but if you have large numbers, the exp of a large number might become undefined. For example, on my Matlab, exp(x)=inf if x>=710

John d’Errico, in his implementation of fminsearchbnd, uses a more robust transformation. This post on the Matlab forum might also be useful.

Not suggesting to change this, just something to be aware of.

robertdkirkby · June 27, 2025, 10:32am

I kind of just fudge these internally. If your x goes above something very big I just do a min(x,C) where C is something large, before doing the operations like ‘take exponential’; so effectively x is capped at C, and thus exp(x) is capped at exp(C) hence avoiding the Infs. There is of course a loss but it is pretty minuscule as is simply implicitly puts that the max value that can be taken by the variable that was constrained to be positive is capped at about 5x10^21 (I cannot think of an economic application where you need a parameter value larger than this, would tend to be a sign the model is not well constructed; even a model of global GDP measured in cents won’t get close to this ).

I have no objection to an improved implementation, but it is simply what I came up with at the time and works for everything I have tried so far.

You can see exactly how the handling of the constraints are done in the first part of the ‘objective function’ used for calibration/estimation. Currently it is lines 6-32ish (looks like C is 50 which gives the exp(C)=5x10^21 described above). This is the ‘convert \hat{\theta} back into \theta’ part described in my previous post:

github.com/vfitoolkit/VFIToolkit-matlab

Estimation/ObjectiveFn/CalibrateLifeCycleModel_objectivefn.m

master

function Obj=CalibrateLifeCycleModel_objectivefn(calibparamsvec, CalibParamNames,n_d,n_a,n_z,N_j,d_grid, a_grid, z_gridvals_J, pi_z_J, ReturnFn, ReturnFnParamNames, Parameters, DiscountFactorParamNames, jequaloneDist,AgeWeightParamNames, ParametrizeParamsFn, FnsToEvaluate, FnsToEvaluateParamNames, usingallstats, usinglcp,targetmomentvec, allstatmomentnames, acsmomentnames, allstatcummomentsizes, acscummomentsizes, AllStats_whichstats, ACStats_whichstats, calibparamsvecindex, calibomitparams_counter, calibomitparamsmatrix, caliboptions, vfoptions,simoptions)
% Note: Inputs are CalibParamNames,TargetMoments, and then everything
% needed to be able to run ValueFnIter, StationaryDist, AllStats and
% LifeCycleProfiles. Lastly there is caliboptions.

% Do any transformations of parameters before we say what they are
penalty=zeros(length(calibparamsvec),1); % Used to apply penalty to objective function when parameters try to leave restricted ranges
for pp=1:length(CalibParamNames)
    if caliboptions.constrainpositive(pp)==1 % Forcing this parameter to be positive
        temp=calibparamsvec(calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1));
        penalty((calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1)))=abs(temp/50).*(temp<-51); % 1 if out of range [Note: 51, rather than 50, so penalty only hits once genuinely out of range]
        % Constrain parameter to be positive (be working with log(parameter) and then always take exp() before inputting to model)
        calibparamsvec(calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1))=exp(calibparamsvec(calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1)));
    elseif caliboptions.constrain0to1(pp)==1
        temp=calibparamsvec(calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1));
        penalty((calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1)))=abs(temp/50).*((temp>51)+(temp<-51)); % 1 if out of range [Note: 51, rather than 50, so penalty only hits once genuinely out of range]
        % Constrain parameter to be 0 to 1 (be working with x=log(p/(1-p)), where p is parameter) then always take 1/(1+exp(-x)) before inputting to model
        calibparamsvec(calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1))=1/(1+exp(-calibparamsvec(calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1))));
        % Note: This does not include the endpoints of 0 and 1 as 1/(1+exp(-x)) maps from the Real line into the open interval (0,1)
        %       R is not compact, and [0,1] is compact, so cannot have a continuous bijection (one-to-one and onto) function from R into [0,1].

This file has been truncated. show original

MichaelTiemann · June 1, 2026, 3:08am

I’m having trouble with fminalgo=5 and constrainpositive. I suspect there’s a missing translation from \theta to \hat\theta somewhere. I’m sending a calibrated price of 3.2832 and my evaluation function is seeing a price of 26.661. I never saw this behavior with fminalgo=1 hence the blame.

aledinola · June 1, 2026, 9:48pm

According to my AI, there are a few bugs.

Proofreading notes: parameter-constraint transformations in VFI Toolkit calibration

Focus: parameter transformations used to handle constraints in the calibration routines.

Relevant files:

SubCodes/ParameterConstraints/ParameterConstraints_TransformParamsToUnconstrained.m
SubCodes/ParameterConstraints/ParameterConstraints_TransformParamsToOriginal.m
SubCodes/ParameterConstraints/ParameterConstraints_PType_TransformParamsToUnconstrained.m
Estimation/Calibration/CalibrateLifeCycleModel.m
Estimation/Calibration/CalibrateLifeCycleModel_PType.m
Estimation/Calibration/CalibrateOLGModel_PType.m
Estimation/Calibration/CalibrateBIHAModel_PType.m
Estimation/Calibration/CalibrateInfHorzAgentModel_PType.m

1. Serious bug: PType constrained parameters are transformed back using the non-PType routine

In PType calibration, the initial vector is transformed using the PType-specific routine:

[calibparamsvec0,caliboptions] = ...
    ParameterConstraints_PType_TransformParamsToUnconstrained( ...
        calibparamsvec0, calibparamsvecindex, CalibParamNames, ...
        nCalibParamsFinder, caliboptions, 1);

This is correct, because PType calibration expands original parameter names into PType-specific blocks.

However, at the end of the PType calibration wrappers, the code transforms back using the non-PType routine:

[calibparamsvec,penalty] = ...
    ParameterConstraints_TransformParamsToOriginal( ...
        calibparamsvec, calibparamsvecindex, CalibParamNames, caliboptions);

This is likely wrong.

The non-PType inverse routine loops over:

for pp = 1:length(CalibParamNames)

But in PType calibration, the relevant number of blocks is:

nCalibParams = size(nCalibParamsFinder,1);

These are not the same object. CalibParamNames contains original parameter names, while nCalibParamsFinder indexes the PType-expanded parameter blocks.

Therefore, if a constrained parameter is PType-specific, some entries may remain in log/logit space instead of being transformed back to the original model scale.

Suggested fix

Add a PType-specific inverse transform:

function [calibparamsvec,penalty] = ...
    ParameterConstraints_PType_TransformParamsToOriginal( ...
        calibparamsvec, calibparamsvecindex, CalibParamNames, ...
        nCalibParamsFinder, caliboptions)

nCalibParams = size(nCalibParamsFinder,1);

penalty = zeros(length(calibparamsvec),1);

for pp = 1:nCalibParams

    idx = calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1);

    if caliboptions.constrainpositive(pp)==1

        temp = calibparamsvec(idx);

        penalty(idx) = abs(temp/50).*(temp < -51);

        calibparamsvec(idx) = exp(temp);

    elseif caliboptions.constrain0to1(pp)==1

        temp = calibparamsvec(idx);

        penalty(idx) = abs(temp/50).*((temp > 51) + (temp < -51));

        calibparamsvec(idx) = 1./(1+exp(-temp));

        calibparamsvec(idx) = calibparamsvec(idx).*(temp > -50);
        calibparamsvec(idx) = calibparamsvec(idx).*(1-(temp > 50)) + (temp > 50);

    end

    if caliboptions.constrainAtoB(pp)==1

        A = caliboptions.constrainAtoBlimits(pp,1);
        B = caliboptions.constrainAtoBlimits(pp,2);

        calibparamsvec(idx) = A + (B-A)*calibparamsvec(idx);

    end

end

if sum(penalty)>0
    penalty = 1/prod(1./penalty(penalty>0));
else
    penalty = 0;
end

end

Then replace, in PType wrappers:

[calibparamsvec,penalty] = ...
    ParameterConstraints_TransformParamsToOriginal( ...
        calibparamsvec, calibparamsvecindex, CalibParamNames, caliboptions);

with:

[calibparamsvec,penalty] = ...
    ParameterConstraints_PType_TransformParamsToOriginal( ...
        calibparamsvec, calibparamsvecindex, CalibParamNames, ...
        nCalibParamsFinder, caliboptions);

Apply this to:

CalibrateLifeCycleModel_PType.m
CalibrateOLGModel_PType.m
CalibrateBIHAModel_PType.m
CalibrateInfHorzAgentModel_PType.m

2. Serious bug in `CalibrateLifeCycleModel.m`: extra `usingallstats` under `fminalgo==8`

In the non-lsqnonlin branch, the objective function receives:

..., usingallstats, usinglcp, usingcustomstats, ...

But in the fminalgo==8 branch, it receives:

..., usingallstats, usingallstats, usinglcp, usingcustomstats, ...

There is an extra usingallstats.

Since fminalgo==8 is the default, this is important.

Suggested fix

Change:

..., FnsToEvaluateParamNames, usingallstats, usingallstats, usinglcp, usingcustomstats, targetmomentvec, ...

to:

..., FnsToEvaluateParamNames, usingallstats, usinglcp, usingcustomstats, targetmomentvec, ...

3. Add checks for positive constraints

Current logic:

calibparamsvec(idx) = max(log(calibparamsvec(idx)),-49.99);

This handles p=0, but not p<0. For p<0, log(p) is complex.

Suggested fix

idx = calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1);
p = calibparamsvec(idx);

if any(p < 0)
    error('Initial guess for positive-constrained parameter is negative.');
end

calibparamsvec(idx) = max(log(p), -49.99);

Apply this in:

ParameterConstraints_TransformParamsToUnconstrained.m
ParameterConstraints_PType_TransformParamsToUnconstrained.m

4. Add checks for 0-to-1 constraints

Current logic:

calibparamsvec(idx) = min(49.99, max(-49.99, ...
    log(calibparamsvec(idx)./(1-calibparamsvec(idx)))));

This handles exactly p=0 and p=1, but not p<0 or p>1.

Suggested fix

idx = calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1);
p = calibparamsvec(idx);

if any(p < 0 | p > 1)
    error('Initial guess for 0-to-1 constrained parameter is outside [0,1].');
end

calibparamsvec(idx) = min(49.99, max(-49.99, log(p./(1-p))));

Apply this in:

ParameterConstraints_TransformParamsToUnconstrained.m
ParameterConstraints_PType_TransformParamsToUnconstrained.m

5. Add checks for A-to-B constraints

Current logic first maps:

p01 = (p-A)/(B-A);

then applies the 0-to-1 logit transform.

This is valid only if:

B > A
initial p is in [A,B]

Suggested fix

idx = calibparamsvecindex(pp)+1:calibparamsvecindex(pp+1);

A = caliboptions.constrainAtoBlimits(pp,1);
B = caliboptions.constrainAtoBlimits(pp,2);

if B <= A
    error('constrainAtoB upper bound must be greater than lower bound.');
end

p = calibparamsvec(idx);

if any(p < A | p > B)
    error('Initial guess for A-to-B constrained parameter is outside [A,B].');
end

calibparamsvec(idx) = (p-A)/(B-A);

Apply this in:

ParameterConstraints_TransformParamsToUnconstrained.m
ParameterConstraints_PType_TransformParamsToUnconstrained.m

6. Penalty is not used in the `lsqnonlin` residual vector

The inverse transform computes a penalty when transformed parameters move too far beyond the artificial cutoffs.

For scalar objectives, the penalty is used. But for lsqnonlin, the objective returns a vector of residuals with caliboptions.vectoroutput==2. The penalty is not appended to that vector.

Suggested fix

In the vectoroutput==2 branch of each objective function, after constructing Obj, add:

if penalty > 0
    Obj = [Obj; sqrt(penalty)];
end

This is natural because lsqnonlin minimizes the sum of squared residuals.

Apply this to relevant objective functions:

CalibrateLifeCycleModel_objectivefn.m
CalibrateOLGModel_*_objectivefn.m
CalibrateBIHAModel_*_objectivefn.m
CalibrateInfHorzAgentModel_objectivefn.m
PType objective functions

7. The inverse transformation formula is conceptually fine

The non-PType inverse transform uses:

p = exp(x)

for positive constraints,

p = 1./(1+exp(-x))

for 0-to-1 constraints, and

p = A + (B-A)*p01

for A-to-B constraints.

These formulas are correct.

The main problems are:

PType wrappers call the non-PType inverse transform.
Initial guesses are not checked before applying log or logit.
lsqnonlin does not use the penalty in the residual vector.
CalibrateLifeCycleModel.m has a duplicated usingallstats argument under fminalgo==8.

Summary of recommended changes

Highest priority:

Add ParameterConstraints_PType_TransformParamsToOriginal.m.
Use it in all PType calibration wrappers.
Fix duplicated usingallstats in CalibrateLifeCycleModel.m under fminalgo==8.

Medium priority:

Add explicit initial-guess checks for:
- positive constraints: parameter must be non-negative
- 0-to-1 constraints: parameter must be in [0,1]
- A-to-B constraints: parameter must be in [A,B]
- A-to-B limits: require B>A

Lower priority:

Add a penalty residual in the vectoroutput==2 / lsqnonlin branch.

aledinola · June 1, 2026, 9:51pm

@robertdkirkby Are you sure this function works?

In particular in line 40 z_grid and pi_z are undefined. So if caliboptions.calibrateshocks==1, the function will fail.
In line 71,
rng(estimoptions.rngindex)
will fail since estimoptions is not defined in the function nor passed as input argument (I think).

robertdkirkby · June 1, 2026, 9:52pm

I am sure that fminalgo=5 does not handle parameter constraints correctly, if that is what you mean

I knew why as soon as Michael mentioned it without checking the code (I was not aware it was broken, but as soon as Michael said it was obvious why without even checking). Will fix later this week. I want to give fminalgo=5 (shooting for stationary general eqm) its own function so that I can fix all the heteroagentstationaryeqm commands in one go and make the shooting easier to maintain going forward.

aledinola · June 1, 2026, 9:55pm

Ah sorry my post was incomplete. The bug that I spotted is related to something else.

robertdkirkby · June 1, 2026, 9:57pm

Serious bug: PType constrained parameters are transformed back using the non-PType routine

as I recall this is not a bug. this is a non-issue, since the ‘transform back’ is just a vector before and after anyway. Unless someone wants to demonstrate this is actually a bug I am going to ignore it.

Serious bug in CalibrateLifeCycleModel.m: extra usingallstats under fminalgo==8

This is a real error, now fixed.

Penalty is not used in the lsqnonlin residual vector

Minor issue, but reasonable and easy to fix.

3, 4 & 5

are “if constrain 0 to 1 but parameter is -3, throw error telling user”, and same for “if constrain A to B…” and “if constrain positive”

Fixed all these in

github.com/vfitoolkit/VFIToolkit-matlab

Estimation/Calibration: fix usingallstats arg shift + validate constraint inputs + lsqnonlin penalty residual

committed 03:21AM - 02 Jun 26 UTC

robertdkirkby

+107 -7

- CalibrateLifeCycleModel.m: remove the duplicate `usingallstats` in the fmina…lgo==8 (lsqnonlin, the default) branch of the objective-function lambda. The extra arg was shifting every subsequent positional input by one and misaligning targetmomentvec/momentnames/whichstats/FnsToEvaluate_* in CalibrateLifeCycleModel_objectivefn. - ParameterConstraints_(PType_)TransformParamsToUnconstrained.m: hoist pp_index, read p_val=calibparamsvec(pp_index) once per block, and add input-validation errors before each transform. Positive: error if p_val<0 (was producing complex log). 0-to-1: error if p_val outside [0,1] (was producing complex log or divide-by-zero). A-to-B: error if B<=A, and error if p_val outside [A,B]. The error fprintf reports which parameter triggered the failure (using CalibParamNames{nCalibParamsFinder(pp,1)} in the PType variant). Transform algebra itself unchanged. - 12 Calibrate*_objectivefn.m (LifeCycle/InfHorzAgent/OLG/BIHA x plain/ PType x Joint/Nested where applicable): in the vectoroutput==2 (lsqnonlin residual-vector) branch, after the existing gather(), append sqrt(penalty) when penalty>0 else 0. The penalty was previously computed but discarded in this branch. Appending unconditionally (with the else-zero) keeps the residual length constant across iterations, which lsqnonlin requires for a stable Jacobian shape. Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>

aledinola · June 1, 2026, 10:13pm

I was reviewing the toolkit folder VFIToolkit-matlab/Estimation at master · vfitoolkit/VFIToolkit-matlab · GitHub and I have found two sets of bugs (I think).

First, there are two bugs in Estimation/ObjectiveFn:

In line 40 z_grid and pi_z are undefined. So if caliboptions.calibrateshocks==1, the function will fail.
In line 71, rng(estimoptions.rngindex)
will fail since estimoptions is not defined in the function nor passed as input argument.

Second, there are some likely bugs flagged in my longer post. It is not easy for me to summarize further the post: you can have a look at Sections 1 (PType constrained parameters are transformed back using the non-PType routine) and 2 (bug CalibrateLifeCycleModel.m: extra usingallstats under fminalgo==8) which explain the bugs. The rest are imprecisions flagged by GPT, so you can ignore them.

Hope this is helpful!

robertdkirkby · June 2, 2026, 3:49am

Cleaned these up in

github.com/vfitoolkit/VFIToolkit-matlab

Calibrate*_objectivefn: fix calibrateshocks==1 codepath (undefined z_grid/pi_z + estimoptions + z_grid_J typo)

committed 03:46AM - 02 Jun 26 UTC

robertdkirkby

+6 -6

The calibrateshocks==1 branch in three non-PType FHorz objective functions was u…nreachable: it referenced variables that aren't in scope, so any attempt to recalibrate shocks would error out immediately. - CalibrateLifeCycleModel_objectivefn / CalibrateOLGModel_Joint_objectivefn / CalibrateOLGModel_Nested_objectivefn: line 40 (43 for OLG_Nested) was passing `z_grid, pi_z` to ExogShockSetup_FHorz, but those names don't exist in the function signature — only `z_gridvals_J, pi_z_J` are received. Switch to feeding the joint-grid forms back in (which ExogShockSetup_FHorz accepts as an input shape), matching the existing PType siblings' pattern. - CalibrateLifeCycleModel_objectivefn line 71: `rng(estimoptions.rngindex)` referenced `estimoptions`, which is neither a function input nor defined anywhere. Rename to `caliboptions.rngindex` (and the comment above it), consistent with this objective only ever receiving `caliboptions`. - CalibrateLifeCycleModel_objectivefn line 72: `SimPanelValues_FHorz_Case1` was being passed `z_grid_J` (undefined) — typo for `z_gridvals_J`. Fixed. Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>

You are never supposed to use the sim panel approach, was just there for an test/paper showing the agent dist approach is better.

robertdkirkby · June 10, 2026, 4:14am

Hopefully should be fixed.

Created a subfn to handle the fminalgo=5, and then modified it to handle parameter constraints.

github.com/vfitoolkit/VFIToolkit-matlab

fminalgo5 shooting into shared subcodes, fix parameter constraint handling

committed 04:13AM - 10 Jun 26 UTC

robertdkirkby

+123 -266

Extract the inline fminalgo==5 shooting loop from all six HeteroAgentStationaryE…qm commands into StationaryGeneralEqm_subcode_fminalgo5 (non-PType: Case1_FHorz, InfHorz) and StationaryGeneralEqm_subcode_fminalgo5_PType (InfHorz_PType, Case1_FHorz_PType, FHorz_PType2L, MixHorz_PType). Previously the shooting loop ignored heteroagentoptions.constrainpositive/ constrain0to1/constrainAtoB: it fed original-space prices to an objective expecting unconstrained ones, and the final answer got transformed a second time. Now the objective is fed the unconstrained vector while the shooting update and convergence are done in original space, and the subcode returns the unconstrained vector so the caller's TransformParamsToOriginal round-trips correctly. PType subcode also fixes the Parameters update to use GEpriceindexes ranges so ptype-dependent prices get all their values. Note: no guard yet against the shooting update leaving the constrained region (e.g. constrainpositive price going negative -> log() goes complex). Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>

Constrain General Eqm Parameters

Proofreading notes: parameter-constraint transformations in VFI Toolkit calibration

1. Serious bug: PType constrained parameters are transformed back using the non-PType routine

Suggested fix

2. Serious bug in CalibrateLifeCycleModel.m: extra usingallstats under fminalgo==8

Suggested fix

3. Add checks for positive constraints

Suggested fix

4. Add checks for 0-to-1 constraints

Suggested fix

5. Add checks for A-to-B constraints

Suggested fix

6. Penalty is not used in the lsqnonlin residual vector

Suggested fix

7. The inverse transformation formula is conceptually fine

Summary of recommended changes

2. Serious bug in `CalibrateLifeCycleModel.m`: extra `usingallstats` under `fminalgo==8`

6. Penalty is not used in the `lsqnonlin` residual vector