Journal of Applied Mathematics
Volume 2018 (2018), Article ID 2983138, 9 pages
https://doi.org/10.1155/2018/2983138

Research Article

A Dynamic Model of PI3K/AKT Pathways in Acute Myeloid Leukemia

Yudi Ari Adi,1,2 Fajar Adi-Kusumo,² Lina Aryati,² and Mardiah S. Hardianti³

¹Department of Mathematics, Faculty of Mathematics and Natural Sciences, Ahmad Dahlan University, Yogyakarta 55166, Indonesia
²Department of Mathematics, Faculty of Mathematics and Natural Sciences, Universitas Gadjah Mada, Yogyakarta 55281, Indonesia
³Department of Internal Medicine, Faculty of Medicine, Universitas Gadjah Mada, Yogyakarta 55281, Indonesia

Correspondence should be addressed to Yudi Ari Adi; yudi.adi@math.uad.ac.id

Received 28 June 2018; Revised 10 October 2018; Accepted 5 November 2018; Published 13 November 2018

Academic Editor: Oluwole D. Makinde

Copyright © 2018 Yudi Ari Adi et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Acute myeloid leukemia (AML) is a malignant hematopoietic disorder characterized by uncontrolled proliferation of immature myeloid cells. In the AML cases, the phosphoinositide 3-kinases (PI3K)/AKT signaling pathways are frequently activated and strongly contribute to proliferation and survival of these cells. In this paper, a mathematical model of the PI3K/AKT signaling pathways in AML is constructed to study the dynamics of the proteins in these pathways. The model is a 5-dimensional system of the first-order ODE which describes the interaction of the proteins in AML. The interactions between those components are assumed to follow biochemical reactions, which are modelled by Hill’s equation. From the numerical simulations, there are three potential components targets in PI3K/AKT pathways to therapy in the treatment of AML patient.

1. Introduction

Acute myeloid leukemia (AML) is a hematological malignancy originating in the bone marrow. It is characterized by the infiltration of the bone marrow, blood, and other tissues by proliferative, abnormally differentiated, and sporadic poorly differentiated cells of the hematopoietic system [1, 2]. The AML is the most common malignancy of hematological system, illustrated by the accumulation of acquired somatic genetic alterations in hematopoietic progenitor cells. This alteration modifies the normal mechanisms of cells proliferation, self-renewal, and differentiation [3, 4]. Based on the validated cytogenetics and molecular abnormalities, the National Comprehensive Cancer Network (NCCN) classifies patients into three risk categories, which are better risk, intermediate risk, and poor risk. Patients with the NPM1 mutation in the absence of FLT3-ITD and CEPBA mutations are classified as favorable risk and patients with FLT3-ITD mutated CN-AML and with P53 mutations are classified as poor risk [5]. Untreated AML patient results in fatal infection, bleeding, or organ infiltration within 1 year of diagnosis but often within weeks to months [6]. The standard therapeutic strategies in a patient with AML are chemotherapy, irradiation, and hematopoietic stem cell transplantation (HSCT) [1, 5, 7–10]. The main objective of those treatments is inducing remission and preventing the relapse [11]. In recent years, despite the potential gain of HSCT, the posttransplantation outcome remains dismal, especially those with high-risk category [10]. Currently, the development of the new therapies has been challenging to further improve the clinical outcome of AML, such as cytotoxic agent, small molecule inhibitor, and targeted therapies [12].

During the last decade, the PI3K/AKT signaling pathway has been studied extensively in human diseases. This pathway plays a significant role in a number of cellular functions, including differentiation, apoptosis, and cell cycle progression [8]. Aberrant PI3K/AKT activation is reported in 50-80% of AML cases [6]. The PI3K/AKT/mTOR network is activated in AML cells through a variety of mechanisms including upstream oncogenes such as FLT3-ITD, KIT, NRAS, and KRAS or autocrine/paracrine growth factors such as VEGF and IGF-1. It can also be activated by altered expression of p110d or phosphorylation of PTEN of pathway components and microenvironmental signals including chemokines and adhesion molecules [13, 14]. The AKT activation is associated with significantly elevated levels of phosphorylation FOXO3a in AML blast cells, supressing its normal function in induction apoptosis and cell cycle regulation [4, 15, 16]. Normally, FOXO3a transcriptionally activates several genes as the target. The FOXO3a binds to the promoter of apoptosis-inducing genes, such as Bim, FasL, and TRAIL, and to the promoter of cell cycle inhibitors, such as p27 and p21. The FOXO3a also activates the autophagy genes Gabarapl1, ATG12, and so forth [7, 9]. Researchers show that phosphorylation of FOXO3a is an adverse prognostic factor in AML associated with increased proliferation and overall survival [9, 16].

