{"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":"gpu","dataSources":[{"sourceId":87793,"databundleVersionId":11403143,"sourceType":"competition"}],"dockerImageVersionId":30919,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"!pip install Bio","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:13:35.788553Z","iopub.execute_input":"2025-03-16T19:13:35.788849Z","iopub.status.idle":"2025-03-16T19:13:42.886177Z","shell.execute_reply.started":"2025-03-16T19:13:35.788827Z","shell.execute_reply":"2025-03-16T19:13:42.885037Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport os\nimport matplotlib.pyplot as plt\nfrom Bio import SeqIO\nfrom Bio.Seq import Seq\nfrom Bio.PDB import PDBParser, PDBIO, Select\nimport seaborn as sns\nfrom collections import defaultdict\nfrom sklearn.decomposition import PCA\nfrom tqdm import tqdm","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:13:45.950921Z","iopub.execute_input":"2025-03-16T19:13:45.951269Z","iopub.status.idle":"2025-03-16T19:13:47.111267Z","shell.execute_reply.started":"2025-03-16T19:13:45.951240Z","shell.execute_reply":"2025-03-16T19:13:47.110591Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# 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')\nvalidation_sequences = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/validation_sequences.csv')\nvalidation_labels = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/validation_labels.csv')\ntest_sequences = pd.read_csv('/kaggle/input/stanford-rna-3d-folding/test_sequences.csv')\n\n# Basic data exploration\nprint(f\"Train sequences: {train_sequences.shape}\")\nprint(f\"Train labels: {train_labels.shape}\")\nprint(f\"Validation sequences: {validation_sequences.shape}\")\nprint(f\"Validation labels: {validation_labels.shape}\")\nprint(f\"Test sequences: {test_sequences.shape}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:15:36.901887Z","iopub.execute_input":"2025-03-16T19:15:36.902227Z","iopub.status.idle":"2025-03-16T19:15:37.157410Z","shell.execute_reply.started":"2025-03-16T19:15:36.902201Z","shell.execute_reply":"2025-03-16T19:15:37.156794Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Check for sequence lengths\ntrain_sequences['seq_length'] = train_sequences['sequence'].str.len()\nvalidation_sequences['seq_length'] = validation_sequences['sequence'].str.len()\ntest_sequences['seq_length'] = test_sequences['sequence'].str.len()\n\nplt.figure(figsize=(10, 6))\nsns.histplot(train_sequences['seq_length'], kde=True, label='Train')\nsns.histplot(validation_sequences['seq_length'], kde=True, label='Validation')\nsns.histplot(test_sequences['seq_length'], kde=True, label='Test')\nplt.title('RNA Sequence Length Distribution')\nplt.xlabel('Sequence Length')\nplt.ylabel('Count')\nplt.legend()\nplt.savefig('sequence_length_distribution.png')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:15:46.466478Z","iopub.execute_input":"2025-03-16T19:15:46.466784Z","iopub.status.idle":"2025-03-16T19:15:47.437741Z","shell.execute_reply.started":"2025-03-16T19:15:46.466761Z","shell.execute_reply":"2025-03-16T19:15:47.437008Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Check for nucleotide distribution\ndef count_nucleotides(sequence):\n    return {\n        'A': sequence.count('A'),\n        'C': sequence.count('C'),\n        'G': sequence.count('G'),\n        'U': sequence.count('U'),\n        'Other': len(sequence) - sequence.count('A') - sequence.count('C') - \n                 sequence.count('G') - sequence.count('U')\n    }\n\ntrain_nucleotides = train_sequences['sequence'].apply(count_nucleotides).apply(pd.Series)\ntrain_nucleotides_norm = train_nucleotides.div(train_sequences['seq_length'], axis=0)\n\nplt.figure(figsize=(10, 6))\nsns.boxplot(data=train_nucleotides_norm[['A', 'C', 'G', 'U']])\nplt.title('Nucleotide Distribution in Training Data')\nplt.ylabel('Proportion')\nplt.savefig('nucleotide_distribution.png')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:15:58.415561Z","iopub.execute_input":"2025-03-16T19:15:58.415891Z","iopub.status.idle":"2025-03-16T19:15:58.925412Z","shell.execute_reply.started":"2025-03-16T19:15:58.415868Z","shell.execute_reply":"2025-03-16T19:15:58.924726Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"\n# Check for multiple conformations in validation data\nval_conformation_counts = {}\nfor col in validation_labels.columns:\n    if col.startswith('x_') and col != 'x_1':\n        val_conformation_counts[col] = (~validation_labels[col].isna()).sum()\n\nplt.figure(figsize=(10, 6))\nplt.bar(val_conformation_counts.keys(), val_conformation_counts.values())\nplt.title('Number of Residues with Multiple Conformations in Validation Data')\nplt.ylabel('Count')\nplt.savefig('multiple_conformations.png')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:16:06.619262Z","iopub.execute_input":"2025-03-16T19:16:06.619546Z","iopub.status.idle":"2025-03-16T19:16:07.095296Z","shell.execute_reply.started":"2025-03-16T19:16:06.619524Z","shell.execute_reply":"2025-03-16T19:16:07.094503Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Temporal distribution analysis\ntrain_sequences['temporal_cutoff'] = pd.to_datetime(train_sequences['temporal_cutoff'])\nplt.figure(figsize=(12, 6))\nsns.histplot(train_sequences['temporal_cutoff'], kde=True, bins=50)\nplt.title('Temporal Distribution of Training Data')\nplt.xlabel('Publication Date')\nplt.ylabel('Count')\nplt.xticks(rotation=45)\nplt.tight_layout()\nplt.savefig('temporal_distribution.png')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:16:16.146739Z","iopub.execute_input":"2025-03-16T19:16:16.147105Z","iopub.status.idle":"2025-03-16T19:16:16.697601Z","shell.execute_reply.started":"2025-03-16T19:16:16.147076Z","shell.execute_reply":"2025-03-16T19:16:16.696750Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Create a proper temporal split\ncasp15_cutoff = pd.to_datetime('2022-05-27')\ntrain_filtered = train_sequences[train_sequences['temporal_cutoff'] < casp15_cutoff]\nprint(f\"Training data before CASP15 cutoff: {train_filtered.shape}\")\n\n# Function to extract MSA information\ndef process_msa(target_id, msa_dir='/kaggle/input/stanford-rna-3d-folding/MSA/'):\n    \"\"\"Process MSA file for a given target_id and return basic statistics.\"\"\"\n    try:\n        msa_file = f\"{msa_dir}{target_id}.MSA.fasta\"\n        if not os.path.exists(msa_file):\n            return {'seq_count': 0, 'coverage': 0, 'exists': False}\n        \n        sequences = list(SeqIO.parse(msa_file, \"fasta\"))\n        if not sequences:\n            return {'seq_count': 0, 'coverage': 0, 'exists': True}\n        \n        seq_count = len(sequences)\n        \n        # Calculate coverage (proportion of non-gap positions)\n        alignments = [str(seq.seq) for seq in sequences]\n        alignment_length = len(alignments[0])\n        coverage = sum(1 for i in range(alignment_length) \n                      if any(ali[i] != '-' for ali in alignments)) / alignment_length\n        \n        return {'seq_count': seq_count, 'coverage': coverage, 'exists': True}\n    except Exception as e:\n        print(f\"Error processing MSA for {target_id}: {e}\")\n        return {'seq_count': 0, 'coverage': 0, 'exists': False, 'error': str(e)}\n\n# Sample MSA processing on a few training examples\nsample_targets = train_sequences['target_id'].sample(min(10, len(train_sequences))).tolist()\nmsa_stats = {target: process_msa(target) for target in sample_targets}\n\n\nprint(\"\\nMSA Statistics for Sample Targets:\")\nfor target, stats in msa_stats.items():\n    print(f\"{target}: {stats}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:16:33.332550Z","iopub.execute_input":"2025-03-16T19:16:33.332856Z","iopub.status.idle":"2025-03-16T19:16:33.438685Z","shell.execute_reply.started":"2025-03-16T19:16:33.332832Z","shell.execute_reply":"2025-03-16T19:16:33.437814Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Process coordinates to create a standard format\ndef extract_coordinates(df, num_structures=1):\n    \"\"\"Extract coordinates from the dataframe and return as a numpy array.\"\"\"\n    coordinates = []\n    for i in range(1, num_structures + 1):\n        if f'x_{i}' in df.columns:\n            struct_coords = df[[f'x_{i}', f'y_{i}', f'z_{i}']].values\n            if not np.isnan(struct_coords).any():\n                coordinates.append(struct_coords)\n    return np.array(coordinates)\n\n# Process a sample target to demonstrate coordinate extraction\nsample_target_id = train_sequences['target_id'].iloc[0]\nsample_target_residues = train_labels[train_labels['ID'].str.startswith(f\"{sample_target_id}_\")]\nsample_coords = extract_coordinates(sample_target_residues)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:16:43.846526Z","iopub.execute_input":"2025-03-16T19:16:43.846847Z","iopub.status.idle":"2025-03-16T19:16:43.884812Z","shell.execute_reply.started":"2025-03-16T19:16:43.846820Z","shell.execute_reply":"2025-03-16T19:16:43.883840Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# Visualize 3D structure of a sample target\nif len(sample_coords) > 0:\n    fig = plt.figure(figsize=(10, 8))\n    ax = fig.add_subplot(111, projection='3d')\n    ax.plot(sample_coords[0][:, 0], sample_coords[0][:, 1], sample_coords[0][:, 2], 'o-')\n    ax.set_title(f'3D Structure of {sample_target_id}')\n    ax.set_xlabel('X')\n    ax.set_ylabel('Y')\n    ax.set_zlabel('Z')\n    plt.savefig('sample_structure_3d.png')\n\n# Deduplicate sequences and merge structures\nsequence_to_structures = defaultdict(list)\nfor idx, row in train_sequences.iterrows():\n    target_id = row['target_id']\n    sequence = row['sequence']\n    residues = train_labels[train_labels['ID'].str.startswith(f\"{target_id}_\")]\n    coords = extract_coordinates(residues)\n    if len(coords) > 0:\n        sequence_to_structures[sequence].append((target_id, coords))\n\n# Show duplicates\nduplicate_sequences = [seq for seq, structs in sequence_to_structures.items() if len(structs) > 1]\nprint(f\"\\nNumber of duplicate sequences: {len(duplicate_sequences)}\")\nif duplicate_sequences:\n    print(f\"Example duplicate: {duplicate_sequences[0]}\")\n    print(f\"Number of structures: {len(sequence_to_structures[duplicate_sequences[0]])}\")\n\n# Function to prepare data for model input\ndef prepare_model_input(sequences_df, labels_df, msa_dir='/kaggle/input/stanford-rna-3d-folding/MSA/'):\n    \"\"\"Prepare features and labels for model training.\"\"\"\n    data = []\n    \n    for idx, row in tqdm(sequences_df.iterrows(), total=len(sequences_df)):\n        target_id = row['target_id']\n        sequence = row['sequence']\n        \n        # Skip sequences with non-standard nucleotides\n        if any(c not in 'ACGU' for c in sequence):\n            continue\n        \n        # One-hot encode sequence\n        seq_onehot = []\n        for c in sequence:\n            if c == 'A':\n                seq_onehot.append([1, 0, 0, 0])\n            elif c == 'C':\n                seq_onehot.append([0, 1, 0, 0])\n            elif c == 'G':\n                seq_onehot.append([0, 0, 1, 0])\n            elif c == 'U':\n                seq_onehot.append([0, 0, 0, 1])\n            else:\n                seq_onehot.append([0.25, 0.25, 0.25, 0.25])\n        \n        # Extract coordinates\n        residues = labels_df[labels_df['ID'].str.startswith(f\"{target_id}_\")]\n        coords = extract_coordinates(residues)\n        \n        # Process MSA (basic features)\n        msa_stats = process_msa(target_id, msa_dir)\n        \n        # Store data\n        data.append({\n            'target_id': target_id,\n            'sequence': sequence,\n            'sequence_onehot': np.array(seq_onehot),\n            'coordinates': coords,\n            'msa_stats': msa_stats\n        })\n    \n    return data\n\n# Prepare a small sample of data\nsample_size = min(10, len(train_sequences))\nsample_data = prepare_model_input(\n    train_sequences.head(sample_size),\n    train_labels\n)\n\nprint(f\"\\nSample data prepared: {len(sample_data)} entries\")\nif sample_data:\n    print(f\"Example entry: target_id={sample_data[0]['target_id']}, \" +\n          f\"sequence length={len(sample_data[0]['sequence'])}, \" +\n          f\"coordinates shape={sample_data[0]['coordinates'].shape if len(sample_data[0]['coordinates']) > 0 else 'None'}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-16T19:17:09.058809Z","iopub.execute_input":"2025-03-16T19:17:09.059227Z","iopub.status.idle":"2025-03-16T19:17:34.518711Z","shell.execute_reply.started":"2025-03-16T19:17:09.059193Z","shell.execute_reply":"2025-03-16T19:17:34.517904Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}