Forest Area, CO 2 Emission, and COVID-19 Case-Fatality Rate: A Worldwide Ecological Study Using Spatial Regression Analysis

: Spatial analysis is essential to understand the spreading of the COVID-19 pandemic. Due to numerous factors of multi-disciplines involved, the current pandemic is yet fully known. Hence, the current study aimed to expand the knowledge on the pandemic by exploring the roles of forests and CO 2 emission in the COVID-19 case-fatality rate (CFR) at the global level. Data were captured on the forest coverage rate and CO 2 emission per capita from 237 countries. Meanwhile, extra demographic and socioeconomic variables were also included to adjust for potential confounding. Associations between the forest coverage rate and CO 2 emission per capita and the COVID-19 CFR were assessed using spatial regression analysis, and the results were further stratiﬁed by country income levels. Although no distinct association between the COVID-19 CFR and forest coverage rate or CO 2 emission per capita was found worldwide, we found that a 10% increase in forest coverage rates was associated with a 2.37‰ (95%CI: 3.12, 1.62) decrease in COVID-19 CFRs in low-income countries; and a 10% increase in CO 2 emission per capita was associated with a 0.94‰ (95%CI: 1.46, 0.42) decrease in COVID-19 CFRs in low-middle-income countries. Since a strong correlation was observed between the CO 2 emission per capita and GDP per capita (r = 0.89), we replaced CO 2 emission with GDP and obtained similar results. Our ﬁndings suggest a higher forest coverage may be a protective factor in low-income countries, which may be related to their low urbanization levels and high forest accessibilities. On the other hand, CO 2 can be a surrogate of GDP, which may be a critical factor likely to decrease the COVID-19 CFR in lower-middle-income countries.


Introduction
At the end of 2019, a newly emerging virus, later known as SARS-CoV-2 or COVID-19, was spotted in China, which has unexpectedly started a global pandemic [1,2]. The novel virus is a sister clade to the previous SARS-CoV [3], less fatal [4] but more contagious [5]. As of February 2022, nearly 400 million confirmed cases and more than five million deaths were recorded globally [6], indicating the great impact of the pandemic on global public health.
In the past two years, numerous studies have been conducted on COVID-19, and those with spatial analysis may play an essential role in understanding the pandemic [7][8][9] because spatial spreading is a primary characteristic of all infectious diseases. Spatial methods have contributed to uncovering the relationships between the COVID-19 casefatality rate (CFR) and sociodemographic or environmental factors, such as unemployment rate [10,11], commuting flows [12], and environmental pollution (e.g., NO 2 and particulate matter) [13,14]. Although many factors have been identified, the COVID-19 pandemic is still full of unknowns because various factors were involved, and their joint influence may The data on the forest coverage rate by country were obtained from the World Bank [32]. The data providers include the Food Agriculture Organization, electronic files, and websites. Moreover, the data on the forest area per capita were also obtained from the Worldometer websites [33]. The data providers include the Food and Agriculture Organization of the United Nations, the United Nations Population Division, the United Nations Statistics Division, and the World Bank. The data on the CO 2 emission by country were obtained from the Emission Database for Global Atmospheric Research, European Commission [34]. Following the latest update released in October 2021, CO 2 emission data are now available for each country for the time period between 1970 and 2020.
The COVID-19 data included in this study were obtained from the Our World in Data website [35], a collection of the COVID-19 data maintained by the organization Our World in Data and updated daily. The dataset includes country-level daily data on confirmed cases, deaths, testing, and other variables of potential interest. There are 67 variables from 237 countries and territories in the dataset by 6 December 2021.
According to the ISO 3166-1 alpha-2 code for spatial modeling, the datasets were linked together using the global geospatial vector database. Variables with missing values in more than 70% of countries were excluded.
In total, 27 variables were included in the final analysis of the presented study. They are: country, country's geographic variables/information (shape, location, longitude, and latitude of the centroid), forest coverage rate (in percentage of land area), forest area per capita, CO 2 emission per capita, total reported COVID-19 cases, total confirmed COVID-19 deaths, population, population density, gross domestic product (GDP) per capita, cardiovascular diseases (CVD) death rate, diabetes prevalence, stringency index of government's response to the COVID-19 pandemic (a composite measure based on 9 response indicators, including school closures, workplace closures, and travel bans, rescaled to a value from 0 to 100 with 100 = strictest response) [36], total tests for COVID-19, total vaccinations for COVID-19, life expectancy, median age of the population, percentage of population aged 65 years or older, percentage of population living in extreme poverty, percentage of female smokers, percentage of male smokers, hospital beds, percentage of population with basic handwashing facilities, and the human development index.
The study used data from publicly available worldwide databases, so it was not required and not possible to involve patients or the public in the design, conduct, reporting, or dissemination plans of our research. Hence, ethical approval is not applicable.
The study was conducted under the Guidelines for Accurate and Transparent Health Estimates Reporting (GATHER) statement [37].

