A mathematical model of food intake

1.   Introduction

"Obesity, anorexia nervosa, cachexia, and starvation are conditions that have a profound medical, social and economic impact on our lives" ^[1]. Overweight and obesity are defined as abnormal or excessive fat accumulation and are typically quantified by values of Body Mass Index (BMI). According to the World Health Organization (WHO), people whose BMI exceeds 25 are categorized as overweight while those whose BMI is over 30 are categorized as obese: "...in 2016, more than 1.9 billion adults aged 18 years and older were overweight. Of these, over 650 million adults were obese. Overall, about $13\%$ of the world's adult population ($11\%$ of men and $15\%$ of women) were obese in 2016. In 2016, $39\%$ of adults aged 18 years and over ($39\%$ of men and $40\%$ of women) were overweight. The worldwide prevalence of obesity nearly tripled between 1975 and 2016. Overweight and obesity are linked to more deaths worldwide than underweight. Globally there are more people who are obese than underweight: this occurs in every region except parts of sub-Saharan Africa and Asia" ^[2].

According to the WHO, there has been a worldwide increase in the availability of high-calorie foods coupled with a widespread decrease of physical activity for much of the world population. These concurrent factors have not been counteracted, so far, with effective educational and social policies ^[2]. We are therefore witnessing a world-wide surge of overweight and obesity of epidemic proportions, even in countries (like in South-East Asia or the Mediterranean) where traditional food culture and necessary physical labor used to maintain the population leaner in the past.

Obesity is a major risk factor for non-communicable diseases such as cardiovascular diseases (mainly heart disease and stroke) ^[3,4], which were the leading cause of death in 2019 ^[5,6] diabetes; musculoskeletal disorders (especially osteoarthritis – a highly disabling degenerative disease of the joints); and even for some types of cancer (including endometrial, breast, ovarian, prostate, liver, gallbladder, kidney, and colon) ^[7,8]. In fact, "... the incidence of obesity and its co-morbidities has increased at a rapid rate over the past two decades..." ^[1]. Of particular interest is the association of obesity with Type 2 Diabetes Mellitus (T2DM) ^[9,10,11], and, also mediated by this, the association of obesity with an increased risk of death ^[12,13]. During the COVID‐19 pandemic in 2019-2020, it was found that the rate of ICU admissions was associated with severe obesity ^[14]. All the evidence emphasizes the danger of obesity in our society.

Eating (food intake) plays an obvious role in the dynamics of BMI. Besides frank eating disorders (on a phychological or physiological basis), also cultural habits and culturally accepted behaviors may result in overeating and hence cause overweight, obesity and insulin resistance. Although the main cause of insulin resistance and T2DM is the skeletal muscle insulin sensitivity ^{[15,16,17,18,19,20]}, the possible effect of dietary intake on T2DM development occurs, if it does, through increase in fat mass and accompanying decrease in insulin sensitivity. Different theories exist linking causally overeating and insulin resistance, viceversa, or both ^{[19,21,22,23]}.

In the quest for a better understanding of the development of metabolic disorders, a mathematical model of food intake is needed. The importance of having a mathematical representation of a meal derives first of all from the general need to use quantitative methods to explain relevant pathophysiological issues ^[24]. Furthermore, modeling the dynamics of food intake would help us characterize its determinants and would provide a hitherto missing entry point into mathematical models of energy metabolism and of the development and management of diabetes mellitus. The purpose of the present work is therefore to build and calibrate such a mathematical model of food intake in man.

In this work, while appetite is in general the desire to eat, hunger is defined as the metabolic need to eat. This means that hunger is one component of appetite, determined by the need to have food in order to survive. On the other hand "...satiety is the feeling of fullness after eating that suppresses the urge to eat for a period of time after a meal." ^[25]. In order to introduce the criteria that have been followed in the design of the model, we begin with analyzing the two processes, meal initiation and meal termination. Meal initiation is a process which drives people to start a meal, while meal termination is the process which drives people to stop the meal. The causes which lead to meal initiation and termination are categorized into two main types, internal and external factors ^[26]. Internal factors are the factors related to the physiological control system (such as hormones and neurotransmitters) that helps us to have the right amount of food for energy balance. The interaction among internal physiological signals results in the feeling of appetite or in the feeling of satiety. When the environment is poor in food availability, this system predominates, the internal factors are the main drivers of appetite. However, when the environment is abundant, external factors tend to dominate the internal ones. External factors are defined as external signals, which drive or interrupt the eating process, such as physical activity ^[27], habits, food appearance, etc..

Regarding the internal factors, the mechanisms which drive hunger are part of the homeostatic system of human energy balance and blood glucose concentration (glycemia), whereas the pleasure after eating is regulated not only by (possible) metabolic factors, but also by the hedonic system of reward through food: people often eat in excess to their metabolic needs to satisfy their desire for palatable foods. The interaction between the two systems, homeorhetic and hedonic, is considered to play a major role in controlling eating behavior ^[28]. In contrast to homeostasis, homeorhesis is viewed as the ability of human to evolve our internal body system over time under different environments ^[29,30].

