Interrelationship between daily COVID-19 cases and average temperature as well as relative humidity in Germany

COVID-19 pandemic continues to obstruct social lives and the world economy other than questioning the healthcare capacity of many countries. Weather components recently came to notice as the northern hemisphere was hit by escalated incidence in winter. This study investigated the association between COVID-19 cases and two components, average temperature and relative humidity, in the 16 states of Germany. Three main approaches were carried out in this study, namely temporal correlation, spatial auto-correlation, and clustering-integrated panel regression. It is claimed that the daily COVID-19 cases correlate negatively with the average temperature and positively with the average relative humidity. To extract the spatial auto-correlation, both global Moran’s \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathscr {I}}$$\end{document}I and global Geary’s \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\mathscr {C}}$$\end{document}C were used whereby no significant difference in the results was observed. It is evident that randomness overwhelms the spatial pattern in all the states for most of the observations, except in recent observations where either local clusters or dispersion occurred. This is further supported by Moran’s scatter plot, where states’ dynamics to and fro cold and hot spots are identified, rendering a traveling-related early warning system. A random-effects model was used in the sense of case-weather regression including incidence clustering. Our task is to perceive which ranges of the incidence that are well predicted by the existing weather components rather than seeing which ranges of the weather components predicting the incidence. The proposed clustering-integrated model associated with optimal barriers articulates the data well whereby weather components outperform lag incidence cases in the prediction. Practical implications based on marginal effects follow posterior to model diagnostics.


COVID-19 and weather situation in Germany.
According to the official 2018 census, the German states considerably vary in population, with North Rhine-Westphalia and Bremen having the highest and lowest population size of about 17,932,651 and 682,986, respectively, out of the total population size of 83,019,213. The states also have varied economic capacities in business, industries, tourism, and education, which affect their population size. For instance, the largely populated states like Bavaria and Baden-Württemberg have booming economy and offer plenty of employment opportunities due to the situation of renowned business centers and industries, whereas low-populated states e.g. Bremen are laid behind (see 64,65 ). Apparently, the number of cases and fatalities relatively depends on the population size. For instance, based on the report from Robert Koch Institute (RKI) on December 16, 2020, the largest populated state shared the highest 7-day incidence cases, and the smallest populated state shared the lowest. Given that the cases are population-driven, the dataset used for this study includes the daily confirmed COVID-19 cases for all the states from the official website of RKI 66 , which was later normal-ized per 100,000 inhabitants using the 2018 population census, see Fig. 1. This dataset spans the time window from March 01, 2020 to December 15, 2020. The normalization was intentional toward making the number of cases comparable across the states so as to allow for appropriate comparison with weather components that do not depend on the population (see similar treatments in 59,67,68 ). Here, the daily cases were defined as the difference of the confirmed cases since the earliest time of the report. As for the accompanying weather components, temperature and relative humidity data were retrieved from climate environment open data 69 . Time series of average temperature and relative humidity were obtained using the records of three weather stations Berlin-Marzahn (Berlin), München-Stadt (Bavaria) and Stuttgart-Schnarrenberg (Baden-Württemberg). This choice was motivated by data availability and the fact that the weather pattern throughout Germany is more or less the same, except in the alps where a negligible percentage of humans live. Average temperature ranges from − 0.766 to 27.13, and average relative humidity ranges from 39.38 to 93.53%. It seems the two weather components have a negative correlation showing equivalence between low temperature and high relative humidity or vice versa. Moreover, looking at the plot of cases by month in Fig. 1 in comparison with the weather components in Fig. 2, it can be seen that cases are generally higher in colder season and considerably reduce during the hot season.
In addition to the reported incidence, the spatial movement of the largest outbreak over the 16 states is also worth investigating. As depicted in Fig. 3, several stages in the timeline can be identified according to the dominance shown by different states. In the first three weeks in March, the largest incidence mainly altered between Hamburg and Baden-Württemberg. Bavaria and Saarland replaced them in the next three weeks. Bavaria hold a local election on March 15, and in the next day, a state of emergency was declared for 14 days with mobility restrictions 70 . Moreover, it is the first state to declare curfew that was imposed on March 20 71 . Saarland neighboring with badly affected French region Grand Est also incurred the same situation at midnight on the same day 72 . Lack of protective clothing and closure of medical practices were also reported from Bavaria 73 . Thus, Bavaria owed the largest incidence from time to time, even after the first few weeks. Outbreaks in initial reception facilities also contributed to the increase of cases in Bavaria. The largest incidence in May and in the first two weeks of June was dominated by Bremen. It was followed by Berlin and North Rhine-Westphalia until the end   Correlation studies. Referred studies in "Introduction" illustrate how meteorological factors correlated with the transmission of COVID-19. Highly transmissible disease like COVID-19 requires pathogens to remain active outside of the host body and relative humidity and temperature affect the virus's survival in the environment 44,77 . Another study engineering a SARS-CoV-2 isolate came across the fact that the virus can survive at least 28 days at ambient temperature 20 • C and 50% relative humidity on non-porous surfaces and is sensible to the variation of the weather components 78 . Therefore, it is considered noteworthy to examine the interrelationship between COVID-19 cases and meteorological factors. Many statistical methods have been used in earlier studies. According to the recent review in 61 , applicable methods other than descriptive analysis are Pearson correlation coefficient, linear, and non-linear regression, LOESS, two-way ANOVA, etc. Wavelet coherency analysis was used in 50 . This study used the Spearman-rank correlation so as to evaluate both the linear and monotonic relationship between two tested covariates. Additionally, auto-correlation between reported COVID-19 cases was also done by piling the spatiotemporal data into one time series, considering that normalized data vary in relatively small numbers. Lags up to 7 days from presently were selected. Therefore, every    82 . The spatial pattern is commonly seen to lie between three extreme cases: locally clustered, random, and locally dispersed. Locally clustered refers to the situation where neighboring states are similar in the level of daily new cases, under which spatial dependency rules out the spatial pattern. Locally dispersed refers to the inverse spatial dependency where neighboring states are dissimilar. Something in between is then referred to as random. Representation of these spatial patterns can be understood with the aid of a chessboard. If the spatial profile of daily cases in all states resembles the chessboard, then the spatial pattern is completely locally dispersed. If all the black cells would have gathered in one spot, then the spatial pattern is completely locally clustered. The random spatial pattern is then recognized from the way the black and white cells locate randomly on the board. This is extreme binary stratification that could never occur in the realism of epidemics, from which the corresponding global measure rarely reaches its bounds. Let us suppose that time is fixed and the daily cases from all states are reported as C = (c 1 , . . . , c S ) ⊤ with mean c . The other main ingredient in spatial auto-correlation is the spatial weight matrix W = (w ij ) , which measures the degree of contiguity among all the states. This study used the binary adjacency matrix, where w ij is 1 in case i and j share a common border or 0 in case otherwise (including diagonal entries). This definition is commonly used in the literature (referred to as "queen case") in contrast to distance-based proximity measure where central locations play a significant role as well as a definition of being a "center" is required to define the distances. Let us write Z = (z 1 , . . . , z S ) ⊤ := C −c and define |W| := i,j w ij . The global Moran's I and Geary's C statistic are given by respectively. According to the formulas, the global Moran's I represents the standardized spatial autocovariance by the variance of the data, while the global Geary's C replaces the autocovariance by the sum of the squared differences in all data values. Both formulas then differ in sensitivity controlled by the autocovariance. In terms of stability against uncertainty in the data, Wijaya et al. in 68 describe how Geary's C tends to vary less significantly than Moran's I when data are perturbed using noise of any kind. The current study presented Geary's C only for the sake of comparison. A measurement 0 < I → 1 (similarly 1 > C → 0 ) indicates the direction toward locally structured spatial pattern; I = 0 (or C = 1 ) random spatial pattern; and 0 > I → −1 (or 1 < C → 2 ) locally dispersed spatial pattern. Statistical inference is usually done under a total randomization assumption to have a decision outcome based on the values of the statistics 83 . The p-value is generated after normalization using the expected values E(I ) = −1/(S − 1) , E(C ) = 1 and variances V(I ) , V(C ) reported in the original studies 79,80 . The null hypothesis is that there is no spatial auto-correlation of the daily cases on the observed S states, meaning that I ≃ E(I ) and C ≃ E(C ) . Therefore, a p-value smaller than a predefined significance level α rejects the null hypothesis whereby either a locally structured or a locally dispersed spatial pattern occurs.
In contrast to the global measures, Moran's scatter plot measures the extent to which a state is considered a "hot spot" or "cold spot" or something in between 83 . It reports the coordinates (Z/σ C , WZ/σ C ) for all states, with σ C = Z ⊤ Z/S denoting the standard deviation of C. As a row-standardized weight matrix is utilized, i.e., |W| = S , the pooled estimator of the regressing linear line for these coordinates passing through the origin is given by (0, I ) . In the present context, a hot spot is defined as a state with a large number of daily cases surrounded by those with large numbers of cases (high-high). In the 2-dimensional Euclidean space, the coordinates of hot spots locate in the upper-right quadrant Q1. A cold spot, on the contrary, defines a state with a small number of cases surrounded by those with small numbers of cases (low-low). The coordinates of cold spots gather in the lower-left quadrant Q3. Other than these, local dispersion may occur falling into the following categories: a state with a small number of cases surrounded by those with large numbers (low-high) in the upper-left quadrant Q2, and a state with a large number of cases surrounded by those with small numbers (high-low) in the lower-right quadrant Q4. From the practical point of view, being a hot spot or cold spot may only rely on the health care capacity to ameliorate the disease burdens without imposing further restrictions to travel around neighboring states, except for those who travel across the border between scattered hot spots and cold spots. A state in a high-low or low-high spatial pattern, however, requires more restriction in traveling to neighboring states as the disease may diffuse (in case of high-low) or be absorbed (in case of low-high).
Simple case-weather relation. Let i and j denote the state and time index where i ∈ {1, . . . , S = 16} and j ∈ {1, . . . , N} . Our approach to modeling daily COVID-19 cases in all states in Germany was based on directly relating collected entities. These include presently (lag-0) reported cases C := (c ij ) , cases reported on the past seven days (lag-1, . . . , lag-7) from presently C −1 := (c i,j−1 ), . . . , C −7 := (c i,j−7 ) , average air temperature T := 1 S ⊗ (t j ) , and lag average relative humidity H := 1 S ⊗ (h j−25 ) corresponding to the cross-correlation result in Fig. 4. The notations 1 S and ⊗ denote the column vector of size S whose entries are 1 and the Kronecker product between two matrices, respectively. The final size of our observations is the entire time window length minus the maximal autoregressive lag, which is N := 290 − 7 = 283 (i.e. from March 8 until December 15, 2020). Let us denote β 0 as the intercept, β ind := (β 1 , . . . , β S−1 ) as the individual-specific effects (cut down by www.nature.com/scientificreports/ one term to avoid linear dependence with the intercept), β −i (for i = 1, . . . , 7 ) as the marginal effects of the lag incidence cases, β T as the marginal effect of the temperature, β H as the marginal effect of the relative humidity, and ε = (ε ij ) as the idiosyncratic error. The direct relationship among these covariates intends to not only skip additional transformations but also return direct marginal effects represented by the coefficients of the corresponding explanatory variables. This reads as which folds The indicator parameters σ (i) take binary values and will serve to drop certain variables in the model specification (by value 0), whenever necessary. This model represents, perhaps, the simplest panel regression model in the following sense. The marginal effects of the lag incidence cases and those of the weather components could have been raised to matrices like in vector autoregression with exogenous variables (VAR-X) models 84 . Besides appending too many parameters (entries of the endogeneous matrices), which may lead to overfitting, VAR-X models also require all the explanatory variables to be covariance stationary (see 85 for details), which is rarely the case for disease and weather data in the subtropics. As the only random spatial pattern was observed from the incidence data for almost all observations, no essential state-crossing marginal effects were expected. State-dependent marginal effects for the weather components were also not considered due to data aggregation and limitation, also to the intention to have unified marginal effects that work on the national level. Moreover, all lags smaller than the optimal values for the weather components were not considered for complexity reduction. For the reason of having straight-forward marginal effects, prior transformations were not applied to any of the variables. Despite its simplicity, the model (1) treats omitted variable bias by including individual-specific effects. These are the simplest terms assuming that the omitted variables only have constant effects on the daily  Figure 4. Spearman-rank correlation coefficients between daily cases from all states in Germany with the average temperature (above) and average humidity (below) on a moving window of 290 observations. Averaging throughout the states obtains the minimum of − 0.5223 (temperature) and maximum of 0.4194 (humidity) corresponding to the lags 0 and 25, respectively. www.nature.com/scientificreports/ COVID-19 cases in all the states. After all, the present study draws forth an outlook for compiling temperature and relative humidity data from all eligible stations as well as data of other confounding factors (e.g. other weather components, human mobility, employment opportunities, mapping of manufactures or public gatherings, etc) that not only add more explanatory variables but also clear up the heteroscedasticity issue.
Model including incidence clustering. Previous studies based their investigation on asking which ranges of weather components correctly predict incidence cases. This study asks a slightly different question: which ranges of incidence cases are correctly predicted by the existing values of the weather components. The values that fail to predict certain incidence cases due to insignificance would deem dropping. In 68 , this clustering strategy was designed to eliminate the weather dependency on the zero incidence cases, handling the zero-inflation problem appropriately. In the context of COVID-19, some extreme cases might have never been related to weather, for example superspreading events [10][11][12][13] and indoor aerosol transmission 8,9 . The basic aim of the clustering is then to correctly place the role of weather where it should have never predicted such events. The use of a transient function to replace this functionality was inapplicable to us, for which bias may arise from the functional choice and its related extension strategy for prediction. The clustering idea departs from stratifying the incidence data into M clusters (� k ) M k=1 separated by barriers θ := (θ k ) M−1 k=1 . In the closed forms, the clusters are given by Here, the incidence data were classified into three clusters ( M = 3 ) on the basis of practicality to call for lower, middle, and upper cluster. In principle, the specification is not bound to such a small number as fitting would be better with more explanatory variables. However, questions regarding complexity and practical interpretations might arise when using a large number of clusters. On the present choice, when for instance T (2) has to be dropped due to insignificance, this simply means that the average temperature fails to predict incidence cases in the range defined by the middle cluster 2 . This model then allows the lone cases to be "unexplained by temperature".
The fact that T k and H k change with the lower and upper barrier θ = (θ l , θ u ) , so does the pooled estimator β =β(θ) where β = (β 0 , β ind , β −1 , . . . , β −7 , β 1 T , . . . , β 3 H ) . Our aim is to find the optimal barriers such that the squared error between data C = (c ij ) and the model approximate C[β](θ) achieves its minimum. Mathematically, the preceding statement translates to the following problem The pooled estimator β follows from the straightforward formula in terms of matrix inverse and multiplication involving explanatory and response variable.

Results
Case-weather cross-correlation and case-specific auto-correlation. Figure 4 represents the correlation coefficients on a moving window of 290 observations with time lags from 0 to 30 days for each state. Notice that the reported daily COVID-19 cases correlated negatively with the average temperature and positively with the average relative humidity. The magnitude of the correlation coefficient with average temperature shows decreasing trends with lag for all the states. With no lag introduced, the correlations are negative and significant for all the states (p-values from 6.27 × 10 −34 to 1.17 × 10 −15 ). Averaging the correlation coefficients throughout the states, the minimum of − 0.5223 was obtained. This negative correlation is comparable up to certain ranges of minimum, maximum and average temperature to the studies in Brazil (with both average ranging from 20.9 to 27 • C and maximum temperature from 23.1 to 34.2 • C in 57 and with average temperature ranging from 16.8 to 27.4 • C in 86 ) as well as the data in 127 countries (with average temperature from − 17.8 to 42.9 • C in 87 ). In New York 88 , the correlation was positive and insignificant for average and minimum temperature but positive and insignificant for the maximum temperature. In Oslo, Norway 89 , the correlation was negative and insignificant for all maximum, minimum, and average temperature with 14 days time lag, but positive and significant correlation was obtained for normal temperature with 0, 5, 6, and 14 days lag. The temperature in Oslo ranged from − 7.5 to 21.9 • C during the study period. COVID-19 cases in Russian Federation exhibited positive significant correlation with minimum (− 17.78 • C to 8.89 • C), maximum (0.56 • C to 27.2 • C) and average temperature (− 2.78 • C to 16.1 • C) 46 .
As far as relative humidity is concerned, it can be observed from Fig. 2 that its average varies from 39.38 to 93.53%. The best lag was found 25 days with the correlation coefficient value of 0.4194 from averaging throughout  50,57 showed that the correlation was positive but not significant with minimum and maximum average humidity. Data from 127 countries 87 led to the conclusion that the relative humidity was correlated negatively and insignificantly with daily new cases. Table 1 shows the case-specific auto-correlations. Generally, Pearson is higher than Spearman-rank correlation coefficient. In addition, both Pearson and Spearman-rank correlation coefficient are significant with minimum 0.78 (p-values ≃ 0 ). From the column of lag-0, the auto-correlation generally swings from a large value at lag-1, then minima at either lag-3 or lag-4, to another large value at lag-7. The same behavior can be observed from the columns lag-1 until lag-3 where decrement rules out the first 4 lags and minima were found at either lag 3 or 4 days from the time series. This finding will set a basis for those in the panel regression models, as seen shortly.  Fig. 5. Although the spatial pattern of the daily cases in all the states changes around with time, it is evident that randomness overwhelms the pattern for most of the time. The progression of p-values (especially below α ) indicates that, generally, no significant difference between Moran's I and Geary's C was observed except on the duration from November until mid of December where Geary's C shows more locally clustered spatial pattern.
The Moran's scatter plot for all the states in Germany was determined for all observations, see Fig. 6. For the sake of serial presentation, indexing the coordinates based on the quadrants is more favorable than plotting them. Overall, the results suggest that all the states show randomness with time in to which spatial pattern they belong. If one solely focuses on the recent observations (November 1 to December 15, 2020), then the following states have the tendency to occupy the following quadrants: Baden-Württemberg, Bavaria, Hesse, Thuringia (Q1); Brandenburg, Rhineland-Palatinate, Saxony-Anhalt (Q2); Hamburg, Mecklenburg-Vorpommern, Lower Saxony, Schleswig-Holstein (Q3); Berlin, Bremen, North Rhine-Westphalia, Saxony (Q4).

Panel regression models.
Variable choices for model specification were investigated. The criteria are based on not only fit and complexity (information-type criterion) but also insignificance, negative marginal effects, and multicollinearity driven by certain variables. For the fit and complexity, a minimal value of Bayesian Information Criterion BIC = −2 log(L) + log(N) · k 91 was sought. The first term of this criterion expresses maximization over the likelihood function L generated from our model and the second term includes the observation size N as well as the number of parameters k. Unlike Akaike Information Criterion (AIC) 92 that would have replaced log(N) by 2, BIC penalizes the number of parameters much more, especially for large observation sizes. Our study aims to drop certain variables toward cutting down BIC and amending insignificance as well as multicollinearity. The standard t-test was used for the significance test. Checking for multicollinearity follows from computing the Inverse Variance Inflation Factor (1/VIF) values for all explanatory variables except the constant. A 1/VIF measures one minus the coefficient of determination derived from an OLS-regression whereby the variable under test serves as the response while the others as the explanatory variables. In this sense, 1/VIF of a value smaller than the rule of thumb 0.1 shows multicollinearity driven by the tested variable 93 . In addition, the p-value of the F-statistic is monitored, which measures if the overall variables are simultaneously significant; of which smaller than α = 0.05 indicates that they are. Not only can the model be designated to be better than just a constant, but multicollinearity can also be diagnosed. Johnston in 94 hinted the existence of multicollinearity as some p-values from t-tests are large while that from F-test is radically small, which agrees to the analytical inves- www.nature.com/scientificreports/ tigation in 95,96 . Besides these aspects, if certain marginal effects would be consistent with our auto-correlation study were also checked. From Table 1, it is seen how cases in the past 7 days positively predict present cases with the least auto-correlations found from cases from the past 3 and 4 days. This led to dropping negative marginal effects corresponding to lag incidence cases that may occur due to a certain model specification.
To deal with the model including incidence clustering (2), the computation of optimal lower and upper barrier (θ l , θ u ) as in (3) is necessary. The characteristic functions embedded in the objective function make the optimization problem non-smooth. The brute-force computations of the objective function in the upper-left triangle of the 50 × 50 grid in the domain [min i,j c ij , max i,j c ij ] 2 and a PSO algorithm 68 were put in comparison. From Fig. 7, PSO outperforms the brute-force computations in locating the optimal barriers that minimize the objective function, also in terms of computation time.
According to Table 2, the BIC value for the simple model (1) is relatively large, exacerbated by large degrees of freedom. The model including incidence clustering (2) gives the least BIC value due to a minimal likelihood function. Additionally, the insignificance of the entire individual-specific effects for both models was spotted. The rationale behind this can be connected to the fact that the entire profile of global and local spatial autocorrelation as well as the largest outbreak ("COVID-19 and weather situation in Germany" and "Spatial pattern") show randomness for almost all observations. Therefore, no state was worth constant recruitment (weighting) for its neighborhood to show a consistent spatial pattern throughout the observations.
Post-estimation diagnostics for all the models including those investigated during model specification were performed. Additional to the models including lag incidence cases and weather components, this study considered the models where either of these entities is present. The fitting results are presented in Table 3. For straightforward marginal effects and computation of optimal barries, the pooled estimator was considered subject to its inefficiency. The test was conducted via the comparison between fixed-effects and random-effects estimator and that between random-effects and pooled estimator. To the former, the two estimators were compared using Durbin-Wu-Hausman test 97,98 , where the fixed-effects estimator is assumed to be consistent, and the randomeffects estimator is efficient and assumed to follow a normal distribution. The null hypothesis suggests that the random-effects estimator is a consistent estimator regardless of the size of the data. According to Table 3, the p-value corresponding to the statistic greater than α = 0.05 indicates that the random-effects estimator is equally   99 , i.e., the model under the random-effects estimator returns zero variance in the state-dependent errors. Apparently, no panel effect was observed for all models except for those that include only weather components, in which case either random-effects or fixed-effects estimator is preferable. The inefficiency of the presented pooled, random-effects, and fixed-effects estimator is confirmed as serial correlation in all the state-dependent errors occurred. Wooldridge test 100 showed this. Therefore, a caveat remains for all models that their standard deviations of the coefficients are smaller and R 2 's are larger than they should be. After all, the pooled estimator is always consistent, even for a relatively small data size. As final practical remarks from the models, all the lag incidence cases give the waving effects in terms of lag where the cases 5 days and Table 3. Fitting results and diagnostics for the models (1) and (2) www.nature.com/scientificreports/ 7 days from presently predict the present cases the least and the most, respectively. Keeping the lag incidence cases, the weather components from model (1) give a consistent prediction with that from the cross-correlation study. Together with clustering, the marginal effects of weather were corrected for model (2). It was observed that temperature fails to predict cases in the upper cluster while relative humidity fails to cases in the middle cluster. Temperature seems to give a larger positive marginal effect for the middle cluster while relative humidity a negative smaller marginal effect for the lower cluster. As far as predictive performance is concerned, several findings can be highlighted. As the larger models exhibit no more issues with insignificance and multicollinearity, neither do the smaller models. For the model variant (1), the smaller models gain R 2 ≈ 0.8544 , BIC ≈ 23,972.15 (only lag incidence cases) and R 2 ≈ 0.3694 , BIC ≈ 30585.35 (only weather components), respectively. Meanwhile the model including only weather components shows the poorest performance; its BIC value is also radically larger than that of the model including only lag incidence cases. For the model (1), the impact of weather is rather small, as the decrease of temperature from a reference value e.g. T ≈ 20 • C to T ≈ 10 • C (i.e. by 50% ) is associated to the increase of COVID-19 cases for all states from e.g. C ≈ 20 by (|β T |10/20) · 100% ≈ 1.475% . When the lag incidence cases were dropped, the increase changes to (0.5681 · 10/20) · 100% ≈ 28% . Moreover, the increase of relative humidity from 60 to 80% (by 33%) is associated to the increase of the cases from C ≈ 20 by 2.46% (with lag incidence cases) and 22.56% (without lag incidence cases). The overall impression indicates the superiority of the model with only lag incidence cases when one designates fit to significantly matter than the number of parameters. For the model including incidence clustering (2), a different profile was obtained when only using non-dropped weather components: R 2 ≈ 0.7948 , BIC ≈ 25517.61 . Here, a significant improvement under incidence clustering becomes evident. Surprisingly, the model including the entire weather components even outperforms that including only lag incidence cases by fit and complexity: R 2 ≈ 0.8692 , BIC ≈ 23494.94 . All marginal effects corresponding to the temperature matrices are negative, and those corresponding to the relative humidity matrices are positive. It was observed that the temperature returns the smallest marginal effect on the COVID-19 cases in the middle cluster and relative humidity in the lower cluster. Besides the significance of the marginal effects, even no multicollinearity was observed. Apart from this, when the predictive ability is evaluated by R 2 and BIC amending multicollinearity and inconsistent predictors, it is still argued that combining lag incidence cases and weather components serve as the best models as presented in Table 3. The corresponding graphical fitting can be seen in Fig. 8.

