sbiofit

Выполните регрессию нелинейного метода наименьших квадратов

Statistics and Machine Learning Toolbox™, Optimization Toolbox™ и Global Optimization Toolbox рекомендуются для этой функции.

Описание

пример

fitResults = sbiofit(sm,grpData,responseMap,estiminfo) оценочные параметры модели SimBiology sm использование регрессии нелинейного метода наименьших квадратов.

grpData isa groupedData object определение данных, чтобы соответствовать. responseMap задает отображение между компонентами модели и данными об ответе в grpData. estimatedInfo EstimatedInfo object это задает предполагаемые параметры в модели sm. fitResults isa OptimResults object или NLINResults object или вектор из этих объектов.

sbiofit использует первую доступную функцию оценки среди следующего: lsqnonlin (Optimization Toolbox), nlinfit (Statistics and Machine Learning Toolbox), или fminsearch.

По умолчанию, каждая группа в grpData является подходящим отдельно, приводя к специфичным для группы оценкам параметра. Если модель содержит активные дозы и варианты, они применяются перед симуляцией.

пример

fitResults = sbiofit(sm,grpData,responseMap,estiminfo,dosing) использует информацию о дозах, заданную матрицей объектов дозы SimBiology dosing вместо того, чтобы использовать активные дозы модели sm если существует кто-либо.

пример

fitResults = sbiofit(sm,grpData,responseMap,estiminfo,dosing,functionName) использует функцию оценки, заданную functionName. Если заданная функция недоступна, предупреждение выдано, и первая доступная функция по умолчанию используется.

пример

fitResults = sbiofit(sm,grpData,responseMap,estiminfo,dosing,functionName,options) использует дополнительные опции, заданные options для функционального functionName.

пример

fitResults = sbiofit(sm,grpData,responseMap,estiminfo,dosing,functionName,options,variants) применяет различные объекты, заданные как variants вместо того, чтобы использовать любые активные варианты модели.

пример

fitResults = sbiofit(_,Name,Value) дополнительные опции использования заданы одним или несколькими Name,Value парные аргументы.

пример

[fitResults,simdata] = sbiofit(_) также возвращается, вектор из SimData возражает simdata использование любого из входных параметров в предыдущих синтаксисах.

Примечание

  • sbiofit объединяет sbionlinfit и sbioparamestim функции оценки. Использование sbiofit выполнять регрессию нелинейного метода наименьших квадратов.

  • sbiofit симулирует модель с помощью a SimFunction object, который автоматически ускоряет симуляции по умолчанию. Следовательно не необходимо запуститься sbioaccelerate прежде чем вы вызовете sbiofit.

Примеры

свернуть все

Фон

В этом примере показано, как соответствовать PK индивидуума профилируют данные к модели с одним отсеком и оценить фармакокинетические параметры.

Предположим, что вы имеете данные о концентрации плазмы препарата от индивидуума и хотите оценить объем центрального отсека и разрешения. Примите, что концентрация препарата по сравнению с профилем времени следует за моноэкспоненциальным снижением Ct=C0e-ket, где Ct концентрация препарата во время t, C0 начальная концентрация, и ke константа скорости устранения, которая зависит от разрешения и объема центрального отсека ke=Cl/V.

Синтетические данные в этом примере были сгенерированы с помощью следующей модели, параметров и дозы:

  • Модель с одним отсеком с болюсным введением и устранением первого порядка

  • Объем центрального отсека (Central) = 1,70 литра

  • Параметр разрешения (Cl_Central) = 0,55 литра/час

  • Постоянная ошибочная модель

  • Доза шарика 10 мг

Загрузите данные и визуализируйте

Данные хранятся как таблица с переменными Time и Conc это представляет ход времени плазменной концентрации индивидуума после внутривенного администрирования шарика, измеренного в 13 различных моментах времени. Переменные модули для Time и Conc час и миллиграмм/литр, соответственно.

clear all
load('data15.mat')
plot(data.Time,data.Conc,'b+')
xlabel('Time (hour)');
ylabel('Drug Concentration (milligram/liter)');

Figure contains an axes. The axes contains an object of type line.

Преобразуйте в groupedData Формат

Преобразуйте набор данных в groupedData объект, который является необходимым форматом данных для подходящего функционального sbiofit для дальнейшего использования. groupedData объект также позволяет вам установить независимую переменную и имена переменных группы (если они существуют). Установите модули Time и Conc переменные. Модули являются дополнительными и только необходимыми для UnitConversion покажите, который автоматически преобразует соответствие с физическими количествами в одну сопоставимую модульную систему.

gData = groupedData(data);
gData.Properties.VariableUnits = {'hour','milligram/liter'};
gData.Properties
ans = struct with fields:
                Description: ''
                   UserData: []
             DimensionNames: {'Row'  'Variables'}
              VariableNames: {'Time'  'Conc'}
       VariableDescriptions: {}
              VariableUnits: {'hour'  'milligram/liter'}
         VariableContinuity: []
                   RowNames: {}
           CustomProperties: [1x1 matlab.tabular.CustomProperties]
          GroupVariableName: ''
    IndependentVariableName: 'Time'

groupedData автоматически определите имя IndependentVariableName свойство к Time переменная данных.

Создайте модель с одним отсеком

Пользуйтесь встроенной библиотекой PK, чтобы создать модель с одним отсеком с болюсным введением и устранением первого порядка, где уровень устранения зависит от разрешения и объема центрального отсека. Используйте configset объект включить модульное преобразование.

pkmd                    = PKModelDesign;
pkc1                    = addCompartment(pkmd,'Central');
pkc1.DosingType         = 'Bolus';
pkc1.EliminationType    = 'linear-clearance';
pkc1.HasResponseVariable = true;
model                   = construct(pkmd);
configset               = getconfigset(model);
configset.CompileOptions.UnitConversion = true;

Для получения дополнительной информации при создании разделенных на отсеки моделей PK с помощью SimBiology® встроенная библиотека, смотрите, Создают Фармакокинетические Модели.

Задайте дозирование

Задайте одну дозу шарика 10 миллиграммов, данных во время = 0. Для получения дополнительной информации при подготовке различных расписаний дозирования, смотрите Дозы в Моделях SimBiology.

dose                = sbiodose('dose');
dose.TargetName     = 'Drug_Central';
dose.StartTime      = 0;
dose.Amount         = 10;
dose.AmountUnits    = 'milligram';
dose.TimeUnits      = 'hour';

Сопоставьте данные об ответе с соответствующим компонентом модели

Данные содержат данные о концентрации препарата, хранимые в Conc переменная. Эти данные соответствуют Drug_Central разновидности в модели. Поэтому сопоставьте данные с Drug_Central можно следующим образом.

responseMap = {'Drug_Central = Conc'};

Задайте параметры, чтобы оценить

Параметры, чтобы поместиться в эту модель являются объемом центрального (Центрального) отсека и уровень раскрываемости преступлений (Cl_Central). В этом случае задайте логарифмическое преобразование для этих биологических параметров, поскольку они ограничиваются быть положительными. estimatedInfo объект позволяет вам указать, что параметр преобразовывает, начальные значения и границы параметра в случае необходимости.

paramsToEstimate    = {'log(Central)','log(Cl_Central)'};
estimatedParams     = estimatedInfo(paramsToEstimate,'InitialValue',[1 1],'Bounds',[1 5;0.5 2]);

Оцените параметры

Теперь, когда вы задали модель с одним отсеком, данные, чтобы соответствовать, сопоставленные данные об ответе, параметры, чтобы оценить, и дозирование, используют sbiofit оценить параметры. Функция оценки по умолчанию, что sbiofit использование изменится, в зависимости от которого тулбоксы доступны. Чтобы видеть, какая функция использовалась во время подбора кривой, проверяйте EstimationFunction свойство соответствующего объекта результатов.

fitConst = sbiofit(model,gData,responseMap,estimatedParams,dose);

Отобразите предполагаемые параметры и постройте результаты

Заметьте, что оценки параметра не были далеки от истинных значений (1.70 и 0.55), которые использовались, чтобы сгенерировать данные. Можно также попробовать различные ошибочные модели, чтобы видеть, могли ли они далее улучшить оценки параметра.

fitConst.ParameterEstimates
ans=2×4 table
         Name         Estimate    StandardError      Bounds  
    ______________    ________    _____________    __________

    {'Central'   }     1.6993       0.034821         1      5
    {'Cl_Central'}    0.53358        0.01968       0.5      2

s.Labels.XLabel     = 'Time (hour)';
s.Labels.YLabel     = 'Concentration (milligram/liter)';
plot(fitConst,'AxesStyle',s);

Figure contains an axes. The axes with title Group One group contains 2 objects of type line. These objects represent OBS1 (Conc), PRED1 (Central.Drug_Central).

Используйте различные ошибочные модели

Попробуйте три других поддерживаемых ошибочных модели (пропорциональный, комбинация постоянных и пропорциональных ошибочных моделей и экспоненциал).

fitProp = sbiofit(model,gData,responseMap,estimatedParams,dose,...
                      'ErrorModel','proportional');
fitExp  = sbiofit(model,gData,responseMap,estimatedParams,dose,...
                      'ErrorModel','exponential');
fitComb = sbiofit(model,gData,responseMap,estimatedParams,dose,...
                      'ErrorModel','combined');

Используйте веса вместо ошибочной модели

Можно задать веса как числовую матрицу, где количество столбцов соответствует количеству ответов. Установка всех весов к 1 эквивалентна постоянной ошибочной модели.

weightsNumeric = ones(size(gData.Conc));
fitWeightsNumeric = sbiofit(model,gData,responseMap,estimatedParams,dose,'Weights',weightsNumeric);

В качестве альтернативы можно использовать указатель на функцию, который принимает вектор из предсказанных значений отклика и возвращает вектор из весов. В этом примере используйте указатель на функцию, который эквивалентен пропорциональной ошибочной модели.

weightsFunction = @(y) 1./y.^2;
fitWeightsFunction = sbiofit(model,gData,responseMap,estimatedParams,dose,'Weights',weightsFunction);

Сравните информационные критерии выбора модели

Сравните логарифмическую правдоподобность, AIC и значения BIC каждой модели, чтобы видеть, какая ошибочная модель лучше всего подбирает данным. Большее значение вероятности указывает, что соответствующая модель подбирает модель лучше. Для AIC и BIC, меньшие значения лучше.

allResults = [fitConst,fitWeightsNumeric,fitWeightsFunction,fitProp,fitExp,fitComb];
errorModelNames = {'constant error model','equal weights','proportional weights', ...
                   'proportional error model','exponential error model',...
                   'combined error model'};
LogLikelihood = [allResults.LogLikelihood]';
AIC = [allResults.AIC]';
BIC = [allResults.BIC]';
t = table(LogLikelihood,AIC,BIC);
t.Properties.RowNames = errorModelNames;
t
t=6×3 table
                                LogLikelihood      AIC        BIC  
                                _____________    _______    _______

    constant error model            3.9866       -3.9732    -2.8433
    equal weights                   3.9866       -3.9732    -2.8433
    proportional weights           -3.8472        11.694     12.824
    proportional error model       -3.8257        11.651     12.781
    exponential error model         1.1984        1.6032     2.7331
    combined error model            3.9163       -3.8326    -2.7027

