WICHTIG: Der Betrieb von goMatlab.de wird privat finanziert fortgesetzt. - Mehr Infos...

Mein MATLAB Forum - goMatlab.de

Mein MATLAB Forum

 
Gast > Registrieren       Autologin?   

Partner:




Forum
      Option
[Erweitert]
  • Diese Seite per Mail weiterempfehlen
     


Gehe zu:  
Neues Thema eröffnen Neue Antwort erstellen

Global Search mit fmincon

 

maki
Forum-Anfänger

Forum-Anfänger


Beiträge: 32
Anmeldedatum: 23.03.20
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 28.05.2020, 11:36     Titel: Global Search mit fmincon
  Antworten mit Zitat      
Hi,
das Skript
Code:
gs = GlobalSearch('Display','iter');
opts=optimset ('MaxTime',1800,'UseParallel',true)
problem= createOptimProblem('fmincon','objective',@Simulation_Time_Courses_BI130,'x0',[0.02,0.02,0.001],...
    'lb',[0.015,0.015,0.00025],...
    'ub',[0.03,0.03,0.004],'options',opts);

x=run(gs,problem)
 


greift auf die Funktion

%% 2. Ebene
%% ODE + Sum of squares
function Sum_of_squares=Simulation_Time_Courses_BI130_myfor(ySpeicher,Exp_Data_VCD,Exp_Data_TCD)
global Data
Data= xlsread('BI130_Average_Daten.xlsx','D3:S30');
Exp_Data = Data (1:28,1:4);
Exp_Data= round(Exp_Data,2);
%% Start values from Experimental Data
c0_global(1) = Data(1,2);  
c0_global(2) = Data(1,4);  
c0_global(3) = Data(1,6);  
c0_global(4) = Data(1,8);  
c0_global(5) = Data(1,10);
c0_global(6) = Data(1,12);
%% Constants
Constants(1) = 0.01;      
Constants(2) = 4*10^-12;  
Constants(3) = 0.4;        
Constants(4) = 1.6;        
Constants(5) = 0;          
Constants(6) = 0.5;        
Constants(7) = 0.03;          
Constants(8) = 0.03;      
Constants(9) = 0.014;      
Constants(10) =0.1;        
Constants(11) = 0.45;    
Constants(12) = 0.085;  
Constants(13) = 0.4;      
Constants(14) = 0.19;      
Constants(15) = 0.3;      
%% Parameters
Parameters(1) =0.055;      
Parameters(2) =0.002;      
Parameters(3) =0.03;      
%% Initial Values
Xv0   = [0.4139*10^6 ,0.3004*10^6    ,0.4169*10^6    ,0.3813*10^6    ,0.3918*10^6    ,0.3914*10^6    ,0.3779*10^6   ,0.4804*10^6    ,0.6151*10^6  ,0.7395*10^6];
Xt0   = [0.4351*10^6 ,0.3372*10^6    ,0.4384*10^6    ,0.3955*10^6    ,0.4067*10^6    ,0.4074*10^6    ,0.3853*10^6   ,0.4876*10^6    ,0.63*10^6  ,0.7521*10^6];
cGlc0 = [6.8      ,6.8         ,6.8         ,6.8         ,6.8         ,6.8         ,6.8         ,6.8         ,6.8        ,6.8];    
cLac0 = [0.3      ,0.3         ,0.3         ,0.3         ,0.3         ,0.3         ,0.3         ,0.3         ,0.3        ,0.3];      
cGln0 = [0.75     ,0.75        ,0.75        ,0.75        ,0.75        ,0.75        ,0.75        ,0.75        ,0.75       ,0.75];      
cLs0  = [1        ,1           ,1           ,1           ,1           ,1           ,1           ,1           ,1          ,1];    
cAmm0 = [1        ,1           ,1           ,1           ,1           ,1           ,1           ,1           ,1          ,1];    
v0    = [0.25     ,0.5         ,1           ,2           ,2           ,10          ,50          ,80          ,400        ,2000];
%%  Time Cell
p =cell (1,10);
p{1}= [0 70.8];
p{2}= [71.28 143.52];
p{3}= [144 216.24];
p{4}= [217.2 287.04];
p{5}= [288.24 358.32];
p{6}= [358.42 429.22];
p{7}= [429.32 474.776 545.192];
p{8}= [545.292 567.991 615.391 684.631];
p{9}= [684.731 707.675 754.876 825.052];
p{10}= [825.152 847.448 894.348 965.532 1059.852];

