Model informed quantification of the feed-forward stimulation of growth hormone by growth hormone-releasing hormone.

AIMS: Growth hormone (GH) secretion is pulsatile and secretion varies highly between individuals. To understand and ultimately predict GH secretion, it is important to first delineate and quantify the interaction and variability in the biological processes underlying stimulated GH secretion. This study reports on the development of a population nonlinear mixed effects model for GH stimulation, incorporating individual GH kinetics and the stimulation of GH by GH-releasing hormone (GHRH).
METHODS: Literature data on the systemic circulation, the median eminence, and the anterior pituitary were included as system parameters in the model. Population parameters were estimated on data from 8 healthy normal weight and 16 obese women who received a 33 μg recombinant human GH dose. The next day, a bolus injection of 100 μg GHRH was given to stimulate GH secretion.
RESULTS: The GH kinetics were best described with the addition of 2 distribution compartments with a bodyweight dependent clearance (increasing linearly from 24.7 L/h for a 60-kg subject to 32.1 L/h for a 100-kg subject). The model described the data adequately with high parameter precision and significant interindividual variability on the GH clearance and distribution volume. Additionally, high variability in the amount of secreted GH, driven by GHRH receptor activation, was identified (coefficient of variation = 90%).
CONCLUSION: The stimulation of GH by GHRH was quantified and significant interindividual variability was identified on multiple parameters. This model sets the stage for further development of by inclusion of additional physiological components to quantify GH secretion in humans.

What is already known about this subject

Growth hormone (GH) secretion is pulsatile and the number of bursts, the amplitude and the burst interval varies greatly between individuals.

Limited information is available on stimulated GH secretion by GH‐releasing hormone while accounting for the individual GH kinetics in humans.

Existing mechanistic models describing the GH system do not include variance on the parameters or quantify the level of interindividual variability, limiting predictive performance for new study designs.

What this study adds

A population nonlinear mixed effects model of stimulated GH secretion, developed on literature information and experimental data, incorporating the individual GH kinetics and the stimulatory effects of GH‐releasing hormone was developed.

The level of variability on GH clearance, distribution volume and the GH‐releasing hormone secretion parameters were quantified, with a significant covariate relationship of bodyweight on GH clearance.

INTRODUCTION

Growth hormone (GH) plays an important role in many biological processes such as growth and cell reproduction, but also in several diseases such as pituitary adenomas. GH is secreted by somatotrophs in the anterior pituitary and can induce a variety of actions, among which the secretion of insulin‐like growth factor‐1 (IGF‐1), responsible for a large part of the growth promoting effects of GH.1, 2, 3 The regulation of GH is a complex interaction of stimulatory, e.g. GH‐releasing hormone (GHRH) and ghrelin, and inhibitory hormones, e.g. somatostatin.1, 4 At the arcuate nucleus in the hypothalamus, neurons secrete pulses of GHRH into the median eminence.5, 6 Via the hypophyseal portal system, GHRH binds to its receptor at the cell surface of somatotrophs in the anterior pituitary, triggering the intracellular signalling cascade, which releases GH into the systemic circulation.7 In response to the increase in GH concentrations in the blood, somatostatin is released from the periventricular region of the hypothalamus to inhibit further GH secretion.8 The interaction between these signals, all having short half‐lives, on a somatotroph results in a highly pulsatile and variable GH profile in plasma.

This variability is further increased in patients with acromegaly, who commonly have a pituitary adenoma that causes GH hypersecretion. The GH concentration–time profiles of these individuals are more stochastic, with higher basal concentrations and higher bursts of GH secretion.9 When assessing treatment effectiveness, sparsely sampled GH or IGF‐1 data are commonly used to inform decision making despite its known risk of misinterpretation due to the large variability in this patient population.10 The use of modelling and simulation techniques can improve interpretation of these results and provide additional information on treatment decisions by incorporating (patho‐)physiological information to describe complex biological systems, especially in endocrinology.11, 12, 13 The simulation of 24‐hour endogenous pulsatile GH profiles could eventually provide important insight to correctly assess a patient's GH secretion, study the effectiveness of treatment, or to inform on the optimal design of new clinical. However, to be able to simulate realistic and patient specific GH profiles, the feed‐forward stimulation of GH by GHRH and the level of variability in the healthy biological system needs to be quantified first.

