<jupyter_start><jupyter_text>Mount Google Drive in colab<jupyter_code>from google.colab import drive
drive.mount("/content/drive")<jupyter_output>Mounted at /content/drive<jupyter_text>🧾Description: We advance the modeling capability of electron particle precipitation from the magnetosphere to the ionosphere through a new database and use of deep learning (DL) tools to gain utility from those data. We have compiled, curated, analyzed, and made available a new and more capable database of particle precipitation data that includes 51 satellite years of Defense Meteorological Satellite Program (DMSP) observations temporally aligned with solar wind and geomagnetic activity data. The new total electron energy flux particle precipitation nowcast model, a neural network called PrecipNet, takes advantage of increased expressive power afforded by ML approaches to appropriately utilize diverse information from the solar wind and geomagnetic activity and, importantly, their time histories. With a more capable representation of the organizing parameters and the target electron energy flux observations, PrecipNet achieves a >50% reduction in errors from a current state‐of‐the‐art model oval variation, assessment, tracking, intensity, and online nowcasting (OVATION Prime), better captures the dynamic changes of the auroral flux, and provides evidence that it can capably reconstruct mesoscale phenomena. We create and apply a new framework for space weather model evaluation that culminates previous guidance from across the solar‐terrestrial research community.<jupyter_code>import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import os 
#
import warnings 
warnings.filterwarnings('ignore')
warnings.filterwarnings('ignore', category=DeprecationWarning)
path_dir = "/content/drive/MyDrive/TMLC/DMSP Particle Precipitate Flux Prediction"
file_name = "AI_Ready_DMSP_Data.csv"
file_path = os.path.join(path_dir,file_name)
print(file_path)
df = pd.read_csv(file_path)

