Travel Time Analysis Project Series 08: Modeling — Feature Selection

Image by Gordon Johnson from Pixabay

[English / 中文] If you are not a member you can read the story here.

This is the eighth article in the Travel Time Analysis Project series. The goal of this project is to leverage data related to accidents, road construction, holidays, and other relevant factors to conduct time series analysis and forecasting. Our aim is to identify the most suitable neural network model for predicting travel times on specific segments of Taiwan’s national highways in designated directions. If you’re interested in learning more about the project’s origins and initial concepts, please refer to the first article in the series: Travel Time Analysis Project Series 01: Overview and Introduction.

Why Is Feature Selection Necessary?

In theory, we can feed all available variables directly into a deep learning model. These models are designed to learn patterns automatically, and in domains like image or speech recognition, this “throw everything in” strategy can work well. However, it also comes with some serious downsides.

When the number of features grows large, training becomes more expensive — both in time and computational resources. It also increases the risk of overfitting, especially when the amount of data isn’t enough to support that many features. On top of that, many variables might be highly correlated or redundant, making it harder for the model to focus on what truly matters.

By applying feature selection before training, we can filter out less informative variables, reduce the dimensionality of the data, and speed up the training process. This not only helps the model generalize better to new data, but also improves interpretability — giving us clearer insights into what drives the predictions. Even in deep learning workflows, thoughtful feature selection remains a valuable and often necessary step.

Feature Selection Methods

Feature selection methods generally fall into three categories: Filter, Wrapper, and Embedded approaches.

Filter Method

• Performed independently of any model, this method ranks and selects features based on their statistical relationship with the target variable.
• Fast and model-agnostic, it is suitable for initial screening. However, some domain knowledge is often needed to make informed decisions.
• Common techniques include:
• Correlation coefficients: such as Pearson or Spearman correlation.
• Univariate selection: such as Chi-square test or ANOVA F-test.
• Variance thresholding: removes features with very low variance (i.e., nearly constant).
• Mutual information: quantifies the amount of shared information between variables.

Wrapper Method

• This method evaluates different combinations of features by training a model on each subset and selecting the best-performing group.
• It takes feature interactions into account and often yields better results.
• However, it is computationally expensive — especially when dealing with many features — as each feature subset requires model training and evaluation. The process can be time-consuming and resource-intensive.
• Common techniques include:
• Recursive Feature Elimination (RFE): iteratively removes the least important features based on model performance.
• Forward/Backward Selection: starts with an empty set or full set of features and adds or removes one feature at a time, evaluating model performance at each step.

Embedded Method

• Feature selection is built into the model training process itself, using the model’s internal weights or structure to assess feature importance.
• It strikes a good balance between efficiency and accuracy, making it suitable for high-dimensional datasets.
• Common techniques include:
• L1 regularization (Lasso Regression): automatically drives unimportant feature weights to zero.
• Tree-based models (e.g., Random Forest, XGBoost): inherently provide feature importance metrics.
• Elastic Net: combines L1 and L2 regularization to balance feature selection and model stability.

In this project, I ultimately chose the Embedded Method, implemented via Random Forest. One of the key advantages is that feature importance scores derived from Random Forests are interpretable, helping us understand which gantry pairs most strongly influence travel time.

Implementation

In this stage, we once again leveraged the previously developed hwttp toolkit. Since we needed to extract and inspect all features, hwttp was especially useful thanks to its built-in feature combinations. The following code is organized into three parts:

• Data Loading: We load the cleaned dataset and integrate additional information such as accidents, congestion time by segment, and holiday indicators based on defined relationships.
• Feature Composition: The hwttp.hwtoolkit module includes prebuilt methods for joining and combining variables.
• Feature Importance and Selection: Using Random Forest to compute feature importance, which is the core focus of this analysis.

import os
import sys
import time
from pathlib import Path
from os import listdir
from os.path import isfile, join

