Estimation of the Basic Reproduction Numbers of the Subtypes H5N1, H5N8, and H5N6 During the Highly Pathogenic Avian Influenza Epidemic Spread Between Farms

It is important to understand pathogen transmissibility in a population to establish an effective disease prevention policy. The basic reproduction number (R0) is an epidemiologic parameter for understanding the characterization of disease and its dynamics in a population. We aimed to estimate the R0 of the highly pathogenic avian influenza (HPAI) subtypes H5N1, H5N8, and H5N6, which were associated with nine outbreaks in Korea between 2003 and 2018, to understand the epidemic transmission of each subtype. According to HPAI outbreak reports of the Animal and Plant Quarantine Agency, we estimated the generation time by calculating the time of infection between confirmed HPAI-positive farms. We constructed exponential growth and maximum likelihood (ML) models to estimate the basic reproduction number, which assumes the number of secondary cases infected by the index case. The Kruskal-Wallis test was used to analyze the epidemic statistics between subtypes. The estimated generation time of H5N1, H5N8, and H5N6 were 4.80 days [95% confidence interval (CI) 4.23–5.38] days, 7.58 (95% CI 6.63–8.46), and 5.09 days (95% CI 4.44–5.74), respectively. A pairwise comparison showed that the generation time of H5N8 was significantly longer than that of the subtype H5N1 (P = 0.04). Based on the ML model, R0 was estimated as 1.69 (95% CI 1.48–2.39) for subtype H5N1, 1.60 (95%CI 0.97–2.23) for subtype H5N8, and 1.49 (95%CI 0.94–2.04) for subtype H5N6. We concluded that R0 estimates may be associated with the poultry product system, climate, species specificity based on the HPAI virus subtype, and prevention policy. This study provides an insight on the transmission and dynamics patterns of various subtypes of HPAI occurring worldwide. Furthermore, the results are useful as scientific evidence for establishing a disease control policy.


INTRODUCTION
Highly pathogenic avian influenza (HPAI) is a highly contagious viral disease that infects domestic poultry and wild birds (1). The HPAI virus can cause an epidemic that may spread rapidly, has a high mortality rate among domestic birds, and devastates the poultry industry (2). Outbreaks of distinct subtypes of HPAI, including H5N1, H5N8, and H5N6, are continually reported worldwide (3)(4)(5), and this global HPAI virus dissemination is caused by migratory wild birds (6). The HPAI crisis appears to be a great threat to not only animal health but also public health worldwide. Furthermore, the World Health Organization reported 860 human infection cases of avian influenza A subtype H5N1 (7) after the first human case of HPAI subtype H5N1, which was reported in Hong Kong in 1997 (8).
In South Korea, outbreaks of three different subtypes of HPAI occurred between 2003 and 2018. The first outbreak of H5N1 occurred from December 2003 to February 2004 and had a high mortality rate at poultry farms, especially among chickens (9). Since then, outbreaks of H5N1 have occurred in 2006, 2008, and 2010 (10)(11)(12). The novel HPAI subtype, H5N8, was first reported in January 2014 at South Korean poultry farms (13). Genetic analyses of viruses isolated from wild birds and poultry farms showed that migratory birds could be responsible for the first wave of H5N8 outbreaks between January and May 2014 (14). After the first wave, two waves of subtype H5N8 occurred during September 2014 to June 2015 and during September 2015 to November 2015 (15). It was reported that these sporadic outbreaks were caused by viruses reintroduced into Korea by migratory waterfowl (16). In November 2016, a novel genotype of H5N6 that was first detected in wild birds in Korea and HPAI infectious cases was reported at poultry farms (17). Another novel H5N8 virus co-circulated with H5N6 virus during the outbreaks in 2016, from February to June 2017 (18). In November 2017, the novel H5N6 virus was detected at a broiler duck farm and in wild mallards, with infection spreading to poultry farms (19).
The main strategies used to prevent and control HPAI outbreaks are based on the prohibition of movement, preemptive culling, and vaccinations in infected areas (20). Therefore, it is important to understand pathogen transmissibility in a population to establish an effective disease prevention policy. The basic reproduction number (R 0 ) is one of the important epidemiologic parameters necessary to understand the characterization of disease and the dynamics in a population (21). R 0 is generally defined as the average number of secondary cases caused by one infectious individual during the entire infectious period in an uninfected population (22). If each infected individual infects more than one other individual, on an average, at any time point, then the epidemic will be sustainable (23). Various methods are used to estimate the reproduction number (24)(25)(26), and these have been implemented in the R program (27) and Excel (28) as ready-made procedures.
Reproduction number estimation has been used to understand HPAI epidemic characteristics and to provide insight regarding control measures for epidemics. These farmto-farm reproduction number estimations were targeted to the HPAI subtype H5N1 and were conducted in Nigeria (29), Romania (30), Thailand (31), Bangladesh (32), India (33), Italy, Canada, and the Netherlands (34). In Korea, there was a mathematical modeling study of the reproduction number for HPAI from 2016 to 2017, but this was limited to the local reproduction number and did not include all epidemics from South Korea (35). We aimed to estimate the serial interval and R 0 of HPAI subtypes H5N1, H5N8, and H5N6, which were associated with nine outbreaks from 2003 to 2018 in Korea, and demonstrate the characterization of each subtype by analyzing HPAI characteristics, including the epidemic days, number of farms, species distribution, serial interval, and R 0 . It is expected that the results of the present study will become a foundation for demonstrating the disease dynamics of each HPAI subtype and its characteristics, as well as for establishing effective HPAI control, not only for traditional HPAI subtype H5N1 but also the emerging subtypes H5N8 and H5N6.

