Analytical solution of l-i SEIR model–Comparison of l-i SEIR model with conventional SEIR model in simulation of epidemic curves

The Susceptible-Exposed-Infectious-Recovered (SEIR) epidemic model has been commonly used to analyze the spread of infectious diseases. This 4-compartment (S, E, I and R) model uses an approximation of temporal homogeneity of individuals in these compartments to calculate the transfer rates of the individuals from compartment E to I to R. Although this SEIR model has been generally adopted, the calculation errors caused by temporal homogeneity approximation have not been quantitatively examined. In this study, a 4-compartment l-i SEIR model considering temporal heterogeneity was developed from a previous epidemic model (Liu X., Results Phys. 2021; 20:103712), and a closed-form solution of the l-i SEIR model was derived. Here, l represents the latent period and i represents the infectious period. Comparing l-i SEIR model with the conventional SEIR model, we are able to examine how individuals move through each corresponding compartment in the two SEIR models to find what information may be missed by the conventional SEIR model and what calculation errors may be introduced by using the temporal homogeneity approximation. Simulations showed that l-i SEIR model could generate propagated curves of infectious cases under the condition of l>i. Similar propagated epidemic curves were reported in literature, but the conventional SEIR model could not generate propagated curves under the same conditions. The theoretical analysis showed that the conventional SEIR model overestimates or underestimates the rate at which individuals move from compartment E to I to R in the rising or falling phase of the number of infectious individuals, respectively. Increasing the rate of change in the number of infectious individuals leads to larger calculation errors in the conventional SEIR model. Simulations from the two SEIR models with assumed parameters or with reported daily COVID-19 cases in the United States and in New York further confirmed the conclusions of the theoretical analysis.

Introduction SIR (Susceptible-Infectious-Recovered) and SIR-derived epidemic mathematical models (such as SEIR model; Susceptible-Exposed-Infectious-Recovered) are commonly used in the analysis of transmission of infectious diseases [1][2][3][4][5][6]. These models have been playing an important role in formulating the proper social interventions to slow down the spread of COVID-19 [7][8][9][10][11][12][13][14][15]. The SEIR model is a 4-compartment model. In the conventional SEIR model, it assumes that a susceptible individual in the compartment S may enter compartment E after the individual is exposed to an infectious source; an exposed individual in compartment E may enter compartment I at a certain rate and become an infectious individual; and an infectious individual in compartment I may enter compartment R at a certain rate and become a recovered individual. Using the four variables S(t), E(t), I(t) and R(t) to represent the numbers of individuals in the four compartments at time t respectively, the relationship between the four variables is defined by the following differential equations [2]: In Eq (1), N is the total population; and the coefficients β, σ and γ are constants. One of the underlying assumptions for the above equations is: the change in E(t), I(t) or R(t) at time t is determined by a fraction of E(t) and I(t) (σE(t) and γI(t)) at time t. This means that after a susceptible individual in compartment S becomes an exposed individual in compartment E at time t, this newly exposed individual is able to enter compartment I at the same rate σE(t) as all other pre-existing individuals in compartment E, without chronological order. However, the reality is different: after a susceptible individual is exposed to the infectious source and becomes an exposed individual in compartment E, the newly exposed individual cannot become an infectious individual and enter compartment I until the exposed individual passes a latent period (l days on average). Likewise, an infectious individual in compartment I can only recover after an infectious period (i days on average). Thus, a l-i chronological order is an inherent relationship that exists among the 4 variables, whereas the conventional SEIR model uses an approximation of temporal homogeneity to calculate the transfer rate of infected individuals moving from compartment E to I and from I to R. Some people may ask, what calculation errors will the conventional SEIR model make when simulating the epidemic data? In this study, we developed a 4-compartment l-i SEIR model from the 3-compartment l-i AIR model published previously [16,17]. The l-i SEIR model, like the l-i AIR model, has considered temporal heterogeneity of individuals in the model. The analytic solution of the l-i SEIR model in closed-form was derived. The 4 compartments in the l-i SEIR model have the same meaning as the 4 compartments in the conventional SEIR model. Both the l-i SEIR model and the conventional SEIR model have the same series of 4 compartments from the compartment S to E to I to R. Thus, we are able to examine how individuals move through each corresponding compartment in the two SEIR models to find what information the conventional SEIR model may miss and what calculation errors in the conventional SEIR model may result from using the temporal homogeneity approximation.

