SimBiology

Segui

Hi, I wrote the code below for estimation of KD from the binding assay for three different cell lines: i.e. bxpc3, colo205, T84 which have different number of antigen molecule on their surface. The code properly work for estimation of kD spearately for each cell line. How can I estimate the kD as one value for different cell lines at the same time? I attached my excel file (exp.xlsx) in which each cell line data has been stored in a spearted sheet. Thanks %% Setup Model case_study = 'T84'; % 'bxpc3' 'colo205' 'T84' model_input_file='exp.xlsx' switch case_study case 'bxpc3' exp_data=readtable(model_input_file,'Sheet','bxpc3'); data.dose = exp_data.dose; data.ro = exp_data.RO_exp; input1=[0.1 0.4 1.2 11.1 33.3 100 300]; RPC_Ag_Tumor = 219202; case 'colo205' exp_data=readtable(model_input_file,'Sheet','colo205'); data.dose = exp_data.dose; data.ro = exp_data.RO_exp; input1=[0.1 0.4 1.2 3.7 11.1 33.3 100]; RPC_Ag_Tumor = 662642; case 'T84' exp_data=readtable(model_input_file,'Sheet','T84'); data.dose = exp_data.dose; data.ro = exp_data.RO_exp; input1=[0.1 0.4 1.2 3.7 11.1 33.3 100 300]; RPC_Ag_Tumor = 169163; end data = structfun( @rmmissing , data , 'UniformOutput' , false); setup; % (I build the model here) %% Data Import gData = groupedData(exp_data); gData.Properties.VariableUnits = {'','second','nanomole','dimensionless'}; gData.Properties.GroupVariableName = 'ID'; gData.Properties.IndependentVariableName = 'Time'; sbiotrellis(gData,'ID','Time',{'RO_exp'},'Marker','+','LineStyle','none'); % Extract data columns into separate variables for the weight expression % evaluation. headings = exp_data.Properties.VariableNames; for i=1:length(headings) next = headings{i}; eval([next ' = exp_data. ' next ';']); end %% Create the data dose. d2 = sbiodose('dose'); d2.TargetName = 'Ab'; d2.LagParameterName = ''; d2 = createDoses(gData, 'dose', '', d2); % Build table of doses. d2 = num2cell(d2); dosesForFit = d2; %% Fitting options % Define response information. responseMap = {'RO = RO_exp'}; % Define objects being estimated and their initial estimates. estimatedInfoObj = estimatedInfo({'kD_Ab'}); % estimatedInfoObj = estimatedInfo({'k_max', 'EC50'}); estimatedInfoObj(1).Bounds = [1E-1 1E2]; % estimatedInfoObj(2).Bounds = [1E-6 1E-4]; % Define Algorithm options. options = optimoptions('particleswarm'); options.Display= 'iter'; % Define fit problem. f = fitproblem('FitFunction', 'sbiofit'); f.Data = gData; f.Model = model; f.Estimated = estimatedInfoObj; f.ErrorModel = 'exponential'; f.ResponseMap = responseMap; f.Weights = []; f.Doses = dosesForFit; f.FunctionName = 'particleswarm'; f.Options = options; f.ProgressPlot = true; f.UseParallel = false; f.Pooled = true; %% Run the model with optimum parameters fitResults = f.fit; plot(fitResults) ci = sbiopredictionci(fitResults); plot(ci) ciParam = sbioparameterci(fitResults); ci2table(ciParam) plot(ciParam) fitResults.ParameterEstimates %update parameters for i=1:length(estimatedInfoObj) xopt = sbioselect(model,'Name',fitResults.EstimatedParameterNames{i}); xopt.value = fitResults.ParameterEstimates.Estimate(i); end params = {'Ab0'}; obs = {'RO'}; input=dose(~isnan(dose)); RO_exp = RO_exp(~isnan(RO_exp)); sfxn = createSimFunction(model, params, obs, d1,'UseParallel',false,'AutoAccelerate',false); sfxn.accelerate; doseTable = getTable(d1); simData=sfxn(input,configsetObj.StopTime,doseTable); %% Plot fontsize = 18; linesize = 2; markersize = 10; f = figure('Position', [14 10 900 600]); colors = ["k", "#D95319", "#A2142F", "#77AC30", "#7E2F8E", "#00CED1", "#FF1493"]; % Define an array of colors set(f, 'Visible', 'on'); set(gca,'fontsize', fontsize); hold on; box on; grid on; set(gca, 'XScale', 'log'); %set(gca, 'YScale', 'log'); for i = 1:size(input,1) if i < size(input,1) line([input1(i) input1(i+1)], [simData(i).Data(end) simData(i+1).Data(end)], 'Color', colors(3),'LineWidth', 2); end % plot(input1(i), simData(i).Data(end), 'O', 'MarkerEdgeColor', colors(3), 'MarkerFaceColor', colors(3), 'markersize', markersize); plot(input1(i), RO_exp(i), 'O', 'MarkerEdgeColor', colors(4), 'MarkerFaceColor', colors(4), 'markersize', markersize); end xlabel('Dose (ug/mL)', 'fontweight', 'bold', 'Fontsize', fontsize, 'Interpreter', 'none'); ylabel('RO (%)', 'fontweight', 'bold', 'Fontsize', fontsize); legend('Simulation','Experiment' , 'NumColumns', 1, 'Location', 'best', 'FontSize', fontsize) saveas(gca, fullfile(save_locations{2}, 'RO'), 'png'); close(f);

