Transmissibility of hand, foot, and mouth disease in 97 counties of China

Hand, foot, and mouth disease (HFMD) is a serious disease burden in the Asia–Pacific region, including China. This study calculated the transmissibility of HFMD at county levels in Jiangsu Province, China, analyzed the differences of transmissibility and explored the possible influencing factors of its transmissibility. We built a mathematical model for seasonal characteristics of HFMD, estimated the effective reproduction number (Reff), and compared the incidence rate and transmissibility in different counties using non-parametric tests, rapid cluster analysis and rank-sum ratio in 97 counties in Jiangsu Province from 2015 to 2020. The average daily incidence rate was between 0 and 4 per 100,000 people in Jiangsu Province from 2015–2020. The Quartile of Reff in Jiangsu Province from 2015 to 2020 was 1.54 (0.49, 2.50). Rugao District and Jianhu District had the highest transmissibility according to the rank-sum ratio. Reff generally decreased in 2017 and increased in 2018 in most counties, and the median level of Reff was the lowest in 2017 (P < 0.05). The transmissibility was different in 97 counties in Jiangsu Province. The reasons for the differences may be related to the climate, demographic characteristics, virus subtypes, vaccination, hygiene and other infectious diseases.

Fitting results of SEIAR model of HFMD in Jiangsu Province from 2015 to 2020. The fitting results of the daily incidence rate over time in the 97 counties in Jiangsu Province from 2015 to 2020 are shown in Fig. 2, 3,4. The correlation analysis between the fitting value and the actual reported value showed that the mean of the coefficient of correlation R 2 was 0.50 ± 0. 15, showing that the model was fitted well (Table S1).

Transmissibility of HFMD in Jiangsu Province from 2015 to 2020. The Quartile of R eff in Jiangsu
Province from 2015 to 2020 was 1.54 (0. 49, 2.50), the 95% reference range was less than 5.88, and the highest R eff could reach 20,000 times the lowest. R eff showed a periodic change in the unit of year, and there was at least one R eff peak in the adjacent years, with the peak greater than 1.0. The median of R eff for each year from 2015 to 2019 was different (χ 2 = 21.283, P = 0.000), and that the median of R eff in 2017 was smaller than that in other to Southern Jiangsu, including Changzhou city No. 1-6 (Jintang district, Suyang county-level-city, Tianning district, Wujin district, Xinbei district, Zhonglou district, respectively), Nanjing city No. 7-17 (Gaochun district, Gulou district, Jianye district, Jiangning district, Liuhe district, Pukou district, Qixia district, Qinhuai district, Sushui district, Xuanwu district, Yuhuatai district, respectively), Suzhou city No.18-27 (Changshu county-levelcity, Gongyeyuan district, Gusu district, Huqiu district, Kunshan county-level-city, Taichng county-level-city, Wujiang district, Wuzhong district, Xiangcheng district, Zhangjiagang county-level-city, respectively), Wuxi city No. 28-34 (Binhu district, Huishan district, Jiangyin county-level-city, Liangxi district, Xishan district, Xin district, Yixing county-level-city, respectively) and Zhenjiang city No. 35-40 (Dantu district, Danyang countylevel-city, Jingkou district, Jurong county-level-city, Runzhou district, Yangzhong county-level-city, respectively). www.nature.com/scientificreports/ years (P < 0.05). The median of R eff in Southern Jiangsu was the smallest among the three regions (P < 0.05), and that the median R eff of Changzhou City was lowest compared with other four cities in Southern Jiangsu (P < 0.05) (Fig. 5). The median of R eff of Yancheng City was the highest compared with other four cities in Northern Jiangsu (P < 0.05) (Fig. 6). According to the periodic change of R eff in different counties, 97 counties could be divided into five types as the following: (1) R eff was basically at a high level and remained above 1.0 with a cyclical change, represented by Huqiu Comparison of the transmissibility in the 97 counties. We compared the transmissibility of 97 counties using RSR. According to the RSR distribution table (Table 1), we constructed the RSR and probit regression equations which could be obtained as RSR = −0.261 + 0.151 × Probit (F = 1813.37, P = 0.000), and the RSR of each district was calculated and used to classify the transmissibility into six levels, as showed in Table 2. Trans-  Sensitivity analysis. We selected a period of incidence data to fit curve with 10 values of parameter κ in the range of 0-1. The result showed that the fitting values had a high degree of coincidence, which indicated that the SEIAR model was not sensitive to this study (Fig. 8).

