SimBiology

フォロー

I have a reaction network model built using simbiology. This includes reactions for synthesis and degradation of various species. So, generally I run steadystate and update the model state before running any simulations. Now, I would like to estimate some of the parameters in the model using sbiofit. My questions is whether sbiofit equilibrates model state for each parameter sample during fitting? Thanks and regards.

I want to estimate the parameters of my own model. The model contains several reactions, most of which obey Uni-Uni Michaelis-Menten kinetic law. There are also additional two reactions that are actually catalyzed by the same enzyme, and thus must have the same parameter values. % Add custom-made kinetic law ratelawexpression1 = 'kcat*e0*(S1*S2/Km1*Km2)/((1+S1/Km1)*(1+S2/Km2))'; mm_1 = sbioabstractkineticlaw('CustomMadeLaw',ratelawexpression1); set(mm_1,'SpeciesVariables',{'S1','S2'}); set(mm_1,'ParameterVariables',{'kcat','e0','Km1','Km2'}); sbioaddtolibrary(mm_1); % Add second reaction to the model r2 = addreaction(model,'B1 + C -> D'); k2 = addkineticlaw(r2,'CustomMadeLaw'); p2a = addparameter(k2,'kcat',110); p2b = addparameter(k2,'e0',0.6); p2c = addparameter(k2,'Km_B1',0.01); p2d = addparameter(k2,'Km_C',0.01); set(k2,'ParameterVariableNames',{'kcat','e0','Km_B1','Km_C'}); set(k2,'SpeciesVariableNames',{'B1','C'}); set(r2,'Name','r2'); % Add third reaction to the model r3 = addreaction(model,'B2 + C -> F'); k3 = addkineticlaw(r3,'CustomMadeLaw'); p3a = addparameter(k3,'kcat',110); p3b = addparameter(k3,'e0',0.6); p3c = addparameter(k3,'Km_B2',0.01); p3d = addparameter(k3,'Km_C',0.01); set(k3,'ParameterVariableNames',{'kcat','e0','Km_B2','Km_C'}); set(k3,'SpeciesVariableNames',{'B2','C'}); set(r3,'Name','r3'); From what I understand when using sbiofit, when performing parameter estimation, the kcat, e0 and Km_C parameters of both r2 and r3 are considered different. But I would like to configure the parameter estimation settings so that those parameters are considered the same (r2.kcat should be considered the same parameter as r3.kcat, r2.e0 is r3.e0, r2.Km_C is r3.Km_C, etc.). How should I write the model or how do I change the settings of sbiofit so that I can do the previously-mentioned task?

Hello everyone, I'm doing fits in Simbiology in order to estimate parameters (Isqnonlin, constant error). The problem is that I got unusual result once. I checked everything (iniatial values, programm and so on) and tried again. It hasn't changed. So I tried also on my friend's PC. The Matlab version and the project are the same. The results on my friend's PC are different and look better. What is the reason for such outcomes? What should I do to get the same results?

Hello, I have a working PK/PD model implemented in SimBiology that relates the kinetics of a drug to it's effect on the growth of a tumor. I'm looking to perform the following task. The goal is to compute the dose at which a minimal efficacy signal is attained. Specifically, I'm looking to estimate the dose at which, four weeks post drug application, the tumor volume is 80% of the volume that is attained when the tumor grows without any applied drug. In other words, the dose at which there is a 20% reduction in the maximum volume attained after four weeks. I was wondering if there is a way to do this using SimBiology? Many thanks for your time in advance.

Suppose I have a therapeutic d that modulates an endogenous substance e. I've a pd model of e (for the sake of argument, a constant synthesis rate, and a linear clearance). I've a pk model for d. I have some sort of model that has e interact and bind with d, and the bound complex has a different clearance rate than e alone. My dataset has therapeutic doses applied to subjects at time zero. But I have undosed data for t<0 - perhaps a dataset of a couple months for e leading into the dosing period. The earlier data helps me set the parameters pertaining to e. But the fitting programs won't accept any data earlier than t=0. Is this correct? Is my workaround to bias all times so that all time values are positive, and to tell SB that the dose is delayed by this bias? Or can I tell SB to accept values of time less than zero?