%% For Loop
n=10;    
t=0;
y=[0,0,0,0,0,0,0,0];

for k = 1:n
 p{k}= round (p{k},2);  
   c0 = [ Xv0(k); Xt0(k); cGlc0(k); cLac0(k);  cGln0(k); cLs0(k); cAmm0(k); v0(k)];
    [tSpeicher,ySpeicher]=ode45(@Model_BI130_myfor,p{k},c0,'',Parameters, Constants);
   
    ySpeicher =ySpeicher(ismember(tSpeicher,p{k}),:);  
    tSpeicher =tSpeicher(ismember(tSpeicher,p{k}));
   
     t = vertcat(t,tSpeicher);
     y= vertcat(y,ySpeicher);

Exp_Data_VCD =Exp_Data(ismember(tSpeicher(:,1),p{k}),2)*1E6;  
Exp_Data_TCD =Exp_Data(ismember(tSpeicher(:,1),p{k}),4)*1E6;    
     
     %% Sum of squares
Magnitude(1)=1E6;
Magnitude(2)=1;
%% Sum of squares
Sum_of_squares1= sum((ySpeicher(:,1))-Exp_Data_VCD).^2/(mean(Exp_Data_VCD)).^2 %VCD
   
Sum_of_squares2= sum((ySpeicher(:,2))-Exp_Data_TCD).^2/(mean(Exp_Data_TCD)).^2 %TCD
 
Sum_of_squares= Sum_of_squares1+ Sum_of_squares2
end
end

 


welches als Unterfunktion

Code:
%% 3. Ebene
%% Kinetiken und das eigentliche Model
function [dcdt] = Model_BI130_myfor(t,c,Parameters,Constants)
%% Parameters()
my_max   = Parameters(1);
my_d_max = Parameters(2);
my_d_min = Parameters(3);
%% constants()
Klys            = Constants(1);
qAmm_uptake_max = Constants(2);
YAmm_Gln        = Constants(3);
YLac_Glc        = Constants(4);
Fsample         = Constants(5);
kAmm            = Constants(6);
KsGlc           = Constants(7);
KsGln           = Constants(8);
qLs_max         = Constants(9);
KsLs            = Constants(10);
qGlc_max        = Constants(11);
qGln_max        = Constants(12);
qLac_uptake_max = Constants(13);
kGlc            = Constants(14);
kGln            = Constants(15);
%% Bennenung der Variablen
Xv   = c(1); %viable cell density
Xt   = c(2); %total cell density
cGlc = c(3); % Glucose concentration
cLac = c(4); % Lactate concentration
cGln = c(5); % Glutamin concentration
cLs  = c(6); % Limiting substrate concentration
cAmm = c(7); % Ammonia concentration
v    = c(8); % Volume
%% Parameter Berechnung
my = my_max  * (cLs/(cLs + KsLs));                                                            
my_d = my_d_min + my_d_max *(cLs/(cLs + KsLs));                                              
qGlc = qGlc_max * (cGlc/(cGlc + kGlc)) * (my/(my + my_max)+0.5);                              
qGln = qGln_max * (cGln/(cGln + kGln));                                                    
qLac = YLac_Glc * qGlc * (cGlc/cLac) - qLac_uptake_max * ((my_max - my)/my_max);            
qAmm = YAmm_Gln * qGln * cGln/cAmm - kAmm * qAmm_uptake_max * ((my_max - my)/my_max);              
qLs  = qLs_max * (cLs/(cLs+KsLs));                                                            
%% Differential Equations
dXvdt   = Xv * (my - my_d);              
dXtdt   = Xv * my - Klys * (Xt - Xv);      
dcGlcdt = -1.0 * Xv * qGlc;              
dcLacdt  = Xv * qLac;                
dcLsdt  = -1.0 * Xv * qLs;              
dcGlndt = -1.0 * Xv * qGln;              
dcAmmdt = Xv * qAmm;              
dvdt    = Fsample;                      


