matlab-fit-simbiology-model

# Fit SimBiology Models

Estimate parameters from data using `fitproblem`, fit population models
with NLME, generate virtual patients, and compute NCA metrics.

## When to Use

- "fit", "estimate", "calibrate" model parameters
- Parameter estimation from experimental/observed data
- Population PK/PD, NLME, mixed effects, inter-individual variability
- Virtual patients, virtual cohorts
- NCA, AUC, Cmax, Tmax, half-life, clearance
- Keywords: "fit", "estimate", "calibrate", "population", "NCA", "AUC"

## When NOT to Use

- Model construction or diagram (use `matlab-build-simbiology-model`)
- Simulation without fitting (use `matlab-simulate-simbiology-model`)
- Sensitivity analysis (use `matlab-simulate-simbiology-model`)

## Must-Follow Rules

### 1. Use `fitproblem` for parameter estimation

Always use `fitproblem` instead of calling `sbiofit` or `sbiofitmixed`
directly. `fitproblem` provides a unified, declarative interface:
```matlab
prob = fitproblem;
prob.Model = model;
prob.Data = data;
prob.ResponseMap = "Species = DataColumn";
prob.Estimated = estimatedInfo({'param'}, 'Bounds', [lo hi]);
results = fit(prob);
```
Do NOT call `sbiofit(model, data, ...)` or `sbiofitmixed(model, data, ...)`
directly — their positional argument signatures are error-prone.

### 2. Fitting requires `groupedData`, NOT a plain table

Always wrap data:
```matlab
data = groupedData(table(...));
data.Properties.IndependentVariableName = 'Time';
```

### 3. `ResponseMap` maps model outputs to data columns

Format is always `"ModelOutput = DataColumnName"`:

```matlab
% Single compartment — use species name on the left
prob.ResponseMap = "Drug = DrugConc";

% Multi-compartment — use qualified name to disambiguate
prob.ResponseMap = "Central.Drug = DrugConc";

% When species name matches data column name, still use the = format
prob.ResponseMap = "Drug = Drug";
```

Use the unqualified species name unless the same species name exists
in multiple compartments (then qualify with `Compartment.Species`).

### 4. Always set bounds

Prevent non-physical values (negative rates, etc.):
```matlab
estimParams = estimatedInfo({'ke','ka'}, ...
    'InitialValue', [0.2, 1.0], ...
    'Bounds', [0.01 1; 0.1 5]);
```

### 5. Use log transform for rate constants

Parameters spanning orders of magnitude (clearances, rate constants)
benefit from log-transform estimation. Set `.Transform` after creation:
```matlab
ei = estimatedInfo({'ke','ka'}, 'InitialValue', [0.1, 0.5], 'Bounds', [0.01 1; 0.1 5]);
ei(1).Transform = 'log';
ei(2).Transform = 'log';
```

Alternative: use `'log(param)'` name syntax (equivalent result):
```matlab
ei = estimatedInfo({'log(ke)','log(ka)'}, 'InitialValue', [0.1, 0.5], 'Bounds', [0.01 1; 0.1 5]);
```

**Important:** `InitialValue` and `Bounds` are always in the
**untransformed** (natural) domain. Do NOT pass `log(value)`.

Available transforms: `'log'`, `'logit'`, `'probit'`

Do NOT pass `'Transform'` as a name-value pair to the `estimatedInfo`
constructor — it errors. Always set the `.Transform` property after.

### 6. Error models for population fitting

Choose the error model that matches the noise structure:
- `'constant'` — absolute noise uniform
- `'proportional'` — noise scales with magnitude (most PK data)
- `'combined'` — both constant and proportional
- `'exponential'` — log-normal residual

### 7. NCA requires `sbioncaoptions` object

Do not use name-value pairs. Column names are camelCase.
EVDose column uses `NaN` for non-dose rows.

## Decision Table

| Scenario | Approach |
|----------|----------|
| Single subject or pooled fit | `fitproblem` with `FitFunction="sbiofit"` |
| Individual fits per subject | `fitproblem` with `Pooled=false` |
| Population NLME (IIV, random effects) | `fitproblem` with `FitFunction="sbiofitmixed"` |
| Model-independent PK metrics | `sbionca` |

## `fitproblem` Workflow (Preferred)

Use `fitproblem` for all parameter estimation. It provides a unified,
declarative interface that replaces direct calls to `sbiofit`/`sbiofitmixed`:

```matlab
% 1. Prepare data
data = groupedData(table(tSample, yData, 'VariableNames', {'Time','Drug'}));
data.Properties.IndependentVariableName = 'Time';

% 2. Define parameters with bounds
estimParams = estimatedInfo({'ke','ka'}, ...
    'InitialValue', [0.2, 1.0], ...
    'Bounds', [0.01 1; 0.1 5]);

% 3. Build the fit problem
prob = fitproblem;
prob.Model = model;
prob.Data = data;
prob.ResponseMap = "Drug = Drug";
prob.Estimated = estimParams;
prob.Doses = dose;                  % optional
prob.FunctionName = 'scattersearch';
prob.ProgressPlot = true;           % show live progress

% 4. Fit
results = fit(prob);

% 5. Inspect
disp(results.ParameterEstimates);
plot(results);
```

### Key `fitproblem` properties

| Property | Purpose |
|----------|---------|
| `Model` | The SimBiology model object |
| `Data` | `groupedData` table |
| `Estimated` | `estimatedInfo` object (**not** `EstimatedParameters`) |
| `ResponseMap` | Maps model species to data columns |
| `Doses` | Dose object(s) (**not** `Dose`) |
| `FitFunction` | `"sbiofit"` (default) or `"sbiofitmixed"` |
| `FunctionName` | Algorithm: `'scattersearch'`, `'nlinfit'`, `'fminsearch'`, `'lsqnonlin'`, `'particleswarm'` |
| `ProgressPlot` | `true` to show live fitting progress |
| `UseParallel` | `true` for parallel evaluation |
| `Pooled` | `true`/`false`/`"auto"` (sbiofit only) |
| `ErrorModel` | `"constant"`, `"proportional"`, `"combined"`, `"exponential"` |
| `Variants` | Variants to apply during fitting |

**Common property name mistakes:** `prob.Estimated` (not `EstimatedParameters`),
`prob.Doses` (not `Dose`), `prob.FunctionName` (not `Algorithm` or `Method`).

### Estimation algorithms

| Method | Use Case |
|--------|----------|
| `'scattersearch'` | Built-in global search, no extra toolbox — **start here** |
| `'nlinfit'` | Default local; smooth problems |
| `'lsqnonlin'` | Bounded least squares (Optimization Toolbox) |
| `'fminsearch'` | Derivative-free, simple problems |