This paper reports a mathematical study analyzing the daily number of COVID-19 infections during the sixth wave (January to March 2022) in 13 Japanese prefectures. Using the Gompertz distribution model, which is commonly employed in biological studies, we estimate the peak position and height of the infection time series data. We apply a linear regression approach for parameter estimation and data fitting, focusing on the relationship between the peak time and the regression time window. The results show that the Gompertz model effectively estimates the height and position of the peak based on the provided data for Japanese prefectures.

Coronaviruses are a diverse family of RNA viruses that infect humans, other mammals and avian species, including livestock. Infections caused by human and animal coronaviruses primarily result in respiratory and enteric diseases. Human coronaviruses have long been known to circulate in the population, typically causes seasonal and generally mild infections. In contrast, SARS-CoV, MERS-CoV and SARS-CoV-2 are highly pathogenic [1, 2]. Consequently, research on coronavirus epidemics has become a critical area of study.

The logistic distribution has been a widely used and important model for describing disease progression because of its effectiveness in characterizing many epidemics [3]. Originally, Verhulst proposed it to model human population growth. Interestingly, more than a decade earlier, Gompertz introduced a model for ani- mal population growth [4]. Although the Gompertz model has been used far less frequently than the logistic model, several studies highlight its usefulness [5, 6].

Time series data on COVID-19 from many countries and regions often exhibit a right-skewed distribution. Similarly, the density function of the Gompertz distribution has a right-skewed form. Consequently, several studies have applied the Gompertz function to model COVID-19 data [7, 8]. In our research, we also utilized the Gompertz distribution to analyze COVID-19 infections and deaths [9, 10].

In Japan, the surge of COVID- 19 cases from January to March 2022 has been commonly referred to as the sixth wave of the outbreak. This paper investigates the daily number of cases reported in thirteen prefectures during this period. In a previous study [11], we applied a linear regression method to estimate model parameters and fit the reported death data from foreign countries. Similarly, this paper employs a linear regression model within the framework of the Gompertz distribution. Our analysis focuses on the peak properties of COVID-19 time series data. To examine the effects of the distance between the peak date and the regression time window, we perform regression analyses via two-time windows. Addition- ally, as a special case, we conduct calculations with three-time windows for Tochigi Prefecture.

We reported preliminary results at symposium [12].

Dataset

Japan is divided into 47 prefectures with a total population of approximately 126 million people. The three major metropolitan areas are centered around Tokyo, Osaka, and Nagoya. This study focuses on 13 prefectures; Tokyo, Kanagawa and Saitama from the Tokyo metropolitan area; Osaka, Hyogo and Wakayama from the Osaka metropolitan area; and Aichi, Gifu and Shizuoka from the Nagoya metropolitan area. Additionally, Hokkaido, Hiroshima, Fukuoka and Tochigi are selected from other regions. Hokkaido, Hiroshima and Fukuoka are among the prefectures with large populations. Hokkaido is included because it is geographically isolated from the other 46 prefectures. Tochigi is chosen because of the necessity of carefully examining the relation- ship between theoretical predictions and observed data in this region. Together, the 13 selected prefectures represent approximately 60% of Japan’s total population.

We downloaded the database of the NHK website (2022); https://www3.nhk.or.jp/news/special/coronavirus/data/” nhk news covid19 prefectures daily data.csv”, downloaded on April 4, 2022.

The database contains daily numbers of COVID-19 infections and deaths from January 2020 until April 2022. We use the daily number of infections from September 1, 2021, to April 3, 2022.

The cumulative distribution function of the Gompertz distribution is defined as

F (t) = exp {−e−y(t)}, y(t) = a(t−b).    (1)

The probability density function of the Gompertz distribution is given by

f(t) = ae−y(t) F (t).                 (2)

Then we have

f(t)/F (t) = ae−y(t).  (3) 

In the analysis, we use this relation.

