Das mathematische Modell in Alverhag y Martin ist ein nichtlineares dynamisches System, das aus vier gruppierten Teilsystemen9 besteht. Die Subsysteme sind kompartimentale Darstellungen des menschlichen Körpers, wobei jedes Kompartiment ein Organ oder Gewebe darstellt, in dem ein wichtiger Prozess des Massenaustauschs durchgeführt wird. Die Kompartimente sind durch den Blutfluss miteinander verbunden. Dann quantifiziert jedes der Subsysteme mittels einer Massenbilanz in den Kompartimenten die Konzentration eines gelösten Stoffes (d. H. Glucose, Insulin, Glucagon oder Inkretine). Eine ausführliche Erläuterung des Systems und seiner Nomenklatur finden Sie in den ergänzenden Informationen.
Das System ist ein Satz von 28 ODEn, die aus nichtlinearen stetigen Funktionen bestehen. Daraus folgt, dass die Lösung des Systems (x(t)) in einer Domäne \(\mathbb {D}\) existiert, solange die Anfangsbedingungen in \(\mathbb {D}\) sind. Als methodischer Ansatz in dieser Arbeit wird die Lösung des Systems aus einer Zustandsraumtheorie als Vektor dargestellt:
wobei \(x(t)=(x_1(t), x_2(t), \ldots , x_{28 }(t))\in \mathbb {D}\subset \mathbb {R}^{28}\) ist halbfein positiv, was bedeutet, dass es zur Menge \(\mathbb {R}^{28}_{+}\). Verwenden der Zustandsdefinition in Gl. (1) ist das System definiert als:
wobei das Vektorfeld \(F(x(t);\pi , \eta ): \rightarrow \mathbb {R}^{28}\) die zeitliche Entwicklung von x(t) ab der Anfangsbedingung (\(x_{0}\)) in der Anfangszeit (\(x_{ t_{0})\), und \(\pi \in \Pi \subset \mathbb {R}^{46}\) enthält die Parameter in den Funktionen, die hämodynamische Prozesse darstellen, während \(\eta \in \)H\(\subset \mathbb {R}^{67}\) enthält die Parameter in den Funktionen, die die Stoffwechselraten des Systems darstellen. Die Parameterwerte des Systems in Gl. (2) finden Sie in den ergänzenden Informationen.
Modellsimulation und Initialisierung
Das mathematische Modell in Gl. (2) simuliert erfolgreich die Blutzuckerdynamik eines gesunden menschlichen Körpers nach intravenöser Glukoseinfusion und oraler Glukoseeinnahme9. Für das Obige wird eine Eingabe in das System betrachtet, die Folgendes enthält: (i) eine kontinuierliche intravenöse Glucose-Infusionsrate (\(r_{IVG}\)), die als Insulinrate in mg\(\cdot \) (dL \ (\cdot \) min.) in das System eingeführt wird)\(^{-1}\), und (ii) eine orale Glukoseaufnahme (\(OGC_{0}\)), die in mg in das System eingeführt wird und mit dem Magenentleerungsprozess verbunden ist (siehe die ergänzenden Informationen). Die Ausgabe des Systems (y) wird als \(x_6=G_{PV}\) und \(x_{14}=I_{PV}\) betrachtet, deren Bedeutung die Glucose- bzw. Die zeitliche Entwicklung von y wird verwendet, um die Modellsimulation mit klinischen Daten zu vergleichen, bei denen die Glukose- und Insulinkonzentrationen während eines Tests aus einer Blutprobe des Unterarms des Patienten entnommen werden. Für alle Simulationen ist das Modell in Gl. (2) wurde numerisch gelöst, indem ein variabler Schritt in der Funktion ode45 (Dormand-Prince) von MATLAB 18 verwendet wurde. Die Simulationszeit wurde als Zeitdauer der klinischen Studie definiert.
Für die Modell-Initialisierung wurden der Basalzustand \(x^B\) und \(x_0\) aus den gelösten Konzentrationen im nüchternen Zustand der Patienten berechnet. Der Zustand \(x ^ B\) wird als mittlere Nüchternglukose- und Insulinkonzentration aus den über mehrere Tage gesammelten Blutproben bestimmt, dies ist \(x_6 ^ B\) bzw. \(x_{14} ^ B \). Der Zustand \ (x_0\) wird als Nüchternglukose- und Insulinkonzentration aus einer Blutprobe zum Zeitpunkt Null des klinischen Tests bestimmt; dies ist \(x_6 (0)\) bzw. \ (x_{14} (0) \). Mathematisch hat der nüchterne Zustand eine physiologische Entsprechung mit dem stationären Zustand des Systems (\(x ^ *\)) in Gl. (2), das ist:
Da dann die interstitiellen, arteriellen und venösen Konzentrationen im Steady State gleich sind, werden die peripheren vaskulären Daten für \(x^B\) und \(x_0\) aus den arteriellen oder venösen Daten berechnet. Die restlichen 26 Komponenten von \(x^B\) und \(x_0\) werden aus der Lösung des Equ. (3).
Stoffwechselraten des Modells
Die in den Zusatzinformationen beschriebenen Subsysteme sind durch die Funktionen gekoppelt, die die Stoffwechselraten der Glucose, Insulin, Glucagon und Inkretine darstellen. Diese Stoffwechselraten werden mathematisch als konstante oder lineare Funktionen der Massenakkumulation in den Kompartimenten modelliert; oder multiplikative Funktionen der metabolischen Basalrate. Insbesondere sind die Stoffwechselraten in den Glucose- und Glucagon-Subsystemen multiplikative Funktionen mit der folgenden allgemeinen Form:
wobei \(r^B\) den Basalwert der Stoffwechselrate r darstellt und jedes M der isolierte Effekt der normalisierten Konzentration von Glucose (\(M^ G\)), Insulin (\(M^I\)) und Glucagon (\(M^{\ }\)) der normalisierten Stoffwechselrate (\(r^N=r/r^ B\)). Das Obige impliziert, dass \(M ^ G = M ^ I = M ^ \ Gamma = 1 \) wenn Glukose, Insulin und Glucagon basal sind, daher \(r = r ^ B \). Um die charakteristischen sigmoidalen Nichtlinearitäten biologischer Datenkorrelationen darzustellen, mit Ausnahme der isolierten Effekte, die Zustände des Systems in Gl. (2) (dh \(M_{HGP} ^ I\) und \(M_{HGU}^I \)), alle isolierten Effekte sind hyperbolische Tangentenfunktionen einer normalisierten Komponente des Zustands, dies ist:
wobei \(x^N_i=x_i/x_i^B\) für \ (i\in {\{1, 2, \ldots 28\}}\) und \(\eta _{j_1}, \eta _{j_2}, \ldots ,\eta _{j_4} \in H\) mit \(j_1, j_2 \ldots j_4 \in \mathbb {N} \le 67\) sind dimensionslose Parameter. Eine Liste mit den Nennwerten der Parameter \(\eta \) finden Sie in den ergänzenden Informationen. Unter Verwendung dieser Werte wird das System in Gl. (2) simuliert die Blutzuckerdynamik nach einer intravenösen Glukoseinfusion oder einer oralen Glukoseaufnahme in einem gesunden menschlichen Körper9. Für die mathematische Modellierung der Blutzuckerdynamik von T2DM muss die Pathophysiologie von T2DM emuliert werden, indem der Wert der Parameter der Funktionen geändert wird, die die Stoffwechselraten darstellen, die für die charakteristische Hyperglykämie verantwortlich sind. Das Obige wird in „Kurvenanpassung“ beschrieben.
Kurvenanpassung
Seit Jahrzehnten haben verschiedene Studien die Stoffwechselprobleme identifiziert, die mit dem Fortschreiten von T2DM bei gesunden Menschen verbunden sind19,20. Es wurde festgestellt, dass diese Probleme mit dem Stoffwechsel von Fetten und Kohlenhydraten zusammenhängen.19,20. Der Stoffwechsel dieses letzteren ist Gegenstand dieser Arbeit.
Hauptsächlich ist die Pathophysiologie des T2DM gekennzeichnet durch19: (i) Insulinresistenz, definiert als beeinträchtigte Wirkung von Insulin auf die Glukoseaufnahme durch periphere Gewebe, (ii) übermäßige Glukoseproduktion in der Leber aufgrund beschleunigter Glukoneogenese und (iii) \ (\ beta \) -Zell-Dysfunktion, dargestellt durch eine beeinträchtigte Pankreasinsulinfreisetzung. Dann werden die mathematischen Funktionen des Systems in Gl. (2) modellierung Die oben genannten Stoffwechselraten sind: die Wirkung von Insulin auf die periphere Glukoseaufnahme (d. H. \(M_{PGU} ^ I\)), die Wirkung von Glucose, Insulin und Glucagon auf die Glukoseproduktion in der Leber (d. H., \(M_{HGP}^{G}\), \(M_{HGP}^{I_{\infty }}\) bzw. \(M_{HGP}^{\Gamma _0}\)) und die Pankreasinsulinfreisetzung (d. h. \(r_{PIR}\)). Da eine geringe Variation der Parameter der nachstehend genannten Stoffwechselraten zu einer Variation der gelösten Stoffkonzentrationen im Modell führt, wird in den folgenden Abschnitten die Terminologie der Sensitivitätsanalyse von Khalil übernommen.21. Daher werden die obigen Stoffwechselraten als empfindliche Stoffwechselraten bezeichnet.
Im Folgenden wurden die sensitiven Stoffwechselraten so ausgewählt, dass sie zu den klinischen Daten von T2DM-Patienten passten. Explizit wird die Anpassung von \(r_{PIR}\) durch mehrere klinische Tests unterstützt, bei denen eine Abnahme der ersten Phase der Pankreasinsulinfreisetzung bei Patienten mit T2DM festgestellt wird22,23,24. Das Obige steht im Einklang mit dem frühen Vorschlag, eine teilweise Beeinträchtigung der Insulinfreisetzung aus dem labilen Kompartiment zu induzieren, um die erste Phase der Insulinfreisetzung bei T2DM-Patienten zu verringern25. Aus diesem Grund wurden die Funktionen, die die erste Phase der Insulinfreisetzung darstellen (X und \(P_\infty \)), und die Zeitvariation der Menge an labilem Insulin, die zur Freisetzung bereit ist, durch eine Sensitivitätsanalyse wie in Khalil21 untersucht, um die Parameter auszuwählen, die einen wesentlichen Beitrag zur Empfindlichkeit auf Lösung \(x(t;\eta,\pi _0) \) . Die ausgewählten Parameter wurden aus den klinischen Daten von T2DM-Patienten identifiziert. Die übrigen Parameter blieben unverändert.
Statischer und dynamischer Anpassungsansatz
Um das Parameteranpassungsproblem zu lösen, sind zwei Dinge erforderlich:
-
Eine Reihe klinischer Daten bei T2DM-Patienten.
-
Eine mathematische Methode, um solche Daten an die Funktion anzupassen, die die sensitiven Stoffwechselraten darstellt.
Der Satz klinischer Daten, die für die Analyse der isolierten Wirkungen verwendet wurden, wurde aus ausgewählten klinischen Tests von T2DM-Patienten erhalten. Die Bedingungen jedes der ausgewählten Artikel stimmen mit denen überein, die ursprünglich für die mathematische Modellierung in Ref.10. Diese Bedingungen sind in Tabelle 1 zusammengestellt. In den ausgewählten Artikeln wurden die klinischen Daten aus einer Menge \ (n_p \) von Personen ohne andere signifikante Anamnese als T2DM entnommen. Für die Kurvenanpassung verwendeten wir jedoch den gemeldeten Mittelwert der Gewebe- / Organantwort auf lokale Änderungen der gelösten Konzentration der \ (n_p \) Probanden. Ursprünglich, um die metabolische Rate \ (r_ {PIR}\) mathematisch zu modellieren, erhielt Grodsky Daten aus einer abgestuften Glucose-Schrittantwort mit der isolierten perfundierten Bauchspeicheldrüse in Ratten25. Da es unmöglich ist, diese Daten vom Menschen zu erhalten, wurden die ausgewählten Parameter dieser metabolischen Rate unter Verwendung klinischer Daten aus einem Input–Output-Ansatz des Systems in Gl. (2). Die Daten wurden einem OGTT in DeFronzo et al.26, wo die Plasmaglukose und die Insulinreaktion auf die orale Einnahme bei neun T2DM-Probanden nach dem Verzehr von 1 g / kg Körpergewicht oraler Glukose gemessen wurden.
Die mathematische Methode zur Anpassung der Funktionen an klinische Daten sind die kleinsten Quadrate (LSM). Im Allgemeinen liegt das LSM darin, dass die folgende Beziehung erfüllt ist27:
wobei z und \(\bar{y}\) Vektoren sind, die n Beobachtungen enthalten, und \(\theta \in \mathbb {R}^{p \times 1}\) ein Vektor von p unbekannten Parametern die empfindliche metabolische Rate. Um \ (\theta \) zu schätzen, werden die n Werte von g für alle z berechnet. Dann ist \(\hat{\theta }\) die Schätzung des Parametervektors entsprechend \(\theta \), der die Restsumme der Quadrate einer Zielfunktion \(Q(\theta )\) über einen Teil des Parametervektors \(\theta \ge 0 \subset \Theta \) minimiert. Die isolierten Effekte der sensitiven metabolischen Raten wurden durch einen statischen Ansatz des LSM an klinische Daten angepasst. Danach wurde ein dynamischer Ansatz des LMS verwendet, um die Parameter der Funktion \(r_{PIR}\) zu identifizieren. Im Folgenden werden beide Ansätze beschrieben.
Im statischen Ansatz werden die unbekannten Parameter aus dem Eq. (5) werden als \ (\theta = ^ T\) gruppiert. Der Vektor \(\hat{\theta }\) wird mit einem iterativen Prozess unter Verwendung der folgenden Zielfunktion geschätzt:
wobei \(y_k\) klinische Daten von der Mittelwert der normalisierten metabolischen Rate bei T2DM-Patienten ist sein Basalwert in Ref.9, und \(z_k\) sind die klinischen Daten des Mittelwerts der normalisierten gelösten Konzentration aus dem Unterarm. Die Minimierung der Zielfunktion in Gl. (7) wurde mit der Funktion lsqcurvefit der optimization Toolbox von MATLAB18 numerisch gelöst. Der iterative Algorithmus, der verwendet wurde, um \ (\hat{\theta }\) zu finden, wurde in Li28 ‚trust-region reflective‘ vorgeschlagen. Nach dem Anpassen werden (\(z_k\),\(y_k\)) grafisch mit den angepassten isolierten Effektfunktionen verglichen. Dann wurden die Werte der Parameter in \(\theta \) durch die Werte in \(\hat{\theta }\) ersetzt.
Im dynamischen Ansatz wurden die ausgewählten Parameter aus \(r_{PIR}\) als \(\theta = ^T\) mit \(l_1, l_2, \ldots l_6 \in \mathbb {N} \le 67\) gruppiert. Der Vektor \(\hat{\theta }\) wurde mit einem iterativen Prozess unter Verwendung der folgenden Zielfunktion geschätzt:
wobei \(y_{1k}\) und \(y_{2k}\) die klinischen Daten sind, die aus dem Mittelwert der Glucose- bzw. Insulinkonzentrationen zum Zeitpunkt \(z_k\) gewonnen wurden, die Gewichte \(w_1\) und \(w_2\) sind der Mittelwert der Basalglucose- bzw. Insulinkonzentrationen; und \(f_1=x_6(z_k,\theta )\), \(f_2=x_{14}(z_k,\theta )\) wurden aus der Modellsimulation erhalten. Die obigen klinischen Daten wurden von DeFronzo et al.26. Das LSM-Problem in Gl. (8) wurde mit der Funktion fmincon der optimization Toolbox von MATLAB18 mit dem iterativen Algorithmus ‚interior-point‘ numerisch gelöst. Nach der Identifizierung der Parameter von \(r_{PIR}\) wurden die Werte in \(\theta \) (aus dem statischen und dynamischen Ansatz) durch \(\hat{\theta }\) ersetzt, um die Pathophysiologie von T2DM nachzuahmen. Im Folgenden wird das resultierende Modell als T2DM-Modell bezeichnet.
Vergleich des T2DM-Modells mit klinischen Daten
Das T2DM-Modell wurde numerisch zum Vergleich mit einem klinischen Test in T2DM simuliert, bei dem die Blutzuckerdynamik nach verschiedenen Stimuli beobachtet wird. In Anbetracht dessen, dass der Weg des Glukoseeintritts in den Körper eine wesentliche Rolle spielt Gesamtglukosehomöostase26wurde das T2DM-Modell für den folgenden Test simuliert: (i) ein programmierter abgestufter intravenöser Glucose-Infusionstest (PGTT) zur Berücksichtigung der schnellen Reaktion der intravenösen Infusionen und (ii) ein OGTT unter Berücksichtigung einer Dosis von 50 g Glucose (50 g-OGTT) und einer Dosis von 75 g Glucose (75 g-OGTT) zur Berücksichtigung von Blutzuckeränderungen aufgrund des Magenentleerungsprozesses und der Auswirkungen des Inkretins.
Die klinischen Daten, die zum Vergleich des DMT2-Modells mit einem PGIGI-Test verwendet wurden, wurden von Carperntier et al.29. In diesem Test wurde die Glucose bei insgesamt 7 Probanden mit DMT2 intravenös verabreicht (d.h. \(n_p = 7\)). Mathematisch gesehen ist dies, dass die Glukose durch \(r_ {IVG}\) während \ (OGC_0 = 0\) zugeführt wurde. Die Dauer des Tests betrug 270 min verteilt wie folgt: Eine basale Probenahmeperiode wurde als \(r_{IVG}=0\) von 0 bis 30 min betrachtet, danach wurden die Schritte der intravenösen Glucoseinfusion als \(r_{IVG} eingeführt}=1, 2, 3, 4, 6\), und 8 mg (dL min)\(^{-1}\) für einen Zeitraum von jeweils 40 min. Die Bedingungen für die Modellsimulation waren \(G_{PV}^B= G_{PV}(0) = 157,5\) mg dL\(^{-1}\) und \(I_ {PV}^B=I_{PV} (0) = 13,02\) mU \ (\h_{L} ^{-1}\).
Die klinischen Daten, die zum Vergleich des DMT2-Modells mit einem OGTT verwendet wurden, wurden von Firth et al.30, und Mari et al.31. In diesem Test wurden 50 und 75 g orale Glucose von insgesamt 13 bzw. 46 Probanden mit DMT2 konsumiert (d. h. \ (n_p = 13 \) oder \ (n_p = 46\)). Mathematisch gesehen ist dies, dass die Glukose durch \(OGC_0\) während \ (r_{IVG} = 0\) geliefert wurde. Für die OGTT betrug die Dauer der Simulation 180 min. Die Bedingungen für die Modellsimulation waren \(OGC_0 =\) 50.000 mg, \(G_{PV} ^B = G_{PV} (0) = 185\) mg dL \ (^{-1}\) und \(I_{PV} ^B = I_{PV} (0) = 14\) mU L \ (^{-1}\) für die 50 g-OGTT. Ferner waren die Bedingungen für die Modellsimulation \(OGC_0 = \) 75.000 mg, \(G_{PV} ^ B = G_{PV} (0) = 176\) mg dL \ (^ {-1}\) und \ (I_{PV} ^ B = I_{PV} (0) = 11,2\) mU L \ (^{-1}\) für die 75 g-OGTT.
Der Unterschied zwischen den klinischen Daten und der Modellsimulation wurde mit dem folgenden statistischen Ausdruck quantifiziert:
wobei \(S_e=\sum _{s=1}^n(x_6(t_s)-G t_s)) ^ 2\), und G ist die Glucosekonzentration, die den T2DM-Patienten zum Zeitpunkt \(t_s \) entnommen wurde. Alle klinischen Tests unterschieden sich von denen, die für die Parameteranpassung verwendet wurden.