Documentation Center

  • Trial Software
  • Product Updates

sbiofitmixed

Fit nonlinear mixed-effects model

Syntax

  • fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo) example
  • fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing) example
  • fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing,functionName) example
  • fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing,functionName,opt) example
  • fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing,functionName,opt,variants) example
  • fitResults = sbiofitmixed(_,'UseParallel',tf) example
  • [fitResults,simDataI,simDataP] = sbiofitmixed(_)

Description

    Note:   This function requires nlmefit and nlmefitsa in Statistics Toolbox™.

example

fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo) performs nonlinear mixed-effects estimation using the SimBiology® model sm and returns a NLMEResults object fitResults.

grpData is a groupedData object specifying the data to fit. responseMap is a string or cell array of strings that maps model components to response data in grpData. covEstiminfo is a CovariateModel object or an array of estimatedInfo objects that defines the parameters to be estimated.

If the model contains active doses and variants, they are applied during the simulation.

example

fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing) uses the dosing information specified by a matrix of SimBiology dose objects dosing instead of using the active doses of the model sm if there are any.

example

fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing,functionName) uses the estimation function specified by functionName that must be either 'nlmefit' or 'nlmefitsa'.

example

fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing,functionName,opt) uses the additional options specified by opt for the estimation function functionName.

example

fitResults = sbiofitmixed(sm,grpData,responseMap,covEstiminfo,dosing,functionName,opt,variants) applies variant objects specified as variants instead of using any active variants of the model.

example

fitResults = sbiofitmixed(_,'UseParallel',tf) provides an option to estimate parameters in parallel if Parallel Computing Toolbox™ is available.

[fitResults,simDataI,simDataP] = sbiofitmixed(_) returns a vector of results objects fitResults, vector of simulation results simDataI using individual-specific parameter estimates, and vector of simulation results simDataP using population parameter estimates.

    Note:   sbiofitmixed unifies sbionlmefit and sbionlmefitsa estimation functions. Use sbiofitmixed to perform nonlinear mixed-effects modeling and estimation.

Examples

expand all

Fit a One-Compartment PK Model to the Phenobarbital Data

This example uses data collected on 59 preterm infants given phenobarbital during the first 16 days after birth. Each infant received an initial dose followed by one or more sustaining doses by intravenous bolus administration. A total of between 1 and 6 concentration measurements were obtained from each infant at times other than dose times, for a total of 155 measurements. Infant weights and APGAR scores (a measure of newborn health) were also recorded. Data is described in [1], a study funded by the NIH/NIBIB grant P41-EB01975.

Load the data.

load pheno.mat ds

Convert the dataset to a groupedData object, a container for holding tabular data that is divided into groups. It can automatically identify commonly used variable names as the grouping variable or independent (time) variable. Display the properties of the data and confirm that GroupVariableName and IndependentVariableName are correctly identified as 'ID' and 'TIME', respectively.

data = groupedData(ds);
data.Properties
ans = 

                Description: ''
       VariableDescriptions: {}
              VariableUnits: {}
             DimensionNames: {'Observations'  'Variables'}
                   UserData: []
                   RowNames: {}
              VariableNames: {'ID'  'TIME'  'DOSE'  'WEIGHT'  'APGAR'  'CONC'}
          GroupVariableName: 'ID'
    IndependentVariableName: 'TIME'

Create a simple one-compartment PK model with bolus dosing and linear clearance to fit such data. Use the PKModelDesign object to construct the model. Each compartment is defined by a name, dosing type, a clearance type, and whether or not the dosing requires a lag parameter. After constructing the model, you can also get a PKModelMap object map that lists the names of species and parameters in the model that are most relevant for fitting.

pkmd = PKModelDesign;
pkmd.addCompartment('Central','DosingType','Bolus','EliminationType','linear-clearance','HasResponseVariable',true,'HasLag',false);
[onecomp, map] = pkmd.construct;