Let U (t) be the cumulative number of infections on the t-th day. Since reported data of daily count u(t) often fluctuate significantly, we use the 7-day moving average m(t) as the daily count. We set t = 1 at the maximum value of m(t).

For Gompertz model analysis, three parameters must be estimated: the total number N, the shape parameter a, and the position parameter b. The parameter N can be eliminated via the relations

U(t) = NF (t),          m(t) = Nf (t). (4)

Our method uses the value M (t), which is defined as

M(t) = m(t)/U (t).                  (5)

From eq. (4), the value M (t) can be estimated as

M(t) ≈ {Nf (t)}/ {NF (t)}.                (6)

Thus, we have

M (t) ≈ ae−a(t−b).                      (7)

This allows us to estimate M (t) from the reported daily numbers.

The final task is to estimatethree parameters. By applying a logarithmic transformation, we define L(t) as

L(t) = − lnM(t) ≈ − ln a + a (t − b). (8)

Thus, L(t) can be approximated by a linear function of t, allowing the parameters a and b to be estimated via the linear regression method. This study performs the regression analysisvia a time window of 12 elements, denoted as W[ts,te=ts+11]. The next step is to estimate the total number

N. We use the average ratio over the interval [ts, ts+11]

N= Ave {U (t)/F (t)}.               (9)

Details of the estimation process can be found in a previouspaper [10].

Tokyo metropolitan area

First, we analyzed the daily infection numbers in Tokyo, Kanagawa, and Saitama during the sixth outbreak. These prefectures collectively account for 24.5% of Japan’s population. Day 1 corresponds to February 5 for Tokyo, February 7 for Kanagawa, and February 8 for Saitama.

Figure 1: Results for Tokyo, Kanagawa, and Saitama. Upper panels, plots of the daily number of infections. Lower panels, plots of linear regression lines. The red solid line represents W [-11,0], and the blue dashed line represents W [-15, -4]. Closed circles indicate reported data.

Table 1 presents the estimated values for peak height and peak position using two-time windows for regression analysis: W [−11, 0] and W [−15, −4]. Note: The column labeled m represents the maximum number of reported cases, me is the estimated value of m, and b is the estimated peak position corresponding to t = 1

Table 1: Results for Tokyo, Kanagawa, and Saitama.

Figure 1 illustrates the reported data and the Gompertz model estimation for Tokyo (left), Kanagawa (center), and Saitama (right). For Tokyo and Kanagawa, the Gompertz model curves using W [−11, 0] align well with the reported data, whereas the results using W [−15, −4] show some deviation. Conversely, for Saitama, the analysis using W [−15, −4] provides a better fit to the reported data than that using W [−11, 0].

Osaka metropolitan area

This subsection presents the results of Gompertz model analysis for Osaka, Hyogo, and Wakayama. These prefectures collectively account for 12.2% of Japan’s population. Day 1 corresponds to February 8 for Osaka, February 7 for Hyogo, and February 5 for Wakayama.

Table 2 lists the numerical results of the regression analysis using W [−11, 0] and W [−15, −4].

Note: The column labeled m represents the maximum number of reported cases, me is the estimated value of m, and b is the estimated peak position corresponding to t = 1.

Table 2: Results for Osaka, Hyogo, and Wakayama.

Figure 2 illustrates the daily numbers of reported cases alongside the theoretical estimations for Osaka (left), Hyogo (center), and Wakayama (right). For Osaka, the theoretical estimation reproduces the reported data very well. Notably, the theoretical results obtained via W [−11, 0] and W [−15, −4] are almost identical in their estimates of the reported data. For Hyogo, the results of Gompertz model analysis using W [−11, 0] and W [−15, −4] demonstrate a reasonable fit to the data. Similar to Osaka, the two theoretical curves closely resemble each other. For Wakayama, the Gompertz model estimation using W [−11, 0] reproduces the daily number of infections very accurately. However, the theoretical curve generated with W [−15, −4] overestimates the reported data in the decreasing phase.