We are often asked to help with parameter estimation problems. This discussion aims to provide guidance for parameter estimation, as well as troubleshooting some of the more common failure modes. We welcome your thoughts and advice in the comments below. *Guidance and Best Practices:* 1. *Make sure your data is formatted correctly.* Your data should have: * a time column (defined as independent variable) that is monotonically increasing within every grouping variable, * one or more concentration columns (dependent variable), * one or more dose columns (with associated rate, if applicable) if you want your model to be perturbed by doses, * a column with a grouping variable is optional. *Note:* the dose column should only have entries at time points where a dose is administered. At time points where the dose is not administered, there should be no entry. When importing your data, MATLAB/SimBiology will replace empty cells with NaNs. Similarly, the concentration column should only have entries where measurements have been acquired and should be left empty otherwise. <</matlabcentral/discussions/uploaded_files/54/Table_example.PNG>> If you import your data first in MATLAB, you can manipulate your data into the right format using datatype ' <https://www.mathworks.com/help/matlab/tables.html table> ' and its methods such as <https://www.mathworks.com/help/matlab/ref/sortrows.html sortrows> , <https://www.mathworks.com/help/matlab/ref/table.join.html join> , <https://www.mathworks.com/help/matlab/ref/innerjoin.html innerjoin> , <https://www.mathworks.com/help/matlab/ref/outerjoin.html outerjoin> , <https://www.mathworks.com/help/matlab/ref/stack.html stack> and <https://www.mathworks.com/help/matlab/ref/unstack.html unstack> . You can then add the data to SimBiology by using the 'Import data from MATLAB workspace' functionality. 2. *Visually inspect data and model response.* Create a simulation task in the SimBiology desktop where you plot your data ( <https://www.mathworks.com/help/simbio/ug/data.html#bveve2m%20 plot external data> ), together with your model response. You can create sliders for the parameters you are trying to estimate (or use <https://www.mathworks.com/help/simbio/ug/configuring-tasks.html#bu80xfr group simulation>). You can then see whether, by varying these parameter values, you can bring the model response in line with your data, while at the same time giving you good initial estimates for those parameters. This plot can also indicate whether units might cause a discrepancy between your simulations and data, and/or whether doses administered to the model are configured correctly and result in a model response. 3. *Determine sensitivity of your model response to model parameters.* The previous section can be considered a manual sensitivity analysis. There is also a more systematic way of performing such an analysis: a global or <https://www.mathworks.com/help/simbio/ug/calculating-sensitivities.html local sensitivity analysis> can be used to determine how sensitive your responses are to the parameters you are trying to estimate. If a model is not sensitive to a parameter, the parameter’s value may change significantly but this does not lead to a significant change in the model response. As a result, the value of the objective function is not sensitive to changes in that parameter value, hindering estimating the parameter’s value effectively. 4. *Choose an optimization algorithm.* SimBiology supports a <https://www.mathworks.com/help/simbio/ug/supported-methods-for-parameter-estimation.html range of optimization algorithms>, depending on the toolboxes you have installed. As a default, we would recommend using lsqnonlin if you have access to the <https://www.mathworks.com/help/optim/index.html Optimization Toolbox>. See troubleshooting below for more considerations choosing an appropriate optimization algorithm. 5. *Map your data to your model components:* Make sure the columns for your dependent variable(s) and dose(s) are <https://www.mathworks.com/help/simbio/ref/sbiofit.html#bt_0umy-1-responseMap mapped> to the corresponding component(s) in your model. 6. *Start small:* bring the estimation task down to the smallest meaningful objective. If you want to estimate 10 parameters, try to start with estimating one or two instead. This will make troubleshooting easier. Once your estimation is set up properly with a few parameters, you can increase the number of parameters. *Troubleshooting* 1. *Are you trying to estimate a parameter that is governed by a rule?* You can’t estimate parameters that are the subject of a rule (initial/repeated assignment, algebraic rule, rate rule), as the rule would supersede the value of the parameter you are trying to estimate. See <https://www.mathworks.com/matlabcentral/answers/428147-infinite-standarderror-with-sbiofit this> topic. 2. *Is the optimization using the correct initial conditions and parameter values?* Check whether - for the fit task - the parameter values and initial conditions that are used for the model, make sense. You can do this by passing the relevant dose(s) and variant(s) to the <https://www.mathworks.com/help/simbio/ref/model.getequations.html getequations> function. In the SimBiology App, you can look at your equation view (When you have your model open, in the Model Tab, click Open -> Equations). Subsequently, - in the Model tab - click "Show Tasks" and select your fit task and inspect the initial conditions for your parameters and species. A typical example of this is when you do a dosing species but ka (the absorption rate) is set to zero. In that case, your dose will not transfer into the model and you will not see a model response. 3. *Are your units consistent between your data and your model?* You can use <https://www.mathworks.com/help/simbio/ref/unitconversion.html unit conversion> to automatically achieve this. 4. *Have you checked your solver tolerances?* The <https://www.mathworks.com/help/simulink/ug/types-of-solvers.html#f11-44943 absolute and relative tolerance> of your solver determine how accurate your model simulation is. If a state in your model is on the order of 1e-9 but your tolerances only allow you to calculate this state with an accuracy down to 1e-8, your state will practically represent a random error around 1e-8. _This is especially relevant if your data is on an order that is lower than your solver tolerances._ In that case, your objective function will only pick up the solver error, rather than the true model response and will not be able to effectively estimate parameters. When you plot your data and model response together and by using a log-scale on the y-axis (right-click on your Live plot, select Properties, select Axes Properties, select Log scale under “Y-axis”) you can also see whether your ODE solver tolerances are sufficiently small to accurately compute model responses at the order of magnitude of your data. A give-away that this is not the case is when your model response appears to randomly vary as it bottoms-out around absolute solver tolerance. <</matlabcentral/discussions/uploaded_files/56/poor_tolerances.png>> _ Tolerances are too low to simulate at the order of magnitude of the data. Absolute Tolerance: 0.001, Relative Tolerance: 0.01_ <</matlabcentral/discussions/uploaded_files/57/good_tolerances.png>> _ Sufficiently high tolerances. Absolute Tolerance: 1e-8, Relative Tolerance: 1e-5_ 5. *Have you checked the tolerances and stopping criteria of your optimization algorithm?* The goal for your optimization should be that it terminates because it meets the imposed tolerances rather than because it exceeds the maximum number of iterations. Optimization algorithms terminate the estimation based on <https://www.mathworks.com/help/optim/ug/tolerances-and-stopping-criteria.html tolerances and stopping criteria>. An example of a tolerance here is that you specify the precision with which you want to estimate a certain parameter, e.g. Cl with a precision down to 0.1 ml/hour. If these tolerances and stopping criteria are not set properly, your optimization could terminate early (leading to loss of precision in the estimation) or late (leading to unnecessarily long optimization compute times). 6. *Have you considered structural and practical identifiability of your parameters?* In your model, there might exist values for two (or more) parameters that result in a very similar model response. When estimating these parameters, the objective function will be very similar for these two parameters, resulting in the optimization algorithm not being able to find a unique set of parameter estimates. This effect is sometimes called aliasing and is a structural identifiability problem. An example would be if you have parallel enzymatic (Km, Vm) and linear clearance (Cl) routes. Practical identifiability occurs when there is not enough data available to sufficiently constrain the parameters you are estimating. An example is estimating the intercompartmental clearance (Q12), when you only have data on the central compartment of a two-compartment model. Another example would be that your data does not capture the process you are trying to estimate, e.g. you don’t have data on the absorption phase but are trying to estimate the absorption constant (Ka). 7. *Have you considered trying another optimization algorithm?* SimBiology <https://www.mathworks.com/help/simbio/ug/supported-methods-for-parameter-estimation.html supports a range of optimization algorithms>, depending on the toolboxes you have installed. There is no single answer as to which algorithm you should use but some general guidelines can help in selecting the best algorithm. * *Non-linear regression:* If your aim is to estimate parameters estimates for each group in your dataset (unpooled) or for all groups (pooled), you can use non-linear regression estimation methods. The optimization algorithms can be broken down into local and global optimization algorithms. You can use a _local optimization algorithm_ when you have good initial estimates for the parameters you are trying to estimate. Each of the local optimization functions has a different default optimization algorithm: <https://www.mathworks.com/help/matlab/ref/fminsearch.html fminsearch> (Nelder-Mead/downhill simplex search method), <https://www.mathworks.com/help/optim/ug/fmincon.html fmincon> (interior-point), <https://www.mathworks.com/help/optim/ug/fminunc.html fminunc> (quasi-newton), <https://www.mathworks.com/help/stats/nlinfit.html nlinfit> (Levenberg-Marquardt), <https://www.mathworks.com/help/optim/ug/lsqcurvefit.html lsqcurvefit>, <https://www.mathworks.com/help/optim/ug/lsqnonlin.html lsqnonlin> (both trust-region-reflective algorithm). As a default, we would recommend using lsqnonlin if you have access to the <https://www.mathworks.com/help/optim/index.html Optimization Toolbox>. Note that all but the fminsearch algorithm are gradient based. If a gradient based algorithm fails to find suitable estimates, you can try fminsearch and see whether that improves the optimization. All local optimization algorithms can get “stuck” in a local minimum of the objective function and might therefore fail to reach the true minimum. _Global optimization algorithms_ are developed to find the absolute minimum of the objective function. You can use global optimization algorithms when your fitting task results in different parameter estimates when repeated with different initial values (in other words, your optimization is getting stuck in local minima). You are more likely to encounter this as you increase the number of parameters you are estimating, as you increase the parameter space you are exploring (in other words, the bounds you are imposing on your estimates) and when you have poor initial estimates (in other words, your initial estimates are potentially very far from the estimates that correspond with the minimum of the objective function). A disadvantage of global optimization algorithms is that these algorithms are much more computationally expensive – they often take significantly more time to converge than the local optimization methods do. When using global optimization methods, we recommend using SimBiology’s built-in <https://www.mathworks.com/help/simbio/ref/sbiofit.html#bvlnfq8 scattersearch> algorithm, combined with lsqnonlin as a local solver. If you have access to the <https://www.mathworks.com/help/gads/index.html Global Optimization Toolbox>, you can try the functions <https://www.mathworks.com/help/gads/ga.html ga> (genetic algorithm), <https://www.mathworks.com/help/gads/patternsearch.html patternsearch> and <https://www.mathworks.com/help/gads/particleswarm.html particleswarm>. Note that some of the global optimization algorithms, including scattersearch, lend themselves well to be accelerated using parallel or distributed computing. * *Estimate category-specific parameters:* If you want to estimate category-specific parameters for multiple subjects, e.g. you have 10 male and 10 female subjects in your dataset and you want to estimate a separate clearance value for each gender while all other parameters will be gender-independent, you can also use non-linear regression. Please refer to <https://www.mathworks.com/help/simbio/ug/estimate-category-specific-PK-parameters-in-a-hierarchical-model.html this> example in the documentation. * *Non-linear mixed effects*: If your data represents a population of individuals where you think there could be significant inter-individual variability you can use <https://www.mathworks.com/help/simbio/ug/what-is-nonlinear-mixed-effects-modeling.html mixed effects modeling> to estimate the fixed and random effects present in your population, while also understanding covariance between different parameters you are trying to estimate. When performing mixed effects estimation, it is advisable to perform fixed effects estimation in order to obtain reasonable initial estimates for the mixed effects estimation. SimBiology supports two estimation functions: <https://www.mathworks.com/help/stats/nlmefit.html nlmefit> (LME, RELME, FO or FOCE algorithms), and <https://www.mathworks.com/help/stats/nlmefitsa.html nlmefitsa> (Stochastic Approximation Expectation-Maximization). Sometimes, these solvers might seem to struggle to converge. In that case, it is worthwhile determining whether your (objective) function tolerance is set too low and increasing the tolerance somewhat. 8. *Does your optimization get stuck?* Sometimes, the optimization algorithm can get stuck at a certain iteration. For a particular iteration, the parameter values that model is simulated with as part of the optimization process, can cause the model to be in a state where the ODE solver needs to take very small time-steps to achieve the tolerances (e.g. very rapid changes of model responses). Solutions can include: changing your initial estimates, imposing lower and upper bounds on the parameters you are trying to estimate, selecting to a different solver, easing solver tolerances (only where possible, see also “Visually inspect data and model response”). 9. *Are you using the proportional error model?* The objective function for the proportional <https://www.mathworks.com/help/simbio/ref/sbiofit.html#bual8r7 error model> contains a term where your response data is part of the denominator. As response variables get close to zero or are exactly zero, this effectively means the objective function contains one or more terms that divide by zero, causing errors or at least very slow iterations of your optimization algorithm. You can try to change the error model to constant or combined to circumvent this problem. Alternatively, you can define separate error models for each response: proportional for those responses that don’t have measurements that contain values close to zero and a constant error model for those responses that do. Guidance and troubleshooting for parameter estimation in SimBiology Excellent points and a very well-written article which I am sure many would find useful! simbiology fitting parameter estimation error non-linear regression mixed effects optimization

