Global Warming Coming? ARIMA & Prophet

Project Proposal

This project aim to study the time series of the Earth's temperature in order to address the issue of global warming. The study of the earth's temperature is a time series depending on many physical variables, and each of these physical variables is a function of time. The following project we want to use the ARIMA (Autoregressive integrated moving average) method to study the behavior of the Earth's temperature.

Libraries

import numpy as np import pandas as pd import matplotlib.pyplot as plt %matplotlib inline import plotly.offline as py import plotly.graph_objs as go import plotly.tools as tls import seaborn as sns import time import warnings

## library import import pandas as pd from datetime import datetime from matplotlib import pyplot as plt

# read temp data temp = pd.read_csv("GlobalTemperatures.csv", \ dtype = {'Data' : str}) global_temp_country = pd.read_csv('GlobalLandTemperaturesByCountry.csv')

global_temp_country_clear = global_temp_country[~global_temp_country['Country'].isin( ['Denmark', 'Antarctica', 'France', 'Europe', 'Netherlands', 'United Kingdom', 'Africa', 'South America'])]

global_temp_country_clear = global_temp_country_clear.replace( ['Denmark (Europe)', 'France (Europe)', 'Netherlands (Europe)', 'United Kingdom (Europe)'], ['Denmark', 'France', 'Netherlands', 'United Kingdom'])

# reset date format global_temp_country_clear['Date'] = global_temp_country_clear["dt"].apply(lambda x : datetime.strptime(x, '%Y-%m-%d')) global_temp_country_clear = global_temp_country_clear.drop(columns=['dt']) temp_country = global_temp_country_clear.dropna() temp_country

split_date ='1952-01-01' temp_beforeIndus = temp_country[temp_country['Date'] < split_date] temp_AfterIndus = temp_country[temp_country['Date'] >= split_date]

After Discovering its hard to distinguish the temperature difference between heat plot of two periods, we decide to directly plot the temperature difference.

countries = np.unique(temp_country['Country']) mean_temp_before = [] mean_temp_after = [] mean_temp = [] for country in countries: difference = temp_AfterIndus[temp_AfterIndus['Country'] == country]['AverageTemperature'].mean() - temp_beforeIndus[temp_beforeIndus['Country'] == country]['AverageTemperature'].mean() mean_temp.append(difference) data = [ dict( type = 'choropleth', locations = countries, z = mean_temp, locationmode = 'country names', text = countries, marker = dict( line = dict(color = 'rgb(0,0,0)', width = 1)), colorbar = dict(autotick = True, tickprefix = '', title = 'Average\nTemperature,\n°C') ) ] layout = dict( title = 'Average land temperature Difference Before& After Industry Revolution', geo = dict( showframe = False, showocean = True, oceancolor = 'rgb(0,255,255)', projection = dict( type = 'orthographic', rotation = dict( lon = 60, lat = 10), ), lonaxis = dict( showgrid = True, gridcolor = 'rgb(102, 102, 102)' ), lataxis = dict( showgrid = True, gridcolor = 'rgb(102, 102, 102)' ) ), ) fig = dict(data=data, layout=layout) py.iplot(fig, validate=False, filename='worldmap')

# reset date format temp['Date'] = temp['dt'].apply(lambda x : datetime.strptime(x, '%Y-%m-%d')) # drop the old date(back the old day, I still hang on in the street) temp = temp.drop(columns=['dt']) temp['year'] = pd.DatetimeIndex(temp['Date']).year temp['month'] = pd.DatetimeIndex(temp['Date']).month # ser datetime as index temp = temp.set_index('Date')

# general view of the data set temp.head(10)

# general view of the data set print("Types of values in temp : ") print(temp.dtypes) print("Shape of values in temp : ") print(temp.shape) temp.describe()

## A General view of temp change (land and ocean) every year since 1769 to 2015 # get the year average temperature temp_avg = temp.groupby('year').mean()

# plot a general view of year average tempature bewteen 1750 to 2015 ax = temp_avg['LandAverageTemperature'].plot() plt.grid() plt.ylabel("year avg temp") plt.title('Yearly land average temperature')

Yearly land average temperature show a high volatility before 1950. But since about 1950, there is an obvious increasing trend, which is also corresponded to the trend of carbon dioxide in the figure before. From this figure, we can know that since the industrial revolution, there is a reasonable evidence of the global warming.

# some basic operations temp_avg['LandAverageTemperatureUncertainty_max'] = temp_avg['LandAverageTemperature'] + temp_avg['LandAverageTemperatureUncertainty'] temp_avg['LandAverageTemperatureUncertainty_min'] = temp_avg['LandAverageTemperature'] - temp_avg['LandAverageTemperatureUncertainty']