Figure 2: Results for Osaka, Hyogo, and Wakayama. The upper panels are plots of the daily number of infections. The lower panels are plots of linear regression lines. The red solid line represents W [-11,0], and the blue dashed line represents W [-15, -4]. Closed circles indicate reported data.

Nagoya metropolitan area

This subsection presents the results of the Gompertz model analysis for Aichi, Gifu, and Shizuoka. The population of these prefectures’ accounts for 10.6% of Japan’s total population. Day 1 corresponds to February 16 for Aichi, February 16 for Gifu, and February 7 for Shizuoka.

Table 3 lists the results of the Gompertz model analysis using W [−11, 0] and W [−15, −4].

Note: The column labeled m represents the maximum number of reported cases, me is the estimated value of m, and b is the estimated peak position corresponding to t = 1.

Table 3: Results for Aichi, Gifu, and Shizuoka.

Figure 3 shows the reported data and the Gompertz model estimation for Aichi (left), Gifu (center), and

Shizuoka (right). The analysis using W [−11, 0] provides good estimates for the peak height (m) and peak position (b) in all three prefectures. For Aichi, the peak at t = 1 has a bumpy structure. In the region of t ≤ −7, the analysis with W [−15, −4] offers a better fit to the reported data than W [−11, 0]. For Gifu, the infection curve exhibits two peaks; the first peak 

occurs at t = −7, and the second peak occurs at t = 1, showing a bumpy structure similar to Aichi. The analysis using W [−15, −4] appears to capture the first peak, whereas W [−11, 0] captures the second peak. For Shizuoka, Table 3 indicates that the analysis using W [−11, 0] provides reasonable estimates for the peak height (m) and the position of t = 1. However, the fitting quality of the analysis using W [-15, -4] is relatively poor, as the model underestimates the peak height.

Figure 3: Results for Aichi, Gifu, and Shizuoka. The upper panels are plots of the daily number of infections. The lower panels are plots of linear regression lines. The red solid line represents W [-11,0], and the blue dashed line represents W [-15, -4]. Closed circles indicate reported data.

Figure 4: Results for Hokkaido, Hiroshima, and Fukuoka. The upper panels are plots of the daily number of infections. The lower panels are plots of linear regression lines. The red solid line represents W [-11,0], and the blue dashed line represents W [-15, -4]. Closed circles indicate reported data.

Other regions

This subsection presents the results of Gompertz model analysis for Hokkaido, Hiroshima, Fukuoka, and Tochigi. Hokkaido is located on Hokkaido Island, which is isolated from mainland Japan, whereas Fukuoka is located on Kyushu Island, the second largest island in Japan. Hiroshima and Tochigi are located on mainland Japan. The combined population of these four prefectures accounts for 12.0% of Japan’s total population. Day 1 corresponds to February 8 for Hokkaido, January 23 for Hiroshima, February 5 for Fukuoka, and February 17 for Tochigi.

Table 4 presents the results of the Gompertz model analysis for Hokkaido, Hiroshima, and Fukuoka using W [−11, 0] and W [−15, −4].

Note: The column labeled m represents the maximum number of reported cases, me is the estimated value of m, and b is the estimated peak position corresponding to t = 1.

Table 4: Results for Hokkaido, Hiroshima, and Fukuoka.

In Figure 4, graphs of the reported data and theoretical results are plot- ted for Hokkaido (left), Hiroshima (center), and Fukuoka (right). For Hokkaido and Hiroshima, the model calculations using W [−11, 0] closely reproduce the reported data. How- ever, the calculation with W [−15, −4] for Hokkaido overestimates the reported infection curve, whereas for Hiroshima, it underestimates the reported curve. For Fukuoka, the analysis using both W [−11, 0] and W [−15, −4] provides a reasonable fit to the reported curve.

Table 5 lists the Gompertz model estimates using W [−11, 0], W [−15, −4], and W [−19, −8] for Tochigi.