Several attempts have been made to capture the hypothalamic and pituitary regulation of GH in physiologically based (systems) models, both in animals and humans.14, 15, 16, 17, 18, 19, 20 Unfortunately, these mechanistic models did not include information on the variability between individuals or the variance in the model parameters, thereby assuming a single typical pulsatile concentration‐time profile of endogenous GH secretion for all individuals. This use of a typical GH profile complicates the ability to judge whether a model is suitable to describe actual clinical data and can therefore not be used for (clinical trial) simulations. Existing empirical models of recombinant human GH (rhGH) do provide information on the interindividual variability in the kinetic parameters, such as clearance and the volume of distribution, but they are of limited use when the model needs to be expanded with physiological information on additional pathways in a different set of individuals, for which only the point estimate of the model can be used and the origin of variability in the response cannot be retrieved.

As a first step towards the mechanistic description of endogenous GH secretion, this study focuses on the quantification of the individual GH kinetics—after administration of rhGH—and the feed‐forward stimulatory properties of GHRH in healthy and obese women using a middle‐out estimation approach. Therefore, data from experiments performed in a cross‐over design, informing on different parts of the biological system in the same individual, were integrated in a population nonlinear mixed effects (NLME) model with the structural model based on physiological information.

METHODS

Subjects and data

The data for this analysis were obtained from a clinical study that has been reported previously.21, 22 The study was approved by the Medical Ethics Committee (METC) of Leiden University Medical Center, Leiden, the Netherlands. All subjects signed an informed consent form prior to the start of the study and the study was executed conforming to Declaration of Helsinki standards.

In this clinical study, 2 experiments were performed in healthy normal weight (n = 8) and lower (n = 8) or upper (n = 8) body obese women (median age = 37 years; interquartile range = [33–42 years], weight = 81.7 kg [69.8–95 kg], body mass index = 30.3 kg/m2 [25.1–43.2 kg/m2]). All subjects participated in 2 experiments: (1) the administration of rhGH; and (2) a GH stimulation test with GHRH. At the start of the rhGH experiment, a 2.5‐hour somatostatin infusion (50 μg/m2/h) was started to inhibit endogenous GH release and plasma samples were obtained every 10 minutes for GH analysis. One hour after the start of the somatostatin infusion, 33 μg of 22‐kDa rhGH was administered as a 5‐minute intravenous infusion, after which the plasma sampling was intensified to every 5 minutes for the first hour after dosing. Then, plasma sampling was reduced to every 10 minutes until the end of the experiment (t = 1.5 h after dosing). The day after experiment 1, a GH stimulation test was performed. Three plasma samples were taken at 10‐minute intervals after which the subjects received a fixed intravenous bolus dose of 100 μg GHRH. Samples were taken at 10‐minute intervals up to 3.5 h after dosing. For modelling purposes, the GHRH dose was converted to nmol by the molecular weight.

The obese women in the study additionally followed a weight loss diet and returned when they had lost 50% of their excess weight, after which the rhGH and stimulated GH occasions were repeated.21 Due to the time between visits and the loss of excess weight, the data of both visits were treated as originating from different individuals.

Serum 22‐kDa GH was measured by immunofluorometric assay (Delphia hGH kit, coefficients of variation ranging from 1.6 to 8.4%) with a lower limit of quantification (LLOQ) of 0.01 ng/mL. Data below the LLOQ were fixed to the LLOQ for analysis.

Additionally, in‐house data on the rhGH experiment in 15 healthy, normal weight female volunteers, subjected to an identical experimental procedure, were added to the dataset.

Model development

To distinguish between the different sources of variability within this population, model development was performed using the middle‐out approach,23 in which literature information on the system (system parameters) was combined with the estimation of population parameters, driven by the available data. This method was applied by first quantifying the individual GH kinetics after administration of rhGH, followed by adding the stimulatory properties of GHRH on GH secretion.

The model structure was based on physiological information of the pituitary and the somatotropic axis. The main physiological components that were included were the systemic circulation, the median eminence and the infundibular stalk, the anterior pituitary and the GHRH receptor at the extracellular surface of the somatotrophs. Unidirectional blood flow from the median eminence to the anterior pituitary and the systemic circulation was assumed, and system parameters were extracted from literature.

Population parameters were estimated using a population NLME modelling approach. A sensitivity analysis was performed on parameters that suffered from numerical instability. Inter‐individual variability on population parameters, assuming a ln‐normal distribution, was tested for a significant (P < .01) improvement in model fit. The investigated residual error structures were additive, proportional or a combination of additive and proportional.