import numpy as np
import pandas as pd
from tqdm import tqdm
import matplotlib.pyplot as plt

from sklearn.ensemble import RandomForestRegressor
from sklearn.metrics import mean_squared_error
from sklearn.model_selection import TimeSeriesSplit

sys.path.append(os.path.abspath('..'))
import hwttp.hwtoolkit as tk

Load Data

data_paths = {'etag_5n_loc': '../data/cleaned/etag_5n_loc.csv',
              'section_info': '../data/cleaned/section_info.csv',
              'hw5_m04a_df': '../data/cleaned/hw5_m04a.csv',
              'congestion_table': '../data/cleaned/congestion_table.csv',
              'calendar_event': '../data/cleaned/calendar_event.csv',
              'road_build_event': '../data/cleaned/202301_10_road_build_event.xlsx',
              'traffic_accident_data': '../data/cleaned/202301_10_traffic_accident_data.xlsx'
             }

rs = tk.hw_df_resource(data_paths)
rs.load_raw_environment_info()
rs.load_raw_event_info()
rs.generate_mile_location_info()

# remove irrelevant holidays
drop_event_list = ['婦女節', '愚人節', '復活節', '地球日', '感恩節', 'Black Friday', 'Cyber Monday', '冬至']
rs.calendar_event = rs.calendar_event[~rs.calendar_event.event_name.isin(drop_event_list)].copy()

# Travel time data
hw5_15watt = pd.read_csv("../data/features/hw5_15watt.csv")
hw5_15watt = hw5_15watt[~hw5_15watt.gf_gt.isin(['03F0150N-03F0140N', 
                                                '03F0201S-03A0041N', 
                                                '03F0201S-03F0217S'])].copy()
hw5_15watt = hw5_15watt[hw5_15watt.TimeStamp<'2024-01-01 00:00:00'].copy()
hw5_15watt['TimeStamp'] = pd.to_datetime(hw5_15watt['TimeStamp'])

At this point, the core dataset contains 280,320 rows and 6 columns. No additional derived variables have been appended yet.

<class 'pandas.core.frame.DataFrame'>
Int64Index: 280320 entries, 131328 to 472799
Data columns (total 6 columns):
 #   Column                 Non-Null Count   Dtype         
---  ------                 --------------   -----         
 0   gf_gt                  280320 non-null  object        
 1   GantryFrom             280320 non-null  object        
 2   GantryTo               280320 non-null  object        
 3   TimeStamp              280320 non-null  datetime64[ns]
 4   WeightedAvgTravelTime  280320 non-null  float64       
 5   TotalTraffic           280320 non-null  float64   

Combinations

In the initial experiments, I opted for a wrapper method — Recursive Feature Elimination (RFE) — to train and validate the model. However, after feeding in a large number of features, it became clear that the method was impractical; each training iteration took far too long to complete.

To address this, I defined function-based handling to systematically expand the feature set. At this stage, all available feature types are included. Below are the simple letter codes I used to represent each type, which together form the full set of features used to estimate travel time between gantries:

• c: Congestion duration for the corresponding road segment
• h: Holiday indicators
• t: Traffic accident information
• r: Road construction data
• p: Travel time from the downstream segment in the previous five minutes

# function of appending features to source dataframe 
def add_c(df):
    output_df = tk.add_congestion_condition(df, rs.congestion_table, rs.milelocation_info_df).copy()
    return output_df

def add_h(df):
    output_df = tk.add_calendar_event(df, rs.calendar_event).copy()
    return output_df

def add_t(df):
    *_, output_df = tk.add_traffic_event(df, rs.traffic_accident_data, rs.milelocation_info_df)
    return output_df

def add_r(df):
    *_, output_df = tk.add_road_build_event(df, rs.road_build_event, rs.milelocation_info_df)
    return output_df

def add_p(df):
    output_df = tk.add_ds_5prev_traveltime(df)
    return output_df