Hello, I have a model built in simbiology that I have confirmed to be stable and running with a set of nominal parameters (10+ parameters). Now I am trying to fit one of these parameters from the model with *sbiofit* using the instructions <https://www.mathworks.com/help/simbio/ref/sbiofit.html#bujy6r8 here>, but am running into some issues. Unfortunately I cannot share the model publicly here at the moment, which I realize makes it a bit difficult to help, but I thought I'd ask anyways to see if I'm doing anything that pops out as immediately incorrect to someone. The general structure of my code is as follows: T = table(time,conc); % time and conc are equally sized vectors of time-concentration data points gData = groupedData(T); responseMap = {'compartment.species= conc'}; paramsToEstimate = {'log(p1)'}; % keep p1 positive p1_i = 0.003; p1_limits = [0.0001 3]; estimatedParams = estimatedInfo(paramsToEstimate,'InitialValue',[p1_i],'Bounds',[p1_limits]); fitConst = sbiofit(model,gData,responseMap,estimatedParams,dose); s.Labels.XLabel = 'Time (hour)'; s.Labels.YLabel = 'Amount (nanomole)'; plot(fitConst,'AxesStyle',s); This code runs successfully, but the estimate never deviates from my initial value guess specified for p1 (0.003) and the StandardError that gets spit out is consistently infinity ('Inf'). I realize this may be difficult to debug without my exact model posted here, but I was wondering if anyone has any pointers as to what might be causing the fitter to behave like this? I've searched around and haven't been able to find anyone who has had this problem on these forums or the rest of the web as far as I can tell. Thanks in advance!

