<jupyter_start><jupyter_text>Sustainable Mobility: classification of electric scooter rides  Get Data (Colab only)<jupyter_code>!git clone https://github.com/DPaletti/mida_acv
!mv mida_acv/data .
!yes|rm -r mida_acv<jupyter_output>Cloning into 'mida_acv'...
remote: Enumerating objects: 249, done.[K
remote: Counting objects: 100% (54/54), done.[K
remote: Compressing objects: 100% (49/49), done.[K
remote: Total 249 (delta 6), reused 46 (delta 2), pack-reused 195[K
Receiving objects: 100% (249/249), 250.21 MiB | 26.30 MiB/s, done.
Resolving deltas: 100% (23/23), done.
Checking out files: 100% (162/162), done.<jupyter_text>Dependencies (Colab only)<jupyter_code>!pip install tsfresh
!pip install rdp
!pip install plotly
!pip install joblib
!pip install sktime
# After installing tsfresh runtime needs to be restarted
exit()<jupyter_output><empty_output><jupyter_text>Imports<jupyter_code>from typing import Tuple, List, Dict, Optional
from tempfile import mkdtemp
from pathlib import Path
import multiprocessing
import pandas as pd
import tsfresh as ts
import sklearn as sk
import scipy as sp
import scipy.fft as fft
import scipy.fftpack as fftpack
import scipy.integrate as integrate
import numpy as np
import plotly as plt
import rdp
import plotly.graph_objects as go
import joblib
import sklearn.pipeline
import sklearn.preprocessing
import sklearn.decomposition
import sklearn.linear_model
import scipy.signal
import tsfresh.feature_extraction
from itertools import repeat
from IPython.core.interactiveshell import InteractiveShell
InteractiveShell.ast_node_interactivity = 'all'
plt.io.renderers.default = 'iframe'<jupyter_output><empty_output><jupyter_text>Definitions  Read data<jupyter_code>def read_data(data_path: str):
  X = pd.DataFrame()
  y_weight=pd.Series(dtype=np.float64)
  y_passengers=pd.Series(dtype=np.int64)
  y_weight_deck=pd.Series(dtype=np.float64)
  y_passengers_deck=pd.Series(dtype=np.int64)
  y_weight_stem=pd.Series(dtype=np.float64)
  y_passengers_stem=pd.Series(dtype=np.int64)
  curr_X = pd.DataFrame()
  i: int = 0
  for placement in {"deck", "stem"}:
      for driver_number in {"single", "double"}:
          for ds in Path(data_path).joinpath(placement, driver_number).iterdir():
              curr_X = pd.read_csv(str(ds))
              curr_X = curr_X.assign(id=i)
              curr_X = curr_X.assign(full_id=placement + "-" + driver_number + "-" + curr_X["Driver"][0])
              curr_X = curr_X.assign(IsDeck=0 if placement == "deck" else 1)
              y_weight.at[i] = curr_X["Weight"][0]
              y_passengers.at[i] = 0 if driver_number=="single" else 1
              if placement == "deck":
                y_weight_deck.at[i] = curr_X["Weight"][0]
                y_passengers_deck.at[i] = 0 if driver_number=="single" else 1
              else:
                y_weight_stem.at[i] = curr_X["Weight"][0]
                y_passengers_stem.at[i] = 0 if driver_number=="single" else 1
                

              curr_X = curr_X.drop(
                  ["Unnamed: 0", "Driver", "Placement", "Weight"], 
                  axis=1
              )
              X = X.append(curr_X)
              i += 1
  to_weight_class =(lambda x: 0 if x<70 else (1 if (x>=70 and x<90) else (2 if (x>=90 and x<110) else 3)))
  return (X.fillna(0),
          X.fillna(0).groupby("id").filter(lambda group: group["full_id"][0].split("-")[0] == "deck").drop("IsDeck", axis=1), 
          X.fillna(0).groupby("id").filter(lambda group: group["full_id"][0].split("-")[0] == "stem").drop("IsDeck", axis=1), 
          y_weight.fillna(0).map(to_weight_class), 
          y_passengers.fillna(0),
          y_weight_deck.fillna(0).map(to_weight_class),
          y_passengers_deck.fillna(0),
          y_weight_stem.fillna(0).map(to_weight_class),
          y_passengers_stem.fillna(0))<jupyter_output><empty_output><jupyter_text>Align signals<jupyter_code>def align_signal(signal_1: np.array, signal_2: np.array, col_name:str):
  signal_1 = (signal_1 - np.mean(signal_1))/np.std(signal_1)
  signal_2 = (signal_2 - np.mean(signal_2))/np.std(signal_2)
  shift = np.argmax(np.convolve(signal_1[::-1],signal_2,mode='valid'))
  return np.append([0]*shift, signal_2[shift: ]), col_name

def align(X: pd.DataFrame):
  out_X = pd.DataFrame(columns=X.columns)
  columns_to_ignore = ["full_id", "IsDeck", "Timestamp", "id", "Latitude", "Longitude", "Speed", "Confidence"]
  for ts in [x for _, x in X.groupby("id").filter(lambda group: group["full_id"][0].split("-")[0]=="deck").groupby("id")]:
    full_id = ts["full_id"][0].split("-")
    for ts_to_align in [x for _, x in X.groupby("id").filter(lambda group: group["full_id"][0].split("-")[0]=="stem").groupby("id")]:
      full_id_to_align = ts_to_align["full_id"][0].split("-")
      if full_id[1:] == full_id_to_align[1:]:
        ts_signals = ts.drop(columns_to_ignore, axis=1)
        ts_to_align_signals = ts_to_align.drop(columns_to_ignore, axis=1)
        with multiprocessing.Pool(multiprocessing.cpu_count()) as p:
          aligned_signals = (p.starmap(align_signal, 
                            [(ts_item[1], ts_to_align_item[1], ts_to_align_item[0])
                            for ts_item, ts_to_align_item 
                            in zip(ts_signals.iteritems(), ts_to_align_signals.iteritems())]))
        temp_ts = pd.DataFrame(columns=ts_to_align.columns)
        temp_ts[columns_to_ignore] = ts_to_align[columns_to_ignore]
        
        for aligned_signal, col in aligned_signals:
          temp_ts[col] = aligned_signal
        out_X = out_X.append(ts)
        out_X = out_X.append(temp_ts)
  return out_X,<jupyter_output><empty_output><jupyter_text>Path simplification<jupyter_code>def get_path(X: pd.DataFrame) -> np.array:
    out = []
    for index, row in X.iterrows():
        out.append((row["Latitude"], row["Longitude"]))
    return np.array(out)

def simplify_path(X: pd.DataFrame,epsilon: float = 1e-6):
    if epsilon <= 0:
      return X
    out_df: pd.DataFrame = pd.DataFrame()
    for df in [x for _, x in X.groupby(["id"])]:
      print("Simplifying: " + df["full_id"][0])
      df = df[rdp.rdp(get_path(df), epsilon=epsilon, return_mask=True)]
      out_df = out_df.append(df)
    return out_df<jupyter_output><empty_output><jupyter_text>Windowing<jupyter_code>def window(X: pd.DataFrame,
           y_weight: pd.Series,
           rolling_direction: int=1,
           min_timeshift: int=0,
           max_timeshift: Optional[int] = None):

    X_rolled = ts.utilities.dataframe_functions.roll_time_series(
            X,
            column_id="id",
            column_sort="Timestamp",
            column_kind=None,
            min_timeshift=min_timeshift,
            max_timeshift=max_timeshift,
            rolling_direction=rolling_direction,
            n_jobs=multiprocessing.cpu_count(),
        )

    y_weight_out = np.empty([0])
    y_passengers_out = np.empty([0])

    for df in [x for _, x in X_rolled.groupby(["id"])]:
        ident = df["id"].values[0][0]
        y_weight_out = np.append(y_weight_out, y_weight[ident])
        y_passengers_out = np.append(y_passengers_out, y_passenger[ident])
    return X_rolled, y_weight_out, y_passengers_out<jupyter_output><empty_output><jupyter_text>Feature Extraction<jupyter_code>def extract_features(X: pd.DataFrame):
  return ts.extract_features(
          X,
          column_id="id",
          column_sort="Timestamp",
          n_jobs=multiprocessing.cpu_count(),
          default_fc_parameters=ts.feature_extraction.MinimalFCParameters()
      )<jupyter_output><empty_output><jupyter_text>Feature Extraction in the frequency domain<jupyter_code>def compute_fft(X: pd.DataFrame):
  X = X.drop(["Timestamp", "Latitude", "Longitude", "Confidence", "id", "full_id"], axis=1)
  out = pd.DataFrame(columns=[i + "_fft" for i in X.columns])
  for column in X.columns:
    out[column+"_fft"] =fft.fft(X[column].values)
  return out

def compute_fft_area_full(X: pd.DataFrame):
  X = X.drop(["Ax_fft", "Ay_fft", "Az_fft", "Jerk_x_fft"], axis=1)
  out = {}
  # TODO care for "id" column and compute windows
  for column in X.columns:
    out[column+"_area"] = integrate.trapz(np.abs(X[column].values))
  return out

def compute_fft_area_split(X: pd.DataFrame):
  X = X.drop(["Speed_fft", "Gx_fft", "Gy_fft", "Gz_fft", "Jerk_y_fft", "Jerk_z_fft"], axis=1)
  out = {}
  half_sample_frequency = 52
  series_length = len(X["Ax_fft"])
  split_frequencies = {"Ax_fft": 4,
                       "Ay_fft": 5,
                       "Az_fft": 10,
                       "Jerk_x_fft": 15}
  frequencies = fftpack.fftfreq(series_length, d=1/104)
  for column in X.columns:
    lower_split = (X[column].values)[np.abs(frequencies) <= split_frequencies[column]]
    upper_split = (X[column].values)[np.abs(frequencies) >  split_frequencies[column]]
    out[column+"_area_lower"] = integrate.trapz(np.abs(lower_split))
    out[column+"_area_upper"] = integrate.trapz(np.abs(upper_split))

  return out

# TODO FFT on the windows<jupyter_output><empty_output><jupyter_text>Feature Selection<jupyter_code>def select_features(X: pd.DataFrame, y: np.array):
  return ts.select_features(X, y)<jupyter_output><empty_output><jupyter_text>Visualization<jupyter_code>def visualize_signal(
    df: pd.DataFrame, signal: str, unit: str, label, *args: Tuple[pd.DataFrame, str], frequency=False
):
    fig = go.Figure(
        go.Scatter(
            mode="markers+lines",
            x=df["Timestamp"] if not frequency else fftpack.fftfreq(len(df[signal if signal not in {"A", "G", "Jerk_"} else signal + "x"]), d=1/104),
            y=df[signal if signal not in {"A", "G", "Jerk_"} else signal + "x"],
            marker={"size": 3},
            name=label,
        )
    )

    fig.update_layout(
        xaxis={"title": "Time (s)"},
        yaxis={"title": signal.replace("_", "") + " (" + unit + ")"},
        width=1920,
        height=1080,
        font=dict(size=18),
        template="plotly_white",
    )

    if signal not in {"A", "G", "Jerk_"}:
        for _df, _label in args:
            fig.add_trace(
                go.Scatter(
                    mode="markers+lines",
                    x=_df["Timestamp"],
                    y=_df[signal],
                    marker={"size": 3},
                    name=_label,
                )
            )
    else:
        dfs: List[pd.DataFrame] = []
        labels: List[str] = []
        for t in args:
            dfs.append(t[0])
            labels.append(t[1])
        for _df, _label, axis in zip(dfs, labels, ["y", "z"]):
            fig.add_trace(
                go.Scatter(
                    mode="markers+lines",
                    x=_df["Timestamp"],
                    y=_df["Speed" if signal == "Speed" else signal + axis],
                    marker={"size": 3},
                    name=_label,
                )
            )
    return fig<jupyter_output><empty_output><jupyter_text>Analysis  Data loading<jupyter_code>X, X_deck, X_stem, y_weight, y_passenger, y_weight_deck, y_passengers_deck, y_weight_stem, y_passengers_stem = read_data("./data/simplified_datasets_1e-06")
#X_aligned = align(X)<jupyter_output><empty_output><jupyter_text>Path simplification<jupyter_code>X_deck_simplified = simplify_path(X_deck)<jupyter_output><empty_output><jupyter_text>Windowing<jupyter_code>X_deck_windowed, y_weight_deck_windowed, y_passengers_deck_windowed = window(X_deck, y_weight_deck, y_passengers_deck, rolling_direction=1, min_timeshift=20,max_timeshift=20)
X_stem_windowed, y_weight_stem_windowed, y_passengers_stem_windowed = window(X_stem, y_weight_stem, y_passengers_stem, rolling_direction=1, min_timeshift=20,max_timeshift=20)<jupyter_output>Rolling: 100%|██████████| 38/38 [00:02<00:00, 15.67it/s]
Rolling: 100%|██████████| 39/39 [00:02<00:00, 14.95it/s]<jupyter_text>Feature Extraction<jupyter_code>deck_features = extract_features(X_deck_windowed.drop("full_id", axis=1))
stem_features = extract_features(X_stem_windowed.drop("full_id", axis=1))
deck_features.to_csv("features_deck_minimal.csv")
stem_features.to_csv("features_stem_minimal.csv")
weight_stem_features = select_features(stem_features, y_weight_stem_windowed)
weight_stem_features.to_csv("selected_features_stem_windowed_20_weight.csv")
passengers_stem_features = select_features(stem_features, y_passengers_stem_windowed)
passengers_stem_features.to_csv("selected_features_stem_windowed_20_passengers.csv")
weight_deck_features = select_features(deck_features, y_weight_deck_windowed)
weight_deck_features.to_csv("selected_features_deck_windowed_20_weight.csv")
passengers_deck_features = select_features(deck_features, y_passengers_deck_windowed)
passengers_deck_features.to_csv("selected_features_deck_windowed_20_passengers.csv")<jupyter_output><empty_output><jupyter_text>Pipeline<jupyter_code>cachedir = mkdtemp()
memory = joblib.Memory(location=cachedir, verbose=10)

pipeline = sklearn.pipeline.Pipeline([#("align", sk.preprocessing.FunctionTransformer(align))
                                 ("windowing", skt.forecasting.model_selection.SlidingWindowSplitter),
                                 ("simplify_paths", sk.preprocessing.FunctionTransformer(simplify_path)),
                                 ("extract_features", sk.preprocessing.FunctionTransformer(extract_features)),
                                 ("pca", sk.decomposition.PCA()),
                                 #("select_features", sk.preprocessing.FunctionTransformer(select_features)),
                                 ("logistic_regression", sk.linear_model.LogisticRegression(max_iter=10000, tol=0.1))],
                                memory=memory)<jupyter_output><empty_output><jupyter_text>Cross Validation<jupyter_code>param_grid = {
    "window__min_timeshift": [0, 10, 100, 1000, 10000],
    "window__max_timeshift": [None, 10, 100, 1000, 10000],
    "simplify_paths": [0, 1e-6, 1e-9, 1e-12, 1e-15],
    "pca__n_components": [5, 15, 30, 45, 64],
    "logistic__C": np.logspace(-4, 4, 4),
}
search = sk.model_selection.GridSearchCV(pipeline, param_grid, n_jobs=-1)
search.fit(pd.DataFrame, y_passenger)<jupyter_output><empty_output><jupyter_text>Fit  Weight Fit<jupyter_code>search.fit(X, y_weight)<jupyter_output><empty_output><jupyter_text>Passenger Fit<jupyter_code>search.fit(X, y_passenger)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13437107357252548, 'fraction_numerical': 0.02191819070221974, 'mean_word_length': 3.596406801411614, 'pattern_counts': {'":': 15, '<': 82, '<?xml version=': 0, '>': 79, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '19', 'dir': 'jupyter-structured-clean-dedup', 'id': '175357', 'license': 'MIT', 'n_tokens_mistral': 5337, 'n_tokens_neox': 5025, 'n_words': 851, 'path': 'mida_acv.ipynb', 'repo_name': 'DPaletti/mida_acv', 'size': '14326'}