[dcdt]=[dXvdt; dXtdt; dcGlcdt; dcLacdt; dcGlndt;dcLsdt; dcAmmdt; dvdt];
end


hat. Laut Profiler verbringt das Skript die meiste Zeit in der ode45. Wie kann ich die Zeit reduzieren?

Viele Grüße,
Martin
Private Nachricht senden Benutzer-Profile anzeigen


Harald
Forum-Meister

Forum-Meister


Beiträge: 24.484
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 28.05.2020, 12:20     Titel:
  Antworten mit Zitat      
Hallo,

wirklich in ode45 oder in den Funktionsauswertungen von Model_BI130_myfor, die daraus aufgerufen werden?

So oder so würde ich an der Stelle mal das Generieren von MEX-Dateien mit MATLAB Coder versuchen.
https://www.mathworks.com/help/code.....er-project-interface.html
Unter "Generate C-Code" ist der entscheidende Punkt dann "Instead of generating a C static library, you can choose to generate a MEX function".

Grüße,
Harald
_________________

1.) Ask MATLAB Documentation
2.) Search gomatlab.de, google.de or MATLAB Answers
3.) Ask Technical Support of MathWorks
4.) Go mad, your problem is unsolvable ;)
Private Nachricht senden Benutzer-Profile anzeigen
 
maki
Themenstarter

Forum-Anfänger

Forum-Anfänger


Beiträge: 32
Anmeldedatum: 23.03.20
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 28.05.2020, 13:22     Titel:
  Antworten mit Zitat      
Hallo Harald,

danke für den Tipp, werde ich versuchen. Der Code lief von Mittwoch auf Donnerstag über Nacht und war immer noch nicht fertig.

Noch eine generelle Frage zum Code:

Wie funktioniert die Kommunikation zwischen GlobalSearch und der Funktion? Woher weiß GlobalSearch, dass es mir die 3 Werte "Parameter()" optimieren soll?

Viele Grüße,
Martin
Private Nachricht senden Benutzer-Profile anzeigen
 
Harald
Forum-Meister

Forum-Meister


Beiträge: 24.484
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 28.05.2020, 13:34     Titel:
  Antworten mit Zitat      
Hallo,

Zitat:
Woher weiß GlobalSearch, dass es mir die 3 Werte "Parameter()" optimieren soll?

Das weiß GlobalSearch gar nicht. Durch das Festlegen von x0 weiß GlobalSearch lediglich, dass 3 Werte optimiert werden sollen.
Da du Parameters in der Funktion Simulation_Time_Courses_BI130_myfor festlegst, bleibt Parameters immer gleich. Die Eingabeargumente dieser Funktion werden nicht verwendet, und das zeigt der Editor auch durch das orange Hinterlegen an.

Sinnvollerweise muss die Funktionsschnittstelle so aussehen:
Code:
function Sum_of_squares=Simulation_Time_Courses_BI130_myfor(Parameters)
 

und die Zeilen 30-32 müssen auskommentiert werden.

Insgesamt, gerade wenn ich mir Data ansehe (ist einerseits global, wird aber andererseits belegt), kann ich nur empfehlen, sich mal seeeehr gründlich mit dem Thema Übergabe von Argumenten auseinanderzusetzen.

Grundsätzlich auch: GlobalSearch macht mit dieser Zielfunktion nichts sinnvolles, weil fmincon schon nichts sinnvolles damit macht. Erst mal schauen, dass fmincon vernünftig läuft (d.h. insbesondere vom Startwert weggeht und einen besseren Zielfunktionswert al für den Startwert liefert). Davor ist es nicht sinnvoll, GlobalSearch zu verwenden.

Grüße,
Harald
_________________

1.) Ask MATLAB Documentation
2.) Search gomatlab.de, google.de or MATLAB Answers
3.) Ask Technical Support of MathWorks
4.) Go mad, your problem is unsolvable ;)
Private Nachricht senden Benutzer-Profile anzeigen
 