# plot a general view of year average tempature bewteen 1750 to 2015 with the Uncertainty fig, ax = plt.subplots() ax.plot(temp_avg['LandAverageTemperature'], color='blue') ax.plot(temp_avg['LandAverageTemperatureUncertainty_max'], color='red', label = 'MaxUncertanity') ax.plot(temp_avg['LandAverageTemperatureUncertainty_min'], color='green', label = 'MinUncertanity') ax.set_ylabel('year avg temp') ax.set_title('Yearly land average temperature with uncertanity') ax.legend()

This figure shows with the improvement of detecting techonology, the uncertainty of land average temperature shows a obvious increase.

# get the month average temperature temp_mon_avg = temp.groupby('month').mean()

# plot a general view of month average tempature between 1750 to 2015 ax = temp_mon_avg['LandAverageTemperature'].plot() plt.grid() plt.ylabel("month avg temp") plt.title('Monthly land average temperature')

This figure shows the seansonal trend of global average temperature. The highest temperature happens on July, and the lowest average temperature occoured in January.

pip install statsmodels

# Decompose data from statsmodels.tsa.seasonal import seasonal_decompose temp_de = temp_avg[['LandAverageTemperature']] # drop na temp_de = temp_de.dropna() # result = seasonal_decompose(temp_de, model='multiplicative', period = 12)

result.plot() plt.xlabel('year') plt.show()

Through the seasonal_decompose package, the trend and seasonal figure prove our conclusion before. 1. the increasing trend of global average temperature. 2. the seansonal trend of land temperature.

import warnings warnings.filterwarnings("ignore")

temp = temp.dropna()

Drop NA value for mahcine learning model training.

Before modeling, we first check whether it is stationary. Stationary means mean, variance, covariance is constant over periods. From the seasonal plot, we can observe the temperature is roughly the same in each period.

from statsmodels.graphics.tsaplots import plot_acf

plot_acf(temp_avg[['LandAverageTemperature']])

The acf plot shows the correlation between the same series with a time lag. When we have one lag the correlation factor is 1, two lag is 0.74...... As we can see the correlation slowly decay, so it is not stationary. It should ideally decay immediately to 0 or negative, showing that the data is stationary.

avg_diff = temp['LandAverageTemperature'].diff(periods = 1)

avg_diff = avg_diff[1:] avg_diff.head()

plot_acf(avg_diff)

The correlation decay to 0, so it is good stationary.

avg_diff.plot()

We can see there is no trend here. It is stationary. This data is good enough to move forward to modeling.

#split_date = '2011-01-01' train_temp = temp[:'2013-1-1'] test_temp = temp['2013-1-1':] predictions = []

pip install pmdarima

Find the best parameter p d q of Arima model

from pmdarima import auto_arima stepwise_fit = auto_arima(train_temp['LandAverageTemperature'], trace=True, suppress_warnings=True)

Print out every combination of the parameter. AIC is the score of every model, and we want to minimize the AIC.

pip install statsmodels

ARIMA(5,1,5)(0,0,0)[0] intercept : AIC=1760.698, Time=11.88 sec

We find the best parameter is 5 1 5

We can see there is no seasonality (0,0,0) in dataset, we can use this Arima model.

import warnings warnings.filterwarnings("ignore")

Train the Model

from statsmodels.tsa.arima.model import ARIMA model=ARIMA(train_temp['LandAverageTemperature'],order=(5,1,5)) #fit model model=model.fit() #model.summary()

Make prediction on start - end

start = len(train_temp) - 1 end = len(train_temp) + len(test_temp) - 2 pred = model.predict(start = start, end = end, type = 'level') #print(pred)

we need to keep date since we want to make plot in future

Plot Result

pred.plot(legend=True) test_temp['LandAverageTemperature'].plot(legend=True)

We plot the prediction and actual data in the same graph, we can see they are pretty close.

Calculate the Error

from sklearn import metrics from sklearn.metrics import mean_squared_error mse = mean_squared_error(test_temp['LandAverageTemperature'], pred) print(mse)

print(metrics.mean_absolute_error(test_temp['LandAverageTemperature'], pred))

print(np.sqrt(metrics.mean_squared_error(test_temp['LandAverageTemperature'], pred)))

The error is pretty small, so the model did good prediction

start = len(train_temp) end = len(train_temp) + len(test_temp) + 10 pred = model.predict(start = start, end = end, type = 'level') #print(pred)

Predict for future date

Now we want to use the model trained to predict for the future.