Discussion
In this study, the seasonally adjusted SRIAR model was used to study the transmissibility of HFMD among the 97 counties in Jiangsu Province, to provide suggestions for local Centers for Disease Control and Prevention, community in Jiangsu Province and other areas with a similar transmissibility of HFMD.
Analysis of the different incidence rates and transmissibility in three regions. The incidence rate of Southern Jiangsu was higher than Northern Jiangsu and peaked in two seasons per year, which was consistent with earlier studies of HFMD in Jiangsu Province 21 . In this study, we found that some counties of Northern Jiangsu had one seasonal peak, two seasonal peak or more than two seasonal peaks. The average R eff of HFMD in Jiangsu Province from 2015 to 2020 was 1.54, which was similar to the research results of studies on foreign and most domestic provinces and regions; while the R eff was lower than that of Shenzhen, Guangdong Province 18 . We found that the R eff in Southern Jiangsu was less than that in Northern Jiangsu, which was contrary to the incidence rate of the regions. We considered that the reasons for the different incidence rates, seasons and transmissibility in three regions were as following: (1) The climate zones of the regions were inconsistent may be influence the spread of HFDM virus and the epidemic virus serotypes. The different pathogens would compete with each other, resulting in the change in seasonal peaks and transmissibility. Liu et al. found that the incidence rate of HFMD in Jiangsu was proportional to the average temperature and rainfall but was negatively correlated with the days of rainfall (≥ 0.1 mm), low temperature, high temperature, and sunshine duration 21 .
Moderately warm environment was found to promote the spread of the HFMD virus 22 . Our team previously used mathematical models to analyze the interaction of the main pathogens of Changsha HFMD and found that different pathogens would compete with each other, resulting in changing transmissibility and an impact on morbidity 23 . In recent years, studies found different disease prototypes in different urban areas of Jiangsu Province and determined that different dominant strains appear alternately 24,25 . We considered that Southern Jiangsu was warmer than Northern Jiangsu, which may have caused the incidence rate to be higher, and the winter in the northern region was too cold to prevent the spread of HFMD. The difference viruses between Northern Jiangsu and Southern Jiangsu, may be affected by the growth environment, and the transmissibility be affected by the interaction between viruses. (2) The demographic characteristics and the hygiene levels of the regions were different. The south of Jiangsu Province was a densely populated area, which had higher birth rates. Studies found that most of the infected patients in this area were infants and children under five years of age. While the   [26][27][28] . The highest development of medical and health undertakings was located in Southern Jiangsu (Nanjing, Suzhou, and Wuxi), and the lowest health environment was mainly in Northern Jiangsu (Lianyun, Yancheng, and Suqian). Although the population base was larger in Southern Jiangsu which could have enabled greater opportunities for population contact and communication, economically developed areas had higher medical and health levels, higher health awareness, and a stronger ability to block the transmission of HFMD, making the transmissibility of Southern Jiangsu relatively weaker than that of Northern Jiangsu. (3) There was a small range of outbreaks. In Northern Jiangsu, for cities with low incidence rates, such as Yancheng (Binhai), the small number of incidences induced data instability in the transmission dynamics model analysis. The uncertainty of the results increased, coupled with the region with a long-term low incidence rate, poor health level, and insufficient experience in HFMD outbreaks, resulting in a strong transmissibility during the outbreak. Therefore, we thought that there may have been small ranges of outbreaks in Northern Jiangsu, some disadvantages of the kinetic model itself, and a weaker ability to cope with the outbreak, making the R eff in the region higher. (4) Emerging contaminants influenced the antibody   29 . A study showed that cord serum long-chain PFASs concentrations significantly correlated with low antibody levels of CA16 and EV71 at three months 30 . A high concentration of PFOS was found in the Yangtze Estuary sediments (72.9-536.7 ng/g) 31 . Southern Jiangsu is near the Yangtze Estuary, the immunity of infants with HFMD may be reduced by fluoride pollutants, leading to a higher probability of incidence. Based on the aforementioned analysis, we suggested that Southern Jiangsu should pay more attention to a wide range of public health publicity during the seasons of HFMD onset. For the central and Northern Jiangsu areas with strong transmissibility of HFMD, improving health services, strengthening investments in healthcare and implementing protective measures that are more helpful in reducing its prevalence are necessary. A study showed that enteroviruses were detected on the surface of environmental items in hospitals and community public playgrounds in areas with a high incidence of HFMD in Jiangsu Province 32 , thus strengthening the need for disinfection during high-incidence periods. Communities are areas prone to cross-infection and should be taken seriously.   33,34 . Because we found that the incidence rate was also the lowest in 2017, the immunity provided by the vaccine against EV71 and the publicity of vaccination reduced the number of susceptible and infected people, thereby reducing the actual transmissibility of HFMD. It is interesting that in many counties, the incidence rates and transmissibility suddenly increased in 2018, and the peak height could have been higher than that in 2015 and 2016. We analyzed the different incidence rates and transmissibility in various years based on following aspects: (1) From the perspective of climate change, although the incidence and spread of HFMD are related to climate factors 21 , according to some meteorological studies, the temperature and rainfall in 2017 and 2018 were not abnormal compared with those in other years 35,36 . Therefore, compared to 2017, the average incidence rate in 2018 increased significantly, which may not have been related to climatic factors. (2) (3) Regarding the impact of other infectious diseases, in this study, we found that the average daily incidence rate of HFMD in the first half of 2020 was 10 times lower than that before. This indicated that the protective measures against coronavirus disease (COVID-19), such as school closures, business discontinuation, frequent hand washing and wearing of masks, and maintaining social distance, have affected the prevalence of HFMD to some extent. Other research also showed that the incidence rate of HFMD was affected by road passenger volume and population mobility during the school terms and Spring Festival. The combined effect was more significant than that of meteorological factors on the epidemic  www.nature.com/scientificreports/ of HFMD every year in Jiangsu province was in spring, summer, and autumn, with better air quality. We thought that the impact of air quality on HFMD in different years in Jiangsu is not clear. Therefore, we suggest that based on the classification of the different transmissibility described by results, specific counties should be selected to monitor the subtypes of HFMD, and an HFMD vaccine for different subtypes should be developed to cope with the changes in epidemic pathogens.