На основе информационных критериев постоянная ошибочная модель (или равные веса) соответствует данным лучше всего, поскольку это имеет самое большое значение логарифмической правдоподобности и самый маленький AIC и BIC.

Отобразите предполагаемые значения параметров

Покажите предполагаемые значения параметров каждой модели.

Estimated_Central       = zeros(6,1);
Estimated_Cl_Central    = zeros(6,1);
t2 = table(Estimated_Central,Estimated_Cl_Central);
t2.Properties.RowNames = errorModelNames;
for i = 1:height(t2)
    t2{i,1} = allResults(i).ParameterEstimates.Estimate(1);
    t2{i,2} = allResults(i).ParameterEstimates.Estimate(2);
end
t2
t2=6×2 table
                                Estimated_Central    Estimated_Cl_Central
                                _________________    ____________________

    constant error model             1.6993                0.53358       
    equal weights                    1.6993                0.53358       
    proportional weights             1.9045                0.51734       
    proportional error model         1.8777                0.51147       
    exponential error model          1.7872                0.51701       
    combined error model             1.7008                0.53271       

Заключение

Этот пример показал, как оценить параметры PK, а именно, объем центрального отсека и параметр разрешения индивидуума, путем подбора кривой данным о профиле PK к модели с одним отсеком. Вы сравнили информационные критерии каждой модели и оценили, что значения параметров различных ошибочных моделей видели, какая модель лучше всего объяснила данные. Финал соответствовал результатам, предложенным обоих, которых постоянные и объединенные ошибочные модели предоставили, самые близкие оценки к значениям параметров раньше генерировали данные. Однако постоянная ошибочная модель является лучшей моделью, как обозначено логарифмической правдоподобностью, AIC и критериями информации о BIC.

Предположим, что у вас есть данные о концентрации плазмы препарата от трех индивидуумов, которых вы хотите использовать, чтобы оценить соответствующие фармакокинетические параметры, а именно, объем центрального и периферийного отсека (Central, Peripheral), уровень раскрываемости преступлений (Cl_Central), и межразделенное на отсеки разрешение (Q12). Примите, что концентрация препарата по сравнению с профилем времени следует за снижением biexponential Ct=Aeat+Bebt, где Ct является концентрацией препарата во время t, и a и b являются наклонами для соответствующих экспоненциальных снижений.

Синтетический набор данных содержит данные о концентрации плазмы препарата, измеренные и в центральных и в периферийных отсеках. Данные были сгенерированы с помощью модели 2D отсека с капельным внутривенным введением и устранением первого порядка. Эти параметры использовались для каждого индивидуума.

 CentralPeripheralQ12Cl_Central
Индивидуум 11.900.680.240.57
Индивидуум 22.106.050.360.95
Индивидуум 31.704.210.460.95

Данные хранятся как таблица с переменными IDВремя, CentralConc, и PeripheralConc. Это представляет ход времени плазменных концентраций, измеренных в восьми различных моментах времени и для центральных и для периферийных отсеков после капельного внутривенного введения.

load('data10_32R.mat')

Преобразуйте набор данных в groupedData объект, который является необходимым форматом данных для подходящей функции sbiofit для дальнейшего использования. groupedData объект также позволяет вам установить независимую переменную и имена переменных группы (если они существуют). Установите модули IDВремя, CentralConc, и PeripheralConc переменные. Модули являются дополнительными и только необходимыми для UnitConversion покажите, который автоматически преобразует соответствие с физическими количествами в одну сопоставимую модульную систему.

gData = groupedData(data);
gData.Properties.VariableUnits = {'','hour','milligram/liter','milligram/liter'};
gData.Properties
ans = 

  struct with fields:

                Description: ''
                   UserData: []
             DimensionNames: {'Row'  'Variables'}
              VariableNames: {'ID'  'Time'  'CentralConc'  'PeripheralConc'}
       VariableDescriptions: {}
              VariableUnits: {1x4 cell}
         VariableContinuity: []
                   RowNames: {}
           CustomProperties: [1x1 matlab.tabular.CustomProperties]
          GroupVariableName: 'ID'
    IndependentVariableName: 'Time'

Создайте график решетки, который показывает профили PK трех индивидуумов.

sbiotrellis(gData,'ID','Time',{'CentralConc','PeripheralConc'},...
            'Marker','+','LineStyle','none');

Пользуйтесь встроенной библиотекой PK, чтобы создать модель 2D отсека с дозированием вливания и устранением первого порядка, где уровень устранения зависит от разрешения и объема центрального отсека. Используйте объект configset включить модульное преобразование.

pkmd                 = PKModelDesign;
pkc1                 = addCompartment(pkmd,'Central');
pkc1.DosingType      = 'Infusion';
pkc1.EliminationType = 'linear-clearance';
pkc1.HasResponseVariable = true;
pkc2                 = addCompartment(pkmd,'Peripheral');
model                = construct(pkmd);
configset            = getconfigset(model);
configset.CompileOptions.UnitConversion = true;

Примите, что каждый индивидуум получает капельное внутривенное введение во время = 0 с общим объемом вливания 100 мг на уровне 50 мг/час. Для получения дополнительной информации при подготовке различных стратегий дозирования, смотрите Дозы в Моделях SimBiology.

dose             = sbiodose('dose','TargetName','Drug_Central');
dose.StartTime   = 0;
dose.Amount      = 100;
dose.Rate        = 50;
dose.AmountUnits = 'milligram';
dose.TimeUnits   = 'hour';
dose.RateUnits   = 'milligram/hour';

Данные содержат измеренные плазменные концентрации в центральных и периферийных отсеках. Сопоставьте эти переменные с соответствующими разновидностями модели, которые являются Drug_Central и Drug_Peripheral.

responseMap = {'Drug_Central = CentralConc','Drug_Peripheral = PeripheralConc'};

Параметры, чтобы оценить в этой модели являются объемами центральных и периферийных отсеков (Central и Peripheral), межразделенное на отсеки разрешение Q12, и уровень раскрываемости преступлений Cl_Central. В этом случае задайте логарифмическое преобразование для Central и Peripheral поскольку они ограничиваются быть положительными. estimatedInfo объект позволяет вам указать, что параметр преобразовывает, начальные значения и (дополнительные) границы параметра.

paramsToEstimate   = {'log(Central)','log(Peripheral)','Q12','Cl_Central'};
estimatedParam     = estimatedInfo(paramsToEstimate,'InitialValue',[1 1 1 1]);

Подбирайте модель ко всем данным, объединенным вместе, то есть, оцените один набор параметров для всех индивидуумов. Метод оценки по умолчанию это sbiofit использование изменится, в зависимости от которого тулбоксы доступны. Видеть который функция оценки sbiofit используемый для подбора кривой, проверяйте EstimationFunction свойство соответствующего объекта результатов.

pooledFit = sbiofit(model,gData,responseMap,estimatedParam,dose,'Pooled',true)
pooledFit = 

  OptimResults with properties:

                   ExitFlag: 3
                     Output: [1x1 struct]
                  GroupName: []
                       Beta: [4x3 table]
         ParameterEstimates: [4x3 table]
                          J: [24x4x2 double]
                       COVB: [4x4 double]
           CovarianceMatrix: [4x4 double]
                          R: [24x2 double]
                        MSE: 6.6220
                        SSE: 291.3688
                    Weights: []
              LogLikelihood: -111.3904
                        AIC: 230.7808
                        BIC: 238.2656
                        DFE: 44
             DependentFiles: {1x3 cell}
    EstimatedParameterNames: {'Central'  'Peripheral'  'Q12'  'Cl_Central'}
             ErrorModelInfo: [1x3 table]
         EstimationFunction: 'lsqnonlin'

Постройте подходящие результаты по сравнению с исходными данными. Несмотря на то, что три отдельных графика были сгенерированы, данные были адаптированы с помощью того же набора параметров (то есть, у всех трех индивидуумов была та же подходящая линия).

plot(pooledFit);

Оцените один набор параметров для каждого индивидуума и смотрите, существует ли какое-либо улучшение оценок параметра. В этом примере с тех пор существует три индивидуума, три набора параметров оцениваются.

unpooledFit =  sbiofit(model,gData,responseMap,estimatedParam,dose,'Pooled',false);

Постройте подходящие результаты по сравнению с исходными данными. Каждый индивидуум был адаптирован по-другому (то есть, каждая подходящая линия уникальна для каждого индивидуума), и каждая линия, казалось, соответствовала хорошо к отдельным данным.

plot(unpooledFit);

Отобразите подходящие результаты первого индивидуума. MSE был ниже, чем та из объединенной подгонки. Это также верно для других двух индивидуумов.

unpooledFit(1)
ans = 

  OptimResults with properties:

                   ExitFlag: 3
                     Output: [1x1 struct]
                  GroupName: 1
                       Beta: [4x3 table]
         ParameterEstimates: [4x3 table]
                          J: [8x4x2 double]
                       COVB: [4x4 double]
           CovarianceMatrix: [4x4 double]
                          R: [8x2 double]
                        MSE: 2.1380
                        SSE: 25.6559
                    Weights: []
              LogLikelihood: -26.4805
                        AIC: 60.9610
                        BIC: 64.0514
                        DFE: 12
             DependentFiles: {1x3 cell}
    EstimatedParameterNames: {'Central'  'Peripheral'  'Q12'  'Cl_Central'}
             ErrorModelInfo: [1x3 table]
         EstimationFunction: 'lsqnonlin'

Сгенерируйте график остаточных значений в зависимости от времени, чтобы сравнить объединенные и необъединенные результаты подгонки. Фигура указывает, что необъединенные подходящие остаточные значения меньше, чем те из объединенной подгонки как ожидалось. В дополнение к сравнению остаточных значений другие строгие критерии могут использоваться, чтобы сравнить подходящие результаты.

t = [gData.Time;gData.Time];
res_pooled = vertcat(pooledFit.R);
res_pooled = res_pooled(:);
res_unpooled = vertcat(unpooledFit.R);
res_unpooled = res_unpooled(:);
plot(t,res_pooled,'o','MarkerFaceColor','w','markerEdgeColor','b')
hold on
plot(t,res_unpooled,'o','MarkerFaceColor','b','markerEdgeColor','b')
refl = refline(0,0); % A reference line representing a zero residual
title('Residuals versus Time');
xlabel('Time');
ylabel('Residuals');
legend({'Pooled','Unpooled'});

Этот пример показал, как выполнить объединенное и необъединенное использование оценок sbiofit. Как проиллюстрировано, необъединенная подгонка составляет изменения из-за определенных предметов в исследовании, и, в этом случае, модель соответствует лучше к данным. Однако объединенная подгонка возвращает параметры всего населения. Если вы хотите оценить параметры всего населения при рассмотрении отдельных изменений использовать sbiofitmixed.

В этом примере показано, как оценить специфичный для категории (такой как молодой по сравнению со старым, штекерным по сравнению с розеткой), отдельно-специфичные, и параметры всего населения с помощью PK профилируют данные от нескольких индивидуумов.

Фон