Theory and methods
l-i SEIR model. The l-i SEIR model is developed from a recently published l-i AIR model based on the following three assumptions: 1. During the spread of an infectious disease, the decreasing rate of the number of susceptible individuals on day n, (S n-1 -S n ), equals to the number of infectious individuals (I n ) times the ratio of S to N (S n-1 /N) and multiplies a rate coefficient β n (b n I n S nÀ 1 N ). This assumption is the same as the one described in the classic SIR model. Here, N is the number of susceptible people right before the infectious disease spread out. If all people in the population are susceptible to the infectious agents before the infectious disease spreads out, N equals to the total population P. However, if a portion of people has immunity to the infectious disease before the infectious disease spreads out, N is smaller than P. It should be emphasized that the impact of other interventions including social distancing, face mask requirements and quarantine measures on the transmission rate of infectious diseases is factored into the coefficient β n . Because the policies and requirements of these interventions may change during the transmission process of an infectious disease, β n is a time-dependent coefficient in the l-i SEIR model.
2. When an individual in compartment S enters compartment E as a newly exposed individual, this individual will stay in compartment E for a latent period (l days on average) before entering compartment I as a newly infectious person. The newly infectious person will stay in compartment I for an infectious period (i days on average) before entering compartment R. The total time length of the latent period and the infectious period is c (c = l+i).
3. Only a fraction α of individuals (I n ) in compartment I are confirmed and reported as daily new cases on day n (y n ). In our calculations in this paper, we will use the 7-day average of y n , which is written as � y n .
I n À I nÀ 1 ¼ ðS nÀ lÀ 1 À S nÀ l Þ À ðS nÀ cÀ 1 À S nÀ c Þ ð2CÞ Regarding the initial conditions of the above equations, we assume: a. (a) For all n<0 before the infectious disease starts to spread, E n = I n = R n = 0 and S n = N.
b. (b) At n = 0, the first person in compartment S is exposed to the infectious source and enters compartment E. Thus, we have S 0 = N-1, E 0 = 1, and I 0 = R 0 = 0.
Furthermore, it should be noted that the impact of different interventions on the transmission rate of infectious diseases is achieved by changing the coefficient (β n in Eq (2A)) of the transmission rate of infectious diseases. When we use l-i SEIR model to analyze COVID-19 cases in the US, β n can be considered a composite coefficient contributed from all states in the US. Thus, each state, when an intervention measure is performed, can make a certain change to the value of β n in the US more or less depending on population of the state, type of the intervention measures, and other related factors. The time-dependent value of β n in the US was determined by fitting αI n (or y n ) in Eq (2E) to the reported daily new cases (� y n ) in the US. Similarly, the time-dependent value of β n in the state of New York (NY) was also determined by fitting y n in Eq (2E) to the reported daily new cases (� y n ) in NY. To compare the l-i SEIR model with the conventional SEIR model in the simulations of daily COVID-19 cases, we assume that the rate coefficient β in Eqs (1A) and (1B) is a time dependent coefficient that can be written as β(t). Furthermore, the following equation, which is similar to Eq (2E) in the l-i SEIR model, should be added into Eq (1) for determining β(t) from the reported daily COVID-19 cases by fitting αI(t) to � y n .
In Eq (1E), t n represents the time on the nth day. Analytical solution of Eqs (2A)-(2D). From Eqs (2A) and (2B) we can derive: By observing how the red subscript "n-1" changes to "n-2" in the related terms during the above derivation, we can easily deduce the final terms in Eq (3A) when the subscript is reduced to 0. From the initial conditions, we know that E 0 = 1 and S 0 -S -l = (N-1)-N = -1 in Eq (3A). Thus, we can simplify the final expression for E n .
Similar to the derivation of expression for E n , we have: I n ¼ I nÀ 1 þ ðS nÀ lÀ 1 À S nÀ l Þ À ðS nÀ cÀ 1 À S nÀ c Þ ¼ I nÀ 1 þ ðS nÀ lÀ 1 À S nÀ cÀ 1 Þ À ðS nÀ l À S nÀ c Þ ¼ ½I nÀ 2 þ ðS nÀ lÀ 2 À S nÀ cÀ 2 Þ À ðS nÀ lÀ 1 À S nÀ cÀ 1 Þ� þ ðS nÀ lÀ 1 À S nÀ cÀ 1 Þ À ðS nÀ l À S nÀ c Þ ¼ I nÀ 2 þ ðS nÀ lÀ 2 À S nÀ cÀ 2 Þ À ðS nÀ l À S nÀ c Þ ¼ I 0 þ ðS À l À S À c Þ À ðS nÀ l À S nÀ c Þ From Eq (2D) we have: Substitution of Eq (3B) into Eq (2A) gives: Thus, the final solution of S n is the following formula: Finally, we obtain the following analytical solutions of the variables S n , E n , I n and R n from Eqs (3A)-(3D): The above closed-form solution, Eq (4), can be easily validated by comparing the calculated R n from the two l-i models when they have the same set of parameters l, i, N, and β n . From Eq (4A), we can calculate the value of S n on any day n if l and c are known and β n is a constant. Then, values of E n , I n and R n can be calculated from Eqs (4B)-(4D) after S n is determined. If β n varies with time, β n on any day n can be determined by fitting αI n to � y n based on Eq (2E) assuming that α has been determined as described previously [18]. Before starting the fitting process, we need to plot the calculated y n , y(t n ) and the reported � y n into one graph in the Excel program, and thus, we can visualize how y n or y(t n ) changes with respect to � y n while we regulate the rate coefficient β n or β(t), respectively, during the fitting process. In this case, we need Eqs (4A), (4C) and (2E) to determine S n , I n and β n , and then find E n and R n from Eqs (4B) and (4D).