maki
Themenstarter

Forum-Anfänger

Forum-Anfänger


Beiträge: 32
Anmeldedatum: 23.03.20
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 03.06.2020, 11:54     Titel:
  Antworten mit Zitat      
Hallo,

ich habe den Code stark vereinfacht und einen anderen Optimierer eingebaut. Wie sage ich jetzt dem optimierer, dass Parameters verändert werden soll?

Hier der Code:

Code:

Xv0   = [0.4139*10^6 ,0.3004*10^6    ,0.4169*10^6    ,0.3813*10^6    ,0.3918*10^6    ,0.3914*10^6    ,0.3779*10^6   ,0.4804*10^6    ,0.6151*10^6  ,0.7395*10^6];% ,0.3*10^6];
Xt0   = [0.4351*10^6 ,0.3372*10^6    ,0.4384*10^6    ,0.3955*10^6    ,0.4067*10^6    ,0.4074*10^6    ,0.3853*10^6   ,0.4876*10^6    ,0.63*10^6  ,0.7521*10^6];% ,0.3*10^6];


p =cell (1,10);
p{1}= [0.24 70.8];
p{2}= [71.28 143.52];
p{3}= [144 216.24];
p{4}= [217.2 287.04];
p{5}= [288.24 358.32];
p{6}= [358.42 429.22];
p{7}= [429.32 474.776 545.192];
p{8}= [545.292 567.991 615.391 684.631];
p{9}= [684.731 707.675 754.876 825.052];
p{10}= [825.152 847.448 894.348 965.532 1059.852];


Parameters0(1)=0.03;
Parameters0(2)=0.0015;

Data= xlsread('BI130_Average_Daten.xlsx','D3:S30');
Exp_Data = Data (1:28,1:4);
Exp_Data= round(Exp_Data,3);

n = 10;
t=0;
y=[413900,435100];
   
for k = 1:n
c0 = [ Xv0(k); Xt0(k)];
p{k}= round (p{k},3);
[tSpeicher,ySpeicher] = ode45(@(t,y) Model_Hass(t,y,Parameters0), p{k}, c0);
    ySpeicher =ySpeicher(ismember(tSpeicher,p{k}),:);  
    tSpeicher =tSpeicher(ismember(tSpeicher,p{k}));

 t = vertcat(t,tSpeicher);
 y= vertcat(y,ySpeicher);
 
Magnitude(1)=1E6;
Exp_Data_VCD_= Exp_Data (1:28,2);
Exp_Data_TCD_= Exp_Data (1:28,4);

Exp_Data_VCD_ =Exp_Data_VCD_(ismember(Exp_Data(:,1),p{k}),:)*1E6;  
Exp_Data_TCD_ =Exp_Data_TCD_(ismember(Exp_Data(:,1),p{k}))*1E6;
 
  options = optimset ('MaxFunEvals', 1000,'PlotFcns',@optimplotfval);
[Estimated_Parameters]= fminsearch (@Objective_Hass_myjederunde,Parameters0,options,ySpeicher,Exp_Data_VCD_,Exp_Data_TCD_);
disp(Estimated_Parameters)


mit dem Model für die ODE

function dcdt = Model_Hass(t,c,Parameters0)


 my  = Parameters0(1);
 my_d= Parameters0(2);

Xv  = c(1);
Xd  = c(2);

dXvdt(1)= (my-my_d)*Xv;  
dXtdt(2)=  my_d*Xv + Xv*(my-my_d);    

   
dcdt=[dXvdt(1);dXtdt(2)];
end


und den Sum of squares als Funktion
function Sum_of_squares = Objective_Hass_myjederunde(Parameters,options,ySpeicher,Exp_Data_VCD_,Exp_Data_TCD_)


Sum_of_squares1= sum((ySpeicher(:,1))-Exp_Data_VCD_).^2/(mean(Exp_Data_VCD_)).^2 %VCD  
Sum_of_squares2= sum((ySpeicher(:,2))-Exp_Data_TCD_).^2/(mean(Exp_Data_TCD_)).^2 %TCD
Sum_of_squares= Sum_of_squares1+ Sum_of_squares2