# put function into the dictionary
element_dict = {'_c': add_c, 
                '_h': add_h, 
                '_t': add_t, 
                '_r': add_r, 
                '_p': add_p
               }

# combinations for gantry pair
combinations = {'05F0001N-03F0150N': ['_c', '_h', '_t', '_r', '_p'],
                '05F0001N-03F0201S': ['_c', '_h', '_t', '_r', '_p'],
                '05F0055N-05F0001N': ['_c', '_h', '_t', '_r', '_p'],
                '05F0287N-05F0055N': ['_c', '_h', '_t', '_r', '_p'],
                '05F0309N-05F0287N': ['_c', '_h', '_t', '_r', '_p'],
                '05F0438N-05F0309N': ['_c', '_h', '_t', '_r', '_p'],
                '05F0438N-05FR143N': ['_c', '_h', '_t', '_r', '_p'],
                '05F0528N-05F0438N': ['_c', '_h', '_t', '_r', '_p'],
               }

Due to the way the data structure was designed, the process of attaching features to the travel time records wasn’t as smooth as expected. Specifically, upstream and downstream-related features had to be merged first before other types could be appended properly.

Each feature-adding step can be executed by calling the corresponding function from the previously defined element_dict, making the process modular and easy to manage.

# dict to store results
results = dict()

# implement all the combinations
for gantry_pair, combo in tqdm(combinations.items()):
    # initialization
    result_df = hw5_15watt.copy()
    output_string = gantry_pair+'_b'
    print(gantry_pair, combo)
    combo_list = combo.copy()

    # bad handling, need to join _p feature first, or multiple columns will appeared.
    if '_p' in combo_list:
        combo_list.remove('_p')
        func = element_dict['_p']
        result_df = func(result_df)
        output_string += '_p'

    print('combo_list = ', combo_list)

    # join rest features from the combinations
    for func_name in combo_list:
        func = element_dict[func_name]
        result_df = func(result_df)
        output_string += func_name

    # at last, limit to specific gantry pair
    result_df = result_df[result_df['gf_gt']==gantry_pair].copy()

    # save to the results dict
    results[output_string] = result_df.copy()

This is one of the completed data tables, containing a total of 63 columns.

<class 'pandas.core.frame.DataFrame'>
Int64Index: 35040 entries, 0 to 35039
Data columns (total 63 columns):
 #   Column                           Non-Null Count  Dtype         