Case-Fatality Rate (CFR)
The CFR of COVID-19 was defined as the number of total deaths due to COVID-19 divided by the number of confirmed COVID-19 cases by 6 December 2020 and multiplied by 1000. The CFR was selected for our study, as it may reflect disease severity and the efficiency of treatment and healthcare response, and strain [38]. Furthermore, the CFR is not constant, which can vary between populations and timeframes, and may depend on the interplay between the causative agent of the disease, the host, the environment, available treatments, and the quality of patient care. For instance, it can increase if the healthcare system collapses due to a sudden increase in cases [39]. The CFR in the current study is a crude rate without adjusting for age and sex due to no more detailed age/sex data available worldwide.

Statistical Analyses
Descriptive statistics (e.g., mean, standard deviation, median, minimum, maximum, and interquartile range (IQR)) were calculated for the variables.
Linear correlations between the variables were examined with Pearson's correlation coefficient [40]. Spatial autocorrelation (or dependence) is the relationship between spatial proximity in observational units and similarity among their values. Positive spatial autocorrelation refers to situations in which the nearer the observational units mean more similar values (and vice versa for its negative counterpart) [41]. This feature violates the assumption of independence among observations to be included in regression analyses. Spatial autocorrelation among the COVID-19 CFRs of the investigated countries was tested using a multivariate linear regression model and the Moran's I test [42]. The autocorrelation was visualized with the aid of the Matérn correlation coefficient [43].
Matérn correlation models were employed to check the relationship between the COVID-19 CFR and the ecological and socioeconomic variables. The latitude and longitude of the centroid of the countries were included as random effects in the Matérn correlation models [44].
Spearman's correlation was employed to probe for general patterns of associations within the core variables (forest, CO 2 emission, and COVID-19 CFR). Skewed variables were log-transformed before linear correlation or regression analysis. The multiple imputation method was employed to handle the missing values in the data (the missing values were assumed to be missing at random). A total of ten copies of the data were created, and each of the complete datasets was analyzed independently. Estimates of parameters of interest were then averaged across the ten copies to create a single estimate according to Rubin's rule [45].
Given the associations might differ by economic level, we also stratified the analyses by income levels of the countries defined by the World Bank [46].
The associations of the investigated ecological and socioeconomic variables with the COVID-19 CFR were presented as changes in the COVID-19 CFR per 10% increase in the ecological and socioeconomic variables. The corresponding 95% confidence intervals (CI) were also provided. An association with a p-value less than 0.05 was considered statistically significant.
All the analyses were performed in R V.4.1.2 (the R Foundation for Statistical Computing, Vienna, Austria) using the package spaMM [47], and in Python V.3.7 (Python Software Foundation) using the packages geopandas and geoplot [48].