One of the internal factors which plays an important role in food intake is ghrelin. Ghrelin is a peptide hormone synthesized primarily in the stomach submucosa when the stomach is not distended by food ^[31]. There is some evidence showing that ghrelin regulates hunger in both short-term meal initiation and long-term energy homeostasis ^[32,31]. Ghrelin levels increase before meals ^[33], and, in turn, ghrelin acts on the hypothalamus by stimulating the production of NeuroPeptideY (NPY) and Agouti-related protein. "Agouti-related protein (AgRP) is a key orexigenic neuropeptide expressed in the hypothalamic arcuate nucleus and a marker for neurons conveying hormonal signals of hunger to the brain" ^[34], with the effect of increasing appetite. There is strong evidence supporting the fact that ghrelin increases food craving and reward-driven eating behaviors ^[35]. This evidence is however not uniform: experiments in mice subjected to ghrelin deficiency did not show significant changes in eating behavior under ad libitum standard rodent chow diet (CHD) ^{[36,37,38,39,40,41,42,43]}. However, reduction in rodent weight under ad libitum CHD was found in ghrelin-receptor-knockout experiments ^{[36,44,45,46,47,48]}. When both ghrelin and ghrelin receptors were suppressed, body weight in rodents was significantly reduced, although no significant change in food intake could be demonstrated ^[36,42]. While the evidence at this moment cannot therefore be considered conclusive, there are grounds to hypothesize that the ghrelin-signaling pathway is important for the control of food intake and body weight. For the purpose of the present model one could suppose ghrelin secretion to be inversely related to stomach volume, neglecting the intermediate role of NPY and thereby making appetite depend directly on circulating ghrelin levels.

While leptin also appears to be an important internal factor influencing energy balance and adipose storage, it is slow-acting and, contrary to ghrelin, its variations may not translate directly into short-term, hour-by-hour variations in food intake ^[49].

An important role in feeding behavior is also played by the components of the glycemia control system. Hyperglycemia (high blood-glucose concentration) from eating high-sugar diets is effective in reducing food intake, even if it is not clear how much of this effect is mediated by the inhibition of ghrelin production due to rising glycemia ^[33]. There is evidence linking the satiety induced by hyperglycemia with its effect in increasing insulin serum concentrations: rising insulinemia in fact inhibits dopamine signaling and interrupts feeding ^[31,50,51]. The role of hypo- and hyper-glycemia themselves has been repeatedly discussed both in generic health literature ^{[52,53,54,55,56,57,58,59]}, in scholarly articles and in books ^{[60,61,62,63,64,65,66,67,68,69]}. It is here maintained that hunger is a symptom of hypoglycemia. There is also some evidence that hyperglycemia by itself may suppress hunger ^[70,71,72]. Even though other contributions ^[73] point to the ineffectiveness of glycemia to regulate hunger, on the whole it would seem reasonable to include glycemia in a food intake model as a possible determinant of appetite, very possibly as a mere marker of more complex influences mediated by the appearance of food in the gastrointestinal tract. Since hunger due to hypoglycemia is not by itself suppressed by stomach filling (e.g. by drinking water or assuming a large amount of calorie-poor fiber), for the present work we may represent the effects of glycemia and ghrelin on appetite independently, even though we are aware that the degree of their interaction is currently under investigation. The ingestion of carbohydrates (in terms of their glucose equivalent ^[74]) and the absorption of glucose into the circulation should be considered as determinants of (the decrease of) appetite through the increase in glycemia they produce.

Gastric emptying is clearly an internal factor to be considered in the modeling of appetite and food intake, because it influences the degree of distension in the stomach (hence the secretion of ghrelin) and because it determines the rate at which glucose-rich or glucose-sparing food enters the absorbing bowel, determining hyperglycemia. Yokrattanasak et al. ^[75] provide a thorough review of the mathematical analysis of gastric emptying and introduce a stochastic pulsatile model, which may improve standard, deterministic exponential-decay representations. Of particular interest in the present context is the recent work by Phillips et al. ^[76], who measured gastric emptying rates in diabetic patients, and by Cao et al. ^[77] who investigate the dose-dependent effect of ghrelin on gastric emptying in rats and the related mechanism of action. While complex, stochastic models of gastric emptying are indeed possible, in the present context we have chosen to use a very simple representation (linear decay), already used by much of existing literature ^[78,79].

The hedonic system can also have a substantial impact on eating control. Several studies of the feeding-reward system have been conducted, some of them in the context of drug addiction, where however the findings may not be directly applicable to the normal homeostatic system, given the context ^{[28,80,81,82,83,84]}. Nevertheless food and drug rewards are both regulated by some common neural transmitters such as Gamma-Aminobutyric Acid (GABA) or Dopamine (DA) ^[28]. Opioid receptor occupation was found to have an impact on both feeding and drug rewarding ^[28,85]. Dopamine enhances reward-associated feeding behavior ^[31]. This hormone can trigger reward-induced overeating and even induce compulsive eating, with loss of control over food intake (food-addiction-like behavior), notwithstanding the rational knowledge of the negative consequences of such behavior ^{[86,87,88,89]}. Here again the evidence is not uniform and some authors found that a decrease in dopamine signaling did not inhibit eating in obese animals ^[31,90,91]. Of importance in connection with the present work is the finding that dopamine signaling and the secretion of dopamine can be activated by ghrelin ^[31,45], with the result that feeding is promoted. However, in conditions of lacking orexin, NPY and AgRP, or under vagal blockade, ghrelin fails to stimulate feeding ^{[31,92,93,94,95,96]}. Given this conflicting evidence, it does not seem appropriate to explicity account for dopamine in the construction of a food intake model. However, the neurobiological mechanisms of food and drug addiction may well be different ^[97,98].

Apart from internal factors, external factors can also have an impact on our eating behavior. It should be noticed that a large role in our daily consumption of food is simply played by habits, defined as repeating patterns of behavior over time. In fact, habits, one of the external factors, may influence our normal eating choices more than appetite does. Not only is habit determined by and determines the actual availability of food (cooking times, lunch-break times etc.), but eating habits influence our perceptions and memories of foods and vice versa. Past perceptions and feelings associated with eating have a strong impact on our next meal and can downright inhibit individual's appetite ^[99].

