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

Этот пример использует данные, собранные по 59 недоношенным детям, данным фенобарбитал в течение первых 16 дней после рождения [1]. ds является таблицей, содержащей разовые концентрацией данные о профиле и ковариационную информацию для каждого младенца (или группа).

load pheno.mat ds

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

Преобразуйте данные в groupedData формат для оценки параметра.

data = groupedData(ds);

Отобразите первые несколько строк data.

data(1:5,:)

ans=5×6 table
    ID    TIME    DOSE    WEIGHT    APGAR    CONC
    __    ____    ____    ______    _____    ____

    1        0     25      1.4        7       NaN
    1        2    NaN      1.4        7      17.3
    1     12.5    3.5      1.4        7       NaN
    1     24.5    3.5      1.4        7       NaN
    1       37    3.5      1.4        7       NaN

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

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

t = sbiotrellis(data, 'ID', 'TIME', 'CONC', 'marker', 'o',...
       'markerfacecolor', [.7 .7 .7], 'markeredgecolor', 'r', ...
       'linestyle', 'none');
t.plottitle = 'Concentration versus Time';

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

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

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

Сопоставьте образцовые разновидности с данными об ответе.

responseMap = 'Drug_Central = CONC';

Задайте предполагаемые параметры

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

и

где $$\theta$$, $$eta$$, и $$\varphi$$ фиксированные эффекты, случайные эффекты и оцененные значения параметров соответственно, вычисленный для каждого младенца (группа) $$i$$. Некоторые произвольные первоначальные оценки для V (объем центрального отсека) и Статья (уровень раскрываемости преступлений) используются здесь в отсутствие лучших эмпирических данных.

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

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

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

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

Соответствуйте модели

Используйте sbiofitmixed, чтобы соответствовать модели с одним отсеком к данным.

nlmeResults = sbiofitmixed(onecomp,data,responseMap,estimatedParams,doses,'nlmefit');

Визуализация результатов

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

fig1 = plot(nlmeResults,'ParameterType','individual');
% Resize the figure.
fig1.hFig.Position(:) = [100 100 1200 800];

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

Предположим, что вы хотите предсказать, как младенцы 1 и 2 ответили бы под различными объемами дозирования. Можно предсказать их ответы можно следующим образом.

Создайте новые объекты дозы с новыми суммами дозы.

dose1 = doses(1);
dose1.Amount = dose1.Amount*2;
dose2 = doses(2);
dose2.Amount = dose2.Amount*1.5;

Используйте функцию predict, чтобы оценить подобранную модель с помощью новых данных о дозировании. Если вы хотите прогнозы ответа в конкретные времена, обеспечьте новый выходной временной вектор. Используйте опцию 'ParameterType', чтобы задать человека или параметры генеральной совокупности, чтобы использовать. По умолчанию predict использует параметры генеральной совокупности, когда вы задаете выходные времена.

timeVec = [0:25:400];
newResults = predict(nlmeResults,timeVec,[dose1;dose2],'ParameterType','population');

Визуализируйте предсказанные ответы при наложении экспериментальных данных для младенцев 1 и 2.

subplot(2,1,1)
plot(data.TIME(data.ID == 1),data.CONC(data.ID == 1),'bo')
hold on
plot(newResults(1).Time,newResults(1).Data,'b')
hold off
ylabel('Concentration')
legend('Observation(CONC)','Prediction')
subplot(2,1,2)
plot(data.TIME(data.ID == 2),data.CONC(data.ID == 2),'rx')
hold on
plot(newResults(2).Time,newResults(2).Data,'r')
hold off
legend('Observation(CONC)','Prediction')
ylabel('Concentration')
xlabel('Time')

Создайте ковариационную модель для ковариационных зависимостей

Предположим, что существует корреляция между объемом и весом, и возможно счетом APGAR и объемом. Рассмотрите эффект веса путем моделирования двух из этих ковариационных зависимостей: объем центральных (Central) и уровень раскрываемости преступлений (Cl_Central) меняется в зависимости от веса. Модель становится

и

Используйте объект CovariateModel задать ковариационные зависимости. Для получения дополнительной информации смотрите, Задают Ковариационную Модель.

covModel = CovariateModel;
covModel.Expression = ({'Central = exp(theta1 + theta2*WEIGHT + eta1)',...
                        'Cl_Central = exp(theta3 + theta4*WEIGHT + eta2)'});
% Use |constructDefaultInitialEstimate| to create an |initialEstimates|
% struct.
initialEstimates = covModel.constructDefaultFixedEffectValues;

Используйте свойство FixedEffectNames отобразить thetas (зафиксированные эффекты) заданный в модели.

covModel.FixedEffectNames

