fit.identifiability¶

Parameter identifiability by profile likelihood.

A fit gives the parameters which describe the data best. Identifiability asks how well the data determines every one of them: a parameter is identifiable if the data constrains it to a finite interval, and non-identifiable if it can be changed without the model describing the data worse.

The profile likelihood answers this by scanning every parameter around the optimum. The scanned parameter is fixed at a value and all other parameters are optimized again, so the profile is the best cost the model reaches with the parameter at that value (Raue et al. 2009):

PL(θi) = min over θj≠i of cost(θ)

The cost of sbmlsim is the cost of scipy.optimize.least_squares, i.e., half the sum of the squared weighted residuals, cost = 0.5 * Σ r². With residuals which are standardized by the errors of the data, 2 * cost is the negative log-likelihood up to a constant, and the profile is compared with the threshold of the likelihood ratio test (Raue et al. 2009, Kreutz et al. 2013):

cost_threshold = cost_min + chi2.ppf(alpha, df) / 2

with df = 1 for the pointwise confidence intervals of single parameters, i.e., a threshold of 1.92 above the minimal cost at 95% confidence, and df = #parameters for simultaneous intervals. The values of a parameter at which the profile crosses the threshold are the bounds of its confidence interval. The intervals are invariant under a transformation of the parameters and may be asymmetric, which is where the intervals of the Fisher information matrix fail for the non-linear models of systems biology (Wieland et al. 2021).

The shape of the profile classifies the parameter (Raue et al. 2009):

• identifiable: the profile crosses the threshold on both sides of the optimum, the confidence interval is finite,
• practically non-identifiable: the profile has a minimum but stays below the threshold on one or both sides, i.e., up to the bound of the parameter. The data does not determine the parameter towards small and/or large values,
• structurally non-identifiable: the profile is flat over the whole scanned range, the parameter is compensated by the other parameters and the data carries no information about it.

The scans run in logarithmic parameter space, as the optimization does. The step along a profile is adaptive (Schälte et al. 2023): a step which raises the cost by more than a fraction of the threshold is halved and repeated, a step which raises it by little is enlarged, so the profile is resolved where it changes. The other parameters start from the previous point of the profile (Simpson & Maclaren 2023) and their paths are stored, so a parameter which is coupled to the scanned one is seen in its path (Maiwald et al. 2016). A scan stops when the profile crosses the threshold, at the bound of the parameter or after max_points; a bound which is reached below the threshold is an open interval (Borisov & Metelkin 2020). A scan which finds a lower cost than the optimum reports it: the fit did not converge, and the threshold is taken relative to the lowest cost of all profiles.

The scans are independent, so they run in the worker pool of the fit runner, two scans per parameter. A parallel analysis needs the same if __name__ == "__main__": guard as a parallel fit.

label

is_identifiable

ProfileSettings(
    alpha=0.95,
    degrees_of_freedom=1,
    initial_step=0.1,
    min_step=0.01,
    max_step=1.0,
    step_factor=2.0,
    max_cost_fraction=0.2,
    max_points=50,
    flatness=0.05,
    reoptimize=True,
    optimizer_kwargs=(lambda: {"diff_step": 0.05})(),
)

delta

threshold(cost)

to_dict()

from_dict(d)

ParameterProfile(
    pid,
    values,
    costs,
    paths,
    converged,
    index_optimum,
    ci_lower=None,
    ci_upper=None,
    identifiability=None,
)

value_optimum

cost_optimum

cost_min

side(direction)

crossing(threshold, direction)

rise(direction)

evaluate(threshold, flatness_cost)

to_dict()

from_dict(d)

IdentifiabilityResult(
    opid,
    parameter_set,
    parameters,
    settings,
    fit_settings,
    cost,
    profiles,
    duration=0.0,
)

pids

cost_min

threshold

flatness_cost

n_identifiable

better_optimum

evaluate()

parameter(pid)

summary_df()

report(path=None, print_output=False)

to_dict()

from_dict(d)

to_json(path=None)

from_json(json_info)

cost_threshold(cost, alpha=0.95, df=1)

profile_likelihood(
    problem,
    settings,
    parameter_set,
    profile_settings=None,
    pids=None,
    n_cores=1,
    serial=False,
    show_progress=True,
)

plot_profiles(result, path=None, ncols=3)

plot_profile(result, pid, path=None)

plot_all(result, output_dir, image_format='svg')