We devote the remaining part of introduction section to reviewing the studies of mathematical modeling of the systems related to the dynamics of food intake.

There has been an increased attention, over the past twenty years, to the possibility of modeling mathematically the development of metabolic disorders such as type 2 diabetes mellitus, also with the intent of allowing quantitative forecasts on the effectiveness of alternative therapies ^{[100,101,102,103,104,105,106,107,108,109,110,111]}. All of the mathematical models for the long-term development of diabetes, which have appeared so far, depend on the specification of the long-term development of insulin resistance as an input or forcing function to the model. Among the readily observable characteristics in a wide population (where e.g. IVGTT's or OGTT's are not routinely done), BMI seems to correlate with insulin resistance ^[112,113]. The correlation varies significantly across individual lifespans and across different populations, and very small increments in physical activity and/or trivial reductions in BMI can lead to large improvements in insulin sensitivity ^[114].

Many more papers have appeared in the literature about mathematical models of the glucose insulin system: besides classical short-term glucose-insulin models ^[115,116] and more recent derivations including gastrointestinal absorption ^[117], Liu ^[118,119] worked on a regulatory system of blood glucose by taking into account the dynamics of receptors at the molecular level; Wu ^[120] proposed a two-compartmental model in order to account for the oscillatory behavior of the glucose-insulin system; Lombarte ^[121] proposed a model including, besides blood glucose and insulin concentrations also intestinal glucose in rats. Other models include lipids in the description of the metabolic response, like in Pratt ^[122], where experimental data for the ingestion of three mixed composition meals over a 24-h period is used to study the response to a variety of mixed meals in fed and fasted conditions. The role of other hormones in regulating glucose metabolism has also been studied, in particular regarding incretins, gastrointestinal hormones that enhance insulin secretion by the pancreas in response to a meal: Girard ^[123] and Holst ^[124,125] review the concept of incretins and describe the biological effects of GIP and GLP-1 in normal subjects and in patients with type 2 diabetes mellitus; Kazakos ^[126] and Holst ^[127] even argument in favor of a possible use of incretin therapy in obesity and diabetes; Masroor et al. presented a model for the homeostasis of glucose through the regulating hormones glucagon and insulin which could be used for studies on glucagon receptor targeted therapy as well as on the artificial pancreas ^[128,129]; finally, a mathematical model was proposed to analyze the impact of pulsatile insulin release on post-prandial glucose and lipid control ^[130].

Once glucose (or other glucose-sparing lipidic or protein substrate) has been absorbed from the gastrointestinal tract, another sub-model should address the resulting glucose concentrations in plasma. The connection between the mathematical modeling of glycemic control and, more generally, of energy metabolism in humans, both in real-life and in artificial or therapeutic conditions, connects directly with the mathematical modeling of food intake and is a major motivation of the present work. Among many authors who have extensively published on glycemic control, our group has also been active, over the past 20 years, on the mathematical modeling of different aspects of the glucose-insulin system: on short-term glycemic control after perturbations ^{[131,132,133,134,135,136,137,138,139,141,142,143,144,145,168]}, on pulsatile insulin secretion ^{[146,147,148]}, on long-term diabetes development and population diabetes incidence ^{[100,102,107,149]}, besides the already mentioned work on stochastic gastric emptying ^[75].

A second step towards the representation of the control of food intake is the mathematical modeling of the digestion and absorption of foodstuffs. Several reviews exist of gastrointestinal tract physiology, including the physiology of the absorption of metabolic substrates and other nutrients ^{[150,151,152,153]}. Strathe ^[154] described a mathematical model for digestion and absorption of nutrients in growing pigs, considering dietary protein, endogenous protein, amino acids, non-amino-acid- and non-protein nitrogen, lipids, fatty acids, starch, sugars and dietary fibre. Salinari et al. ^[155] dealt with glucose absorption in the human. Taghipoor ^[156] developed a model of digestion in the small intestine focusing on dietary fiber. Several other meal absorption models have been proposed ^{[157,158,159]}. Again, within the framework of a relatively complex representation of the factors influencing food intake, the submodel concerned with gastrointestinal absorption of nutrients could be kept relatively simple.

Whereas, as touched upon above, a large body of mathematical modeling work exists in the areas of gastric emptying, food digestion, absorption and metabolism, very little has been done to date, as far as we know, about the actual act of eating. In Barnwell ^[160], a mathematical model of two decision-making circuits of feeding and satiety was proposed, considering the interaction of two neural centers (lateral hypothalamus and ventromedial hypothalamus), with feeding being terminated by feedbacks such as increased glycemia, temperature, or stomach distension. Boston et al. ^[161] modeled the rate of eating during a meal as a normal distribution, and then, using parametric deconvolution, investigated differences in eating patterns between healthy individuals and subjects with night eating syndrome. Cameron et al. ^[162] investigated the possibility of detecting meals and estimating meal sizes so that patients would be free from data input in order to calculate insulin dosages: it is of interest that in their work they use a Bayesian paradigm with a uniform prior for each single meal. In the work by Balakrishnan et al. ^[163], a personalized model of exercise, meal and insulin interventions for type 1 diabetic children and adolescents is proposed, to be used to predict adverse events such as exercise-induced hypoglycemia. Finally Jacquier et al. ^[164] propose a mathematical model to describe body weight, fat mass, fat-free mass, energy expenditure and food intake dynamics in rats, where food intake is regulated by available food and hunger, and where hunger and satiety are triggered by ghrelin and leptin plasma concentrations, respectively. Finally, Murillo et al. ^[165] developed a societal model describing the effect of transmissible cultural and psychological factors in shaping eating behavior at the population level, without however addressing the actual mechanics of feeding.

