Abstract
The increase in the country’s population attracted the establishment of more schools, both public and private schools, to cater for the increasing number of students. However, there have been dynamics of students’ population both in public and private schools through transfer from one category of school to the other, through completion of the learning period, and through dropout due to unknown reasons which have subjected both the public and private schools to compete in order to maintain a good number of students. In this work, a modified Lotka–Volterra model of schools and nonenrolled entities population in the education system is studied. Private schools and nonenrolled entities play the role of a predator in public schools. Again, public schools and nonenrolled entities play the role of predators in private schools. Holling type II functional responses have been integrated in the analysis of the Lotka–Volterra model. The equilibrium points are established and their stability are determined using the Routh–Hurwitz criterion and eigenvalue method. Global stability has been done for the positive equilibrium point. Bifurcation is also done around the positive equilibrium point. Finally, a graphical illustration of various parameter is derived to show their effect on schools when they are varied. It is revealed that the increase in parameters , , and greatly affects the schools population as they are the ones leading to predation in school. Therefore, proper strategies should be developed to focus on reducing the mentioned parameters to avoid leading schools’ population to extinct.
1. Introduction
Education stakeholders and strategists across countries work tirelessly to ensure that their respective country’s education system runs smoothly. Prey-predator models are usually used to study the interaction of animals and plants. This study applies the concept of a prey-predator model to study students’ population dynamics in schools. Prey-predator models can show different behaviors such as steady states, oscillations, and bifurcations depending on the parameters [1]. The study of the dynamics of prey-predator systems is one of the most dominant subjects in mathematical ecology due to its universal existence and importance [1]. The schools category population is considered in this study. The first category is the public schools which are owned and controlled by the government in their day-to-day operations. Teachers in such schools are government employees. The second category is private schools. United Nations Educational Scientific Cultural Organization (UNESCO) defines them as any school that is not operated by the government but is controlled and managed, whether for profit or not, by a private body like a community, foundation, faith-based organization, nongovernmental organizations (NGOs), private proprietor, or private enterprise [2]. This study also considers students who are out of school (nonenrolled) but are of school-going age as the third category and termed it as the nonenrolled category.
The first two categories of schools face challenges like peer influence among students and parents causing transfer from one category to the other, population predation between categories, and competition for students [3]. These factors create population dynamics in the respective schools. The school population is an important parameter for economic development. Predicting the future population is a core factor for development planning and decision-making [4]. Quality of education has decreased in public schools owing to overcrowded classrooms attracting more private schools [5]. Private schools tend to concentrate in an area where public schools perform poorly [5]. In the same way Italian mathematician Volterra in the model of fishing argued that, when fishing was good, the number of fishermen increased. After sometime, fish declined and the number of fishermen also declined. After some time, the cycle repeated [6]. However, the development of private schools increases the competition in the schooling system [7]. This theoretically leads to efficiency gain in terms of both cost and quality as private and public schools compete to attract students.
Due to prevailing situations in the country, where the country’s activities are geared towards achieving the Sustainable Development Goals (SDGs) which in addition is aimed at the realization of the vision 2030, the countries activities are thus focused on a better education system. This initiated the transfer of students to different categories of schools to acquire better education for self-reliance [3, 8]. Generally, the education sector shapes the nation through good leadership which can result to better development. For this matter, it is a secondary requirement to many people throughout the world. From the literature review, a range of social statistical techniques have been used to explain school population dynamics [2, 3, 5, 7–9]. But till now, very little research has been carried out on the interaction of public schools, private schools, and the nonenrolled entity population. Thus, it is important that approaches to modelling school population dynamics be considered based on past enrolment data incorporating the effects of competition and predation.
In this paper, a modified prey-predator model for three categories, namely; public schools, private schools, and the nonenrolled entity has been formed. So the uniqueness in this paper is that both public and private schools act as a predator as well as a prey. Nonenrolled entity is purely a predator. In this work, our main goal is to study the dynamical behavior of these three categories of schools, which is very important for the economic development [4].
Usually, Lotka–Volterra model comprises two equations, but this study presents the Lotka–Volterra model with three equations hence an improved one [10]. Holling extended the system of Lotka–Volterra and came up with the idea of functional response [11]. Thus, Holling type II functional response is used since our model adopts a logistic growth as proposed by [11]. Holling argued that when the prey is scattered in an area, the predator has to spend more of their time searching them as opposed to when they are concentrated in an area, where the predator will only have to spend most of their time handling the prey. Holling type II functional response with logistic growth is thus given by [1].
Lotka–Volterra model can be used to show the competitive relationship between organisms, companies, markets, or organizations [12–19]. Thus, this study considers it as also much applicable in schools’ population competition. The main purpose of modelling population interactions is to understand what causes such fluctuations [6]. This is only possible when the parameters of the model are varied [1, 19]. This study is helpful to the country for future population prediction [20]. This can help the country to be prepared to cater for the future changes that are likely to happen in schools. However, parameter estimation in this model was challenging. This is due to difficulty in getting predation values as well as the nonenrolled entity values as this data is sensitive and lies with the government. These are the challenges faced in modelling.
For the analysis of the model, the feasible region of the model is established and boundedness of the model is determined. The local stability of the equilibrium points determined with the help of a Wolfram Mathematica software is studied using the Routh–Hurwitz criterion [1, 21]. For the global stability, this study adopted the Lyapunov function. Finally, numerical simulation is done to verify the stabilities which cannot be determined analytically.
2. Model Formulation and Assumptions
In this work, the coexistence of the public schools, private schools, and nonenrolled entity is discussed. To formulate the model, we let , , and denotes the populations of public schools, private schools, and nonenrolled entity, respectively, at any given time t and they are subject to non-negative initial conditions , , and . The parameters and represents the recruitment rate in the categories and the environmental carrying capacity, respectively.
Figure 1 shows a modified Lotka–Volterra model with Holling type II is formulated where public, private, and nonenrolled entities assume the role of prey and predator depending on the scenario. In this model, the predator (private schools) population grows by enrolment and transfer from public schools and again, the predator (nonenrolled) gains population as a result of decreased enrolment in schools and dropouts. The total recruitment into the said categories is denoted by which is subdivided into the following categories as the recruitment in public schools, as the recruitment in the private schools, and () for the nonenrolled entities. Since public and private schools consider their carrying capacity during enrolment, their model adopts a logistic growth [22], while the nonenrolled entities assume an exponential growth. The population dynamics in the categories is due to the predation effect known as Holling type II response [11], exit rate, and the competition as represented by , and , respectively. To develop an equation for a category, for example, public school, we take the enrolment which is carried out by considering the carrying capacity, then we subtract the exit, transfer to private, and drop out due to unknown factor, less predation from nonenrolled, and private schools, and then we add the population that is gained by the public when predating private schools. By this, we arrive at equation (1). The exit rate is a fraction of the study period which is a quarter for the secondary schools. For the nonenrolled entities, the study considers their exit rate as the period which a person cannot come back to school to study. Anybody who has not joined the two categories of schools is considered to be in the nonenrolled category even if they joined tertiary institutions.
Moreover, since a model is just a representation of the reality, it cannot capture all the features of the school setup. Therefore, this study considers the following assumptions in order to formulate the model:(i)The study does not include special schools and group of schools in the country(ii)Dropout rate due to other factors other the ones discussed in the model is insignificant(iii)The rate of student retention in a class will be insignificant(iv)Public and private schools cannot predate the nonenrolled entities
From the above assumptions and illustrations, a mathematical model of the three categories has been developed in Figure 1.
3. Model Equations
with initial conditions , , and .
4. Feasible Region of the Solution
Consider equation (1), by inspection method the term and . Thus, . It follows, . On integrating, . . Exponential of a negative is always positive.
Consider equation (2), by inspection method, the term , and is .
Thus, , i.e., the exponential of a negative is positive.
Similarly, for equation (3), we have .
5. Boundedness of the Solution
For the system to be mathematically meaningful, it is necessary to show that its state variables are positive and bounded for all time [1]. That is, the solution of the system with a positive initial value will remain positive for all time .
On simplifying the above,
Thus,
6. Steady States
To study the stability of our proposed model, the equilibrium points of the system need to be determined. The possible equilibrium to be considered are Case 1: No schools exist; , , The trivial equilibrium points are Case 2: Only public schools exist; , , There are two equilibrium points say where Case 3: Only private schools exist; , , . The equilibrium points are , where and . Case 4: No enrolment in schools; , , . The equilibrium points are and when , where . Case 5: No private schools only; , , . Three equilibrium points are obtained , and , where where , , , , , , , , , , , , and . Case 6: Full enrolment; , , . Three equilibrium points are obtained; , , and , where where , , , , , , , , , , , , , and . Case 7: No public schools only; , , . There are three equilibrium points, namely, and , where, and which corresponds to the trivial equilibrium point. where , , , , , , , and . Case 8: Coexistence; , , .
There are many equilibrium points are , , , , , and . This corresponds to a trivial equilibrium point.where , , , , , , , , , ,, , , , , , , , , , , , , , , , , , , , , , , , , , , , , and .
Given that , , , , , , , and .
7. Stability Analysis
7.1. Local Stability
The local stability of the system is determined by the nature of the eigenvalues of the variational matrix, that is, for negative eigenvalues imply a stable otherwise unstable. In the case where the matrix is large such that the nature of eigenvalues cannot be directly determined, we will apply Routh–Hurwitz criteria to determine stability.
The system has eight solutions.
.
The variational matrix of the system is given bywhere
Theorem 1. The trivial equilibrium point, , is stable if the eigenvalues of the variational matrix are all less than zero. Otherwise, the system will be unstable.
Proof. The variational matrix of (16) evaluated at is given byBy letting and , eigenvalues are .
The equilibrium point is stable if and only if , , and . Otherwise, the system will be unstable.
Theorem 2. Equilibrium point, , at the point , , , by Routh–Hurwitz criterion is locally asymptotically stable if and . Otherwise, the system will be unstable.
Proof. corresponds to a trivial equilibrium point.The variational matrix is given bywhereThe characteristic equation is given by, , where , and .
Theorem 3. Equilibrium point, , at the point , , is stable if and only if the eigenvalues of the variational matrix are less than zero.
Proof. The variational matrix is given bywhere , , , , , .
The eigenvalues are , and .
The eigenvalues are stable if , , and .
Theorem 4. Equilibrium point, , for point , , is stable if and only if , and otherwise the system will be unstable.
Proof. The variational matrix is given byBy choosing and , the eigenvalues areThe eigenvalues will be negative if and only if , and for it to be stable, otherwise unstable.
Theorem 5. Equilibrium point, , for case , , by the Routh–Hurwitz criterion is stable if and otherwise unstable.
Proof. The variational matrix of the equilibrium point at , , is given byThe characteristic equation is given by , where
Theorem 6. Equilibrium point, , for the case , , is stable if and , otherwise the system will be unstable.
Proof. The variational matrix for the equilibrium point at , , is given bywhereThe characteristic equation is given by , where , , and .
Theorem 7. Equilibrium point, , for , , , by Routh–Hurwitz is stable if and , otherwise the system will be unstable.
Proof. The variational matrix for the equilibrium point at , , is given bywhereThe resulting characteristic equation is given by , where , , and .
Theorem 8. Equilibrium point, , for the case , , is stable if and , otherwise the system will be unstable.
Proof. The variational matrix for the equilibrium point at , , is given bywhereThe resulting characteristic equation is given bywhere , , and .
7.2. Global Stability
Linear stability analysis tells us how a system behaves near an equilibrium point. It does not however tell us anything about what happens farther away from equilibrium [23]. For higher dimensions, we consider using a Lyapunov function. This paper considers Lyapunov proposed by [1].
Theorem 9. Letbe the Lyapunov function, where are to be taken properly such that the above function satisfies the Lyapunov conditions.
, where and . The time derivative of V is implies that the system is just stable and implies that is globally asymptotically stable.
Proof. Getting a time derivative of (36),by substituting the values of results toBy factorizing ,At special cases, , , and implying that , , and are also equilibrium points, then substituting results toFurther simplifying, is the same as implying . Thus . Similarly, and . Since and , then they can be chosen such that , and that the equilibrium is globally asymptotically stable.
8. Bifurcation Analysis
In the Lotka–Volterra system, many parameters are used to describe and formulate the system. But if parameters used in the model are changed, other types of dynamical behavior may occur and the critical parameter values at which such transitions happen are called bifurcation points [1]. This study addresses the transition factors of the prey-predator interactions. This transition with respect to small changes is called Hopf bifurcation. Hopf bifurcation occurs at the point where the system has nonhyperbolic equilibrium with purely imaginary eigenvalues [1].
This study considers as the bifurcation parameters and as the critical values of the bifurcation parameters.
Theorem 10. The positive equilibrium over real domain is considered. Let be the continuous differential function of . Then, . Since is the positive root of . Thus, Hopf bifurcation of the interior equilibrium occurs at if and only if and .
Proof. By the first condition above that , the characteristic equation of the variational matrix in Theorem 8 can be written as .
The roots of the equation are .
There exist an interval of since is continuous of all its roots.
The general form of complex roots is given byWe verify the transversality condition .
By substituting the general form of the complex roots into the characteristic equation and getting the derivative results toBy implicit differentiation of equation (43), we obtainBy letting , , , and .
Thus, we haveFor at , therefore .
By comparing the values of using equations (45) and (46)., only if . Thus, the transversality condition is satisfied and hence Hopf bifurcation occurs at .
9. Numerical Simulation
In this section, the numerical analysis is done with the help of Wolfram Mathematica by considering the parameters values of the model cited in Table 1. Due to the difficulty in parameter estimation, we relied heavily on previously published work related to this study. The initial values of the variables used to generate graphs in Figures 2–6 are , [24] in page 10 and since the data for the nonenrolled was unavailable, the study considered the school-going age and further extracted the data from the 2019 national census statistics data. From the numerical analysis, we have found that all the eight equilibrium points are unstable. Since our model lies in the real domain, we will not verify for imaginary equilibrium points.
(a)
(b)
(c)
For Theorem 1, the eigenvalues are , , and ; hence, the system is unstable
Theorem 2, the equilibrium point is . Now, we have , , and and by Routh–Hurwitz, the system is unstable at that equilibrium point.
Theorem 3, the equilibrium point , the eigenvalues are , , and ; hence, Theorem 3 is unstable since it has positive eigenvalues.
Theorem 4, the eigenvalues are , , and , thus confirms that Theorem 4 is unstable.
Theorem 5, The equilibrium points are and which are unstable since , , and .
and ; we have , , .
Theorem 6, the equilibrium points are and . We have , , and .
and . We have , , and . From the above-given theory, it confirms that the system is unstable for Theorem 6.
Theorem 7, the equilibrium points are and . By Routh–Hurwitz, we have , , and . implies the system is unstable.
and . We have , , and . Since , the equilibrium points are unstable.
Theorem 8, the equilibrium points are . By Routh–Hurwitz criterion, we have , , and . Since , the equilibrium points are unstable.
. We have , , and .
. We have , , and . . We have , , and .
. We have , , and .
. We have , , and . Hence, Theorem 8 is unstable.
10. Discussion of the Results
Just like disease, predation help in population regulation [15]. But in this study, predation from nonenrolled entity needs to be avoided at all cost to save the schools population by avoiding school dropouts from learners. In this study, the system of ODEs formed had eight equilibrium points. From the numerical simulation, all the eight equilibrium points were unstable. This implies that the system cannot tend to one of the states of having no schools’ population and no nonenrolled entity. In Figure 2, the parameter and are varied from 0.0015 to 1.3 and 0.0008 to 1, respectively. In both cases as the value of the parameter increase the population of the public school and private school decreases. When the value of is 1, the population of private school becomes extinct. This implies these parameters are the ones responsible for the decrease in schools’ population. In Figure 4 as and increases the population of the nonenrolled entity increases. This is because and are the parameters that lead to predation of the schools population as seen in Figure 2. These parameters further cause a decrease in the private population as seen in Figure 3. In Figure 5, public and private schools’ population are almost directly proportional to each other while for public and the nonenrolled, increases in the nonenrolled causes a decrease in the public population. This is due to predation of the nonenrolled as confirmed by Figure 4. The same situation happens in when private schools are plotted against the nonenrolled entity as seen in Figure 5. Figure 6 shows the variations of schools’ population in the long run which was plotted for six years, public and private schools will be directly proportional while for public and the nonenrolled, public school populations increase steadily faster but it reaches a time when the nonenrolled overwhelm the situation and thus public population remains constant. Figure 7 shows in the absence of the nonenrolled the behavior of public and private schools in the long run. Figure 8 show a bifurcation diagram with respect to . The bifurcation diagram confirms the theoretical results of Theorem 8. The system of schools and the nonenrolled entity tend to coexistence when as in Theorem 8 illustrates the coexistence of the system.
11. Conclusion
In this paper, we have studied a model of competition of students’ population growth. The first and second categories of the population are the public schools’ population and private schools’ population, respectively, which grows according to the logistic growth equation. The third category is the nonenrolled entity which grows exponentially. In this work, eight equilibrium points were determined and then studied for the local stability and global stability of the system. The local stability of the equilibrium points has been studied using the eigenvalue method and the Routh–Hurwitz criterion. A numerical simulation has been done to verify the local stability. From the numerical simulation, it was found that all the eight equilibrium points are unstable. However, the system was found to be globally asymptotically stable using the Lyapunov function. Finally, from the results, this study suggests that the government of Kenya through the ministry of education and any concerned body should consider implementing strategies that are aimed at reducing the parameters , , and . This will save the schools population from becoming to an extinct. This is important because education is the backbone of a country for the development purposes as the nation requires people who are wise and can help to solve the problems of the country. By employing strategies that mitigate these parameters, it will also help to prevent future occurrences of the same. Bifurcation was carried out for and which show coexistence when greater than 0.2.
In this work, we have developed the Lotka–Volterra model for the system of ODEs to asses to effects of predation, school dropout, and competition among schools with Holling type II response. We used MATLAB software to generate graphs. However, concerning our future studies, we intend to do real data fitting to the model. Another direction that is also of interest is to do optimal control theory and stochastic modelling.
Data Availability
The materials and data used to support the study have been stated in the body and are included in the reference section of this paper.
Conflicts of Interest
The authors declare that there are no conflicts of interest.
Acknowledgments
The authors acknowledge the University of Embu for the insightful comments.