Descriptive Characteristics of the Variables
In total, 264,861,480 reported COVID-19 cases and 5,255,849 confirmed deaths related to COVID-19 between 1 January 2020 and 6 December 2021, from 194 countries were included in the study. The descriptive statistics of the variables are shown in Table 1. Regarding the outcome and main factors of interests, we found a strong correlation between the COVID-19 CFR and CO 2 emission per capita (spearman ρ = −0.37, p < 0.001). However, we did not observe a distinct linear correlation between the COVID-19 CFR and forest coverage rate (spearman ρ = 0.01, p = 0.859) or forest area per capita (spearman ρ = 0.07, p = 0.314) ( Figure 1).
In the linear correlation analysis for log-transformed variables (Supplemental Figure  S1), high linear correlations were found between the forest coverage rate and forest area per capita (r = 0.70, p < 0.001); between CO 2 emission per capita and GDP per capita (r = 0.89, p < 0.001), the median age of the population (r = 0.76, p < 0.001), percentage of the population living in extreme poverty (r = −0.83, p < 0.001), percentage of the population with basic handwriting facilities (r = 0.80, p > 0.001), and the human development index (r = 0.87, p < 0.001); and between life expectancy and the percentage of the population aged 65 years or older (r = 0.85, p < 0.001), therefore, forest area per capita, GDP per capita, median age of the population, percentage of population aged 65 years or older, percentage of the population living in extreme poverty, percentage of the population with basic handwriting facilities, and the human development index were excluded in later regression analysis. Regarding the outcome and main factors of interests, we found a strong correlation between the COVID-19 CFR and CO2 emission per capita (spearman ρ = −0.37, p < 0.001).
However, we did not observe a distinct linear correlation between the COVID-19 CFR and forest coverage rate (spearman ρ = 0.01, p = 0.859) or forest area per capita (spearman ρ = 0.07, p = 0.314) ( Figure 1). In the linear correlation analysis for log-transformed variables (Supplemental Figure  S1), high linear correlations were found between the forest coverage rate and forest area per capita (r = 0.70, p < 0.001); between CO2 emission per capita and GDP per capita (r = 0.89, p < 0.001), the median age of the population (r = 0.76, p < 0.001), percentage of the population living in extreme poverty (r = −0.83, p < 0.001), percentage of the population with basic handwriting facilities (r = 0.80, p > 0.001), and the human development index (r = 0.87, p < 0.001); and between life expectancy and the percentage of the population aged 65 years or older (r = 0.85, p < 0.001), therefore, forest area per capita, GDP per capita, median age of the population, percentage of population aged 65 years or older, percentage of the population living in extreme poverty, percentage of the population with basic handwriting facilities, and the human development index were excluded in later regression analysis.

Worldwide Distributions of Forest, CO2 Emission, and COVID-19 CFR
Distributions of forest coverage rate and area per capita, CO2 emission per capita, and COVID-19 CFRs of the countries with data available are shown in Figure 2. The highest values for the forest coverage rate and area per capita were found in North America, northern Europe, northern South America, and southern Africa; for CO2 emissions per capita, in North America, Western Europe, the Middle East, and Australia; for the COVID-19 CFR, in southern North America, western South America, northern Africa, and eastern Asia.

Worldwide Distributions of Forest, CO 2 Emission, and COVID-19 CFR
Distributions of forest coverage rate and area per capita, CO 2 emission per capita, and COVID-19 CFRs of the countries with data available are shown in Figure 2. The highest values for the forest coverage rate and area per capita were found in North America, northern Europe, northern South America, and southern Africa; for CO 2 emissions per capita, in North America, Western Europe, the Middle East, and Australia; for the COVID-19 CFR, in southern North America, western South America, northern Africa, and eastern Asia.

Spatial Autocorrelation of the COVID-19 CFR
The residuals of the COVID-19 CFR from the common (non-spatial) multivariate linear regression model show apparent spatial dependence among the countries/territories ( Figure 3). The p-value from Moran's I test for the spatial autocorrelation of the residuals is 0.021.

Spatial Autocorrelation of the COVID-19 CFR
The residuals of the COVID-19 CFR from the common (non-spatial) multivariate linear regression model show apparent spatial dependence among the countries/territories (Figure 3). The p-value from Moran's I test for the spatial autocorrelation of the residuals is 0.021.

