Since the years 1850–1900, the rapid and relentless phenomenon of global warming has severely and negatively altered Earth’s environment and its inhabitants. NASA describes global warming as “the long-term heating of Earth’s surface observed since the pre-industrial period (between 1850 and 1900) due to human activities, primarily fossil fuel burning, which increases heat-trapping greenhouse gas levels in Earth’s atmosphere” ( NASA 2023). Carbon dioxide (CO ₂) is naturally present in the atmosphere as part of the Earth’s carbon cycle (the natural circulation of carbon among the atmosphere, oceans, soil, plants, and animals) ( EPA 2023). Unlike oxygen or nitrogen (which make up most of our atmosphere), greenhouse gases absorb heat radiating from the Earth’s surface and re-release it in all directions–including back toward the surface ( Lindsey 2022). In order to preserve Earth’s future, action must be taken now. The first step is tracking how atmospheric carbon dioxide emissions affect global warming rates.

The correlation between global temperature and carbon dioxide levels has been investigated by many researchers in the literature ( Florides and Christodoulides 2009; Macedo and Madaleno 2023; Hereher 2016; Palmer et al. 2007; Woodward and Gray 1993).

In the article, Global warming and carbon dioxide through science, ( Florides and Christodoulides 2009) revealed that there is no way to definitively say whether or not CO ₂ directly impacts global warming and temperature increases. This study used three independent sets of data (collected from ice cores and chemistry) to perform an analysis. Through a specific regression analysis of their data, they found that the data stating that there is a correlation between CO ₂ concentration and temperature relies heavily on specific choices of data. From this conclusion, they found that for both chemistry data and ice core data, “one cannot be positive that indeed such a correlation exists”. Through further research into the topic, there is evidence that CO ₂ is not inherently harmful to the environment, in addition to the fact that there is very little evidence linking CO ₂ levels and temperature. In the end, they concluded, “It is our view that it is not yet sufficient, let alone rigorous, evidence that anthropogenic CO ₂ increase is indeed the main factor contributing towards the global warming of the 20 th century”.

Evidence from a Maximum Entropy Approach ( Macedo and Madaleno 2023) uses a statistical approach based on maximum entropy to conduct a study that supports the results from different techniques that found that CO ₂ does, in fact, impact the increase in global temperature. This study, along with a recent review of emerging literature, led them to the conclusion that the impact of CO ₂ on our atmosphere is trending toward a detrimental conclusion. They further state that if people continue emitting CO ₂ into the atmosphere, it will have a detrimental impact on human, plant, and animal health in the long run.

A time series-related analysis was used in Hereher (2016) to study the global temperature trends of land surface temperatures and climate changes based on the temperature data in Egypt. This study reports serious climate change in Egypt by detecting variability of land surface temperature (LST) over the last decade at selected locations, with varied geomorphological characteristics and human stressors. The time series for the land surface temperatures were acquired from satellite images from 2003 to 2014, totaling 276 images. The analysis suggests that the LST in Egypt increases by 1.54 °C/decade. The variation of LST depends on latitude, geology, topography, and surface albedo. The dataset used in this analysis is named MODIS, from the NASA Land Processes Distributed Active Archive Center website. This analysis found that the time-series MODIS LST data proved sufficient for the short-term monitoring of land surface temperature variations in Egypt. It is noted that geology, topography, and surface albedo have significant impacts on the LST of Egypt. Further, it was revealed that between the years of 2003-2014, the LST of Egypt increased by 0.3–1.06 °C/decade, and for urban areas, the excess LST is higher at 1.54 °C/decade.

According to a new isothermal analysis by Palmer et al. (2007) in order to produce a more accurate depiction of the underlying warming. They proposed a new analysis of millions of ocean temperature profiles intended to filter out local dynamical changes. According to the authors’ comments, it was a difficult analysis since oceans do not warm uniformly across the globe. They present decadal-scale analyses of the ocean’s thermal state relative to a fixed isotherm. This new diagnostic is less prone to the influence of dynamical processes at both high and low frequencies, and the results present a more globally uniform picture of ocean warming. The limitation of the isothermal analysis is that high-latitude oceans were not included.

Using the Autoregressive Moving Average (ARMA) model, a study of time series data on global warming and trends is presented in Woodward and Gray (1993). They concluded that atmospheric greenhouse gases will affect the continued projection of the warming trend of global temperatures. Further, they stated that there is no conclusive evidence that this trend will continue. This is due to the difference between data on long-term trends and data on random trends with the length of the temperature series.