Note: The column labeled m represents the maximum number of reported cases, me is the estimated value of m, and b is the estimated peak position corresponding to t = 1.

Table 5: Results for Tochigi.

Figure 5 shows the daily reported cases alongside the theoretical results for Tochigi. A double-peaked epi-demic curve is observed. The first (lower) peak occurs at t = −11, with a peak height of 819. Table 5 shows that the Gompertz model estimate me using W [−19, −8] provides a bet- ter fit than those using W [−11, 0] and W [−15, −4]. The estimates b = −12 and me = 796 using W [−15, −4] seem to correspond to the first peak.

In the period of the sixth wave outbreak, sub-emergency measures were issued for many prefectures. In Tokyo, sub-emergency measures were applied from January 21 till March 21, or −14 ≤ t ≤ 45. For Osaka, sub-emergency measures were issued in the period −11 ≤ t ≤ 27, and for Aichi, −25 ≤ t ≤ 34. For Hiroshima, sub-emergency measures were applied January 9 till March 6, or −15 ≤ t ≤ 41. We have checked the linearity of L(t) for 14 prefectures, and have investigated the relation- ship to the sub-emergency measures. However, we can find no systematic relationship between them.

Notably, our approach closely resembles the K-value model pro- posed by Nakano and Ikeda [7]. The K-value is an indicator of the spread rate and is defined as

After some modifications

K(t) = {U (t) − U (t − 7)}/U (t),

we have a relation using a seven-day moving average

Thus, we can show

K(t) ≈ 7M (t),

where M (t) is defined in eq. (5) of the Gompertz model.

The Gompertz function has been shown to fit the growth curve of chicken heart [13]. Here, we present his derivation of the Gompertz model. Let ϕ (U) be the growth rate defined

Figure 5: Daily numbers of Tochigi. The red solid line represents W [-11,0], the blue dashed line represents W [-15, -4], and the black dashed line represents W [-19, -8] Closed circles indicate reported data.

Figure 6: Logarithm of growth rate for Osaka

As

He has shown that ϕ can be approximately given by a decreasing exponential function of time

ϕ (U) = Ae−at,                     (11)

with constants A and a. Figure 6 shows the logarithm of the growth rate of infections in Osaka, where the daily count U (t +1) −U (t) is replaced by the moving average m (t + 1).

The growth rate ϕ may be approximated by

We integrate the equation

and obtain the solution with the form of the Gompertz function

where B = A/a. The constant C depends on the initial condition.

We consider the application of a flexible growth function proposed by Richards [14]. Zou et al. reported successful results when a logistic growth model was used to describe the COVID-19 outbreak in China [15]. Richards’ function can describe a wide range of distributions, including logistic and Gompertz distributions, as sub-models. The cumulative distribution functions of the logistic and Gompertz models can be transformed to the function z(t)

z(t) = α +βρt,                   (13)

where α, β and ρ are parameters. We are now analyzing the COVID- 19 outbreak data using a non-linear regression method.

We present the analysis of infection counts for 13 prefectures during the sixth wave of the COVID-19 outbreaks in Japan. This study encompasses approximately 60% of the Japanese population. Following previous research, we modeled the daily number of infections via the Gompertz distribution. Using the linear regression approach with the time window of W [-11,0], we accurately estimate the peak height and peak date of the reported infection data. However, the estimation quality for W [-15, -4] varies across prefectures. Notably, our analysis of Tochigi reveals that the time window W [-19, -8] provides better estimates of peak characteristics than W [-11,0]. This finding indicates that the estimation using the time window of the early phase must be carefully examined.

Overall, our analysis demonstrates that modeling epidemic curves within the framework of the Gompertz function holds great promise for studying COVID-19 dynamics in Japanese prefectures.

Data availability

The data generated and analyzed during the current study are available from the corresponding author upon reasonable request.

Author contributions

The present study was conducted equally by the authors.

Competing interests

The authors declare no competing interests.