{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":102335,"databundleVersionId":12518947,"sourceType":"competition"}],"dockerImageVersionId":31090,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import os, json, joblib, numpy as np, pandas as pd\nimport random\nfrom pathlib import Path\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nfrom sklearn.model_selection import train_test_split, StratifiedGroupKFold\nfrom sklearn.preprocessing import StandardScaler, LabelEncoder\nfrom sklearn.utils.class_weight import compute_class_weight\nfrom sklearn.metrics import classification_report, confusion_matrix\n\nfrom tensorflow.keras.utils import Sequence, to_categorical, pad_sequences\nfrom tensorflow.keras.models import Model, load_model\nfrom tensorflow.keras.layers import (\n    Input, Conv1D, BatchNormalization, Activation, add, MaxPooling1D, Dropout,\n    Bidirectional, LSTM, GlobalAveragePooling1D, Dense, Multiply, Reshape,\n    Lambda, Concatenate, GRU, GaussianNoise\n)\nfrom tensorflow.keras.regularizers import l2\nfrom tensorflow.keras.optimizers import Adam\nfrom tensorflow.keras.callbacks import EarlyStopping, ModelCheckpoint, ReduceLROnPlateau\nfrom tensorflow.keras import backend as K\nimport tensorflow as tf\nimport polars as pl\nfrom scipy.spatial.transform import Rotation as R\n\ndef seed_everything(seed):\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    random.seed(seed)\n    np.random.seed(seed)\n    tf.random.set_seed(seed)\n    tf.experimental.numpy.random.seed(seed)\n    os.environ['TF_CUDNN_DETERMINISTIC'] = '1'\n    os.environ['TF_DETERMINISTIC_OPS'] = '1'\n\nseed_everything(seed=42)\n\n# Training Configuration\nTRAIN = True\nRAW_DIR = Path(\"/kaggle/input/cmi-detect-behavior-with-sensor-data\")  # Update this to your data directory\nEXPORT_DIR = Path(\"./models\")\nEXPORT_DIR.mkdir(exist_ok=True)\n\nBATCH_SIZE = 128\nPAD_PERCENTILE = 95\nLR_INIT = 5e-4\nWD = 3e-3\nMIXUP_ALPHA = 0.4\nEPOCHS = 100\nPATIENCE = 40\nN_FOLDS = 5\n\nprint(\"▶ Training mode activated · tensorflow\", tf.__version__)\n\n#=============================================================================\n# Utility Functions (same as inference)\n#=============================================================================\n\ndef time_sum(x):\n    return K.sum(x, axis=1)\n\ndef squeeze_last_axis(x):\n    return tf.squeeze(x, axis=-1)\n\ndef expand_last_axis(x):\n    return tf.expand_dims(x, axis=-1)\n\ndef se_block(x, reduction=8):\n    ch = x.shape[-1]\n    se = GlobalAveragePooling1D()(x)\n    se = Dense(ch // reduction, activation='relu')(se)\n    se = Dense(ch, activation='sigmoid')(se)\n    se = Reshape((1, ch))(se)\n    return Multiply()([x, se])\n\ndef residual_se_cnn_block(x, filters, kernel_size, pool_size=2, drop=0.3, wd=1e-4):\n    shortcut = x\n    for _ in range(2):\n        x = Conv1D(filters, kernel_size, padding='same', use_bias=False,\n                   kernel_regularizer=l2(wd))(x)\n        x = BatchNormalization()(x)\n        x = Activation('relu')(x)\n        x = se_block(x)\n    \n    if shortcut.shape[-1] != filters:\n        shortcut = Conv1D(filters, 1, padding='same', use_bias=False,\n                         kernel_regularizer=l2(wd))(shortcut)\n        shortcut = BatchNormalization()(shortcut)\n    \n    x = add([x, shortcut])\n    x = Activation('relu')(x)\n    x = MaxPooling1D(pool_size)(x)\n    x = Dropout(drop)(x)\n    return x\n\ndef attention_layer(inputs):\n    score = Dense(1, activation='tanh')(inputs)\n    score = Lambda(squeeze_last_axis)(score)\n    weights = Activation('softmax')(score)\n    weights = Lambda(expand_last_axis)(weights)\n    context = Multiply()([inputs, weights])\n    context = Lambda(time_sum)(context)\n    return context\n\nclass MixupGenerator(Sequence):\n    def __init__(self, X, y, batch_size, alpha=0.2):\n        self.X, self.y = X, y\n        self.batch = batch_size\n        self.alpha = alpha\n        self.indices = np.arange(len(X))\n    \n    def __len__(self):\n        return int(np.ceil(len(self.X) / self.batch))\n    \n    def __getitem__(self, i):\n        idx = self.indices[i*self.batch:(i+1)*self.batch]\n        Xb, yb = self.X[idx], self.y[idx]\n        lam = np.random.beta(self.alpha, self.alpha)\n        perm = np.random.permutation(len(Xb))\n        X_mix = lam * Xb + (1-lam) * Xb[perm]\n        y_mix = lam * yb + (1-lam) * yb[perm]\n        return X_mix, y_mix\n    \n    def on_epoch_end(self):\n        np.random.shuffle(self.indices)\n\ndef remove_gravity_from_acc(acc_data, rot_data):\n    if isinstance(acc_data, pd.DataFrame):\n        acc_values = acc_data[['acc_x', 'acc_y', 'acc_z']].values\n    else:\n        acc_values = acc_data\n    \n    if isinstance(rot_data, pd.DataFrame):\n        quat_values = rot_data[['rot_x', 'rot_y', 'rot_z', 'rot_w']].values\n    else:\n        quat_values = rot_data\n    \n    num_samples = acc_values.shape[0]\n    linear_accel = np.zeros_like(acc_values)\n    gravity_world = np.array([0, 0, 9.81])\n    \n    for i in range(num_samples):\n        if np.all(np.isnan(quat_values[i])) or np.all(np.isclose(quat_values[i], 0)):\n            linear_accel[i, :] = acc_values[i, :]\n            continue\n        \n        try:\n            rotation = R.from_quat(quat_values[i])\n            gravity_sensor_frame = rotation.apply(gravity_world, inverse=True)\n            linear_accel[i, :] = acc_values[i, :] - gravity_sensor_frame\n        except ValueError:\n            linear_accel[i, :] = acc_values[i, :]\n    \n    return linear_accel\n\ndef calculate_angular_velocity_from_quat(rot_data, time_delta=1/200):\n    if isinstance(rot_data, pd.DataFrame):\n        quat_values = rot_data[['rot_x', 'rot_y', 'rot_z', 'rot_w']].values\n    else:\n        quat_values = rot_data\n    \n    num_samples = quat_values.shape[0]\n    angular_vel = np.zeros((num_samples, 3))\n    \n    for i in range(num_samples - 1):\n        q_t = quat_values[i]\n        q_t_plus_dt = quat_values[i+1]\n        \n        if np.all(np.isnan(q_t)) or np.all(np.isclose(q_t, 0)) or \\\n           np.all(np.isnan(q_t_plus_dt)) or np.all(np.isclose(q_t_plus_dt, 0)):\n            continue\n        \n        try:\n            rot_t = R.from_quat(q_t)\n            rot_t_plus_dt = R.from_quat(q_t_plus_dt)\n            delta_rot = rot_t.inv() * rot_t_plus_dt\n            angular_vel[i, :] = delta_rot.as_rotvec() / time_delta\n        except ValueError:\n            pass\n    \n    return angular_vel\n\ndef calculate_angular_distance(rot_data):\n    if isinstance(rot_data, pd.DataFrame):\n        quat_values = rot_data[['rot_x', 'rot_y', 'rot_z', 'rot_w']].values\n    else:\n        quat_values = rot_data\n    \n    num_samples = quat_values.shape[0]\n    angular_dist = np.zeros(num_samples)\n    \n    for i in range(num_samples - 1):\n        q1 = quat_values[i]\n        q2 = quat_values[i+1]\n        \n        if np.all(np.isnan(q1)) or np.all(np.isclose(q1, 0)) or \\\n           np.all(np.isnan(q2)) or np.all(np.isclose(q2, 0)):\n            angular_dist[i] = 0\n            continue\n        try:\n            r1 = R.from_quat(q1)\n            r2 = R.from_quat(q2)\n            relative_rotation = r1.inv() * r2\n            angle = np.linalg.norm(relative_rotation.as_rotvec())\n            angular_dist[i] = angle\n        except ValueError:\n            angular_dist[i] = 0\n    \n    return angular_dist\n\ndef build_two_branch_model(pad_len, imu_dim, tof_dim, n_classes, wd=1e-4):\n    inp = Input(shape=(pad_len, imu_dim+tof_dim))\n    imu = Lambda(lambda t: t[:, :, :imu_dim])(inp)\n    tof = Lambda(lambda t: t[:, :, imu_dim:])(inp)\n    \n    # IMU deep branch\n    x1 = residual_se_cnn_block(imu, 64, 3, drop=0.1, wd=wd)\n    x1 = residual_se_cnn_block(x1, 128, 5, drop=0.1, wd=wd)\n    \n    # TOF/Thermal lighter branch\n    x2 = Conv1D(64, 3, padding='same', use_bias=False, kernel_regularizer=l2(wd))(tof)\n    x2 = BatchNormalization()(x2); x2 = Activation('relu')(x2)\n    x2 = MaxPooling1D(2)(x2); x2 = Dropout(0.2)(x2)\n    x2 = Conv1D(128, 3, padding='same', use_bias=False, kernel_regularizer=l2(wd))(x2)\n    x2 = BatchNormalization()(x2); x2 = Activation('relu')(x2)\n    x2 = MaxPooling1D(2)(x2); x2 = Dropout(0.2)(x2)\n    \n    merged = Concatenate()([x1, x2])\n    \n    xa = Bidirectional(LSTM(128, return_sequences=True, kernel_regularizer=l2(wd)))(merged)\n    xb = Bidirectional(GRU(128, return_sequences=True, kernel_regularizer=l2(wd)))(merged)\n    xc = GaussianNoise(0.09)(merged)\n    xc = Dense(16, activation='elu')(xc)\n    \n    x = Concatenate()([xa, xb, xc])\n    x = Dropout(0.4)(x)\n    x = attention_layer(x)\n    \n    for units, drop in [(256, 0.5), (128, 0.3)]:\n        x = Dense(units, use_bias=False, kernel_regularizer=l2(wd))(x)\n        x = BatchNormalization()(x); x = Activation('relu')(x)\n        x = Dropout(drop)(x)\n    \n    out = Dense(n_classes, activation='softmax', kernel_regularizer=l2(wd))(x)\n    return Model(inp, out)\n\n#=============================================================================\n# Data Loading and Feature Engineering\n#=============================================================================\n\ndef load_and_preprocess_data():\n    \"\"\"Load and preprocess the training data\"\"\"\n    print(\"Loading training data...\")\n    \n    # Load main training data\n    train_df = pd.read_csv(RAW_DIR / \"train.csv\")\n    train_demographics = pd.read_csv(RAW_DIR / \"train_demographics.csv\")\n    \n    print(f\"Training data shape: {train_df.shape}\")\n    print(f\"Unique sequences: {train_df['sequence_id'].nunique()}\")\n    print(f\"Gesture distribution:\\n{train_df['gesture'].value_counts()}\")\n    \n    return train_df, train_demographics\n\ndef engineer_features(df_seq):\n    \"\"\"Apply feature engineering to a sequence\"\"\"\n    # Remove gravity from accelerometer data\n    linear_accel = remove_gravity_from_acc(df_seq, df_seq)\n    df_seq['linear_acc_x'] = linear_accel[:, 0]\n    df_seq['linear_acc_y'] = linear_accel[:, 1] \n    df_seq['linear_acc_z'] = linear_accel[:, 2]\n    df_seq['linear_acc_mag'] = np.sqrt(df_seq['linear_acc_x']**2 + \n                                      df_seq['linear_acc_y']**2 + \n                                      df_seq['linear_acc_z']**2)\n    df_seq['linear_acc_mag_jerk'] = df_seq['linear_acc_mag'].diff().fillna(0)\n    \n    # Calculate angular velocity\n    angular_vel = calculate_angular_velocity_from_quat(df_seq)\n    df_seq['angular_vel_x'] = angular_vel[:, 0]\n    df_seq['angular_vel_y'] = angular_vel[:, 1]\n    df_seq['angular_vel_z'] = angular_vel[:, 2]\n    df_seq['angular_distance'] = calculate_angular_distance(df_seq)\n    \n    # TOF sensor aggregations\n    for i in range(1, 6):\n        pixel_cols = [f\"tof_{i}_v{p}\" for p in range(64)]\n        tof_data = df_seq[pixel_cols].replace(-1, np.nan)\n        df_seq[f'tof_{i}_mean'] = tof_data.mean(axis=1)\n        df_seq[f'tof_{i}_std'] = tof_data.std(axis=1)\n        df_seq[f'tof_{i}_min'] = tof_data.min(axis=1)\n        df_seq[f'tof_{i}_max'] = tof_data.max(axis=1)\n    \n    return df_seq\n\ndef prepare_dataset(train_df, train_demographics):\n    \"\"\"Prepare the complete dataset for training\"\"\"\n    print(\"Preparing dataset with feature engineering...\")\n    \n    # Define feature columns (same as inference)\n    base_features = ['acc_x', 'acc_y', 'acc_z', 'rot_w', 'rot_x', 'rot_y', 'rot_z']\n    \n    # Thermal features\n    thermal_features = [f'thm_{i}' for i in range(1, 6)]\n    \n    # Engineered features\n    engineered_features = [\n        'linear_acc_x', 'linear_acc_y', 'linear_acc_z', 'linear_acc_mag', \n        'linear_acc_mag_jerk', 'angular_vel_x', 'angular_vel_y', 'angular_vel_z', \n        'angular_distance'\n    ]\n    \n    # TOF aggregated features\n    tof_agg_features = []\n    for i in range(1, 6):\n        tof_agg_features.extend([f'tof_{i}_mean', f'tof_{i}_std', f'tof_{i}_min', f'tof_{i}_max'])\n    \n    # TOF pixel features\n    tof_pixel_features = []\n    for i in range(1, 6):\n        for p in range(64):\n            tof_pixel_features.append(f'tof_{i}_v{p}')\n    \n    feature_cols = base_features + thermal_features + engineered_features + tof_agg_features + tof_pixel_features\n    \n    # Process sequences\n    sequences = []\n    labels = []\n    subjects = []\n    sequence_ids = []\n    \n    unique_seq_ids = train_df['sequence_id'].unique()\n    total_sequences = len(unique_seq_ids)\n    \n    for idx, seq_id in enumerate(unique_seq_ids):\n        if idx % 100 == 0:\n            print(f\"Processing sequence {idx+1}/{total_sequences}: {seq_id}\")\n            \n        seq_data = train_df[train_df['sequence_id'] == seq_id].copy()\n        \n        # Apply feature engineering\n        seq_data = engineer_features(seq_data)\n        \n        # Extract features and fill missing values\n        seq_features = seq_data[feature_cols].ffill().bfill().fillna(0).values.astype('float32')\n        \n        sequences.append(seq_features)\n        labels.append(seq_data['gesture'].iloc[0])\n        subjects.append(seq_data['subject'].iloc[0])\n        sequence_ids.append(seq_id)\n    \n    # Calculate padding length\n    sequence_lengths = [len(seq) for seq in sequences]\n    pad_len = int(np.percentile(sequence_lengths, PAD_PERCENTILE))\n    \n    print(f\"Sequence length statistics:\")\n    print(f\"  Min: {min(sequence_lengths)}\")\n    print(f\"  Max: {max(sequence_lengths)}\")\n    print(f\"  Mean: {np.mean(sequence_lengths):.1f}\")\n    print(f\"  {PAD_PERCENTILE}th percentile (pad_len): {pad_len}\")\n    \n    # Pad sequences\n    X = pad_sequences(sequences, maxlen=pad_len, padding='post', \n                     truncating='post', dtype='float32')\n    \n    # Encode labels\n    label_encoder = LabelEncoder()\n    y_encoded = label_encoder.fit_transform(labels)\n    y_categorical = to_categorical(y_encoded)\n    \n    # Fit scaler\n    print(\"Fitting scaler on all training data...\")\n    all_data = np.vstack([seq.reshape(-1, seq.shape[-1]) for seq in X])\n    scaler = StandardScaler()\n    scaler.fit(all_data)\n    \n    # Scale the data\n    X_scaled = np.zeros_like(X)\n    for i in range(len(X)):\n        X_scaled[i] = scaler.transform(X[i])\n    \n    print(f\"Final dataset shape: {X_scaled.shape}\")\n    print(f\"Feature dimensions: {len(feature_cols)}\")\n    print(f\"Number of classes: {len(label_encoder.classes_)}\")\n    \n    return {\n        'X': X_scaled,\n        'y': y_categorical,\n        'y_encoded': y_encoded,\n        'subjects': np.array(subjects),\n        'sequence_ids': np.array(sequence_ids),\n        'labels': labels,\n        'feature_cols': feature_cols,\n        'pad_len': pad_len,\n        'scaler': scaler,\n        'label_encoder': label_encoder,\n        'gesture_classes': label_encoder.classes_\n    }\n\n#=============================================================================\n# Training Pipeline\n#=============================================================================\n\ndef train_model(data_dict, fold_idx, train_idx, val_idx, model_name_prefix=\"model\"):\n    \"\"\"Train a single model for one fold\"\"\"\n    print(f\"\\n{'='*60}\")\n    print(f\"Training {model_name_prefix} - Fold {fold_idx + 1}\")\n    print(f\"{'='*60}\")\n    \n    X, y = data_dict['X'], data_dict['y']\n    pad_len = data_dict['pad_len']\n    n_classes = len(data_dict['gesture_classes'])\n    \n    # Split data\n    X_train, X_val = X[train_idx], X[val_idx]\n    y_train, y_val = y[train_idx], y[val_idx]\n    \n    print(f\"Train samples: {len(X_train)}, Validation samples: {len(X_val)}\")\n    \n    # Calculate class weights\n    y_train_encoded = data_dict['y_encoded'][train_idx]\n    class_weights = compute_class_weight(\n        'balanced', \n        classes=np.unique(y_train_encoded), \n        y=y_train_encoded\n    )\n    class_weight_dict = dict(enumerate(class_weights))\n    \n    # Determine feature dimensions\n    imu_dim = 7  # acc_x, acc_y, acc_z, rot_w, rot_x, rot_y, rot_z\n    total_dim = X.shape[-1]\n    tof_dim = total_dim - imu_dim\n    \n    # Build model\n    model = build_two_branch_model(pad_len, imu_dim, tof_dim, n_classes, wd=WD)\n    \n    # Compile model\n    model.compile(\n        optimizer=Adam(learning_rate=LR_INIT, weight_decay=WD),\n        loss='categorical_crossentropy',\n        metrics=['accuracy']\n    )\n    \n    # Callbacks\n    model_path = EXPORT_DIR / f\"{model_name_prefix}_{fold_idx}.h5\"\n    callbacks = [\n        EarlyStopping(\n            monitor='val_loss',\n            patience=PATIENCE,\n            restore_best_weights=True,\n            verbose=1\n        ),\n        ModelCheckpoint(\n            model_path,\n            monitor='val_loss',\n            save_best_only=True,\n            verbose=1\n        ),\n        ReduceLROnPlateau(\n            monitor='val_loss',\n            factor=0.5,\n            patience=15,\n            min_lr=1e-7,\n            verbose=1\n        )\n    ]\n    \n    # Create mixup generator\n    mixup_gen = MixupGenerator(X_train, y_train, BATCH_SIZE, MIXUP_ALPHA)\n    \n    # Train model\n    history = model.fit(\n        mixup_gen,\n        epochs=EPOCHS,\n        validation_data=(X_val, y_val),\n        callbacks=callbacks,\n        class_weight=class_weight_dict,\n        verbose=1\n    )\n    \n    # Load best weights and evaluate\n    model.load_weights(model_path)\n    val_pred = model.predict(X_val, verbose=0)\n    val_pred_classes = np.argmax(val_pred, axis=1)\n    val_true_classes = np.argmax(y_val, axis=1)\n    \n    # Print evaluation metrics\n    print(f\"\\nFold {fold_idx + 1} Results:\")\n    print(f\"Validation Accuracy: {np.mean(val_pred_classes == val_true_classes):.4f}\")\n    \n    return history, model\n\ndef main():\n    \"\"\"Main training pipeline\"\"\"\n    print(\"🚀 Starting BFRB Classification Training Pipeline\")\n    \n    # Load and prepare data\n    train_df, train_demographics = load_and_preprocess_data()\n    data_dict = prepare_dataset(train_df, train_demographics)\n    \n    # Save preprocessing artifacts\n    print(\"\\nSaving preprocessing artifacts...\")\n    np.save(EXPORT_DIR / \"feature_cols.npy\", data_dict['feature_cols'])\n    np.save(EXPORT_DIR / \"sequence_maxlen.npy\", data_dict['pad_len'])\n    np.save(EXPORT_DIR / \"gesture_classes.npy\", data_dict['gesture_classes'])\n    joblib.dump(data_dict['scaler'], EXPORT_DIR / \"scaler.pkl\")\n    \n    # Setup cross-validation\n    sgkf = StratifiedGroupKFold(n_splits=N_FOLDS, shuffle=True, random_state=42)\n    \n    # Train models\n    fold_results = []\n    \n    for fold_idx, (train_idx, val_idx) in enumerate(\n        sgkf.split(data_dict['X'], data_dict['y_encoded'], data_dict['subjects'])\n    ):\n        # Train first model type (D-111)\n        history1, model1 = train_model(\n            data_dict, fold_idx, train_idx, val_idx, \"D-111\"\n        )\n        \n        # Train second model type (v0629) \n        history2, model2 = train_model(\n            data_dict, fold_idx, train_idx, val_idx, \"v0629\"\n        )\n        \n        fold_results.append({\n            'fold': fold_idx,\n            'history1': history1,\n            'history2': history2\n        })\n    \n    print(f\"\\n🎉 Training completed! Models saved in {EXPORT_DIR}\")\n    print(\"Preprocessing artifacts saved:\")\n    print(\"  - feature_cols.npy\")\n    print(\"  - sequence_maxlen.npy\") \n    print(\"  - gesture_classes.npy\")\n    print(\"  - scaler.pkl\")\n    print(f\"  - {N_FOLDS}x2 model files (.h5)\")\n\nif __name__ == \"__main__\":\n    main()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-08-15T03:39:26.652204Z","iopub.execute_input":"2025-08-15T03:39:26.652979Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}