Global temperature anomalies and carbon dioxide data were obtained from ( NASA 2023) global climate change website. The temperature anomalies (in Celsius) and carbon dioxide (parts per million) were used from 1960 to 2022. The temperature anomalies are used instead of the absolute temperature data, as they accurately represent the temperature variability over larger areas. Further, they give a frame of reference that allows more meaningful comparisons between locations and more accurate calculations of temperature trends ( NOAA 2023).

This section provides an overview of the definitions and test formulas that we use in this article to examine the time series models. Time series analysis can be used to examine the outcomes of either a planned or unplanned intervention as well as to better understand the underlying naturalistic process and the pattern of change over time ( Velicer and Molenaar 2012). In time series analysis, regularity is essential, and if time series data exhibits irregular behavior over time, establishing meaningful conclusions about the underlying causes will be challenging. This regular behavior is explained by the concept of stationarity. In general, to satisfy the conditions of stationarity, the time series must satisfy some strong conditions, and those conditions are too strong for most of the applications. Therefore, most time series seen in real life are analyzed using a weaker version of stationarity.

The weakly stationary time series should satisfy the following conditions ( Shumway and Stoffer 2017). (1)

The expected value (mean) of the time series is constant and does not depend on time.

If two time series are taken into account simultaneously, such as x t and y t , they are said to be jointly stationary if the individual time series are stationary and the cross-covariance function solely depends on the lag value.

To evaluate the model’s goodness of fit, we utilize Akaike’s Information Criterion (AIC), the bias-corrected Akaike’s Information Criterion (AICc), and Bayesian Information Criterion (BIC). Those criteria measure the goodness of fit of the models by balancing the number of parameters of the model and the error of the fit. If k is the number of parameters of the model, n is the number of observations, and SSE is the sum squared error of the fit, then AIC is given by AIC = ln SSE / n + n + 2 k n . (1)

The AICc is given by AICc = ln SSE / n + n + k n − k − 2 . (2)

With the Bayesian corrected term, the BIC takes the form BIC = ln SSE / n + k ln n n . (3)

More details about those information criteria can be found in Akaike (1974, 1973, 1969) and Schwarz (1978).

Detrending and differencing are the two fundamental methods for transforming a non-stationary time series into a stationary one. The independence of the residuals can be visually demonstrated by using the residual plots. The autocorrelation plots and cross-correlation plots can be used the explain the nature of the correlation at different lag values. If μ x , t and μ y , t are the means of the time series x t and y t respectively, the autocovariance function of x t is defined by cov x t + h x t = E x t + h − μ x , t + h x t − μ x , t , (4) where h represents the time shift or lag value. The cross-covariance function for x t and y t is given by cov x t + h y t = E x t + h − μ x , t y t − μ y , t (5)

The autocorrelation and cross-correlation functions are normalized versions of the above formulas. For a given time series x t , and y t the autocorrelation (ACF) and cross-correlation (CCF) functions are given by ACF x s x t = cov x s x t cov x s x s cov x t x t , (6) and CCF x s y t = cov x s y t cov x s x s cov y t y t , (7) respectively.

Because they depend on the locations of time points s and t , the autocovariance and cross-covariance functions can change during the course of the series. If the autocovariance function depends on the separation h = | t − s | rather than the points where the time series are situated, we are able to analyze sample time series data when there is only one series available.

We use the global temperature anomalies from 1960 to 2022, obtained from NASA (2023). As demonstrated in Figure 1 most of the planet is warming (yellow, orange, and red). Only a few locations, most of them in the southern hemisphere oceans, cooled over this time period. According to the fact highlighted in Lindsey and Dahlman (2023) earth’s temperature has risen by an average of 0.14° Fahrenheit per decade since 1880. The temperature anomaly variability (in Celsius), including pre-industrial time, is shown in Figure 2. The trend is shown by the red dashed line while smoothing splines with the smoothing parameter value spar is 0.5, is demonstrated by the blue spline.

The trend and irregular behavior of the temperature series require a modification to make it stationary before additional analysis can be done to uncover its other properties. The stationary time series are easy to study using established principles created in the time series literature because of their predictable long-term behavior.

We investigate detrending and differencing techniques to transform the original time series into stationary time series.

We analyze the basic global temperature trend model and temperature and carbon dioxide models. We compare two models using AIC, AICc, BIC, R-squared value, and sum squared error (SSE). The main objective of these models is to further analyze the effect of carbon dioxide on global temperatures. The basic trend model for temperature (model 1): T = β 0 + β 1 t + w t (19) where T is the global temperature anomaly (°C) and t is the time in years ( 1960 ≤ t ≤ 2022 ), w t is the white noise and is assumed to be normally distributed with mean zero and fixed variance (say) σ 2 .