Multiple variables (covariates) were explored to resolve the unexplained variability in the population: age, height, bodyweight, body mass index and body surface area. Covariate relationships were explored by visualization of the posthoc Bayesian estimates vs the covariates, and judged based on their Pearson correlation coefficient. Potential relationships were tested in the model using linear or power relationships and judged on a significant improvement in model fit (P < .01) with a backward elimination step (P < .01).

GH kinetics

Structural model development explored the addition of 1 or 2 single adjusting compartments (SACs), originating from the systemic circulation in order to encompass the distribution tissues that are in fast or slow equilibrium with the plasma compartment.24 The improvement in model fit after the inclusion of a baseline secretion parameter (GH baseline) of endogenous GH, mimicking endogenous GH release not fully blocked by the administered somatostatin, was investigated.

GH stimulation

After development of the GH kinetics section of the model, the estimated population parameters for the GH kinetics were fixed, thereby linking the kinetics of GH to the individual response of stimulation by GHRH. The GHRH receptor activation rate (k ), the receptor inactivation rate (k ), and the amount of GH released by a bound GHRH receptor (GHRH‐stimulated secretion) were estimated. The inclusion of a transit compartment, causing a delay between the activation of the GHRH receptor and the release of GH in the anterior pituitary, was explored. A baseline secretion parameter, independent of the exogenous GHRH, as was tested in the GH kinetics model, was tested for significance. Additionally, it was explored whether there was a significant improvement in the model fit after the inclusion of the individual posthoc Bayesian estimates for the GH kinetics compared to the population parameter estimates.

Local sensitivity of system parameters

System parameters that were included as point estimates in the models were subjected to a local sensitivity analysis. These parameters were increased or decreased by a factor of 2 or 10, and the impact on the area under the curve (AUC) for GH in plasma was evaluated up to 4 hours after GHRH administration by simulation of 500 concentration–time profiles while including interindividual variability and residual error. To assess whether the estimated population parameters could account for a possible bias in the system parameters, the increase and decrease of the system parameters were combined with a re‐estimation step of the population parameters. The resulting parameter estimates were then included in a new simulation (n = 500) and judged on a bias with the original model.

Model evaluation and internal validation

Models were evaluated on basis of the objective function value (OFV; −2 log‐likelihood), numerical evaluation, and goodness‐of‐fit (GOF) graphs.25 A significant (P < .01) improvement in the model fit was based on a drop in the OFV of minimally 6.64 points after addition of 1 additional degree of freedom in a nested model. The relative standard errors (RSEs) of population parameters and η‐shrinkage were considered acceptable when below 50% and 30%, respectively. The condition number, the ratio of the highest to lowest eigenvalue, was used to identify model overparameterization and should remain under the value of 1000. GOF graphs included the population (PRED) and individual model predictions (IPRED) vs observations, which should show a homogenous scatter around the line of unity. The conditional weighted residuals with interaction (CWRESI) vs PRED and time after dose should have the majority of the data between the [−2,2] interval and be homogenously distributed around 0. The GOF plots were checked for outliers and structural model misspecifications.

Models were internally validated using a nonparametric bootstrap with 1000 samples in order to compare the mean bootstrap results and the 95% confidence interval with the model parameter estimates, after resampling from the original dataset. Models were further internally validated using a confidence interval visual predictive check (VPC), which was created separately for the GH part of the model (data from experiment 1) and the feed‐forward stimulation model (data from experiment 2). The median and 80% prediction interval of the simulated model, with their corresponding 95% confidence intervals, were obtained from 500 simulations of the original dataset. The 80% distribution of the data was compared to the simulated intervals and judged on structural bias and the ability to correctly capture the variability in the data.

Software

Data transformation and graphical analysis was performed in R (V3.5.1).26 NLME modelling was performed in NONMEM V7.327 in conjunction with Perl‐speaks‐NONMEM V4.6.0.28

Nomenclature of targets and ligands

Key protein targets and ligands in this article are hyperlinked to corresponding entries in http://www.guidetopharmacology.org, the common portal for data from the IUPHAR/BPS Guide to PHARMACOLOGY.

RESULTS

System parameters

