\n\n1. [Glimpse of Data](#PREPARATION)\n * [Import packages](#IMPORT)\n2. [Reading S5p data Whether and No2](#READS5P)\n3. [Visuals](#V1)\n4. [Anomaly Analysis](#OUTLIER)\n 4.1 [Gaussian](#OUTLIER1)\n 4.2 [Isolation Forest](#OUTLIER2)\n 4.3 [One Class SVM](#OUTLIER3)\n\n5. [Prediction using LSTM with Python](#LSTM1)\n 5.1 [Get the root mean squared error (RMSE)](#LSTM1)\n\n6. [Arima with Python](#AR)\n 6.1. [Rolling Forecast ARIMA Model](#AR2)\n\n7. [Time series prediction using Prophet in Python](#PRO)\n 7.1. [Forecast quality evaluation](#PRO1)\n 7.2. [Incorporating the Effects of Weather Condition](#PRO2)\n 7.3. [Forecast quality evaluation](#PRO3)\n 7.4. [Save The Model](#PRO4)\n8. [Prediction of NO2 density for each primary_fuel throughout the year](#M1)\n 8.1. [Forecast quality evaluation for Power Plant over the year](#M2)\n 8.3. [Outlier Analysis of Power Plant - Coal over the year](#M3)\n \n9. [About the data](#ABOUTDATA) \n10. [Ending note](#END) \n\n "},{"metadata":{},"cell_type":"markdown","source":"## In this notebook we investigated the presence of NO2 concentration in air, considering its constant increase over days, years. Owing to accurate future air quality estimates, the need for detecting the anomalously high increase in the concentration of pollutants cannot be adjourned. This study is helpful in educating the government for decision making and people about spatiotemporal, geographical, and economic conditions responsible for anomalously high NO2 concentrations in air. In this work, we modeling the solution and analyze the impacts of air pollution for each region in Porto Rico for each primary_fuel in the year."},{"metadata":{},"cell_type":"markdown","source":"
Glimpse of Data - Power Plants
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{},"cell_type":"markdown","source":"#
Imports
\n\n[Back to Table of Contents](#top)\n\nWe are using a stack: ``numpy``, ``pandas``, ``sklearn``, ``matplotlib``, ``rasterio``, ``plotly``."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"import os, random\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\npd.set_option('max_columns', 200)","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from sklearn.covariance import EllipticEnvelope\nfrom sklearn.ensemble import IsolationForest\nfrom sklearn.svm import OneClassSVM\nfrom sklearn import preprocessing","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"# !pip install plotly\n# !pip install fbprophet","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"from keras.models import Sequential\nfrom keras.layers import Dense, LSTM\n\nimport datetime as dt\nfrom statsmodels.tsa.arima_model import ARIMA\n\nfrom fbprophet import Prophet\nfrom fbprophet.plot import plot_plotly\n\nimport plotly.offline as py\nfrom matplotlib import pyplot\npy.init_notebook_mode()\n","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"%%HTML\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Reading S5p data Whether and No2
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n \nYou can verify the data with more details here: [Prepare Data for Modeling](https://www.kaggle.com/caesarlupum/ds4g-go-to-the-green-future#-Satellite-Information)"},{"metadata":{"_cell_guid":"","_uuid":"","trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"import pandas as pd\nno2_weather = pd.read_csv(\"../input/s5p-data-csv/no2_weather.csv\")\ns5p_no2_pictures_df = pd.read_csv(\"../input/s5p-data-csv/s5p_no2_pictures_df.csv\")\nweather_pictures_df = pd.read_csv(\"../input/s5p-data-csv/weather_pictures_df.csv\")","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"import rasterio as rio\ndef split_column_into_new_columns(dataframe,column_to_split,new_column_one,begin_column_one,end_column_one):\n for i in range(0, len(dataframe)):\n dataframe.loc[i, new_column_one] = dataframe.loc[i, column_to_split][begin_column_one:end_column_one]\n return dataframe","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Power plants on Puerto Rico"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"power_plants = pd.read_csv('/kaggle/input/ds4g-environmental-insights-explorer/eie_data/gppd/gppd_120_pr.csv')\npower_plants = split_column_into_new_columns(power_plants,'.geo','latitude',50,66)\npower_plants = split_column_into_new_columns(power_plants,'.geo','longitude',31,48)\npower_plants['latitude'] = power_plants['latitude'].astype(float)\na = np.array(power_plants['latitude'].values.tolist()) # 18 instead of 8\npower_plants['latitude'] = np.where(a < 10, a+10, a).tolist() \n\npower_plants_df = power_plants.sort_values('capacity_mw',ascending=False).reset_index()\npower_plants_df['img_idx_lt']=(((18.6-power_plants_df.latitude)*148/(18.6-17.9))).astype(int)\npower_plants_df['img_idx_lg']=((67.3+power_plants_df.longitude.astype(float))*475/(67.3-65.2)).astype(int)\npower_plants_df['plant']=power_plants_df.name.str[:3]+power_plants_df.name.str[-1]+'_'+power_plants_df.primary_fuel\npower_plants=power_plants_df[['name','latitude','longitude','primary_fuel','capacity_mw','img_idx_lt','img_idx_lg','plant']]\npower_plants","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"data.info()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
drop nan values
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"data.dropna(axis=0, inplace=True)\nprint(data.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"data = data.set_index('start_date')\ndata.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"weather_pictures_df.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"weather_pictures_df.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"no2_weather['start_date'] = pd.to_datetime(no2_weather['start_date'])\nno2_weather['no2_emission_sum'] = (no2_weather['no2_emission_sum'] - 32) * 5/9\n# plot the data\nno2_weather.plot(x='start_date', y='no2_emission_sum')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Anomaly Analysis
\n\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{},"cell_type":"markdown","source":"Owing to accurate future air quality estimates, need for detecting the anomalously high increase in concentration of pollutants cannot be adjourned. The presence of NO2 concentration in air is investigated in this notebook, considering its constant increase over years as well as its inevitable health risks. Furthermore, spatiotemporal segments with anomalously high NO2 concentrations for Porto Rico.\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Suppress warnings \nimport warnings\nwarnings.filterwarnings(\"ignore\", category=DeprecationWarning)\nwarnings.filterwarnings(\"ignore\", category=UserWarning)\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\nfrom IPython.display import HTML\n\nHTML('')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Anomaly Detection | Developing And Evaluating An Anomaly Detection System"},{"metadata":{},"cell_type":"markdown","source":"
Feature engineering
\nExtracting some features\n"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"# An estimation of anomly population of the dataset (necessary for several algorithm)\noutliers_fraction = 0.01\n\nno2_weather['day'] = no2_weather['start_date'].dt.day\n# the day of the week (Monday=0, Sunday=6) and if it's a week end day or week day.\nno2_weather['DayOfTheWeek'] = no2_weather['start_date'].dt.dayofweek","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Suppress warnings \nimport warnings\nwarnings.filterwarnings(\"ignore\", category=DeprecationWarning)\nwarnings.filterwarnings(\"ignore\", category=UserWarning)\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\nfrom IPython.display import HTML\n\nHTML('')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Anomaly Detection Problem | Motivation — [ Machine Learning | Andrew Ng ]"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Take useful feature and standardize them \ndata_SVM = no2_weather[['time_epoch','DayOfTheWeek','day','no2_emission_sum','temperature_2m_above_ground_mean','specific_humidity_2m_above_ground_mean','relative_humidity_2m_above_ground_mean','u_component_of_wind_10m_above_ground_mean','v_component_of_wind_10m_above_ground_mean','total_precipitation_surface_mean']]\ndata_SVM.dropna(axis=0, inplace=True)\nprint(data_SVM.shape)\n\n\nmin_max_scaler = preprocessing.StandardScaler()\nnp_scaled = min_max_scaler.fit_transform(data_SVM)\n# train one class SVM \nmodel = OneClassSVM(nu=0.95 * outliers_fraction) #nu=0.95 * outliers_fraction + 0.05\ndata_SVM = pd.DataFrame(np_scaled)\nmodel.fit(data_SVM)\n# add the data to the main \nno2_weather['anomaly26'] = pd.Series(model.predict(data_SVM))\nno2_weather['anomaly26'] = no2_weather['anomaly26'].map( {1: 0, -1: 1} )\nprint(no2_weather['anomaly26'].value_counts())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"\n
"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"timeLags = np.arange(1,10*48*7)\nautoCorr = [no2_weather.no2_emission_sum.autocorr(lag=dt) for dt in timeLags]\nplt.figure(figsize=(19,8))\nplt.plot(1.0/(48*7)*timeLags, autoCorr)\nplt.xlabel('time lag [weeks]')\nplt.ylabel('correlation coeff', fontsize=12)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"AutoCorrelation 10 weeks depth\n\nThe NO2 demand seems to be driven by a weekly trend: on certain days of the week, is higher than the others. We simply prove this computing autocorrelation."},{"metadata":{},"cell_type":"markdown","source":"
Prediction using LSTM with Python
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Suppress warnings \nimport warnings\nwarnings.filterwarnings(\"ignore\", category=DeprecationWarning)\nwarnings.filterwarnings(\"ignore\", category=UserWarning)\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\nfrom IPython.display import HTML\n\nHTML('')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Recurrent Networks can be improved to remember long range dependencies by using whats called a Long-Short Term Memory (LSTM) Cell. Let's build one using just numpy! I'll go over the cell components as well as the forward and backward pass logic\n "},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"data_ = no2_weather.loc[:,['start_date','no2_emission_sum']] \ndata_['start_date'] = pd.to_datetime(no2_weather['start_date'])\ndata_.set_index('start_date', inplace=True)\ndata_ = data_.resample(\"1D\").sum() # day sum","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"#Create a new dataframe with only the 'no2_emission_sum column\ndata_2 = data_.filter(['no2_emission_sum'])\n#Convert the dataframe to a numpy array\ndataset = data_2.values\n#Get the number of rows to train the model on\ntraining_data_len = int(np.ceil( len(dataset) * .8 ))\ntraining_data_len","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"#Scale the data\nfrom sklearn.preprocessing import MinMaxScaler\nscaler = MinMaxScaler(feature_range=(0,1))\nscaled_data = scaler.fit_transform(dataset)\n# scaled_data","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Create the training data set andC create the scaled training data set"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train_data = scaled_data[0:int(training_data_len), :]\n#Split the data into x_train and y_train data sets\nx_train = []\ny_train = []\n\nfor i in range(129, len(train_data)):\n x_train.append(train_data[i-129:i, 0])\n y_train.append(train_data[i, 0])\n# if i<= 61:\n# print(x_train)\n# print(y_train)\n# print()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"# Convert the x_train and y_train to numpy arrays \nx_train, y_train = np.array(x_train), np.array(y_train)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"#Reshape the data\nx_train = np.reshape(x_train, (x_train.shape[0], x_train.shape[1], 1))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Build the LSTM model"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"model = Sequential()\nmodel.add(LSTM(50, return_sequences=True, input_shape= (x_train.shape[1], 1)))\nmodel.add(LSTM(50, return_sequences= False))\nmodel.add(Dense(25))\nmodel.add(Dense(1))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Compile and Train the model"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"\nmodel.compile(optimizer='adam', loss='mean_squared_error')\nmodel.fit(x_train, y_train, batch_size=1, epochs=1)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"#Create the testing data set\ntest_data = scaled_data[training_data_len - 129: , :]\n#Create the data sets x_test and y_test\nx_test = []\ny_test = dataset[training_data_len:, :]\nfor i in range(129, len(test_data)):\n x_test.append(test_data[i-129:i, 0])\n \n# Convert the data to a numpy array\nx_test = np.array(x_test)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"\n
Get the root mean squared error (RMSE)
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\n"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"# Reshape the data\nx_test = np.reshape(x_test, (x_test.shape[0], x_test.shape[1], 1 ))\n# Get the models predicted values \npredictions = model.predict(x_test)\npredictions = scaler.inverse_transform(predictions)\n# Get the root mean squared error (RMSE)\nrmse = np.sqrt(np.mean(((predictions - y_test) ** 2)))\nrmse","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Visualize the data\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Plot the data\ntrain = data_[:training_data_len]\nvalid = data_[training_data_len:]\nvalid['Predictions'] = predictions","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Visualize the data\nplt.figure(figsize=(16,8))\nplt.title('Model')\nplt.xlabel('start_date', fontsize=18)\nplt.ylabel('NO2 Emission Sum', fontsize=18)\nplt.plot(train['no2_emission_sum'])\nplt.plot(valid[['no2_emission_sum', 'Predictions']])\nplt.legend(['Train', 'Val', 'Predictions'], loc='lower right')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"\n\n
30 days NO2 average
\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"f, ax = plt.subplots(figsize=(14,8))\npd.plotting.register_matplotlib_converters() # Add this \ndata_.plot(ax=ax, color='C0')\ndata_.rolling(window=30, center=True).mean().plot(ax=ax, ls='-', lw=3, color='C3')\nax.grid(ls=':')\nax.legend(['daily values','30 days No2 average'], frameon=False, fontsize=14)\n\n[l.set_fontsize(13) for l in ax.xaxis.get_ticklabels()]\n[l.set_fontsize(13) for l in ax.yaxis.get_ticklabels()]\nax.set_xlabel('date', fontsize=15)\nax.set_ylabel('N02 values', fontsize=15);\n# ax.axvline('2018', color='0.8', lw=8, zorder=-1)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"high NO2 mean in September"},{"metadata":{},"cell_type":"markdown","source":"
ARIMA with Python
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n "},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Suppress warnings \nimport warnings\nwarnings.filterwarnings(\"ignore\", category=DeprecationWarning)\nwarnings.filterwarnings(\"ignore\", category=UserWarning)\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\nfrom IPython.display import HTML\n\nHTML('')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The Analysis of Time Series"},{"metadata":{},"cell_type":"markdown","source":"Printing a summary of the fit model.\nThis summarizes the coefficient values used as well as the skill of the fit on the on the in-sample observations."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"from statsmodels.tsa.arima_model import ARIMA\nfrom sklearn.metrics import mean_squared_error\n\narimaM = ARIMA(data, order=(5,1,0))\narimaMfit = arimaM.fit(disp=0)\nprint(arimaMfit.summary())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We get a line plot of the residual errors, suggesting that there may still be some trend information not captured by the model and \n\nwe get a density plot of the residual error values, suggesting the errors are Gaussian, but may not be centered on zero."},{"metadata":{},"cell_type":"markdown","source":"
# plot residual errors
\n\n "},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"errors = pd.DataFrame(arimaMfit.resid)\nerrors.plot()\npyplot.show()\nerrors.plot(kind='kde')\npyplot.show()\nprint(errors.describe())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The distribution of the residual errors is displayed. \nThe results show that indeed there is a bias in the prediction (a non-zero mean in the residuals)."},{"metadata":{},"cell_type":"markdown","source":"\n
Rolling Forecast ARIMA Model
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"X = data.values\nsize = int(len(X) * 0.70)\nlimitCount = 40\ntrain, test = X[0:size], X[size:size+limitCount]\nhistory = [x for x in train]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can also calculate a final mean squared error score (MSE) and (RMSLE) for the predictions, providing a point of comparison for other ARIMA configurations."},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"def rmsle(y, y_pred):\n return np.sqrt(mean_squared_error(y, y_pred))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"pred = []\nfor t in range(len(test)):\n model = ARIMA(history, order=(5,1,0))\n model_fit = model.fit(disp=0)\n output = model_fit.forecast()\n yhat = output[0]\n pred.append(yhat)\n obs = test[t]\n history.append(obs)\n print('pred=%f, exp=%f' % (yhat, obs))\nerror = mean_squared_error(test, pred)\nerror2 = rmsle(pred,test)\n\nprint('Mean Squared Error: %.3f' % error)\nprint('RMSLE: %.3f' % error)\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"A line plot is created showing the expected values (blue) compared to the rolling forecast predictions (red). \nWe can see the values show some trend and are in the correct scale"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"pyplot.plot(test)\npyplot.plot(pred, color='red')\npyplot.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Time series prediction using Prophet in Python
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# Suppress warnings \nimport warnings\nwarnings.filterwarnings(\"ignore\", category=DeprecationWarning)\nwarnings.filterwarnings(\"ignore\", category=UserWarning)\nwarnings.filterwarnings(\"ignore\", category=FutureWarning)\nfrom IPython.display import HTML\n\nHTML('')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Delivered by Sean Taylor (Facebook) at the 2018 New York R Conference at Work-Bench on April April 20 and 21"},{"metadata":{},"cell_type":"markdown","source":"
Advantages of using Prophet
\n- Accommodates seasonality with multiple periods\n- Prophet is resilient to missing values\n- Best way to handle outliers in Prophet is to remove them\n- Fitting of the model is fast\n- Intuitive hyper parameters which are easy to tune"},{"metadata":{},"cell_type":"markdown","source":"Define Prophet dataset"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"X= no2_weather[['start_date','no2_emission_sum','temperature_2m_above_ground_mean','specific_humidity_2m_above_ground_mean','relative_humidity_2m_above_ground_mean',\n 'u_component_of_wind_10m_above_ground_mean','v_component_of_wind_10m_above_ground_mean','total_precipitation_surface_mean']]\ny=no2_weather['no2_emission_sum']","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Creating the data set for Prophet"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train_dataset= pd.DataFrame()\ntrain_dataset['ds'] = pd.to_datetime(no2_weather[\"start_date\"])\ntrain_dataset['y']=y\ntrain_dataset.head(2)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"prophet_basic = Prophet()\nprophet_basic.fit(train_dataset)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Predicting the values for the future
\n\nFor predicting the values using Prophet, we need to create a dataframe with ds(datetime stamp) containing the dates for which we want to make the predictions.\nWe use make_future_dataframe() to which we specify the number of days to extend into the future. By default it includes dates from the history"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"future= prophet_basic.make_future_dataframe(periods=30)\nfuture.tail()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Our prediction contains historical dates with 30 days."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"forecast=prophet_basic.predict(future)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
\nWe can plot the trend and seasonality, components of the forecast."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"fig1 = prophet_basic.plot_components(forecast)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The trend shows that the days as Monday, Weednesday have high values and the last month of year have high values. "},{"metadata":{},"cell_type":"markdown","source":"
Forecast quality evaluation
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\n\nLet's evaluate the quality of the algorithm by calculating the error metrics for the last 30 days that we predicted. For this, we will need the observations yi and the corresponding predicted values y^i .\n\nLet's look into the object forecast that the library created for us:"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"from fbprophet.diagnostics import cross_validation, performance_metrics\ndf_cv = cross_validation(prophet_basic, horizon='30 days')\ndf_p = performance_metrics(df_cv)\ndf_p.head(30)","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from fbprophet.plot import plot_cross_validation_metric\nfig_mape = plot_cross_validation_metric(df_cv, metric='mape')\nfig_rmse = plot_cross_validation_metric(df_cv, metric='rmse')\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It shows that 10 days forecast results in around 10% error."},{"metadata":{},"cell_type":"markdown","source":"#### Adding ChangePoints to Prophet\nChangepoints are the datetime points where the time series have abrupt changes in the trajectory."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"from fbprophet.plot import add_changepoints_to_plot\nfig = prophet_basic.plot(forecast)\na = add_changepoints_to_plot(fig.gca(), prophet_basic, forecast)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can view the dates where the chagepoints occurred"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"prophet_basic.changepoints[:10]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can change the inferred changepoint range by setting the changepoint_range"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"pro_change= Prophet(changepoint_range=0.9)\nforecast = pro_change.fit(train_dataset).predict(future)\nfig= pro_change.plot(forecast);\na = add_changepoints_to_plot(fig.gca(), pro_change, forecast)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Decreasing the changepoint_prior_scale to 0.001 to make the trend less flexible"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"pro_change= Prophet(n_changepoints=20, yearly_seasonality=True, changepoint_prior_scale=0.08)\nforecast = pro_change.fit(train_dataset).predict(future)\nfig= pro_change.plot(forecast);\na = add_changepoints_to_plot(fig.gca(), pro_change, forecast)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Incorporating the effects of weather condition
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\n\n\n> Now we add as extra regressors in the fbprophet model"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train_dataset['temperature_2m_above_ground_mean'] = X['temperature_2m_above_ground_mean']\ntrain_dataset['specific_humidity_2m_above_ground_mean'] = X['specific_humidity_2m_above_ground_mean']\ntrain_dataset['relative_humidity_2m_above_ground_mean'] = X['relative_humidity_2m_above_ground_mean']\ntrain_dataset['u_component_of_wind_10m_above_ground_mean'] = X['u_component_of_wind_10m_above_ground_mean']\ntrain_dataset['v_component_of_wind_10m_above_ground_mean'] = X['v_component_of_wind_10m_above_ground_mean']\ntrain_dataset['total_precipitation_surface_mean'] = X['total_precipitation_surface_mean']\n\ntrain_X= train_dataset[:200]\ntest_X= train_dataset[200:]\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"#Additional Regressor\npro_regressor= Prophet()\npro_regressor.add_regressor('temperature_2m_above_ground_mean')\npro_regressor.add_regressor('specific_humidity_2m_above_ground_mean')\npro_regressor.add_regressor('relative_humidity_2m_above_ground_mean')\npro_regressor.add_regressor('u_component_of_wind_10m_above_ground_mean')\npro_regressor.add_regressor('v_component_of_wind_10m_above_ground_mean')\npro_regressor.add_regressor('total_precipitation_surface_mean')\n#Fitting the data\npro_regressor.fit(train_X)\nfuture_data = pro_regressor.make_future_dataframe(periods=30) # 30 days\n#forecast the data for Test data\nforecast_data = pro_regressor.predict(test_X)\npro_regressor.plot(forecast_data);","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Predicted data is the blue shaded region at the end."},{"metadata":{},"cell_type":"markdown","source":"
Forecast quality evaluation
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\n\nLet's evaluate the quality of the algorithm by calculating the error metrics for the last 30 days that we predicted. For this, we will need the observations yi and the corresponding predicted values y^i .\n"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"df_cv_reg = cross_validation(pro_regressor, horizon='30 days')\ndf_p_reg = performance_metrics(df_cv_reg)\ndf_p_reg.head(30)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The RMSE for 30 days its 0.984912."},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"fig_mape_reg = plot_cross_validation_metric(df_cv_reg, metric='mape')\nfig_mape_reg = plot_cross_validation_metric(df_cv_reg, metric='rmse')\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It shows that 10 days forecast results in around ~7% error."},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from fbprophet.plot import plot_plotly, add_changepoints_to_plot\nimport plotly.offline as py\n\nfig_d_reg = plot_plotly(pro_regressor, forecast_data)\n\npy.iplot(fig_d_reg) \n\nfig_d_reg = pro_regressor.plot(forecast_data,xlabel='Date',ylabel='N02 values')\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"- ds — forecast date\n- yhat — forecast value for the given date\n- yhat_lower — lower forecast boundary for the given date\n- yhat_uppet — upper forecast boundary for the given date\nCalling plot function for Prophet model displays how the model was trained according to training data (black points — training data, blue line — forecast value, light blue area — forecast boundaries):"},{"metadata":{},"cell_type":"markdown","source":"## We can see that the Best Approach is Prophet."},{"metadata":{},"cell_type":"markdown","source":"
Save The Model
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\nThe model should be re-trained when new data becomes available. There is no point to re-train model, if data is not changed. Save model instead and use it again, when user wants to call predict function. Use pickle functionality for that:\n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"import pickle\nwith open('forecast_model_No2.pckl', 'wb') as fout:\n pickle.dump(pro_regressor, fout)\nwith open('forecast_model_No2.pckl', 'rb') as fin:\n m2 = pickle.load(fin)\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Prediction of NO2 density for each primary_fuel throughout the year
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)"},{"metadata":{},"cell_type":"markdown","source":"
Regional NO2 Density
"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from datetime import datetime\nfiles=[]\nfor dirname, _, filenames in os.walk('/kaggle/input/ds4g-environmental-insights-explorer/eie_data/s5p_no2'):\n for filename in filenames:\n files.append(os.path.join(dirname, filename))\n\n# read all the absorbing aerosol index data into one list of arrays\nno2_first_day=[]\nno2_first_key=[]\nno2_arr=[]\nband=0 # NO2_column_number_density\nfor i in range(0,len(files)):\n no2_first_day.append(datetime.strptime(files[i][76:91], '%Y%m%dT%H%M%S').date())\n no2_first_key.append(datetime.strptime(files[i][76:91], '%Y%m%dT%H%M%S').toordinal()+1) # correction of + 1 day in order to sync on climate data\n no2_arr.append(rio.open(files[i]).read(band+1))\n\n","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from pandas.plotting import register_matplotlib_converters\nregister_matplotlib_converters()\n\na=[]\na_pos=[]\nfor i in range(0,len(no2_arr)): \n a.append(np.nanmean(no2_arr[i]))\n a_pos.append(np.nanmean(np.clip(no2_arr[i],0,10000)))\n \nno2_rgn=pd.DataFrame({ 'start_date': no2_first_day,'no2_rgn' : a_pos, 'key_date' : no2_first_key })\nno2_rgn=no2_rgn.sort_values('start_date')\nno2_rgn=no2_rgn.reset_index()","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"# read only the NO2 index arrays with a nan-percentage <5% into one list of arrays for calculation of local NO2 data\nfiles=[]\nfor dirname, _, filenames in os.walk('/kaggle/input/ds4g-environmental-insights-explorer/eie_data/s5p_no2'):\n for filename in filenames:\n files.append(os.path.join(dirname, filename))\n\nno2_first_day=[]\nno2_first_key=[]\nno2_arr=[]\nband=0 # NO2_column_number_density\nfor i in range(0,len(files)):\n a=rio.open(files[i]).read(band+1)\n if pd.isnull(a).sum().sum() < 3515:\n no2_first_day.append(datetime.strptime(files[i][76:91], '%Y%m%dT%H%M%S').date())\n no2_first_key.append(datetime.strptime(files[i][76:91], '%Y%m%dT%H%M%S').toordinal()+1) # correction of + 1 day in order to sync on climate data\n no2_arr.append(np.clip(a,0,10000)) # clip negative values to zero","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Forecast quality evaluation for Power Plant over the year
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\nLet's evaluate the quality of the algorithm by calculating the error metrics for the last 30 days that we predicted. For this, we will need the observations yi and the corresponding predicted values y^i \n"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"df_cv_reg = cross_validation(pro_regressor, horizon='30 days')\ndf_p_reg = performance_metrics(df_cv_reg)\ndf_p_reg.head(30)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The RMSE for 30 days its 0.1.778483e-11."},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"fig_mape_reg = plot_cross_validation_metric(df_cv_reg, metric='mape')\nfig_mape_reg = plot_cross_validation_metric(df_cv_reg, metric='rmse')","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from fbprophet.plot import plot_plotly, add_changepoints_to_plot\nimport plotly.offline as py\n\nfig_d_reg = plot_plotly(pro_regressor, forecast_data)\n\npy.iplot(fig_d_reg) \n\nfig_d_reg = pro_regressor.plot(forecast_data,xlabel='Date',ylabel='no2 Density by Location of Plant values')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Save The Model for primary_fuel"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"import pickle\nwith open('forecast_model_No2Density.pckl', 'wb') as fout:\n pickle.dump(pro_regressor, fout)\nwith open('forecast_model_No2Density.pckl', 'rb') as fin:\n m2 = pickle.load(fin)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Another Example - Verification A.E.S.Corp "},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"power_plants_df[['name','primary_fuel','plant']][power_plants_df['primary_fuel']=='Coal']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"train_dataset['A.E.S._Corp._mean'] = X['A.E.S._Corp._mean']\n# Additional Regressor\npro_regressor= Prophet()\npro_regressor.add_regressor('A.E.S._Corp._mean')\ntrain_X= train_dataset[:200]\ntest_X= train_dataset[200:]\n\n#Fitting the data\npro_regressor.fit(train_X)\nfuture_data = pro_regressor.make_future_dataframe(periods=30) # 30 days\n#forecast the data for Test data\nforecast_data = pro_regressor.predict(test_X)\npro_regressor.plot(forecast_data);","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Forecast quality evaluation for region\nLet's evaluate the quality of the algorithm by calculating the error metrics for the last 30 days that we predicted. For this, we will need the observations yi and the corresponding predicted values y^i \n"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"df_cv_reg = cross_validation(pro_regressor, horizon='30 days')\ndf_p_reg = performance_metrics(df_cv_reg)\ndf_p_reg.head(30)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"fig_mape_reg = plot_cross_validation_metric(df_cv_reg, metric='mape')\nfig_mape_reg = plot_cross_validation_metric(df_cv_reg, metric='rmse')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"
Outlier Analysis of Power Plant - Coal over the year
\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n \n\n\nIdentifying outliers for A.E.S._Corp "},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"month_p_fuel = X[['start_date', 'no2_density_locationofplant', 'A.E.S._Corp._mean']].copy()\nmonth_p_fuel['date'] = pd.to_datetime(X[\"start_date\"])\nmonth_p_fuel['date'] = month_p_fuel['date'].dt.month\nmonth_p_fuel = month_p_fuel.groupby(['date','A.E.S._Corp._mean']).sum()\nmonth_p_fuel\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"month_p_fuel_agg = month_p_fuel.groupby(['date', 'A.E.S._Corp._mean']).agg(['sum'])\nmonth_p_fuel_agg = month_p_fuel_agg.reset_index()\nlevel_0 = month_p_fuel_agg.columns.droplevel(0)\nlevel_1 = month_p_fuel_agg.columns.droplevel(1)\nlevel_0 = ['' if x == '' else '-' + x for x in level_0]\n\nmonth_p_fuel_agg.columns = level_1 + level_0\nmonth_p_fuel_agg.rename_axis(None, axis=1)\n# month_p_fuel_agg.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"from plotly import tools, subplots\nimport plotly.offline as py\npy.init_notebook_mode(connected=True)\nimport plotly.graph_objs as go\nimport plotly.express as px\n\nfig_total = px.line(month_p_fuel_agg, x='date', y='no2_density_locationofplant-sum', color='A.E.S._Corp._mean', render_mode='svg')\nfig_total.update_layout(title='Total NO2 aspect in A.E.S._Corp._mean - 12 months')\nfig_total.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The sum, facetted for A.E.S._Corp._mean aspect, shows some aberrant values, for example in month 4, April. August have the least N02 value. "},{"metadata":{},"cell_type":"markdown","source":"### Outlier Analysis of Power Plant - Coal day of week"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"dayweek_p_fuel = X[['start_date', 'no2_density_locationofplant', 'A.E.S._Corp._mean']].copy()\ndayweek_p_fuel['dateofweek'] = pd.to_datetime(X[\"start_date\"])\ndayweek_p_fuel['dateofweek'] = dayweek_p_fuel['dateofweek'].dt.dayofweek\ndayweek_p_fuel = dayweek_p_fuel.groupby(['dateofweek','A.E.S._Corp._mean']).sum()\n\ndayofweek_p_fuel_agg = dayweek_p_fuel.groupby(['dateofweek', 'A.E.S._Corp._mean']).agg(['sum'])\ndayofweek_p_fuel_agg = dayofweek_p_fuel_agg.reset_index()\nlevel_0 = dayofweek_p_fuel_agg.columns.droplevel(0)\nlevel_1 = dayofweek_p_fuel_agg.columns.droplevel(1)\nlevel_0 = ['' if x == '' else '-' + x for x in level_0]\n\ndayofweek_p_fuel_agg.columns = level_1 + level_0\ndayofweek_p_fuel_agg.rename_axis(None, axis=1)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Outlier Analysys - day of the week with Monday=0, Sunday=6"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"from plotly import tools, subplots\nimport plotly.offline as py\npy.init_notebook_mode(connected=True)\nimport plotly.graph_objs as go\nimport plotly.express as px\n\nfig_total = px.line(dayofweek_p_fuel_agg, x='dateofweek', y='no2_density_locationofplant-sum', color='A.E.S._Corp._mean', render_mode='svg')\nfig_total.update_layout(title='Total NO2 aspect in A.E.S._Corp._mean - day of week')\nfig_total.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The sum, faceted for A.E.S.Corp._ mean aspect shows some aberrant values. For example, in general, Thursday (3) has the bigger aberrant value. Tuesday (1) has the least N02 value."},{"metadata":{},"cell_type":"markdown","source":"## Our modeling investigate regions, primary_fuel model and we can decompose emission factors between plants,over time and identifying anomaly events."},{"metadata":{},"cell_type":"markdown","source":" \nInspired by: [Exploratory Data Analysis and Factor Model](https://www.kaggle.com/ragnar123/exploratory-data-analysis-and-factor-model-idea),\n[Modelling of emissions of power plants](https://www.kaggle.com/tiurii/ds4g-modelling-of-emissions-of-power-plants)\n\nsource: [Survey](https://www.tandfonline.com/doi/full/10.1080/10962247.2019.1577314?scroll=top&needAccess=true),[Arima python](https://machinelearningmastery.com/arima-for-time-series-forecasting-with-python/), [anomaly detection](https://github.com/Vicam/Unsupervised_Anomaly_Detection/blob/master/Anomaly%20detection%2C%20different%20methods%20on%20a%20simple%20example.ipynb), [Intro LSTM](https://colah.github.io/posts/2015-08-Understanding-LSTMs/),[LSTM](https://www.kaggle.com/faressayah/stock-market-analysis-prediction-using-lstm), [prophet facebook](https://towardsdatascience.com/time-series-prediction-using-prophet-in-python-35d65f626236), [forecast in python](https://towardsdatascience.com/forecasting-in-python-with-facebook-prophet-29810eb57e66)"},{"metadata":{},"cell_type":"markdown","source":" \n\n#
About the data
\n\n\n[Back to Table of Contents](#top)\n\n[General Findinds](#theend)\n\n\n[Global Power Plant database ](https://developers.google.com/earth-engine/datasets/catalog/WRI_GPPD_power_plants) by WRI\n> Description\nThe Global Power Plant Database is a comprehensive, open source database of power plants around the world. It centralizes power plant data to make it easier to navigate, compare and draw insights for one’s own analysis. The database covers approximately 30,000 power plants from 164 countries and includes thermal plants (e.g. coal, gas, oil, nuclear, biomass, waste, geothermal) and renewables (e.g. hydro, wind, solar). Each power plant is geolocated and entries contain information on plant capacity, generation, ownership, and fuel type. It will be continuously updated as data becomes available.\n\n[Sentinel 5P OFFL NO2](https://developers.google.com/earth-engine/datasets/catalog/COPERNICUS_S5P_OFFL_L3_NO2) by [EU/ESA/Copernicus](https://sentinel.esa.int/web/sentinel/user-guides/sentinel-5p-tropomi/document-library)\n> Sentinel-5 Precursor\nSentinel-5 Precursor is a satellite launched on 13 October 2017 by the European Space Agency to monitor air pollution. The onboard sensor is frequently referred to as Tropomi (TROPOspheric Monitoring Instrument). The OFFL/NO2 is a dataset that provides offline high-resolution imagery of **NO2 concentration**.\n\n[Global Forecast System 384-Hour Predicted Atmosphere Data](https://developers.google.com/earth-engine/datasets/catalog/NOAA_GFS0P25) by NOAA/NCEP/EMC\n> The Global Forecast System (GFS) is a weather forecast model produced by the National Centers for Environmental Prediction (NCEP). The GFS dataset consists of selected model outputs (described below) as gridded forecast variables. The 384-hour forecasts, with 3-hour forecast interval, are made at 6-hour temporal resolution (i.e. updated four times daily). Use the 'creation_time' and 'forecast_time' properties to select data of interest.\n\n[Global Land Data Assimilation System](https://developers.google.com/earth-engine/datasets/catalog/NASA_GLDAS_V021_NOAH_G025_T3H) by NASA\n> Global Land Data Assimilation System (GLDAS) ingests satellite and ground-based observational data products. Using advanced land surface modeling and data assimilation techniques, it generates optimal fields of land surface states and fluxes. his dataset provided by NASA ingest satellite.\n\nParticipants may also consider using other public datasets related to trade commodities for fuel types, total fuel consumed, and/or data from the [US Energy Information Agency (EIA)](https://www.eia.gov/state/data.php?sid=RQ#CarbonDioxideEmissions)."},{"metadata":{},"cell_type":"markdown","source":"\n\n#
Don't hesitate to give your suggestions in the comment section.
\n\n\n#
Remember the upvote button is next to the fork button, and it's free too! ;)