Analysis of the comprehensive comparison results of transmissibility in the 97 counties. Jianhu
District in Northern Jiangsu and Rugao District in Central Jiangsu had the strongest comprehensive evaluation of transmissibility, though the trends of transmissibility were different. Jianhu District maintained a high transmissibility from 2015 to 2016 and had a downward trend from 2017 to 2019. Based on previous research, we found that in Yancheng City, where Jianhu is located, HFMD was highly prevalent among infants and that the higher the level of maternal antibodies against to EV71, the stronger the protection for infants 53 . Therefore, the implementation of vaccine immunization in Jianhu District provided a certain protective effect. However, it is still necessary to further strengthen the propaganda and education and detect whether there is a new virus subtype epidemic. The transmissibility in Rugao District had an uptrend trend from 2017 to 2020, with research showing that the HFMD in Rugao District had been more serious in recent years, and the incidence rate and incidence ratio were the highest fpr class C infectious diseases 54 , with critically ill infected patients in 2015 to 2020. There were inappropriate nursing practices and poor health environments in the rural areas in this district, and the number of vaccination individuals was low. We need to focus on improving the health environment, strengthening public health marketing and healthcare education, improving the awareness of epidemic prevention, and improving epidemic situation monitoring, especially the analysis and monitoring of virus subtypes of severe patients.