To this day, no formal, usable, quantitative model of the mechanisms whereby we start and stop our food intake is available for helping us understand normal and pathological meal patterns and the determinants of nutrition-related disorders. Given the background just discussed, we believe there is much scope for a realistic mathematical model of food intake. Such a model could be used as entry point in full-scale simulations of food absorption and metabolism, both in health and disease, and could help us understand the interrelationships of the several components involved. The present work intends therefore to formulate a realistic model of the initiation and cessation of food intake.

2.   Materials and methods

2.1.   The model

The model is a dynamical system composed of ordinary differential equations and algebraic equations, evolving over time $t$ in minutes. It is a continuous-time Markov Chain model, where the feeding state switches between fasting and eating with variable probability intensities, determined by prevailing hormonal and metabolic conditions. If the state of the system is suitably expanded, then a Markov approach can be used, since past events are reflected in current states.

Figure 1 shows a block-diagram schema of the relevant quantities in the model and of their interrelations. Red pointed arrows signify positive (stimulatory) effects, blue blocked arrows signify negative (inhibitory) effects.

In the following we list, define and discuss in turn the variables included in the model.

In this version of the model we assume for simplicity that physical activity (translating into an increased apparent linear elimination rate constant for glucose from the circulation) is higher during the day and lower at night. $E$ (/min) is therefore an indicator function of physical activity:

An important determinant of food intake is the actual availability of food, often determined by habitual activities like interrupting work, cooking food etc. $H$ (#) is an a priori, habits-determined switch-to-intake indicator at given times. We assume a truncated Gaussian-shaped intensity of given peak weight over the intervals corresponding to each one of four meals, plus a baseline snacking weight.

$L$ (pM) is the ghrelin plasma concentration. Ghrelin concentration decays linearly due to non-specific elimination from the bloodstream, consistent with general polypeptide catabolism, while secretion is actually limited exponentially by stomach filling:

A simple single-compartment model for glycemia is used in the present work as an example of how to close the food intake - glycemia loop (it can be replaced in applications by whatever other glucose metabolism models the Investigator deems appropriate). Denoting as $G$ (mM) the glucose plasma concentration, we can write

$A$ (#) is appetite (i.e. the desire for food determined by hormonal concentrations, metabolite concentrations, and stomach volume), in arbitrary units.

$k_{ji}$ (/min) is the probability intensity of switching from the fasting ("jeune", j) state to the eating ("intake", i) state, combining the contributions of appetite and of habitual food availabilty and consumption. Notice that the order of suffixes reflects the usual notation in Continuous-Time Markov Chain literature ($k_{ji}$ indicating switch from state $j$ to state $i$) and goes in the reverse direction with respect to usual compartmental modeling notation ($k_{BA}$ indicating transfer to compartment $B$ from compartment $A$).

$k_{ij}$ (/min) is the probability intensity of switching from the eating ("intake", i) state to the fasting ("jeune", j) state. We assume a baseline constant probability intensity to stop eating, associated with increased probability of stopping eating as the stomach fills up and as glycemia increases.

$u$ (#) denotes a Uniformly distributed random number, resampled at each evaluation time.

$\chi_i$ (#) is the indicator function of the feeding or intake state, switching with probabilities determined by the Markov Chain intensities.

Notice that the switching probabilities of the Markov process $X(t)$ are given by:

therefore, the indicator function of the "intake" state, $\chi_i$, evolves over time as:

In the numerical implementation of the model, care must be taken that the $t_{\Delta}$ term multiplying the appropriate intensities corresponds exactly to the integration time interval at the beginning of which the algebraic variables (among which the intake indicator) are computed.

$Q$ (g) is the quantity of food ingested over the time interval $t_{\Delta}$:

$S$ (g) is the stomach food content, modelled to exhibit linear emptying and incremented by the current rate of food intake.

2.2.   Parameter definitions and values

A few of the parameters in the model could be calibrated from information retrieved in the published literature, while others were calibrated in order to obtain model forecasts replicating, in the average, observed food intake patterns as reported in the 2015-2016 NHANES survey ^[166].

In the following we discuss the rationale for attributing specific values to some of the model parameters. Is should be noted that the parameter values referring to ghrelin are taken from the work of Cummings ^[167], whose studied sample consisted of 9 females out of 10 subjects: the results of simulations from our model might therefore be more appropriate for women rather than men. We take the molecular weight of ghrelin at 3371 Dalton.

- $t_{\Delta}$ : interval between successive observation times ($min$)

It is set up to be 2 minutes. This value has to be explicitly declared in the model implementation since the actual probability of starting or stopping eating depends on both the corresponding intensity and the length of time over which such intensity applies.

- $E^{TL}$ : time lower bound of activity period ($min$)

- $E^{TU}$ : time upper bound of activity period ($min$)

- $t_{M_m}^L$ : time lower bound of $m^{th}$ meal ($min$)

- $t_{M_m}^U$ : time upper bound of $m^{th}$ meal ($min$)

lower and upper bound for mealtimes determine when habits make it possible to actually access food for that meal, and are set to be rather wide, also taking into account the large amount of food intake outside canonical meal times observed in NHANES. It should be noticed that meal time-windows may overlap.

- $t_{M_m}^{\mu}$ : mean time of $m^{th}$ meal ($min$)

- $t_{M_m}^{\sigma}$ : standard deviation of time for $m^{th}$ meal ($min$)

- $w_{M_m}^{peak}$ : peak switch-to-intake weight of $m^{th}$ meal ($\#$)

- $w_{snacks}$ : switch-to-intake weight of snacks ($\#$)

Mean and standard deviation of food intake timing for the $m^{th}$ meal reflect the relative frequency of food intake within the permissible time window for each meal, with the weight expressing the relative importance of that meal within the overall daily food intake.

- $k_{XL}$ : rate of decrease in ghrelin level (elimination of ghrelin from the bloodstream) ($/min$). According to previous work by Toghaw et al. ^[140], the rate of elimination of ghrelin level is $0.02 \; /min$.

- $k_{XG}$ : glucose apparent linear elimination rate ($/min$). From ^[140], glucose clearance rate is $0.089 \; L/min$. Since $CL = 0.089$ and $V_G = 12.4\, L$ (defined in the section of determined parameters below),

- $k_{XGE}$ : additional glucose tissue uptake rate during exercise $/min$

We have hypothesized that a moderate level of aerobic exercise adds on the average, throughout the day, another 50% of tissue glucose uptake from the bloodstream, hence $k_{XGE} = 0.0036 /min$

- $\rho_{GS}$ : glucose bioavailability from food, the fraction of glucose in food available for absorption ($\#$). We assume that a fraction $\rho_{GS} = 0.9$ of the ingested glucose is actually absorbed into the bloodstream. As indicated by Beaugerie et al., the mean percentage of absorbed oligosaccharides and polysaccharides with degree of polymerization $> 2$ was $90\% \pm 3\%$ ^[169]. Similarly, the percentage of unabsorbed carbohydrates from breakfast was found to be $12.6\% \pm 10.6 \%$, yielding a percentage of absorbed carbohydrates of $87.4\% \pm 10.6 \%$ ^[170]. In animal studies, the apparent intestinal absorption of glucose was seen to be even lower: Atkinson et al. injected 36 $mg \; ^{14} C$-labelled glucose into the intestine in dogs and found that approx $70\%$ of the injected $^{14} C$ was absorbed as glucose ^[171].

- $\eta_{G}$ : glucose contribution from food (grams glucose per gram food) ($g/g$)

We assume an indicative value for $\eta_{G} = 0.2 g/g$. In order to derive this value, we considered a sample meal consisting of a hamburger and a beer. A McDonald hamburger (Big Mac) ^[172] whose weight is 99 grams contains 46 grams of carbohydrate. A 12-ounce (or 340-gram) bottle of microbrewed beer ^[173] contains 12 grams carbohydrates and 12 grams of sugars. Adding up the weight of hamburger bun and beer we come to about 439 grams, with approximately 0.159 g of carbohydrate per gram food. Of these, approximately $50\%$ ^[174] are directly absorbed as glucose, but not only the rest of absorbable carbohydrate is eventually absorbed as glucose, there is also some glucose-sparing effect of the ingested lipids (Randle's cycle ^[175]) and proteins. We therefore set the overall glucose contribution (direct supply plus sparing effect) of ingested food at approximately 0.2 g/g.

- $W$ : Body Weight ($kg$)

We consider here the world average weight of an adult ^[176], equal to $62 \; kg$. Notice that while this value may (substantially) differ from the average weight of the NHANES respondents, it is used in the model only for the purpose of calculating distribution volumes (which would scale up with increasing body size) and thereby metabolite and hormonal concentrations, which conversely should remain approximately the same (reflecting the equilibrium between scaled up production and scaled up elimination).

- $A^{max}$ : appetite maximum ($\#$)

We assume the maximum of appetite to be 300 (arbitrary units).

- $L_{A50}$ : ghrelin plasma concentration when appetite is half of its maximum ($pM$)

We assume ghrelin to affect appetite according to a non-linear, saturable curve, with ghrelin concentration at $50\%$ appetite equal to $120 \, pM$.

- $\lambda_{AG}$ : a factor representing the effect of decreasing glucose concentration on appetite ($/mM$)

The factor is assumed to be $0.3 \, /mM$. In other words, with this value of $\lambda_{AG}$ we assume that appetite will be reduced by half with a glucose increment of $log \, 2\lambda_{AG} = 2.31 mM$

- $\rho^{ji}_{HA}$ : conversion factor from habits and appetite weight to probability intensity of switching from fasting (j) to eating (i) ($/min$)

- $w_{A}$ : relative weight of Appetite (versus Habits) in inducing eating behavior ($\#$)

It is assumed to be 0.1. This low value reflects the fact that habits largely determine the intake of food at regular times, while appetite plays a role in a fewer occasions, as mentioned above.

- $k_{ij0}$ : baseline probability intensity of stopping eating ($/min$)

- $\rho^{ij}_S$ : conversion factor from stomach contents to probability intensity of switching from eating (i) to fasting (j) ($/min/g$)

This value has been calibrated against NHANES data to be around 0.0001.

- $\rho^{ij}_G$ : conversion factor from glycemia to probability intensity of switching from eating (i) to fasting (j) ($/min/mM$)

This value, like the previous one, has been calibrated against NHANES data to be around 0.0001. Notice that the implication of this and the above value is that stomach filling plays a larger role in stopping food intake than hyperglycemia, since $\rho^{ij}_S$ multiplies the grams of food in the stomach (in the hundreds or thousands) while $\rho^{ij}_G$ multiplies millimolars of glycemia (a few units or at most tens of units).

- $k_S$ : rate of food intake when eating ($g/min$)

The average rate of food intake has been calibrated to be around $16.5 g/min$.

- $t_{1/2}$ : half emptying time of stomach ($min$) According to Rocco ^[79], the half-emptying time $t_{1/2}$ of the stomach is approximately 89 minutes. However, the emptying time of the stomach depends substantially from the kind of food ingested ^{[177,178,179]}, and other sources, among which Locatelli ^[180], indicate shorter emptying times for different subjects. We have therefore assumed an overall half-time for stomach emptying of approximately 45 minutes.

- $L^{max}$ : ghrelin level at maximum ($pM$)

- $S_{50}$ : Stomach content at $50\%$ ghrelin secretion. ($g$)

We assume Stomach content to drive ghrelin secretion according to a non-linear, saturable curve, with stomach content at $50\%$ ghrelin secretion equal to $150 \; g$.

- $L_{ss}$ : ghrelin level at steady state ($pM$)

We assume that the plasma ghrelin level reaches steady state at midnight. Since the plasma ghrelin level at midnight is $600 \; pg / mL \; = (600)(1000)/3371 = 178 \; pmol/L$ ^[167], $L_{ss} = 178\, pM$.

- $G_{ss}$ : Glucose plasma concentration at steady state ($mM$)

From ^[181], normal glucose fasting level is between 3.9 to 5.5 mM. We set indicatively $G_{ss} = 5 \; mM$.

- $S_0$ : initial value of stomach food content ($g$)

It is usual to assume that at midnight the stomach is empty and thus we set $S_0 = 0 g$.

Some of the model parameters are not free to vary, but are determined by steady state considerations or other conditions and depend on the values of the other (non-determined) parameters in the model. In the following the determined parameters are defined.

The maximum rate of ghrelin production is computed from Steady-State considerations on ghrelin:

The rate of decrease in ghrelin secretion due to the amount of food in the stomach is computed from the stomach content at 50% ghrelin secretion, assuming an exponential decay:

$L_0$ is defined as initial value of ghrelin plasma concentration ($pmol / L$) which is set to be the ghrelin level at steady state $L_{SS}$ (or the plasma ghrelin level at midnight) equal to $178 \; pmol/l$.

$V_G$ is glucose distribution volume ($L$). From Toghaw et al. ^[140], the glucose distribution volume is approximately $0.2 \; L \, / kgBW$ (liters per kilogram Body Weight). For the adult world average 62 kg person, direct computation yields

$G_0$ denotes initial value of glucose plasma concentration ($mM$). It is assumed to equal glucose plasma concentration at steady state $G_{ss}$ which is $5\, mM$.

The assumed constant rate of glucose entry into plasma (from hepatic glucose output) is computed from Steady-State considerations at fasting:

The stomach emptying rate is computed from the stomach emptying half-time:

In this work, we adapt our proposed model to the available dietary data of the years 2015-2016 from NHANES ^[166]. More specifically, the data used in the present work is extracted from the data file DR1IFF. For each respondent, the weight (in grams) of food intake is considered based on the time when food intake starts. Respondent sequence number, time of intake, and weight of the intake are stored in the original NHANES variables $SEQN$, $DR1\_020$, and $DR1IGRMS$, respectively.

Tables 1 and 2 summarize respectively the model variables and the model parameters.

3.   Results

The model has been implemented in a mixed R/C++ environment where R-3.6.1 is used as a user interface and for producing graphics ^[182], and where an underlying computational engine in compiled C++ (Microsoft Visual Studio Community Edition 2019) is called to perform the actual model calculations.

Figure 2 portrays the time course, over 24 hours, of the model's main state variables. Meals are initiated, continued and interrupted as a consequence of glycemia, ghrelin dynamics and food availability or eating habits, without resorting to the action of any external trigger.

Figure 3 shows six sample days produced by the model using parameter values identical with those used for Figure 2: it can be appreciated that the meal pattern can vary considerably in both meal size and meal timing. The model may thus reproduce diverse patterns occurring in real life, such as prevalence of food intake at breakfast (panel a) or dinner (panel b), few large meals (panel c) or irregular food intake throughout the day (panel d), excessive food intake (panel e) or near-fasting (panel f).

The following figures show the adaptation of the model to the NHANES 2015-2016 survey data ^[166]. The distribution of food intake (grams) from 8505 individuals polled from the NHANES 2015-2016 database over the 24 hours of the day, by hourly intervals, is shown in Figure 4.

Correspondingly, Figure 5 shows the distribution of food intake (grams) averaged over 400 samples of 8505 virtual subjects each generated with the model, throughout the 24 hours and by hourly intervals. It can be appreciated that the overall pattern closely replicates the NHANES data.

When considering as typical meal periods the three intervals 6:00 A.M to 09:00 AM (breakfast), 12:00 to P.M. and 03:00 PM (lunch), 6:00 P.M. to 09:00 P.M. (dinner) as well as the agammaegated food intake outside of these three meals (snacks), we observe an average food intake recorded by NHANES of 480,495,540 and 1207 grams respectively (notice the large food intake outside of regular meal hours); the corresponding averages over 400 samples of 8505 model-generated virtual subjects are 490,445,554 and 1161 grams respectively.

The model can be used to obtain predictions about the likely effect of potential maneuvers directed at modifying the pattern of food intake, within the framework of scenarios of interest. In order to exemplify, in the following we study four possible scenarios where couples of parameter values are altered from their baseline values (baseline values are in every case those reported in Table 2). It should be kept in mind that in each case the actual time-course of the model variables depends both on the current parameter values and on the driving stochastic process defined in Eq. 2.8, incorporated in Eq. 2.10. In order to visualize the effect of changing parameter values, in each of the following figures the stochastic process is identical in the four panels (thereby representing for each of the four panels the same accidental combination of life events and subjective feelings and impulses of the virtual subject under consideration). In the first panel the baseline parameter values are used, and in each of the other three panels first one then the other then both of a couple of parameters are modified to show the results of different system regulations or different responses to the same external and internal environmental influences. In order to interpret the results, they should be compared with an average intake, from NHANES 2015-2016, of 2722 grams daily, which, under the assumption of a glucose-equivalent content of food equal to 0.2 grams/gram, correspond to 544 grams of glucose in the form of direct intake or of sparing effect.

We first assess what is the effect, in a normal everyday meal pattern situation, of changing the type of food ingested (carbohydrate overload, Figure 6). We simulate both decreasing the stomach emptying time (e.g. by consuming high-carbohydrate fast food, desserts and sugar-loaded drinks) and increasing the total carbohydrate content of the ingested food. The framework scenario we consider may be likened to a western-style workday pattern, with a large breakfast, a relatively smaller lunch and an early, substantial dinner.

We next investigate the potential to reduce carbohydrate intake by changing food composition in a scenario of habitual evening snacking. In this case, we hypothesize either an increase in the stomach emptying time (e.g. by adding vegetables or switching from carbohydrates to fat and protein in the diet), or a substantial reduction in the caloric density in the diet obtained by drinking an amount of plain water equal to the amount of the other food ingested.

As an alternative strategy to the reduction in food and carbohydrate intake with respect to the variation in diet composition, we explore the effects of paying more attention to the actual feeling of hunger (mindfulness) together with slowing down by 25% the rate of food intake, in a context of a large midday meal, such as used to be traditionally consumed in Mediterranean countries (Figure 8).

Lastly, we study a scenario of regular and controlled food intake, such as could be observed in an individual carefully following a diet, in order to see what changes might be produced by modifying the switching intensities (fasting to feeding or feeding to fast). Figure 9 portrays the reciprocal relationship of the tendency to initiate food intake with the tendency to stop eating. The probability to stop a meal is decreased by hedonic factors such as the palatability of the available food, as well as by psychological factors such as stress-related compensatory food intake. Conversely, the probability of initiating a meal also reflects practical factors such as the need to actually prepare the food and sit at the table before actual eating may start: the attendant lag between internal, appetite signals and the actual initiation of food intake may however be irrelevant if meal initiation is largely determined by habits.

4.   Discussion

Nutrient intake is obviously essential to the survival of the individual and of the species, and evolution has endowed humans with the ability of actually enjoying eating. The interrelationship of eating with psychological or emotional states is complex. Macht et al. ^[183] assess emotional states experienced in everyday life and examine the subjective motivation to eat associated with these emotional states. Things become even more complex in overweight individuals who are dieting: for example, Canetti et al. ^[184] summarize findings on the reciprocal interactions between emotions and food intake, concentrating on the psychological and emotional consequences of losing weight and dieting. Pecina et al. ^[185] examine the separate mechanisms of hedonic impact (food 'liking') and incentive motivation (food 'wanting'), raising the possibility that divergent neural control of hedonic and motivational processes may be associated with excessive eating and obesity. Complete reviews of the endocrinological controllers of food intake are offered by Fulton et al. ^[186] and Sam et al. ^[187]. It is clear that the subtlety and the complexity of the many endocrine and psychological factors involved in eating ^{[31,188,189,190,191,192]}, with the feeling of satiety following food intake ^{[186,193,194,195,196]} or with the influence of visceral sensations on the desire to eat ^[197] cannot be captured by a relatively simple model such as the one proposed here. On the other hand, metabolic research, in particular research on the mathematical modeling of energy metabolism, glucose-insulin control and diabetes development ^{[1,102,103,107,118,119,120,131,132,133,134,135,136,137,138,139,140,141,142,143,144,145,146,147,149,155,157,158,159,198,199,200,201]} would greatly benefit from a reasonable quantitative model of food intake, that could drive simulations of large groups of virtual patients subject to alternative hypothesized interventions.

An important basic premise, which needs to be explicitly declared, concerns the time-scale of interest. In the present work we are restricting ourselves to dynamics taking place over hours to (a few) days. It would indeed be of great interest to extend this analysis to body mass changes occurring over several years, particularly in connection with the current "epidemics" of obesity and diabetes in western and westernizing countries. This would implicate the consideration of other hormonal influences, notably leptin, as well as culturally determined lifestyle factors that interact with and possibly override hormonal metabolic control ^[202]. This is however a topic that will need to be addressed in depth in future work.

The main physiologic determinants of appetite seem to be a low blood concentration of glucose and a high blood concentration of the hormone ghrelin, presumably related to stomach distension, besides other diverse hormonal influences ^{[28,203,204,205]}. The dynamics of glycemia, stomach volume and ghrelin concentration are explicitly represented in our model.

One major determinant of feeding behavior, which we do consider in our model, is the actual availability of food, which we represent with an arbitrary variable, waxing at habitual mealtimes and waning in-between. This food availability of course comprises both "hard" determinant factors (for example the opening and closing hours of restaurant facilities in communal situations such as large companies, hospitals, army barracks etc.), as well as soft influencing factors (such as habitual office lunch break times). Peer pressure is in fact very important in inducing alimentation patterns, even to the extent that obesity appears to be significantly related with culturally spread behavior ^[206].

The tendency to continue eating once a meal has started or, respectively, to continue doing other things before a meal has actually begun ^[186,207] are also explicitly represented in our model through a continuous-time Markov stochastic process, whose intensities (to switch states between fasting and feeding) are function of the current values of the other state variables. The tendency to continue eating and individual increases in the probability of starting eating, given the same appetite level, can be related to the 'liking' element of feeding behavior as discussed by several authors ^{[185,196,207,208,209,210,211,212]}, whereas the "appetite" state variable of the present model directly translates the 'wanting' component described by these same authors (see also Finlayson et al. for a review ^[213]). It should be noted in this context that the appearance of an unusually palatable dish (or, conversely, of an unusually unappealing dish) is translated, in the current model, in extreme values of the driving random process. Moreover, the possible existence of metabolic inertia has been neglected in the present formulation of the model, in order not to expand unduly the state space while retaining markovianity.

Recent literature also points to the role of nutrients arriving in the distal ileum as a possible cause of the suppression of appetite ^[214,215]. By its very nature, given the length of the ileal tube and the necessary transit time, the "ileal brake" effect would act in suppressing appetite over the long term, rather than acutely controlling the food ingestion during a single meal. Further, while exposure of the ileum to lipids, carbohydrates and proteins may decrease gut motility and gastric emptying, the evidence that such effects may translate into long-term suppression of food intake is not clear. For both these reasons the ileal brake effect has not been represented in the currently proposed model.

It was found that ghrelin peaks before meals and declines after meals ^[167] and that ghrelin injections stimulate food intake rapidly and transiently, primarily by increasing appetite ^{[32,33,216,217]} and gastric emptying rate ^{[216,217,218,219,220]}. While ghrelin is the only known peripheral signal that has appetite-stimulating (orexigenic) effects, appetite-suppressing (anorexigenic) effects are determined by several other substances besides leptin, such as glucose, insulin, cholecystokinin, pancreatic polypeptide (PP) and Peptide YY (PYY) ^{[221,222,223]}. In summary, as the stomach empties it produces more ghrelin; increased circulating ghrelin induces firing of NPY and AgRP neurons, which in turn results in increasing appetite. After eating, leptin levels increase and inhibit the firing of NPY and AgRP neurons, determining a feeling of satiety. Since in the short term the effects of an increase in leptin are similar to the effects of a decrease in ghrelin and to the effects of an increase in glycemia, and since leptin has a major role as a long-term controller of body mass ^[224], for simplicity leptin was not explicitly represented in the present short-term model.

Short-term feeding and physical activity are not necessarily related to long-term patterns, nor are short-term behaviors related to long-term outcomes such as overweight and obesity. A short-term model does not determine long-term outcomes, especially given the nonlinear dynamics between short-term behaviors and long-term (chronic/habitual) patterns and outcomes. However, a short-term model would be extremely useful to provide input for a more general representation of metabolic events subsequent to feeding over days or weeks.

The model presented here captures several of the most important known factors determining food intake in man: stomach filling and emptying, glycemia, plasma ghrelin concentration, habit, probability to switch between feeding and fasting as determined by the practical arrangements of real life as well as by random events. The model appears to be consistent with the main ideas emerging from the rather extensive literature survey conducted for the present work. The simulations show time-courses of food intake, which are in fact consistent with observed, variable human behavior.

The model parameters have been calibrated by comparison with the 2015-2016 data from NHANES ^[225]. NHANES is a survey program designed to observe the health and nutritional status of adults and children in the US. This initiative has been conducted since the 1970s, and is being administered by the US federal government, a reputable institution. This survey is unique both for its scale and because it integrates both interviews and physical examinations. The source data are collected over a large-scale campaign of observation across the country for the past 50 years to produce a snapshot of the average conditions of the US population. The data from NHANES has been selected to calibrate our model and the data qualifies as merely indicative, rather than quantitatively precise. Some calibrated parameter values may thus more closely represent the eating behavior of the North American population: for example, it is of interest that, according to the data, a very large amount of food intake happens outside of regular meal hours. NHANES data may not be completely quantitatively reliable, since it is based on self-reporting, which may result in under-reported values ^{[226,227,228]}. However, NHANES is the largest publicly available data collection and, with the above caveats, is suitable for the calibration of our model. There may be a lack of physiological plausibility in NHANES over-reporters, respondents who mentioned a caloric intake twice as large (or more) as daily energy requirements. It should be noticed that in the present work, however, intake is considered in grams rather than calories.

This first mathematical formulation of a plausible model of food intake in man presents several limitations. The model deals with the relationship between ghrelin and glycemia in a rather summary fashion, by assuming independence and summation of their effects on appetite: such relationship is likely to be more complex and new versions of the model will have to incorporate new experimental results in order to better represent it. The myriad psychological factors modulating appetite and feeding behavior are also dealt with in a rather cursory way: in this case it is to be hoped that, by making parameter values vary randomly, following appropriate distributions, the effects of the concurrence of these many influences may be captured, to generate plausible populations of virtual subjects for stochastic simulation purposes.

Notwithstanding the above limitations, we hope that the model presented here can be used as a realistic driving input function for other mathematical models, representing the cascade of absorption, metabolism and control of glycemia over a period of several hours to a few days.

Conclusions

A plausible mathematical model of food intake, incorporating appetite control by glycemic levels, ghrelin, habits and the probability intensity of switching between feeding and fasting has been developed, for use as input function in realistic simulations of virtual subjects under free-living conditions.

Acknowledgments

This research project is supported by Mahidol University and is partially supported by the Center of Excellence in Mathematics, the Commission on Higher Education, Thailand.