---  ------                           --------------  -----         
 0   gf_gt                            35040 non-null  object        
 1   GantryFrom                       35040 non-null  object        
 2   GantryTo                         35040 non-null  object        
 3   TimeStamp                        35040 non-null  datetime64[ns]
 4   WeightedAvgTravelTime            35040 non-null  float64       
 5   TotalTraffic                     35040 non-null  float64       
 6   ds_prev_1_WATT                   35040 non-null  float64       
 7   ds_prev_2_WATT                   35040 non-null  float64       
 8   ds_prev_3_WATT                   35040 non-null  float64       
 9   ds_prev_4_WATT                   35040 non-null  float64       
 10  ds_prev_5_WATT                   35040 non-null  float64       
 11  congestion_syndrome              35040 non-null  int64         
 12  holiday_continue                 35040 non-null  int64         
 13  holiday_length                   35040 non-null  float64       
 14  dayofweek                        35040 non-null  int64         
 15  holiday_name                     4416 non-null   object        
 16  holiday_name_七夕情人節               35040 non-null  int64         
 17  holiday_name_中元節                 35040 non-null  int64         
 18  holiday_name_中秋節                 35040 non-null  int64         
 19  holiday_name_二二八紀念日              35040 non-null  int64         
 20  holiday_name_元宵節                 35040 non-null  int64         
 21  holiday_name_兒童與清明節              35040 non-null  int64         
 22  holiday_name_勞動節                 35040 non-null  int64         
 23  holiday_name_國慶日                 35040 non-null  int64         
 24  holiday_name_教師節                 35040 non-null  int64         
 25  holiday_name_母親節                 35040 non-null  int64         
 26  holiday_name_父親節                 35040 non-null  int64         
 27  holiday_name_白色情人節               35040 non-null  int64         
 28  holiday_name_端午節                 35040 non-null  int64         
 29  holiday_name_聖誕節                 35040 non-null  int64         
 30  holiday_name_萬聖節                 35040 non-null  int64         
 31  holiday_name_西洋情人節               35040 non-null  int64         
 32  holiday_name_跨年元旦                35040 non-null  int64         
 33  holiday_name_農曆新年                35040 non-null  int64         
 34  accident_mileage                 35040 non-null  float64       
 35  event_occurrence                 35040 non-null  object        
 36  event_exclusion                  35040 non-null  object        
 37  handling_minutes                 35040 non-null  float64       
 38  accident_type                    35040 non-null  int64         
 39  death_count                      35040 non-null  float64       
 40  injuries_count                   35040 non-null  float64       
 41  inner_shoulder_flag              35040 non-null  float64       
 42  inner_lane_flag                  35040 non-null  float64       
 43  middle_inner_lane_flag           35040 non-null  float64       
 44  middle_lane_flag                 35040 non-null  float64       
 45  middle_outer_lane_flag           35040 non-null  float64       
 46  outer_lane_flag                  35040 non-null  float64       
 47  outer_shoulder_flag              35040 non-null  float64       
 48  ramp_flag                        35040 non-null  float64       
 49  overturn_accident_flag           35040 non-null  float64       
 50  construction_accident_flag       35040 non-null  float64       
 51  hazardous_material_vehicle_flag  35040 non-null  float64       
 52  on_fire_vehicle_flag             35040 non-null  float64       
 53  smoking_vehicle_flag             35040 non-null  float64       
 54  mainlane_disruption_flag         35040 non-null  float64       
 55  accident_vehicle_count           35040 non-null  float64       
 56  light_truck_count                35040 non-null  float64       
 57  passenger_car_count              35040 non-null  float64       
 58  bus_count                        35040 non-null  float64       
 59  heavy_truck_count                35040 non-null  float64       
 60  road_build                       35040 non-null  int64         
 61  total_block_count                35040 non-null  int64         
 62  road_block_count                 35040 non-null  int64         
dtypes: datetime64[ns](1), float64(31), int64(25), object(6)
memory usage: 17.1+ MB

Make sure to save the data in Parquet format — this will speed up data loading in the later stages before running model cross-validation.

# fast save, next time no need to wait for long processing
for df_name, df in results.items():
    object_columns = df.select_dtypes(include='object').columns
    df[object_columns] = df[object_columns].astype('string')
    df.to_parquet(f'../data/features/all_features_by_intergrantry/{df_name}.parquet')

Feature Importance | Feature Selection

The feature importance calculation method used in this study is particularly interesting. By applying TimeSeriesSplit, the dataset is incrementally divided into multiple time-based folds. Each fold serves as the basis for feature importance computation, effectively simulating a scenario of cumulative time progression. This design helps capture the influence of features that emerge over long-term temporal changes. Feature importance scores from each fold are aggregated and averaged, resulting in a stable estimate that reflects long-term patterns.

The core model used is Random Forest. As an ensemble of multiple decision trees, Random Forest remains robust even under high data variability, and its resulting feature importance scores are generally reliable.

In each loop, data is filtered according to different gantry pairs, and the rolling_window_feature_selection function is called with the model and data to perform feature selection. The results are then saved into the results_fs dictionary for later use.

