SimBiology parameter estimation for a series of single dose-response data

1 view (last 30 days)
I would like to estimate set of parameters from a series of experiments with different doses and single measurement points. I.e. there is one concentration measurement for each dose, all administered at time=0.
(It is similar to the a recent discussion but I could not find solution to my problem.)
I get reasonably looking parameter estimates but when reentering the doses in the model the reults don't quite match. I attach the working sbproj and script in PS. Would appreciate some comments what is missing.
Best, M
PS:
So far my code reads:
%% multiDOSE.m
clear all
DOSE_TIME = [0,0,0,0,0,0]'; % day
DOSE_AMT = [1,2,3,4,5,6]'; % mg
DOSEdata = groupedData(table(DOSE_TIME, DOSE_AMT));
DOSEdata.Properties.IndependentVariableName = 'DOSE_TIME';
DOSEdata.Properties.VariableUnits = {'day','milligram'};
DOSEdata.Properties
sampleDose = sbiodose('sample', 'TargetName', 'CELL.DRUG');
doseArray = createDoses(DOSEdata,'DOSE_AMT','',sampleDose);
TIME = [10,10,10,10,10,10]'; % day
CONC = [14,44,88,105,192,254]'; % ng/mL
CONCdata = groupedData(table(TIME, CONC));
CONCdata.Properties.VariableUnits = {'day','milligram/milliliter'};
s = sbioloadproject("multiDOSE.sbproj");
modelObj = s.m1;
paramsToEstimate = {'log(kdeg)','log(kin)'};
paramToEst = estimatedInfo(paramsToEstimate,'InitialValue',[1,1]);
opt = optimset('PlotFcns',@optimplotfval,'MaxIter',150);
result1 = sbiofit(modelObj, CONCdata, 'Species1 = CONC', paramToEst, ...
doseArray, 'fminsearch', opt,'Pooled', true);
estValues1 = result1.ParameterEstimates

Accepted Answer

Arthur Goldsipe
Arthur Goldsipe on 3 Aug 2022
I think the problem is that you need to add a grouping variable to your DOSEdata and your CONCdata. As the code is currently written, you are doing a single simulation with 6 different doses at time 0, for a cumultive dose of sum(1:6), or 21 milligram. And then you're simultaneously trying to fit all 6 measurements at time 10 to the result of a single simulation at time 10. (We allow different repeated measurements at the same time in order to support technical or biological replicates.)
I'll attach my updated version of your script. But the main changes are that (1) I added the following group identifier line near the top:
ID = (1:6)';
and (2) I updated the groupedData setup as follows:
DOSEdata = groupedData(table(ID, DOSE_TIME, DOSE_AMT));
DOSEdata.Properties.VariableUnits = {'', 'day','milligram'};
CONCdata = groupedData(table(ID, TIME, CONC));
CONCdata.Properties.VariableUnits = {'', 'day','milligram/milliliter'};
If you like, you could create a single groupedData that contains both your dosing information and the measurements. You just need to use NaN values as appropriate to indicate no dose or no measurement at the specific time. For example:
ID = [1 1 2 2]';
TIME = [0 10 0 10]';
DOSE = [1 nan 2 nan]';
CONC = [nan 14 nan 44]';
data = groupedData(table(ID, TIME, DOSE, CONC))
data = 4×4 groupedData ID TIME DOSE CONC __ ____ ____ ____ 1 0 1 NaN 1 10 NaN 14 2 0 2 NaN 2 10 NaN 44
  3 Comments
emjey
emjey on 5 Aug 2022
Thanks Arthur, it works!
The second solution looks very efficient. Does that mean the sbiofit would read like this:
result1 = sbiofit(modelObj, data, 'Species1 = CONC', paramToEst, ...
data, 'fminsearch', opt,'Pooled', true);
with 'data' replacing 'doseArray' and 'CONCdata' with NONMEM-like keywords 'DOSE' and 'CONC'? Do you need to specify the specify the Properties as well?
Best, Maciej
Arthur Goldsipe
Arthur Goldsipe on 5 Aug 2022
No, you still need to pass doseArray in to sbiofit. But you would create it from the same data you pass in to sbiofit. So, the code might look something like this:
data = groupedData(table(ID, TIME, DOSE, CONC));
doseArray = createDoses(data, 'DOSE_AMT', '', sampleDose)
result1 = sbiofit(modelObj, data, 'Species1 = CONC', paramToEst, ...
doseArray, 'fminsearch', opt, 'Poosed', true);
The main reason I designed things this way was for additional flexibility. Your dosing information can be embedded in the same data set as your experimental measurements, but it doesn't have to be.

Sign in to comment.

More Answers (0)

Products


Release

R2022a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!