pred.plot(legend=True) test_temp['LandAverageTemperature'].plot(legend=True)

start = len(train_temp) - 1 end = len(train_temp) + len(test_temp) + 480 pred = model.predict(start = start, end = end, type = 'level')

We use this model to predict for the next 40 years

pred_df = pred.to_frame()

pred_df['Date'] = pred_df.index

pred_df['year'] = pd.DatetimeIndex(pred_df['Date']).year pred_avg = pred_df.groupby('year').mean()

ax = pred_avg['predicted_mean'].plot() plt.ylabel("year avg temp") plt.title('Yearly land average temperature')

The prediction model did good prediction up to year 2017. After year 2017, the time lag is too long and model doesn't work well. So in the future, we would try to figure out how to predict for long period of time.

Based on our current prediction, in year 2017, the average temperature will increase from 9.560 to 9.574.

For comparing the predicting accuracy better, we employed another time series model - Prophet.

We build the prophet model first, and let porphet do the same job as the ARIMA did before.

''' library import ''' from prophet import Prophet import pandas as pd from datetime import datetime from matplotlib import pyplot as plt from sklearn.metrics import mean_squared_error # read the data temp = pd.read_csv("GlobalTemperatures.csv", \ dtype = {'Data' : str}) # reset date format temp['Date'] = temp['dt'].apply(lambda x : datetime.strptime(x, '%Y-%m-%d')) # drop the old date(back the old day, I still hang on in the street) temp = temp.drop(columns=['dt']) temp['year'] = pd.DatetimeIndex(temp['Date']).year temp['month'] = pd.DatetimeIndex(temp['Date']).month # ser datetime as index temp = temp.set_index('Date') # set ds and y for prophet temp['ds'] = temp.index temp['y'] = temp['LandAverageTemperature'] temp_new = temp[['ds', 'y']] # set train and test data set train_temp = temp_new[:'2013-1-1'] test_temp = temp_new['2013-1-1':] # build prophet model model = Prophet(changepoint_prior_scale=20) model.fit(train_temp) # build the future dataset for Prophet future = model.make_future_dataframe(periods=30, freq='M') fax = model.predict(future) y_true = temp['y']['2013-1-1':].values y_pred = fax['yhat'][-36:].values mse = mean_squared_error(y_true, y_pred) print(mse)

The MSE of Prophet is 0.1141353608465227

We used Prophet tot predict the temperature between Jan, 2013 to Dec, 2015 (Same as the work in ARIMA). Compare the predict data and the ture data, Prophet hold a accuracy as 0.1141353608465227 MSE. It's reasonable to believe that under such a big size data set, the prediction of Prophet got a good result.

Compare this value with the accuracy in ARIMA (MSE with 0.15716738194270644), Prophet got a bit higher accuracy.

plt.plot(y_true, label='Actual') plt.plot(y_pred, label='Predicted') plt.legend() plt.show() future2 = model.make_future_dataframe(periods=47, freq='M') fax2 = model.predict(future2) y_true = temp['y']['2013-1-1':].values y_pred = fax2['yhat'][-48:].values plt.plot(y_true, label='Actual') plt.plot(y_pred, label='Predicted') plt.legend() plt.show()

We built the time series trend of global avg temperature in the actual data and the Prophet predict data.

From the first figure, it can be observed that the Porphet could capture the chaning trend of global temperature well. There is no obvious gap between the predict and the true value.

From the second figure, we additionally forecast the global land average temperature in 2016. The seasonal trend of temperaturein 2016 also be captured well.

In this project, we made the heatmap of the global land avg temperature. From that, we found that there was an obvious increase of global temperature since the industrial revolution. And we also made the time series plot and the trend plot of land avg temperature to prove this point. After data processing, we employed the ARIMA model and Prphet model to predict the temperature from 2013 to 2015. Through comparing the true value, both model achieved a good result with MSE lower than 0.2. The accuracy in Prophet is a bit higher, but both result were reasonable and acceptable. Additionally, we used the two models to forecast the temperature in 2016. The results are exciting because both two models could capture the seasonal trend well in 2016.

Furthermore, we employed the ARIMA to forecast the global avg temperature until 2055. ARIMA captured a basic increasing trend before about 2016, but after that, ARIAM cannot capture the changing trend anymore. Which means our model may not have a good ability for long term forecasting. This is an unresolved problem and may need our further research.

And also for futher research, maybe combine different machine learning models may help us to get a better predicting accuracy.

.css-15w88e5{color:var(--chakra-colors-fg-neutral-primary);font-weight:inherit;letter-spacing:-0.09px;}Project Proposal