The system parameters that were obtained from literature on the volume of the systemic circulation, the pituitary blood flow and the kinetics of GHRH are depicted in Table 1. The total volume of the median eminence and the infundibular stalk compartment was approximated using a truncated cone formula, based on magnetic resonance imaging results.33 For GHRH, the distribution volume was assumed to be equal to the volume of the extracellular fluid with fast distribution kinetics due to the short half‐life of GHRH, mimicking 1‐compartment distribution kinetics.36

Table 1

Overview of physiological parameters identified in literature

Parameter	Value	Reference
Plasma volume (mL)	(40.5*HT) + (8.4*WT) ‐ 4811	29
Red cell volume (mL)	(16.4*HT) + (5.7*WT) ‐ 1649	29
Pituitary volume (mm3)	506.8	30
Anterior pituitary volume (mm3)	72.5% of total pituitary volume	31
	367
Ophthalmic artery flow (mL/min)	4% of total cerebral flow	32
	28.68
Qpituitary (mL/min)	11.47
Median eminence + infundibulum volume (mm3)	36	33
Half‐life GHRH (min)	6.8	34
Volume of distribution GHRH (extracellular fluid, L)	0.247*WT	35

HT, height in CM; WT, bodyweight in kg.

The blood flow kinetics of hypophyseal arteries in humans was essential in the development of this model, but was to our knowledge not reported in the literature. However, this blood flow to the pituitary could be correlated with the blood flow in the ophthalmic artery (OA), since the OA branches of the internal carotid artery at the same location as the superior hypophyseal artery (SHA; C6). The OA blood flow has been quantified as 4% of the total cerebral flow (28.68 mL/min) with an artery diameter of ~1.25 mm.32, 37 The SHA, with a diameter of ~0.5 mm, further splits into the inferior hypophyseal artery, of which the remainder passes through the median eminence and the infundibular stalk to the anterior pituitary.38 As such, it was assumed that the blood flow to the anterior pituitary was solely based on the difference in diameter between the OA and the SHA, resulting in an approximation of the blood flow of 11.47 mL/min. As plasma concentrations were quantified, the blood flow was corrected for by the ratio between plasma and red cell volume, based on the weight and height of a subject (Table 1).29

The GHRH receptor activation and inactivation rates (k and k ) were estimated as population parameters in the model using the data from the GH stimulation test. Since the GHRH receptor density on somatotrophs was unknown, no mass transfer between the anterior pituitary compartment and the GHRH receptor could be included. Therefore, the total amount of available GHRH receptors was fixed to 1.

The schematic structure of the model for GH and GHRH is depicted in Figure 1. The ordinary differential equations of the model, with the parameterization of all rate constants and compartment initializations, and the NONMEM model control stream is available in Supporting Information 1.

Figure 1

Model structure with the growth hormone‐releasing hormone (GHRH) and growth hormone (GH) kinetics. ECF, extracellular fluid; SAC, single adjusting compartment; QSAC, blood flow between plasma and SAC; Qpituitary, pituitary (blood) flow; kel‐GHRH, elimination rate constant for GHRH; CLGH, clearance of GH; kact, rate constant for receptor activation; kinact, rate constant for receptor inactivation

A total of 11 observations were below the LLOQ, all before the administration of rhGH. Exploratory analysis of the data at the start of the somatostatin infusion (t = −1 h) showed an initial wide scatter followed by a clear reduction in plasma GH concentrations until rhGH administration at t = 0 h (Figure 2).

Figure 2

Visual predictive check of plasma growth hormone concentrations relative to the time of recombinant human growth hormone (rhGH) administration (vertical dashed line). The black dots represent the observations, with their 10–50–90% distribution (black dashed lines). The shaded areas represent the 95% confidence intervals of the median prediction (red) and the 10th and the 90th prediction interval (blue)

The initial model structure overestimated the observations and the addition of SACs was explored. The addition of 2 peripheral SACs in the model provided a significant improvement in the model fit compared to the additional of a single SAC. However, the intercompartmental flow (QSAC‐fast) was high and could not be accurately estimated, indicating a fast equilibrium between the systemic circulation and this distribution compartment. This value was therefore fixed to 10 000 L/h to allow a quick/instantaneous equilibrium, which was required to accurately describe the data. A sensitivity analysis on the QSAC‐fast (range: 1000–5*105 L/h) did not result in significant changes in the population parameter estimates or the OFV. The inclusion of a 0‐order baseline GH secretion (GH baseline) in the pituitary compartment, resembling limited endogenous GH secretion despite somatostatin blockage, improved the model fit in the low concentration ranges and significantly reduced the OFV. Significant interindividual variability was included on: GH clearance (CL ), GH baseline and the distribution volume of the SAC compartment in slow equilibrium (V ).