Describe the experimentally measured response by mapping the appropriate model component to the response variable. In other words, indicate which species in the model corresponds to which response variable in the data. The PKModelMap property Observed indicates that the relevant species in the model is Drug_Central, which represents the drug concentration in the system. The relevant data variable is CONC, which you visualized previously.

map.Observed
responseMap = 'Drug_Central = CONC';
ans = 

    'Drug_Central'

The parameters to estimate in this model are the volume of the central compartment Central and the clearance rate Cl_Central. The PKModelMap property Estimated lists these relevant parameters. The underlying algorithm of sbiofit assumes parameters are normally distributed, but this assumption may not be true for biological parameters that are constrained to be positive, such as volume and clearance. Specify a log transform for the estimated parameters so that the transformed parameters follow a normal distribution. Use an estimatedInfo object to define such transforms and initial values (optional).

map.Estimated
estimatedParams = estimatedInfo({'log(Central)','log(Cl_Central)'},'InitialValue',[1 1]);
ans = 

    'Central'
    'Cl_Central'

Each infant received a different schedule of dosing. The amount of drug is listed in the data variable DOSE. To specify these dosing during fitting, create dose objects from the data. These objects use the property TargetName to specify which species in the model receives the dose. In this example, the target species is Drug_Central, as listed by the PKModelMap property Dosed.

map.Dosed
ans = 

    'Drug_Central'

Create a sample dose with this target name and then use the createDoses method of groupedData object data to generate doses for each infant based on the dosing data DOSE.

sampleDose = sbiodose('sample','TargetName','Drug_Central');
doses = createDoses(data,'DOSE','',sampleDose);

Fit the model.

[nlmeResults,simI,simP] = sbiofitmixed(onecomp,data,responseMap,estimatedParams,doses,'nlmefit');

Visualize the fitted results using individual-specific parameter estimates.

plot(nlmeResults,'ParameterType','individual');

Visualize the fitted results using population parameter estimates.

plot(nlmeResults,'ParameterType','population');

Input Arguments

expand all

sm — SimBiology modelSimBiology model object

SimBiology model, specified as a SimBiology model object. The active configset object of the model contains solver settings for simulation. Any active doses and variants are applied to the model during simulation unless specified otherwise using the dosing and variants input arguments, respectively.

grpData — Data to fitgroupedData object

Data to fit, specified as a groupedData object.

The name of the time variable must be defined in the IndependentVariableName property of grpData. For instance, if the time variable name is 'TIME', then specify it as follows.

grpData.Properties.IndependentVariableName = 'TIME';

grpData must at least two groups, and the name of grouping variable name must be defined in the GroupVariableName property of grpData. For example, if the grouping variable name is 'GROUP', then specify it as follows.

grpData.Properties.GroupVariableName = 'GROUP';

A group usually refers to a set of measurements that represent a single time course, often corresponding to a particular individual or experimental condition.

responseMap — Mapping information of model components to response datastring | cell array of strings

Mapping information of model components to response data in grpData, specified as a string or cell array of strings.

Each string is an equation-like expression, similar to assignment rules in SimBiology. It contains the name (or qualified name) of a quantity (species, compartment, or parameter) in the model sm, followed by the character '=' and the name of a variable in grpData. For clarity, white spaces are allowed between names and '='.

For example, if you have the concentration data 'CONC' in grpData that corresponds to a model species 'Drug_Central', you can specify the mapping information as follows.

responseMap = 'Drug_Central = CONC';

To unambiguously name a species, use the qualified name, which includes the name of the compartment. To name a reaction-scoped parameter, use the reaction name to qualify the parameter. If the name is not a valid MATLAB® variable name, surround it by square brackets such as [reaction 1].[parameter 1].

An error is issued if any (qualified) name matches two components of the same type. However, you can use a (qualified) name that matches two components of different types, and the function first finds the species with the given name, followed by compartments and then parameters.

covEstiminfo — Estimated parametersvector of estimatedInfo objects | CovariateModel object