df.head()
print(df.columns.tolist())<jupyter_output>['Datetimes', 'SC_AACGM_LAT', 'SC_AACGM_LTIME', 'ELE_TOTAL_ENERGY_FLUX', 'ELE_TOTAL_ENERGY_FLUX_STD', 'ELE_AVG_ENERGY', 'ELE_AVG_ENERGY_STD', 'ID_SC', 'sin_ut', 'cos_ut', 'sin_doy', 'cos_doy', 'sin_SC_AACGM_LTIME', 'cos_SC_AACGM_LTIME', 'F107', 'Bz', 'By', 'Bx', 'AE', 'AL', 'AU', 'SymH', 'PC', 'vsw', 'vx', 'psw', 'borovsky', 'newell', 'F107_6hr', 'Bz_6hr', 'By_6hr', 'Bx_6hr', 'AE_6hr', 'AL_6hr', 'AU_6hr', 'SymH_6hr', 'PC_6hr', 'vsw_6hr', 'vx_6hr', 'psw_6hr', 'borovsky_6hr', 'newell_6hr', 'F107_5hr', 'Bz_5hr', 'By_5hr', 'Bx_5hr', 'AE_5hr', 'AL_5hr', 'AU_5hr', 'SymH_5hr', 'PC_5hr', 'vsw_5hr', 'vx_5hr', 'psw_5hr', 'borovsky_5hr', 'newell_5hr', 'F107_3hr', 'Bz_3hr', 'By_3hr', 'Bx_3hr', 'AE_3hr', 'AL_3hr', 'AU_3hr', 'SymH_3hr', 'PC_3hr', 'vsw_3hr', 'vx_3hr', 'psw_3hr', 'borovsky_3hr', 'newell_3hr', 'F107_1hr', 'Bz_1hr', 'By_1hr', 'Bx_1hr', 'AE_1hr', 'AL_1hr', 'AU_1hr', 'SymH_1hr', 'PC_1hr', 'vsw_1hr', 'vx_1hr', 'psw_1hr', 'borovsky_1hr', 'newell_1hr', 'F107_45min', 'Bz_45min', 'By_45min', '[...]<jupyter_text>Summary Statistics<jupyter_code>df.shape
df.info()
df['Datetimes'].dtypes<jupyter_output><empty_output><jupyter_text>Check Null Values<jupyter_code>df.isnull().sum()[df.isnull().sum() > 0]
df['Datetimes'] = pd.to_datetime(df['Datetimes'])
df.info()<jupyter_output><class 'pandas.core.frame.DataFrame'>
RangeIndex: 1945887 entries, 0 to 1945886
Columns: 154 entries, Datetimes to newell_5min
dtypes: datetime64[ns](1), float64(153)
memory usage: 2.2 GB<jupyter_text>Target - ELE_TOTAL_ENERGY_FLUX<jupyter_code>df['ELE_TOTAL_ENERGY_FLUX'].describe()
sns.displot(data=df,x='ELE_TOTAL_ENERGY_FLUX',kind="kde")
df['Year'] = df['Datetimes'].dt.year
df['Month'] = df['Datetimes'].dt.month
df['Day'] = df['Datetimes'].dt.day
sns.lineplot(x='Year',y='ELE_TOTAL_ENERGY_FLUX',data=df)
sns.lineplot(x='Month',y='ELE_TOTAL_ENERGY_FLUX',data=df)
sns.lineplot(x='Day',y='ELE_TOTAL_ENERGY_FLUX',data=df)
new_df = df.drop('Datetimes',axis=1)
new_df.shape<jupyter_output><empty_output><jupyter_text>https://github.com/rmcgranaghan/precipNet/blob/master/Final__Data_Read_And_Prepare.ipynb<jupyter_code># Choose columns to drop
cols_to_drop = [c for c in new_df.columns if ('1min' in c) | ('3min' in c) | ('4min' in c) | ('_5min'in c) | ('15min' in c)| ('5hr' in c) | ('F107_' in c) | ('vx' in c) | ('AE' in c) | ('newell' in c) | ('borovsky' in c) | ('STD' in c) | ('AVG' in c) | ('SC_AACGM_LTIME'==c)]
new_df = new_df.drop(columns=cols_to_drop)
new_df.shape
print('number of input features total = {}'.format(new_df.shape[1]))<jupyter_output>number of input features total = 76<jupyter_text>Reference :https://www.researchgate.net/publication/350946002_Toward_a_Next_Generation_Particle_Precipitation_Model_Mesoscale_Prediction_Through_Machine_Learning_a_Case_Study_and_Framework_for_Progress Two final points about data preparation. * First, we remove observations that represent unrealistically large fluxes, deemed to be those exceeding the 99.995th percentile (equal to 7.3673E+13 [eV/cm2/s/ster] or ∼118erg/cm2/sr/s), which was empirically determined to separate values that were nonphysical. * Second, we combine the northern and southern hemispheres. The sun-synchronous orbits of the DMSP spacecraft restricts their MLT coverage, necessitating the combination of northern and southern hemisphere observations to com-plete the local time coverage<jupyter_code># Define the threshold above which to remove energy flux values
eflux_threshold_quantile = 99.995
eflux_threshold = .99995

# Define the conversion factor from eV to erg
erg_eV_factor = (1.60218e-12)<jupyter_output><empty_output><jupyter_text>Remove Outliers<jupyter_code># Remove outliers

print(new_df.shape)
threshold_drop_val = np.percentile(new_df['ELE_TOTAL_ENERGY_FLUX'],eflux_threshold_quantile)
print('... dropping rows where total electron energy flux is greater than {:.4E} [eV/cm2/s/ster]'.format(threshold_drop_val))
print('--------> number dropped = {}'.format(len( np.argwhere(new_df['ELE_TOTAL_ENERGY_FLUX'].values > threshold_drop_val) )))
#
new_df = new_df[ (new_df['ELE_TOTAL_ENERGY_FLUX'].values < threshold_drop_val) ]
print(new_df.shape)<jupyter_output>(1945887, 76)
... dropping rows where total electron energy flux is greater than 7.3673E+13 [eV/cm2/s/ster]
--------> number dropped = 98
(1945789, 76)<jupyter_text>Separate features and target<jupyter_code>X = new_df.drop('ELE_TOTAL_ENERGY_FLUX',axis=1)
Y = new_df['ELE_TOTAL_ENERGY_FLUX']<jupyter_output><empty_output><jupyter_text>Feature-engine: Drop Correlated Features<jupyter_code>!pip install --quiet feature-engine
import feature_engine
from feature_engine.selection import DropCorrelatedFeatures
tr = DropCorrelatedFeatures(variables=None,method='pearson',threshold=0.8)
drop_features= tr.fit_transform(X)
drop_features.shape
cols_to_drop = drop_features.columns.tolist()
cols_to_drop<jupyter_output><empty_output><jupyter_text>Final Features<jupyter_code>final_features = X.drop(cols_to_drop,axis=1)<jupyter_output><empty_output><jupyter_text>Train Test Split<jupyter_code>from sklearn.model_selection import train_test_split
X_train,X_test,y_train,y_test = train_test_split(final_features,Y)
X_train.shape<jupyter_output><empty_output><jupyter_text>Normalize the data<jupyter_code>from sklearn.preprocessing import StandardScaler
sc = StandardScaler()

X_train_scaled = sc.fit_transform(X_train)
X_test_scaled = sc.transform(X_test)
train_scaled = pd.DataFrame(X_train_scaled,columns=X_train.columns)
test_scaled = pd.DataFrame(X_test_scaled,columns=X_test.columns)
train_scaled.head()
test_scaled.head()<jupyter_output><empty_output><jupyter_text>Save the standarascaler trained data<jupyter_code>from joblib import dump
sc_filename = "standardscaler.bin"
sc_path = os.path.join(path_dir,sc_filename)
print(sc_path)
dump(sc,sc_path)<jupyter_output>/content/drive/MyDrive/TMLC/DMSP Particle Precipitate Flux Prediction/standardscaler.bin<jupyter_text>Import Keras packages<jupyter_code>import tensorflow as tf
from tensorflow import keras
from tensorflow.keras.models import Sequential
from tensorflow.keras.layers import Dense,Dropout
from tensorflow.keras.callbacks import ReduceLROnPlateau
callbacks = ReduceLROnPlateau(monitor='val_loss',
                              factor=0.1,
                              patience=5,
                              verbose=1,
                              mode="auto",
                              min_lr=0.001)<jupyter_output><empty_output><jupyter_text>Create Model<jupyter_code>def create_model():
  model = Sequential()
  model.add(Dense(256, input_dim=33, activation= 'relu' ))
  model.add(Dense(256, activation= 'relu'))
  model.add(Dropout(0.3)),
  model.add(Dense(256,  activation= 'relu' ))
  model.add(Dropout(0.3)),
  model.add(Dense(1,  activation= 'linear' ))
  #compile the model 
  model.compile(loss='mean_absolute_error', optimizer='adam', metrics=['mean_absolute_error'])
  return model
model = create_model()
model.summary()<jupyter_output>Model: "sequential_2"
_________________________________________________________________
 Layer (type)                Output Shape              Param #   
=================================================================
 dense_8 (Dense)             (None, 256)               8704      
                                                                 
 dense_9 (Dense)             (None, 256)               65792     
                                                                 
 dropout_4 (Dropout)         (None, 256)               0         
                                                                 
 dense_10 (Dense)            (None, 256)               65792     
                                                                 
 dropout_5 (Dropout)         (None, 256)               0         
                                                                 
 dense_11 (Dense)            (None, 1)                 257       
                                                      [...]<jupyter_text>Train the Model<jupyter_code>history = model.fit(X_train_scaled, y_train, validation_data = (X_test_scaled,y_test),epochs=15, batch_size=40)
 plt.plot(history.history['loss']) 
plt.plot(history.history['val_loss']) 
plt.title('Model loss') 
plt.ylabel('Loss') 
plt.xlabel('Epoch') 
plt.legend(['Train', 'Test'], loc='upper left') 
plt.show()
from sklearn.metrics import mean_squared_error
y_pred = model.predict(X_test_scaled, batch_size=1, verbose=1)<jupyter_output>266006/486448 [===============>..............] - ETA: 6:56
{'fraction_non_alphanumeric': 0.12527202195098874, 'fraction_numerical': 0.027911817579714258, 'mean_word_length': 4.69811320754717, 'pattern_counts': {'":': 0, '<': 63, '<?xml version=': 0, '>': 68, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 12, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '17', 'dir': 'jupyter-structured-clean-dedup', 'id': '48081', 'license': 'MIT', 'n_tokens_mistral': 3845, 'n_tokens_neox': 3508, 'n_words': 773, 'path': 'DMSP_Particle_Precipitate_Flux_Prediction.ipynb', 'repo_name': 'plaban1981/TMLC', 'size': '10569'}