{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.12.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":118765,"databundleVersionId":15231210,"sourceType":"competition"},{"sourceId":14604295,"sourceType":"datasetVersion","datasetId":9328538}],"dockerImageVersionId":31259,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"!pip install --no-index /kaggle/input/datasets/kami1976/biopython-cp312/biopython-1.86-cp312-cp312-manylinux2014_x86_64.manylinux_2_17_x86_64.manylinux_2_28_x86_64.whl","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2026-02-14T21:42:05.433446Z","iopub.execute_input":"2026-02-14T21:42:05.434492Z","iopub.status.idle":"2026-02-14T21:42:10.757613Z","shell.execute_reply.started":"2026-02-14T21:42:05.434439Z","shell.execute_reply":"2026-02-14T21:42:10.756306Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport time\nimport warnings\nwarnings.filterwarnings('ignore')\n\n# ============================================================\n# DATA LOADING\n# ============================================================\nDATA_PATH = '/kaggle/input/stanford-rna-3d-folding-2/'\n\ntrain_seqs = pd.read_csv(DATA_PATH + 'train_sequences.csv')\ntest_seqs  = pd.read_csv(DATA_PATH + 'test_sequences.csv')\ntrain_labels = pd.read_csv(DATA_PATH + 'train_labels.csv')\n\n# Also load validation data to expand template pool\nval_seqs   = pd.read_csv(DATA_PATH + 'validation_sequences.csv')\nval_labels = pd.read_csv(DATA_PATH + 'validation_labels.csv')\n\n# Merge train + validation into one pool\ntrain_seqs = pd.concat([train_seqs, val_seqs], ignore_index=True)\ntrain_labels = pd.concat([train_labels, val_labels], ignore_index=True)\n\n# Build coordinate dictionary\ndef process_labels(labels_df):\n    coords_dict = {}\n    prefixes = labels_df['ID'].str.rsplit('_', n=1).str[0]\n    for id_prefix, group in labels_df.groupby(prefixes):\n        coords_dict[id_prefix] = group.sort_values('resid')[['x_1', 'y_1', 'z_1']].values\n    return coords_dict\n\ntrain_coords_dict = process_labels(train_labels)\n\n# Setup aligner (same as original code)\nfrom Bio.Align import PairwiseAligner\n\naligner = PairwiseAligner()\naligner.mode = 'global'\naligner.match_score = 2\naligner.mismatch_score = -1.5\naligner.open_gap_score   = -8\naligner.extend_gap_score = -0.4\naligner.query_left_open_gap_score  = -8\naligner.query_left_extend_gap_score = -0.4\naligner.query_right_open_gap_score = -8\naligner.query_right_extend_gap_score = -0.4\naligner.target_left_open_gap_score = -8\naligner.target_left_extend_gap_score = -0.4\naligner.target_right_open_gap_score = -8\naligner.target_right_extend_gap_score = -0.4\n\n# ============================================================\n# PART 1: Length filter analysis (fast, no alignment needed)\n# ============================================================\nprint(\"=\" * 60)\nprint(\"PART 1: Length Filter Analysis (30% threshold)\")\nprint(\"=\" * 60)\n\nvalid_train = train_seqs[train_seqs['target_id'].isin(train_coords_dict.keys())].copy()\n\ntest_lengths  = test_seqs['sequence'].str.len().values\ntrain_lengths = valid_train['sequence'].str.len().values\n\nprint(f\"\\nTotal test sequences:  {len(test_seqs)}\")\nprint(f\"Total train sequences (train+val): {len(train_seqs)}\")\nprint(f\"Train with coords:     {len(valid_train)}\")\n\nlength_pass_counts = []\nfor t_len in test_lengths:\n    ratio = np.abs(train_lengths - t_len) / np.maximum(train_lengths, t_len)\n    n_pass = np.sum(ratio <= 0.3)\n    length_pass_counts.append(n_pass)\n\nlength_pass_counts = np.array(length_pass_counts)\n\nprint(f\"\\n--- Templates passing 30% length filter per test sequence ---\")\nprint(f\"Min:    {length_pass_counts.min()}\")\nprint(f\"Max:    {length_pass_counts.max()}\")\nprint(f\"Mean:   {length_pass_counts.mean():.1f}\")\nprint(f\"Median: {np.median(length_pass_counts):.1f}\")\nprint(f\"\\nTest seqs with 0 templates passing length filter: \"\n      f\"{np.sum(length_pass_counts == 0)} / {len(test_seqs)}\")\nprint(f\"Test seqs with < 5 templates:  {np.sum(length_pass_counts < 5)}\")\nprint(f\"Test seqs with < 10 templates: {np.sum(length_pass_counts < 10)}\")\nprint(f\"Test seqs with < 30 templates: {np.sum(length_pass_counts < 30)}\")\n\nprint(f\"\\n--- Length distributions ---\")\nprint(f\"Test  lengths: min={test_lengths.min()}, max={test_lengths.max()}, \"\n      f\"mean={test_lengths.mean():.0f}, median={np.median(test_lengths):.0f}\")\nprint(f\"Train lengths: min={train_lengths.min()}, max={train_lengths.max()}, \"\n      f\"mean={train_lengths.mean():.0f}, median={np.median(train_lengths):.0f}\")\n\n# ============================================================\n# PART 2: Actual alignment score analysis (slower)\n# ============================================================\nprint(\"\\n\" + \"=\" * 60)\nprint(\"PART 2: Best Alignment Score per Test Sequence\")\nprint(\"=\" * 60)\n\nbest_scores = []\nbest_template_ids = []\nn_candidates_per_test = []\nstart = time.time()\n\nfor idx, row in test_seqs.iterrows():\n    tid = row['target_id']\n    query_seq = row['sequence']\n    query_len = len(query_seq)\n    \n    best_score = -1.0\n    best_tid = None\n    n_cand = 0\n    \n    for _, trow in valid_train.iterrows():\n        train_tid = trow['target_id']\n        train_seq = trow['sequence']\n        train_len = len(train_seq)\n        \n        # Skip self-matches (test targets may exist in train/val)\n        if train_tid == tid:\n            continue\n        \n        if abs(train_len - query_len) / max(train_len, query_len) > 0.3:\n            continue\n        \n        n_cand += 1\n        raw_score = aligner.score(query_seq, train_seq)\n        norm_score = raw_score / (2 * min(query_len, train_len))\n        \n        if norm_score > best_score:\n            best_score = norm_score\n            best_tid = train_tid\n    \n    best_scores.append(best_score)\n    best_template_ids.append(best_tid)\n    n_candidates_per_test.append(n_cand)\n    \n    elapsed = time.time() - start\n    if idx % 5 == 0:\n        print(f\"  [{idx+1}/{len(test_seqs)}] {tid}: \"\n              f\"best_score={best_score:.4f}, candidates={n_cand} | {elapsed:.1f}s\")\n\nbest_scores = np.array(best_scores)\n\n# ============================================================\n# PART 3: Categorize match quality\n# ============================================================\nprint(\"\\n\" + \"=\" * 60)\nprint(\"PART 3: Match Quality Categories\")\nprint(\"=\" * 60)\n\nGOOD_THRESHOLD = 0.7\nMEDIOCRE_THRESHOLD = 0.4\nPOOR_THRESHOLD = 0.2\n\nn_good     = np.sum(best_scores >= GOOD_THRESHOLD)\nn_mediocre = np.sum((best_scores >= MEDIOCRE_THRESHOLD) & (best_scores < GOOD_THRESHOLD))\nn_poor     = np.sum((best_scores >= POOR_THRESHOLD) & (best_scores < MEDIOCRE_THRESHOLD))\nn_bad      = np.sum((best_scores < POOR_THRESHOLD) & (best_scores >= 0))\nn_none     = np.sum(best_scores < 0)\n\ntotal = len(best_scores)\nprint(f\"\\nScore thresholds: good >= {GOOD_THRESHOLD}, mediocre >= {MEDIOCRE_THRESHOLD}, poor >= {POOR_THRESHOLD}\")\nprint(f\"\\n{'Category':<20} {'Count':>6} {'Percent':>8}\")\nprint(\"-\" * 36)\nprint(f\"{'Good (>= 0.7)':<20} {n_good:>6} {100*n_good/total:>7.1f}%\")\nprint(f\"{'Mediocre (0.4-0.7)':<20} {n_mediocre:>6} {100*n_mediocre/total:>7.1f}%\")\nprint(f\"{'Poor (0.2-0.4)':<20} {n_poor:>6} {100*n_poor/total:>7.1f}%\")\nprint(f\"{'Bad (< 0.2)':<20} {n_bad:>6} {100*n_bad/total:>7.1f}%\")\nprint(f\"{'No match at all':<20} {n_none:>6} {100*n_none/total:>7.1f}%\")\n\nprint(f\"\\n--- Score distribution ---\")\nvalid_scores = best_scores[best_scores >= 0]\nif len(valid_scores) > 0:\n    print(f\"Min:    {valid_scores.min():.4f}\")\n    print(f\"Max:    {valid_scores.max():.4f}\")\n    print(f\"Mean:   {valid_scores.mean():.4f}\")\n    print(f\"Median: {np.median(valid_scores):.4f}\")\n    print(f\"Std:    {valid_scores.std():.4f}\")\n\n# ============================================================\n# PART 4: Per-sequence detail table\n# ============================================================\nprint(\"\\n\" + \"=\" * 60)\nprint(\"PART 4: Per-Sequence Details\")\nprint(\"=\" * 60)\n\nresults_df = pd.DataFrame({\n    'target_id': test_seqs['target_id'],\n    'seq_length': test_seqs['sequence'].str.len(),\n    'n_candidates': n_candidates_per_test,\n    'best_score': best_scores,\n    'best_template': best_template_ids,\n    'quality': pd.cut(best_scores, \n                       bins=[-np.inf, 0, POOR_THRESHOLD, MEDIOCRE_THRESHOLD, GOOD_THRESHOLD, np.inf],\n                       labels=['no_match', 'bad', 'poor', 'mediocre', 'good'])\n})\n\nprint(results_df.to_string(index=False))\n\nresults_df.to_csv('template_match_analysis.csv', index=False)\nprint(\"\\nSaved to template_match_analysis.csv\")\n\n# ============================================================\n# PART 5: Plots\n# ============================================================\nfig, axes = plt.subplots(2, 2, figsize=(14, 10))\n\nax = axes[0, 0]\nax.hist(valid_scores, bins=30, edgecolor='black', alpha=0.7, color='steelblue')\nax.axvline(GOOD_THRESHOLD, color='green', linestyle='--', label=f'Good ({GOOD_THRESHOLD})')\nax.axvline(MEDIOCRE_THRESHOLD, color='orange', linestyle='--', label=f'Mediocre ({MEDIOCRE_THRESHOLD})')\nax.axvline(POOR_THRESHOLD, color='red', linestyle='--', label=f'Poor ({POOR_THRESHOLD})')\nax.set_xlabel('Best Normalized Alignment Score')\nax.set_ylabel('Count')\nax.set_title('Distribution of Best Template Match Scores')\nax.legend()\n\nax = axes[0, 1]\nax.scatter(results_df['seq_length'], results_df['best_score'], \n           alpha=0.6, s=40, c='steelblue', edgecolor='black', linewidth=0.5)\nax.axhline(GOOD_THRESHOLD, color='green', linestyle='--', alpha=0.5)\nax.axhline(MEDIOCRE_THRESHOLD, color='orange', linestyle='--', alpha=0.5)\nax.set_xlabel('Test Sequence Length')\nax.set_ylabel('Best Alignment Score')\nax.set_title('Match Quality vs Sequence Length')\n\nax = axes[1, 0]\nax.hist(length_pass_counts, bins=30, edgecolor='black', alpha=0.7, color='coral')\nax.set_xlabel('Number of Training Templates Passing Length Filter')\nax.set_ylabel('Count')\nax.set_title('Template Pool Size per Test Sequence')\n\nax = axes[1, 1]\nsizes = [n_good, n_mediocre, n_poor, n_bad, n_none]\nlabels_pie = ['Good', 'Mediocre', 'Poor', 'Bad', 'No match']\ncolors = ['#2ecc71', '#f39c12', '#e74c3c', '#8e44ad', '#95a5a6']\nnonzero = [(s, l, c) for s, l, c in zip(sizes, labels_pie, colors) if s > 0]\nif nonzero:\n    ax.pie([x[0] for x in nonzero], labels=[x[1] for x in nonzero],\n           colors=[x[2] for x in nonzero], autopct='%1.1f%%', startangle=90)\nax.set_title('Match Quality Distribution')\n\nplt.tight_layout()\nplt.savefig('template_match_analysis.png', dpi=150, bbox_inches='tight')\nplt.show()\nprint(\"\\nPlot saved to template_match_analysis.png\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-02-14T21:52:42.214373Z","iopub.execute_input":"2026-02-14T21:52:42.215524Z","iopub.status.idle":"2026-02-14T21:55:45.504906Z","shell.execute_reply.started":"2026-02-14T21:52:42.215427Z","shell.execute_reply":"2026-02-14T21:55:45.503394Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}