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

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

load pheno.mat ds

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

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

data = groupedData(ds);

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

data(1:5,:)

ans =

  5×6 groupedData

    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 = 4×1 cell
    {'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.