Spatial Autocorrelation of the COVID-19 CFR
The residuals of the COVID-19 CFR from the common (non-spatial) multivariate linear regression model show apparent spatial dependence among the countries/territories (Figure 3). The p-value from Moran's I test for the spatial autocorrelation of the residuals is 0.021. The estimated spatial autocorrelation coefficient of the COVID-19 CFR between two locations against their distance is shown in Supplemental Figure S2, with a strength parameter ν = 10.33 and a decay parameter ρ = 0.89. Basically, the autocorrelation coefficient  Figure S2).

Association of Ecological and Socioeconomic Variables with COVID-19 CFR
After controlling for the spatial dependence, the multivariable-adjusted associations of the investigated ecological variables (i.e., forest coverage rate and CO 2 emission per capita) and socioeconomic variables are presented as changes in the COVID-19 CFR per 10% increase in the variables in Table 2.
Worldwide, no statistically significant associations of forest coverage rate or CO 2 emission per capita with COVID-19 CFRs were found. However, the stratified results by economic level of the countries indicate that per 10% increase in forest coverage rate was associated with a 2.37‰ (95%CI: 1.62-3.12‰) decrease in the COVID-19 CFR in low-income countries, and per 10% increase in CO 2 emission per capita was associated with a 0.94‰ (95%CI: 0.42-1.46‰) decrease in the COVID-19 CFR in lower-middle-income countries ( Figure 4).

Discussion
This study examined the association between COVID-19 CFRs and forest or CO2. Although no general patterns were observed at the global level, some notable results were observed in specific countries. Our findings imply the indicative functions of forests and CO2 emission in the global pandemic. We could not provide solid mechanisms for explaining the observed phenomenon due to our study design. Nevertheless, we have tried to propose some potential clues or reasons to help understand the findings, which are presented in the following sections.

Forest Coverage and COVID-19 CFR
The current study included forest coverage rate and forest area per capita to indicate forest resources. The global analysis did not reveal a statistically significant association between the COVID-19 CFR and any of the forest variables. This result is inconsistent with a previous study where a dose-response association was spotted between green space and a reduced risk of COVID-19 mortality among counties in the USA from January to July 2020 [23]. However, another study reported that greenness was not associated with

Discussion
This study examined the association between COVID-19 CFRs and forest or CO 2 . Although no general patterns were observed at the global level, some notable results were observed in specific countries. Our findings imply the indicative functions of forests and CO 2 emission in the global pandemic. We could not provide solid mechanisms for explaining the observed phenomenon due to our study design. Nevertheless, we have tried to propose some potential clues or reasons to help understand the findings, which are presented in the following sections.