Предположим, что вы имеете данные о концентрации плазмы препарата от 30 индивидуумов и хотите оценить фармакокинетические параметры, а именно, объемы центрального и периферийного отсека, разрешения и межразделенного на отсеки разрешения. Примите, что концентрация препарата по сравнению с профилем времени следует за снижением biexponential Ct=Ae-at+Be-bt, где Ct концентрация препарата во время t, и a и b наклоны для соответствующих экспоненциальных снижений.

Загрузка данных

Эти синтетические данные содержат ход времени плазменных концентраций 30 индивидуумов после дозы шарика (100 мг), измеренных в разное время и для центральных и для периферийных отсеков. Это также содержит категориальные переменные, а именно, Sex и Age.

clear
load('sd5_302RAgeSex.mat')

Преобразуйте в groupedData Формат

Преобразуйте набор данных в groupedData объект, который является необходимым форматом данных для подходящего функционального sbiofit. groupedData объект также позволяет вам независимую переменную набора и имена переменных группы (если они существуют). Установите модули IDВремя, CentralConc, PeripheralConc, Age, и Sex переменные. Модули являются дополнительными и только необходимыми для UnitConversion покажите, который автоматически преобразует соответствие с физическими количествами в одну сопоставимую модульную систему.

gData = groupedData(data);
gData.Properties.VariableUnits = {'','hour','milligram/liter','milligram/liter','',''};
gData.Properties
ans = struct with fields:
                Description: ''
                   UserData: []
             DimensionNames: {'Row'  'Variables'}
              VariableNames: {1x6 cell}
       VariableDescriptions: {}
              VariableUnits: {1x6 cell}
         VariableContinuity: []
                   RowNames: {}
           CustomProperties: [1x1 matlab.tabular.CustomProperties]
          GroupVariableName: 'ID'
    IndependentVariableName: 'Time'

IndependentVariableName и GroupVariableName свойства были автоматически установлены в Time и ID переменные данных.

Визуализируйте данные

Отобразите данные об ответе для каждого индивидуума.

t = sbiotrellis(gData,'ID','Time',{'CentralConc','PeripheralConc'},...
                'Marker','+','LineStyle','none');
% Resize the figure.
t.hFig.Position(:) = [100 100 1280 800];

Figure contains 30 axes. Axes 1 with title ID 1 contains 2 objects of type line. These objects represent CentralConc, PeripheralConc. Axes 2 with title ID 6 contains 2 objects of type line. Axes 3 with title ID 11 contains 2 objects of type line. Axes 4 with title ID 16 contains 2 objects of type line. Axes 5 with title ID 21 contains 2 objects of type line. Axes 6 with title ID 26 contains 2 objects of type line. Axes 7 with title ID 2 contains 2 objects of type line. Axes 8 with title ID 7 contains 2 objects of type line. Axes 9 with title ID 12 contains 2 objects of type line. Axes 10 with title ID 17 contains 2 objects of type line. Axes 11 with title ID 22 contains 2 objects of type line. Axes 12 with title ID 27 contains 2 objects of type line. Axes 13 with title ID 3 contains 2 objects of type line. Axes 14 with title ID 8 contains 2 objects of type line. Axes 15 with title ID 13 contains 2 objects of type line. Axes 16 with title ID 18 contains 2 objects of type line. Axes 17 with title ID 23 contains 2 objects of type line. Axes 18 with title ID 28 contains 2 objects of type line. Axes 19 with title ID 4 contains 2 objects of type line. Axes 20 with title ID 9 contains 2 objects of type line. Axes 21 with title ID 14 contains 2 objects of type line. Axes 22 with title ID 19 contains 2 objects of type line. Axes 23 with title ID 24 contains 2 objects of type line. Axes 24 with title ID 29 contains 2 objects of type line. Axes 25 with title ID 5 contains 2 objects of type line. Axes 26 with title ID 10 contains 2 objects of type line. Axes 27 with title ID 15 contains 2 objects of type line. Axes 28 with title ID 20 contains 2 objects of type line. Axes 29 with title ID 25 contains 2 objects of type line. Axes 30 with title ID 30 contains 2 objects of type line.

Настройте модель 2D отсека

Пользуйтесь встроенной библиотекой PK, чтобы создать модель 2D отсека с дозированием вливания и устранением первого порядка, где уровень устранения зависит от разрешения и объема центрального отсека. Используйте configset объект включить модульное преобразование.

pkmd                                    = PKModelDesign;
pkc1                                    = addCompartment(pkmd,'Central');
pkc1.DosingType                         = 'Bolus';
pkc1.EliminationType                    = 'linear-clearance';
pkc1.HasResponseVariable                = true;
pkc2                                    = addCompartment(pkmd,'Peripheral');
model                                   = construct(pkmd);
configset                               = getconfigset(model);
configset.CompileOptions.UnitConversion = true;

Для получения дополнительной информации при создании разделенных на отсеки моделей PK с помощью SimBiology® встроенная библиотека, смотрите, Создают Фармакокинетические Модели.

Задайте дозирование

Примите, что каждый индивидуум получает дозу шарика 100 мг во время = 0. Для получения дополнительной информации при подготовке различных стратегий дозирования, смотрите Дозы в Моделях SimBiology.

dose             = sbiodose('dose','TargetName','Drug_Central');
dose.StartTime   = 0;
dose.Amount      = 100;
dose.AmountUnits = 'milligram';
dose.TimeUnits   = 'hour';

Сопоставьте данные об ответе с соответствующими компонентами модели

Данные содержат измеренную плазменную концентрацию в центральных и периферийных отсеках. Сопоставьте эти переменные с соответствующими компонентами модели, которые являются Drug_Central и Drug_Peripheral.

responseMap = {'Drug_Central = CentralConc','Drug_Peripheral = PeripheralConc'};

Задайте параметры, чтобы оценить

Задайте объемы центральных и периферийных отсеков Central и Peripheral, межразделенное на отсеки разрешение Q12, и разрешение Cl_Central как параметры, чтобы оценить. estimatedInfo объект позволяет вам опционально указать, что параметр преобразовывает, начальные значения и границы параметра. Начиная с обоих Central и Peripheral ограничиваются быть положительным, задать логарифмическое преобразование для каждого параметра.

paramsToEstimate    = {'log(Central)', 'log(Peripheral)', 'Q12', 'Cl_Central'};
estimatedParam      = estimatedInfo(paramsToEstimate,'InitialValue',[1 1 1 1]);

Оцените отдельно-специфичные параметры

Оцените один набор параметров для каждого индивидуума путем установки 'Pooled' аргумент пары "имя-значение" false.

unpooledFit =  sbiofit(model,gData,responseMap,estimatedParam,dose,'Pooled',false);

Отображение результатов

Постройте подходящие результаты по сравнению с исходными данными для каждого индивидуума (группа).

t = plot(unpooledFit);
% Resize the figure.
t.hFig.Position(:) = [100 100 1280 800];

Figure contains 30 axes. Axes 1 with title Group 1 contains 4 objects of type line. These objects represent OBS1 (CentralConc), OBS2 (PeripheralConc), PRED1 (Central.Drug_Central), PRED2 (Peripheral.Drug_Peripheral). Axes 2 with title Group 6 contains 4 objects of type line. Axes 3 with title Group 11 contains 4 objects of type line. Axes 4 with title Group 16 contains 4 objects of type line. Axes 5 with title Group 21 contains 4 objects of type line. Axes 6 with title Group 26 contains 4 objects of type line. Axes 7 with title Group 2 contains 4 objects of type line. Axes 8 with title Group 7 contains 4 objects of type line. Axes 9 with title Group 12 contains 4 objects of type line. Axes 10 with title Group 17 contains 4 objects of type line. Axes 11 with title Group 22 contains 4 objects of type line. Axes 12 with title Group 27 contains 4 objects of type line. Axes 13 with title Group 3 contains 4 objects of type line. Axes 14 with title Group 8 contains 4 objects of type line. Axes 15 with title Group 13 contains 4 objects of type line. Axes 16 with title Group 18 contains 4 objects of type line. Axes 17 with title Group 23 contains 4 objects of type line. Axes 18 with title Group 28 contains 4 objects of type line. Axes 19 with title Group 4 contains 4 objects of type line. Axes 20 with title Group 9 contains 4 objects of type line. Axes 21 with title Group 14 contains 4 objects of type line. Axes 22 with title Group 19 contains 4 objects of type line. Axes 23 with title Group 24 contains 4 objects of type line. Axes 24 with title Group 29 contains 4 objects of type line. Axes 25 with title Group 5 contains 4 objects of type line. Axes 26 with title Group 10 contains 4 objects of type line. Axes 27 with title Group 15 contains 4 objects of type line. Axes 28 with title Group 20 contains 4 objects of type line. Axes 29 with title Group 25 contains 4 objects of type line. Axes 30 with title Group 30 contains 4 objects of type line.

Для необъединенной подгонки, sbiofit всегда возвращается, каждый заканчивается объект для каждого индивидуума.

Исследуйте оценки параметра на зависимости от категории

Исследуйте необъединенные оценки, чтобы видеть, существуют ли какие-либо специфичные для категории параметры, то есть, если некоторые параметры связаны с одной или несколькими категориями. Если существуют какие-либо зависимости от категории, может быть возможно уменьшать количество степеней свободы путем оценки только специфичных для категории значений для тех параметров.

Сначала извлеките ID и значения категории для каждого ID

catParamValues = unique(gData(:,{'ID','Sex','Age'}));

Добавьте переменные в таблицу, содержащую оценку каждого параметра.

allParamValues            = vertcat(unpooledFit.ParameterEstimates);
catParamValues.Central    = allParamValues.Estimate(strcmp(allParamValues.Name, 'Central'));
catParamValues.Peripheral = allParamValues.Estimate(strcmp(allParamValues.Name, 'Peripheral'));
catParamValues.Q12        = allParamValues.Estimate(strcmp(allParamValues.Name, 'Q12'));
catParamValues.Cl_Central = allParamValues.Estimate(strcmp(allParamValues.Name, 'Cl_Central'));

Постройте оценки каждого параметра для каждой категории. gscatter требует Statistics and Machine Learning Toolbox™. Если у вас нет его, используйте другие альтернативные функции построения графика, такие как plot.

h           = figure;
ylabels     = {'Central','Peripheral','Cl\_Central','Q12'};
plotNumber  = 1;
for i = 1:4
    thisParam = estimatedParam(i).Name;
    
    % Plot for Sex category
    subplot(4,2,plotNumber);
    plotNumber  = plotNumber + 1;
    gscatter(double(catParamValues.Sex), catParamValues.(thisParam), catParamValues.Sex);
    ax          = gca;
    ax.XTick    = [];
    ylabel(ylabels(i));
    legend('Location','bestoutside')
    % Plot for Age category
    subplot(4,2,plotNumber);
    plotNumber  = plotNumber + 1;
    gscatter(double(catParamValues.Age), catParamValues.(thisParam), catParamValues.Age);
    ax          = gca;
    ax.XTick    = [];
    ylabel(ylabels(i));
    legend('Location','bestoutside')
end
% Resize the figure.
h.Position(:) = [100 100 1280 800];