Estimated parameters, specified as a vector of estimatedInfo objects or a CovariateModel object that defines the estimated parameters in the model sm, their initial estimates (optional), and their relationship to group-specific covariates included in grpData (optional). If this is a vector of estimatedInfo objects, then no covariates are used, and all parameters are estimated with group-specific random effects.

You can also specify parameter transformations if necessary. Supported transforms are log, logit, and probit. For details, see estimatedInfo object and CovariateModel object.

dosing — Dosing information[] | 2-D matrix of dose objects

Dosing information, specified as an empty array or 2-D matrix of dose objects (ScheduleDose object or RepeatDose object). If empty, no doses are applied during simulation, even if the model has active doses. If not empty, the matrix must have a single row or one row per group in the data grpData. If it has a single row, the same doses are applied to all groups during simulation. If it has multiple rows, each row is applied to a separate group, in the same order as the groups appear in grpData.

Multiple columns are allowed so that you can apply multiple dose objects to each group. Each column of doses must reference the same components in the model sm. Specifically, all doses (except for empty doses) in a column must have the same values for TargetName, DurationParameterName, and LagParameterName. If some groups require more doses than others, then fill in the matrix with dummy doses that are either default doses or empty doses.

A default dose has default values for all properties, except for the Name property. An empty dose has a dose amount of 0, thus having no effect on the model. Create a default dose as follows.

d1 = sbiodose('d1');

In addition to manually constructing dose objects, if grpData has dosing information, you can use the createDoses method of groupedData object grpData to construct doses.

functionName — Estimation function namestring

Estimation function name, specified as a string. Choices are 'nlmefit' or 'nlmefitsa'.

opt — Options specific to estimation functionstruct

Options specific to the estimation function, specified as a structure. The structure contains fields and default values that are the name-value pair arguments accepted by nlmefit and nlmefitsa, except the following that are not supported.

  • 'FEConstDesign'

  • 'FEGroupDesign'

  • 'FEObsDesign'

  • 'FEParamsSelect'

  • 'ParamTransform'

  • 'REConstDesign'

  • 'REGroupDesign'

  • 'REObsDesign'

  • 'Vectorization'

'REParamsSelect' is only supported when covEstiminfo is a vector of estimatedInfo objects.

Use the statset function only to set the 'Options' field of the opt structure as follows.

opt.Options = statset('Display','iter','TolX',1e-3,'TolFun',1e-3);

For other supported name-value pair arguments (see nlmefit and nlmefitsa), set them as follows.

opt.ErrorModel = 'proportional';
opt.ApproximationType = 'LME';

variants — Variants[] | vector of variant objects

Variants, specified as an empty array or vector of variant objects. If empty, no variants are used even if the model has active variants.

tf'UseParallel' optiontrue | false

'UseParallel' option, specified as true or false. If true, and Parallel Computing Toolbox is available, the function performs parameter estimation in parallel.

Output Arguments

expand all

fitResults — Estimation resultsNLMEResults object

Estimation results, returned as an NLMEResults object.

simDataI — Simulation resultsvector of SimData objects

Simulation results, returned as a vector of SimData objects representing simulation results for each group (or individual) using fixed-effect and random-effect estimates (individual-specific parameter estimates).

The states reported in simDataI are the states that were included in the responseMap input argument as well as any other states listed in the StatesToLog property of the runtime options (RuntimeOptions) of the SimBiology model sm.

simDataP — Simulation resultsvector of SimData objects

Simulation results, returned as a vector of SimData objects representing simulation results for each group (or individual) using only fixed-effect estimates (population parameter estimates).

The states reported in simDataP are the states that were included in the responseMap input argument as well as any other states listed in the StatesToLog property of the runtime options (RuntimeOptions) of the SimBiology model sm.

References

[1] Grasela Jr, T.H., Donn, S.M. (1985) Neonatal population pharmacokinetics of phenobarbital derived from routine clinical data. Dev Pharmacol Ther. 8(6), 374–83.

See Also

| | | | | | |

Was this topic helpful?