Forest Coverage and COVID-19 CFR
The current study included forest coverage rate and forest area per capita to indicate forest resources. The global analysis did not reveal a statistically significant association between the COVID-19 CFR and any of the forest variables. This result is inconsistent with a previous study where a dose-response association was spotted between green space and a reduced risk of COVID-19 mortality among counties in the USA from January to July 2020 [23]. However, another study reported that greenness was not associated with COVID-19 mortality in all counties in the USA [49]. These differences may come from different measures for the "natural environment" (leaf area index for the previous study), statistical methods, and spatiotemporal frameworks employed for analysis. The null effect we observed indicates that the role of and measures for "exposure to nature" should be considered cautiously. Currently, global urbanization has gathered most of the population into cities [50]. However, most forest areas are outside cities, especially in highly urbanized countries. Since distance strongly affects green space utilization [51,52], urban green spaces within or around cities may be the most popular and accessible type of natural environment for urban residents during the pandemic [53]. Thus, the forest variables calculated based on national forest area may not necessarily indicate the actual utilization or exposure to natural environments among residents in urbanized countries.
In the analysis stratified by income, we found that both the forest coverage rate and forest area per capita (Supplemental Table S1) were associated with decreased COVID-19 CFRs in low-income countries. Previous studies have speculated on some mechanisms to understand the decreased COVID-19 CFR, such as pollution reduction, physical activity promotion, and immunoregulation induced by microorganisms in natural environments [22,23,54]. However, these mechanisms did not consider low-income coun-tries, and they cannot explain the discrepancies among countries with different income levels. As discussed above, the accessibility of forests may be a problem in urbanized countries. In contrast, those low-income countries are less urbanized due to being highly reliant on agriculture [55], and more residents are distributed in rural areas, such as living in small villages and small towns. Therefore, the residents in low-income countries may live closer to natural environments than residents in urbanized countries. According to Pichlerová, M., D. Önkal, A. Bartlett, J. Výbošt'ok, and V. Pichler [17], forest coverage rate and settlement size are both factors of forest accessibility. A higher forest coverage and a smaller settlement size both indicate more forest visits. Therefore, for residents scattered in rural areas, a higher forest coverage rate may imply a greater chance of forest visits. In addition to physically entering forests, forests around neighborhoods may also benefit them remotely. For example, a recent study in Italy found that the percentage of COVID-19 deaths was lower in greener areas when compared to the regions with low forest coverage [56]. The phenomenon, explained by the author, might partially result from nondeciduous Mediterranean plants, which emit volatile immunomodulatory and antiviral compounds. These volatile biogenic compounds are likely to benefit local residents and reduce COVID-19 CFRs [18]. According to the classification of the World Bank, most of the low-income countries are in tropical or subtropical Africa, where the climate is mild, and evergreen forests play a dominant role. Therefore, residents in low-income countries may have more chances to be affected by forests in their residential neighborhoods.

CO 2 Emission and COVID-19 CFR
We did not observe a distinct association between CO 2 emission per capita and COVID-19 CFRs worldwide, which is inconsistent with a recent study, where data from 119 countries or regions indicated a nonlinear relationship between the CFR and CO 2 emission [27]. Some points may explain the inconsistency, such as the differences in countries or regions included, dates of data access, and the measures for CO 2 emission (per capita vs. per country). Nevertheless, our results may, to a certain extent, challenge the view that CO 2 emission is a critical factor likely to increase total coronavirus cases and death rates [27].
Aside from being a common environmental indicator [57], CO 2 emission may somewhat reflect economic conditions, as economic growth mostly depends on the consumption of energy resources [58]. For example, many studies have suggested the potential relationship between CO 2 emissions and GDP [59]. A recent study on COVID-19 has underlined the method of estimating the changes in CO 2 emission with the aid of GDP when detailed energy activity data are unavailable [60]. In the current study, we also observed a strong positive correlation between the GDP per capita and CO 2 emission (r = 0.89), implying that CO 2 emission may be a surrogate of GDP. For this reason, we conducted a sensitivity analysis, which includes GDP (per capita) instead of CO 2 emission in the model. We did not observe a distinct association between the GDP per capita and COVID-19 CFR in the global analysis either (Supplemental Table S2), suggesting that GDP might not be a critical factor associated with the COVID-19 CFR worldwide. This finding appears counter-intuitive, as GDP is associated with financial resources for health protection [29]. However, we should realize that pressures on medical and public health services due to COVID-19 varied among countries. For well-developed countries, such as Japan, healthcare services were affordable, and the medical system has remained functional throughout [61]. Therefore, other factors may play dominant roles in well-developed countries, such as age structure [62,63], which may mask the role of GDP in fighting the COVID-19 pandemic globally.
In the analysis stratified by income level, we found that CO 2 emission was associated with a decreased COVID-19 CFR in lower-middle-income countries, and the association remained statistically significant in the sensitivity analysis when CO 2 emission was replaced with GDP per capita (Supplemental Tables S2 and S3). This result is not surprising, given the reasons mentioned above. Compared to upper-middle and high-income countries, where medical and infrastructure resources are usually sufficient [64], lower-middle-income countries are more vulnerable to the pandemic [65], especially developing countries with poor economic conditions [66]. Therefore, a better economic condition may be more critical for lower-middle-income countries to fight against the pandemic.
In addition, the association between GDP and the COVID-19 CFR in low-income countries was ambiguous when the model was loaded with different forest variables (Supplemental Tables S2 and S3). We cannot conclude if GDP is an essential factor for low-income countries. Theoretically, they may face a similar crisis to lower-middle-income countries because resource shortage has also been reported as a common issue in lowincome countries [64]. Further studies are warranted to check their differences.