Figure contains 8 axes. Axes 1 contains 2 objects of type line. These objects represent Female, Male. Axes 2 contains 2 objects of type line. These objects represent Old, Young. Axes 3 contains 2 objects of type line. These objects represent Female, Male. Axes 4 contains 2 objects of type line. These objects represent Old, Young. Axes 5 contains 2 objects of type line. These objects represent Female, Male. Axes 6 contains 2 objects of type line. These objects represent Old, Young. Axes 7 contains 2 objects of type line. These objects represent Female, Male. Axes 8 contains 2 objects of type line. These objects represent Old, Young.

На основе графика кажется, что молодые индивидуумы склонны иметь более высокие объемы центральных и периферийных отсеков (Central, Peripheral) чем старые индивидуумы (то есть, объемы, кажется, являются возрастными). Кроме того, штекеры имеют тенденцию иметь более высокий уровень раскрываемости преступлений (Cl_Central) чем розетки, но противоположное для параметра Q12 (то есть, разрешение и Q12, кажется, сексуально-специфичны).

Оцените специфичные для категории параметры

Используйте 'CategoryVariableName' свойство estimatedInfo объект задать, который категория использовать во время подбора кривой. Используйте 'Sex' как группа к пригодному для разрешения Cl_Central и Q12 параметры. Используйте 'Age' как группа к пригодному для Central и Peripheral параметры.

estimatedParam(1).CategoryVariableName = 'Age';
estimatedParam(2).CategoryVariableName = 'Age';
estimatedParam(3).CategoryVariableName = 'Sex';
estimatedParam(4).CategoryVariableName = 'Sex';
categoryFit = sbiofit(model,gData,responseMap,estimatedParam,dose)
categoryFit = 
  OptimResults with properties:

                   ExitFlag: 3
                     Output: [1x1 struct]
                  GroupName: []
                       Beta: [8x5 table]
         ParameterEstimates: [120x6 table]
                          J: [240x8x2 double]
                       COVB: [8x8 double]
           CovarianceMatrix: [8x8 double]
                          R: [240x2 double]
                        MSE: 0.4362
                        SSE: 205.8690
                    Weights: []
              LogLikelihood: -477.9195
                        AIC: 971.8390
                        BIC: 1.0052e+03
                        DFE: 472
             DependentFiles: {1x3 cell}
    EstimatedParameterNames: {'Central'  'Peripheral'  'Q12'  'Cl_Central'}
             ErrorModelInfo: [1x3 table]
         EstimationFunction: 'lsqnonlin'

При подборе кривой по категориям (или группа), sbiofit всегда возвращается, каждый заканчивается объект, не один для каждого уровня категории. Это вызвано тем, что и штекерные и розеточные индивидуумы считаются частью той же оптимизации с помощью тех же ошибочных параметров модели и параметров ошибок, так же для молодых и старых индивидуумов.

Постройте результаты

Постройте специфичные для категории предполагаемые результаты.

t = plot(categoryFit);
% Resize the figure.
t.hFig.Position(:) = [100 100 1280 800];

Figure contains 30 axes. Axes 1 with title Group 1 contains 4 objects of type line. These objects represent OBS1 (CentralConc), OBS2 (PeripheralConc), PRED1 (Central.Drug_Central), PRED2 (Peripheral.Drug_Peripheral). Axes 2 with title Group 6 contains 4 objects of type line. Axes 3 with title Group 11 contains 4 objects of type line. Axes 4 with title Group 16 contains 4 objects of type line. Axes 5 with title Group 21 contains 4 objects of type line. Axes 6 with title Group 26 contains 4 objects of type line. Axes 7 with title Group 2 contains 4 objects of type line. Axes 8 with title Group 7 contains 4 objects of type line. Axes 9 with title Group 12 contains 4 objects of type line. Axes 10 with title Group 17 contains 4 objects of type line. Axes 11 with title Group 22 contains 4 objects of type line. Axes 12 with title Group 27 contains 4 objects of type line. Axes 13 with title Group 3 contains 4 objects of type line. Axes 14 with title Group 8 contains 4 objects of type line. Axes 15 with title Group 13 contains 4 objects of type line. Axes 16 with title Group 18 contains 4 objects of type line. Axes 17 with title Group 23 contains 4 objects of type line. Axes 18 with title Group 28 contains 4 objects of type line. Axes 19 with title Group 4 contains 4 objects of type line. Axes 20 with title Group 9 contains 4 objects of type line. Axes 21 with title Group 14 contains 4 objects of type line. Axes 22 with title Group 19 contains 4 objects of type line. Axes 23 with title Group 24 contains 4 objects of type line. Axes 24 with title Group 29 contains 4 objects of type line. Axes 25 with title Group 5 contains 4 objects of type line. Axes 26 with title Group 10 contains 4 objects of type line. Axes 27 with title Group 15 contains 4 objects of type line. Axes 28 with title Group 20 contains 4 objects of type line. Axes 29 with title Group 25 contains 4 objects of type line. Axes 30 with title Group 30 contains 4 objects of type line.

Для Cl_Central и Q12 параметры, все штекеры имели те же оценки, и так же для розеток. Для Central и Peripheral параметры, у всех молодых индивидуумов были те же оценки, и так же для старых индивидуумов.

Оцените параметры всего населения

Чтобы лучше сравнить результаты, подбирайте модель ко всем данным, объединенным вместе, то есть, оцените один набор параметров для всех индивидуумов путем установки 'Pooled' аргумент пары "имя-значение" true. Предупреждающее сообщение говорит вам, что эта опция проигнорирует любую специфичную для категории информацию (если они будут существовать).

pooledFit = sbiofit(model,gData,responseMap,estimatedParam,dose,'Pooled',true);
Warning: CategoryVariableName property of the estimatedInfo object is ignored when using the 'Pooled' option.

Постройте результаты

Постройте подходящие результаты по сравнению с исходными данными. Несмотря на то, что отдельный график был сгенерирован для каждого индивидуума, данные были адаптированы с помощью того же набора параметров (то есть, у всех индивидуумов была та же подходящая линия).

t = plot(pooledFit);
% Resize the figure.
t.hFig.Position(:) = [100 100 1280 800];

Figure contains 30 axes. Axes 1 with title Group 1 contains 4 objects of type line. These objects represent OBS1 (CentralConc), OBS2 (PeripheralConc), PRED1 (Central.Drug_Central), PRED2 (Peripheral.Drug_Peripheral). Axes 2 with title Group 6 contains 4 objects of type line. Axes 3 with title Group 11 contains 4 objects of type line. Axes 4 with title Group 16 contains 4 objects of type line. Axes 5 with title Group 21 contains 4 objects of type line. Axes 6 with title Group 26 contains 4 objects of type line. Axes 7 with title Group 2 contains 4 objects of type line. Axes 8 with title Group 7 contains 4 objects of type line. Axes 9 with title Group 12 contains 4 objects of type line. Axes 10 with title Group 17 contains 4 objects of type line. Axes 11 with title Group 22 contains 4 objects of type line. Axes 12 with title Group 27 contains 4 objects of type line. Axes 13 with title Group 3 contains 4 objects of type line. Axes 14 with title Group 8 contains 4 objects of type line. Axes 15 with title Group 13 contains 4 objects of type line. Axes 16 with title Group 18 contains 4 objects of type line. Axes 17 with title Group 23 contains 4 objects of type line. Axes 18 with title Group 28 contains 4 objects of type line. Axes 19 with title Group 4 contains 4 objects of type line. Axes 20 with title Group 9 contains 4 objects of type line. Axes 21 with title Group 14 contains 4 objects of type line. Axes 22 with title Group 19 contains 4 objects of type line. Axes 23 with title Group 24 contains 4 objects of type line. Axes 24 with title Group 29 contains 4 objects of type line. Axes 25 with title Group 5 contains 4 objects of type line. Axes 26 with title Group 10 contains 4 objects of type line. Axes 27 with title Group 15 contains 4 objects of type line. Axes 28 with title Group 20 contains 4 objects of type line. Axes 29 with title Group 25 contains 4 objects of type line. Axes 30 with title Group 30 contains 4 objects of type line.

Сравните остаточные значения

Сравните остаточные значения CentralConc и PeripheralConc ответы для каждой подгонки.

t = gData.Time;
allResid(:,:,1) = pooledFit.R;
allResid(:,:,2) = categoryFit.R;
allResid(:,:,3) = vertcat(unpooledFit.R);

h = figure;
responseList = {'CentralConc', 'PeripheralConc'};
for i = 1:2
    subplot(2,1,i);
    oneResid = squeeze(allResid(:,i,:));
    plot(t,oneResid,'o');
    refline(0,0); % A reference line representing a zero residual
    title(sprintf('Residuals (%s)', responseList{i}));
    xlabel('Time');
    ylabel('Residuals');
    legend({'Pooled','Category-Specific','Unpooled'});
end
% Resize the figure.
h.Position(:) = [100 100 1280 800];

Figure contains 2 axes. Axes 1 with title Residuals (CentralConc) contains 4 objects of type line. These objects represent Pooled, Category-Specific, Unpooled. Axes 2 with title Residuals (PeripheralConc) contains 4 objects of type line. These objects represent Pooled, Category-Specific, Unpooled.

Как показано в графике, необъединенная подгонка произвела лучшую подгонку к данным, когда это соответствует данным к каждому индивидууму. Это ожидалось, поскольку это использовало большую часть количества степеней свободы. Подгонка категории уменьшала количество степеней свободы подгонкой данных к двум категориям (пол и возраст). В результате остаточные значения были больше, чем необъединенная подгонка, но еще меньшими, чем подгонка населения, которая оценила всего один набор параметров для всех индивидуумов. Подходящая категорией сила быть хорошим компромиссом между необъединенным и объединенным подбором кривой при условии, что любая иерархическая модель существует в ваших данных.

Этот пример использует дрожжи гетеротримерные данные модели белка G и экспериментальные данные, о которых сообщают [1]. Для получения дополнительной информации о модели, смотрите раздел Background в Сканировании Параметра, Оценке Параметра и Анализе чувствительности в Дрожжах Гетеротримерный Цикл Белка G.

Загрузите модель белка G.

sbioloadproject gprotein

Храните экспериментальные данные, содержащие курс времени для части активного белка G.

time = [0 10 30 60 110 210 300 450 600]';
GaFracExpt = [0 0.35 0.4 0.36 0.39 0.33 0.24 0.17 0.2]';

Создайте groupedData основанный на объектах на экспериментальных данных.

tbl = table(time,GaFracExpt);
grpData = groupedData(tbl);

Сопоставьте соответствующий компонент модели с экспериментальными данными. Другими словами, укажите, которому разновидность в модели соответствует который переменная отклика в данных. В этом примере сопоставьте параметр модели GaFrac к переменной GaFracExpt экспериментальных данных от grpData.

responseMap = 'GaFrac = GaFracExpt';

Используйте estimatedInfo объект задать параметр модели kGd в качестве параметра быть оцененным.

estimatedParam = estimatedInfo('kGd');

Выполните оценку параметра.

fitResult = sbiofit(m1,grpData,responseMap,estimatedParam);

Просмотрите предполагаемое значение параметров kGd.

fitResult.ParameterEstimates
ans=1×3 table
     Name      Estimate    StandardError
    _______    ________    _____________

    {'kGd'}    0.11307      3.4439e-05  