def rolling_window_feature_selection(X, y, model, n_splits=5):
    tscv = TimeSeriesSplit(n_splits=n_splits)
    feature_importances = np.zeros(X.shape[1])

    for train_index, test_index in tscv.split(X):
        X_train, X_test = X.iloc[train_index], X.iloc[test_index]
        y_train, y_test = y.iloc[train_index], y.iloc[test_index]

        model.fit(X_train, y_train)
        y_pred = model.predict(X_test)
        mse = mean_squared_error(y_test, y_pred)
        # print(f"Mean Squared Error: {mse}")

        # accumulate feature importance
        feature_importances += model.feature_importances_

    # calculate average feature importance
    feature_importances /= n_splits
    return feature_importances

# store feature selection's results
results_fs = dict()

for item in tqdm(results.items()):
    cur_df_name = item[0]
    cur_df = item[1]

    if cur_df_name in results_fs.keys():
        continue

    # separate features and target variables
    X = cur_df.drop(columns=['WeightedAvgTravelTime', 'gf_gt', 'GantryFrom', 'GantryTo', 'TimeStamp', 'TotalTraffic'])
    y = cur_df['WeightedAvgTravelTime']

    # define base model
    base_model = RandomForestRegressor(n_estimators=100, random_state=42)

    # calculate feature importance
    feature_importances = rolling_window_feature_selection(X, y, base_model, n_splits=180)
    feature_importances_df = pd.DataFrame({'feature': X.columns,
                                           'importance': feature_importances
                                          }).sort_values(by='importance', ascending=False)

    # display most important feature
    selected_features = feature_importances_df['feature'].head(10).values
    print(f"df = {cur_df_name} Selected features: {selected_features}")

    results_fs[cur_df_name] = feature_importances_df.copy()

print('feature selection complete')

Lastly, we take a quick look at the feature importance distribution for a few selected gantry pairs. Here, we present only three plots as reference.

for i in results_fs.keys():
    print(i)
    display(results_fs[i].query('importance > 0'))

image by author

In the following experiments, I will use 0.00005 as the threshold for feature selection. This value is based on practical experience. While it is theoretically possible to retain all features for training, doing so would significantly increase training time without necessarily improving model performance. Setting a reasonable threshold helps preserve sufficient information while keeping the model efficient.

For the pair 05F0001N-03F0150N, the day of the week shows the strongest impact on travel time, followed by the length of consecutive holidays, and then effects from specific holidays. For 05F0287N-05F0055N, the most significant factor is the travel time from downstream segments, followed by day-of-week and some minor congestion-related features. In contrast, 05F0438N-05FR143N appears more influenced by road construction and day of the week.

In summary, feature importance varies significantly across different road segments. By using feature importance as a basis for selection, we can ensure that the resulting models are more aligned with the specific travel time characteristics of each target segment.

(Remember to save all feature importance results for use in the upcoming multivariate modeling phase.)

To improve the efficiency of upcoming multivariate model development, we brought back the features we carefully constructed earlier and used a Random Forest model — an embedded method — to evaluate feature importance. Now, we’re entering the part of this series that excites me the most, and also the final stage of the modeling process: Multivariate Modeling.

Note: This article was initially written in Chinese, then translated by ChatGPT and adjusted by myself. Slight differences in meaning may exist.

Travel Time Analysis Project Series

• Overview and Introduction
• Data Crawling (資料爬蟲)
• Data Cleaning (資料清洗)
• Data Exploration (資料探索)
• Preliminary EDA (初步探索性數據分析)
• EDA for National Freeway №5 (國道五號的探索性數據分析)
• Time Series Analysis and Modeling (時間序列分析與建模)
• Baseline Modeling (基線模型建模)
• Univariate Modeling (單變量建模)
• Feature Selection (特徵篩選)
• Multivariate Modeling (多變量建模)
• Cross-Validation and Performance Analysis (交叉驗證與成果分析)
• Closing Remarks

If you enjoyed this article, you can:

• 👏 Clap a few times — no need to hit all 50
• ✉️Follow me on Medium to stay updated and keep me motivated
• 😜 Connect with me on LinkedIn — I’d love to exchange ideas