Visual exploration of the covariate correlations suggested a relationship between an individual's CL and bodyweight. The inclusion of a linear covariate relationship (Figure 3) gave the highest reduction in OFV (−18 points) and reduced the η‐distribution correlation plot to a homogenous scatter around 0, superior to the inclusion of a power relationship. The coefficient of variation of CLGH was reduced from 19.8 to 16.5%, indicating that bodyweight accounted for a small but significant explanation of the variability on CLGH. No other covariate relationships were identified.

Figure 3

Posthoc Bayesian estimates of growth hormone (GH) clearance vs the bodyweight of subjects. Dashed line indicates the included linear covariate relationship on GH clearance. Red, healthy normal weight; green, lower body obese; blue, upper body obese

The GOF plots showed an adequate model fit, indicated by the close scatter around the line of unity in the PRED and IPRED vs observations (Figure S1A, B). The wider distribution in the lower concentration region of the PRED originated from the pre‐rhGH dose observations. The homogenous scatter around 0 for the CWRESI indicates no structural model misspecifications over the concentration range or vs time after dose (Figure S1C, D). The condition number was low (6.89) and the RSEs of population parameters were all below 30%. The final parameter estimates for the GH kinetics are given in Table 2. The VPC (Figure 2) shows that the model was able to capture the general trend and the variability in the population. However, 1 subject showed consistently high values outside of the 80% prediction interval. The exclusion of this individual did not change the parameter estimates, except for a nonsignificant decrease in the slope of the covariate relationship from 0.185 to 0.157 L/h/kg.

Table 2

Parameter estimates of the GH kinetics

Parameter	Estimate [RSE%] (CV%)	Shrinkage (%)	Bootstrap mean [95%‐CI]
CLGH‐slope (L/h/kg)	0.185 [6.1]	‐	0.188 [0.081–0.285]
CLGH‐intercept (L/h)	26.5 [3.65]	‐	26.58 [23.90–30.09]
VSAC‐Fast (L)	1.17 [29.4]	‐	1.26 [0.27–2.3]
QSAC‐Fast (L/h)	10 000a	‐	‐
VSAC‐Slow (L)	2.29 [6.63]	‐	2.27 [1.95–2.58]
QSAC‐Slow (L/h)	12.1 [11.4]	‐	11.84 [8.95–14.44]
GH baseline (μg/h)	1.04 [9.49]	‐	1.03 [0.75–1.36]
ω2 CLGH	0.0268 (16.5)	2.41	0.0250 [0.013–0.041]
ω2 VSAC‐Slow	0.0714 (27.2)	17.8	0.069 [0.029–0.122]
ω2 GH baseline	0.701 (101)	8.07	0.71 [0.36–1.18]
σ2 proportional error	0.0415	7.64	0.0417 [0.029–0.057]

CLGH equation, CLGH‐intercept + CLGH‐slope* (bodyweight – 70); RSE, relative standard error; CV%, coefficient of variation; CI, confidence interval; GH, growth hormone; SAC, single adjusting compartment.

Indicate fixed parameter.

Simultaneous estimation of k , k , and the GHRH‐stimulated secretion parameters resulted in numerical instability of the model. Due to the short half‐life of GHRH and the quick onset of the GH‐releasing effects of GHRH, the k needed to be fixed in order to improve the numeric stability. Based on a sensitivity analysis on k , a value of 100/nmol GHRH*h was implemented in the model, which resulted in the lowest OFV. The addition of a transit compartment between the activated receptor and the GHRH‐stimulated secretion lowered the OFV by 17 points but increased the RSE's from 2% to 62% for k and from 20 to 57% for the GHRH‐stimulated secretion, and was therefore excluded from the structural model. The inclusion of the posthoc Bayesian estimates of the individual GH kinetics significantly improved the model fit (ΔOFV = −84). Significant interindividual variability was identified on the GHRH‐stimulated secretion, which lowered the OFV by 1557 points (coefficient of variation = 89.7%). No covariates were identified that were explanatory of this variability.