Limitations
There are several notable limitations to the current study. First, to indicate the typical natural environment worldwide, we used forest coverage rate and forest area per capita in the analysis. As discussed above, these measures may not necessarily indicate a potential exposure to nature. There is not a satisfactory measure to indicate the actual exposure to natural environments at a worldwide level. Even if the results were limited to the residential areas, common measures such as the normalized difference vegetation index or enhanced vegetation index may not explore the properties of natural environments, as well as their qualities, accessibility, and availability [67]. Future studies may refer to the data of population mobility, such as location reports of smartphone applications, and incorporate it with green space indicators to quantify actual green space usage.
Second, there is no globally accepted definition of a COVID-19-related death. Thus, reporting bias cannot be ruled out. Moreover, given the fact that medical standards can vary with regions, it is not clear if some low-income countries can effectively report COVID-19-related deaths. Thus, their CFR values may be underestimated.
Third, though we tried to access the latest data worldwide, there are some discrepancies in their timeframes. Due to the impact of COVID-19, the patterns of CO 2 emissions may be changed with time and region. Therefore, the identified relationships may need re-examination in further work.

Conclusions
This study aimed to examine the roles of forest and CO 2 emission on COVID-19 CFRs while controlling known confounders. We conducted a spatial analysis that included worldwide-level data obtained from public access databases. Although we did not observe a distinct association between the COVID-19 CFR and forest coverage rate or forest area per capita worldwide, we found that the forest coverage rate and forest area per capita were associated with lower COVID-19 CFRs in low-income countries, which is preliminarily speculated to be related to their lower urbanized levels.
On the other hand, we did not find a distinct association between CO 2 emissions per capita and COVID-19 CFRs at the worldwide level either. Nevertheless, we found that CO 2 emissions per capita were associated with lower CFRs in lower-middle-income countries. Notably, we assumed that CO 2 emission might be a surrogate of GDP as they showed a strong correlation. Similar results were obtained when we replaced the CO 2 emission per capita with GDP per capita in the regression. These findings suggest that economic conditions may be a critical factor for lower-middle-income countries to fight against the current pandemic.
Generally, our findings suggest that forest and CO 2 emissions are two potential ecological factors that may help to monitor and predict the global pandemic. Given the limitation of the public access data, further research is warranted to uncover the mechanisms behind the associations.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/f13050736/s1, Figure S1: Pairwise scatter plots and Pearson's correlation coefficients of the independent variables; Figure S2: Strength and decay of the spatial autocorrelation between pair of locations: Table S1: Changes in COVID-19 CFR (‰) per 10% increase in the ecological and socioeconomic variables (including forest area per capita instead of forest cover); Table S2: Changes in COVID-19 CFR (‰) per 10% increase in the ecological and socioeconomic variables (including GDP per capita instead of CO 2 emission per capita); Table S3: Changes in COVID-19 CFR (‰) per 10% increase in the ecological and socioeconomic variables (including GDP per capita instead of CO 2 emission per capita; and forest area per capita instead of forest cover).
Author Contributions: Conceptualization, Y.C.; methodology, Y.C.; data curation, Y.C.; writing-original draft preparation, Y.C. and H.L.; writing-review and editing, Y.C. and H.L.; supervision, Y.C. and G.Z.; project administration, Y.C. All authors have read and agreed to the published version of the manuscript.
Funding: This research received no external funding.

Conflicts of Interest:
The authors declare no conflict of interest.