end


Mir ist bewusst, dass ich Parameters0() keinen festen Wert vorgeben darf. Laut Syntax von fminsearch dient es als Startpunkt.
1.)Ist die Übergabe von Parameters0 bei fminsearch richtig? Wenn ich die Syntax richtig verstehe, müsste es dem Optimierer jetzt klar sein, was optimiert werden soll.
2.)
Zitat:
Objective_Hass_myjederunde(Parameters,options,ySpeicher,Exp_Data_VCD_,Exp_Data_TCD_)
hier ist mir bewusst, dass Parameters0 und options nicht verwendet werden. Entferne ich diese jedoch aus der Klammer, hat die Funktion zu viele Input Argumente. Woran liegt das?

3.) Allerdings stoße ich gerade auf den Fehler
Code:
Index in position 2 exceeds array bounds (must not exceed 1).

Error in Objective_Hass_myjederunde (line 5)
Sum_of_squares2= sum((ySpeicher(:,2))-Exp_Data_TCD_).^2/(mean(Exp_Data_TCD_)).^2 %TCD

Error in fminsearch (line 200)
fv(:,1) = funfcn(x,varargin{:});

Error in Call_Model_Hass_myfor1 (line 49)
[Estimated_Parameters]= fminsearch (@Objective_Hass_myjederunde,Parameters0,options,ySpeicher,Exp_Data_VCD_,Exp_Data_TCD_);
und mir ist nicht ganz bewusst, wieso.

Kann mir jemand bei den 3 Fragen weiterhelfen?

Viele Grüße,
Martin
Private Nachricht senden Benutzer-Profile anzeigen
 
Harald
Forum-Meister

Forum-Meister


Beiträge: 24.484
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 03.06.2020, 13:25     Titel:
  Antworten mit Zitat      
Hallo,

1. Du verwendest die nicht dokumentierte Syntax, feste Variablen als zusätzliche Parameter zu übergeben. Davon würde ich stark abraten, eben weil nicht dokumentiert und daher auch nicht mehr offiziell unterstützt. Stattdessen siehe hier
https://www.mathworks.com/help/opti.....ing-extra-parameters.html
am besten mit anonymous function handles.

2. Deine Zielfunktion arbeitet in der Tat nicht mit den Parametern. Bei Objective_Hass_myjederunde ist kein Zusammenhang zu der gelösten DGL vorhanden. Die DGL soll ja wahrscheinlich für jeden Parametersatz neu gelöst werden. Dann muss das in die Zielfunktion statt vorab berechnet zu werden. Ein Beispiel dazu findest du hier:
https://www.mathworks.com/help/opti.....ntialEquationODEExample-3

3. Löst sich vermutlich von selbst, wenn 1. und 2. umgesetzt sind. Falls nicht: Debuggen. Der Code erwartet anscheinend ySpeicher als Matrix mit 2 Spalten. Es ist aber wohl nur ein Vektor.

Grüße,
Harald
_________________

1.) Ask MATLAB Documentation
2.) Search gomatlab.de, google.de or MATLAB Answers
3.) Ask Technical Support of MathWorks
4.) Go mad, your problem is unsolvable ;)
Private Nachricht senden Benutzer-Profile anzeigen
 
Neues Thema eröffnen Neue Antwort erstellen



Einstellungen und Berechtigungen
Beiträge der letzten Zeit anzeigen:

Du kannst Beiträge in dieses Forum schreiben.
Du kannst auf Beiträge in diesem Forum antworten.
Du kannst deine Beiträge in diesem Forum nicht bearbeiten.
Du kannst deine Beiträge in diesem Forum nicht löschen.
Du kannst an Umfragen in diesem Forum nicht mitmachen.
Du kannst Dateien in diesem Forum posten
Du kannst Dateien in diesem Forum herunterladen
.





 Impressum  | Nutzungsbedingungen  | Datenschutz | FAQ | goMatlab RSS Button RSS

Hosted by:


Copyright © 2007 - 2024 goMatlab.de | Dies ist keine offizielle Website der Firma The Mathworks

MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.