Differences in the definition of rate change of variables E, I and R between l-i SEIR model and conventional SEIR model
Both the conventional SEIR model and the l-i SEIR model use the average latent period (l) and the average infectious period (i) to determine the rate of change in model variables E, I and R. The main difference between the two SEIR models is about how to use l (or σ = 1/l) and i (or γ = 1/i) to calculate the rate of change in E, I and R. This difference can be clearly seen by comparing Eq (1) with Eq (2). For example, both Eq (1B) and Eq (2B) are used to determine the rate of change of E(t) or E n with respect to time t or day n, respectively. Here, E(t) or E n is the number of people in the latent period. Both Eq (1B) and Eq (2B) have the following form: Here, v E is the rate of change in E(t) or E n , which is written as dE(t)/dt in Eq (1B) or (E n −E n −1 ) in Eq (2B). The term v E_in is the number of people who enter compartment E per unit time, which is written as bSðtÞIðtÞ N in Eq (1B) and ( b n S nÀ 1 I n N ) in Eq (2B). So far, there is no difference between Eq (1B) and Eq (2B) except that the former uses the differential form while the latter uses the difference form. The main difference between the two equations appears in the term v E_out . Based on the conventional SEIR model, if the number of people in compartment E is E (t) and the average latent period for each person is l days, then the number of people who leave compartment E per day can be defined as The method using this definition to find v E_out sounds reasonable and is easy to perform in the calculation program. However, this definition ignores temporal heterogeneity from the initial stage of infection to the recovery stage of the infectious disease among the individuals (or samples) in compartment E and I. Is there any other reasonable and easy-to-implement definition, which has taken the temporal heterogeneity of samples into account, for calculating v E_out ? To answer this question, let us make some assumptions (Fig 1). (1) There are currently 256 people in compartment E at the beginning of the morning of day 5 since the infectious disease started and spread. (2) The average latent period for each person is 4 days. (3) Among the 256 people in compartment E, 24 were infected (entered compartment E) on day 1, 56 on day 2, 80 on day 3 and 96 on day 4. Because we have known that the latent period l is 4 days and 24 people were infected on day 1, all these 24 people will transfer from compartment E (latent period) to compartment I (infectious period) on day 5 after they passed their 4 days of latent period. Thus, based on Eq (2B) of the l-i SEIR model, v E_out on day 5 is 24/day. It is easy to understand that v E_out on day 6 is 56/day because 56 infected people entered compartment E on day 2 (or 4 days ago). Simply speaking, the l-i SEIR model considers this temporal heterogeneity of samples and applies the principle of "first-in, first-out" to find v E_out . In comparison, the conventional SEIR model uses approximation of temporal homogeneity to find the rate of changes in E(t), I(t) and R(t) in (Eqs (1B)-(1D)) without considering the chronological order of samples in each compartment. The calculated v E_out from Eq (1B) of the conventional SEIR model is 64/day (256/4) on day 5. Among the 64 people, 24 (96/4), 20 (80/4), 16 (56/4) and 6 (24/4) people come from those who entered compartment E on day 4, 3, 2 and 1, respectively. It can be seen that v E_out calculated from Eq (1B) of the conventional SEIR model includes infected people in compartment E for all time frames from day 1 to day 4 because any person in compartment E has the same probability of moving from compartment E to compartment I, regardless of the chronological order of the person entering compartment E. In the above example, the number of daily new infections in compartment E increases with time from day 1 to day 4. In this rising phase of infection cases, the v E_out on day 5 calculated from the conventional SEIR model Eq (1B) is 64/day. However, only 24 people in compartment E have passed their latent period on day 5, so the v E_out determined from the l-i SEIR model Eq (2B) is 24/day. Thus, the conventional SEIR model overestimates v E_out in the rising phase of infection cases. Clearly, if the rising rate of infection cases is greater, then the calculation error of v E_out from the conventional SEIR model will be greater. Vice versa, it is not difficult to understand that if the number of daily new infections in compartment E decreases with time, the v E_out calculated from the conventional SEIR model Eq (1B) will be smaller than the v E_out determined from the l-i SEIR model Eq (2B). Therefore, the conventional SEIR model underestimates v E_out in the falling phase of infection cases. Also, if the falling rate of infection cases is greater, then the calculation error of v E_out from the conventional SEIR model will be greater. In the same way, we can discuss the differences in the definition of the rate change of variables I and R between l-i SEIR model and the conventional SEIR model.