Limitations.
Owing to shortcomings in the data, this study had some limitations. In this model, factors that may affect the transmissibility, such as age and gender, were not included, which may have impacted the results. The actual data of possible influencing factors, such as climate characteristics, virus types, population data, were  www.nature.com/scientificreports/ also not collected for a correlation analysis with transmissibility of HFMD in various counties. We did not consider that the infected individuals were also affected by other bacterial or viruses, which affected the course of HFMD. Additionally, we did not consider that these infectious individuals could be reinfected with HFMD after recovering from HFMD.

Conclusion
1. The epidemic situation of HFMD in Jiangsu Province in 2015-2019 was more severe than in 2009-2013. The impact of COVID-19 was related to a reduction in the HFMD incidence rate in Jiangsu Province in 2020. 2. The outbreaks and transmissibility of HFMD in Jiangsu had regional and seasonal characteristics. The higher the incidence rate in the three regions, the lower the transmissibility. The peak period of the epidemic changed from season to season. 3. The differences in the incidence rates and transmissibility of HFMD in Jiangsu Province were related to the climate, population, virus subtypes, vaccination, and other infectious diseases; the difference in virus subtypes may have been the most important factor. 4. Rugao District in Central Jiangsu and Jianhu District in Northern Jiangsu had the strongest transmissibility of HFMD among these 97 counties in Jiangsu Province. The vaccination rate should be increased in Jianhu District, and health marketing, health conditions, and virus subtype monitoring should be reinforced in Rugao District. 5. The transmissibility was similar between some cities or regions, suggesting that the representative areas should be selected for virus subtype surveillance according to the characteristics of transmissibility in Jiangsu Province.