The model parameter estimates of the GH stimulation part of the model are reported in Table 3. The confidence interval VPC is shown in Figure 4, which shows that the model is able to adequately describe the trend of the data and is able to capture the existing variability in the GH stimulation test. The GOF figures are shown in Figure S2. Population model predictions and observations show a high level of variability, indicative of the variability in GH secretion (Figure S2A). A small bias can be observed in the CWRESI over time, indicating that there is an underprediction of the GH concentrations around 2 hours after dosing (Figure S2C). This bias may originate from endogenous GH secretion or from a delay in the peak GH concentrations, which could not be estimated in the current model. A single outlier was identified (CWRESI = 7.5) which was due to the release of endogenous GH. The majority of the data (>95%) remains within the acceptance criteria of [−2,2]. The model was estimated using FOCE+I, ADVAN13, TOL = 7, NSIG = 3 and SIGL = 6.

Table 3

Parameter estimates of the growth hormone stimulation part of the model

Parameter	Estimate [RSE%] (CV%)	Shrinkage (%)	Bootstrap mean [95%‐CI]
kact (/nmol*h)	100a	‐	‐
kinact (/h)	1.46 [4.43]	‐	1.46 [1.32–1.57]
GHRH‐stimulated secretion (mg/h)	99.8 [14.9]	‐	101.1 [71.0–137.3]
ω2 GHRH‐stimulated secretion	0.59 (89.7)	< 0.01	0.58 [0.32–0.87]
σ2 proportional error	0.225	3	0.224 [0.178–0.280]

RSE, relative standard error; CV%, coefficient of variation; CI, confidence interval; GHRH, growth hormone‐releasing hormone.

Indicates fixed parameter.

Figure 4

Visual predictive check of plasma growth hormone (GH) concentrations relative to the time of GH‐releasing hormone administration (vertical dashed line). The black dots represent the observations, with their 10–50–90% distribution (black dashed lines). The shaded areas represent the 95% confidence intervals of the median prediction (red) and the 10th‐ and the 90th prediction interval (blue)

Local sensitivity analysis of system parameters

Using the estimated model parameters, a 2‐ and 10‐fold change in the half‐life of GHRH, the k , the pituitary volume or the distribution volume of GHRH changed the AUC of the GH stimulation by GHRH (Figure S3A/C). No impact in the median eminence volume or the ophthalmic artery flow on the AUC were observed. When the k and the GHRH‐stimulated secretion were re‐estimated with the altered systems parameters, all scenario's showed equal AUC's compared to the base model with only a minimal deviation when the half‐life of GHRH was changed (Figure S3B/D). These results indicate that the k and the GHRH‐stimulated secretion parameters can account for a bias in the included system parameters in the model.

DISCUSSION

The model described here was able to characterize the GH kinetics and the response to the administration of rhGH and the stimulation of GH by GHRH in normal weight and obese women. It allowed for quantification of the level of variability on the CL , which was partly explained by differences in bodyweight, on the V , and on the GHRH‐stimulated secretion. The study design affected different parts of the same biological system in the same subjects in 2 experiments, which allowed for a data‐driven approach for the estimation of the parameters in the model with high parameter accuracy (RSE < 30%).

Previous research indicated that the volume of distribution of GH was around 43.8 mL/kg (range 35.4–57.4 mL/kg) in healthy subjects and 54.9 mL/kg (39.9–84.4 mL/kg) in acromegaly patients.39 This would result in a volume of distribution of 3.1 and 3.8 L, respectively, for a 70‐kg person. For a typical subject in our study (70 kg, 1.7 m) the distribution volume was higher but in the same order of magnitude (6.1 L), which may be due to the different study population in our study or differences in statistical methods used in the calculation of the distribution volume. In our population, no significant differences between the groups (normal weight, upper body obese, lower body obese) in the volume of distribution was identified, indicating that increased weight was of little influence on the volume of distribution in this population or that the weight and height dependent system parameters already fully accounted for this variability. The effect of bodyweight/obesity on GH clearance was not significant after a noncompartmental analysis previously reported on the same data.22 The use of a more advanced population NLME analysis in this study was able to better characterize this significant linear covariate relationship, in which GH clearance would increase linearly from 24.7 L/h for a 60‐kg subject to 32.1 L/h for a 100‐kg subject.

