Selection of the study area
Shanghai is situated at 31°12′N north latitude and 121°30′E east longitude in the east part of the alluvial plain of the Yangtze Delta, adjacent to the Yangtze River estuarine and the East China Sea. It has four distinct seasons and abundant precipitation with a subtropical monsoon climate. The mean annual temperature in Shanghai is approximately 17 ℃, and the mean annual precipitation is more than 1100 mm, with 53% occurring between June and September. The study area is situated in the center of Shanghai, China, with a total area of 37.37 km2, measuring 6.15 km from east to west and 11.93 km from south to north (Fig. 1). In 2018, the study region included 14 sub-districts with a resident population of 1,057,700 .
The monthly total precipitation and monthly mean maximum and minimum temperatures were calculated based on data from the China Meteorological Administration . Weather variables were recorded at Xujiahui, which is located approximately 2 km from the study area.
The abundance of Ae. albopictus was analyzed by using the MOT (Tian®, Kaiqi Co. Ltd, Shanghai, China) and LT methods. The MOT trap  consists of a transparent cylindrical plastic jar (100 mm high, 70 mm diameter, 66 mm internal diameter) with a concave bottom (20 mm inward) and a black top cover with three conical openings of 100 mm diameter. When used as a collection container, a white filter paper of 70 mm diameter, used as an egg deposition substrate, was placed inside the bottom of the MOT, and 25 ml of dechlorination water was poured into the jar to keep the paper moist but not submerged. MOTs were placed outdoors on grasslands, kept away from direct sunlight, rain, and wind at ground level by a skilled technician, and maintained unchanged until the end of the study.
The MOTs were placed once a month in 2018 between April and November. To enhance the data in peak period, in 2019, the frequency was increased to once a week in week 20, week 23, week 25, weeks 27–39 and weeks 41–46. MOTs that were removed, emptied, or interfered for any reasons were excluded from further analysis. After four days, each MOT was collected to the laboratory. Species identification was performed with a stereomicroscope, and the MOT positivity index (MPI) was calculated as: MPI = number of the Aedes-positive MOTs / total number of MOTs retrieved × 100%. In 2018, the center of the study area was under construction and could not be placed with MOTs, so we chose the rest area and divided it into 128 grids with a side length of approximately 500 m as surveillance units. We found a residential area with vegetative coverage as the monitoring point in each unit and put eight MOTs around a center point (Fig. 1). The MPIs of the eight MOTs were used to represent the Ae. albopictus density at that point. In 2019, we added five surveillance units in the middle of the study area, which lead to a total of 133 surveillance units in 2019.
The LTs baited with carbon dioxide. On the third Wednesday of each week from May to November, the Centers for Disease Control and Prevention set two LTs in every sub-district throughout Shanghai as part of a city-wide mosquito surveillance program. LTs were usually collected from 4 pm to 10 pm. Contents of LTs were sent to the laboratory for species identification, and only female mosquitos from the traps were collected for data analyses.
The locations of MOTs and LTs were georeferenced using GPS CHCNAV X360H (WGS 1984 coordinate system) and later projected to the Shanghai local coordinate system. The georeferenced positions of MOTs and LTs in 2018 and 2019 are presented in Fig. 1.
The monitoring data were assigned into MOT2018, MOT2019, LT2018 and LT2019 groups and analyzed. The total collections of each LT and the mean MPI of each unit were calculated. Geostatistical analyses were conducted using ArcGIS 10.3 Spatial Statistics Tools (ESRI, China), and data in Microsoft Excel 2019 were imported into ArcGIS. The area of the minimum enclosing rectangle around the input features was used for cluster analysis based on Euclidean distance.
The average nearest neighbor (the mean distance from each feature to its nearest neighboring feature) analysis was then performed to calculate the nearest neighbor index. The index compared the mean distance to the null hypothesis states in which traps were randomly distributed, with index < 1 indicating a trend of clustering, while index > 1 indicating dispersion or competition.
We evaluated whether the mosquito abundances were spatially autocorrelated by calculating the incremental spatial autocorrelation (Global Moran’s I at multiple distances). The Global Moran’s I was tested using the permutation procedure based on feature locations and attribute values (total collections of each LT and mean MPI of each unit) against the null hypothesis (the absence of spatial autocorrelation). This analysis identifies the spatial patterns in the study area but did not indicate where such clusters occur, which could be determined by Local Moran’s I. The local Moran's I analysis was performed to assess the presence of hot spots with statistically significant clusters, cold spots, and spatial outliers.
The indicator Kriging  was used for geostatistical analysis. The Kriging and co-kriging techniques were performed using the ArcGIS 10.3 Geostatistical Analyst extension, including exploratory statistical analysis, variogram modeling, and producing a probability or standard error of indicator. The Kriging interpolation method was used to quantify the spatial structure of the data and predict species abundance at unsampled locations on the basis of the mean MPI of each unit. MPI was made into a binary (0 or 1) variable, where 1 means MPI value above the threshold, and 0 means MPI value below the threshold. The resulting interpolation map shows the probabilities of exceeding the threshold. The co-kriging uses information on several variables. We used several thresholds for MPI and then the binary data on each threshold (covariable) to predict the threshold of primary interest by co-kriging. The prediction surface with the lowest error output and standard error surface was applied and subtracted from the study area. A leave-one-out cross-validation method was used to determine whether the Kriging interpolation provided reliable estimates of the indicator at unsampled locations. The criteria used for accurate prediction in the cross-validation were requested to be the following: Root mean square standardized (RMSS) approximately 1, mean standardized (MS) approximately 0, and root mean square (RMS) approximately the average standard error (ASE).