Data Collection
The epidemic data of HPAI outbreaks in Korea were collected by the Animal and Plant Quarantine Agency (APQA) in Gimcheon, Korea ( Table 1). In Korea, three HPAI subtypes occurred from 2003 to 2017, HPAI subtype H5N1 occurred in a total of 214 poultry farms, H5N8 occurred in 469 farms, and H5N6 occurred in total 362 farms. The livestock owner (including the manager) or veterinarian who found an animal with clinical signs and suspected HPAI was required to report the case to the APQA according to the Prevention of Contagious Animal Disease Act. Cloacal, fecal, and blood samples were collected from sick or dead poultry in reported poultry farms, and HPAI virus was confirmed using reverse-transcriptase polymerase chain reaction at the Avian Influenza Research and Diagnosis Department of the APQA. If the suspected farm was confirmed as HPAI-positive and deemed an infected premise (IP), then depopulation of farms with infected poultry and depopulation of all poultry farms in the protection zone were conducted. If a depopulated farm was found to be positive, then it was defined as a positive premise (PP) (36). Both IP and PP were considered cases in this study. The epidemic curve of these HPAI cases was depicted using the "incidence" package in R (37) to illustrate the weekly reported number of poultry farms in the International Organization for Standardization (ISO) week date system (37) (Figure 1). In Korea, there were no poultry farms infected with two HPAI subtypes simultaneously, and each farm only had one subtype in each outbreak.
Based on the APQA epidemiology reports, the HPAI outbreaks were classified as waves when the period between cases was longer than 1 month (38). As a result of this classification, outbreaks of the subtype H5N8, which occurred in 2014 and 2016, were classified as four and two waves, respectively. Four outbreaks, including the H5N1 outbreak in 2003, outbreak in 2006, the fourth wave of the H5N8 outbreak in 2014, and the second wave of the H5N8 outbreak in 2016, were excluded from the analysis because the samples were too small to calculate R 0 .