Calculation methods and programs
Simulations and calculations of model variables S, E, I, R and y were based on the analytical solution of the l-i SEIR model Eq (4) with Eq (2E) (or Eq (2)) and the conventional SEIR model Eq (1). The parameters l, i and coefficient α in l-i SEIR model equations for the simulations of the COVID-19 transmission in the United States and NY were determined with the methods described in previous studies [17,18]. All calculation programs were written in Excel [19,20]. COVID-19 data in the US and NY obtained from Worldometer website (https:// www.worldometers.info/coronavirus/country/us/) and Wikipedia website (COVID-19 pandemic in New York (state)-Wikipedia), respectively.

Results and discussion
The conventional SEIR model can't generate propagated epidemic curves when l>i like the l-i SEIR model In Fig 2, we compared the number of infectious individuals I n calculated from the l-i SEIR model with the one calculated from the conventional SEIR model assuming l>i (Fig 2). It can be seen that the l-i SEIR model generated a propagated epidemic curve for the number of infectious individuals I n as l>i by assuming l = 1/σ = 8, i = 1/γ = 2, β n = β(t) = 1 and N = 3.3x10 8 (the solid line in Fig 2A). However, under the same conditions (l = 1/σ = 8, i = 1/γ = 2 and N = 3.3x10 8 ), the conventional SEIR model generated an epidemic curve, which increased in a near-exponential form, for the number of infectious individuals I(t) (the dashed line in Fig  2A). The one in Fig 2B shows the daily measles cases (an example of propagated epidemic curves) reported in Aberdeen, South Dakota, USA from October 15, 1970 to January 16, 1971 [21].
From Fig 2A, we can see that if l>i, then the I n and I(t) curves, representing the number of infectious individuals, calculated from the two models, are largely different. The l-i SEIR model can generate a propagated epidemic curve (solid line in Fig 2A) that is similar to the reported propagated epidemic curve in Fig 2B. In contrast, the conventional SEIR model  Fig 2A) does not generate a propagated epidemic curve even if the average latent period (1/σ = 8) is greater than the average infectious period (1/γ = 2). Why the l-i SEIR model can but the conventional SEIR model cannot simulate the propagated curve when latent period is longer than the infectious period? If we assume that the latent period l = 10 days, the infectious period i = 2 days, and that 10 infected people are all currently in the infectious period (compartment I), and they lead to 20 new people being infected during the 2 days of infectious period, then we will only see a peak with 10 cases of infectious people in the current infectious period. The reason is that the newly infected 20 people are still in the latent period, and we can only see the next peak with the 20 new cases after the 10 days of latent period. Keeping transmission in this way, we will observe a propagated epidemic curve. This result indicates that the chronological order of latent period-infectious period should be applied to epidemic models for assuring the chronological relationship among the model variables.

Calculation errors generated from the conventional SEIR model increase with the rate of change in the number of infection cases
The above theoretical analysis on the conventional SEIR model shows that if the rising rate of infection cases is greater, the calculation error of v E_out from the conventional SEIR model will be greater. To test this outcome of theoretical analysis, we examined how the differences in the model variables S, E, I and R between the two SEIR models varied with the rate of change in the number of infection cases by regulating parameters and transmission coefficients in the two SEIR models. To change the rate of infection cases in the simulations, we set β n = β(t) = 1, N = 3.3x10 8 in Fig 3A-3D Fig 3D. In Fig 3E & 3F, we set l = 1/σ = 3, i = 1/γ = 7 with varied β n and β(t). Under these conditions, both SEIR models generated sigmoidal curves for model variables S and R, and bell-shaped curves for model variables E and I. The curve S(t) (blue dashed line) is very close to the curve S n (blue solid line) as l = 1/σ = 2, i = 1/γ = 8 ( Fig  3A), and the peak positions of E(t) and I(t) are very close to the peak positions of E n and I n respectively. However, it should be pointed out that the time course of change in R(t) (purple dashed line) is significantly different from the time course of change in R n (purple solid line): R(t) increases earlier and faster than R n in the rising phase of I(t) approximately before the I(t) peak and increases slower than R n in the falling phase of I(t) approximately after the I(t) peak. This result indicates that due to ignoring chronological order, the conventional SEIR model does overestimate (or underestimate) the transfer rate from compartment E to I and from I to R in the rising (or falling) phase of I(t).
In comparison to those parameters used for simulating the curves in Fig 3A, if we change parameters l (= 1/σ) to 3 and i (= 1/γ) to 7, we can see larger differences between the model variables (S n , E n , In, Rn) in l-i SEIR model and their corresponding variables (S(t), E(t), I(t), R (t)) in SEIR model (Fig 3B). However, these apparently large differences can be partially reduced by simply moving the initial date (the day on which the first person is infected) in the conventional SEIR model forward or backward by a few days. As shown in Fig 3C, the differences in model variables (S, E, I and R) between the two SEIR models become similar to those shown in Fig 3A after postponing the initial date of the conventional SEIR model by 3 days. Carefully comparing Fig 3C with 3A, we can know that the peak height, rising rate before the peak, and falling rate after the peak of I(t) and I n in Fig 3C are smaller than those of I(t) and I n in Fig 3A. Correspondingly, R(t) is closer to R n in Fig 3C than in Fig 3A. If assuming l = 1/σ = 4, i = 1/γ = 5 and postponing the initial date of the conventional SEIR model by 6 days, the shapes of the simulated curves S, E, I and R in Fig 3D are similar to what we have seen in Fig  3A and 3C, but the peak height, rising rate before the peak, and falling rate after the peak of I