Предположим, что вы хотите построить результаты симуляции модели с помощью предполагаемого значения параметров. Можно или повторно выполнить sbiofit функционируйте и задайте, чтобы возвратить дополнительный второй выходной аргумент, который содержит результаты симуляции, или используйте fitted метод, чтобы получить результаты, не повторно выполняя sbiofit.

[yfit,paramEstim] = fitted(fitResult);

Постройте результаты симуляции.

sbioplot(yfit);

Figure contains an axes. The axes with title States versus Time contains 7 objects of type line. These objects represent G, Gd, Ga, RL, R, Gbg, GaFrac.

В этом примере показано, как оценить задержку, прежде чем доза шарика была введена и длительность дозы с помощью модели с одним отсеком.

Загрузите набор выборочных данных.

load lagDurationData.mat

Отобразите данные на графике.

plot(data.Time,data.Conc,'x')
xlabel('Time (hour)')
ylabel('Conc (milligram/liter)')

Figure contains an axes. The axes contains an object of type line.

Преобразуйте в groupedData.

gData = groupedData(data);
gData.Properties.VariableUnits = {'hour','milligram/liter'};

Создайте модель с одним отсеком.

pkmd                    = PKModelDesign;
pkc1                    = addCompartment(pkmd,'Central');
pkc1.DosingType         = 'Bolus';
pkc1.EliminationType    = 'linear-clearance';
pkc1.HasResponseVariable = true;
model                   = construct(pkmd);
configset               = getconfigset(model);
configset.CompileOptions.UnitConversion = true;

Добавьте два параметра, которые представляют задержку и длительность дозы. Параметр задержки задает задержку, прежде чем доза будет введена. Параметр длительности задает отрезок времени, который требуется, чтобы ввести дозу.

lagP = addparameter(model,'lagP');
lagP.ValueUnits = 'hour';
durP = addparameter(model,'durP');
durP.ValueUnits = 'hour';

Создайте объект дозы. Установите LagParameterName и DurationParameterName свойства дозы к именам параметров задержки и длительности, соответственно. Установите суммарную дозу на 10 миллиграммов, которая была суммой, используемой, чтобы сгенерировать данные.

dose                = sbiodose('dose');
dose.TargetName     = 'Drug_Central';
dose.StartTime      = 0;
dose.Amount         = 10;
dose.AmountUnits    = 'milligram';
dose.TimeUnits      = 'hour';
dose.LagParameterName = 'lagP';
dose.DurationParameterName = 'durP';

Сопоставьте разновидности модели с соответствующими данными.

responseMap = {'Drug_Central = Conc'};

Задайте параметры задержки и длительности как параметры, чтобы оценить. Логарифмическое преобразование параметры. Инициализируйте их к 2 и установите верхнюю границу и нижнюю границу.

paramsToEstimate    = {'log(lagP)','log(durP)'};
estimatedParams     = estimatedInfo(paramsToEstimate,'InitialValue',2,'Bounds',[1 5]);

Выполните оценку параметра.

fitResults = sbiofit(model,gData,responseMap,estimatedParams,dose,'fminsearch')
fitResults = 
  OptimResults with properties:

                   ExitFlag: 1
                     Output: [1x1 struct]
                  GroupName: One group
                       Beta: [2x4 table]
         ParameterEstimates: [2x4 table]
                          J: [11x2 double]
                       COVB: [2x2 double]
           CovarianceMatrix: [2x2 double]
                          R: [11x1 double]
                        MSE: 0.0024
                        SSE: 0.0213
                    Weights: []
              LogLikelihood: 18.7511
                        AIC: -33.5023
                        BIC: -32.7065
                        DFE: 9
             DependentFiles: {1x2 cell}
    EstimatedParameterNames: {'lagP'  'durP'}
             ErrorModelInfo: [1x3 table]
         EstimationFunction: 'fminsearch'

Отобразите результат.

fitResults.ParameterEstimates
ans=2×4 table
      Name      Estimate    StandardError    Bounds
    ________    ________    _____________    ______

    {'lagP'}     1.986        0.0051568      1    5
    {'durP'}     1.527         0.012956      1    5

plot(fitResults)

Figure contains an axes. The axes with title Group One group contains 2 objects of type line. These objects represent OBS1 (Conc), PRED1 (Central.Drug_Central).

Входные параметры

свернуть все

Модель SimBiology в виде SimBiology model object. Активное configset object из модели содержит настройки решателя для симуляции. Любые активные дозы и варианты применяются к модели в процессе моделирования, если не задано в противном случае использование dosing и variants входные параметры, соответственно.

Данные, чтобы соответствовать в виде a groupedData object.

Имя переменной времени должно быть задано в IndependentVariableName свойство grpData. Например, если именем переменной времени является 'TIME', затем задайте его можно следующим образом.

grpData.Properties.IndependentVariableName = 'TIME';

Если данные содержат больше чем одну группу измерений, имя сгруппированной переменной должно быть задано в GroupVariableName свойство grpData. Например, если именем сгруппированной переменной является 'GROUP', затем задайте его можно следующим образом.

grpData.Properties.GroupVariableName = 'GROUP';
Группа обычно обращается к набору измерений, которые представляют один курс времени, часто соответствуя конкретному индивидууму или экспериментальному условию.

Примечание