In SimBiology I have constructed a simple 2 Compartment Model to describe drug distribution over a set time course. The Model uses 6 parameters. I have a set of 11 individual data sets for the drug, however they are taken from different disease stages (early & late) which can be quantified as the "Time stage of Disease". Previous literature has stated a relationship between the stage of the disease and drug accumulation. However the specific relationship is not known. How do I incorporate the stage of disease as a covariate into my model? I'm unclear how to specify a relationship between the stage and the parameter, most literature says to use a NLME model but again I'm unclear about the implementation of this

I have seen some similar questions about this but they haven't quite answered what I need to do. For my set-up I have a simple compartmental model which describes drug distribution, the volumes of these compartments is dependent on a parameter ( _Bodyweight_) I am trying to parameterise the model using several datasets which are pooled to provide a single parameter set However each dataset comes from a different experiment where the bodyweights are different and so I cannot use a single set bodyweight for all as it bias' my parameter estimations Does anyone know how I could set a parameter estimation so that for each dataset a different value for the parameter _Bodyweight_ is used?

Hello Matlab community, I have created a schedule dosage to be supplied to central compartment in the model. I have successfully used it in the simulation task. However, I am having trouble figuring out how to use that dosage for "fit parameter" task. The dataset I wish to fit has a independent variable (time) and a dependent variable (concentration). I have tried creating a separate column in dataset and provide dosage information there. However, for the time points dosages are given, I do not have the measured concentration at the same time points. This creates a missing data in the dataset for the concentration. Can experts please help me out with the situation. thanks, appreciate that. Best chhavi