Fig 3. Examination of differences between curves S(t), E(t), I(t) and R(t) of the conventional SEIR model and their corresponding curves S n , E n , I n and R n of the l-i SEIR model at different rates of change in the number of infection cases. Simulations in (A)-(D)
were performed by assuming parameters β n = β(t) = 1 and N = 3.3x10 8 , and (A) l = 1/σ = 2, i = 1/γ = 8, (B) l = 1/σ = 3, i = 1/γ = 7, (C) parameters same as those in (B), but the initial date (the day on which the first person was infected) in the conventional SEIR model was postponed by 3 days, and (D) l = 1/σ = 4, i = 1/γ = 5, the initial date in the conventional SEIR model was postponed by 6 days. (E) Simulations were performed by assuming parameters N = 3.3x10 8 , l = 1/σ = 3, i = 1/γ = 7, and downregulating both β n and β(t) to slow down I(t) and matching I(t) and I n to each other. (F) Parameters are the same as those used in (E), but I(t) and I n were further slowed down by reducing β n and β(t).
https://doi.org/10.1371/journal.pone.0287196.g003 (t) and I n in Fig 3D further decrease compared to those in Fig 3A and 3C. Correspondingly, R (t) is closer to R n in Fig 3D than in both Fig 3A and 3C. When using the same parameters as those used in Fig 3D, but decreasing the transmission rate β, we can further reduce I(t) and I n and their rising rates as shown in Fig 3E & 3F. Correspondingly, R(t) is further close to R n in Fig 3E & 3F comparing to that in Fig 3D. To quantify these results from Fig 3A-3F, we made the following two definitions: (1) The normalized maximal calculation error in R(t), which is defined as (R(t)-R n ) max / (R n ) max . Here, (R(t)-R n ) max is the maximal difference between R(t) and R n , and (R n ) max is the maximal R n within the time range considered. The normalized maximal calculation error will be used to check the closeness between R(t) and R n calculated from the two SEIR models.
(2) The relative rising rate, which is roughly inversely proportional to the time required to go from the bottom of I(t) (at~5% of I(t)) to the top of I(t). For example, the bottom to top time of I(t) in Fig 3A is 10 days, so the relative rising rate of I(t) is 1/10 or 0.1 per day.
The factor (R(t)-R n ) in the formula given in the above definition (1) is used to obtain the difference between R(t) of the conventional SEIR model ignoring chronological order of infected individuals in the model and R n of the l-i SEIR model taking into account of the chronological order of individuals in the model. Thus, this formula is appropriate to examine the calculation errors caused by ignoring the chronological order of individuals in the conventional SEIR model. Because the chronological order of infected individuals also exists in the transmission process of infectious diseases in the real world, the calculation errors of the conventional SEIR model resulted from ignoring the chronological order of infected individuals also exist in the analysis of epidemic data in the real world.
The relative rising rates of I(t) and the normalized maximal calculation errors of R(t) obtained from Fig 3 are listed in Table 1.