Serial Interval and Generation Time
A serial interval is the time between successive cases in a chain of transmission, estimated from the interval between clinical onsets in patients (25). We estimated the serial interval of HPAI as the time between the reported date of the first farm with infected cases and secondary farm with infected cases. This estimation was based on the investigation of the epidemic pathway of HPAI transmission, which shows the epidemiologic relationship between the infector and infectee. According to the APQA investigations, HPAI transmission could be possible through wild migratory birds, wild animals, farm owners, managers, staff, vehicles related to the poultry industry, and airborne transmission from nearby infected farms. The epidemic transmission pathway investigation was conducted by an APQA epidemiologic investigator visiting and interviewing the places suspected to be associated with the infected farms, including animal facilities such as hatcheries, feed factories, and live bird markets. The APQA investigated vehicles, people, livestock, and their products that entered an infected farm from 21 days prior to infection and estimated the disease transmissions.
In addition to investigating via interview, the APQA used geographic information to identify HPAI viral transmission by vehicles. In Korea, vehicles related to the poultry industry transporting poultry, poultry products, medicines, feed, and feces must be registered with the Korea Animal Health Integrated System (KAHIS; http://www.kahis.go.kr). The movements of livestock-related vehicles are reported to the KAHIS, making it possible to track the movement of vehicles, people, livestock, and animal products.
Through these interviews and vehicle information, the disease transmission pathway via transportation and human movement was identified. If a clear epidemiologic link to the infected farm could not be found through interviews and movement tracking, then we hypothesized that the farm might have been infected with HPAI by wild migratory birds or wild animals. We then excluded infection thought to be caused by wild birds or wild animals during the estimation of the serial interval because it is not possible to observe the serial interval of virus transmission from wild birds and animals.
The generation time is the modeling term describing the time duration from the onset of transmittable infection in a primary case to the onset of infection in a secondary case infected from the primary case. We defined the generation time as the difference between suspected infection days of the primary farm and secondary farm, which was measured through epidemiologic investigation (Figure 2). The suspected infection day was estimated according to the day reported by the farm owner after clinical symptoms were found in the poultry and the period between the infection and latent period of each HPAI subtype in the poultry species. We estimated the suspected infection date from the day the clinical symptoms were reported by subtracting the periods between infection and clinical symptoms. For H5N1, the periods between infection and clinical symptoms were assumed to be 2 days for chickens (9), 4 days for ducks (39), and 3.8 days for other poultry species (9). For H5N8, the periods were 3.2 days for chickens (40), 8.0 days for ducks (15), and 2.0 days for other species (41). For H5N6, the periods were 2.6 days for chickens (42), 4.6 days for ducks (43), and 3.0 days for other species (43).
Based on the generation time between case farms, we calculated the discretized generation time distribution using a function (est.GT) in the R0 package (27). Discretization is performed on the grid [0, 0.5), [0.5, 1.5), [1.5, 2.5), etc. . . where the unit is time interval of days (27). Time-to-event data were assumed to follow a parametric distribution with a probability density function (PDF). The distribution of generation time is expressed in the form of parametric distribution such as "gamma, " "lognormal, " or "Weibull, " using maximum likelihood. The mean and standard deviation of generation time is provided in the desired time units. The calculated distribution of the generation time in each subtype and outbreaks is depicted in Figure 3.

Model Assumption and Data Analysis
The study model is based on the susceptible-infected-removed (SIR) compartmental model (44), which divides poultry farms into compartment. A susceptible farm (S) becomes infectious (I) through contact with the possible disease and is then removed (R) by depopulation. The dynamics of an epidemic can be described as the equation given below when N is the sum of S, I, and R.
In this model, β is a parameter, which controls how much the disease can be transmitted through the exposure of HPAI virus, and γ is a parameter, which expresses how many poultry farms can be removed in a specific period. In this model, the average number of secondary infections caused by an infected host, R 0 , equals β/ γ (45). We constructed exponential growth (EG) and maximum likelihood (ML) models to estimate early reproduction numbers using the R0 package (27) in R (version 3.3.0). The EG model assumes that the initial reproduction ratio can be associated with the EG rate during the early epidemic phase (24). The formula is R 0 = 1/M (-r), where r denotes the initial EG rate and M stands for the moment generating function of generation time distribution. In the initial EG model, a period from day 1 to day 14 of the epidemics was chosen when the outbreak's growth was exponential. The 14-day interval was selected based on Korea's standstill policy (38). When HPAI outbreaks are reported in South Korea, a standstill policy is implemented for the vehicles, to reduce the spread of HPAI. This policy is intended to minimize the contact between vehicles and suppress the HPAI dissemination. We determined that this intervention affects the basic reproduction number of HPAI. Therefore, we specified the exponential growth for the first 14 days of the epidemic wave in each case.
A function (est.R0.EG) in the R0 package was used (27). We used a sensitivity test in EG to select the period during which growth is exponential as optimized time windows. We used the "sensitivity analysis" function to compute the deviance R-squared statistic over a range of periods.
The ML estimation model assumes that the number of secondary cases caused by an index case is Poisson-distributed with the expected value R 0 (25). The log-likelihood (LL) of R 0 was defined as LL This model assume that the number of new cases at indexing time t as N= {N t }, t = 0,. . . T, a generation time distribution w, and µ t which represent the total number of cases produced by the earlier case N t . The likelihood must be calculated on a period of exponential, and the deviance R-squared measure may be used to select the best period that maximized the likelihood. In this study, the range was set as 0.01-50, in which the maximum must be searched. A function (est.R0.ML) in the R0 package was used (27). The goodness of fit of each model was calculated using the chi-square goodness of fit test in R.
The Kruskal-Wallis test was used to determine the statistical differences in epidemics between subtypes (46). The epidemic days, number of farms, cases per day, poultry species distribution of farms, generation time, and R 0 estimated by EG and ML of the three subtypes H5N1, H5N8, and H5N6 were analyzed. The significance level was α = 0.05. These statistical analyses were performed using SPSS 22.0 (IBM, Armonk, NY, USA).

HPAI Epidemic in Korea
We investigated 12 HPAI outbreaks of three subtypes, H5N1, H5N8, and H5N6, that occurred from 2003 to 2018 in Korea. Table 1 presents a summary of the epidemic data, including the period of outbreaks and the number of infected farms that were investigated. The weekly epidemic curves of HPAI outbreaks are shown in Figure 1 based on the ISO 8601 week date system. The H5N1 HPAI outbreaks (except for the 2008 outbreaks) began between November and February, when the lowest temperature drops below 0 • C ( Figure 1A). Regarding H5N8 in 2014, the second and third waves recurred in September 2015 and 2016, respectively ( Figure 1B). However, the second wave of H5N6 in 2016 occurred in June (Figure 1C). The longest outbreak was the second wave of H5N8 in 2014, which occurred over 260 days. The shortest outbreak was the fourth wave of H5N8 in 2014, which occurred over 14 days. The outbreaks with the most cases (340 poultry farms) and cases per day (3.579 cases per day) were the H5N6 outbreaks in 2016. Regarding H5N8 in 2014, more than 72% of the occurrences were in ducks; however, there was no apparent species specificity for subtypes H5N1 and H5N6.

Serial Interval and Basic Reproduction Number of HPAI in Korea
We selected nine outbreaks with sufficient number of premises to calculate R 0 and analyzed the generation time and initial R 0 using the EG and ML methods ( Table 2). Generation time distributions are illustrated by each HPAI subtype as the PDF in Figure 3. Generation time of H5N1 were estimated between 4.58 and 5.24 days (Figure 3A), generation time of H5N8 were estimated to have 6 days or more (6.01-8.23 days) (Figure 3B), and generation time of H5N6 were estimated between 5.02 and 5.91 days (Figure 3C). R 0 was estimated as 1.65-2.20 for subtype H5N1, 0.03-1.56 for subtype H5N8, and 1.03-1.24 for subtype H5N6 using EG methods. Using ML methods, R 0 was estimated as 1.68-1.95 for subtype H5N1, 1.03-1.83 for subtype H5N8, and 1.37-1.60 for subtype H5N6.
Most of the R 0 in the EG and ML methods were similar, except for the second and third waves of H5N8 in 2014. The R value obtained by the EG method was <1 for the second and third waves of H5N8 in 2014. To select the optimal time windows, sensitivity results of the time windows and R 0 were used ( Table 2). Optimized time windows selected by sensitivity tests accounted for 69.14% of the outbreak periods, on an average, and the optimal R 0 values in optimized time windows were <1 for the subtypes H5N1 and H5N6 outbreaks.

Epidemic Statistics Between Subtypes
The average values of the number of epidemic days, infected poultry farms, species distribution, and infected farms per day for the three subtypes of nine selected outbreaks were determined ( Table 3). The average numbers of epidemic days were 86.0 for H5N1, 108.0 for H5N6, and 143.8 for H5N8. The average numbers of farms were 69.0 for H5N1, 107.8 for H5N8, and 181.0 for H5N6. Regarding the species distribution, subtype H5N8 was more highly distributed among duck farms (74.2%) than other subtypes (37.7% for H5N1 and 42.5% for H5N6).
The Kruskal-Wallis H test showed a statistically significant difference in mean generation time among the different subtypes [χ 2 (2) = 6.444; p = 0.040], with mean rank scores of 2.33 for subtype H5N1, 7.50 for H5N8, and 4.00 for H5N6. The pairwise comparison showed that the mean H5N8 generation time (7.58 days) was significantly longer than the H5N1 generation time (4.80 days) (P = 0.03) ( Table 3). There were no significant differences among subtypes in epidemic days, number of farms, cases per day, species distributions, or reproduction number.

DISCUSSION
HPAI outbreaks occur continually worldwide and have become a major threat to animal and human public health. In South Korea, eight outbreaks with multiple waves of infections occurred between 2003 and 2018; these involved three different HPAI subtypes, H5N1, H5N8, and H5N6, and massively damaged the poultry industry. Therefore, it is important to understand the HPAI transmissibility at poultry farms to control outbreaks by establishing an effective prevention policy. An effective tool for understanding disease characteristics is the R 0 , which is generally defined as the average number of secondary cases caused by one infected individual (21). Therefore, we investigated the transmission dynamics of the HPAI subtypes H5N1, H5N8, and H5N6 by estimating the generation time and R 0 . To the best of our knowledge, no previous study has attempted to estimate R 0 of various HPAI subtypes and perform comparative analyses among them. This could be the first study to investigate the disease transmission dynamics of HPAI subtypes H5N1, H5N8, and H5N6, which are emerging worldwide.  The R 0 of HPAI H5N1 in Korea estimated in this study was between 1.68 and 1.95, according to the ML method ( Table 1). The R 0 of subtype H5N1 has previously been estimated in countries such as Italy (1. Despite being the same subtype of HPAI, the estimated R 0 subtype H5N1 varied across countries. We assumed that several factors, such as geographic distribution of poultry farms, mixed farming systems, poultry product supply system, and climate, were associated with this difference. We believe that unique characteristics of the poultry industry in Korea and climatic differences are the major causes for these observed differences. We speculate that the estimated R0 may be related to characteristics of the Korean poultry industry, such as the coexistence of large-scale commercial farms and small family farms. Among the Organization for Economic Cooperation and Development (OECD) countries, Korea has the lowest availability of arable land per capita (0.03 hectare in 2016) (48). This land scarcity is an important factor leading to high stocking densities (49). A previous study suggested that farms with large flocks and the presence of a neighboring farm within 500 m were risk factors of HPAI at Korean broiler duck farms (50). This high stocking and local density of large-scale poultry farms could increase the likelihood of massive infections when HPAI outbreaks occur in Korea.
Small family poultry farms also represent a biosecurity risk during HPAI outbreaks. Most of these small farms sell live poultry to local markets without going through slaughterhouses; this could be a pathway for the spread of HPAI viruses. Additionally, there was an obvious lack of information regarding the official statistics of poultry farms too small to be defined as agricultural holders in Korea (51). This includes establishments with <0.1 hectares of land or with sales of agricultural products per year or value of agricultural animals less than KRW 1.2 million (USD 1,090).
Secondly, we hypothesize that climate factors during the epidemic period may affect R 0 in these countries. Climate factors could affect HPAI transmission and persistence by altering bird migration, virus shedding between hosts, and virus survival outside the host (52). Climate change is considered to influence the wild bird species composition and their migration cycle, and these changes will affect the transmission intensity of disease (53). Furthermore, temperature and humidity could be related to viral persistence in the host and environment. An influenza virus transmission experiment using a guinea pig model suggested that relative low humidity and cold temperature were favorable for spreading influenza (54). Liu et al. (55) showed that the environmental temperature decreased shortly before HPAI H5N1 outbreaks in domestic poultry in Eurasia between 2005 and 2006. Additionally, AI viral infectivity remained at lower temperatures (<17 • C) during an in vivo test (56). Therefore, it is assumed that our estimated R 0 in Korea is higher than the R 0 in Thailand and Bangladesh, where the average annual temperatures and humidity are higher. Based on these results, we assumed that the climate factors were closely related to the R 0 estimated in several countries in terms of virus transmission and survivability.
In 2016, two novel HPAI subtypes, H5N6 and H5N8, occurred simultaneously. HPAI H5N6 occurred from November 2016 to February 2017, whereas subtype H5N8 occurred from February to April 2016; the first wave and second wave occurred in June. Although these two subtypes occurred simultaneously, both were novel viruses. The genetic clade analysis suggested that Korean H5N6 viruses are novel reassortments of multiple virus subtypes, and it is difficult for H5N6 virus reassortment to occur during outbreaks that could increase the possibility of viral subtype mutation (5). Additionally, an infection experiment involving wild mandarin ducks demonstrated a difference in viral shedding and viral tropism in H5N8 and H5N6 viruses within the same clade of 2.3.4.4 H5 HPAI viruses (57). Based on these findings, both subtypes were independent of each other, and the virus infectivity could also be different; therefore, different R 0 was expected.
However, our estimated initial R 0 value in 2016 suggested a similarity between the reproduction number represented in subtypes H5N8 (1.70) and H5N6 (1.60) ( Table 2). Apart from the difference in transmissibility of each virus subtype, the level of transmission between farms in the field may be similar between the two subtypes. However, this presumes that the values of R 0 of the two subtypes were similarly calculated because the biosecurity policy implemented during the outbreaks was identical. The basic reproductive number is affected by the rate of contacts in the host population, the probability of infection being transmitted during contact, and the duration of infectiousness (58). Therefore, it can be estimated that the R 0 of two different subtypes were similar due to the reduction of the poultry population through preemptive culling and the reduction of contact between farms because of the standstill (59).
The quarantine against HPAI in Korea has changed over 14 years after the first HPAI epidemic in 2003. The HPAI prevention policy changed dramatically, especially before and after H5N8 epidemics in 2014. Before the outbreaks, Korea Animal Health Integrated System (KAHIS) was established in 2013 to monitor livestock vehicle movement. In this system, all poultry-related vehicles must be registered with KAHIS and equipped with a global positioning system mandatorily (60). Also during the epidemics, the preemptive depopulation of the protective zone was changed from a radius of 500 m−3 km, and inspections were conducted more than once before releasing poultry and poultry products (36), The influence of these quarantine policy can also be seen in the changes in the R 0 values of each wave of subtype H5N8 that occurred between 2014 and 2016. For H5N8 in 2014, the initial R 0 of each wave showed a tendency to decrease as the outbreak progressed gradually ( Table 2). This would indicate that the effectiveness of control measures for HPAI were increasing while the waves were passing.
In the Kruskal-Wallis model, H5N1 and H5N8 subtypes showed statistically significant differences in generation time (P = 0.03) ( Table 3). However, there were no significant differences in the epidemic characteristics of the subtypes. There was also no statistical significance in the R 0 obtained through the EG and ML models. This generation time difference in the two subtypes might be associated with subtype pathogenicity in the poultry species. The spread of H5N1 viruses in the field was quickly controlled as a result of the rapid diagnosis of the infections due to the high pathogenicity of these viruses in poultry. In contrast, subtypes H5N6 and H5N8 clustered as clade 2.3.4. H5NX viruses are usually mild in ducks, leading to delayed diagnosis of infections and persistent spread in the wild (61). Therefore, the H5N8 subtype could possibly spread the HPAI virus over a relatively longer period than the H5N1 subtype which could be driven by sub-clinical spread in ducks.
In conclusion, this study showed the characterization of each subtype by analyzing the HPAI characteristics, including the epidemics, number of farms, species distribution, generation time, and R 0 of HPAI subtypes H5N1, H5N8, and H5N6, which were associated with nine outbreaks in Korea between 2003 and 2018. R 0 , which is estimated by the generation time, index case, and secondary cases, is essential for identifying the characteristics of HPAI. In particular, our findings suggest that the estimated R 0 might be influenced by the HPAI subtype and might be associated with the seasonal aspects during the early stage, species specificity by virus subtype, and prevention policy. We believe that the results of the present study are helpful for demonstrating the disease dynamics of each HPAI subtype and its characteristics and, thus greatly assist in better disease control strategies. It could be possible to establish systematic quarantine policies to reduce the socio-economic losses caused by HPAI, Especially differences observed between countries with different poultry raising systems and climatic conditions. This study provided insight regarding HPAI transmission of the traditional subtype H5N1 and newly emerging subtypes H5N8 and H5N6. Further research on the basic reproduction numbers of the HPAI subtypes occurring worldwide is required to understand the global dynamics of HPAI transmission.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author/s.

AUTHOR CONTRIBUTIONS
W-HK designed the study, investigated data collection, reviewed the data, performed data analysis, and participated in manuscript preparation. SC supervised the project, administrated the project, acquired funds, and participated in the manuscript review.