Methods
Ethics declarations. This efort of disease control was part of CDC's routine responsibility in Jiangsu Province, China. Therefore, institutional review and informed consent were not required for this study. All data analyzed were anonymized, and does not contain any personal privacy or identity information, so the ethics approval documents may be exempted.
Data sources. Jiangsu Province is located in the eastern coastal area of the Chinese mainland. From Southern Jiangsu to Northern Jiangsu, the climate transition from subtropical zone to warm temperate zone, with mild climate and moderate rainfall. HFMD is a class C legal infectious disease in China. Doctors must report these cases with suspected HFMD, including suspected cases, clinical cases, and experimental diagnostic cases within 24 h to the network direct reporting system for the monitoring information of statutory infectious diseases. The data of the HFMD cases used in this study were obtained from the China Information System for Disease Control and Prevention, including the number of cases, deaths reported daily and date of onset. The case types included clinical diagnosis and laboratory diagnosis cases. Demographic information was obtained from the statistical yearbook of Jiangsu Province, including the number of permanent residents at the end of the year, birth rate, and mortality rate. According to the statistical yearbook of Jiangsu Province, Jiangsu Province is divided into three regions. The three regions are bounded by the Huaihe River and irrigation canal, with a subtropical humid monsoon climate in the south and a warm temperate humid and semi-humid monsoon climate in the north. The size of three regions is similar, but there are obvious differences in social and economic development: the regional economy is the highest in Southern Jiangsu, followed by Central Jiangsu and then Northern Jiangsu 13 .
Case definition. The diagnosis of HFMD was carried out according to the guideline issued by the National Health and Family Planning Commission of the People's Republic of China 14 . Patients with HFMD, whether probable or confirmed, were classified as severe if they had any neurological complications (aseptic meningitis, encephalitis, encephalomyelitis, acute flaccid paralysis, and autonomic nervous system dysregulation), or cardiopulmonary complications (pulmonary edema, pulmonary haemorrhage, and cardiorespiratory failure), or both; otherwise, patients were categorized as mild cases.
According to the diagnostic criteria (2018 version) of HFMD 14 , the confirmed cases were confirmed by via enzyme-linked immunosorbent assays (ELISA), reverse transcription polymerase chain reaction (RT-PCR), real-time PCR, or virus isolation.
The transmission models of HFMD. According to the epidemiological feature of HFMD and our previous studies 6,7,9 , the SEIAR model could be used for the simulation in the model, the population was divided into susceptible individuals (S), exposed individuals (E), infectious individuals (I), asymptomatic individuals (A) and recovery individuals (R). The model diagram is shown in Fig. S1.
The differential equations of the model were used to describe the dynamic changes of each state. The corresponding model equations were as follows: www.nature.com/scientificreports/ 1. The model assumed that HFMD could not propagate vertically, and that all of the infectious individuals were infected through contact. Then we set birth rate (br), the natural mortality rate (dr), and the mortality rate of the infectious individuals (f). The mortality rate of all kinds of people in the disease spectrum was low, and the mortality rate of population attributable to HFMD was even lower, we set the mortality rate of the whole population as the sum of the mortality of the whole population and the mortality of HFMD. 2. Transmission of HFMD occurred via person-person, and the transmissibility between infectious individual and asymptomatic one may be different. So, the k was defined as the relative transmissibility rate of asymptomatic to symptomatic individuals. At the same time, we assumed the S would be potentially infectious as long as they are in contact with infectious individuals or asymptomatic individuals, and the coefficient of the infection rate was set as β. 3. Infectious individuals (I) and asymptomatic individuals (A) came from the susceptible individuals, so we considered that there was a certain proportion of exposed individuals pE (0⩽p⩽1) transformed into I after incubation, another part of exposed individuals (1-p) E were transformed into A after incubation as well. At a certain time (t), we set the speed of the development speed from the E to I pathway as ω (0⩽ω⩽1), and the development speed of E to A pathway as ω′. So the proportional coefficient of E to I was set as pω, and E to A was set as (1-p) ω′. 4. In our model, I and A may move to R, and the speed of recovering was in direct proportion to the number of individuals. The proportional coefficients were γ and γ′ respectively. 5. When I and A moved to R, we assumed that the infectious individuals recovered from the virus type they were diagnosed with, and that these recovered individuals had permanent immunity against this virus: thus, R was set as the end of the model. 1. The birth rate (br) and death rate (dr) were collected from 97 counties' statistical yearbooks in Jiangsu Province. 2. Studies showed that the proportion of dominant infection ranges were 19-47% 2,15,16 , selecting the median value 44.23%, therefore p = 0.4423. 3. The ranges of the incubation period (1/ω) were 3-7 days 2,17,18 , selecting the median value 5 days, therefore ω = 0.2. The latent period was set to 5 days, therefore ωʹ = 0.2. 4. The duration of symptomatic infection was 2 weeks 8,18 , therefore, the rate of disease removal γ = 0.0714.
The duration of asymptomatic infection ranged from 2 to 4 weeks 15,16 , Median of 3 weeks was chosen as the disease removal rate of asymptomatic patients, therefore, γ' = 0.0476. 5. The mortality of symptomatic infection ranged from 0.0001 to 0.0005 19,20 , selecting the median value 0.0003.
Parameter β is estimated by curve fitting. 6. There is no clear data or references to support the parameter κ, which is still uncertain. Therefore, in this study, we assumed κ = 1 for calculation, and sensitivity analysis was carried out to calculate its impact on the model.
The significance of each variable and parameter of the model is shown in Table 3.

Transmissibility evaluation index.
In this study, the population was not completely susceptible and artificially adopted some prevention and control measures, so we chose the effective reproduction number (R eff ) to calculate transmissibility. The calculation formula was as follows: