Does solar ultraviolet radiation play a role in COVID-19 infection and deaths? An environmental ecological study in Italy.

A significantly stronger impact in mortality and morbidity by COVID-19 has been observed in the northern Italian regions compared to the southern ones. The reasons of this geographical pattern might involve several concurrent factors. The main objective of this work is to investigate whether any correlations exist between the spatial distribution of COVID-19 cases and deaths in the different Italian regions and the amount of solar ultraviolet (UV) radiation at the Earth's surface. To this purpose, in this environmental ecological study a mixed-effect exponential regression was built to explain the incidence of COVID-19 based on the environmental conditions, and demographic and pathophysiologic factors. Observations and estimates of the cumulative solar UV exposure have been included to quantify the amount of radiation available e.g., for pre-vitamin D3 synthesis or SARS-CoV-2 inactivation by sunlight. The analysis shows a significant correlation (p-value <5 × 10-2) between the response variables (death percentage, incidence of infections and positive tests) and biologically effective solar UV radiation, residents in nursing homes per inhabitant (NHR), air temperature, death percentage due to the most frequent comorbidities. Among all factors, the amount of solar UV radiation is the variable contributing the most to the observed correlation, explaining up to 83.2% of the variance of the COVID-19 affected cases per population. While the statistical outcomes of the study do not directly entail a specific cause-effect relationship, our results are consistent with the hypothesis that solar UV radiation impacted on the development of the infection and on its complications, e.g. through the effect of vitamin D on the immune system or virus inactivation by sunlight. The analytical framework used in this study, based on commonly available data, can be easily replicated in other countries and geographical domains to identify possible correlations between exposure to solar UV radiation and the spread of the pandemic.

Introduction

COVID-19 infection has heavily affected the Italian population and, as of May 2020, 232,997 cases were recorded in the country, with very heterogeneous clinical manifestations. Despite unavoidable difficulties of collecting accurate data due to the different ways of recording clinical information by the regional health attendants, two features of the COVID-19 outbreak in Italy can be pointed out. First, the pandemic event proved to be particularly aggressive in older people (78.4% of the dead patients aged between 60 and 89 years), notably in nursing home residents (NHR). Second, the northern administrative regions experienced much more cases (up to nearly 1% of the population, in the Aosta Valley) and deaths (up to 1.6 per 1000 inhabitants, in Lombardy) than the southern ones (only 0.06% cases, and 0.05 deaths per 1000 inhabitants, in Calabria). The latter characteristics can be hardly explained by the activation of the lockdown regulations in the most affected regions, starting from 23rd February 2020, and by the subsequent lockdown in the entire country, since 9th March 2020. Indeed, recent analyses (La Rosa et al., 2020) detected SARS-CoV-2 RNA in the sewage of the main cities in northern Italy as early as in December 2019 and January 2020. Given that the infectious impact dose of viable SARS-CoV-2 required to cause infection in another person may be easily reached, for instance, with a 20-minute unprotected conversation (Burke et al., 2020), and considering that travels from north to south of the country were very frequent before the lockdown, it appears unlikely that the containment measures put into effect as of February and March were the only important discriminating factor.

Several studies, described in more depth in Section 2, have already attempted to explain the geographical and temporal patterns of the pandemic both in Italy and worldwide by one or more environmental factors, possibly driving the virulence and the severity of the disease. Among them, solar ultraviolet (UV) radiation at the Earth's surface is especially noteworthy, since it can impact on both spread and mortality of various infections in multiple ways (Section 2.1). These above-described issues have triggered our study, the first in Italy employing actual measurements of UV radiation, with the aim of better understanding such a high rate of infection in the north of the country, and, notably, to test whether the different ambient UV distribution over the Italian territory supports the hypothesis of an active role played by solar UV radiation in the virus transmission and outcomes.

Indeed, Italy, due to its wide latitude extension and its different environmental conditions, is the ideal test area to explore a correlation between the clinical outcomes of the pandemic and the geographical distribution of solar UV irradiance. To this purpose, we exploit the most recent satellite products from the Tropospheric Emission Monitoring Internet Service (TEMIS), validated for the first time globally on a country-wide scale using a novel dataset of ground-based accurate solar UV measurements from several Italian stations at latitudes ranging from 46.6°N in Trentino-South Tyrol, to 35.6 N° in the Lampedusa island (Section 3). The explanatory power of solar UV radiation is additionally compared to the relative importance of other environmental (Section 3.4.1), comorbidity and demographic factors (3.4.2, 3.4.3). The results of our analysis, reported in Section 4, and further discussed in Section 5, may provide a significant step towards supporting the role of solar UV radiation in the COVID-19 disease transmission and clinical outcomes. Furthermore, the statistical setup employed in this study can be applied to other series and geographical domains in an easy and effortless way, being it based on publicly available data, in order to provide, and to motivate deeper and more focused investigations.

Literature review

This section presents some of the main studies about the biological effects of UV radiation with reference to the present pandemic (Section 2.1), and some examples of ecological studies related to COVID-19 and including different environmental variables (Section 2.2). While an exhaustive review is out of the scope of this article and can be found elsewhere (e.g., Dobricic et al., 2020), here we offer a quick overview of the current scientific knowledge with the aim of framing the problem and motivating the present study.

Biological effects of UV radiation and possible connections to the present pandemic

Solar UV radiation may affect the spread of the SARS-CoV-2 virus and the outcomes of COVID-19 in patients through multiple pathways. First of all, high levels of solar radiation, especially in the UV-B region (formally defined in the wavelength range 280–315 nm, but mainly detectable at wavelengths larger than 290 nm at the Earth's surface, due to absorption by the atmosphere), might reduce the dissemination of the virus, notably SARS-CoV-2, by acting as a disinfectant for contaminated nonporous materials, as found by Sagripanti and Lytle (2020). Additionally, Ratnesar-Shumate et al. (2020) detected inactivation of SARS-CoV-2 on surfaces in culture media when exposed to simulated sunlight representative of the summer solstice at 40°N latitude at sea level on a clear day, and mainly attributed this process to UV-B radiation. Herman et al. (2020) used the data from both previous studies to determine the inactivation times from 290 to 315 nm UV-B solar radiation for SARS coronaviruses CoV and CoV-2 using global satellite data. The exposure duration needed to reduce the virus particles by 90% depends on the assumed inactivation effective doses and the type of virus, ranging from ca. 90 min for SARS CoV viruses for midlatitude sites between March and September, to even less for SARS CoV-2.

Second, the exposure to UV radiation favours the synthesis of vitamin D in the body. In fact, exposure to UV radiation determines the photo-conversion of the pro-vitamin D3 (7-dehydrocholesterol) in the skin to pre-vitamin D3 (Norval et al., 2010). Subsequently, it meets two distinct hydroxylations, the first in the liver, where it is converted to 25(OH)D, and the second in the kidney, with synthesis of its active form (1-25OH2D3 or calcitriol). Finally, it binds to specific receptors (VDRs), acting on different tissues with a similar hormonal mechanism.

Among other effects, such as the regulation of calcium and phosphate metabolism (Bischoff-Ferrari et al., 2019) and the possible prevention of several diseases (Maretzke et al., 2020), vitamin D plays a well-established role in the innate and adaptative immune defense of viral and bacterial infections (Calder et al., 2020). In this perspective, a recent review (Grant et al., 2020) has highlighted that vitamin D is able to keep the integrity of the tight junctions and the pulmonary barrier, by its antiviral properties and its possible role in mitigating pneumonia and hyper inflammation (Borella et al., 2014). These properties are thought to originate from the vitamin D ability to modulate gene expression by activating the VDRs, ubiquitously distributed in many target cells, including immune cells, and by promoting the expression of antimicrobic peptides such as cathelicidins and beta-defensins, which also have antiviral and immunomodulatory activity (Tripathi et al., 2015). Indeed, VDRs have been found on both resident immune cells and respiratory epithelial cells (Lange et al., 2009) and a metanalysis, considering interventional studies in more than 11,000 patients, reported that vitamin D supplements reduce by two-thirds the incidence of acute respiratory infections (not induced by COVID-19) in patients with low levels of blood serum 25(OH)D which is a good indicator of how much vitamin D is in the body (Martineau et al., 2019).

Some observational outcomes endorse the idea that vitamin D (and hence an environmental component such as solar UV radiation) may also play a role in the prevention and development of the SARS-CoV-2 pandemic event (Isaia and Medico, 2020).

Notably, hypovitaminosis D could have contributed to the individual susceptibility to the virus. Although direct evidences are not yet well established (Rhodes et al., 2020), this hypothesis is supported by circumstantial evidences that associate the outcomes of COVID-19 to vitamin D levels (Mitchell, 2020): specifically, a correlation between vitamin D status and COVID-19 incidence and death rates has been reported in twenty European countries. Notably in Spain, Italy and Switzerland severe low Vitamin D levels in the aging and vulnerable populations were found (Ilie et al., 2020). The terrestrial amounts of solar UV radiation mediated by vitamin D status could also explain differences in susceptibility of particular ethnic groups: this occurred, for example, in Connecticut (USA), where African-American population known at higher hypovitaminosis D risk (O'Connor et al., 2013), have a higher rate of infection and death by COVID-19, in comparison to lighter-skinned people at a comparable latitude (Laurencin and McClinton, 2020). Moreover, a retrospective, observational analysis of deidentified tests was performed in over 190,000 patients from all US 50 states to determine if circulating 25-hydroxyvitamin D [25(OH)D] levels were associated with severe acute respiratory disease coronavirus 2 (SARS-CoV-2) positivity rates: SARS-CoV-2 positivity was strongly and inversely associated with circulating 25(OH)D levels, and this relationship persisted across latitudes, races/ethnicities, both sexes, and age range (Kaufman et al., 2020). Finally, a very interesting study demonstrated that administration of a high dose of 25-hydroxyvitamin D, a main metabolite of vitamin D endocrine system, significantly reduced the need for Intensive Care Unit treatment of patients requiring hospitalization due to proven COVID-19 (Castillo et al., 2020).

Since more than 90% of the vitamin D in the majority of individuals is effectively produced by solar UV radiation (mainly in the UV-B band) and only a modest supplement can be provided from food (though less than 20%), there are large seasonal differences in its production mainly caused by the variation of solar elevation and different skin exposure to UV radiation (Webb et al., 2010). This is particularly true for Italy, where UV-B irradiance reveals large differences during the year, with winter months characterised by an UV irradiance about tenfold lower than summer months (Calgani et al., 2016); the same study reports, as a consequence, that the seasonal stratification of vitamin D concentrations revealed an evident trend with the minimum mean value recorded in April and a maximum mean value obtained in September. Therefore, blood sampling seasonality should be regarded as an important preanalytical factor in vitamin D assessment and the amount of total vitamin D synthesized during the summer should be high enough to maintain the levels >50 nmol/l throughout the remaining part of the year (Bonelli et al., 2016). Although Italy is located in southern Europe, a large part of the elderly is lacking vitamin D (Lips et al., 2019). In particular, previous research proved that more than 70% of the elderly women in Italy have low levels of vitamin D (Isaia et al., 2003). Unfortunately, no recent 25(OH)D analyses are available on the whole Italian territory, which compelled us to use the exposure to ambient solar UV radiation as an alternative proxy of the vitamin D status in this study.

Ecological studies on environmental factors and COVID-19

Several studies carried out throughout the world have attempted to establish a relationship between one or more environmental factors and the spread (or the consequences) of the present pandemic. The most frequently studied environmental variables air temperature and relative humidity. Iqbal et al. (2020) for example, explored the temporal correlation between daily average temperature and daily new COVID-19 cases in Wuhan, China, in the period between January and March 2020, and found a positive association between both datasets. On the other hand, based on data collected in the same period in the 10 most affected provinces in China, Shahzad et al. (2020a) showed that the relationship between temperature and COVID-19 actually differs for different provinces, the correlation between the temporal series being mostly positive for some areas and negative for others. Sajadi et al. (2020) found that the initial spread of the virus occurred on a narrow latitude belt, roughly located between 30°N and 50°N, characterised by specific conditions of temperature (5–11 °C) and absolute humidity (4–7 g/m3), consistently with the behaviour of a seasonal respiratory virus. In a recent study, Bherwani et al. (2020) included temperature and relative humidity in their epidemiological model for India, which also combined socio-behavioural aspects such as social distancing and lockdown regulations. They found fair correlation with temperature, compatible with disruption of the lipid layer of coronavirus at higher temperatures, but undermined by the most influencing factor of social distancing and confounded by other co-varying environmental factors. No clear role of relative humidity was established.

Air quality is another factor potentially impacting on the spread and impacts of the pandemic (Sciomer et al., 2020). Particulate matter (PM), for example, might act as a carrier of the virus, whose RNA was found in PM samples collected in northern Italy, the epicenter of the pandemic in Europe (Setti et al., 2020), or contribute, together with other atmospheric pollutants, to comorbidity and worsening the prognosis (Fattorini and Regoli, 2020), by triggering inflammatory response at molecular, cellular and organ levels (Mescoli et al., 2020). Fareed et al. (2020) correlated the air quality index and average relative humidity with COVID-19 mortality in Wuhan between January and March 2020, showing a negative relationship between air quality and deaths, and relative humidity and deaths. For the same period, Bianconi et al. (2020) investigated the spatial correlation between PM and demographic factors with the total cases and deaths for the different Italian regions and provinces using adjusted regression models, and found a significant association between exposure to PM, COVID-19 incidence and death. An even wider dataset, encompassing PM2.5, PM10, nitrogen dioxide, carbon monoxide, benzene, sulfur dioxide and ozone was used by Cazzolla Gatti et al. (2020) to train a machine learning algorithm and predict total cases and deaths in the Italian regions and provinces up to June 2020. Air pollution was found to be the main predictor for the effects of SARS-CoV-2, an increase of pollution by 5–10% being responsible of an increase by 19–28% of the cases and by 4–14% of the deaths. Using daily concentrations of carbon monoxide, nitrogen dioxide, sulfur dioxide, PM2.5, and PM10 between January and April 2020 in 5 states in USA and 5 provinces in China, Shakoor et al. (2020) found significant correlation between the reported cases and deaths with pollution. Similarly, Doğan et al. (2020) found positive correlation between pollution and COVID-19 daily cases in New Jersey, USA, in the period between March and July 2020, in addition to positive correlation of the cases with relative humidity and negative correlation with temperature. Shahzad et al. (2020b), correlated daily averages of temperature and PM2.5 with daily new cases of COVID-19 (up to July 2020) in the 4 most populated regions of Spain, showing that cases decrease at higher temperature and lower PM2.5. Generally speaking, although different research methods are used in the scientific literature and various confounding factors are not always taken into account, the major findings are consistent and highlight the contribution of pollutants such as PM2.5 and nitrogen dioxide in enhancing the COVID-19 spread and lethality (Copat et al., 2020).

Fewer studies have taken solar UV radiation into account as a possible influencing factor in the COVID-19 pandemic. Indeed, its important role has already been proved in the case of other diseases, such as the influenza virus in northern Europe between 2010 and 2018 (Ianevski et al., 2019), and the 1918–1919 pandemic influenza in the USA, for which UV radiation was empirically deducted from latitude (Grant and Giovannucci, 2009). For example, Gunthe et al. (2020), found a negative relationship between the number of COVID-19 cases and the UV-Index (a widely used parameter to provide the level of potentially harmful UV radiation to the public) by examining data from more than 80 cities in the world until the beginning of March 2020. Ahmadi et al. (2020) studied the spatial correlation among the infection rate and a number of socio-demographic and environmental factors in Iran in February and March, including solar radiation, and found that low values of wind speed, relative humidity, and solar radiation exposure are positively correlated with the infection rate. Sehra et al. (2020) explored the association between COVID-19 incidence and maximum daily temperature, precipitation, and solar UV radiation in the USA and showed that a higher UV index was associated with a lower incidence rate. Using a more recent dataset, Tang et al. (2020) found significant spatial and temporal negative correlation in the USA between COVID-19 (and other coronaviruses) cases and UV radiation (weekly and monthly radiation weighted by the erythemal, vitamin D and DNA damage action spectra and integrated over a given time period). Finally, Fazzini et al. (2020) examined the temporal relationships between several climate parameters, including solar radiation, and new daily positive swabs for COVID-19 for the Lombardy region in the period from March to April 2020. The multiple linear model resulting from the variables that explain at least 3% of the total variability showed that the spread of the pandemic is favoured by increased sunshine (in contrast to most studies) and high relative humidity, while the contribution of temperature is negligible. To the best of our knowledge, no further studies relating COVID-19 cases and solar radiation, notably in the UV band, in Italy have been conducted so far.

Materials and methods

Study design and chosen investigation period

Our multidisciplinary study involves data of solar UV ambient exposure (i.e., energy of the solar radiation reaching the Earth's surface over a specified time period per unit area) to quantify the amount of radiation necessary to activate pre-vitamin D3 synthesis. These are statistically analyzed, together with other environmental (Section 3.4.1) comorbidity and demographic factors (3.4.2, 3.4.3), to explain the geographical distribution of COVID-19 mortality and infections among the different Italian regions.

Fig. 1 shows the evolution of the pandemic in Italy, in terms of daily new cases and deaths (right vertical axis, logarithmic scale), and the yearly behaviour of the solar UV radiation reaching the surface (left axis) at the northernmost and southernmost measurement stations of the considered dataset (Section 3.4.1), i.e. Renon (46.6°N) and Lampedusa (35.6°N). The lockdown measures are additionally represented by vertical bars: initial shutdown of some activities (e.g., schools) in the northern regions (23rd February 2020, marker 1); closure of schools of every order on the national territory (4th March 2020, marker 2), of retail business and restaurants (9th–11th March 2020, marker 3), and of all other “unnecessary” activities (21th March 2020, marker 4). The “phase 2” started on 4th May 2020 (marker 5), with the progressive easing of the regulations, culminating in returning to the free circulation across the regional boundaries on 3rd June 2020 (marker 6). The period chosen in this study to investigate the evolution of COVID-19 cases and outcomes (response variables, Section 3.3) spans the interval between 25th February and 31st May 2020, i.e. it includes the whole “phase 1” and the beginning of “phase 2”, and encompasses 97% of the total cases in the “first wave” of the pandemic in Italy and 96% of the total deaths. In other words, these data describe the impact of the pandemic on a rather uniform population that was not prepared for an event of this size and strength. We refrained from considering the summer period, since the number of the cases is greatly reduced (thus increasing the statistical noise), and other variables might have altered the evolution of the pandemic. Among them, the mean age of the affected decreased (younger people were affected by COVID-19 during the summer period), better treatments in hospitals and better preparation of the people (frailty populations have been more protected by the social isolation and distancing and by the use of masks). Notably, the mortality observed from summer to autumn is probably altered by these provisions.

Fig. 1

Yearly course of the UV radiation at the surface (left axis) at the Renon (46.6°N, blue line) and Lampedusa (35.6°N, brown line) stations, and evolution of the pandemic in Italy (right vertical axis, logarithmic scale), in terms of daily new cases (continuous line) and deaths (dotted line). The lockdown measures are additionally represented by vertical dashed lines (see main text for further explanations).

Fig. 1 also shows that, although the day-to-day variability of solar UV radiation can be large due to the impact of clouds, the seasonal course (envelope) of the UV exposure is rather smooth. Moreover, since a long exposure interval (several months) is required in order for cholecalciferol to accumulate in the fat tissues of the human body (Section 2.1), the expected effects from vitamin D will be modest and difficult to interpret on the temporal scale of the current duration of the pandemic. Indeed, based on previous research, e.g. in Great Britain (Hyppönen and Power, 2007) and the United States (Kroll et al., 2015), the average 25(OH)D levels in the population are projected to increase by no more than 20–30% in the period February–May considered in the present study. Incidentally, this also justifies the positive correlation (contrasting to most studies) between new daily positive swabs and UV exposure found by Fazzini et al. (2020), in the too short period between March and April 2020. For this reason, no temporal correlation between the response variables and the environmental, demographic and pathophysiologic factors was attempted in this study, which is based instead on spatial correlations.

Statistical analysis

Different univariate models were built as linear regressions of the logarithm base 10 (Sajadi et al., 2020) of the response variables (listed in Section 3.3 and also reported in Table 1 ), as a function of each of several factors aggregated on a regional basis (independent variables, Section 3.4 and Table 2a, Table 2b, Table 2c ). For each independent variable, the parameters of the regression were fitted with the least-squares method and the statistics of the regression analysis were calculated. Multivariate analyses were also performed in order to determine how large is the increase of the explained variance when using multiple predictors, and the logarithm base 10 of every response variable was fitted by means of a multiple linear regression to several independent variables simultaneously. In this case, the independent variables were selected as the subset of variables that resulted as statistically significant in the univariate analysis (significance threshold value set to 0.05). The analysis was performed in the MATLAB R2019b environment.

Table 1

Covid-19 data for each Italian region obtained from January to May 2020. These results were used as response variables in the statistical analysis.

Italian regions (N)a	North Lat.a	Populationb	# of deathsc	# of casesc	# of swabsc	Deaths/pop (%)	Affected/pop (%)	Affected/swabs (%)
Aosta Valley (1)	45.74	125,666	143	1184	15,199	0.114	0.942	7.79
Piedmont (2)	45.04	4,356,406	3867	30,637	319,135	0.089	0.703	9.60
Liguria (3)	44.25	1,550,640	1465	9663	106,421	0.094	0.623	9.08
Lombardy (4)	45.28	10,060,574	16,112	88,968	753,966	0.160	0.884	11.8
Trentino-South Tyrol (5)	46.04	1,072,276	753	7027	154,780	0.070	0.655	4.54
Veneto (6)	45.26	4,905,854	1918	19,152	669,650	0.039	0.390	2.86
Friuli Venezia-Giulia (7)	46.04	1,215,220	333	3273	134,139	0.027	0.269	2.44
Emilia-Romagna (8)	44.30	4,459,477	4114	27,790	325,410	0.092	0.623	8.54
Tuscany (9)	43.46	3,729,641	1041	10,104	251,970	0.028	0.271	4.01
Umbria (10)	43.07	882,015	76	1431	70,493	0.009	0.162	2.03
Marche (11)	43.46	1,525,271	987	6730	103,698	0.065	0.441	6.49
Lazio (12)	41.54	5,879,082	735	7728	255,894	0.013	0.131	3.02
Abruzzo (13)	42.21	1,311,580	405	3222	75,634	0.031	0.246	4.26
Molise (14)	41.34	305,617	22	436	14,631	0.007	0.143	2.98
Campania (15)	40.21	5,801,692	412	4802	201,765	0.007	0.083	2.38
Apulia (16)	41.07	4,029,053	504	4494	118,575	0.013	0.112	3.79
Basilicata (17)	40.38	562,869	27	399	29,776	0.005	0.071	1.34
Calabria (18)	38.54	1,947,131	97	1158	70,182	0.005	0.059	1.65
Sicily (19)	38.07	4,999,891	274	3443	150,349	0.005	0.069	2.29
Sardinia (20)	39.13	1,639,591	130	1356	57,215	0.008	0.083	2.37

Numbers in round brackets are referred to Regions reported in Fig. 3.

Latitude of the capital city.

Italian National Institute of Statistics (ISTAT) – data recorded 31/12/2019.

http://www.salute.gov.it/portale/nuovocoronavirus/dettaglioContenutiNuovoCoronavirus.jsp?area=nuovoCoronavirus&id=5351&lingua=italiano&menu=vuoto.

Table 2a

Environmental factors for each Italian region used as independent variables in the statistical analysis.

Region	Vitamin D UV exposure (MJ/m2) a	PM10 (μg/m3) b	Relative humidity (%) c	Air temperature (°C)
	(Jun–Dec 2019)	(2015–2019)	(Jan–May 2020)	(Jan–May 2020)
Aosta Valley (1)	0.873	18.0	67.5	5.7
Piedmont (2)	0.899	30.6	71.2	9.6
Liguria (3)	0.925	21.8	66.9	11.9
Lombardy (4)	0.900	31.4	70.3	10.0
Trentino-South Tyrol (5)	0.849	18.8	69.1	6.4
Veneto (6)	0.888	32.3	69.7	10.2
Friuli Venezia Giulia (7)	0.870	22.2	67.4	10.1
Emilia-Romagna (8)	0.924	27.9	70.9	10.7
Tuscany (9)	0.957	22.1	69.8	11.4
Umbria (10)	0.993	23.5	71.4	10.5
Marche (11)	0.966	27.0	68.2	11.4
Lazio (12)	1.033	26.5	71.3	12.1
Abruzzo (13)	1.001	24.4	66.7	11.0
Molise (14)	1.036	18.9	69.8	10.1
Campania (15)	1.063	29.7	68.6	12.3
Apulia (16)	1.066	23.7	66.8	12.9
Basilicata (17)	1.073	18.4	67.0	10.4
Calabria (18)	1.119	22.9	62.8	13.0
Sicily (19)	1.173	24.3	64.7	13.8
Sardinia (20)	1.093	22.2	71.3	13.4

Numbers in round brackets are referred to Regions reported in Fig. 3.

Calculated based on the TEMIS estimates.

Downloaded from the website of the European Environment Agency.

Estimated at 2 m above the surface by the COSMO 2l model.

Table 2b

Comorbidity factors for each Italian region used as independent variables in the statistical analysis.

Region	Ischemic heart disease deaths/pop (%)a	Circulatory system diseases deaths/pop (%)a	Cerebrovascular diseases deaths/pop (%)a	Diabetes mellitus deaths/pop (%)b
	(2013)	(2013)	(2013)	(2009)
Aosta Valley (1)	0.048	0.129	0.029	0.051
Piedmont (2)	0.043	0.147	0.043	0.059
Liguria (3)	0.048	0.142	0.035	0.039
Lombardy (4)	0.045	0.130	0.033	0.051
Trentino-South Tyrol (5)	0.048	0.127	0.026	0.074
Veneto (6)	0.044	0.132	0.030	0.059
Friuli Venezia Giulia (7)	0.052	0.139	0.033	0.048
Emilia-Romagna (8)	0.046	0.134	0.031	0.044
Tuscany (9)	0.042	0.137	0.040	0.081
Umbria (10)	0.056	0.143	0.035	0.108
Marche (11)	0.051	0.137	0.035	0.061
Lazio (12)	0.057	0.152	0.034	0.051
Abruzzo (13)	0.060	0.164	0.038	0.042
Molise (14)	0.059	0.170	0.041	0.102
Campania (15)	0.069	0.197	0.051	0.059
Apulia (16)	0.047	0.145	0.031	0.045
Basilicata (17)	0.051	0.162	0.039	0.047
Calabria (18)	0.048	0.175	0.044	0.033
Sicily (19)	0.053	0.179	0.053	0.082
Sardinia (20)	0.041	0.129	0.033	0.059

Numbers in round brackets are referred to Regions reported in Fig. 3.

http://www.cuore.iss.it/indicatori/mortalita.

https://www.epicentro.iss.it/igea/diabete/Mortalita.

Table 2c

Demographic and social factors for each Italian region used as independent variables in the statistical analysis.

Region	Mean age (years)a	# of NHR/population (%)b	Average mortality rate (%)c
	(2019)	(2016)	(2018)
Aosta Valley (1)	45.6	1.07	1.17
Piedmont (2)	46.5	0.99	1.23
Liguria (3)	48.5	0.82	1.43
Lombardy (4)	44.7	1.11	0.99
Trentino-South Tyrol (5)	43.2	1.24	0.88
Veneto (6)	45.1	0.82	1.00
Friuli Venezia Giulia (7)	47.0	1.00	1.19
Emilia-Romagna (8)	45.7	0.77	1.12
Tuscany (9)	46.5	0.56	1.16
Umbria (10)	46.5	0.78	1.14
Marche (11)	46.1	0.59	1.12
Lazio (12)	44.6	0.38	0.97
Abruzzo (13)	45.7	0.54	1.12
Molise (14)	46.3	0.15	1.21
Campania (15)	42.2	0.43	0.92
Apulia (16)	44.2	0.34	0.96
Basilicata (17)	45.3	0.53	1.11
Calabria (18)	44.0	0.33	1.01
Sicily (19)	43.5	0.45	1.04
Sardinia (20)	46.3	0.47	0.99

Numbers in round brackets are referred to Regions reported in Fig. 3.

Italian National Institute of Statistics (ISTAT).

NHR = Nursing Home Residents https://www.cergas.unibocconi.eu/wps/wcm/connect/73024a58-4ae1-4ccf-9280-47b58ea2e448/Cap5OASI_2019.pdf?MOD=AJPERES&CVID=mWPGsIR.

https://www.statista.com/statistics/568083/death-rate-in-italy-by-region/.

Fig. 3

Cumulative ambient vitamin D UV exposures averaged for each Italian administrative region from 1st June to 31st December 2019. The calculations take the average altitude and the geographical distribution of the population into account. The numbers identify the respective Regions, as reported in Table 1, Table 2a, Table 2b, Table 2c.

Response variables

The study was conducted considering the available clinical outcomes of all 20 Italian regions. We employed these official data, aggregated for each region, in order to investigate the COVID-19 outcome variations throughout the national territory. In particular, we selected the following clinical data: (i) the percentage of positive pharyngeal nose swabs, for viral RNA detection, relative to the population of the corresponding region, (ii) the percentage testing positive for COVID-19 relative to the number of swab tests performed regionally, (iii) the COVID-19 death percentage relative to the population of the corresponding region. The data (reported in Table 1) were provided by the different health institutions and centrally registered at the National Institute for Health (Istituto Superiore di Sanità - ISS) according to the Mandatory Reporting System related to Ministerial Decree of 15th December 1990 - “Infectious and Transmissible Diseases Information System” (http://www.salute.gov.it).

Independent variables

Several parameters were considered in the regression as independent variables, described in the next paragraphs in detail.

Environmental data

Italy is a Mediterranean country which extends in a wide latitudinal zone, ranging from about 47 to 36°N. Its topography is very complex, mainly because of the presence of the Apennines extending from the north to the south, and the Alps at the north, the latter surrounding the wide Po basin. The large altitude gradients due to the presence of such high mountains and the latitudinal extent of more than 10° affect the meteorological conditions which vary significantly throughout the country, even within short horizontal distances of a few kilometers, resulting in a correspondingly large variability in the levels and the spectral distribution of the solar UV irradiance that reaches the surface (Meloni et al., 2000). Topography also impacts air quality, by limiting the exchange of air masses and creating permanent atmospheric pollution hotspots, such as the Po basin (actually one of the strongest in Europe, also due to the very populated urban settlements), or, conversely, by favouring winds contributing to the dispersion of pollutants, e.g. close to the coastline. The large variability of the environmental data among the different regions is clearly visible in Table 2a. These factors, described below in more detail, were aggregated for each Italian region as spatial averages, after weighting their original 2-D fields over the respective distribution of the population density (European Commission, Global Human Settlement, GHS-POP product release R2019A).

Solar UV irradiance at the Earth's surface

In the absence of recent clinical data assessing the distribution of 25(OH)D in Italy, we considered the ambient UV exposure spectrally weighted with the pre-vitamin D3 action curve, as defined by the International Commission on Illumination (CIE) (Bouillon et al., 2006), as a proxy to quantify the capability of ambient UV radiation to initiate the synthesis of pre-vitamin D3 in the human skin (the resulting quantity will be shortened hereafter to “vitamin D UV exposure”). Although the solar spectral irradiance is greater in the UV-A region (formally defined between 315 and 400 nm), the CIE action spectrum for pre-vitamin D3 peaks at much shorter wavelengths. This makes the 300–320 nm spectral range (mostly overlapping with the nominal UV-B region, defined in Section 2.1) the most biologically effective one. No radiation at wavelengths shorter than 290 nm reaches the Earth's surface due to absorption by the upper atmosphere. Both estimates of vitamin D UV exposure from the space, enabling uninterrupted monitoring over wide areas, and in-situ measurements at the ground, known for their high accuracy and reliability, were employed. Notably, the daily vitamin D UV dose operational product (v2.0) from TEMIS (Van Geffen et al., 2017) hosted by the Royal Netherlands Meteorological Institute (KNMI http://www.temis.nl/uvradiation/UVarchive.html) was used to calculate the cumulative exposure for each Italian region in the June–December 2019 period. As already mentioned, such a long interval was chosen since cholecalciferol accumulated in the fat tissues of the human body from solar UV exposure lasts for several months (MacLaughlin et al., 1982). The TEMIS product is based on radiative transfer models driven either by space-measured quantities, such as total ozone and cloud cover, or by climatologies, such as for ground albedo (Zempila et al., 2017). To account for the complex topography of Italy, notably in the Alps and Apennines, and avoid underestimation of the vitamin D UV exposure, we further resampled the TEMIS exposures on a finer grid (by rescaling from the original 0.25° × 0.25° to 0.025° × 0.025°) by employing a digital terrain model and an altitude amplification factor of 5%/1000 m (Zempila et al., 2017; Vitt et al., 2020). The resulting quantity was validated, for the first time globally on a country-wide scale, against ground-based UV measurements (spectroradiometers and broadband radiometers) from the nine Italian stations shown in Fig. 2 , selected on the basis of their high quality and traceability to the International System of Units (Diémoz et al., 2011; Fountoulakis et al., 2020). This comparison reveals large and statistically significant correlation (Pearson's coefficient = 0.88, p = 1.7 × 10−3) between TEMIS and the ground-based instruments and relatively low differences (within 10%, comparable to the expanded uncertainties of most of the used instruments and within the uncertainty of satellite products). Overall, the dependence of UV radiation on the latitude of the Italian stations, also found in previous research (Siani et al., 2013; Fountoulakis et al., 2020b), is well captured by TEMIS (Fig. 3 ).

Fig. 2

Cumulative ambient vitamin D UV exposure, in mega Joule per square metre (1 MJ = 106 J), for the period 1st June–31st December 2019 from ground-based measurements (horizontal axis) and TEMIS v2.0 (corrected for the surface altitude, vertical axis), for the nine selected Italian stations. These latter are listed in order of decreasing latitude, from the north to the south.

Air temperature and relative humidity

To account for possible influences of air temperature and relative humidity on the pandemic and its effects (Section 2.2), the time interval considered for these two variables is the same as the clinical data collection period, i.e. February–May 2020. For a better spatial coverage, the air temperature and relative humidity at 2 m above the surface were obtained for the whole Italian territory from the analysis of a numerical prediction model, COSMO 2I. Owing to its high resolution (about 2.2 km grid step) and its accurate physical scheme (non-hydrostatic, fully compressible), COSMO 2I allows a good representativeness of the meteorological fields even on complex terrains (Diémoz et al., 2019). A detailed description of COSMO can be found elsewhere (Baldauf et al., 2011).

Aerosols

Concentration measurements of airborne particulate matter with a diameter of 10 μm or less (named as PM10) were downloaded from the website of the European Environment Agency (https://www.eea.europa.eu/data-and-maps/dashboards/air-quality-statistics-expert-viewer). These consist in yearly regional averages from more than 500 stations in Italy, which were further averaged on the 2015–2019 period to account for the possible long-term effects of chronic exposure to PM10 on comorbidities (e.g., Fattorini and Regoli, 2020). The geographical pattern of PM10 mostly arises from the environmental, social and economic characteristics of the different Italian regions, the rate of industrialization, and the population density. Since these are permanent or slowly changing features, the resulting geographical distribution is expected not to change even during the shorter period of the pandemic in Italy and could thus account for other possible short-term effects on human health and virus transmission.

Comorbidity factors

For the second group of variables, three chronic pathologies (Table 2b) were arbitrarily chosen as highly prevalent causes of death in Italy in normal conditions, especially in elderly patients (data recorded in 2017): type 2 diabetes mellitus (causing 4% of deaths), cardio-vascular diseases (causing 10.4% of deaths), cerebrovascular diseases (causing 9.2% of deaths). The data for each Italian region, as reported by National Institute for Health (ISS) (http://old.iss.it/statistica/index.php?lang=1&tipo=28), were included into the statistical analysis.

Demographic and social factors

Finally, the third group includes the mean age, NHR per 100 inhabitants and percent mortality rate (demographic factors). These data (Table 2c) were provided by the National Institute for Statistic (ISTAT).

Results

The analysis of the regression on the percentage testing positive for COVID-19 relative to the population of the corresponding region (hereafter named as affected/pop) and on the COVID-19 death percentage relative to the population of the corresponding region (deaths/pop) yielded very similar results (Table 3a, Table 3b ). For both response variables, the univariate analyses of the estimated coefficients showed significant correlations (p-value <5 × 10−2) only with the following factors: the ambient vitamin D UV exposure, air temperature, NHR, the death percentage due to circulatory system diseases, to cerebrovascular diseases and to diabetes mellitus. An inverse relationship was found for all these factors, the only exception being NHR, which showed a direct relationship, i.e. COVID-19 deaths percentage increased with an increase in NHR. The observed relationships will be discussed in the next section.

Table 3a

Output of the univariate regression on the COVID-19 affected cases per population by region SE = standard error. Significant level at 0.05 (The variables with p-value <0,05 are marked in bold) CI = confidence interval.

Variable	Regression coefficient	SE	P-value	% of variation explained	Effect size (%)	95% CI
Vitamin D UV exposure	−3.97	0.42	<5.0 × 10−4	83.3	−8.6	−10.4	−6.8
PM10	0.02	0.02	3.9 × 10−1	4.1	1.1	−1.5	3.7
RH	0.06	0.04	1.3 × 10−1	12.4	10.2	−3.0	25.1
T	−0.14	0.03	7.0 × 10−4	48.1	−3.4	−5.1	−1.7
Mean age	0.10	0.06	1.1 × 10−1	13.5	11.2	−2.7	27.1
NHR/population %	108.72	18.79	<5.0 × 10−4	65.0	1.7	1.1	2.3
Mortality rate	0.11	0.07	1.4 × 10−1	11.8	2.8	−1.0	6.6
Ischemic heart diseases	−22.50	12.70	9.3 × 10−2	14.8	−2.6	−5.5	0.5
Circulatory system diseases	−14.25	3.57	8.0 × 10−4	47.0	−4.8	−7.2	−2.3
Cerebrovascular diseases	−33.45	11.38	8.8 × 10−3	32.4	−2.8	−4.7	−0.8
Diabetes mellitus	−15.73	3.05	<5.0 × 10−4	59.6	−2.1	−3.0	−1.3

RH = relative humidity.

NHR = nursing home residents.

T= air temperature.

Table 3b

Output of the univariate regression on the COVID-19 deaths percentage by region SE = standard error. Significant level at 0.05, (The variables with p-value <0,05 are marked in bold) CI = confidence interval.

Variable	Regression coefficient	SE	P-value	% of variation explained	Effect size (%)	95% CI
Vitamin D UV exposure	−4.77	0.66	<5.0 × 10−4	74.2	−10.3	−13.0	−7.4
PM10	0.03	0.03	2.4 × 10−1	7.5	1.8	−1.3	5.1
RH	0.06	0.05	2.2 × 10−1	8.2	10.5	−6.3	30.4
T	−0.15	0.05	6.8 × 10−3	34.2	−3.6	−6.0	−1.1
Mean age	0.12	0.08	1.6 × 10−1	10.8	12.9	−5.0	34.0
NHR/population %	132.86	25.56	<5.0 × 10−4	60.0	2.1	1.2	2.9
Mortality rate	0.12	0.09	2.1 × 10−1	8.7	3.0	−1.8	8.0
Ischemic heart diseases	−31.35	15.87	6.4 × 10−2	17.8	−3.6	−7.2	0.2
Circulatory system diseases	−17.91	4.59	1.0 × 10−3	45.9	−5.9	−9.0	−2.8
Cerebrovascular diseases	−41.22	14.69	1.2 × 10−2	30.4	−3.4	−5.9	−0.9
Diabetes mellitus	−18.49	4.28	<5.0 × 10−4	51.0	−2.5	−3.7	−1.3

RH = relative humidity.

NHR = nursing home residents.

T= air temperature.

The same tables show that mean age variations for each Italian region did not significantly correlate with affected/pop (p = 1.1 × 10−1) and deaths/pop (p = 1.6 × 10−1). In particular, the mean age among the Italian regions did not show significant differences, thus preventing this factor to be statistically significant when correlated to COVID-19 infections and deaths. In other words, the geographical variation of the mean age in Italy did not drive the geographical distribution of COVID-19 cases in a significant way. We also tested different metrics to describe the age distribution in each region, such as percentage of people older or younger than a specified threshold, and assess their correlation with the response variables, however the overall conclusion does not change, i.e. the population is about the same age throughout Italy.

The effect size of a specific factor is quantified as the percent change of the response variable resulting from a 1% increase of the mean of the factor itself, over the whole Italian territory. The vitamin D UV exposure had the larger effect size, both for affected/pop, −8.6% (95% CI: −10.4% to −6.8%) and death/pop, −10.3% (95% CI: −13.0% to −7.4%). The vitamin D UV exposure had also the largest percentage of variation explained, 83.3% and 74.2% for affected/pop and death/pop, respectively, meaning that the univariate regression based on this factor showed the best goodness of fit (Fig. 4a–b). Interestingly enough, the analysis of the regression coefficient proved that, among the regional demographic factors, only the NHR had a significant effect on the response variables, with p-value <5 × 10−4 for both response variables.

Fig. 4

Percent of affected cases (a) and mortality (b) per region and percent affected cases to the performed swabs per region (c), versus vitamin D UV exposure. Fitted values are derived from linear model of the logarithm of the variable versus vitamin D UV exposure.

To assess whether the inclusion of all significant factors led to a large increase of the explained variance, both for the affected/pop and death/pop response variables, multivariate analyses were performed on the subsets of factors presented above, considering only those with p-value <5 × 10−2. It is important to highlight that the multivariate regression is not intended as a predictive model, since the goal of this regression is just to compare the maximum % of variation explained using all the relevant independent variables to the one from the univariate models. The combined effect of the vitamin D UV exposure, air temperature, NHR, death percentage due to circulatory system diseases, to cerebrovascular diseases and to diabetes mellitus for each region on the response variables was analyzed. In terms of % of variation explained, it is interesting to notice that the combined effect of multiple factors did not substantially improve the goodness of fit of the regressions. In fact, the formulated multivariate regressions had a total explained variation of 78.5% for death/pop and of 85.5% for affected/pop, i.e. only around 5% more than the values obtained with the univariate regressions with the vitamin D UV exposure. Considering all the available factors in the multivariate analyses, i.e. including those with p-value >5 × 10−2 from the univariate analysis, thus likely overfitting the response variables, obviously gave the maximum explained variations achievable, equal to 91.1% and to 95.8% for death/pop and affected/pop, respectively. These values were roughly 15% higher than the variations explained by the vitamin D UV exposure in the univariate regressions.

The analysis of the regression coefficients on the percentage testing positive for COVID-19 relative to the number of swabs performed in the corresponding region (affected/swabs %) showed an overall worse index of goodness of fit (fifth column, % of variation explained, in Table 3c ), compared to the other response variables (fifth column, % of variation explained, in Table 3a, Table 3b).

Table 3c

Output of the univariate regression on the COVID-19 affected cases per swabs by region SE = standard error. Significant level at 0.05 (The variables with p-value <0,05 are marked in bold) CI = confidence interval.

Variable	Regression coefficient	SE	P-value	% of variation explained	Effect size (%)a	95% CI
Vitamin D UV exposure	−1.89	0.53	2.2 × 10−3	41.4	−4.2	−6.6	−1.7
PM10	0.02	0.01	1.9 × 10−1	9.5	1.1	−0.6	2.8
RH	0.03	0.03	2.1 × 10−1	8.7	5.6	−3.3	15.3
T	−0.05	0.03	8.9 × 10−2	15.2	−1.3	−2.8	0.2
Mean age	0.06	0.04	1.9 × 10−1	9.3	6.1	−3.2	16.4
NHR/population %	50.21	17.90	1.2 × 10−2	30.4	0.8	0.2	1.4
Mortality rate	0.07	0.05	1.5 × 10−1	11.4	1.8	−0.7	4.4
Ischemic heart diseases	−14.00	8.69	1.2 × 10−1	12.6	−1.6	−3.7	0.5
Circulatory system diseases	−7.17	2.85	2.1 × 10−2	26.1	−2.4	−4.4	−0.4
Cerebrovascular diseases	−14.26	8.73	1.2 × 10−1	12.9	−1.2	−2.7	0.3
Diabetes mellitus	−7.11	2.77	2.0 × 10−2	26.8	−1.0	−1.8	−0.2

RH = relative humidity.

NHR = nursing home residents.

T= air temperature.

The effect size is the estimated mortality change due to a 1% increase of the relative variable average.

The univariate analysis of the estimated coefficients only showed significant correlations (p-value <5 × 10−2) between the response variable and the vitamin D UV exposure (Fig. 4c), NHR, and the death percentage due to circulatory system diseases and to diabetes mellitus. The relationships between the response variable and the factors were the same as the other response variables presented above. The effect sizes of these factors on affected/swabs were roughly half the values of the effect sizes of the same factors on death/pop and on affected/pop. For this response variable as well, the vitamin D UV exposure had the largest percentage of variation explained, even though its value was as low as 41.4%. This analysis revealed that air temperature, death percentage due to vascular diseases and to diabetes mellitus did not have a significant correlation with affected/swabs, as opposed to affected/pop and death/pop. The multivariate regression had an explained variation of 43.8%, that was almost equal to the value obtained from the univariate regression with the vitamin D UV exposure as the only independent variable.

Discussion

As shown in the previous section, the most effective variables explaining the geographical distribution of infections and mortality from COVID-19 in Italy were vitamin D UV exposure, temperature, NHR, and some chronic diseases (cardio and cerebrovascular diseases, and diabetes mellitus). These are discussed hereafter.

Comparison with previous studies

First of all, the inverse relation between COVID-19 cases and vitamin D UV exposure found in this study is in line with most previous studies. For example, it was shown through regression modeling that, as of May 2020, countries in the Northern hemisphere were manifesting relatively high COVID-19 mortality, with an estimated 4.4% increase in mortality per degree of latitude, north of 28°N (p = 3.1 × 10−2), after adjustment for age of population (Rhodes et al., 2020). Our analysis revealed an even larger increase for Italy, of about 40% per degree of latitude (about 30% for affected/pop and 16% for affected/swabs). Hypovitaminosis D could have, for example, played a role in the northern regions during late winter, when the vitamin D accumulated in the fat tissues after solar UV exposure during summer and autumn is exhausted. Likewise, we could also speculate that increasing solar UV levels in late spring could have led to increased SARS-CoV-2 inactivation, and vitamin D production could have helped the public health measures (physical distancing and lockdown restrictions) in reducing the COVID-19 cases (e.g., the deaths were 919 on 27th March, 333 on 27th April, 117 on 27th May). Also, only few deaths were counted during the next months, beyond the period considered in the present study, despite the increase of cases (8 on 27th June, 5 on 27th July, 5 on 27th August). Our hypothesis is consistent with previous scientific literature (Section 2.1), and with some controlled clinical data reporting the association of vitamin D deficiency (<12 ng/ml) with higher risk of invasive mechanical ventilation and/or death (Radujkovic et al., 2020) and the presence of low 25(OH)D levels, usually used to check the vitamin D status, in COVID-19 patients versus other hospitalized patients (D'Avolio et al., 2020). Moreover, unpublished data collected in 33 unselected COVID-19-positive hospitalized patients from two hospitals of the Turin (northern Italy) area (4 women, 29 men, age between 41 and 87) show low 25(OH)D plasma levels in 100% of the patients, with an average level of 15.4 (±9.7) ng/ml. This value should be compared to the threshold of Vitamin D insufficiency (20 ng/ml) which, according to the Italian guidelines provided by many scientific societies (Nuti et al., 2019), is considered the minimum acceptable level in the general population (while in the elderly it should be increased to 30 ng/ml). Additionally, our data are consistent with the results of previous medical studies in Italy (Section 2.1), showing that part of the population is lacking or has low levels of vitamin D.

Besides, the inverse relationship between the reported COVID-19 cases and air temperature is consistent with the current understanding of the most favorable environmental conditions to the transmission of SARS-CoV-2 (Section 2.2). The direct relationship with the NHR per inhabitant was also expected and would point out the fact that promiscuous environments, such as nursing homes, are exposed to a higher risk than private homes.

Unexpected findings

The inverse relation between the number of deaths and common comorbidities was unexpected and deserves separate discussion. Considering the “past” comorbidity deaths (the latest data available referring to the year 2017) as a proxy for cardiovascular risk, the statistical analysis revealed that both COVID-19 death and incidence are lower in those regions where cardiovascular risk is higher. Therefore, our regionally-averaged data reveal that comorbidities alone are not able to justify the geographical distribution of COVID-19 incidence and death. In support of this important finding, the data recently published by the National Institute for Statistics (ISTAT) and by the Italian National Institute of Health (ISS) (https://www.istat.it/it/files//2020/07/ReportISS_Istat_Cause-di-morte-Covid.pdf) demonstrate that COVID-19 is directly responsible (from January to May 2020) for 89% of the deaths of people testing positive, while only a small fraction (11%) of the deaths is caused by cardiovascular diseases (4.6%), tumors (2.4%), respiratory illness (1%), type 2 diabetes (0.6%), dementia (0.6%), gastrointestinal diseases (0.5%) and other causes (1.3%). However, it should be also considered that a severe comorbidity rather than single illnesses alone could better explain a cause-effect relation with COVID-19, although current data do not allow to confirm or deny this statement. Future epidemiological investigations are necessary to address this question. Comorbidities must of course be taken into account as important factors when different age classes are considered, also given that they are frequently associated with hypovitaminosis D.

Finally, a short comment should be made about the lower variation explained by the regression describing affected/swabs. Although, in principle, this ratio should represent the best proxy of the infections and reduce the bias due to the different number of swabs performed in each region (“the more you search for, the more you will find”), the tests were performed following different regional strategies, in different ways and on different groups (e.g. old people, nursing home residents, doctors or nurses, relatives of infected etc.). Furthermore, one person could be tested several times during the course of the disease. These issues introduce large uncertainties in the infected/swabs ratio.

Sources of uncertainty

More generally, the results obtained by our study were affected by a number of uncertainties. First of all, it is well known that, in general, all observational retrospective ecological studies may suffer from ecological fallacy. Moreover, possible geographical differences in the initial number of infected people, and different regional impacts of the lockdown measures could not be excluded and might have contributed to the observed correlation. Also, the input data, e.g. the number of infections, might not be entirely accurate since, as already mentioned, the tests may not be uniformly distributed within the Italian territory. Likewise, the number of deaths officially reported might have been underestimated due to organizational reasons, although this bias could affect the whole territory in a similar way. Even the independent variables, such as the vitamin D UV exposure, are affected by uncertainties. For example, the ambient UV exposure differs from the UV exposure received by anatomical sites (i.e. personal exposure) which depends upon different postures under several environmental conditions and the duration of exposure (Schmalwieser, 2020). Furthermore, the exact shape of the effective pre-vitamin D3 action spectrum and the conditions under which it should be used are still under debate (Norval et al., 2010). Finally, the pre-vitamin D3 levels are also determined, for the same personal exposure, by the phototype and other personal characteristics, which might not be as uniform in Italy as in other countries (Yeum et al., 2016). However, despite the existing uncertainties, and in the absence of recent countrywide results from screening campaigns of the levels of 25(OH)D in Italy, the spatial variability of the cumulative solar UV exposure obtained from TEMIS is expected to be a sufficiently representative proxy of the spatial variability of the vitamin D adequacy in the population, as already proven by previous studies (Kelly et al., 2016; O'Sullivan et al., 2017; Prodam et al., 2016). Large scale analyses of the levels of 25(OH)D in Italians are however necessary in order to further clarify and definitely prove this assumption, and to discriminate the effect by vitamin D from the enhanced SARS-CoV-2 inactivation by UV radiation.

Conclusions

In this study, we tested the correlation between the presently COVID-19 available clinical data for every Italian region and several environmental, comorbidity and demographic factors, obtained from public sources.

Among all statistically significant factors, including comorbidities and nursing home residents, the vitamin D UV exposure was found to be by far the most effective in describing the geographical distribution of the pandemic in Italy, for all the response variables (infections per inhabitants, deaths per inhabitants and infections per swabs). In particular, the percentage of variation explained by vitamin D UV exposure in univariate regressions ranges from 41.4% for the percentage of infections per number of swabs to 83.3% for the percentage of infections per inhabitants. For multivariate regressions the percentage of variation explained by the most statistically relevant variables (solar UV radiation, air temperature, nursing home residents, comorbidities) ranges from 43.8% to 85.5%, respectively i.e. only slightly larger than the one from the univariate regression.

Although the analysis can be further improved using more refined statistical techniques, such as non-linear models or machine learning algorithms, the main outcomes of the study are not expected to change significantly, also given the relatively low number of response variables (equal to the number of the Italian regions, for each regression).

While the results do not directly imply a specific cause-effect relationship, they support the hypothesis that exposure to solar UV radiation – together with other factors such as social, demographic, and environmental conditions – affects SARS-CoV-2 transmission and its outcomes. The statistical setup and the available data did not allow us to discriminate the protective effect of vitamin D, synthesized after exposure to the sun, and the effect of inactivation of the virus by UV radiation in the outdoor environment. However, supplementation of vitamin D can be recommended to reduce susceptibility to infection, especially in the elderly and frailty population. In fact, among all possible pathways, hypovitaminosis D might be considered one factor, at least “guilty by association” able to contribute to COVID-19 associated diseases and to the pandemic spread. Some practical implications follow. Considering the large hypovitaminosis D occurrence in many countries, adequate exposure to the sun for sufficient vitamin D, balancing the negative and positive effects of the solar UV radiation, should be fostered, in compliance with physical distancing measures. This is particularly true for the elderly, since: (i) in general, they seldom set themselves to the sun and their skin has a lower efficiency in the synthesis of vitamin D compared to young people (MacLaughlin and Holick, 1985); (ii), they do not gladly eat food containing vitamin D, such as fat cheese, butter, mushroom, some fish; (iii) the assumption of food fortified with vitamin D is not widespread in southern European countries, as Italy. Therefore, public health campaigns should be promoted to increase consumption of vitamin D rich food or properly controlled pharmaceutical supplementation, and to encourage more adequate sunlight exposure, especially in countries exposed at the risk of vitamin D insufficiency or deficiency (Lips et al., 2019; Isaia et al., 2003).

Since the statistical framework employed in the present study relies on commonly available data, it can be easily applied to other contexts and geographical domains. Follow-up research could include, for example, an extension of the present analysis to the whole European continent.

Funding

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

CRediT authorship contribution statement

Giancarlo Isaia: Conceptualization, Supervision, Writing - review & editing. Henri Diémoz: Data curation, Methodology, Formal analysis, Writing - review & editing. Francesco Maluta: Formal analysis, Writing - review & editing. Ilias Fountoulakis: Formal analysis, Validation. Daniela Ceccon: Data curation. Alcide di Sarra: Data curation, Writing - review & editing. Stefania Facta: Data curation. Francesca Fedele: Data curation. Giuseppe Lorenzetto: Data curation. Anna Maria Siani: Data curation, Writing - review & editing. Gianluca Isaia: Supervision, Writing - original draft.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.