The trend model (model 1) in Eq. 19 can be improved by introducing CO ₂ as an explanatory variable.

Figure 9 shows the autocorrelation values for the carbon dioxide series at different lag values. Due to the considerable strength of the lag-related correlation in the second model, we adjust the carbon dioxide data for its mean, C ¯ = 358.9675 ( ppm). If C t is the carbon dioxide level at time t , then we have model 2: T = β 0 ′ + β 1 ′ t + β 2 ′ C ′ + w t (20) where C ′ = C t − C ¯ , 1960 ≤ t ≤ 2022 , is the adjusted carbon dioxide level (in ppm). If k is the number of parameters of the model and df is the degree of freedom, Table 2 shows the summary statistics for the two models.

The value k yielding the minimum AIC, AICc, and BIC specifies the best model. We note that the second model is substantially better than the first model. Model 2, which includes CO ₂ accounting for 92.4% of the variability ( R 2 value of model 2), which was 90.5% ( R 2 value of model 1) without that. Also, it gives the best value for AIC and BIC. In addition, we can notice that AIC and AICc are nearly equal. To calculate those values, we used the formulas given in the equations 1, 2, and 3. Without using the aforementioned formulas, one can utilize the regression model summaries to obtain these values, although there might be a few minor variations from the corresponding values shown in Table 2.

Further, the trend model (model 1) can be compared to the model with carbon dioxide (model 2) with the null hypothesis, H 0 : β 3 = 0 . The corresponding F statistic: F = SSE r − SSE / q − r SSE / n − q − 1 (21) where SSE r is the sum squared error of the reduced model, and q and r are the numbers of predictor variables in the full and reduced models, respectively. When q = 2 , r = 1 , and n = 63 F = 0.632 − 0.495 / 1 0.495 / 60 = 16.606 , (22) which exceeds F 1 , 60 0.001 = 11.973 . Hence, model 2 is a better prediction model compared to model 1.

Further for the purpose of predicting the global temperature the prediction model can be given by, T ̂ t = 4.943 − 0.002 t + 0.012 C t − 358.968 (23) where T ̂ t is the estimated temperature (anomaly) at time t . A negative (but very small) weight is present as the time coefficient. But a relatively larger constant ( β 0 ) value mitigates the effect of the negative weight. The positive weight of the carbon dioxide indicates a positive contribution to the global temperature when C t > 358.968 ( ppm). Because of this, whether carbon dioxide has a positive or negative effect on the rise in global temperatures depends on its relative value to the average annual carbon dioxide level (in this study, the average was computed throughout the years from 1960 to 2022).

Analyzing the potential leading and lagging relationships between the global temperature and carbon dioxide series is another intriguing investigation. Leading and lagging relationships might be helpful when one time series is used to predict another. Let T t and C t represent the global temperature and CO ₂ time series respectively. Consider the model of the form T t = β t l C t − l + w t (24) where β t l is a real value that depends on t or l . C t lead T t if l > 0 , and lag T t if l < 0 . If w t is uncorrelated with C t time series, the cross-correlation is given by cov T s C t = cov β t l C s − l + w s C t = cov β t l C s − l C t , t > 0 , s > 0 . (25)

According to the graphical demonstration of the cross-correlation function in Figure 11a for the original data, no time series leads or lags other series since peak shows occur at l = 0 . When viewed in relation to the zero-lag value, the cross-correlation is approximately symmetric. However, those original series are not stationary. We take integration order 2 into consideration when we analyze the behavior of the cross-correlation functions of the (weakly) stationary CO ₂ and global temperature anomaly series. Although the first integrated temperature series is stationary, we use the differencing order 2 ( d = 2 ) to transform both series in order to tackle the dimension issues. According to Figure 11b most dominant cross-correlation for stationary series occurs at l = 0 and l = 1 . Those integrated series are jointly stationary if the cross-correlation or cross-covariance function is a function only on the lag value.

According to the summary statistics criterion in Table 2, the impact of CO ₂ on the global temperature cannot be negligible. The models we’ve examined can be improved by adding additional greenhouse gases and other factors that could have an impact on the global temperature. Although carbon dioxide is not the only greenhouse gas that can affect the temperature, it occupies a significant amount of space compared to other greenhouse gases.