Определение модельной функции нелинейной модели смешанных эффектов как системы ОДУ в Matlab nlmefit - PullRequest
0 голосов
/ 02 сентября 2018

Я изучаю нелинейные модели смешанных эффектов и пытаюсь оценить параметры системы ODE, используя это. Я пробовал пример Matlab для нелинейных смешанных эффектов на с использованием nlmefit

Основываясь на примере на этой странице, я хочу изменить код, чтобы можно было найти параметры модели следующей системы ($ \ lambda, \ rho, \ eta, \ delta, c $).
enter image description here

(В качестве попытки я подхожу к тем же данным, найденным в примере Matlab от load indomethacin. Итак, данные состоят из того, как изменяется концентрация лекарства у 6 субъектов, измеренных за период времени 8 часов, и в моем коде Я предполагаю, что concentration дается (T_u + T_v). Я хочу приспособить модель, решая систему, используя ode45, а не аналитическое решение T_u + T_v.

Это код, который я написал,

clear
close all
load indomethacin

nlme_model = @(PHI,t)virusModel(PHI,time(1:11),[0.05,0,0]);

phi0 = [0.09 0.06212 0.02588 0.08414 0.02];

[phi,PSI,stats] = nlmefit(time,concentration,subject, ...
[],nlme_model,phi0)

PHI = repmat(phi,1,6) + ... % Fixed effects
[b(1,:);b(2,:);b(3,:);b(4,:);b(5,:)]; % Random effects
RES = zeros(11,6); % Residuals
colors = 'rygcbm';
%plotting the results
for I = 1:6

fitted_model = virusModel([PHI(1,I),PHI(2,I),PHI(3,I),PHI(4,I),PHI(5,I)],time(1:11),[0.05,0,0]);
tI = time(subject == I);
cI = concentration(subject == I);
RES(:,I) = cI - fitted_model(tI);
subplot(2,3,I)
scatter(tI,cI,20,colors(I),'filled')
hold on
plot(tI,fitted_model(tI),'Color',colors(I))
axis([0 8 0 3.5])
xlabel('Time (hours)')
ylabel('Concentration (mcg/ml)')
legend(num2str(I),'Subject')
end
figure(2)
normplot(RES(:))  

и функция модели,

function output= virusModel(PHI,thours,initialCond)

   [~,values] = ode45(@virusEquations,thours,initialCond);

    function s=virusEquations(t,y)

    s=zeros(3,1);

    s(1)=PHI(:,1)-PHI(:,2).*y(1)-PHI(:,3).*y(1).*y(3);%differential equaiton for T_u
    s(2)=PHI(:,3).*y(1)*y(3)-PHI(:,4).*y(2);%differential equaiton for T_i
    s(3)=PHI(:,4).*y(2)-PHI(:,5)*y(3);%differential equaiton for V

    end
output=values(:,1)+values(:,2);
end

Я не уверен, что это написано правильно, так как я получаю ошибку как

Error using nlmefit>LMfit (line 1848)
After the initial refinement of the fixed effects with the Levenberg-Marquardt algorithm some columns of the Jacobian are
effectively zero at BETA0, indicating that the model is insensitive to some of its parameters.  That may be because those
parameters are not present in the model, or otherwise do not affect the predicted values.  It may also be due to numerical
underflow in the model function, which can sometimes be avoided by choosing better initial parameter values, or by rescaling or
recentering.  
  1. Может кто-нибудь сообщить мне, правильно ли я вызываю систему ode и правильно ли реализую смешанные эффекты?

  2. Кроме того, я хочу знать, как я могу передать начальные условия Система Оде. Можно ли указать конкретные начальные условия для каждого из субъектов?

...