Hello everyone, I am trying to fit a model of mine on a remote linux machine (Debian GNU/Linux 11 (bullseye)). I'd like to use MATLAB as a command line application and not use a X11 interface. I have not found any tutorial or suggestion in the documentation, as all examples either use the Simbiology app or the 'classical' Desktop version of MATLAB. I have already installed MATLAB on the server with all necessary add-ons, does someone have suggetions on how to proceed? Also considering that I am working on my local pc for the development nof the model and I need to 'transfer' my data to the MATLAB worksopace on the server.

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?

Hi all, With sbiofit (SimBiology) the fit results give me standard errors obtained from the Jacobian calculation. I would like to compute the standard errors with the Hessiana matrix, but I don't know how to get it analytically from my object function and the fit parameters. Thank you in advance, Valentina

Hello everyone, I am trying to plot the results of the confidence intervals for the model predictions. Since I have 17 ResponseNames, my figure is overloaded with images (see attached picture). How do I separate the various plots from the figure? I am using this function: sbioparameterci ciPredPooled = sbiopredictionci(pooledFit); plot(ciPredPooled); Thank you in advance for your help. Best regards, Valentina

Hi, I am currently working with SimBiology and have implemented a compartmental model to follow the Lutatera therapy. I am using the 2021 version of Matlab. I have modified sbiofit to implement the object function in order to perform a fit with Bayes values. Comparing it with another software that already uses Bayes, I could see that the fitted parameters have values that differ little. The problem is that with matlab, the standard errors that sbiofit calculates are too large. At first I thought it was just because of the degrees of freedom, but unfortunately, even when correcting them, I still get too high errors on the fitted parameters. Upon redoing the code, I noticed that the Jacobian does not take into account my Bayes values, nor the standard Bayes error. numjac is the function that calculates the Jacobian. My question is whether it is possible, and how, to incorporate Bayes standard errors into the standard errors of my parameters. In theory, from the fit, the errors of the fitted parameters should be no larger than the Bayes errors. Has anyone ever solved this problem in other functions used for the fit? Or are there functions that already use Bayes that can help me solve my problem? Regards

Hello team, I am using command line sbiofit to fit my simbiology model with my experimental data by parameter estimation. In addition, after I got the fitting results, I know that I could simply plot the fitted results versus the original data by plot(fitResults); However, I would like to access the simulated results of my best fit to carry out NCA, but I could not find the results in the SimBiology.fit.OptimResults. Could you please show me a clue on how to get the piece of the simulated data? Thank you very much. Have a nice day. Best, Jesse