I have a multi-compartmental PK model that I am trying to parameterise utilising two different data-sets. One data-set is a short term data set (0 -> 24hrs) and the second data set is a single data point at 168 hrs. The only difference between the two is that the initial bolus dosage is a different amount (and so I cannot combine the two sets) How can I get SimBiology to run a parameter estimation where the output parameters are based off of both data sets? I am running MATLAB R2016b on Mac & Version 5.5 of SimBiology

HI there, Q1: i'm wondering if there are two distinct ODE models and each representing a biological system that may have some relationship together (for example cell cycle pathway and apoptosis pathway), is it possible to mix them in one model using the parameters and initial values same as the original model (without changing the parameters)? Because there are some mutual elements in these models, the modelers made some assumption to get a correct dynamics. So, when you combine the two model, the dynamic of a particular element (X) will be different (incorrect). Q2: My question is that since the dynamics of element X in each of these single models is correct, is it possible to give the correct time course (dynamic) of X to SIMBIOLOGY in order to have the correct dynamic at last(and Simbiology probably fine-tunes some other parameters)? Q3: Is it possible that this process makes the dynamics of other species in the model incorrect? Q4: And another question is that if it is possible to fix some parameters (known parameters) and want Simbiology to just estimate the unknown parameters (and we may have the dynamics of some of the species but not all)? Q5: I remember that I tried to give one time course (of one specie) to simbilogy but the answer was quite wrong! DO you know why it may happens?