ans = 4x1 cell array
    {'theta1'}
    {'theta3'}
    {'theta2'}
    {'theta4'}

Используйте свойство FixedEffectDescription показать описания соответствующих фиксированных эффектов (thetas) используемый в ковариационном выражении. Например, theta2 является фиксированным эффектом для коварианта веса, который коррелирует с объемом (Central), обозначенный как 'Центральный / ВЕС'.

disp('Fixed Effects Description:');

Fixed Effects Description:

disp(covModel.FixedEffectDescription);

    'Central'
    'Cl_Central'
    'Central/WEIGHT'
    'Cl_Central/WEIGHT'

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

initialEstimates.theta1 = nlmeResults.FixedEffects.Estimate(1);
initialEstimates.theta3 = nlmeResults.FixedEffects.Estimate(2);
covModel.FixedEffectValues = initialEstimates;

Соответствуйте модели

nlmeResults_cov = sbiofitmixed(onecomp,data,responseMap,covModel,doses,'nlmefit');

Отобразите подходящие параметры и ковариации

disp('Estimated Fixed Effects:');

Estimated Fixed Effects:

disp(nlmeResults_cov.FixedEffects);

      Name          Description        Estimate    StandardError
    ________    ___________________    ________    _____________

    'theta1'    'Central'              -0.45664      0.078933   
    'theta3'    'Cl_Central'            -5.9519        0.1177   
    'theta2'    'Central/WEIGHT'        0.52948      0.047342   
    'theta4'    'Cl_Central/WEIGHT'     0.61954      0.071386   

disp('Estimated Covariance Matrix:');

Estimated Covariance Matrix:

disp(nlmeResults_cov.RandomEffectCovarianceMatrix);

              eta1        eta2  
            ________    ________

    eta1    0.046503           0
    eta2           0    0.041609

Визуализация результатов

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

fig2 = plot(nlmeResults_cov,'ParameterType','individual');
% Resize the figure.
fig2.hFig.Position(:) = [100 100 1200 800];

Используйте новые ковариационные данные, чтобы оценить подобранную модель

Предположим, что вы хотите исследовать ответы младенцев 1 и 2 использующих различных ковариационных данных, а именно, WEIGHT. Можно сделать это путем определения новых данных WEIGHT. Переменная ID данных соответствует отдельным младенцам.

newData = data(data.ID == 1 | data.ID == 2,:);
newData.WEIGHT(newData.ID == 1) = 1.3;
newData.WEIGHT(newData.ID == 2) = 1.4;

Моделируйте ответы младенцев 1 и 2 использования новых ковариационных данных.

[newResults_cov, newEstimates] = predict(nlmeResults_cov,newData,[dose1;dose2]);

newEstimates содержит обновленные оценки параметра для каждого человека (младенцы 1 и 2) после того, как модель будет переоценена с помощью новых ковариационных данных.

newEstimates

newEstimates=4×3 table
    Group        Name        Estimate 
    _____    ____________    _________

      1      'Central'          2.5596
      1      'Cl_Central'    0.0065965
      2      'Central'          1.7123
      2      'Cl_Central'    0.0064806

Сравните с ориентировочными стоимостями от исходной подгонки с помощью старых ковариационных данных.

nlmeResults_cov.IndividualParameterEstimates( ...
            nlmeResults_cov.IndividualParameterEstimates.Group == '1' | ...
            nlmeResults_cov.IndividualParameterEstimates.Group == '2',:)

ans=4×3 table
    Group        Name        Estimate 
    _____    ____________    _________

      1      'Central'          2.6988
      1      'Cl_Central'    0.0070181
      2      'Central'          1.8054
      2      'Cl_Central'    0.0068948

Визуализируйте новые результаты симуляции вместе с экспериментальными данными для младенца 1 и 2.

subplot(2,1,1)
plot(data.TIME(data.ID == 1),data.CONC(data.ID == 1),'bo')
hold on
plot(newResults_cov(1).Time,newResults_cov(1).Data,'b')
hold off
ylabel('Concentration')
legend('Observation(CONC)','Prediction','Location','NorthEastOutside')
subplot(2,1,2)
plot(data.TIME(data.ID == 2),data.CONC(data.ID == 2),'rx')
hold on
plot(newResults_cov(2).Time,newResults_cov(2).Data,'r')
hold off
legend('Observation(CONC)','Prediction','Location','NorthEastOutside')
ylabel('Concentration')
xlabel('Time')

Ссылки

[1] Grasela, T. H. Младший, и С. М. Донн. "Неонатальная фармакокинетика генеральной совокупности фенобарбитала выведена от стандартных клинических данных". Фармакол Dev Там 1985:8 (6). 374-83. Краткий обзор PubMed