In the model described here, GH is stimulated by GHRH in an attempt to quantify the feed‐forward mechanism present in GH control. However, existing biological knowledge reports on several feedback mechanisms that also control GH secretion (e.g. GH, somatostatin and IGF‐1),40 which cannot be identified using the currently available data. To quantitate the negative feedback, a different study design would be required that uses repetitive stimulation (e.g. multiple doses of GHRH). An alternative approach to estimate the impact of this feedback component could be the estimation of pulsatile secretion of GHRH underlying an endogenous GH profile, as a new component to this model, thus mimicking the hypothalamic function. Secondly, GH release will increase endogenous somatostatin concentrations that in turn block the release of GH. As such, a high burst of GH prior to the start of the experiments will have increased the endogenous somatostatin concentration, and thereby block part of the stimulatory exogenous GHRH effects during the experiment. To identify such feedback mechanisms, measures of target site concentrations of somatostatin in the anterior pituitary would be required, which is currently not feasible. In the current model results, this mechanism may cause the true variability of the GHRH‐stimulated secretion parameter to be overestimated. Additionally, the impact of the feedback from IGF‐1 was assumed to be limited due to the slow change in response of IGF‐1 to increases in GH concentrations. Despite the omission of several known feedback mechanisms, this model forms the basis for quantifying the feed‐forward relationship between GHRH and GH and understanding where variability in the somatotropic axis resides.

The local sensitivity analysis of the system parameters showed that the re‐estimation of the k and the GHRH‐stimulated secretion resulted in similar GH AUCs. This indicates that these parameters remain empirical in nature, since they are able to correct for a possible bias from system parameters, and therefore cannot be compared directly with experimental data. Furthermore, the system parameters were included as point estimates in the model, without any variability between individuals. Interindividual variability in these values can account for part of the identified variability in the model in which, for example, small changes in the half‐life of GHRH can impact the release of GH (Supporting Information Figure S3). Additionally, the estimated GHRH‐stimulated secretion parameter was estimated as a linear relationship based on a single GHRH dose, the addition of multiple GHRH dosing levels can inform on the shape of this relationship in future research. This model is therefore an approximation of the biological system, which, when available, can be updated with additional information to explain the k and k , or the GHRH receptor density in the human pituitary.

The data for this analysis were obtained from a heterogeneous healthy normal weight and obese female population. The structural model can be used for other populations when taking into account the general challenges in GH research. Many different reference standards, sampling methods and GH assays are used, which limits the comparability between studies.10, 41 For example, to account for a discrepancy between different GH assays, the GHRH‐stimulated secretion parameter can be re‐estimated, but the other parameters should remain in the same order of magnitude. Similarly, the use of a different GHRH analogue, with a different half‐life, may increase the duration of GH stimulation, which can be incorporated in the k parameter.

In conclusion, the presented model was able to capture the interaction between GHRH and GH, and quantify the interindividual variability in GH kinetics and the GH stimulation by GHRH in healthy and obese women. To improve the robustness of the model, more data from a larger population with wider distributions in age and bodyweight in men and women are needed to identify additional covariates that may explain the currently observed variability. A future step would be the expansion of this model by the addition of the hypothalamic control of GHRH to study endogenous pulsatile GH secretion in humans and the interaction with IGF‐1.

COMPETING INTERESTS

There are no competing interests to declare.

CONTRIBUTORS

M.J.V.E.: design and execution of project, writing of manuscript. J.B.: design of project and review of manuscript. P.H.V.D.G.: design of project and review of manuscript. J.S.: design and execution of project, review of manuscript.

Data S1 Supporting Information

Click here for additional data file.

FIGURE S1 Goodness of fit graphs of the recombinant human growth hormone (rhGH) concentrations (A) Population model predictions (PRED) vs observations, (B) individual model predictions (IPRED) vs observations, (C) conditional weighted residuals with interaction (CWRESI) vs time after rhGH administration, (D) CWRESI vs population model predictions.

Click here for additional data file.

FIGURE S2 Goodness of fit graphs of the feed‐forward growth hormone (GH) stimulation model. (A) Population model predictions (PRED) vs observations, (B) individual model predictions (IPRED) vs observations, (C) conditional weighted residuals with interaction (CWRESI) vs time after GH‐releasing hormone administration, (D) CWRESI vs population model predictions.

Click here for additional data file.

FIGURE S3 Simulation based local sensitivity analysis with a 2‐ or 10‐fold change in the system parameters (A, C) and with re‐estimated parameters (B, D). Boxes = 25–75% distribution of data, horizontal line in box = median, whiskers = 1.5× interquartile range, solid dots = data outside the 1.5× interquartile range.

Click here for additional data file.