Hello, I have a question regarding the regression tools in SimBiology. I have an ODE model which I would like to fit for a given dataset using sbiofit or sbiofitmixed. When i want to find some parameter set using sbiofitmixed, the ODE model can return with invalid values (e.g. complex, NaN or Inf) which completely halts the search process. My experience with sbiofit is that the algorithms can recover from these values (as opposed to sbiofitmixed) for the vast majority of cases, but in some occasion the procedure still halts. Even if I set bounds on my parameters, there could be a (hidden) nonlinear relationship which makes the ODE invalid for a specific combination of the fitted parameter values. The following error message can occur using sbiofitmixed: Error using sbiofitmixed (line 249) A simulation errored or returned NaN, Inf, or complex values during parameter estimation. This typically indicates that the estimated parameter values converged to values that are invalid for the model. My question is if there is any way to recover sbiofitmixed from these parameter values? For example defining an event during the ODE simulation which stops the integration. It would be nice, if there was a way to find a different trial point for the current iteration in such cases for the algorithm. Bonus question: Is there any 'best practice' in the literature to handle these situations? For example if I code the algorithm by hand, I could just stop the integration if some state variable exceeds some given value. At these instances I could just set the cost function to some high value (for example 1e5), however I have a gut feeling that this is not the best solution for this issue. Thank you in advance, Ben

Hi, In the example "Fit a One-Compartment Model to an Individual's PK Profile" for the function sbiofit (the example is found in the sbiofit help page), an initial dose of 10 milligram is used: dose = sbiodose('dose'); dose.TargetName = 'Drug_Central'; dose.StartTime = 0; dose.Amount = 10; dose.AmountUnits = 'milligram'; dose.TimeUnits = 'hour'; I just wonder where that value is coming from as it seems I cannot find it anywhere in the example. Should it be also a value that needs to be estimated, or it is a value that needs to be known from experimental/synthetic data? Many thanks, Jacopo

Hello, Is it possible to specify different variants for different groups when using the 'sbiofit' function in code? From the documentation, it appears that 'sbiofit' only supports the specification of a single set of variants to be applied to all of the groups. Thank you, Abed Alnaif

I'm fitting a model to a dataset. I have two variables I'm fitting simultaneously: the therapy, and the quantity of an endogenous protein. Call the model variables Drug and Prot. The corresponding columns in the dataset woudl be DrugObs and ProtObs. Drug and Prot are of different magnitude. I fit the model: [fitcon, simdat] = sbiofit( m1, gmidata, resmap, estpars, doses, 'UseParallel',true); If I use plot(fitcon) I get a trellis plot of data and mode for each group (that is, each subject). This is great! However, since the magnitudes are different I'd like to plot the PRED and OBS of each variable in a different plot. Also, it's useful to look at these things in both log and linear form. So is there any easy way to say plot(fitcon,'PREDVarToPlot','Drug', 'OBSVarToPlot', 'DrugObs'); % and plot(fitcon,'PREDVarToPlot','Prot', 'OBSVarToPlot', 'ProtObs'); % and also plot(fitcon,'PREDVarToPlot','Drug', 'OBSVarToPlot', 'DrugObs',YScale,'log'); % and plot(fitcon,'PREDVarToPlot','Prot', 'OBSVarToPlot', 'ProtObs',YScale,'log'); % I can do this by going into the plot but its tedious: fh = gcf; fhc = fh.Children nfhc = numel(fhc) for i = 1:nfhc test the specific figure child to see if its an axes or a legend if axes get axes children fhcc = fhc(i).Childen % These are lines fhcc(1).Visible = 'off' % This is one way. Any way to delete the line completely? fhcc(3).Visible = 'off' % fhcc([1 3]) are two of the four plotted variables corresponding to etiher drug, or protein. Could be fhcc([2 4]) elseif legend castr = fhc(i).String castr = castr([2 4]) end Another approach would be to allow the second variable to be plotted using a different (right hand) yscale. I've not tried this in script, but as I understand it I'd have to loop through current axes and create a new Child axes with YAxisLocation = 'right'. Not sure then how I would transfer one set of variables (Drug and DrugObs or Prot and ProtObs) to that new axis. Again, it would be ideal if there was an ability to say: plot(fitcon, 'LeftAxis', 'Drug','RightAxis','Prot'); Are there any easy ways to do this, or must I create a script function? Thanks!

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!