Discussion
In this study, lags from the cross-correlation between average temperature and relative humidity were extracted to synthesize suitable variables in the regression models. Additionally, case-specific auto-correlation supports the model specification where lag-3 and lag-4 incidence would rather be insignificant predictors for the present incidence. Spatial auto-correlation using global Moran's I and global Geary's C was investigated in the framework of analyzing the spatial effect in COVID-19 transmission. The global measures indicate random spatial patterns most of the time, except there were either local clusters or dispersion in recent observations from November 1 to December 15, 2020. Moran's scatter plot was then used to disclose the local behavior of the spatial pattern. The result shows that the distribution of the hot spots and cold spots generally changed with time. The random spatial pattern justifies the model specification where the individual-or state-specific effects that would have served to endow specific states with constant weighting factors, were dropped.
In the simple random-effects model, the average temperature and lag relative humidity were shown to affect the incidence significantly, however, the resulting coefficient of determination is comparably much smaller than whenever only lag incidence cases were used; panel effect also raises in the former case. For the reason of placing the correct role of weather in predicting certain ranges of incidence, the weather components were grouped with the aid of a clustering strategy. The new clustering-integrated model accompanied by optimal barriers shows good agreement with the data whereby weather components outperform lag incidence cases in the prediction. On this matter, the fixed-effects estimator was the only presumably consistent estimator that also tackles the panel effect. For all models, it was observed that every explanatory variable competes against the others to be a significant predictor. Therefore, model choice together with its consequences (marginal effects), depend entirely on the decision-maker. Marginal effects can be guidance when a model is chosen a priori. When R 2 and BIC matter a lot, our recommendation is to opt for the clustering-integrated model with lag incidence cases and lag weather components. There it was found that temperature and relative humidity have negative, relatively small marginal effects on the cases in the lower cluster (below 13 cases per 100,000 inhabitants); the temperature has a large positive marginal effect on the cases in the middle cluster (between 13 and 36 cases per 100,000 inhabitants) and no marginal effect on the upper cluster (above 36 cases per 100,000 inhabitants); relative humidity has a large positive marginal effect on the upper cluster but none on the middle cluster. The clustering-integrated model with only weather components is recommended when weather receives more privilege than lag incidence cases. Our result is consistent with the cross-correlation study that temperature has negative marginal effects while relative humidity has positive marginal effects on the incidence in all clusters. The middle cluster receives the smallest marginal effect from temperature and the lower cluster from relative humidity. This hints physical consequences that temperature can only predict incidence cases during hot (summer) and cold season (winter), where cases clearly distinguish against each other from the data, not during transitional seasons (spring and fall). Relative humidity, on the other hand, is less likely to predict sinking cases during the hot season.

Conclusion
This study focused on the interrelationship between two weather components overlapping in many previous studies (average temperature and relative humidity) and COVID-19 incidence in Germany. Cross-correlation, case-specific auto-correlation, and spatial auto-correlation analysis were done to determine suitable variables www.nature.com/scientificreports/ and to explain the negligible panel effect in the panel random-effects models. In addition, the findings from the spatial auto-correlation provide the placement of the 16 states in the four quadrants from Moran's scatter plot and appropriate policy regarding traveling restrictions. The increasing demand for confounding factors to explain various incidence levels has been neutralized by the aid of incidence clustering. This strategy supports the idea of considering only certain hypothetical factors predicting COVID-19 incidence and general regression modeling wherein explanatory variables are limited. This localization of incidence that is correctly predicted by the two weather components has profound implications for public health authorities. The modeling does not only determine the extent of the prediction via marginal effects but also paves the way for precautionary actions amidst upcoming weather.