In a study of mathematical modelling of cancer, the dynamics of autologous immune system in chronic myeiloid leukemia have been constructed by Clapp et al. [17]. For mathematical models in leukemia and lymphoma, we refer the reader to Clapp and Levy [15]. Hovewer, these models are cellular level modelling, which study the interactions between cells and there is no known published model that studies the dynamics of acute myeloid leukemia (AML) in molecular level or biopathway. For modelling in biopathway, the model of PI3K/AKT pathway has been conducted in [18, 19]. The modelling of the interaction protein-protein in cell repair regulation also has been constructed in [20]. However, such models have not specifically studied a specific disease. In a previous work, Adi et al. in [21] studied the mathematical model of PI3K/AKT signaling pathway in AKT phosphorylation. The model did not include the FOXO3a protein that has been known as prognostic factor in AML. This model also does not follow Hill’s equation that describes the substrate-enzyme interaction that has multiple ligand-binding sites. In this paper, we construct a mathematical model of PI3K/AKT pathway in AML which considers FOXO3a, the most important downstream pathway of AKT. With deep understanding of PI3K/AKT signaling pathways, the important parameters that play a role in the development of AML will be identified. Furthermore, a strategy can be determined in the treatment of AML disease through targeted therapy. For the model, the Michaelis Menten kinetics and Hill’s equation in some components of biochemical reactions are used.

2. Model Development

The mathematical model is constructed by extending the AKT phosphorylation model in [21]. The extended model is defined by adding FOXO3a, which is a potential downstream pathway of AKT in the leukemic progression. Figure 1 shows the simplification of the complex network diagram of PI3K/AKT signaling pathways which drives the AML cells. Our model is focused on discussing the activities of five proteins in PI3K/AKT pathways that have been observed to be significant players in AML cells. We assume that the interactions of the protein follow the Michaelis Menten and Hill’s equation. The mutation of the growth regulatory genes such as FLT3-ITD is common in AML cases and results in the activation of PI3K [13]. The activation of PI3K then catalyzes the phosphorylation of phosphatidylinositol bisphosphate (PIP2) which can be phosphorylated at the D3 position of the inositol ring on extracellular stimulation, resulting in the formation of phosphatidylinositol trisphosphate (PIP3). In this model, the formation is grouped as a single process, called the PI3K level, and denoted by . This formation can be reversed by tumor suppressor PTEN that catalyzes PIP3 dephosphorylation into PIP2. The PIP3 dephosphorylation that is assumed follows Hill’s equation with coefficient 4. It is according to the fact that PIP3 has four binding sites of PH domain, that is, with PTEN, SHIP1, InsP4, and AKT [1].

Figure 1: Simplified diagram of PI3K/AKT/FOXO3a pathway.

Inactive AKT binds PIP3 which enables 3-phosphoinositide-dependent kinase-1 (PDK1) to phosphorylate AKT at Thr308. For full activation, AKT is also phosphorylated at Ser473 by mTORC2. PDK1 and mTORC2 are grouped as a single enzyme catalyzing the phosphorylation (activation) of AKT. Activated AKT is regulated by protein phosphatase 2A (PP2A) and pleckstrin homology domain leucine-rich repeat protein phosphatase (PHLPP). The phosphate PP2A preferentially dephosphorylates AKT on the Thr308 site, while PHLPP specifically dephosphorylates AKT on Ser473 site. In this model, the two phosphatases are grouped into a single enzyme catalyzing the dephosphorylation of AKT. According to the fact that AKT has two binding sites, this protein activity is assumed to follow Hill’s equation with coefficient 2. In the next downstream pathways, the activated AKT phosphorylates and inhibits the forkhead transcription factor, FOXO3a. FOXO3a phosphorylation promotes its translocation from the nucleus to the cytoplasm. Phosphorylation by AKTp on Thr 24, Ser 256, and Ser 318 inhibits FOXO3a activities by increasing nuclear export and this in turn increases proliferation. The FOXO3a in the cytoplasm, denoted by FOXO3ap, is the interaction with the 14-3-3 nuclear export protein. This interaction preventing nuclear reimport by concealing nuclear localization signals and promotes the FOXO3ap degradation by the proteasome [7, 16]. In this model, the growth of FOXO3a is assumed to follow the logistic model as well as the translocation and relocation of the cytoplasm and nucleus by phosphorylation and dephosphorylation. According to the fact that FOXO3a has three binding sites, Hill’s equation with coefficient 3 is used. The FOXO3ap enhances the expression and phosphorylation of RTKs which could in turn activate and sustain this pathways [3]. Thus, FOXO3a indirectly interacts and enhances PIP3 activity, resulting in a positive feedback loop. The reactivation or dephosphorylation FOXO3a, which is mediated by PP2A, promotes the relocation to the nucleus. Based on the diagram in Figure 1, a mathematical model is defined as follows:

The variables , , , , and represent the concentration of PIP3, AKT, AKT phosphorylation (AKTp), FOXO3a, and FOXO3a phosphorylation (FOXO3ap), respectively.

In the next section, we do some numerical simulation to understand the dynamics of protein in PI3K/AKT pathways. The numerical simulation will be run in two different situations based on the existence of the FOXO3a translocation from the nucleus to the cytoplasm as a normal cell or AML cell to understand the dynamics of AKT/FOXO pathways. In the normal cell, the activities of AKTp do not induce the translocation of FOXO3a from the nucleus to the cytoplasm. In the AML cell, aberrant PI3K/AKT signaling pathway results in phosphorylation of FOXO3a leading to cytoplasmic mislocalization and consequent degradation of these proteins [9].

The situation of normal and AML cells is distinguished based on differences in some parameter values. First, the constant rate of PIP3 dephosphorylation in AML cells is lower than normal cells as a result of various abnormal mechanisms in the PI3K signal upstream pathway, for example, PTEN deletion [8]. Second, the value of dephosphorylation rate of AKT in AML cells is smaller than the one in the normal cells. This is due to the fact that, in AML, there is a decrease of PIP3 level, a protein phosphatase that plays a role in the dephosphorylation of AKTp [8]. Furthermore, in normal cells, it is assumed that AKTp does not induce the translocation of FOXO3a from the nucleus to the cytoplasm so that the parameter value of phosphorylation of FOXO3a is tending to zero. The last difference is that the rate of dephosphorylation of FOXO3ap in AML is smaller than that in normal cells. This is due to the degradation of FOXO3ap in AML cells in the cytoplasm [7, 9].

3. Results and Discussion

The model equations (1)–(5) are not sufficiently accessible to allow us to conduct the mathematical analysis. Therefore, in this paper, we only provide numerical simulations. In this section, the numerical results of system (1)–(5) are simulated by employed the Runge-Kutta method of order 4 to provide the integration in some cases depending on the parameter values. The parameter values used in the system are based on the clinical data that can be obtained in some medical literature as in Tables 1 and 2.

Table 1: Initial concentration of the molecular component.

Protein
Concentration ()
References

PI3K
0.01 – 0.1
[21]

PIP3
0.7 – 0.8
[21]

AKT
0.01 – 1.0
[20, 21]

PP2A
0.004 – 0.15
[21]

FOXO3a
0.01 – 1.4
Assumed

Table 2: Parameter values and kinetic rates being used.

Parameter
Description
Unit
Value
References

PI3K level

0.01 – 0.1
[21]

Increase activation PIP3 by FOXO3ap

0.0083

Constant rate of PIP3 dephosphorylation by PTEN

0.0006 – 0.21
[21]

PIP3 degradation

0.001 – 0.01
[22]

AKT production rate

0.036 – 0.108
[21, 22]

Constant rate of AKT phosphorylation

1 – 20
[21]

Constant rates of AKTp dephosphorylation by PP2A

0.36 – 13.5
[21, 22]

AKT degradation rate

0.063 – 0.08
[21]

AKTp degradation rate

0.0008 – 0.1
[21, 22]

FOXO3a production rate

0.002 – 0.5

FOXO3a degradation rate by 14-3-3 protein

0.004 – 0.28

Constant rate of FOXO3a phosphorylation

0 – 0.33

Constant rate of FOXO3a dephosphorylation by PP2A

0.000297 – 2.92

FOXO3ap degradation rate

0.033 – 0.125

Michaelis constant of PIP3 dephosphorylation

0.01 – 1
[18, 20]

Michaelis constant of AKT phosphorylation

0.1
[20]

Michaelis constant of AKTp dephosphorylation

0.08 – 0.4
[20, 22]

Michaelis constant of FOXO3a phosphorylation

0.1

Michaelis constant of FOXO3ap dephosphorylation

0.1

Parameters values are assumed to be the same as other transcription factors such as MDM2 in [18, 22].

Tables 1 and 2 show the kinetic rates and the initial concentration levels of the various proteins in PI3K and AKT pathways. For AML cells the parameter values are ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; and For normal cells, similar parameter values are used, except ; ; ; and .

The dynamics concentrations of FOXO3a and FOXO3ap in the AML and the normal cell are given in Figure 2. Figure 2(a) shows that existence of the FOXO3a phosphorylation plays a role in the increasing of FOXO3ap concentration. The presence of inactive FOXO3a in the cytoplasm, FOXO3ap in our model, means that the apoptosis mechanism is not working properly so that there is no cell death. The FOXO3ap also enhances proliferation cell, leading to accumulation of abnormal cells because they do not stop growing when they should. The lifespan of the white blood cell in a myeloid lineage is about 3-12 days [23]. Therefore, there should be apoptosis between 3-12 days characterized by low-level FOXO3ap after that time, which does not occur in the AML cell. It can be seen that FOXO3ap reaches a peak in 100 minutes and then decreases and oscillates to a certain level (see Figure 2(a)).

(a)

(b)

(c)

Figure 2: Dynamic of FOXO3a and FOXO3ap concentration in AML cell (a), dynamic of FOXO3a and FOXO3ap concentration in the normal cell (b), and FOXO3ap concentration in the absence of FOXO3a phosphorylation (c).

The increasing concentration of the FOXO3ap would affect the decreasing concentration of FOXO3a. It is illustrated in Figure 2(a) that the FOXO3a concentration initially increases and peaks within 38 minutes. Moreover, the FOXO3a concentration immediately decreases and oscillates in low concentration with small amplitude. Under the normal condition, the concentration of FOXO3a transcription factor in the nucleus is much higher than those in the cytoplasm, FOXO3ap; see Figure 2(b). The concentration of FOXO3a in the normal cell reaches the maximal level in short time. The increasing of FOXO3a is followed by the slightly increasing FOXO3ap in much lower concentration. It indicates that FOXO3a is not translocated to the cytoplasm. It shows that FOXO3a promotes apoptosis and cell cycle regulation as well. Thus the balancing of cell cycle regulation can be well preserved. We note that if the constant rate of FOXO3a phosphorylation is set to be zero, the concentration of FOXO3ap will be zero (see Figure 2(c)).

Figure 3 shows the oscilation of FOXO3ap and FOXO3a concentration in AML cell. The oscilation of FOXO3ap concentration is given in Figure 3(a) and the oscilation of FOXO3a concentration is given in Figure 3(b). The oscillation indicates that cell cycle and proliferation continue to occur and there is no maturation of these cells. The behavior of FOXO3a can be used to identify the existence of AML disease. In addition, from the simulation, it is known that high level of FOXO3a and low level of FOXO3ap in the normal cell do not oscillate as in AML cell; see Figure 4. It indicates that FOXO3a works properly in cell cycle regulation and apoptosis.

(a)

(b)

Figure 3: The oscillation of the FOXO3ap concentration (a) and FOXO3a concentration (b) in AML cell.

(a)

(b)

Figure 4: The FOXO3ap (a) and FOXO3a (b) concentration in the normal cell.

Next, we will see the effect of phosphorylation of FOXO3a on the other proteins, such as PIP3, AKT, and AKTp. Figures 5 and 6 show the comparison of PIP3, AKT, and AKTp concentration in the normal cell and AML cell. It can be seen that, in the longtime behavior, the concentrations of PIP3 and AKTp in the AML cells are higher than those in the normal cells, while AKT in the AML cells is lower than that in the normal cells. In the AML cells, the concentration of PIP3 is gradually increasing and oscillates at a certain level as a response to FOXO3a phosphorylation. The PIP3 concentration reaches the maximum level and remains at that level for a long time and oscillates with small amplitude; see Figure 5(b). Figure 5(c) illustrates the concentration of PIP3 in the normal cells is at a low level without oscillation.

(a)

(b)

(c)

Figure 5: Dynamics of PIP3 in the normal cell and AML cell (a). Concentration of PIP3 in AML cell oscillates at a certain level (b), while PIP3 in normal cell does not oscillate (c).

(a)

(b)

Figure 6: Dynamics of AKT (a) and AKTp (b) in normal and AML cell.

Figure 6 shows the differences between AKT and AKTp behavior in the normal cell and AKT and AKTp in the AML cell, respectively. Under normal conditions, when the levels of PIP3 decrease, the AKT activity is attenuated by dephosphorylation by phosphatase. Figure 6(a) shows the AKT in AML cell subsequently sustained at lower concentrations than in the normal cell. In the AML cell, as the result of FOXO3a translocation from the nucleus to the cytoplasm, the concentration level of AKT decreases quickly, while AKTp immediately increases and remains at a certain level. Figure 6(b) shows that the AKTp in AML cell is subsequently sustained at a higher concentration than in the normal cell. Figure 7 tested the model by increasing the constant rate of FOXO3a phosphorylation from 1 (Figure 7(a)) to 2 (Figure 7(b)), while keeping all other parameter values the same as in Figure 2. The effect of this increase is that the greater the rate of FOXO3a phosphorylation, the lower concentrations of FOXO3a and FOXO3ap. This condition is due to the fact that the greater value of FOXO3a phosphorylation rate will accelerate the translocation of FOXO3a from the nucleus to the cytoplasm and lead to proteasome degradation. The increase of FOXO3a phosphorylation rate does not extremely affect the dynamics of AKT and AKTp, while the PIP3 concentration becomes slightly lower.

(a)

(b)

Figure 7: Dynamics of PIP3, AKT, AKTp, FOXO3a, and FOXO3ap in the AML cell with increasing rate of FOXO3a phosphorylation from 1.0 (a) to 2.0 (b).

4. Conclusions

As shown in the numerical simulation, the key components in driven AML cell are high levels of PIP3, AKTp, and FOXO3ap, that is, inactive FOXO3a in the cytoplasm. These results suggest that these three components are potential targets for AML therapy, of course with due regard to the other proteins that mediated protein interactions. For example, the parameter values associated with FOXO3a are taken from proteins that have a similar function to FOXO3a, such as MDM2. It could be very useful for determining reasonable ranges for the rate of various biochemical reactions involved. As more medical facts are known about the PI3K/AKT signaling pathways in AML, the model may be needed to be modified. It is possible that other equations representing rates of change of other proteins or other signaling pathways that integrated with PI3K/AKT pathways may have to be added to the system. Mathematical analysis of the model may also be useful in understanding protein interactions in this pathway. In the future studies, we will analyze mathematically the dynamics of the system and study the bifurcation related to the variation of its parameter values.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This work was supported by the Ministry of Research and Higher Education (Kemenristek DIKTI) of Indonesia (Grant no. 109/SP2H/LT/DPRM/2018). Special thanks are due to the Cancer Modelling Team UGM for the discussions during the research.

References

P. Faduola, A. Hakim, J. Mansnérus, A. Imai, and R. O'Neill, “Acute myeloid leukaemia - therapy - past, present and future,” Translational Biomedicine, vol. 4, no. 2, 2013.

G. Mengistu, F. Balcha, and S. Britton, “Clinical presentation of onchocerciasis among indigenous and migrant farmers in Southwest Ethiopia,” East African Medical Journal, vol. 76, no. 11, pp. 635–638, 1999.

Y. Shen, J. Bai, and A. He, “Role of mTOR signaling pathway in acute myeloid leukemia,” International Journal of Clinical and Experimental Medicine, vol. 9, no. 2, pp. 637–647, 2016.

S. Takahashi, “Downstream molecular pathways of FLT3 in the pathogenesis of acute myeloid leukemia: biology and therapeutic implications,” Journal of Hematology & Oncology, vol. 4, article 13, 2011.

NCCN, NCCN Clinical Practice Guidelines in Oncology, 2016.

G. A. Horne, R. Kinstrie, and M. Copland, “Novel drug therapies in myeloid leukemia,” Pharmaceutical patent analyst, vol. 4, no. 3, pp. 187–205, 2015.

R. S. Nho and P. Hergert, “FoxO3a and disease progression,” World Journal of Biological Chemistry, vol. 5, no. 3, pp. 346–354, 2014.

R. Polak and M. Buitenhuis, “The PI3K/PKB signaling module as key regulator of hematopoiesis: Implications for therapeutic strategies in leukemia,” Blood, vol. 119, no. 4, pp. 911–923, 2012.

C. M. Santamaría, M. C. Chillón, R. García-Sanz et al., “High FOXO3a expression is associated with a poorer prognosis in AML with normal cytogenetics,” Leukemia Research, vol. 33, no. 12, pp. 1706–1709, 2009.

G. J. Schiller, P. Tuttle, and P. Desai, “Allogeneic Hematopoietic Stem Cell Transplantation in FLT3-ITD-Positive Acute Myelogenous Leukemia: The Role for FLT3 Tyrosine Kinase Inhibitors Post-Transplantation,” Biology of Blood and Marrow Transplantation, vol. 22, no. 6, pp. 982–990, 2016.

T. Grafone, M. Palmisano, C. Nicci, and S. Storti, “An overview on the role of FLT3-tyrosine kinase receptor in acute myeloid leukemia: Biology and treatment,” Oncology Reviews, vol. 6, no. 1, pp. 64–74, 2012.

T. M. Kadia, F. Ravandi, J. Cortes, and H. Kantarjian, “New drugs in acute myeloid leukemia,” Annals of Oncology, vol. 27, no. 5, pp. 770–778, 2016.

M. Andreeff, Targeted Therapy of Acute Myeloid Leukemia, Springer New York, New York, NY, USA, 2015.

O. Lindblad, E. Cordero, A. Puissant et al., “Aberrant activation of the PI3K/mTOR pathway promotes resistance to sorafenib in AML,” Oncogene, vol. 35, no. 39, pp. 5119–5131, 2016.

G. Clapp and D. Levy, “A review of mathematical models for leukemia and lymphoma,” Drug Discovery Today: Disease Models, vol. 16, pp. 1–6, 2015.

S. M. Kornblau, N. Singh, Y. Qiu, W. Chen, N. Zhang, and K. R. Coombes, “Highly phosphorylated FOXO3A is an adverse prognostic factor in acute myeloid leukemia,” Clinical Cancer Research, vol. 16, no. 6, pp. 1865–1874, 2010.

G. D. Clapp, T. Lepoutre, R. El Cheikh et al., “Implication of the autologous immune system in BCR-ABL transcript variations in chronic myelogenous leukemia patients treated with imatinib,” Cancer Research, vol. 75, no. 19, pp. 4053–4062, 2015.

B. W. Keng and B. D. Aguda, “Akt versus p53 in a network of oncogenes and tumor suppressor genes regulating cell survival and death,” Biophysical Journal, vol. 91, no. 3, pp. 857–865, 2006.

G. Wang, Analysis of Complex Diseases, CRC Press, Taylor & Francis Group, 2014.

F. Adi-Kusumo and A. Wiraya, “Mathematical modeling of the cells repair regulations in nasopharyngeal carcinoma,” Mathematical Biosciences, vol. 277, pp. 108–116, 2016.

Y. A. Adi, F. Adi-Kusumo, L. Aryati, and M. S. Hardianti, “Modelling inhibition of AKT phosphorylation in acute myeloid leukemia,” in Proceedings of the 2016 Conference on Fundamental and Applied Science for Advanced Technology, ConFAST 2016, Indonesia, January 2016.

R. Karabekmez, “Modeling of Cancer Signaling Pathways,” UWSpace, 2013.

B. Young, G. O'Dowd, and P. Woodford, Wheater’s Functional Histology: A Text and Colour Atlas, vol. 64, Churchill Livingstone, Elsevier, 2014.