pre_process.py

# -*- coding: utf-8 -*-
"""Pre-process.ipynb

Automatically generated by Colaboratory.

Original file is located at
    https://colab.research.google.com/drive/1o1ZtxG7s8UWQ1utXv__IxlBVmmS69rFc

# Space Weather

Notebook to reproduce/improve forecasts obtained by paper by Siciliano et at. *Forecasting SYM‐H Index: A Comparison Between Long Short‐Term Memory and Convolutional Neural Networks* ([paper](https://agupubs.onlinelibrary.wiley.com/doi/epdf/10.1029/2020SW002589)).
"""

import pandas as pd
import datetime
import seaborn as sns
import matplotlib.pyplot as plt
import numpy as np
import json
import os
sns.set_style('dark')
plt.rcParams['font.size'] = 16

#from google.colab import drive
#drive.mount('/content/gdrive')

import os
dir_path = os.getcwd()

# Results directory
results_dir = 'ACE'

drive_dir = os.path.join(dir_path, results_dir)
if not os.path.exists(drive_dir):
    os.makedirs(drive_dir)

# {'n_layers': 1, 'n_unit': 68, 'lr': 0.0026336545657049317, 'lb': 25}

"""# Download Database

Database obtabined from: https://omniweb.gsfc.nasa.gov
"""

# 5 minute sampling
#!wget -q https://www.dropbox.com/s/hcblxnxzg1qrjsg/gicDB_5min.csv?dl=1 -O gicDB.csv

"""

# Data pre-processing"""

data = pd.read_csv("gicDB.csv")
data.head()

"""# Just for to check data, remove overflows"""

data = data.replace(99999.9,None)
data = data.replace(9999.99,None)
data = data.replace(9999.9,None)
data = data.replace(999.99,None)

# Get month and day to build proper data field
data.loc[:,'Date'] = pd.to_datetime((data.Year*1000+data.Day).apply(str), format = "%Y%j")
data['month'] = data['Date'].dt.month
data['day'] = data['Date'].dt.day

# Date field
data["time"] = pd.to_datetime(data[['Year','month','day','Hour','Minute']])
data["date"] = pd.to_datetime((data.Year*1000+data.Day).apply(str), format = "%Y%j")
data['month'] = data['date'].dt.month
data['day'] = data['date'].dt.day

data["time"] = pd.to_datetime(data[['Year','month','day','Hour','Minute']])
data = data.sort_values(by='time')

data.interpolate(method ='linear', limit_direction ='forward')

data["SYM/H nT"] = data["SYM/H nT"].astype(np.float64)
data["BX nT (GSE GSM)"] = data["BX nT (GSE GSM)"].astype(np.float64)
data["BY nT (GSM)"] = data["BY nT (GSM)"].astype(np.float64)
data["BZ nT (GSM)"] = data["BZ nT (GSM)"].astype(np.float64)

# Add new variables to dataframe
data["B**2"] = data["BX nT (GSE GSM)"]**2+data["BY nT (GSM)"]**2+data["BZ nT (GSM)"]**2
data["BY**2"] = data["BY nT (GSM)"]**2

data.head()

data.columns

feature_keys = ['Field magnitude average nT', 'BX nT (GSE GSM)', 'BY nT (GSE)', 'BZ nT (GSE)', 'BY nT (GSM)', 'BZ nT (GSM)', 'Speed km/s', 'Vx Velocity km/s', 'Vy Velocity km/s', 'Vz Velocity km/s', 'Proton Density n/cc', 'SYM/H nT', 'B**2', 'BY**2']
titles = feature_keys

colors = [
    "blue",
    "orange",
    "green",
    "red",
    "purple",
    "brown",
    "pink",
    "gray",
    "olive",
    "cyan",
]

date_time_key = "time"

"""# Visualize"""

# Create a few plots with the values for all data
def show_raw_visualization(data):
    time_data = data[date_time_key]
    fig, axes = plt.subplots(
        nrows=7, ncols=2, figsize=(15, 20), dpi=80, facecolor="w", edgecolor="k"
    )
    for i in range(len(feature_keys)):
        key = feature_keys[i]
        c = colors[i % (len(colors))]
        t_data = data[key]
        t_data.index = time_data
        t_data.head()
        ax = t_data.plot(
            ax=axes[i // 2, i % 2],
            color=c,
            title="{} - {}".format(titles[i], key),
            rot=25,
        )
        ax.legend([titles[i]])
    plt.tight_layout()
    plt.savefig('All_Variables.png')
    plt.close(fig)


show_raw_visualization(data)

df = data
non_floats = []
for col in df:
    if df[col].dtypes != "float64":
        non_floats.append(col)
df = df.drop(columns=non_floats)

# Plot correlations between variables
def show_heatmap(df):
    plt.figure(figsize=(15,10))
    plt.matshow(df.corr(), fignum=1)
    plt.xticks(range(df.shape[1]), df.columns, fontsize=14, rotation=90)
    plt.gca().xaxis.tick_bottom()
    plt.yticks(range(df.shape[1]), df.columns, fontsize=14)

    cb = plt.colorbar()
    cb.ax.tick_params(labelsize=14)
    plt.title("Feature Correlation Heatmap", fontsize=14)
    plt.savefig('Heatmap.png')


show_heatmap(df)

"""# Create Pandas for Train, Validation and Test"""

# Data Normalization to make training easier
def normalize(data):
    data_mean = data.mean(axis=0)
    data_std = data.std(axis=0)
    return data_mean, data_std, (data - data_mean) / data_std

# Check when a new Storm starts
interface = [0]

def checkBreaks(data) :

  storms = data['Dataset'].unique()
  for st in storms :

    features = data.loc[data['Dataset'] == st]
    interface.append(interface[-1]+features.shape[0])

  return interface

"""Create train and validation data frames"""

selected_features = ["SYM/H nT", "BZ nT (GSM)", "B**2", "BY**2"]

data_types = ['Train', 'Validation', 'Test']

for dtype in data_types :
  features = data.loc[data['Type'] == dtype][selected_features]
  features.index = data.loc[data['Type'] == dtype][date_time_key]
  mean, std, norm_features = normalize(features.values)
  breaks = checkBreaks(data.loc[data['Type'] == dtype])
  norm_data = pd.DataFrame(norm_features,columns=selected_features)

  # Normalized data
  file_name =os.path.join(drive_dir, dtype+'_norm.gzip')
  norm_data.to_parquet(file_name,compression='gzip')

  # Processed data
  file_name =os.path.join(drive_dir, dtype+'.gzip')
  features.to_parquet(file_name,compression='gzip')

  file_name =os.path.join(drive_dir, dtype+'.json')
  summDict = {"mean":mean.tolist(), "std" : std.tolist(), "breaks" : breaks, "features":selected_features, "time_key" : date_time_key}
  with open(file_name, 'w') as outfile:
    json.dump(summDict, outfile)


# Processed data
file_name =os.path.join(drive_dir, 'all_data.gzip')
data.to_parquet(file_name,compression='gzip')