{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":87793,"databundleVersionId":11403143,"sourceType":"competition"}],"dockerImageVersionId":30918,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# 🧬 Stanford RNA 3D Folding Competition: Structure Prediction Masterclass 🧠\n\n<a href=\"https://www.kaggle.com/competitions/stanford-ribonanza-rna-folding\" target=\"_blank\">\n  <img src=\"https://img.shields.io/badge/Kaggle-Competition-blue?style=for-the-badge&logo=kaggle\" alt=\"Kaggle Competition\">\n</a>\n\n🌟 **Welcome to the Ultimate RNA 3D Structure Prediction Notebook!** 🌟\n\n*Unlock the secrets of RNA folding through graph deep learning. This comprehensive guide combines biological insights with advanced AI techniques to tackle one of molecular biology's most challenging problems.*\n\n---\n\n## 🚀 Notebook Roadmap: From Sequence to Structure\n\n<div style=\"padding: 15px; border: 2px solid #2ecc71; border-radius: 10px; margin: 20px 0;\">\n🔍 **Section 1: Exploratory Data Analysis (EDA)**  \n   - 📌 3D coordinate distribution analysis  \n   - 📌 Sequence pattern visualization  \n   - 📌 Structural thermodynamics insights\n</div>\n\n<div style=\"padding: 15px; border: 2px solid #3498db; border-radius: 10px; margin: 20px 0;\">\n🔧 **Section 2: Data Preprocessing Pipeline**  \n   - 🧬 RNA sequence encoding  \n   - 🎯 Coordinate normalization & outlier detection  \n   - 🧩 Dynamic padding for variable-length sequences\n</div>\n\n<div style=\"padding: 15px; border: 2px solid #9b59b6; border-radius: 10px; margin: 20px 0;\">\n🤖 <strong>Section 3: Graph Neural Network Architecture</strong>  \n</div>\n\n<div style=\"padding: 15px; border: 2px solid #e67e22; border-radius: 10px; margin: 20px 0;\">\n⚙️ <strong>Section 4: Model Training Strategies</strong><br>  \n🎯 Loss: Weighted MAE + Structural Consistency<br>  \n⏱️ Early stopping with 3D validation<br>  \n🔄 Gradient clipping (norm=1.0)  \n</div>\n","metadata":{"_uuid":"ba75c96f-68e8-4eb7-8159-5d762f35366f","_cell_guid":"4cd7fda1-faf2-4318-8c84-a231a515b338","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"markdown","source":"## 1. Data  Loading and Exploration","metadata":{"_uuid":"4974e6d3-499d-4e64-ab6f-eb6b1fd4ba99","_cell_guid":"2e35b03c-2c0b-49fe-945d-e632cd4e7eeb","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"pip install --upgrade pip\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:06:21.553138Z","iopub.execute_input":"2025-03-19T15:06:21.553489Z","iopub.status.idle":"2025-03-19T15:06:30.295103Z","shell.execute_reply.started":"2025-03-19T15:06:21.553457Z","shell.execute_reply":"2025-03-19T15:06:30.293990Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!pip install spektral","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:20:57.759753Z","iopub.execute_input":"2025-03-19T17:20:57.760043Z","iopub.status.idle":"2025-03-19T17:21:02.576147Z","shell.execute_reply.started":"2025-03-19T17:20:57.760021Z","shell.execute_reply":"2025-03-19T17:21:02.575145Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"pip install --upgrade spektral\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:21:02.577135Z","iopub.execute_input":"2025-03-19T17:21:02.577474Z","iopub.status.idle":"2025-03-19T17:21:06.180305Z","shell.execute_reply.started":"2025-03-19T17:21:02.577442Z","shell.execute_reply":"2025-03-19T17:21:06.179237Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import os\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\n# TensorFlow/Keras for deep learning model\nimport tensorflow as tf\nfrom tensorflow.keras.models import Model\nfrom tensorflow.keras.layers import Input, Embedding, Conv1D, BatchNormalization, Dropout\nfrom tensorflow.keras.preprocessing.sequence import pad_sequences\nfrom tensorflow.keras.callbacks import EarlyStopping\n#For GCN MODEL\nfrom spektral.layers import GCNConv\n\n# Set random seed for reproducibility\nnp.random.seed(42)\ntf.random.set_seed(42)","metadata":{"_uuid":"1c8a6414-5cce-4cf9-82dc-74c2b1937448","_cell_guid":"95be3c1b-d91d-4ff2-aeb1-16b50fa1936e","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T17:22:13.333231Z","iopub.execute_input":"2025-03-19T17:22:13.333524Z","iopub.status.idle":"2025-03-19T17:22:25.885899Z","shell.execute_reply.started":"2025-03-19T17:22:13.333499Z","shell.execute_reply":"2025-03-19T17:22:25.885215Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"\nWe load the CSV files provided in the competition:\n- `train_sequences.csv`\n- `train_labels.csv`\n- `validation_sequences.csv` & `validation_labels.csv`\n- `test_sequences.csv`\n- `sample_submission.csv`\n\n","metadata":{"_uuid":"7d1c7eed-82d4-46ea-859c-d7a7490d9852","_cell_guid":"08036a25-432d-4eca-bbfe-961e87c3f83a","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"TRAIN_SEQ_PATH = '/kaggle/input/stanford-rna-3d-folding/train_sequences.csv'\nTRAIN_LABELS_PATH = '/kaggle/input/stanford-rna-3d-folding/train_labels.csv'\nVALID_SEQ_PATH = '/kaggle/input/stanford-rna-3d-folding/validation_sequences.csv'\nVALID_LABELS_PATH = '/kaggle/input/stanford-rna-3d-folding/validation_labels.csv'\nTEST_SEQ_PATH  = '/kaggle/input/stanford-rna-3d-folding/test_sequences.csv'\nSAMPLE_SUB_PATH = '/kaggle/input/stanford-rna-3d-folding/sample_submission.csv'\n\ntrain_sequences = pd.read_csv(TRAIN_SEQ_PATH)\ntrain_labels = pd.read_csv(TRAIN_LABELS_PATH)\nvalid_sequences = pd.read_csv(VALID_SEQ_PATH)\nvalid_labels = pd.read_csv(VALID_LABELS_PATH)\ntest_sequences = pd.read_csv(TEST_SEQ_PATH)\nsample_submission = pd.read_csv(SAMPLE_SUB_PATH)\n\ntrain_labels.fillna(0, inplace=True)\nvalid_labels.fillna(0, inplace=True)\n\nprint(\"Train Sequences Shape:\", train_sequences.shape)\nprint(\"Train Labels Shape:\", train_labels.shape)\nprint(\"Validation Sequences Shape:\", valid_sequences.shape)\nprint(\"Validation Labels Shape:\", valid_labels.shape)\nprint(\"Test Sequences Shape:\", test_sequences.shape)\n\nprint(\"\\nTrain Sequences Head:\")\nprint(train_sequences.head())\nprint(\"\\nTrain Labels Head:\")\nprint(train_labels.head())","metadata":{"_uuid":"bd1b8c77-f907-4a3b-af02-5141a285d2ac","_cell_guid":"9f130e5f-7f0d-424f-8634-f304c96c9437","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T17:22:25.886987Z","iopub.execute_input":"2025-03-19T17:22:25.887593Z","iopub.status.idle":"2025-03-19T17:22:26.320583Z","shell.execute_reply.started":"2025-03-19T17:22:25.887567Z","shell.execute_reply":"2025-03-19T17:22:26.319870Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**histogram of sequence lengths across train, validation, and test sets**\nThe model pads sequences to a maximum length (4298 in the training set), but understanding the distribution can inform padding strategies or model architecture (e.g., handling variable lengths better).","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(10, 6))\nplt.hist(train_sequences['sequence'].str.len(), bins=50, alpha=0.7, label='Train')\nplt.hist(valid_sequences['sequence'].str.len(), bins=50, alpha=0.7, label='Validation')\nplt.hist(test_sequences['sequence'].str.len(), bins=50, alpha=0.7, label='Test')\nplt.xlabel('Sequence Length')\nplt.ylabel('Frequency')\nplt.title('Distribution of RNA Sequence Lengths')\nplt.legend()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:22:26.321990Z","iopub.execute_input":"2025-03-19T17:22:26.322204Z","iopub.status.idle":"2025-03-19T17:22:26.791320Z","shell.execute_reply.started":"2025-03-19T17:22:26.322186Z","shell.execute_reply":"2025-03-19T17:22:26.790508Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Calculating and visualizing the proportion of each nucleotide**\n\nThe frequency of A, C, G, and U might influence folding patterns or model bias.\n","metadata":{}},{"cell_type":"code","source":"from collections import Counter\ntrain_nucleotides = ''.join(train_sequences['sequence'])\ncounts = Counter(train_nucleotides)\nplt.bar(counts.keys(), counts.values())\nplt.title('Nucleotide Composition in Training Sequences')\nplt.xlabel('Nucleotide')\nplt.ylabel('Count')\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:22:26.792228Z","iopub.execute_input":"2025-03-19T17:22:26.792467Z","iopub.status.idle":"2025-03-19T17:22:26.950482Z","shell.execute_reply.started":"2025-03-19T17:22:26.792447Z","shell.execute_reply":"2025-03-19T17:22:26.949845Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**histograms or boxplots of coordinates and compute statistics**\n\nUnderstanding the range and spread of x, y, z coordinates can help assess data scale and whether normalization is needed.","metadata":{}},{"cell_type":"code","source":"coords = np.vstack([train_labels_dict[tid] for tid in train_labels_dict])\nfig, axes = plt.subplots(1, 3, figsize=(15, 5))\nfor i, coord in enumerate(['x_1', 'y_1', 'z_1']):\n    axes[i].hist(coords[:, i], bins=50)\n    axes[i].set_title(f'Distribution of {coord}')\nplt.tight_layout()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:22:26.951261Z","iopub.execute_input":"2025-03-19T17:22:26.951549Z","iopub.status.idle":"2025-03-19T17:22:27.318109Z","shell.execute_reply.started":"2025-03-19T17:22:26.951511Z","shell.execute_reply":"2025-03-19T17:22:27.316891Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Convert to datetime and explore trends over time**","metadata":{}},{"cell_type":"code","source":"train_sequences['temporal_cutoff'] = pd.to_datetime(train_sequences['temporal_cutoff'])\nplt.figure(figsize=(10, 6))\ntrain_sequences['temporal_cutoff'].dt.year.value_counts().sort_index().plot(kind='bar')\nplt.title('Sequences by Year of Temporal Cutoff')\nplt.xlabel('Year')\nplt.ylabel('Number of Sequences')\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:22:32.502184Z","iopub.execute_input":"2025-03-19T17:22:32.502490Z","iopub.status.idle":"2025-03-19T17:22:32.854285Z","shell.execute_reply.started":"2025-03-19T17:22:32.502468Z","shell.execute_reply":"2025-03-19T17:22:32.853166Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 2. Data Preprocessing\n\n### Sequence Encoding\n\nWe map each nucleotide to an integer:\n- A: 1, C: 2, G: 3, U: 4  \nUnknown characters are mapped to 0.","metadata":{"_uuid":"7b3ae69d-4427-4c9f-931e-9ea1176f5886","_cell_guid":"f84ddf2c-b665-420f-a7ec-43cfe711c0a4","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"nucleotide_map = {'A': 1, 'C': 2, 'G': 3, 'U': 4}\n\ndef encode_sequence(seq):\n    \"\"\"Encodes a RNA sequence into a list of integers based on nucleotide_map.\"\"\"\n    return [nucleotide_map.get(ch, 0) for ch in seq]\n\n# Apply encoding to all sequence files\ntrain_sequences['encoded'] = train_sequences['sequence'].apply(encode_sequence)\nvalid_sequences['encoded'] = valid_sequences['sequence'].apply(encode_sequence)\ntest_sequences['encoded'] = test_sequences['sequence'].apply(encode_sequence)","metadata":{"_uuid":"2f38e291-5974-4c43-93a0-fb024c343ba0","_cell_guid":"3c95facf-6580-4e0d-bf92-9cfb9539a551","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T17:22:36.127524Z","iopub.execute_input":"2025-03-19T17:22:36.127864Z","iopub.status.idle":"2025-03-19T17:22:36.148678Z","shell.execute_reply.started":"2025-03-19T17:22:36.127839Z","shell.execute_reply":"2025-03-19T17:22:36.147733Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"** Processing Label Data**\n\nEach row in the labels CSV is for one residue, with an `ID` formatted as `target_id_resid`.\nWe group rows by `target_id` and sort by residue number.\nHere, we use the first structure (x_1, y_1, z_1) as our target coordinates.","metadata":{"_uuid":"fbc2c63d-cb94-493c-9d4f-e0cf078d07c1","_cell_guid":"c3189e08-66bf-4548-8f33-3422b765f88c","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def process_labels(labels_df):\n    \"\"\"\n    Processes a labels DataFrame by grouping rows by target_id.\n    Returns a dictionary mapping target_id to an array of coordinates (seq_len, 3).\n    \"\"\"\n    label_dict = {}\n    for idx, row in labels_df.iterrows():\n        # Split ID into target_id and residue number (assumes format \"targetid_resid\")\n        parts = row['ID'].split('_')\n        target_id = \"_\".join(parts[:-1])\n        resid = int(parts[-1])\n        # Extract the coordinates; they should be numeric (missing values already set to 0)\n        coord = np.array([row['x_1'], row['y_1'], row['z_1']], dtype=np.float32)\n        if target_id not in label_dict:\n            label_dict[target_id] = []\n        label_dict[target_id].append((resid, coord))\n    \n    # Sort residues by resid and stack coordinates\n    for key in label_dict:\n        sorted_coords = sorted(label_dict[key], key=lambda x: x[0])\n        coords = np.stack([c for r, c in sorted_coords])\n        label_dict[key] = coords\n    return label_dict\n\n# Process training and validation labels\ntrain_labels_dict = process_labels(train_labels)\nvalid_labels_dict = process_labels(valid_labels)","metadata":{"_uuid":"1d94dbc2-d120-4c82-9d70-711e6b4a68c1","_cell_guid":"afbd9f74-c685-4e0a-9262-005656f79ca9","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T17:22:36.695418Z","iopub.execute_input":"2025-03-19T17:22:36.695728Z","iopub.status.idle":"2025-03-19T17:22:43.764172Z","shell.execute_reply.started":"2025-03-19T17:22:36.695702Z","shell.execute_reply":"2025-03-19T17:22:43.763457Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom scipy.sparse import coo_matrix\n\nnucleotide_map = {'A': 1, 'C': 2, 'G': 3, 'U': 4}\n\ndef encode_sequence(seq):\n    return np.array([nucleotide_map.get(ch, 0) for ch in seq], dtype=np.int8)\n\ntrain_sequences['encoded'] = train_sequences['sequence'].apply(encode_sequence)\nvalid_sequences['encoded'] = valid_sequences['sequence'].apply(encode_sequence)\ntest_sequences['encoded'] = test_sequences['sequence'].apply(encode_sequence)\n\ndef process_labels(labels_df):\n    label_dict = {}\n    for idx, row in labels_df.iterrows():\n        parts = row['ID'].split('_')\n        target_id = \"_\".join(parts[:-1])\n        resid = int(parts[-1])\n        coord = np.array([row['x_1'], row['y_1'], row['z_1']], dtype=np.float32)\n        if target_id not in label_dict:\n            label_dict[target_id] = []\n        label_dict[target_id].append((resid, coord))\n    for key in label_dict:\n        sorted_coords = sorted(label_dict[key], key=lambda x: x[0])\n        coords = np.stack([c for r, c in sorted_coords])\n        label_dict[key] = coords\n    return label_dict\n\ntrain_labels_dict = process_labels(train_labels)\nvalid_labels_dict = process_labels(valid_labels)\n\ndef create_graph_dataset(sequences_df, labels_dict, max_samples=None, max_len=None):\n    X_list, A_list, y_list, target_ids = [], [], [], []\n    for idx, row in sequences_df.iterrows():\n        tid = row['target_id']\n        if tid in labels_dict:\n            seq = row['encoded']\n            seq_len = min(len(seq), max_len) if max_len else len(seq)  # Cap sequence length\n            X = seq[:seq_len].reshape(-1, 1)  # Truncate if needed\n            edges = [(i, i+1) for i in range(seq_len-1)] + [(i+1, i) for i in range(seq_len-1)]\n            A = np.array(edges, dtype=np.int16).T  # [2, n_edges]\n            y = labels_dict[tid][:seq_len]  # Match labels to truncated sequence\n            X_list.append(X)\n            A_list.append(A)\n            y_list.append(y)\n            target_ids.append(tid)\n            if max_samples and len(X_list) >= max_samples:\n                break\n    return X_list, A_list, y_list, target_ids\n\n# Set max_len globally and pass it to ensure consistency\nmax_len = max(min(len(seq), 122) for seq in train_sequences['encoded'])  # Use train max\nprint(\"Maximum sequence length (train):\", max_len)\n\nX_train, A_train, y_train, train_ids = create_graph_dataset(train_sequences, train_labels_dict, max_samples=100, max_len=max_len)\nX_valid, A_valid, y_valid, valid_ids = create_graph_dataset(valid_sequences, valid_labels_dict, max_samples=100, max_len=max_len)\n\ndef pad_graph(X, A, y, max_len):\n    seq_len = X.shape[0]\n    if seq_len < max_len:\n        X_pad = np.pad(X, ((0, max_len - seq_len), (0, 0)), mode='constant', constant_values=0)\n        A_pad = A.copy()  # Edge list stays as-is\n        y_pad = np.pad(y, ((0, max_len - seq_len), (0, 0)), mode='constant', constant_values=0)\n    else:\n        X_pad, A_pad, y_pad = X, A, y\n    return X_pad, A_pad, y_pad\n\nX_train_pad, A_train_pad, y_train_pad = [], [], []\nfor x, a, y in zip(X_train, A_train, y_train):\n    x_p, a_p, y_p = pad_graph(x, a, y, max_len)\n    X_train_pad.append(x_p)\n    A_train_pad.append(a_p)\n    y_train_pad.append(y_p)\n\nX_valid_pad, A_valid_pad, y_valid_pad = [], [], []\nfor x, a, y in zip(X_valid, A_valid, y_valid):\n    x_p, a_p, y_p = pad_graph(x, a, y, max_len)\n    X_valid_pad.append(x_p)\n    A_valid_pad.append(a_p)\n    y_valid_pad.append(y_p)\n\nX_train_pad = np.array(X_train_pad, dtype=np.int8)\ny_train_pad = np.array(y_train_pad, dtype=np.float32)\nX_valid_pad = np.array(X_valid_pad, dtype=np.int8)\ny_valid_pad = np.array(y_valid_pad, dtype=np.float32)\n\ndef edges_to_sparse_adj(edges, n_nodes):\n    row, col = edges[0], edges[1]\n    data = np.ones(len(row), dtype=np.float32)\n    return coo_matrix((data, (row, col)), shape=(n_nodes, n_nodes)).toarray()\n\nA_train_adj = [edges_to_sparse_adj(a, max_len) for a in A_train_pad]\nA_valid_adj = [edges_to_sparse_adj(a, max_len) for a in A_valid_pad]\nA_train_adj = np.array(A_train_adj, dtype=np.float32)\nA_valid_adj = np.array(A_valid_adj, dtype=np.float32)\n\nprint(\"X_train_pad shape:\", X_train_pad.shape)\nprint(\"A_train_adj shape:\", A_train_adj.shape)\nprint(\"y_train_pad shape:\", y_train_pad.shape)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:22:43.765255Z","iopub.execute_input":"2025-03-19T17:22:43.765493Z","iopub.status.idle":"2025-03-19T17:22:50.786704Z","shell.execute_reply.started":"2025-03-19T17:22:43.765473Z","shell.execute_reply":"2025-03-19T17:22:50.785748Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**Creating Datasets and Padding**\n\nWe match each target sequence with its corresponding coordinate labels.\nThen we pad sequences and coordinate arrays to a uniform length.\n\nPadded positions in coordinates are set to 0.","metadata":{"_uuid":"acfc3ec8-97ac-44cf-a2dc-d1a721b1c8d8","_cell_guid":"1508c42f-9b88-4934-be7d-c75fcc36195c","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"def create_dataset(sequences_df, labels_dict):\n    \"\"\"\n    Creates a dataset from a sequences DataFrame and a labels dictionary.\n    Returns:\n        X: list of encoded sequences,\n        y: list of coordinate arrays,\n        target_ids: list of target ids.\n    \"\"\"\n    X, y, target_ids = [], [], []\n    for idx, row in sequences_df.iterrows():\n        tid = row['target_id']\n        if tid in labels_dict:\n            X.append(row['encoded'])\n            y.append(labels_dict[tid])\n            target_ids.append(tid)\n    return X, y, target_ids\n\n# Create training and validation datasets\nX_train, y_train, train_ids = create_dataset(train_sequences, train_labels_dict)\nX_valid, y_valid, valid_ids = create_dataset(valid_sequences, valid_labels_dict)\n\n# Determine maximum sequence length from training set\nmax_len = max(len(seq) for seq in X_train)\nprint(\"Maximum sequence length (train):\", max_len)\n\n# Pad the sequences (padding value = 0)\nX_train_pad = pad_sequences(X_train, maxlen=max_len, padding='post', value=0)\nX_valid_pad = pad_sequences(X_valid, maxlen=max_len, padding='post', value=0)\n\n# Function to pad coordinate arrays\ndef pad_coordinates(coord_array, max_len):\n    L = coord_array.shape[0]\n    if L < max_len:\n        pad_width = ((0, max_len - L), (0, 0))\n        return np.pad(coord_array, pad_width, mode='constant', constant_values=0)\n    else:\n        return coord_array\n\n# Pad coordinate arrays\ny_train_pad = np.array([pad_coordinates(arr, max_len) for arr in y_train])\ny_valid_pad = np.array([pad_coordinates(arr, max_len) for arr in y_valid])\n\n# Check for any NaN values in the targets\nprint(\"Any NaN in y_train_pad?\", np.isnan(y_train_pad).any())\nprint(\"Any NaN in y_valid_pad?\", np.isnan(y_valid_pad).any())\n\nprint(\"X_train_pad shape:\", X_train_pad.shape)\nprint(\"y_train_pad shape:\", y_train_pad.shape)","metadata":{"_uuid":"f00ec075-2695-461e-a423-b8b251c43acd","_cell_guid":"be76be8c-682f-418d-94d0-4f6bb519b344","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T17:22:50.788393Z","iopub.execute_input":"2025-03-19T17:22:50.788705Z","iopub.status.idle":"2025-03-19T17:22:50.935161Z","shell.execute_reply.started":"2025-03-19T17:22:50.788680Z","shell.execute_reply":"2025-03-19T17:22:50.934265Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Install dependencies (uncomment if running in a fresh environment like Kaggle)\n# !pip install torch transformers tensorflow spektral biopython viennarna\n\nimport numpy as np\nimport pandas as pd\nimport tensorflow as tf\nimport torch\nfrom transformers import AutoTokenizer, AutoModel\nimport RNA  # ViennaRNA\nfrom spektral.layers import GCNConv\nfrom tensorflow.keras import layers, Model\n\n# Set random seeds for reproducibility\nnp.random.seed(42)\ntf.random.set_seed(42)\ntorch.manual_seed(42)\n\n\n# Load data\ntrain_sequences = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/train_sequences.csv')\ntrain_labels = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/train_labels.csv')\ntest_sequences = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/test_sequences.csv')\n\n# Fill NaN values in labels\ntrain_labels.fillna(0, inplace=True)\n\n# Secondary structure prediction\ndef get_secondary_structure(sequence):\n    ss, _ = RNA.fold(sequence)\n    return ss\n\ndef build_adjacency_matrix(sequence, ss):\n    n = len(sequence)\n    adj = np.zeros((n, n))\n    for i in range(n-1):\n        adj[i, i+1] = 1\n        adj[i+1, i] = 1\n    stack = []\n    for i, char in enumerate(ss):\n        if char == '(':\n            stack.append(i)\n        elif char == ')':\n            j = stack.pop()\n            adj[i, j] = 1\n            adj[j, i] = 1\n    return adj\n\n# Apply secondary structure and adjacency matrix\ntrain_sequences['ss'] = train_sequences['sequence'].apply(get_secondary_structure)\ntrain_sequences['adj'] = train_sequences.apply(lambda row: build_adjacency_matrix(row['sequence'], row['ss']), axis=1)\ntest_sequences['ss'] = test_sequences['sequence'].apply(get_secondary_structure)\ntest_sequences['adj'] = test_sequences.apply(lambda row: build_adjacency_matrix(row['sequence'], row['ss']), axis=1)\n\n# Cap max_len to avoid memory issues\nmax_len = min(max(train_sequences['sequence'].str.len().max(), test_sequences['sequence'].str.len().max()), 500)\nprint(f\"Capped max_len: {max_len}\")\n\n\n# Load pre-trained ESM2 model and tokenizer\ntokenizer = AutoTokenizer.from_pretrained(\"facebook/esm2_t12_35M_UR50D\")\npretrained_model = AutoModel.from_pretrained(\"facebook/esm2_t12_35M_UR50D\").to('cuda' if torch.cuda.is_available() else 'cpu')\n\n# Batch processing for embeddings\ndef tokenize_and_embed_batch(sequences, max_len=max_len, batch_size=32):\n    embeddings = []\n    for i in range(0, len(sequences), batch_size):\n        batch = sequences[i:i + batch_size]\n        inputs = tokenizer(batch, padding='max_length', truncation=True, max_length=max_len, return_tensors=\"pt\").to(pretrained_model.device)\n        with torch.no_grad():\n            outputs = pretrained_model(**inputs)\n        embeddings.append(outputs.last_hidden_state.cpu().numpy())\n        torch.cuda.empty_cache()\n    return np.concatenate(embeddings, axis=0)\n\n# Get embeddings\nX_train_embeddings = tokenize_and_embed_batch(train_sequences['sequence'].tolist())\nX_test_embeddings = tokenize_and_embed_batch(test_sequences['sequence'].tolist())\n\n# Prepare adjacency matrices\ndef pad_or_truncate_adj(adj, max_len):\n    n = adj.shape[0]\n    if n > max_len:\n        return adj[:max_len, :max_len]\n    return np.pad(adj, ((0, max_len - n), (0, max_len - n)), 'constant')\n\nA_train = np.array([pad_or_truncate_adj(adj, max_len) for adj in train_sequences['adj']], dtype=np.float32)\nA_test = np.array([pad_or_truncate_adj(adj, max_len) for adj in test_sequences['adj']], dtype=np.float32)\n\n# Prepare labels\ndef pad_or_truncate_coords(coords, max_len):\n    n = coords.shape[0]\n    if n > max_len:\n        return coords[:max_len]\n    return np.pad(coords, ((0, max_len - n), (0, 0)), 'constant')\n\ny_train = []\nfor target_id in train_sequences['target_id']:\n    target_labels = train_labels[train_labels['ID'].str.startswith(target_id)]\n    coords = target_labels[['x_1', 'y_1', 'z_1']].values\n    y_train.append(pad_or_truncate_coords(coords, max_len))\ny_train = np.array(y_train, dtype=np.float32)\n\n# Verify shapes\nprint(f\"X_train_embeddings shape: {X_train_embeddings.shape}\")\nprint(f\"A_train shape: {A_train.shape}\")\nprint(f\"y_train shape: {y_train.shape}\")\n\n\n# Define inputs with explicit float32\nembedding_input = layers.Input(shape=(max_len, 480), dtype=tf.float32, name=\"embeddings\")\nadj_input = layers.Input(shape=(max_len, max_len), dtype=tf.float32, name=\"adj\")\nmask_input = layers.Input(shape=(max_len,), dtype=tf.float32, name=\"mask\")  # No default_value\n\n# GNN layers with explicit mask\ngcn1 = GCNConv(128, activation='relu', use_bias=True)([embedding_input, adj_input], mask=mask_input)\ngcn2 = GCNConv(64, activation='relu', use_bias=True)([gcn1, adj_input], mask=mask_input)\n\n# Coordinate prediction\ncoords = layers.Dense(3, activation='linear', dtype=tf.float32)(gcn2)\n\n# Model\nmodel = Model(inputs=[embedding_input, adj_input, mask_input], outputs=coords)\n\n# Compile\nmodel.compile(optimizer=tf.keras.optimizers.Adam(1e-4),\n              loss='mae',\n              metrics=['mae'])\n\n# Create mask array (all ones for all sequences)\nmask_train = np.ones((len(train_sequences), max_len), dtype=np.float32)\n\n# Train\nhistory = model.fit([X_train_embeddings, A_train, mask_train], y_train,\n                    validation_split=0.2,\n                    epochs=50,\n                    batch_size=8,\n                    callbacks=[tf.keras.callbacks.EarlyStopping(patience=5)])\n\n\n# Predict in batches\ndef predict_in_batches(model, X_emb, A, mask, batch_size=8):\n    preds = []\n    for i in range(0, len(X_emb), batch_size):\n        batch_emb = X_emb[i:i + batch_size]\n        batch_adj = A[i:i + batch_size]\n        batch_mask = mask[i:i + batch_size]\n        preds.append(model.predict([batch_emb, batch_adj, batch_mask], verbose=0))\n    return np.concatenate(preds, axis=0)\n\nmask_test = np.ones((len(test_sequences), max_len), dtype=np.float32)\npredictions = predict_in_batches(model, X_test_embeddings, A_test, mask_test)\n\nsubmission_rows = []\nfor idx, row in test_sequences.iterrows():\n    target_id = row['target_id']\n    seq_len = min(len(row['sequence']), max_len)\n    pred_coords = predictions[idx][:seq_len, :]\n    for i, coords in enumerate(pred_coords):\n        submission_rows.append({\n            'ID': f\"{target_id}_{i+1}\",\n            'resname': row['sequence'][i],\n            'resid': i+1,\n            **{f\"x_{j+1}\": coords[0] for j in range(5)},\n            **{f\"y_{j+1}\": coords[1] for j in range(5)},\n            **{f\"z_{j+1}\": coords[2] for j in range(5)}\n        })\n\nsubmission_df = pd.DataFrame(submission_rows)\nsubmission_df.to_csv(\"submission.csv\", index=False)\nprint(\"Submission file saved as submission.csv\")\nprint(\"Submission DataFrame shape:\", submission_df.shape)\nprint(submission_df.head(10))\n\nmodel.summary()\nprint(\"Training history:\", history.history)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:47:01.519894Z","iopub.execute_input":"2025-03-19T17:47:01.520243Z","iopub.status.idle":"2025-03-19T17:54:16.953865Z","shell.execute_reply.started":"2025-03-19T17:47:01.520220Z","shell.execute_reply":"2025-03-19T17:54:16.952656Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 3.1 Basic CNN Model Training\n\nIn this section, we build a basic CNN-based model.\nThe model uses:\n- An Embedding layer  \n- Two Conv1D blocks (with BatchNormalization and Dropout)  \n- A final Conv1D layer (kernel size 1) to output 3 coordinates per residue","metadata":{"_uuid":"8d832ed7-a808-47b1-9bd1-4ca8aa72b337","_cell_guid":"9616835f-b774-4de4-bb7d-99fffe07dc65","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"# Define hyperparameters for the CNN model\nvocab_size = max(nucleotide_map.values()) + 1  # +1 for padding token 0\nembedding_dim = 16\nnum_filters = 64\nkernel_size = 3\ndrop_rate = 0.2\n\n# Build the CNN model\ninput_seq_cnn = Input(shape=(max_len,), name='input_seq')\nx_cnn = Embedding(input_dim=vocab_size, output_dim=embedding_dim, mask_zero=True, name='embedding')(input_seq_cnn)\n\n# First convolutional block\nx_cnn = Conv1D(filters=num_filters, kernel_size=kernel_size, padding='same', activation='relu', name='conv1')(x_cnn)\nx_cnn = BatchNormalization(name='bn1')(x_cnn)\nx_cnn = Dropout(drop_rate, name='drop1')(x_cnn)\n\n# Second convolutional block\nx_cnn = Conv1D(filters=num_filters, kernel_size=kernel_size, padding='same', activation='relu', name='conv2')(x_cnn)\nx_cnn = BatchNormalization(name='bn2')(x_cnn)\nx_cnn = Dropout(drop_rate, name='drop2')(x_cnn)\n\n# Final convolution to output 3 coordinates per residue (x, y, z)\noutput_coords_cnn = Conv1D(filters=3, kernel_size=1, padding='same', activation='linear', name='predicted_coords')(x_cnn)\n\ncnn_model = Model(inputs=input_seq_cnn, outputs=output_coords_cnn)\ncnn_model.compile(optimizer='adam', loss='mse')\n\ncnn_model.summary()","metadata":{"_uuid":"b5b77c10-dad8-4527-8f28-174c04afefe2","_cell_guid":"eb0e2b5f-347a-45b9-a475-a2a4108dd9f1","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T17:22:50.936267Z","iopub.execute_input":"2025-03-19T17:22:50.936597Z","iopub.status.idle":"2025-03-19T17:22:53.256653Z","shell.execute_reply.started":"2025-03-19T17:22:50.936566Z","shell.execute_reply":"2025-03-19T17:22:53.256005Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 3.2 Model Training\n\nWe train the CNN model using early stopping to monitor the validation loss.\nWith the NaN issues addressed in the data, training should proceed without nan losses.","metadata":{"_uuid":"93c354e3-05d7-4632-89cb-b9c843c37977","_cell_guid":"7744df17-746a-438d-82eb-3fda8874d03c","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"early_stop_cnn = EarlyStopping(monitor='val_loss', patience=5, restore_best_weights=True)\n\nhistory_cnn = cnn_model.fit(X_train_pad, y_train_pad,\n                            validation_data=(X_valid_pad, y_valid_pad),\n                            epochs=50,\n                            batch_size=16,\n                            callbacks=[early_stop_cnn],\n                            verbose=1)\n\n# Plot training and validation loss\nplt.figure(figsize=(8, 5))\nplt.plot(history_cnn.history['loss'], label='Train Loss (CNN)')\nplt.plot(history_cnn.history['val_loss'], label='Val Loss (CNN)')\nplt.xlabel(\"Epoch\")\nplt.ylabel(\"MSE Loss\")\nplt.title(\"CNN Model Training vs. Validation Loss\")\nplt.legend()\nplt.show()","metadata":{"_uuid":"8b6e7d19-9670-4f06-9995-61f156f5588e","_cell_guid":"fbd23461-6579-49ed-adba-4c414b97a74b","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T15:07:01.897574Z","iopub.status.idle":"2025-03-19T15:07:01.897991Z","shell.execute_reply":"2025-03-19T15:07:01.897808Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 4. GCN model building and training","metadata":{}},{"cell_type":"code","source":"import numpy as np\nfrom scipy.sparse import coo_matrix\nimport tensorflow as tf\nfrom tensorflow.keras.models import Model\nfrom tensorflow.keras.layers import Input, Embedding, Dense, Dropout, BatchNormalization, Layer\n\n# Preprocessing (same as before)\nmax_len = 500\nmax_samples = 200\n\nX_train, A_train, y_train, train_ids = create_graph_dataset(train_sequences, train_labels_dict, max_samples=max_samples, max_len=max_len)\nX_valid, A_valid, y_valid, valid_ids = create_graph_dataset(valid_sequences, valid_labels_dict, max_samples=max_samples, max_len=max_len)\n\nX_train_pad, A_train_pad, y_train_pad = [], [], []\nfor x, a, y in zip(X_train, A_train, y_train):\n    x_p, a_p, y_p = pad_graph(x, a, y, max_len)\n    X_train_pad.append(x_p)\n    A_train_pad.append(a_p)\n    y_train_pad.append(y_p)\n\nX_valid_pad, A_valid_pad, y_valid_pad = [], [], []\nfor x, a, y in zip(X_valid, A_valid, y_valid):\n    x_p, a_p, y_p = pad_graph(x, a, y, max_len)\n    X_valid_pad.append(x_p)\n    A_valid_pad.append(a_p)\n    y_valid_pad.append(y_p)\n\nX_train_pad = np.array(X_train_pad, dtype=np.int8)\nA_train_adj = [edges_to_sparse_adj(a, max_len) for a in A_train_pad]\nA_train_adj = np.array(A_train_adj, dtype=np.float32)\ny_train_pad = np.array(y_train_pad, dtype=np.float32)\n\nX_valid_pad = np.array(X_valid_pad, dtype=np.int8)\nA_valid_adj = [edges_to_sparse_adj(a, max_len) for a in A_valid_pad]\nA_valid_adj = np.array(A_valid_adj, dtype=np.float32)\ny_valid_pad = np.array(y_valid_pad, dtype=np.float32)\n\n# Model\nvocab_size = 5\nembedding_dim = 16\ngcn_units = 64\ndrop_rate = 0.2\n\nclass SimpleGCN(Layer):\n    def __init__(self, units, activation=None, **kwargs):\n        super(SimpleGCN, self).__init__(**kwargs)\n        self.units = units\n        self.activation = tf.keras.activations.get(activation)\n\n    def build(self, input_shape):\n        self.kernel = self.add_weight(name=\"kernel\",\n                                      shape=(input_shape[0][-1], self.units),\n                                      initializer=\"glorot_uniform\",\n                                      trainable=True)\n        self.bias = self.add_weight(name=\"bias\",\n                                    shape=(self.units,),\n                                    initializer=\"zeros\",\n                                    trainable=True)\n        super(SimpleGCN, self).build(input_shape)\n\n    def call(self, inputs):\n        x, adj = inputs\n        features = tf.matmul(x, self.kernel)\n        output = tf.matmul(adj, features)\n        output = output + self.bias\n        if self.activation is not None:\n            output = self.activation(output)\n        return output\n\ninput_nodes = Input(shape=(max_len,), name='input_nodes')  # [batch, 500]\ninput_adj = Input(shape=(max_len, max_len), name='input_adj')  # [batch, 500, 500]\n\nx = Embedding(input_dim=vocab_size, output_dim=embedding_dim)(input_nodes)\nx = SimpleGCN(gcn_units, activation='relu')([x, input_adj])\nx = BatchNormalization()(x)\nx = Dropout(drop_rate)(x)\nx = SimpleGCN(gcn_units, activation='relu')([x, input_adj])\nx = BatchNormalization()(x)\nx = Dropout(drop_rate)(x)\noutput_coords = Dense(3, activation='linear', name='predicted_coords')(x)\n\ngcn_model = Model(inputs=[input_nodes, input_adj], outputs=output_coords)\n\ndef masked_mse(y_true, y_pred):\n    mask = tf.cast(tf.reduce_any(tf.not_equal(y_true, 0), axis=-1), tf.float32)\n    mse = tf.reduce_mean(tf.square(y_true - y_pred), axis=-1)  # Manual MSE\n    return tf.reduce_sum(mse * mask) / tf.reduce_sum(mask)\n\ngcn_model.compile(optimizer='adam', loss=masked_mse)\ngcn_model.summary()\n\n# Squeeze X to match model input\nX_train_pad = X_train_pad.squeeze(axis=-1)  # [200, 500]\nX_valid_pad = X_valid_pad.squeeze(axis=-1)  # [200, 500]\n\n# Verify shapes\nprint(\"X_train_pad shape:\", X_train_pad.shape)\nprint(\"A_train_adj shape:\", A_train_adj.shape)\nprint(\"y_train_pad shape:\", y_train_pad.shape)\n\n# Train\nhistory=gcn_model.fit([X_train_pad, A_train_adj], y_train_pad,\n              validation_data=([X_valid_pad, A_valid_adj], y_valid_pad),\n              epochs=25, batch_size=32, verbose=1)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:07:01.899023Z","iopub.status.idle":"2025-03-19T15:07:01.899439Z","shell.execute_reply":"2025-03-19T15:07:01.899266Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import numpy as np\nfrom scipy.sparse import coo_matrix\nimport tensorflow as tf\nfrom tensorflow.keras.models import Model\nfrom tensorflow.keras.layers import Input, Embedding, Dense, Dropout, BatchNormalization, Layer\nfrom tensorflow.keras import regularizers\n\n# Preprocessing (same)\nmax_len = 500\nmax_samples = 200\n\nX_train, A_train, y_train, train_ids = create_graph_dataset(train_sequences, train_labels_dict, max_samples=max_samples, max_len=max_len)\nX_valid, A_valid, y_valid, valid_ids = create_graph_dataset(valid_sequences, valid_labels_dict, max_samples=max_samples, max_len=max_len)\n\nX_train_pad, A_train_pad, y_train_pad = [], [], []\nfor x, a, y in zip(X_train, A_train, y_train):\n    x_p, a_p, y_p = pad_graph(x, a, y, max_len)\n    X_train_pad.append(x_p)\n    A_train_pad.append(a_p)\n    y_train_pad.append(y_p)\n\nX_valid_pad, A_valid_pad, y_valid_pad = [], [], []\nfor x, a, y in zip(X_valid, A_valid, y_valid):\n    x_p, a_p, y_p = pad_graph(x, a, y, max_len)\n    X_valid_pad.append(x_p)\n    A_valid_pad.append(a_p)\n    y_valid_pad.append(y_p)\n\nX_train_pad = np.array(X_train_pad, dtype=np.int8).squeeze(axis=-1)  # [200, 500]\nA_train_adj = np.array([edges_to_sparse_adj(a, max_len) for a in A_train_pad], dtype=np.float32)\ny_train_pad = np.array(y_train_pad, dtype=np.float32)\n\nX_valid_pad = np.array(X_valid_pad, dtype=np.int8).squeeze(axis=-1)  # [200, 500]\nA_valid_adj = np.array([edges_to_sparse_adj(a, max_len) for a in A_valid_pad], dtype=np.float32)\ny_valid_pad = np.array(y_valid_pad, dtype=np.float32)\n\n# Model with regularization\nvocab_size = 5\nembedding_dim = 16\ngcn_units = 64\ndrop_rate = 0.2\n\nclass SimpleGCN(Layer):\n    def __init__(self, units, activation=None, **kwargs):\n        super(SimpleGCN, self).__init__(**kwargs)\n        self.units = units\n        self.activation = tf.keras.activations.get(activation)\n\n    def build(self, input_shape):\n        self.kernel = self.add_weight(name=\"kernel\",\n                                      shape=(input_shape[0][-1], self.units),\n                                      initializer=\"glorot_uniform\",\n                                      regularizer=regularizers.l2(0.01),  # Add L2 reg\n                                      trainable=True)\n        self.bias = self.add_weight(name=\"bias\",\n                                    shape=(self.units,),\n                                    initializer=\"zeros\",\n                                    trainable=True)\n        super(SimpleGCN, self).build(input_shape)\n\n    def call(self, inputs):\n        x, adj = inputs\n        features = tf.matmul(x, self.kernel)\n        output = tf.matmul(adj, features)\n        output = output + self.bias\n        if self.activation is not None:\n            output = self.activation(output)\n        return output\n\ninput_nodes = Input(shape=(max_len,), name='input_nodes')\ninput_adj = Input(shape=(max_len, max_len), name='input_adj')\n\nx = Embedding(input_dim=vocab_size, output_dim=embedding_dim)(input_nodes)\nx = SimpleGCN(gcn_units, activation='relu')([x, input_adj])\nx = BatchNormalization()(x)\nx = Dropout(drop_rate)(x)\nx = SimpleGCN(gcn_units, activation='relu')([x, input_adj])\nx = BatchNormalization()(x)\nx = Dropout(drop_rate)(x)\noutput_coords = Dense(3, activation='linear', name='predicted_coords')(x)\n\ngcn_model = Model(inputs=[input_nodes, input_adj], outputs=output_coords)\n\ndef masked_mse(y_true, y_pred):\n    mask = tf.cast(tf.reduce_any(tf.not_equal(y_true, 0), axis=-1), tf.float32)\n    mse = tf.reduce_mean(tf.square(y_true - y_pred), axis=-1)\n    mask_sum = tf.reduce_sum(mask)\n    return tf.cond(mask_sum > 0, \n                   lambda: tf.reduce_sum(mse * mask) / mask_sum, \n                   lambda: tf.constant(0.0))  # Avoid divide-by-zero\n\ngcn_model.compile(optimizer='adam', loss=masked_mse)\ngcn_model.summary()\n\n# Debug shapes and values\nprint(\"X_train_pad shape:\", X_train_pad.shape)\nprint(\"A_train_adj shape:\", A_train_adj.shape)\nprint(\"y_train_pad shape:\", y_train_pad.shape)\nprint(\"X_valid_pad shape:\", X_valid_pad.shape)\nprint(\"A_valid_adj shape:\", A_valid_adj.shape)\nprint(\"y_valid_pad shape:\", y_valid_pad.shape)\nprint(\"y_valid_pad sample (first 5):\", y_valid_pad[0, :5])\nprint(\"y_valid_pad max/min:\", y_valid_pad.max(), y_valid_pad.min())\n\n# Train\nhistory1=gcn_model.fit([X_train_pad, A_train_adj], y_train_pad,\n              validation_data=([X_valid_pad, A_valid_adj], y_valid_pad),\n              epochs=25, batch_size=32, verbose=1)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:07:01.900279Z","iopub.status.idle":"2025-03-19T15:07:01.900664Z","shell.execute_reply":"2025-03-19T15:07:01.900494Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"plt.plot(history1.history['loss'], label='Train Loss')\nplt.plot(history1.history['val_loss'], label='Validation Loss')\nplt.xlabel('Epoch')\nplt.ylabel('Loss')\nplt.legend()\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:07:01.901827Z","iopub.status.idle":"2025-03-19T15:07:01.902227Z","shell.execute_reply":"2025-03-19T15:07:01.902053Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 4.2 Hybrid Transformer + Graph Neural Network (GNN) with Fine-Tuning","metadata":{}},{"cell_type":"code","source":"!pip install transformers biopython viennarna","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T17:23:45.595278Z","iopub.execute_input":"2025-03-19T17:23:45.595617Z","iopub.status.idle":"2025-03-19T17:23:51.088548Z","shell.execute_reply.started":"2025-03-19T17:23:45.595593Z","shell.execute_reply":"2025-03-19T17:23:51.087443Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom Bio import SeqIO\nimport RNA  # ViennaRNA package\n\n# Load data (as in your notebook)\ntrain_sequences = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/train_sequences.csv')\n\n# Predict secondary structure\ndef get_secondary_structure(sequence):\n    ss, _ = RNA.fold(sequence)  # Dot-bracket notation\n    return ss\n\n# Adjacency matrix from secondary structure\ndef build_adjacency_matrix(sequence, ss):\n    n = len(sequence)\n    adj = np.zeros((n, n))\n    # Covalent bonds\n    for i in range(n-1):\n        adj[i, i+1] = 1\n        adj[i+1, i] = 1\n    # Base pairs\n    stack = []\n    for i, char in enumerate(ss):\n        if char == '(':\n            stack.append(i)\n        elif char == ')':\n            j = stack.pop()\n            adj[i, j] = 1\n            adj[j, i] = 1\n    return adj\n\ntrain_sequences['ss'] = train_sequences['sequence'].apply(get_secondary_structure)\ntrain_sequences['adj'] = train_sequences.apply(lambda row: build_adjacency_matrix(row['sequence'], row['ss']), axis=1)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:09:01.335936Z","iopub.execute_input":"2025-03-19T15:09:01.336273Z","iopub.status.idle":"2025-03-19T15:16:13.516956Z","shell.execute_reply.started":"2025-03-19T15:09:01.336245Z","shell.execute_reply":"2025-03-19T15:16:13.516198Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"from transformers import AutoModel, AutoTokenizer\nfrom spektral.layers import GCNConv\nimport tensorflow as tf\nfrom tensorflow.keras import layers, Model\n\n# Load pre-trained RNA model\ntokenizer = AutoTokenizer.from_pretrained(\"facebook/esm2_t12_35M_UR50D\")  # Example, adjust to RNA model\npretrained_model = AutoModel.from_pretrained(\"facebook/esm2_t12_35M_UR50D\")\n\n# Tokenize sequences\ndef tokenize_sequences(sequences):\n    return tokenizer(sequences, padding=True, truncation=True, return_tensors=\"tf\")\n\nX_train = tokenize_sequences(train_sequences['sequence'].tolist())\n\n# Define hybrid model\nmax_len = 200  # Adjust based on your data\ninput_ids = layers.Input(shape=(max_len,), dtype=tf.int32, name=\"input_ids\")\nadj_input = layers.Input(shape=(max_len, max_len), dtype=tf.float32, name=\"adj\")\n\n# Transformer embeddings\ntransformer_outputs = pretrained_model(input_ids).last_hidden_state  # [batch, max_len, hidden_size]\n\n# GNN layers\ngcn1 = GCNConv(128, activation='relu')([transformer_outputs, adj_input])\ngcn2 = GCNConv(64, activation='relu')([gcn1, adj_input])\n\n# Coordinate prediction\ncoords = layers.Dense(3, activation='linear')(gcn2)  # [batch, max_len, 3]\n\nmodel = Model(inputs=[input_ids, adj_input], outputs=coords)\n\n# Compile\nmodel.compile(optimizer=tf.keras.optimizers.Adam(1e-4),\n              loss='mae',\n              metrics=['mae'])","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:16:27.036724Z","iopub.execute_input":"2025-03-19T15:16:27.037055Z","iopub.status.idle":"2025-03-19T15:16:45.420941Z","shell.execute_reply.started":"2025-03-19T15:16:27.037028Z","shell.execute_reply":"2025-03-19T15:16:45.419628Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Prepare adjacency matrices\nA_train = np.array(train_sequences['adj'].tolist())\ny_train = np.array(train_labels[['x', 'y', 'z']].values)  # Adjust based on your labels\n\n# Train\nhistory = model.fit([X_train['input_ids'], A_train], y_train,\n                    validation_split=0.2,\n                    epochs=50,\n                    batch_size=16,\n                    callbacks=[tf.keras.callbacks.EarlyStopping(patience=5)])","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Test data\nX_test = tokenize_sequences(test_sequences['sequence'].tolist())\nA_test = np.array(test_sequences['adj'].tolist())\npredictions = model.predict([X_test['input_ids'], A_test])\n\n# Build submission (similar to your GCN code)\nsubmission_rows = []\nfor idx, row in test_sequences.iterrows():\n    target_id = row['target_id']\n    pred_coords = predictions[idx][:len(row['sequence']), :]\n    for i, coords in enumerate(pred_coords):\n        submission_rows.append({\n            'ID': f\"{target_id}_{i+1}\",\n            'resname': row['sequence'][i],\n            'resid': i+1,\n            **{f\"x_{j+1}\": coords[0] for j in range(5)},\n            **{f\"y_{j+1}\": coords[1] for j in range(5)},\n            **{f\"z_{j+1}\": coords[2] for j in range(5)}\n        })\n\nsubmission_df = pd.DataFrame(submission_rows)\nsubmission_df.to_csv(\"submission.csv\", index=False)","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 5. Generating Predictions and Submission File\n\nFor each test sequence, we predict the 3D coordinates using our trained CNN model.\n\nThe submission requires 5 sets of coordinates per target. In this baseline, we replicate the same predicted structure 5 times.","metadata":{"_uuid":"fe88f530-9cd2-4dd5-94f2-41e2021c50f9","_cell_guid":"fe7b78f7-a78f-4cda-b597-5cbf25d98f28","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"# Prepare test data: pad sequences to same length as training\nX_test = test_sequences['encoded'].tolist()\nX_test_pad = pad_sequences(X_test, maxlen=max_len, padding='post', value=0)\n\n# Predict coordinates using the trained CNN model\npredictions = cnn_model.predict(X_test_pad)\n\n# Build submission rows. Each row corresponds to a residue from a test target.\nsubmission_rows = []\nfor idx, row in test_sequences.iterrows():\n    target_id = row['target_id']\n    # Get predicted coordinates (shape: [max_len, 3])\n    pred_coords = predictions[idx]\n    # Determine actual sequence length\n    seq_length = len(row['encoded'])\n    pred_coords = pred_coords[:seq_length, :]  # only actual residues\n    \n    # For each residue, create a row in the submission file\n    for i in range(seq_length):\n        coords = pred_coords[i, :]\n        # Replicate the same prediction 5 times for submission format\n        submission_rows.append({\n            'ID': f\"{target_id}_{i+1}\",\n            'resname': row['sequence'][i],\n            'resid': i+1,\n            **{f\"x_{j+1}\": coords[0] for j in range(5)},\n            **{f\"y_{j+1}\": coords[1] for j in range(5)},\n            **{f\"z_{j+1}\": coords[2] for j in range(5)}\n        })\n\nsubmission_df = pd.DataFrame(submission_rows)\nprint(\"Submission DataFrame shape:\", submission_df.shape)\nprint(submission_df.head(10))","metadata":{"_uuid":"d1e938fc-98b7-499c-b05b-b6edd621c9d0","_cell_guid":"fbab4072-6a39-49e3-b6f7-57e64556e9d7","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T15:07:01.903160Z","iopub.status.idle":"2025-03-19T15:07:01.903536Z","shell.execute_reply":"2025-03-19T15:07:01.903398Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 6.generating predictions for GCN model ","metadata":{}},{"cell_type":"code","source":"# Prepare test data\nX_test = test_sequences['encoded'].tolist()\nA_test = [np.zeros((len(seq), len(seq))) for seq in X_test]\nfor i, seq in enumerate(X_test):\n    for j in range(len(seq) - 1):\n        A_test[i][j, j + 1] = 1\n        A_test[i][j + 1, j] = 1\n\nX_test_pad = []\nA_test_pad = []\nfor x, a in zip(X_test, A_test):\n    x_p = np.pad(x, (0, max_len - len(x)), mode='constant', constant_values=0).reshape(-1, 1)\n    a_p = np.pad(a, ((0, max_len - len(a)), (0, max_len - len(a))), mode='constant', constant_values=0)\n    X_test_pad.append(x_p)\n    A_test_pad.append(a_p)\n\nX_test_pad = np.array(X_test_pad)\nA_test_pad = np.array(A_test_pad)\n\n# Predict coordinates\npredictions = gcn_model.predict([X_test_pad, A_test_pad])\n\n# Build submission\nsubmission_rows = []\nfor idx, row in test_sequences.iterrows():\n    target_id = row['target_id']\n    pred_coords = predictions[idx]\n    seq_length = len(row['encoded'])\n    pred_coords = pred_coords[:seq_length, :]\n    \n    for i in range(seq_length):\n        coords = pred_coords[i, :]\n        submission_rows.append({\n            'ID': f\"{target_id}_{i+1}\",\n            'resname': row['sequence'][i],\n            'resid': i+1,\n            **{f\"x_{j+1}\": coords[0] for j in range(5)},\n            **{f\"y_{j+1}\": coords[1] for j in range(5)},\n            **{f\"z_{j+1}\": coords[2] for j in range(5)}\n        })\n\nsubmission_df = pd.DataFrame(submission_rows)\nprint(\"Submission DataFrame shape:\", submission_df.shape)\nprint(submission_df.head(10))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-19T15:07:01.904314Z","iopub.status.idle":"2025-03-19T15:07:01.904633Z","shell.execute_reply":"2025-03-19T15:07:01.904471Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 6. Saving the Submission File\n\nFinally, we save the submission file as `submission.csv`.","metadata":{"_uuid":"7d9bbfd9-5f22-46bb-9a74-86c64c32623e","_cell_guid":"e122805e-51a2-411f-81f8-072ff85e4f5a","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false}}},{"cell_type":"code","source":"submission_df.to_csv(\"submission.csv\", index=False)\nprint(\"Submission file saved as submission.csv\")","metadata":{"_uuid":"2632f190-329c-4373-9632-b03b99d89e90","_cell_guid":"f1b1e66f-6ed2-4fed-b7bc-b5b9d2a306c2","trusted":true,"collapsed":false,"jupyter":{"outputs_hidden":false},"execution":{"iopub.status.busy":"2025-03-19T15:07:01.905509Z","iopub.status.idle":"2025-03-19T15:07:01.905793Z","shell.execute_reply":"2025-03-19T15:07:01.905681Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"","metadata":{}},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}