Comparison of l-i SEIR model with the conventional SEIR model in simulations of COVID-19 transmission in the Unites States and in NY
To compare the two SEIR models in the simulations of epidemic data in the real world, we separately applied the differential equations of the conventional SEIR model (Eq (1)) and the analytical solution (Eq (4) with Eq (2E)) or the recursive equation (Eq (2)) of l-i SEIR model to simulate the number of daily COVID-19 cases in the United States (reported before mid-January 2021) and in NY in 2020 (reported before the end of June, 2020) (Fig 4A-4D). Because many social intervention measures were performed for reducing the rate of COVID-19 Table 1 transmission, the actual rate coefficient (β(t) and β n ) of COVID-19 transmission in the United States within certain time frames may be much smaller than 1. As a result, the rate of COVID-19 transmission (Eqs (1A) and (2A)) and the number of infectious cases (I(t) and I n ) were greatly reduced. Thus, it is interesting to examine whether we would see a significant reduction in calculation errors of the conventional SEIR model compared to the l-i SEIR model. We previously determined the parameters (l, i and α) in l-i AIR model for the simulations of daily COVID-19 cases in the United States [18] and NY [17]. Because l-i SEIR model was developed from l-i AIR model, the parameters l, i and α in l-i SEIR model are the same as those in l-i AIR model when the two l-i models are used in simulating the same set of epidemic data in the real world. In all of the simulations below, we let l = 4, i = 10 [17,18]; and let α = 0.01453 [18] for the simulation of COVID-19 cases in the US and α = 0.01176 [20] for the simulations of COVID-19 cases in NY. The time-dependent β n in l-i SEIR model was determined by regulating the value of β n day after day for fitting the calculated y n or αI n (black solid line) to the reported � y n (red dots) based on Eq (2E) (Fig 4A and 4C). To find suitable sets of parameters for the conventional SEIR model in the simulations of daily COVID-19 cases in the US and NY, we first tried to use the same sets of parameters as those used in l-i SEIR model. For example, we used l = 4, i = 10 and α = 0.01453 in l-i SEIR model for simulating daily COVID-19 cases in the US, and assumed that the first infected person in the US was exposed to COVID-19 on Feb 2, 2020. Therefore, in the conventional SEIR model, we let σ = 1/l = 1/ 4 = 0.25, γ = 1/i = 1/10 = 0.1, β(t n ) = β n , and assumed that the first infected person in the US was exposed to COVID-19 also on Feb 2, 2020. In this way, we could immediately calculate a y (t n ) curve of daily COVID-19 cases using Eq (1) in our Excel program. By visually comparing the calculated y(t n ) curve with the reported daily COVID-19 cases � y n in the US on one graph in the Excel program, we found that to match y(t n ) with � y n , we needed to postpone the date (in initial conditions), on which the first person was exposed to COVID-19 in the US, by 7 days from Feb 2, 2020 in the l-i SEIR model to Feb 9, 2020 in the conventional SEIR model. After this was done, the y(t n ) curve calculated from the conventional SEIR model was similar in shape to the curve of � y n reported from the US. Then by slightly regulating β(t n ) day by day, we were able to match the calculated y(t n ) (black dashed line in Fig 4A) with the reported � y n in the US very well. Using the same procedure, we were able to match the calculated y(t n ) (black dashed line in Fig 4C) with the reported � y n in NY very well too. From Fig 4A and 4C, we can see that the y(t n ) curve (black dashed line) simulated from the conventional SEIR model almost overlaps with the y n curve (black solid line) simulated from the l-i SEIR model.

Normalized Maximal calculation error of R (t) Data Source and parameters (l, i, β) used in simulations
After this fitting process was completed, the coefficient β n and β(t n ) were determined, and the computation program in Excel immediately obtained and showed (N-S n ), E n , I n and R n curves of the l-I SEIR model and (N-S(t)), E(t), I(t) and R(t) curves of the conventional SEIR model (Fig 4B and 4D). It can be seen from Fig 4B and 4D that the I n curve and the I(t) curve almost overlap with each other. However, the curves of R n and (N-S n ) are more or less different from the curves of R(t) and (N-S(t)), respectively. Results in Figs 3 and 4 tell us that although we can force one of 4 variables in the conventional SEIR model to be very close to its corresponding variable in the l-i SEIR model by choosing parameters and regulating the timedependent rate coefficients, significant calculation errors may still be generated in other 3 variables because of ignoring chronological order of individuals in the model.
In Fig 4D, R(t) rises earlier and faster than R n before the I(t) peak, and then slower than  Table 1. Using these data in Table 1, we plotted the normalized maximal calculation errors of R(t) vs the relative rising rate of I(t) (Fig 5). The result shows that the normalized maximal calculation errors of R(t) are linear with the relative rising rate of I(t) with a correlation coefficient r = 0.95.
The difference between (N-S(t)) and (N-S n ) is similar to the difference between R(t) and R n in Fig 4B and 4D. (N-S n ) and (N-S(t)) are the calculated total numbers of COVID-19 infections from l-i SEIR model and from conventional SEIR model, respectively. These total numbers of infections, which include asymptomatic cases and any other infected but undetected cases, on some certain dates in the US and NY were measured by antibody tests or other methods and reported in publications or public websites. The calculated (N-S n ) and (N-S(t)) in the US on September 30, 2020 are 52.5 million and 54.5 million respectively, both of which are close to 52.9 million, the real value of (N-S) in the US, which was estimated on the same date in a recent paper [22]. The calculated (N-S n ) and (N-S(t)) in NY are 2.66 million and 2.90 million on April 23, 2020 respectively, both of which are close to 2.71 million (real value of (N-S) in NY) detected by antibody tests on the same date [23]. However, the conventional SEIR model gave relatively large calculation errors both in the US and in NY. From the above data, it is easy to calculate the relative calculation errors of (N-S(t)) with respect to the real value of (N-S) in the US on September 30, 2020 (E US ), and with respect to the real value of (N-S) in NY on April 23, 2020 (E NY ). The results show that E US is 3% and E NY is 7%. Correspondingly, the relative rising rate of I(t) is 0.0118 day -1 in the US and 0.0455 day -1 in NY (Table 1). These results show that the relative calculation error of (N-S(t)) increases as the relative rising rate of I(t) increases, which is consistent with the result shown in Fig 5. Although the conventional SEIR model made some calculation errors, the normalized maximal calculation error made from the conventional SEIR model in the simulations of transmission of infectious diseases in the real world may be small. This is because the transmission rate of these infectious diseases may be small or may have been greatly slowed down by different social interventions. As shown in Fig 4B and Table 1, in the simulations of COVID-19 transmission in the US, the normalized maximal calculation error of R(t) made by the conventional SEIR model is only~6% (Figs 4B & 5). However, because the normalized maximal calculation error increases if the relative rising rate of infection cases increases as we can see in Fig  In summary, the comparison of the l-i SEIR model with the conventional SEIR model was performed in the simulations of epidemic curves. It was observed that without considering the temporal heterogeneity of infected individuals, the conventional SEIR model couldn't generate propagated epidemic curves when l>i like the l-i SEIR model, whereas the propagated epidemic curve was reported in the real world. Both theoretical analysis and simulations show that due to ignoring chronological order of infected individuals in compartments E and I, the conventional SEIR model overestimates (or underestimates) the transfer rate from compartment E to I and from I to R in the rising (or falling) phase of I(t). The calculation errors of variables in the conventional SEIR model increase with the relative rising rate of infection cases. In the simulation of daily COVID-19 cases in the US, because the rising rate is relative slow (~0.0118/day), the simulated S(t), E(t), I(t), R(t) curves from the conventional SEIR model are close to the simulated curves S n , E n , I n , R n from the l-i SEIR model. However, in the simulation of daily cases of COVID-19 outbreak in NY with a much faster relative rising rate of infection cases (0.0455/day), there are significant differences between R(t) and R n , and between (N-S(t)) and (N-S n ). These results indicate that one needs to use the conventional SEIR model with caution when simulating an epidemic curve with a fast rising rate, such as a highly transmissible infectious disease outbreak without stringent interventions to slow the spread.