sbiofit использует categorical функция, чтобы идентифицировать группы. Если какие-либо значения группы преобразованы в то же значение categorical, затем те наблюдения обработаны как принадлежащий той же группе. Например, если некоторые наблюдения не имеют никакой информации о группе (то есть, пустой символьный вектор ''то categorical преобразует пустые символьные вектора в <undefined>, и эти наблюдения обработаны как одна группа.

Отображение информации компонентов модели к grpDataВ виде вектора символов, строки, вектора строки или массива ячеек из символьных векторов.

Каждый вектор символов или строка являются подобным уравнению выражением, похожим на правила присвоения в SimBiology. Это содержит имя (или полностью определенное имя) количества (разновидности, отсек или параметр) или observable объект в модели sm, сопровождаемый символьным '=' и имя переменной в grpData. Для ясности пробелы позволены между именами и '='.

Например, если у вас есть данные о концентрации 'CONC' в grpData для разновидности 'Drug_Central', можно задать его можно следующим образом.

responseMap = 'Drug_Central = CONC';

Чтобы назвать разновидность однозначно, используйте полностью определенное имя, которое включает имя отсека. Чтобы назвать ограниченный по объему реакцией параметр, используйте имя реакции, чтобы квалифицировать параметр.

Если имя компонента модели или grpData имя переменной не является допустимым именем переменной MATLAB®, окружите его квадратными скобками, такими как:

responseMap = '[Central 1].Drug = [Central 1 Conc]';

Если само имя переменной содержит квадратные скобки, вы не можете использовать его в выражении, чтобы задать информацию об отображении.

Ошибка выпущена, если какое-либо (квалифицированное) имя совпадает с двумя компонентами того же типа. Однако можно использовать (квалифицированное) имя, которое совпадает с двумя компонентами различных типов, и функция сначала находит разновидности с именем, сопровождаемым отсеками и затем параметрами.

Предполагаемые параметры в виде EstimatedInfo object или вектор из estimatedInfo объекты, который задает предполагаемые параметры в модели sm, и другая дополнительная информация, такая как их первоначальные оценки, преобразования, связала ограничения и категории. Поддерживаемыми преобразованиями является log, logit, и probit. Для получения дополнительной информации смотрите Преобразования Параметра.

Можно задать границы для всех методов оценки. Нижняя граница должна быть меньше верхней границы. Для получения дополнительной информации см. boundValues.

При использовании scattersearch, необходимо задать конечные преобразованные границы для каждого предполагаемого параметра.

При использовании fminsearch, nlinfit, или fminunc с границами целевая функция возвращает Inf если границы превышены. Когда вы включаете опции, такие как FunValCheck, оптимизация может ошибка, если границы превышены во время оценки. При использовании nlinfit, это может сообщить о предупреждениях о якобиане, являющемся плохо обусловленным или не бывшем способном оценивать, слишком близок ли конечный результат к границам.

Если вы не задаете Pooled аргумент пары "имя-значение", sbiofit использование CategoryVariableName свойство estiminfo решить, должны ли параметры быть оценены для каждого индивидуума, группы, категории или всех индивидуумов в целом. Используйте Pooled опция, чтобы заменить любой CategoryVariableName значения. Для получения дополнительной информации о CategoryVariableName свойство, смотрите EstimatedInfo object.

Примечание

sbiofit использует categorical функция, чтобы идентифицировать группы или категории. Если какие-либо значения группы преобразованы в то же значение categorical, затем те наблюдения обработаны как принадлежащий той же группе. Например, если некоторые наблюдения не имеют никакой информации о группе (то есть, пустой символьный вектор '' как значение группы), затем categorical преобразует пустые символьные вектора в <undefined>, и эти наблюдения обработаны как одна группа.

Информация о дозах в виде пустого массива или 2D матрица объектов дозы (ScheduleDose object или RepeatDose object).

  • Если пустой, никакие дозы не применяются в процессе моделирования, даже если модель имеет активные дозы.

  • Если не пустой, матрица должна иметь одну строку или одну строку на группу во входных данных. Если это имеет одну строку, те же дозы применяются ко всем группам в процессе моделирования. Если это имеет несколько строк, каждая строка применяется к отдельной группе в том же порядке, как группы появляются во входных данных.

  • Несколько столбцов позволены так, чтобы можно было применить несколько объектов дозы к каждой группе. Все дозы в столбце должны сослаться на те же компоненты в модели. Например, дозы в том же столбце должны иметь ту же цель дозы (TargetName). Если вы параметрируете какое-либо свойство дозы, то всем дозам в рамках столбца нужно было параметрировать это свойство к тому же параметру. Если некоторые группы требуют большего количества доз, чем другие, то заполняют матрицу с (фиктивными) дозами по умолчанию.

  • Доза по умолчанию имеет значения по умолчанию для всех свойств, за исключением Name свойство. Создайте дозу по умолчанию можно следующим образом.

    d1 = sbiodose('d1');

  • В дополнение к ручному построению объектов дозы с помощью sbiodose, если grpData имеет информацию о дозах, можно использовать createDoses метод, чтобы создать дозы.

Имя функции оценки в виде вектора символов или строки. Выбор следующие.

  • 'fminsearch'

  • 'nlinfit' (Statistics and Machine Learning Toolbox требуется.)

  • 'fminunc' (Optimization Toolbox требуется.)

  • 'fmincon' (Optimization Toolbox требуется.)

  • 'lsqcurvefit' (Optimization Toolbox требуется.)

  • 'lsqnonlin' (Optimization Toolbox требуется.)

  • 'patternsearch' (Global Optimization Toolbox требуется.)

  • 'ga' (Global Optimization Toolbox требуется.)

  • 'particleswarm' (Global Optimization Toolbox требуется.)

  • 'scattersearch'

По умолчанию, sbiofit использует первую доступную функцию оценки среди следующего: lsqnonlin (Optimization Toolbox), nlinfit (Statistics and Machine Learning Toolbox), или fminsearch. Тот же приоритет применяется к локальному выбору решателя по умолчанию за scattersearch.

Для сводных данных поддерживаемых методов и подходящих опций, см. Поддерживаемые Методы для Оценки Параметра SimBiology.

Опции, характерные для оценки, функционируют в виде struct или optimoptions объект.

  • statset struct для nlinfit

  • optimset struct для fminsearch

  • optimoptions объект для lsqcurvefit, lsqnonlin, fmincon, fminunc, particleswarm, ga, и patternsearch

  • struct для scattersearch

Смотрите Опции по умолчанию для Алгоритмов Оценки для получения дополнительной информации и опции по умолчанию, сопоставленные с каждой функцией оценки.

Варианты в виде пустого массива или вектора из различных объектов. Если пустой, никакие варианты не используются, даже если модель имеет активные варианты.

Аргументы в виде пар имя-значение

Задайте дополнительные разделенные запятой пары Name,Value аргументы. Name имя аргумента и Value соответствующее значение. Name должен появиться в кавычках. Вы можете задать несколько аргументов в виде пар имен и значений в любом порядке, например: Name1, Value1, ..., NameN, ValueN.

Пример: 'ErrorModel','constant','UseParallel',true задает постоянную ошибочную модель и запускать параллельные симуляции во время оценки параметра.

Ошибочные модели, используемые для оценки в виде вектора символов, строки, представляют в виде строки вектор, массив ячеек из символьных векторов, категориальный вектор или таблицу.

Если это - таблица, это должно содержать одну переменную, которая является вектор-столбцом ошибочных имен модели. Имена могут быть массивом ячеек из символьных векторов, вектором строки или вектором из категориальных переменных. Если таблица имеет больше чем одну строку, то RowNames свойство должно совпадать с именами переменной отклика, заданными в правой стороне responseMap. Если таблица не использует RowNames свойство, n th ошибка сопоставлен с n th ответ.

Если вы задаете только одну ошибочную модель, то sbiofit оценки один набор параметров ошибок для всех ответов.

Если вы задаете несколько ошибочных моделей с помощью категориального вектора, вектора строки, или массив ячеек из символьных векторов, длина векторного массива или массива ячеек должна совпадать с количеством ответов в responseMap аргумент.

Можно задать несколько ошибочных моделей, только если вы используете эти методы: lsqnonlin, lsqcurvefit, fmincon, fminuncfminsearch, patternsearch, ga, и particleswarm.

Существует четыре встроенных ошибочных модели. Каждая модель задает ошибку стандартный средний нуль и переменная (Gaussian) модульного отклонения e, результаты симуляции f и один или два параметра a и b.

  • 'constant': y=f+ae

  • 'proportional': y=f+b|f|e

  • 'combined': y=f+(a+b|f|)e

  • 'exponential': y=fexp(ae)

Примечание

  • Если вы задаете ошибочную модель, вы не можете задать веса за исключением постоянной ошибочной модели.

  • Если вы используете пропорциональную или объединенную ошибочную модель во время подбора кривой данных, стараетесь не задавать точки данных время от времени, где решение (результат симуляции) является нулем, или устойчивое состояние является нулем. В противном случае можно столкнуться с проблемами деления на нуль. Рекомендуется, чтобы вы удалили те точки данных из набора данных. Для получения дополнительной информации на функциях оценки параметра ошибок, смотрите Оценку Наибольшего правдоподобия.

Веса использовали для оценки в виде матрицы действительных положительных весов, где количество столбцов соответствует количеству ответов или указателю на функцию, который принимает вектор из предсказанных значений отклика и возвращает вектор из действительных положительных весов.

Если вы задаете ошибочную модель, вы не можете использовать веса за исключением постоянной ошибочной модели. Если ни один 'ErrorModel' или 'Weights' задан, по умолчанию функция использует постоянную ошибочную модель с равными весами.

Подходящий флаг опции, задавая, соответствовать ли каждому индивидууму (false) или объедините все отдельные данные (true).

Когда true, функция выполняет подхождение всем индивидуумам или группам, одновременно использующим те же оценки параметра и fitResults скалярный объект результатов.

Когда false, функция соответствует каждой группе или индивидууму, отдельно использующему группу - или отдельно-специфичные параметры и fitResults вектор из объектов результатов одним результатом для каждой группы.

Примечание

Используйте эту опцию, чтобы заменить любой CategoryVariableName значения estiminfo.

Отметьте, чтобы включить распараллеливание в виде true или false. Если true и Parallel Computing Toolbox™ доступен, sbiofit поддержки несколько уровней распараллеливания, но только один уровень используется за один раз.

Для необъединенной подгонки ('Pooled' = false) для нескольких групп каждая подгонка запущена параллельно.

Для объединенной подгонки ('Pooled' = true), распараллеливание происходит на уровне решателя. Другими словами, расчеты решателя, такие как оценки целевой функции, запущены параллельно. Для получения дополнительной информации смотрите Несколько Оценок Параметра параллельно.

Отметьте, чтобы использовать чувствительность параметра, чтобы определить градиенты целевой функции в виде true или false. По умолчанию это - true для fmincon, fminunc, lsqnonlin, и lsqcurvefit методы. В противном случае это - false. Если это верно, SimBiology всегда использует sundials решатель, независимо от того, что вы выбрали как SolverType свойство в Configset объект.

SimBiology использует метод приближения комплексного шага, чтобы вычислить чувствительность параметра. Такая расчетная чувствительность может использоваться, чтобы определить градиенты целевой функции во время оценки параметра, чтобы улучшить подбор кривой. Поведение по умолчанию sbiofit должен использовать такую чувствительность, чтобы определить градиенты каждый раз, когда алгоритм основан на градиенте и если модель SimBiology® поддерживает анализ чувствительности. Для получения дополнительной информации о требованиях модели и анализе чувствительности, смотрите Анализ чувствительности в SimBiology.

Отметьте, чтобы показать прогресс оценки параметра в виде true или false. Если true, новая фигура открывается содержащий графики, которые показывают логарифмическую правдоподобность, критерии завершения и оцененные параметры для каждой итерации. Эта опция не поддерживается для nlinfit метод.

Если вы используете функции Optimization Toolbox (fminunc, fmincon, lsqcurvefit, lsqnonlin), рисунок также показывает Оптимальность Первого порядка (Optimization Toolbox) график.

Для необъединенной подгонки каждая линия на графиках представляет индивидуума. Для объединенной подгонки одна строка представляет всех индивидуумов. Линия становится увядшей, когда подгонка завершена. Графики также отслеживают прогресс, когда вы запускаетесь sbiofit (и для объединенных и для необъединенных подгонок) в параллельных использующих удаленных кластерах. Для получения дополнительной информации см. График Прогресса.

Выходные аргументы

свернуть все

Результаты оценки, возвращенные как a OptimResults object или NLINResults object или вектор из этих объектов. Оба объекта результатов являются подклассами LeastSquaresResults object.

Если функция использует nlinfit (Statistics and Machine Learning Toolbox), затем fitResults isa NLINResults object. В противном случае, fitResults OptimResults object.

Когда 'Pooled' установлен в false, функция соответствует каждой группе, отдельно использующей специфичные для группы параметры и fitResults вектор из объектов результатов с, каждый заканчивается объект для каждой группы.

Когда 'Pooled' установлен в true, функция выполняет подхождение всем индивидуумам или группам, одновременно использующим те же оценки параметра и fitResults скалярный объект результатов.

Когда 'Pooled' не используется, и CategoryVariableName значения estiminfo весь <none>, fitResults один объект результата. Это - то же поведение как установка 'Pooled' к true.

Когда 'Pooled' не используется, и CategoryVariableName значения estiminfo весь <GroupVariableName>, fitResults вектор из объектов результатов. Это - то же поведение как установка 'Pooled' к false.

Во всех других случаях, fitResults скалярный объект, содержащий оцененные значения параметров для различных групп или категорий, заданных CategoryVariableName.

Результаты симуляции, возвращенные как вектор из SimData объекты, представляющие результаты симуляции для каждой группы или индивидуума.

Если 'Pooled' опцией является false, затем каждая симуляция использует специфичные для группы значения параметров. Если true, затем все симуляции используют те же значения параметров (всего населения).

О состояниях сообщают в simdata состояния, которые были включены в responseMap входной параметр и любые другие состояния, перечисленные в StatesToLog свойство опций во время выполнения (RuntimeOptions) из модели SimBiology sm.

Больше о

свернуть все

Оценка наибольшего правдоподобия

SimBiology оценивает параметры методом наибольшего правдоподобия. Вместо того, чтобы непосредственно максимизировать функцию правдоподобия, SimBiology создает эквивалентную проблему минимизации. Каждый раз, когда возможно, оценка формулируется как оптимизация метода взвешенных наименьших квадратов (WLS), которая минимизирует сумму квадратов взвешенных остаточных значений. В противном случае оценка формулируется как минимизация отрицания логарифма вероятности (NLL). Формулировка WLS часто сходится лучше, чем формулировка NLL, и SimBiology может использовать в своих интересах специализированные алгоритмы WLS, такие как алгоритм Levenberg-Marquardt, реализованный в lsqnonlin и lsqcurvefit. SimBiology использует WLS, когда существует одна ошибочная модель, которая является постоянной, пропорциональной, или экспоненциальной. SimBiology использует NLL, если у вас есть объединенная ошибочная модель или модель нескольких-ошибок, то есть, модель, имеющая ошибочную модель для каждого ответа.

sbiofit поддерживает различные методы оптимизации и передает в сформулированном WLS или выражении NLL к методу оптимизации, который минимизирует его. Для простоты каждое выражение, показанное ниже, принимает только одну ошибочную модель и один ответ. Если существуют множественные ответы, SimBiology берет сумму выражений, которые соответствуют ошибочным моделям данных ответов.

 Выражение, которое минимизируется
Метод взвешенных наименьших квадратов (WLS)Для постоянной ошибочной модели, iN(yifi)2
Для пропорциональной ошибочной модели, iN(yifi)2fi2/fgm2
Для экспоненциальной ошибочной модели, iN(lnyilnfi)2
Для числовых весов, iN(yifi)2wgm/wi
Отрицательная логарифмическая правдоподобность (NLL)Для объединенной ошибочной модели и модели нескольких-ошибок, iN(yifi)22σi2+iNln2πσi2

Переменные определяются следующим образом.

N

Количество экспериментальных наблюдений

yi

i th экспериментальное наблюдение

fi

Ожидаемое значение i th наблюдение

σi

Стандартное отклонение i th наблюдение.

  • Для постоянной ошибочной модели, σi=a

  • Для пропорциональной ошибочной модели, σi=b|fi|

  • Для объединенной ошибочной модели, σi=a+b|fi|

fgm

fgm=(iN|fi|)1N

wi

Вес i th ожидаемое значение

wgm

wgm=(iNwi)1N

Когда вы используете числовые веса или функцию веса, веса приняты, чтобы быть обратно пропорциональными отклонению ошибки, то есть, σi2=a2wi где a является постоянным параметром ошибок. Если вы используете веса, вы не можете задать ошибочную модель кроме постоянной ошибочной модели.

Различные методы оптимизации имеют различные требования к функции, которая минимизируется. Для некоторых методов оценка параметров модели выполняется независимо от оценки ошибочных параметров модели. Следующая таблица обобщает ошибочные модели и любые отдельные формулы, используемые для оценки ошибочных параметров модели, где a и b являются ошибочными параметрами модели, и e является стандартным средним нулем и переменной (Gaussian) модульного отклонения.

Ошибочная модельФункция оценки параметра ошибок
'constant': yi=fi+aea2=1NiN(yifi)2
'exponential': yi=fiexp(ae)a2=1NiN(lnyilnfi)2
'proportional': yi=fi+b|fi|eb2=1NiN(yififi)2
'combined': yi=fi+(a+b|fi|)eПараметры ошибок включены в минимизацию.
Весаa2=1NiN(yifi)2wi

Примечание

nlinfit только поддержите одну ошибочные модели, не модели нескольких-ошибок, то есть, специфичные для ответа ошибочные модели. Для объединенной ошибочной модели это использует итеративный алгоритм WLS. Для других ошибочных моделей это использует алгоритм WLS, аналогичный описанному ранее. Для получения дополнительной информации смотрите nlinfit (Statistics and Machine Learning Toolbox).

Опции по умолчанию для алгоритмов оценки

Следующая таблица обобщает опции по умолчанию для различных функций оценки.

ФункцияОпции по умолчанию
nlinfit (Statistics and Machine Learning Toolbox)

sbiofit использует структуру опций по умолчанию, сопоставленную с nlinfit, за исключением:

FunValCheck = 'off'
DerivStep = max(eps^(1/3), min(1e-4,SolverOptions.RelativeTolerance)), где SolverOptions свойство соответствует активной модели configset object.
fmincon (Optimization Toolbox)

sbiofit использует структуру опций по умолчанию, сопоставленную с fmincon, за исключением:

Display = 'off'
FunctionTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
OptimalityTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
Algorithm = 'trust-region-reflective' когда 'SensitivityAnalysis' верно, или 'interior-point' когда 'SensitivityAnalysis' является ложным.
FiniteDifferenceStepSize = max(eps^(1/3),min(1e-4,SolverOptions.RelativeTolerance)), где SolverOptions свойство соответствует активной модели configset object.
TypicalX = 1e-6*x0, где x0 массив преобразованных первоначальных оценок.
fminunc (Optimization Toolbox)

sbiofit использует структуру опций по умолчанию, сопоставленную с fminunc, за исключением:

Display = 'off'
FunctionTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
OptimalityTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
Algorithm = 'trust-region' когда 'SensitivityAnalysis' верно, или 'quasi-newton' когда 'SensitivityAnalysis' является ложным.
FiniteDifferenceStepSize = max(eps^(1/3),min(1e-4,SolverOptions.RelativeTolerance)), где SolverOptions свойство соответствует активной модели configset object.
TypicalX = 1e-6*x0, где x0 массив преобразованных первоначальных оценок.
fminsearch

sbiofit использует структуру опций по умолчанию, сопоставленную с fminsearch, за исключением:

Display = 'off'
TolFun = 1e-6*abs(f0), где f0 начальное значение целевой функции.
lsqcurvefit (Optimization Toolbox), lsqnonlin (Optimization Toolbox)

Требует Optimization Toolbox.

sbiofit использует структуру опций по умолчанию, сопоставленную с lsqcurvefit и lsqnonlin, за исключением:

Display = 'off'
FunctionTolerance = 1e-6*norm(f0), где f0 начальное значение целевой функции.
OptimalityTolerance = 1e-6*norm(f0), где f0 начальное значение целевой функции.
FiniteDifferenceStepSize = max(eps^(1/3),min(1e-4,SolverOptions.RelativeTolerance)) , где SolverOptions свойство соответствует активной модели configset object.
TypicalX = 1e-6*x0, где x0 массив преобразованных первоначальных оценок.
patternsearch (Global Optimization Toolbox)

Требует Global Optimization Toolbox.

sbiofit использует объект опций по умолчанию (optimoptions) сопоставленный с patternsearch, за исключением:

Display = 'off'
FunctionTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
MeshTolerance = 1.0e-3
AccelerateMesh = true
ga (Global Optimization Toolbox)

Требует Global Optimization Toolbox.

sbiofit использует объект опций по умолчанию (optimoptions) сопоставленный с ga, за исключением:

Display = 'off'
FunctionTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
MutationFcn = @mutationadaptfeasible
particleswarm (Global Optimization Toolbox)

Требует Global Optimization Toolbox.

sbiofit использует следующие опции по умолчанию для particleswarm алгоритм, за исключением:
Display = 'off'
FunctionTolerance = 1e-6*abs(f0), где f0 начальное значение целевой функции.
InitialSwarmSpan = 2000 or 8; 2000 для предполагаемых параметров без преобразования; 8 для предполагаемых параметров с log, logit, или probit преобразовывания.
scattersearchСмотрите поля точек алгоритм поиска.

Рассейте алгоритм поиска

scattersearch метод реализует глобальный алгоритм оптимизации [2], который обращается к некоторым проблемам оценки параметра динамические модели, такие как сходимость к локальным минимумам.

Обзор алгоритма

Алгоритм выбирает подмножество точек от начального пула точек. В том подмножестве некоторые точки являются лучшими по качеству (то есть, самое низкое значение функции), и некоторые случайным образом выбраны. Алгоритм итеративно оценивает точки и исследует различные направления вокруг различных решений найти лучшие решения. Во время этого шага итерации алгоритм заменяет любое старое решение на новое лучшего качества. Итерации продолжают, пока любому критерию остановки не соответствуют. Это затем запускает локальный решатель на лучшей точке.

Инициализация

Чтобы запустить поля точек поиск, алгоритм сначала решает общее количество необходимых точек (NumInitialPoints). По умолчанию общим количеством является 10*N, где N является количеством предполагаемых параметров. Это выбирает NumInitialPoints точки (строки) от InitialPointMatrix. Если InitialPointMatrix не имеет достаточного количества точек, алгоритм вызывает функцию, определяемую в CreationFcn сгенерировать дополнительные необходимые точки. По умолчанию латинская выборка гиперкуба используется, чтобы сгенерировать эти дополнительные точки. Алгоритм затем выбирает подмножество NumTrialPoints точки от NumInitialPoints 'points'. Часть (FractionInitialBest) из подмножества содержит лучшие точки в терминах качества. Остающиеся точки в подмножестве случайным образом выбраны.

Шаги итерации

Алгоритм выполняет итерации на точках в подмножестве можно следующим образом:

  1. Задайте гиперпрямоугольники вокруг каждой пары точек при помощи относительных качеств (то есть, значения функции) этих точек как мера смещения, чтобы создать эти прямоугольники.

  2. Оцените новое решение в каждом прямоугольнике. Если новое решение превосходит исходное решение по характеристикам, замените оригинал на новый.

  3. Примените движение - вне стратегии к улучшенным решениям и использованию многообещающие направления, чтобы найти лучшие решения.

  4. Запустите локальный поиск в каждом LocalSearchInterval итерация. Используйте LocalSelectBestProbability вероятность, чтобы выбрать лучшую точку как начальную точку для локального поиска. По умолчанию решение случайно, давая равному шансу выбрать лучшую точку или случайную точку от испытательных точек. Если новое решение превосходит старое решение по характеристикам, замените старый на новый.

  5. Замените любую остановленную точку, которая не производит нового решения для поражения после MaxStallTime секунды с другой точкой от начального набора.

  6. Оцените критерий остановки. Прекратите выполнять итерации, если каким-либо критериям соответствуют.

Алгоритм затем запускает локальный решатель на лучшей замеченной точке.

Критерий остановки

Алгоритм выполняет итерации, пока он не достигает останавливающегося критерия.

Остановка опцииОстановка теста
FunctionTolerance и MaxStallIterations

Относительное изменение в лучшем значении целевой функции по последнему MaxStallIterations меньше FunctionTolerance.

MaxIterations

Количество итераций достигает MaxIterations.

OutputFcn

OutputFcn может остановить итерации.

ObjectiveLimit

Лучшее значение целевой функции в итерации меньше чем или равно ObjectiveLimit.

MaxStallTime

Лучшее значение целевой функции не изменилось в последнем MaxStallTime секунды.

MaxTime

Функциональное время выполнения превышает MaxTime секунды.

Опции алгоритма

Вы создаете опции для алгоритма с помощью struct.

ОпцияОписание
CreationFcn

Обработайте к функции, которая создает дополнительные точки, необходимые для алгоритма. Значением по умолчанию является вектор символов 'auto', который использует латинскую выборку гиперкуба.

Функциональная подпись: points = CreationFcn(s,N,lb,ub), где s является общим количеством произведенных точек, N является количеством предполагаемых параметров, lb является нижней границей, и ub является верхней границей. Если какой-либо выход от функции превышает границы, эти результаты являются усеченными к границам.

Display

Level of display возвращен в командную строку.

  • 'off' или 'none' (значение по умолчанию) не отображает вывода.

  • 'iter' дает итеративное отображение.

  • 'final' отображения только окончательный результат.

FractionInitialBest

Числовой скаляр от 0 через 1. Значением по умолчанию является 0.5. Этот номер является частью NumTrialPoints это выбрано как лучшие точки из NumInitialPoints 'points'.

FunctionTolerance

Числовой скаляр от 0 через 1. Значением по умолчанию является 1e-6. Решатель останавливается если относительное изменение в лучшем значении целевой функции по последнему MaxStallIterations меньше FunctionTolerance. Эта опция также используется, чтобы удалить дублирующиеся локальные решения. Смотрите XTolerance для деталей.

InitialPointMatrix

Начальная буква (или частичный) набор точек. M-by-N действительная конечная матрица, где M является числом точек и N, является количеством предполагаемых параметров.

Если M < NumInitialPoints, затем scattersearch создает больше точек так, чтобы общим количеством строк был NumInitialPoints.

Если M > NumInitialPoints, затем scattersearch использует первый NumInitialPoints 'Строки' .

Значением по умолчанию являются преобразованные значения начальной буквы предполагаемых параметров, сохраненных в InitialTransformedValue свойство EstimatedInfo object, то есть, [estiminfo.InitialTransformedValue].

LocalOptions

Опции для локального решателя. Это может быть struct (созданный с optimset или statset (Statistics and Machine Learning Toolbox)) или optimoptions Объект (Optimization Toolbox), в зависимости от локального решателя. Значением по умолчанию является вектор символов 'auto', который использует опции по умолчанию выбранного решателя за некоторыми исключениями. В дополнение к этим исключениям следующие опции ограничивают время, проведенное в локальном решателе, потому что это неоднократно называется:

  • MaxFunEvals (максимальное количество позволенных вычислений функции) = 300

  • MaxIter (максимальное количество позволенных итераций) = 200

LocalSearchInterval

Положительное целое число. Значением по умолчанию является 10. scattersearch алгоритм применяет локальный решатель к одной из испытательных точек после первой итерации и снова каждого LocalSearchInterval итерация.

LocalSelectBestProbability

Числовой скаляр от 0 через 1. Значением по умолчанию является 0.5. Это - вероятность выбора лучшей точки как начальная точка для локального поиска. В других случаях одна из испытательных точек выбрана наугад.

LocalSolver

Вектор символов или строка, задающая имя локального решателя. Поддерживаемыми методами является 'fminsearch', 'lsqnonlin', 'lsqcurvefit', 'fmincon', 'fminunc', 'nlinfit'.

Локальный решатель по умолчанию выбран со следующим приоритетом:

  • Если Optimization Toolbox доступен, решатель lsqnonlin.

  • Если Statistics and Machine Learning Toolbox доступен, решатель nlinfit.

  • В противном случае решатель fminsearch.

MaxIterations

Положительное целое число. Значением по умолчанию является вектор символов 'auto' представление 20*N, где N является количеством предполагаемых параметров.

MaxStallIterations

Положительное целое число. Значением по умолчанию является 50. Решатель останавливается если относительное изменение в лучшем значении целевой функции по последнему MaxStallIterations итерации меньше FunctionTolerance.

MaxStallTime

Положительная скалярная величина. Значением по умолчанию является Inf. Решатель останавливается если MaxStallTime секунды передали начиная с последнего улучшения лучше всего замеченного значения целевой функции. Здесь, время является стенкой, показывают время в противоположность циклам процессора.

MaxTime

Положительная скалярная величина. Значением по умолчанию является Inf. Решатель останавливается если MaxTime секунды передали с начала поиска. Время здесь означает, что стенка показывает время в противоположность циклам процессора.

NumInitialPoints

Положительное целое число, которое является >= NumTrialPoints. Решатель генерирует NumInitialPoints точки прежде, чем выбрать подмножество испытательных точек (NumTrialPoints) для последующих шагов. Значением по умолчанию является вектор символов 'auto', который представляет 10*N, где N является количеством предполагаемых параметров.

NumTrialPoints

Положительное целое число, которое является >= 2 и <= NumInitialPoints. Решатель генерирует NumInitialPoints начальные точки прежде, чем выбрать подмножество испытательных точек (NumTrialPoints) для последующих шагов. Значением по умолчанию является вектор символов 'auto', который представляет первое четное число n для который n2n10*N, где N является количеством предполагаемых параметров.

ObjectiveLimit

Скаляр. Значением по умолчанию является -Inf. Решатель останавливается, если лучшее значение целевой функции в итерации меньше чем или равно ObjectiveLimit.

OutputFcn

Указатель на функцию или cell-массив указателей на функцию. Выходные функции могут считать итеративные данные и остановить решатель. Значением по умолчанию является [].

Подписью выходной функции является stop = myfun(optimValues,state), где:

  • stop логический скаляр. Установите на true остановить решатель.

  • optimValues структура, содержащая информацию об испытательных точках с полями.

    • bestx найденная точка лучшего решения, соответствуя значению функции bestfval.

    • bestfval лучшее (самое низкое) найденное значение целевой функции.

    • iteration номер итерации.

    • medianfval среднее значение целевой функции среди всех текущих испытательных точек.

    • stalliterations количество итераций начиная с последнего изменения в bestfval.

    • trialx матрица текущих испытательных точек. Каждая строка представляет одну точку, и количество строк равно NumTrialPoints.

    • trialfvals вектор из значений целевой функции для испытательных точек. Это - матрица для lsqcurvefit и lsqnonlin методы.

  • state вектор символов, дающий состояние текущей итерации.

    • 'init' – Решатель не начал выполнять итерации. Ваша выходная функция может использовать это состояние, чтобы открыть файлы или настроить структуры данных или графики для последующих итераций.

    • 'iter' – Решатель является продолжением его итераций. Как правило, это состояние - то, где ваша выходная функция выполняет свою работу.

    • 'done' – Решатель достигает останавливающегося критерия. Ваша выходная функция может использовать это состояние, чтобы вымыться, такие как закрытие любых файлов, которые это открыло.

TrialStallLimit

Положительное целое число, со значением по умолчанию 22. Если конкретная испытательная точка не улучшается после TrialStallLimit итерации, это заменяется другой точкой.

UseParallel

Логический флаг, чтобы вычислить целевую функцию параллельно. Значением по умолчанию является false.

XTolerance

Числовой скаляр от 0 через 1. Значением по умолчанию является 1e-6. Эта опция задает, как близко две точки должны быть должны считать их идентичными для создания вектора из локальных решений. Решатель вычисляет расстояние между парой точек с norm, Евклидово расстояние. Если два решения в XTolerance расстояние друг друга и имеет значения целевой функции в FunctionTolerance друг из друга решатель считает их идентичными. Если оба условия не соблюдают, решатель сообщает о решениях как отличных.

Получить отчет каждого потенциального локального минимального, набора XTolerance к 0. Чтобы получить отчет меньшего количества результатов, установите XTolerance к большему значению.

Несколько оценок параметра параллельно

Существует два способа использовать параллельные вычисления для оценки параметра.

Установите 'UseParallel' к истине

Включить распараллеливание для sbiofit, установите пару "имя-значение" 'UseParallel' к true. Функция поддерживает несколько уровней распараллеливания, но только один уровень используется за один раз. Для необъединенного пригодного для нескольких групп (или индивидуумы), каждая подгонка запускается параллельно. Для объединенной подгонки распараллеливание происходит на уровне решателя, если решатель поддерживает его. Таким образом, sbiofit устанавливает параллельную опцию соответствующего метода оценки (решатель) к истине, и вычисления функции возражения выполняются параллельно. Например, градиенты вычисляются параллельно для основанных на градиенте методов. Если вы выполняете объединенную подгонку с несколькими группами, использующими решатель, который не имеет параллельной опции, симуляции запущенный параллельно для каждой группы во время оптимизации (оценка наибольшего правдоподобия).

Использование parfeval или parfor

Можно также вызвать sbiofit внутри a parfor цикл или использование a parfeval внутри a for- цикл, чтобы выполнить несколько оценок параметра параллельно. Рекомендуется, чтобы вы использовали parfeval потому что эти параллельные оценки, запущенные асинхронно. Если одна подгонка производит ошибку, она не влияет на другие подгонки.

При попытке найти глобальный минимум, можно использовать глобальные решатели, такой как particleswarm (Global Optimization Toolbox) или ga (Global Optimization Toolbox) (Global Optimization Toolbox требуется). Однако, если вы хотите задать начальные условия и запустить подгонки параллельно, смотрите следующий пример, который показывает, как использовать обоих parfor и parfeval.

Модель и Setup данных

Загрузите модель белка G.

sbioloadproject gprotein

Храните экспериментальные данные, содержащие курс времени для части активного белка G [1].

time = [0 10 30 60 110 210 300 450 600]';
GaFracExpt = [0 0.35 0.4 0.36 0.39 0.33 0.24 0.17 0.2]';

Создайте a groupedData object на основе экспериментальных данных.

tbl = table(time,GaFracExpt);
grpData = groupedData(tbl);

Сопоставьте соответствующий элемент модели с экспериментальными данными.

responseMap = 'GaFrac = GaFracExpt';

Задайте параметр, чтобы оценить.

paramToEstimate    = {'kGd'};

Сгенерируйте начальные значения параметров для kGd.

rng('default');
iniVal = abs(normrnd(0.01,1,10,1));
fitResultPar = [];

Параллельный Setup пула

Запустите параллельный пул с помощью локального профиля.

poolObj = parpool('local');
Starting parallel pool (parpool) using the 'local' profile ...
Connected to the parallel pool (number of workers: 6).

Используя parfeval (Рекомендуемый)

Во-первых, задайте указатель на функцию, который использует локальную функцию sbiofitpar для оценки. Убедитесь функциональный sbiofitpar задан в конце скрипта.

optimfun = @(x) sbiofitpar(m1,grpData,responseMap,x);

Выполните несколько оценок параметра параллельно через parfeval использование различных начальных значений параметров.

for i=1:length(iniVal)
    f(i) = parfeval(optimfun,1,iniVal(i)); 
end
fitResultPar = fetchOutputs(f);

Обобщите результаты для каждого запуска.

allParValues = vertcat(fitResultPar.ParameterEstimates);
allParValues.LogLikelihood = [fitResultPar.LogLikelihood]';
allParValues.RunNumber = (1:length(iniVal))';
allParValues.Name = categorical(allParValues.Name);
allParValues.InitialValue = iniVal;
% Rearrange the columns.
allParValues = allParValues(:,[5 1 6 2 3 4]);
% Sort rows by LogLikelihood.
sortrows(allParValues,'LogLikelihood')
ans=10×6 table
    RunNumber    Name    InitialValue    Estimate    StandardError    LogLikelihood
    _________    ____    ____________    ________    _____________    _____________

        9        kGd        3.5884         3.022           0.127         -1.2843   
       10        kGd        2.7794         2.779        0.029701         -1.2319   
        3        kGd        2.2488        2.2488        0.096013         -1.0786   
        2        kGd        1.8439         1.844         0.28825        -0.90104   
        6        kGd        1.2977        1.2977        0.011344        -0.48209   
        4        kGd       0.87217       0.65951        0.003583          0.9279   
        1        kGd       0.54767       0.54776       0.0020424          1.5323   
        7        kGd       0.42359       0.42363       0.0024555          2.6097   
        8        kGd       0.35262       0.35291      0.00065289          3.6098   
        5        kGd       0.32877       0.32877      0.00042474          4.0604   

Задайте локальную функцию sbiofitpar это выполняет использование оценки параметра sbiofit.

function fitresult = sbiofitpar(model,grpData,responseMap,initialValue)
estimatedParam = estimatedInfo('kGd');
estimatedParam.InitialValue = initialValue;
fitresult = sbiofit(model,grpData,responseMap,estimatedParam);
end

Используя parfor

В качестве альтернативы можно выполнить несколько оценок параметра параллельно через parfor цикл.

parfor i=1:length(iniVal)
    estimatedParam = estimatedInfo(paramToEstimate,'InitialValue',iniVal(i));
    fitResultTemp = sbiofit(m1,grpData,responseMap,estimatedParam); 
    fitResultPar = [fitResultPar;fitResultTemp];
end

Закройте параллельный пул.

delete(poolObj);

Оценка параметра с гибридными решателями

sbiofit глобальные методы оптимизации поддержек, а именно, ga (Global Optimization Toolbox) и particleswarm (Global Optimization Toolbox) (требуемый Global Optimization Toolbox). Чтобы улучшить результаты оптимизации, эти методы позволяют вам запустить гибридную функцию после того, как глобальный решатель остановится. Гибридная функция использует конечную точку, возвращенную глобальным решателем как его начальная точка. Поддерживаемые гибридные функции:

Убедитесь, что ваша гибридная функция принимает ваши ограничения задач. Таким образом, если ваши параметры ограничены, используйте соответствующую функцию (такой как fmincon или patternsearch) для ограниченной оптимизации. Если не ограниченный, использовать fminunc, fminsearch, или patternsearch. В противном случае, sbiofit выдает ошибку.

Для проиллюстрированного примера смотрите, Выполняют Гибридную Оптимизацию Используя sbiofit.

Ссылки

[1] И, T-M., Kitano, H. и Саймон, M. (2003). Количественная характеристика дрожжей гетеротримерный цикл белка G. PNAS. 100, 10764–10769.

[2] Габор, A., и Банга, J.R. (2015). Устойчивая и эффективная оценка параметра динамические модели биологических систем. Системная биология BMC. 9:74.

Расширенные возможности

Введенный в R2014a