{"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.10.12"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":87793,"databundleVersionId":12276181,"isSourceIdPinned":false,"sourceType":"competition"},{"sourceId":10933126,"sourceType":"datasetVersion","datasetId":6798383}],"dockerImageVersionId":30918,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false},"papermill":{"default_parameters":{},"duration":6119.7991,"end_time":"2025-05-27T22:19:44.394148","environment_variables":{},"exception":null,"input_path":"__notebook__.ipynb","output_path":"__notebook__.ipynb","parameters":{},"start_time":"2025-05-27T20:37:44.595048","version":"2.6.0"}},"nbformat_minor":4,"nbformat":4,"cells":[{"id":"b80f46bf","cell_type":"markdown","source":"# **RNA 3D Structure Prediction Pipeline 🧬**\n\n## **Overview 📜**\n\nThis project implements a comprehensive pipeline for the Stanford RNA 3D Folding competition, aiming to predict the three-dimensional structure of RNA molecules from their nucleotide sequences. The latest version employs an **Adaptive Hybrid Energy Field Ensemble with Dynamical Attractors** approach using **Boltzmann weighting** and **adaptive temperature sampling**, representing a significant advancement in RNA structure prediction.\n\n---\n\n## **Key Components 🧩**\n\n### **1. Data Processing and Management 📊**\n\n**Memory-Optimized Data Loading**\n\n* Chunk-based loading for large datasets\n* Automatic datatype optimization (int8/16/32, float32, category)\n* Threshold-based categorical conversion\n* Special value handling (-1.0e+18)\n\n**Data Exploration and Analysis**\n\n* Sequence length distribution analysis\n* Nucleotide frequency calculation\n* 3D coordinate distribution visualization\n* ID mapping and structure verification\n\n**Feature Engineering**\n\n* One-hot encoding of RNA sequences (A, C, G, U, N)\n* Sequence padding for uniform input dimensions\n* Correlation preservation between coordinates\n* Structure normalization and centralization\n\n---\n\n### **2. Phase Space Conformational Modeling 🌌**\n\n**Conformational Phase Space Construction**\n\n* Generation of diverse base predictions from multiple models\n* Extraction of RNA-specific structural characteristics\n* Dimensionality reduction to create manageable phase space\n* Phase space topology analysis with density estimation\n\n**Dynamical Attractor Identification**\n\n* Clustering to identify attractor basins in phase space\n* Attractor stability assessment through Langevin dynamics simulation\n* Convergence rate calculation and dynamic stability scoring\n* Basin size and shape analysis for each attractor\n\n**Thermodynamic Principles Application**\n\n* Conversion of dynamic scores to pseudo-energies\n* Boltzmann distribution weighting with temperature control\n* Weighted structure selection based on attractor properties\n* Energy landscape modeling with multiple stable states\n\n---\n\n### **3. Reference-Based Modeling Strategy 🧮**\n\n**Multi-Seed Ensemble Approach**\n\n* Expanded set of base models for richer phase space coverage\n* Balanced seed selection for diverse predictions\n* Fixed seeds for reproducibility\n* Weighted ensemble creation based on TM-score performance\n\n**RNA Property-Specific Parameter Optimization**\n\n* GC content-based adjustments for thermodynamic stability\n* Sequence length-specific parameter tuning\n* RNA motif detection (tetraloops, pseudoknots, G-quadruplexes)\n* Size-specific dynamics modeling (persistence length, rigidity)\n\n**Advanced Parameter Search**\n\n* Two-phase parameter search (broad + refined)\n* Noise level and correlation parameter optimization\n* Ablation analysis to identify critical components\n* Grid search with randomized exploration\n\n---\n\n### **4. Structure Generation with Physical Dynamics 🎯**\n\n**Adaptive Temperature Sampling**\n\n* Physics-based sampling with Worm-Like Chain model principles\n* Temperature adjustment based on RNA properties\n* Progressive noise levels representing energy states\n* Boltzmann distribution-based structure generation\n\n**Metastable States Modeling**\n\n* Identification of potential metastable conformations\n* Diverse structure sampling from different energy basins\n* Representation of alternative folding pathways\n* Weighted ensemble based on basin stability\n\n**Dynamical Systems-Based Structure Generation**\n\n* Langevin dynamics simulation for trajectory analysis\n* Structure selection based on attractor proximity\n* Diversity-aware selection to cover conformational space\n* Hybrid energy field ensemble creation\n\n---\n\n### **5. Structure Validation and Refinement ✅**\n\n**Biophysical Validation Checks**\n\n* Bond distance constraints (0.8–7.0 Å)\n* Clash detection with tolerance threshold\n* Invalid bond percentage analysis\n* Structure completeness verification\n\n**Domain-Aware Structure Analysis**\n\n* Natural hinge points identification\n* Junction detection between helices\n* GC/AU content-based structural adjustments\n* Sequence motif recognition for structural features\n\n---\n\n### **6. Evaluation Metrics Implementation 📏**\n\n**Structure Similarity Assessment**\n\n* TM-score calculation with size-dependent scaling\n* Exact TM-score with multiple rotation schemes\n* Enhanced RNA-specific evaluation metrics\n* Structure-specific performance analysis\n\n**Statistical Performance Analysis**\n\n* TM-score distribution visualization\n* Cross-validation for parameter robustness\n* Size-category performance breakdowns\n* Dynamic stability vs. prediction accuracy analysis\n\n---\n\n### **7. Submission Generation and Pipeline Management 🔄**\n\n**Robust Submission Generation**\n\n* Multi-layer fallback system for error handling\n* Ensemble prediction with attractor-based weighting\n* Multiple structure generation per sequence\n* Graceful degradation with simplified approaches\n\n**Pipeline Engineering**\n\n* Exception handling throughout the pipeline\n* File existence and validity checking\n* Progress tracking and reporting\n* Checkpointing intermediate results\n\n---\n\n## **Methodology 🔍**\n\nThe pipeline now implements a sophisticated **Adaptive Hybrid Energy Field Ensemble with Dynamical Attractors** approach:\n\n### **1. Phase Space Construction and Attractor Identification**\n\n* Multiple base models generate diverse predictions to sample the conformational phase space\n* Dimensionality reduction creates a manageable 3D representation of this space\n* Clustering techniques identify potential attractor basins representing stable states\n* Langevin dynamics simulations assess the stability and convergence properties of each attractor\n\n### **2. Thermodynamic-Based Ensemble Creation**\n\n* Dynamic scores are converted to pseudo-energies using statistical mechanics principles\n* Boltzmann weighting (T = 0.2) assigns probabilities to different attractor states\n* Structures are selected from each attractor basin with diversity considerations\n* A weighted ensemble structure is created as the primary model\n\n### **3. RNA-Specific Adaptive Sampling**\n\n* Sequence properties (GC content, length) inform parameter adjustments\n* Adaptive noise levels based on attractor stability (typically 0.17–0.21)\n* RNA motif detection influences thermodynamic parameters\n* Worm-Like Chain model with twist-bend coupling for realistic flexibility\n\n### **4. Multi-Level Fallback System**\n\n* Primary approach uses phase space and attractor dynamics\n* Secondary approach employs metastable state detection and Boltzmann weighting\n* Tertiary approach uses balanced ensemble with parameter diversity\n* Final fallback uses simplified reference model if all else fails\n\n---\n\n## **Technical Implementation Highlights**\n\n* **Expanded Model Ensemble**: Uses 8 base models with diverse seeds (84, 294, 1134, 420, 252, 1001, 3780, 756) to ensure comprehensive conformational space coverage.\n* **Sophisticated Phase Space Analysis**: Implements density estimation, clustering, and topology analysis to identify meaningful attractor basins.\n* **Physical Dynamics Simulation**: Employs Langevin dynamics to evaluate attractor stability and convergence properties through multi-step simulations.\n* **Enhanced Structure Selection**: Uses both proximity to attractors and diversity considerations to select representative structures from the conformational landscape.\n* **Adaptive Parameter Tuning**: Automatically adjusts noise levels (approximately 0.15–0.21) based on attractor stability, sequence properties, and RNA characteristics.\n\n---\n\nThis approach represents a significant advancement in RNA structure prediction by incorporating principles from **statistical mechanics**, **polymer physics**, and **dynamical systems theory** to create a physically informed model of RNA folding.\n\n---","metadata":{"papermill":{"duration":0.025642,"end_time":"2025-05-27T20:37:48.186713","exception":false,"start_time":"2025-05-27T20:37:48.161071","status":"completed"},"tags":[]}},{"id":"629939c7","cell_type":"markdown","source":"## Library Imports 📚🔧","metadata":{"papermill":{"duration":0.021066,"end_time":"2025-05-27T20:37:48.229428","exception":false,"start_time":"2025-05-27T20:37:48.208362","status":"completed"},"tags":[]}},{"id":"f61034f1","cell_type":"code","source":"!pip install /kaggle/input/biopython-1-85/biopython-1.85-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl","metadata":{"papermill":{"duration":7.310722,"end_time":"2025-05-27T20:37:55.562634","exception":false,"start_time":"2025-05-27T20:37:48.251912","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:16:22.380103Z","iopub.execute_input":"2025-10-13T02:16:22.380569Z","iopub.status.idle":"2025-10-13T02:16:30.129840Z","shell.execute_reply.started":"2025-10-13T02:16:22.380533Z","shell.execute_reply":"2025-10-13T02:16:30.128154Z"}},"outputs":[],"execution_count":null},{"id":"92493bdc","cell_type":"code","source":"# =============================================================================\n# 🔧 REPRODUCIBILITY CONFIGURATION\n# =============================================================================\n\n# Master seed for entire pipeline\nMASTER_SEED = 8339\n\n# Environment variables for deterministic operations\nimport os\nos.environ[\"OMP_NUM_THREADS\"] = \"1\"\nos.environ[\"OPENBLAS_NUM_THREADS\"] = \"1\"\nos.environ[\"MKL_NUM_THREADS\"] = \"1\"\nos.environ[\"VECLIB_MAXIMUM_THREADS\"] = \"1\"\nos.environ[\"NUMEXPR_NUM_THREADS\"] = \"1\"\nos.environ['PYTHONHASHSEED'] = str(MASTER_SEED)\nos.environ['TF_DETERMINISTIC_OPS'] = '1'\nos.environ['SKLEARN_SEED'] = str(MASTER_SEED)\n\n# =============================================================================\n# 📦 STANDARD LIBRARY IMPORTS\n# =============================================================================\n\nimport sys\nimport time\nimport gc\nimport random\nimport traceback\nimport warnings\nimport glob\nimport pickle\nimport hashlib\nfrom datetime import datetime\nfrom collections import Counter, defaultdict\nfrom concurrent.futures import ProcessPoolExecutor\n\n# Configure warnings\nwarnings.filterwarnings('ignore')\n\n# =============================================================================\n# 🔢 NUMERICAL COMPUTING\n# =============================================================================\n\nimport numpy as np\nimport pandas as pd\nfrom scipy.spatial.distance import pdist, squareform\nfrom scipy.optimize import minimize\n\n# =============================================================================\n# 🤖 MACHINE LEARNING CORE\n# =============================================================================\n\nfrom sklearn.cluster import DBSCAN\nfrom sklearn.utils import check_random_state\nfrom sklearn.base import clone\nfrom sklearn.decomposition import PCA\nfrom sklearn.model_selection import KFold\nfrom sklearn.utils import resample\n\n# =============================================================================\n# 🔥 DEEP LEARNING (PyTorch)\n# =============================================================================\n\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.data import DataLoader, TensorDataset\n\n# =============================================================================\n# 📊 VISUALIZATION\n# =============================================================================\n\nimport matplotlib.pyplot as plt\nimport matplotlib.colors as mcolors\nfrom mpl_toolkits.mplot3d import Axes3D\nimport seaborn as sns\nfrom matplotlib.gridspec import GridSpec\n\n# =============================================================================\n# 🧬 BIOINFORMATICS AND STRUCTURAL BIOLOGY\n# =============================================================================\n\ntry:\n    from Bio.PDB import MMCIFParser, Select\n    from Bio.PDB.MMCIF2Dict import MMCIF2Dict\n    from Bio.Seq import Seq\n    from Bio.SeqUtils import seq1\n    import Bio.PDB\n    print(\"✅ BioPython libraries loaded successfully\")\nexcept ImportError:\n    print(\"⚠️ BioPython not available - some features may be limited\")\n\n# =============================================================================\n# 🌐 DATA STRUCTURES AND ALGORITHMS\n# =============================================================================\n\ntry:\n    import networkx as nx\n    print(\"✅ NetworkX loaded successfully\")\nexcept ImportError:\n    print(\"⚠️ NetworkX not available - graph analysis features disabled\")\n\n# =============================================================================\n# 🎯 REPRODUCIBILITY SETUP\n# =============================================================================\n\ndef set_global_seed(seed):\n    \"\"\"Set global seed for reproducibility across all libraries\"\"\"\n    random.seed(seed)\n    np.random.seed(seed)\n    torch.manual_seed(seed)\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    \n    # CUDA reproducibility (if available)\n    if torch.cuda.is_available():\n        torch.cuda.manual_seed(seed)\n        torch.cuda.manual_seed_all(seed)\n        torch.backends.cudnn.deterministic = True\n        torch.backends.cudnn.benchmark = False\n\n# Apply global seed\nset_global_seed(MASTER_SEED)\n\n# =============================================================================\n# 🔧 CONDITIONAL IMPORTS AND CONFIGURATION\n# =============================================================================\n\n# TensorFlow Configuration (if available)\ntry:\n    import tensorflow as tf\n    \n    # Set TensorFlow seed for reproducibility\n    tf.random.set_seed(MASTER_SEED)\n    \n    # Configure threads for deterministic operations\n    tf.config.threading.set_inter_op_parallelism_threads(1)\n    tf.config.threading.set_intra_op_parallelism_threads(1)\n    \n    # Enable determinism for ops (available in TF 2.9+)\n    try:\n        tf.config.experimental.enable_op_determinism()\n        print(\"✅ TensorFlow determinism enabled\")\n    except AttributeError:\n        print(\"⚠️ TensorFlow experimental determinism not available in this version\")\n        \nexcept ImportError:\n    print(\"⚠️ TensorFlow not available - PyTorch will be used for deep learning\")\n\n# =============================================================================\n# 📁 DIRECTORY CONFIGURATION\n# =============================================================================\n\n# Kaggle paths configuration\nDATA_DIR = \"/kaggle/input/stanford-rna-3d-folding\"\nOUTPUT_DIR = \"/kaggle/working\"\nWORKING_DIR = \"/kaggle/working\"\n\n# Ensure output directory exists\nos.makedirs(OUTPUT_DIR, exist_ok=True)\n\n# =============================================================================\n# ✅ FINAL REPRODUCIBILITY VERIFICATION\n# =============================================================================\n\n# Force sequential thread model for NumPy (redundant but ensures consistency)\nnp.random.seed(MASTER_SEED)\n\n# Verify PyTorch determinism\ntorch.manual_seed(MASTER_SEED)\nif torch.cuda.is_available():\n    torch.cuda.manual_seed_all(MASTER_SEED)\n\nprint(f\"🎯 Master Seed: {MASTER_SEED}\")\nprint(f\"🐍 Python: {sys.version.split()[0]}\")\nprint(f\"📊 NumPy: {np.__version__}\")\nprint(f\"🐼 Pandas: {pd.__version__}\")\nprint(f\"🔥 PyTorch: {torch.__version__}\")\nprint(f\"🤖 Device: {'CUDA' if torch.cuda.is_available() else 'CPU'}\")\nprint(f\"📁 Data Dir: {DATA_DIR}\")\nprint(f\"💾 Output Dir: {OUTPUT_DIR}\")\nprint(\"✅ All libraries configured with reproducibility settings\")\nprint(\"=\" * 80)","metadata":{"papermill":{"duration":3.968275,"end_time":"2025-05-27T20:37:59.553656","exception":false,"start_time":"2025-05-27T20:37:55.585381","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:16:30.131975Z","iopub.execute_input":"2025-10-13T02:16:30.132482Z","iopub.status.idle":"2025-10-13T02:16:55.065791Z","shell.execute_reply.started":"2025-10-13T02:16:30.132435Z","shell.execute_reply":"2025-10-13T02:16:55.064526Z"}},"outputs":[],"execution_count":null},{"id":"3fb5bd5d","cell_type":"markdown","source":"## 🧬 RNA 3D Structure Prediction and Analysis Pipeline 🔬","metadata":{"papermill":{"duration":0.021798,"end_time":"2025-05-27T20:37:59.597939","exception":false,"start_time":"2025-05-27T20:37:59.576141","status":"completed"},"tags":[]}},{"id":"e74ac42c","cell_type":"code","source":"# Directories and files adjusted for the new competition\nDATA_DIR = os.getenv('DATA_DIR', '/kaggle/input/stanford-rna-3d-folding/')\nmain_files = [\n    \"train_sequences.csv\", \n    \"train_labels.csv\", \n    \"validation_sequences.csv\", \n    \"validation_labels.csv\", \n    \"test_sequences.csv\",\n    \"sample_submission.csv\"\n]\n\nDEFAULT_THRESHOLD = 0.4  # Default threshold after analysis\n\ndef optimize_dataframe(df, inplace=False, category_threshold=DEFAULT_THRESHOLD):\n    \"\"\"\n    Optimizes the DataFrame to save memory.\n    \"\"\"\n    if category_threshold < 0 or category_threshold > 1:\n        raise ValueError(\"category_threshold must be between 0 and 1.\")\n    \n    if not inplace:\n        df = df.copy()\n    \n    for col in df.columns:\n        col_type = df[col].dtype\n        if np.issubdtype(col_type, np.integer):\n            c_min, c_max = df[col].min(), df[col].max()\n            if c_min > np.iinfo(np.int8).min and c_max < np.iinfo(np.int8).max:\n                df[col] = df[col].astype(np.int8)\n            elif c_min > np.iinfo(np.int16).min and c_max < np.iinfo(np.int16).max:\n                df[col] = df[col].astype(np.int16)\n            elif c_min > np.iinfo(np.int32).min and c_max < np.iinfo(np.int32).max:\n                df[col] = df[col].astype(np.int32)\n        elif np.issubdtype(col_type, np.floating):\n            if df[col].min() > np.finfo(np.float32).min and df[col].max() < np.finfo(np.float32).max:\n                df[col] = df[col].astype(np.float32)\n        if col_type == object:\n            unique_vals = len(df[col].unique())\n            if unique_vals / len(df) < category_threshold:\n                df[col] = df[col].astype('category')\n    \n    return df\n\ndef load_main_data(chunksize=50000):\n    \"\"\"\n    Loads the main files.\n    \"\"\"\n    data = {}\n    for file_name in main_files:\n        file_path = os.path.join(DATA_DIR, file_name)\n        if os.path.exists(file_path):\n            chunks = pd.read_csv(file_path, on_bad_lines='skip', low_memory=False, chunksize=chunksize)\n            dataframes = [optimize_dataframe(chunk, category_threshold=DEFAULT_THRESHOLD) for chunk in chunks]\n            data[file_name] = pd.concat(dataframes, ignore_index=True)\n        else:\n            print(f\"File {file_path} not found!\")\n    return data\n\ndef check_data_integrity(original_df, optimized_df):\n    \"\"\"\n    Checks if the optimization did not alter the data.\n    \"\"\"\n    try:\n        pd.testing.assert_frame_equal(original_df, optimized_df, check_like=True)\n        print(\"Integrity check passed: No changes in data after optimization.\")\n    except AssertionError as e:\n        print(f\"Data integrity check failed: {e}\")\n\ndef check_duplicates(df):\n    \"\"\"\n    Checks for duplicates in the DataFrame.\n    \"\"\"\n    duplicates = df[df.duplicated(keep=False)]\n    if not duplicates.empty:\n        print(f\"Warning: Duplicates found in the dataset. Number of duplicates: {duplicates.shape[0]}\")\n        return duplicates\n    else:\n        print(\"No duplicates found.\")\n    return None\n\ndef test_thresholds(df):\n    \"\"\"\n    Tests different thresholds for DataFrame optimization.\n    \"\"\"\n    thresholds = np.linspace(0.1, 0.9, 9)\n    memory_usages = []\n    for threshold in thresholds:\n        optimized_df = optimize_dataframe(df.copy(), category_threshold=threshold)\n        memory_usages.append(optimized_df.memory_usage(deep=True).sum() / 1024**2)\n    return thresholds, memory_usages\n\ndef plot_memory_usage(thresholds, memory_usages):\n    \"\"\"\n    Plots memory usage versus thresholds.\n    \"\"\"\n    plt.figure(figsize=(10, 6))\n    plt.plot(thresholds, memory_usages, marker='o', linestyle='-')\n    plt.title(\"Memory Usage vs. Threshold\")\n    plt.xlabel(\"Threshold\")\n    plt.ylabel(\"Memory Usage (MB)\")\n    plt.grid(True)\n    plt.show()\n\ndef analyze_sequence_data(df_sequences):\n    \"\"\"\n    Analyzes RNA sequence data.\n    \"\"\"\n    # Basic information\n    print(f\"Total sequences: {len(df_sequences)}\")\n    print(f\"Available columns: {df_sequences.columns.tolist()}\")\n    \n    # Sequence analysis\n    if 'sequence' in df_sequences.columns:\n        # Distribution of sequence lengths\n        seq_lengths = df_sequences['sequence'].apply(len)\n        print(f\"\\nSequence length statistics:\")\n        print(f\"Minimum: {seq_lengths.min()}\")\n        print(f\"Maximum: {seq_lengths.max()}\")\n        print(f\"Average: {seq_lengths.mean():.2f}\")\n        \n        # Nucleotide count\n        nucleotides = ['A', 'C', 'G', 'U']\n        nucleotide_counts = {n: df_sequences['sequence'].str.count(n).sum() for n in nucleotides}\n        total_nucleotides = sum(nucleotide_counts.values())\n        \n        print(\"\\nNucleotide distribution:\")\n        for n, count in nucleotide_counts.items():\n            print(f\"{n}: {count} ({count/total_nucleotides*100:.2f}%)\")\n    \n    return df_sequences\n\ndef analyze_label_data(df_labels):\n    \"\"\"\n    Analyzes 3D coordinate data (labels).\n    \"\"\"\n    print(f\"Total entries in labels: {len(df_labels)}\")\n    print(f\"Available columns: {df_labels.columns.tolist()}\")\n    \n    # Analysis of 3D coordinates if available\n    coord_columns = [col for col in df_labels.columns if col.startswith(('x_', 'y_', 'z_'))]\n    if coord_columns:\n        print(f\"\\nCoordinate columns found: {len(coord_columns)}\")\n        \n        # Basic statistics of coordinates\n        for i in range(1, 6):  # For the 5 possible structures\n            x_col = f'x_{i}'\n            y_col = f'y_{i}'\n            z_col = f'z_{i}'\n            \n            if x_col in df_labels.columns and y_col in df_labels.columns and z_col in df_labels.columns:\n                print(f\"\\nStatistics for structure {i}:\")\n                print(f\"X - Mean: {df_labels[x_col].mean():.2f}, Std: {df_labels[x_col].std():.2f}\")\n                print(f\"Y - Mean: {df_labels[y_col].mean():.2f}, Std: {df_labels[y_col].std():.2f}\")\n                print(f\"Z - Mean: {df_labels[z_col].mean():.2f}, Std: {df_labels[z_col].std():.2f}\")\n    \n    return df_labels\n\ndef create_submission_template(test_df, sample_submission_df):\n    \"\"\"\n    Creates a submission template based on test data.\n    \"\"\"\n    # Check if sample_submission.csv is available\n    if sample_submission_df is None:\n        print(\"Sample submission file not found. Creating a new template.\")\n        \n        # Create a new DataFrame for submission\n        submission_df = pd.DataFrame()\n        \n        # Example code to fill the template (adjust as needed)\n        ids = []\n        resnames = []\n        resids = []\n        \n        for _, row in test_df.iterrows():\n            sequence = row['sequence']\n            target_id = row['target_id']\n            \n            for i, nucleotide in enumerate(sequence, 1):\n                ids.append(f\"{target_id}_{i}\")\n                resnames.append(nucleotide)\n                resids.append(i)\n        \n        submission_df['ID'] = ids\n        submission_df['resname'] = resnames\n        submission_df['resid'] = resids\n        \n        # Add coordinate columns (5 structures)\n        for i in range(1, 6):\n            submission_df[f'x_{i}'] = 0.0\n            submission_df[f'y_{i}'] = 0.0\n            submission_df[f'z_{i}'] = 0.0\n    else:\n        submission_df = sample_submission_df.copy()\n        print(\"Submission template created based on the provided example.\")\n    \n    return submission_df\n\ndef main():\n    start_time = time.time()\n    \n    # Load main data\n    print(\"Loading main data...\")\n    main_data = load_main_data()\n    \n    # Check which files were loaded\n    print(\"\\nLoaded files:\")\n    for file_name, df in main_data.items():\n        print(f\"- {file_name}: {df.shape if df is not None else 'Not found'}\")\n    \n    # Analyze training sequence data\n    if \"train_sequences.csv\" in main_data:\n        print(\"\\n===== Training Sequences Analysis =====\")\n        analyze_sequence_data(main_data[\"train_sequences.csv\"])\n    \n    # Analyze training label data\n    if \"train_labels.csv\" in main_data:\n        print(\"\\n===== Training Labels Analysis =====\")\n        analyze_label_data(main_data[\"train_labels.csv\"])\n    \n    # Check for duplicates in training data\n    if \"train_sequences.csv\" in main_data:\n        print(\"\\nChecking for duplicates in training sequences...\")\n        check_duplicates(main_data[\"train_sequences.csv\"])\n    \n    # Create submission template\n    if \"test_sequences.csv\" in main_data:\n        print(\"\\nCreating submission template...\")\n        submission_template = create_submission_template(\n            main_data[\"test_sequences.csv\"],\n            main_data.get(\"sample_submission.csv\")\n        )\n        print(f\"Submission template shape: {submission_template.shape}\")\n        print(f\"First rows of the submission template:\")\n        print(submission_template.head())\n    \n    # Calculate execution time\n    end_time = time.time()\n    print(f\"\\nRuntime: {end_time - start_time:.2f} seconds\")\n    \n    return main_data\n\nif __name__ == '__main__':\n    main_data = main()","metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","papermill":{"duration":0.805926,"end_time":"2025-05-27T20:38:00.426212","exception":false,"start_time":"2025-05-27T20:37:59.620286","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:18:58.475093Z","iopub.execute_input":"2025-10-13T02:18:58.475565Z","iopub.status.idle":"2025-10-13T02:18:59.303646Z","shell.execute_reply.started":"2025-10-13T02:18:58.475534Z","shell.execute_reply":"2025-10-13T02:18:59.302544Z"}},"outputs":[],"execution_count":null},{"id":"dfcc6d5c","cell_type":"markdown","source":"## Directory Explorer & CSV Verification for RNA3D 🗂️🔬","metadata":{"papermill":{"duration":0.022421,"end_time":"2025-05-27T20:38:00.471068","exception":false,"start_time":"2025-05-27T20:38:00.448647","status":"completed"},"tags":[]}},{"id":"b7748b6c","cell_type":"code","source":"# Updated main directory\ndir_main = \"/kaggle/input/stanford-rna-3d-folding/\"\n\n# List all files and directories in the main directory\ntry:\n    all_files = os.listdir(dir_main)\n    print(f\"All files and directories in '{dir_main}':\")\n    \n    for file in all_files:\n        # Check if it's a file or directory\n        full_path = os.path.join(dir_main, file)\n        type_desc = \"directory\" if os.path.isdir(full_path) else \"file\"\n        size = os.path.getsize(full_path) / 1024  # Size in KB\n        print(f\" - {file} ({type_desc}, {size:.2f} KB)\")\n        \n        # If it's a directory, list up to 5 files inside it\n        if os.path.isdir(full_path):\n            try:\n                internal_files = os.listdir(full_path)[:5]  # Limit to 5 files\n                if internal_files:\n                    print(f\"   First files in '{file}':\")\n                    for internal_file in internal_files:\n                        print(f\"    * {internal_file}\")\n                    if len(os.listdir(full_path)) > 5:\n                        print(f\"    * ... and {len(os.listdir(full_path)) - 5} more file(s)\")\n                else:\n                    print(f\"   '{file}' is empty\")\n            except Exception as e:\n                print(f\"   Error listing contents of '{file}': {e}\")\nexcept Exception as e:\n    print(f\"Error listing directory {dir_main}: {e}\")\n\n# Check the structure of the main CSV files\nmain_files = [\n    \"train_sequences.csv\", \n    \"train_labels.csv\", \n    \"validation_sequences.csv\", \n    \"validation_labels.csv\", \n    \"test_sequences.csv\",\n    \"sample_submission.csv\"\n]\nprint(\"\\nChecking main CSV files:\")\n\nfor file in main_files:\n    full_path = os.path.join(dir_main, file)\n    if os.path.exists(full_path):\n        # Get file size\n        size_mb = os.path.getsize(full_path) / (1024 * 1024)  # Size in MB\n        \n        # Read the first lines to check the structure\n        try:\n            import pandas as pd\n            df = pd.read_csv(full_path, nrows=1)\n            print(f\"\\n{file} ({size_mb:.2f} MB):\")\n            print(f\"Columns: {df.columns.tolist()}\")\n            print(f\"Example:\")\n            print(df.head())\n        except Exception as e:\n            print(f\"Error reading {file}: {e}\")\n    else:\n        print(f\"{file} not found.\")","metadata":{"papermill":{"duration":0.250206,"end_time":"2025-05-27T20:38:00.743655","exception":false,"start_time":"2025-05-27T20:38:00.493449","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:19:07.591944Z","iopub.execute_input":"2025-10-13T02:19:07.592423Z","iopub.status.idle":"2025-10-13T02:19:07.871281Z","shell.execute_reply.started":"2025-10-13T02:19:07.592391Z","shell.execute_reply":"2025-10-13T02:19:07.869964Z"}},"outputs":[],"execution_count":null},{"id":"be26b2b0","cell_type":"markdown","source":"## RNA3D Data Checker 🔍🧬","metadata":{"papermill":{"duration":0.025692,"end_time":"2025-05-27T20:38:00.794556","exception":false,"start_time":"2025-05-27T20:38:00.768864","status":"completed"},"tags":[]}},{"id":"f5485cea","cell_type":"code","source":"# Updated main directory\ndir_main = \"/kaggle/input/stanford-rna-3d-folding/\"\n\ndef load_data():\n    \"\"\"\n    Loads the main CSV files from the Stanford RNA 3D Folding competition.\n    Returns a dictionary with DataFrames.\n    \"\"\"\n    main_files = [\n        \"train_sequences.csv\", \n        \"train_labels.csv\", \n        \"validation_sequences.csv\", \n        \"validation_labels.csv\", \n        \"test_sequences.csv\",\n        \"sample_submission.csv\"\n    ]\n    \n    data = {}\n    for file_name in main_files:\n        file_path = os.path.join(dir_main, file_name)\n        if os.path.exists(file_path):\n            try:\n                data[file_name] = pd.read_csv(file_path)\n                print(f\"File {file_name} loaded successfully. Shape: {data[file_name].shape}\")\n            except Exception as e:\n                print(f\"Error loading {file_name}: {e}\")\n        else:\n            print(f\"File {file_name} not found.\")\n            data[file_name] = None\n    \n    return data\n\ndef compare_columns(main_data):\n    \"\"\"\n    Compares columns between different DataFrames.\n    \"\"\"\n    # List all available keys\n    print(\"\\nLoaded files:\")\n    print(list(main_data.keys()))\n    \n    # Compare columns between train_sequences.csv and test_sequences.csv\n    if \"train_sequences.csv\" in main_data and \"test_sequences.csv\" in main_data:\n        train_cols = set(main_data[\"train_sequences.csv\"].columns)\n        test_cols = set(main_data[\"test_sequences.csv\"].columns)\n        \n        print(\"\\nColumns in train_sequences.csv:\")\n        print(list(main_data[\"train_sequences.csv\"].columns))\n        \n        print(\"\\nUnique columns in train_sequences.csv (not present in test_sequences.csv):\")\n        print(train_cols - test_cols)\n        \n        print(\"\\nUnique columns in test_sequences.csv (not present in train_sequences.csv):\")\n        print(test_cols - train_cols)\n    \n    # Compare columns between train_labels.csv and validation_labels.csv\n    if \"train_labels.csv\" in main_data and \"validation_labels.csv\" in main_data:\n        train_label_cols = set(main_data[\"train_labels.csv\"].columns)\n        val_label_cols = set(main_data[\"validation_labels.csv\"].columns)\n        \n        print(\"\\nColumns in train_labels.csv:\")\n        print(list(main_data[\"train_labels.csv\"].columns))\n        \n        print(\"\\nColumns in validation_labels.csv:\")\n        print(list(main_data[\"validation_labels.csv\"].columns))\n        \n        print(\"\\nUnique columns in validation_labels.csv (not present in train_labels.csv):\")\n        print(val_label_cols - train_label_cols)\n    \n    # Compare columns between validation_labels.csv and sample_submission.csv\n    if \"validation_labels.csv\" in main_data and \"sample_submission.csv\" in main_data:\n        val_label_cols = set(main_data[\"validation_labels.csv\"].columns)\n        sample_cols = set(main_data[\"sample_submission.csv\"].columns)\n        \n        print(\"\\nColumns in sample_submission.csv:\")\n        print(list(main_data[\"sample_submission.csv\"].columns))\n        \n        print(\"\\nUnique columns in validation_labels.csv (not present in sample_submission.csv):\")\n        print(val_label_cols - sample_cols)\n        \n        print(\"\\nUnique columns in sample_submission.csv (not present in validation_labels.csv):\")\n        print(sample_cols - val_label_cols)\n\ndef analyze_structure_format(main_data):\n    \"\"\"\n    Analyzes the format of 3D structures (coordinates).\n    \"\"\"\n    if \"validation_labels.csv\" in main_data and main_data[\"validation_labels.csv\"] is not None:\n        df = main_data[\"validation_labels.csv\"]\n        \n        # Find all coordinate columns (x_1, y_1, z_1, etc.)\n        coord_cols = [col for col in df.columns if col.startswith(('x_', 'y_', 'z_'))]\n        \n        # Group by structure\n        structures = {}\n        for col in coord_cols:\n            # Extract structure number (e.g., \"x_1\" -> 1)\n            parts = col.split('_')\n            if len(parts) == 2:\n                struct_num = int(parts[1])\n                coord_type = parts[0]\n                \n                if struct_num not in structures:\n                    structures[struct_num] = []\n                \n                structures[struct_num].append(col)\n        \n        print(\"\\nStructure of the labels file:\")\n        print(f\"Total structures found: {len(structures)}\")\n        \n        # Show details of the first structure\n        if structures:\n            first_struct = min(structures.keys())\n            print(f\"\\nDetails of structure {first_struct}:\")\n            print(f\"Columns: {sorted(structures[first_struct])}\")\n            \n            # Check for missing values\n            for col in structures[first_struct]:\n                missing = df[col].isna().sum()\n                total = len(df)\n                print(f\"{col}: {missing} missing values ({missing/total*100:.2f}%)\")\n            \n            # Check the range of non-missing values for the first structure\n            for col in structures[first_struct]:\n                non_null = df[col][df[col] != -1.0e+18]  # Values that are not -1.0e+18\n                if not non_null.empty:\n                    print(f\"{col} - Range: [{non_null.min():.3f}, {non_null.max():.3f}]\")\n\ndef main():\n    # Load the data\n    main_data = load_data()\n    \n    # Compare columns between different files\n    compare_columns(main_data)\n    \n    # Analyze the format of 3D structures\n    analyze_structure_format(main_data)\n    \n    return main_data\n\nif __name__ == '__main__':\n    main_data = main()","metadata":{"papermill":{"duration":0.367479,"end_time":"2025-05-27T20:38:01.184632","exception":false,"start_time":"2025-05-27T20:38:00.817153","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:19:12.509250Z","iopub.execute_input":"2025-10-13T02:19:12.509641Z","iopub.status.idle":"2025-10-13T02:19:12.861601Z","shell.execute_reply.started":"2025-10-13T02:19:12.509610Z","shell.execute_reply":"2025-10-13T02:19:12.860246Z"}},"outputs":[],"execution_count":null},{"id":"7a3ca756","cell_type":"markdown","source":"## Integrated RNA3D Sequence and Structure Analyzer 🔬🧬","metadata":{"papermill":{"duration":0.023056,"end_time":"2025-05-27T20:38:01.230458","exception":false,"start_time":"2025-05-27T20:38:01.207402","status":"completed"},"tags":[]}},{"id":"7bd1804e","cell_type":"code","source":"# Initialize seed to control randomness\nnp.random.seed(0)\n\n# Directories and files adjusted for the new competition\nDATA_DIR = os.getenv('DATA_DIR', '/kaggle/input/stanford-rna-3d-folding/')\nmain_files = [\n   \"train_sequences.csv\", \n   \"train_labels.csv\", \n   \"validation_sequences.csv\", \n   \"validation_labels.csv\", \n   \"test_sequences.csv\",\n   \"sample_submission.csv\"\n]\n\nDEFAULT_THRESHOLD = 0.4  # Default threshold after analysis\n\ndef optimize_dataframe(df, inplace=False, category_threshold=DEFAULT_THRESHOLD):\n   \"\"\"\n   Optimizes the DataFrame to save memory.\n   \"\"\"\n   if category_threshold < 0 or category_threshold > 1:\n       raise ValueError(\"category_threshold must be between 0 and 1.\")\n   \n   if not inplace:\n       df = df.copy()\n   \n   for col in df.columns:\n       col_type = df[col].dtype\n       if np.issubdtype(col_type, np.integer):\n           c_min, c_max = df[col].min(), df[col].max()\n           if c_min > np.iinfo(np.int8).min and c_max < np.iinfo(np.int8).max:\n               df[col] = df[col].astype(np.int8)\n           elif c_min > np.iinfo(np.int16).min and c_max < np.iinfo(np.int16).max:\n               df[col] = df[col].astype(np.int16)\n           elif c_min > np.iinfo(np.int32).min and c_max < np.iinfo(np.int32).max:\n               df[col] = df[col].astype(np.int32)\n       elif np.issubdtype(col_type, np.floating):\n           # First check if it's not the special value -1.0e+18\n           if df[col].min() > np.finfo(np.float32).min and df[col].max() < np.finfo(np.float32).max:\n               df[col] = df[col].astype(np.float32)\n       if col_type == object:\n           unique_vals = len(df[col].unique())\n           if unique_vals / len(df) < category_threshold:\n               df[col] = df[col].astype('category')\n   \n   return df\n\ndef load_main_data(chunksize=50000):\n   \"\"\"\n   Loads the main files.\n   \"\"\"\n   data = {}\n   for file_name in main_files:\n       file_path = os.path.join(DATA_DIR, file_name)\n       if os.path.exists(file_path):\n           chunks = pd.read_csv(file_path, on_bad_lines='skip', low_memory=False, chunksize=chunksize)\n           dataframes = [optimize_dataframe(chunk, category_threshold=DEFAULT_THRESHOLD) for chunk in chunks]\n           data[file_name] = pd.concat(dataframes, ignore_index=True)\n           print(f\"File {file_name} loaded successfully. Shape: {data[file_name].shape}\")\n       else:\n           print(f\"File {file_path} not found!\")\n           data[file_name] = None\n   return data\n\ndef filter_columns_by_prefix(df, prefix=\"x_\"):\n   \"\"\"\n   Filters and counts the number of columns in a DataFrame based on a provided prefix.\n   \n   :param df: DataFrame where filtering will be applied.\n   :param prefix: Prefix to be used for filtering. Ex: \"x_\", \"y_\", \"z_\".\n   :return: List of filtered columns.\n   \"\"\"\n   filtered_columns = [col for col in df.columns if col.startswith(prefix)]\n   return filtered_columns\n\ndef count_nucleotides(df, column_name='sequence'):\n   \"\"\"\n   Counts the frequency of each nucleotide in a specific column of a DataFrame.\n   \n   :param df: DataFrame containing the sequences.\n   :param column_name: Name of the column containing the sequences. Default is 'sequence'.\n   :return: Counter object with the nucleotide counts.\n   \"\"\"\n   from collections import Counter\n\n   # Check if the column exists in the DataFrame\n   if column_name not in df.columns:\n       raise ValueError(f\"Column '{column_name}' not found in DataFrame.\")\n   \n   # Concatenate all sequences and count nucleotides\n   all_sequences = ''.join(df[column_name].tolist())\n   nucleotide_counts = Counter(all_sequences)\n   \n   return nucleotide_counts\n\ndef get_columns_without_missing_values(df):\n   \"\"\"\n   Returns columns without any missing values in the DataFrame.\n   \n   :param df: DataFrame to be checked.\n   :return: List of columns without missing values.\n   \"\"\"\n   missing_values = df.isnull().sum()\n   return missing_values[missing_values == 0].index.tolist()\n\ndef get_empty_columns(df):\n   \"\"\"\n   Returns columns that are completely empty in the DataFrame.\n   \n   :param df: DataFrame to be checked.\n   :return: List of empty columns.\n   \"\"\"\n   missing_values = df.isnull().sum()\n   return missing_values[missing_values == df.shape[0]].index.tolist()\n\ndef plot_coord_distributions(df_labels, prefix='x_', max_structures=5):\n   \"\"\"\n   Plots the distribution of coordinates (x, y, or z) for up to max_structures structures.\n   \n   :param df_labels: DataFrame containing the coordinates.\n   :param prefix: Prefix of columns to be plotted ('x_', 'y_', or 'z_').\n   :param max_structures: Maximum number of structures to show.\n   \"\"\"\n   # Find coordinate columns with the specified prefix\n   coord_cols = filter_columns_by_prefix(df_labels, prefix)\n   \n   # Limit to the maximum number of structures\n   coord_cols = sorted(coord_cols)[:max_structures]\n   \n   if not coord_cols:\n       print(f\"No column with prefix '{prefix}' found.\")\n       return\n   \n   # Set up the plot\n   fig, axes = plt.subplots(1, len(coord_cols), figsize=(16, 4))\n   if len(coord_cols) == 1:\n       axes = [axes]  # Ensure axes is iterable even with a single subplot\n   \n   # Plot histograms for each column\n   for i, col in enumerate(coord_cols):\n       # Filter special values (-1.0e+18) if present\n       values = df_labels[col]\n       filtered_values = values[values > -1.0e+17]  # Cutoff value to filter -1.0e+18\n       \n       axes[i].hist(filtered_values, bins=30, alpha=0.7)\n       axes[i].set_title(f'Distribution of {col}')\n       axes[i].set_xlabel('Value')\n       axes[i].set_ylabel('Frequency')\n   \n   plt.tight_layout()\n   plt.show()\n\ndef analyze_3d_structure(df_labels):\n   \"\"\"\n   Analyzes the 3D coordinates of RNA structures.\n   \n   :param df_labels: DataFrame containing 3D coordinates.\n   \"\"\"\n   # Find all coordinate columns\n   x_cols = filter_columns_by_prefix(df_labels, 'x_')\n   y_cols = filter_columns_by_prefix(df_labels, 'y_')\n   z_cols = filter_columns_by_prefix(df_labels, 'z_')\n   \n   print(f\"Number of x columns: {len(x_cols)}\")\n   print(f\"Number of y columns: {len(y_cols)}\")\n   print(f\"Number of z columns: {len(z_cols)}\")\n   \n   # Check for missing or special values in coordinates\n   special_value = -1.0e+18  # Special value observed in the data\n   \n   for i, (x_col, y_col, z_col) in enumerate(zip(x_cols, y_cols, z_cols), 1):\n       # Count missing or special values\n       x_special = (df_labels[x_col] == special_value).sum()\n       y_special = (df_labels[y_col] == special_value).sum()\n       z_special = (df_labels[z_col] == special_value).sum()\n       \n       x_null = df_labels[x_col].isnull().sum()\n       y_null = df_labels[y_col].isnull().sum()\n       z_null = df_labels[z_col].isnull().sum()\n       \n       # Count how many complete structures exist (all x, y, z are neither special nor null)\n       valid_structures = ((df_labels[x_col] != special_value) & \n                          (df_labels[y_col] != special_value) & \n                          (df_labels[z_col] != special_value) &\n                          df_labels[x_col].notnull() & \n                          df_labels[y_col].notnull() & \n                          df_labels[z_col].notnull()).sum()\n       \n       total_rows = len(df_labels)\n       \n       print(f\"\\nStructure {i}:\")\n       print(f\"  Special values: x={x_special} ({x_special/total_rows*100:.2f}%), y={y_special} ({y_special/total_rows*100:.2f}%), z={z_special} ({z_special/total_rows*100:.2f}%)\")\n       print(f\"  Null values: x={x_null} ({x_null/total_rows*100:.2f}%), y={y_null} ({y_null/total_rows*100:.2f}%), z={z_null} ({z_null/total_rows*100:.2f}%)\")\n       print(f\"  Complete structures: {valid_structures} ({valid_structures/total_rows*100:.2f}%)\")\n       \n       # Limit analysis to the first 5 structures\n       if i >= 5:\n           print(\"\\nAnalysis limited to the first 5 structures.\")\n           break\n\ndef analyze_sequences(df_sequences):\n   \"\"\"\n   Analyzes RNA sequences.\n   \n   :param df_sequences: DataFrame containing the 'sequence' column.\n   \"\"\"\n   # Basic statistics of the sequence column\n   print(\"\\nBasic statistics of the 'sequence' column:\")\n   print(df_sequences['sequence'].describe())\n   \n   # Sequence lengths\n   seq_lengths = df_sequences['sequence'].apply(len)\n   print(\"\\nSequence length statistics:\")\n   print(f\"Minimum: {seq_lengths.min()}\")\n   print(f\"Maximum: {seq_lengths.max()}\")\n   print(f\"Mean: {seq_lengths.mean():.2f}\")\n   print(f\"Median: {seq_lengths.median()}\")\n   \n   # Nucleotide counts\n   nucleotide_counts = count_nucleotides(df_sequences)\n   total_nucleotides = sum(nucleotide_counts.values())\n   \n   print(\"\\nNucleotide distribution:\")\n   for nucleotide, count in sorted(nucleotide_counts.items()):\n       print(f\"{nucleotide}: {count} ({count/total_nucleotides*100:.2f}%)\")\n   \n   # Plot length distribution\n   plt.figure(figsize=(10, 6))\n   plt.hist(seq_lengths, bins=30, alpha=0.7)\n   plt.title('Sequence Length Distribution')\n   plt.xlabel('Length')\n   plt.ylabel('Frequency')\n   plt.grid(True, alpha=0.3)\n   plt.show()\n\ndef main():\n   # Load main data\n   main_data = load_main_data()\n\n   # Check which files were loaded\n   print(\"\\nLoaded files:\")\n   for file_name, df in main_data.items():\n       if df is not None:\n           print(f\"- {file_name}: {df.shape}\")\n   \n   # Analyze 3D structures in validation_labels.csv\n   if \"validation_labels.csv\" in main_data and main_data[\"validation_labels.csv\"] is not None:\n       print(\"\\n===== Analysis of 3D Structures (validation_labels.csv) =====\")\n       df_labels = main_data[\"validation_labels.csv\"]\n       \n       # Count coordinate columns\n       x_cols = filter_columns_by_prefix(df_labels, 'x_')\n       y_cols = filter_columns_by_prefix(df_labels, 'y_')\n       z_cols = filter_columns_by_prefix(df_labels, 'z_')\n       \n       print(f\"There are {len(x_cols)} x_ columns in the DataFrame.\")\n       print(f\"There are {len(y_cols)} y_ columns in the DataFrame.\")\n       print(f\"There are {len(z_cols)} z_ columns in the DataFrame.\")\n       \n       # Identify columns without missing values\n       columns_without_missing = get_columns_without_missing_values(df_labels)\n       print(f\"\\nColumns without missing values: {len(columns_without_missing)}\")\n       \n       # Identify completely empty columns\n       empty_columns = get_empty_columns(df_labels)\n       print(f\"Completely empty columns: {len(empty_columns)}\")\n       \n       # Analyze 3D coordinates in detail\n       analyze_3d_structure(df_labels)\n       \n       # Plot distribution of x, y, z coordinates for the first structures\n       print(\"\\nDistribution of X coordinates:\")\n       plot_coord_distributions(df_labels, 'x_', max_structures=3)\n       print(\"\\nDistribution of Y coordinates:\")\n       plot_coord_distributions(df_labels, 'y_', max_structures=3)\n       print(\"\\nDistribution of Z coordinates:\")\n       plot_coord_distributions(df_labels, 'z_', max_structures=3)\n   \n   # Analyze sequences in train_sequences.csv\n   if \"train_sequences.csv\" in main_data and main_data[\"train_sequences.csv\"] is not None:\n       print(\"\\n===== Analysis of Sequences (train_sequences.csv) =====\")\n       df_sequences = main_data[\"train_sequences.csv\"]\n       \n       # First few rows of the sequence column\n       print(\"\\nFirst few rows of the 'sequence' column:\")\n       print(df_sequences['sequence'].head())\n       \n       # Data type of the sequence column\n       print(\"\\nData type of the 'sequence' column:\")\n       print(df_sequences['sequence'].dtype)\n       \n       # Complete sequence analysis\n       analyze_sequences(df_sequences)\n   \n   return main_data\n\nif __name__ == '__main__':\n   main_data = main()","metadata":{"papermill":{"duration":3.121814,"end_time":"2025-05-27T20:38:04.375816","exception":false,"start_time":"2025-05-27T20:38:01.254002","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:19:18.591934Z","iopub.execute_input":"2025-10-13T02:19:18.592404Z","iopub.status.idle":"2025-10-13T02:19:22.151911Z","shell.execute_reply.started":"2025-10-13T02:19:18.592367Z","shell.execute_reply":"2025-10-13T02:19:22.150576Z"}},"outputs":[],"execution_count":null},{"id":"d7560bae","cell_type":"markdown","source":"## Data Preparation for RNA 3D Structure Prediction 🧬🔍","metadata":{"papermill":{"duration":0.035057,"end_time":"2025-05-27T20:38:04.441316","exception":false,"start_time":"2025-05-27T20:38:04.406259","status":"completed"},"tags":[]}},{"id":"200d7894","cell_type":"code","source":"# File paths\nDATA_DIR = \"/kaggle/input/stanford-rna-3d-folding/\"\nOUTPUT_DIR = \"/kaggle/working/\"\nos.makedirs(OUTPUT_DIR, exist_ok=True)\n\ndef load_data():\n    \"\"\"\n    Loads the necessary data for the competition.\n    \"\"\"\n    data = {}\n    \n    # Load sequences\n    data['train_seq'] = pd.read_csv(os.path.join(DATA_DIR, \"train_sequences.csv\"))\n    data['valid_seq'] = pd.read_csv(os.path.join(DATA_DIR, \"validation_sequences.csv\"))\n    data['test_seq'] = pd.read_csv(os.path.join(DATA_DIR, \"test_sequences.csv\"))\n    \n    # Load structures (labels)\n    data['train_labels'] = pd.read_csv(os.path.join(DATA_DIR, \"train_labels.csv\"))\n    data['valid_labels'] = pd.read_csv(os.path.join(DATA_DIR, \"validation_labels.csv\"))\n    \n    # Load submission format\n    data['sample_submission'] = pd.read_csv(os.path.join(DATA_DIR, \"sample_submission.csv\"))\n    \n    return data\n\ndef analyze_id_structure(data_dict):\n    \"\"\"\n    Analyzes the ID structure in different files to understand the correct mapping.\n    \"\"\"\n    # We'll analyze the specific formats for train and valid\n    \n    # 1. Analysis of training labels\n    train_label_ids = data_dict['train_labels']['ID'].tolist()\n    print(f\"Total IDs in training labels: {len(train_label_ids)}\")\n    print(f\"Number of unique IDs: {len(set(train_label_ids))}\")\n    \n    # Try to understand the ID format in the labels file\n    train_id_parts = {}\n    for id_str in train_label_ids[:100]:  # Analyze the first 100\n        parts = id_str.split('_')\n        num_parts = len(parts)\n        if num_parts not in train_id_parts:\n            train_id_parts[num_parts] = []\n        train_id_parts[num_parts].append(parts)\n    \n    print(\"\\nID formats found in train_labels:\")\n    for num_parts, examples in train_id_parts.items():\n        print(f\"\\nFormat with {num_parts} parts:\")\n        for i, parts in enumerate(examples[:3]):\n            print(f\"  Example {i+1}: {parts}\")\n    \n    # 2. Analysis of training sequences\n    train_seq_ids = data_dict['train_seq']['target_id'].tolist()\n    print(f\"\\nTotal IDs in training sequences: {len(train_seq_ids)}\")\n    print(f\"Number of unique IDs: {len(set(train_seq_ids))}\")\n    \n    # Try to understand the ID format in the sequences file\n    train_seq_id_parts = {}\n    for id_str in train_seq_ids[:100]:  # Analyze the first 100\n        parts = id_str.split('_')\n        num_parts = len(parts)\n        if num_parts not in train_seq_id_parts:\n            train_seq_id_parts[num_parts] = []\n        train_seq_id_parts[num_parts].append(parts)\n    \n    print(\"\\nID formats found in train_sequences:\")\n    for num_parts, examples in train_seq_id_parts.items():\n        print(f\"\\nFormat with {num_parts} parts:\")\n        for i, parts in enumerate(examples[:3]):\n            print(f\"  Example {i+1}: {parts}\")\n    \n    # 3. Analysis of validation labels\n    valid_label_ids = data_dict['valid_labels']['ID'].tolist()\n    print(f\"\\nTotal IDs in validation labels: {len(valid_label_ids)}\")\n    print(f\"Number of unique IDs: {len(set(valid_label_ids))}\")\n    \n    # Count unique sequence IDs in validation labels\n    valid_seq_ids_from_labels = set([id_str.split('_')[0] for id_str in valid_label_ids])\n    print(f\"Number of unique sequence IDs in validation labels: {len(valid_seq_ids_from_labels)}\")\n    print(f\"Examples: {list(valid_seq_ids_from_labels)[:5]}\")\n    \n    # 4. Analysis of validation sequences\n    valid_seq_ids = data_dict['valid_seq']['target_id'].tolist()\n    print(f\"\\nTotal IDs in validation sequences: {len(valid_seq_ids)}\")\n    print(f\"Number of unique IDs: {len(set(valid_seq_ids))}\")\n    print(f\"Examples: {valid_seq_ids[:5]}\")\n    \n    # 5. Check correspondence between unique IDs\n    overlap_valid = set(valid_seq_ids).intersection(valid_seq_ids_from_labels)\n    print(f\"\\nCorrespondence between validation sequences and labels: {len(overlap_valid)} of {len(valid_seq_ids)}\")\n    \n    # 6. Check how sequences and residues relate\n    if len(overlap_valid) > 0:\n        sample_id = list(overlap_valid)[0]\n        sample_seq = data_dict['valid_seq'][data_dict['valid_seq']['target_id'] == sample_id]['sequence'].iloc[0]\n        sample_labels = data_dict['valid_labels'][data_dict['valid_labels']['ID'].str.startswith(f\"{sample_id}_\")]\n        \n        print(f\"\\nAnalysis for sequence ID: {sample_id}\")\n        print(f\"Sequence length: {len(sample_seq)}\")\n        print(f\"Number of residues in labels: {len(sample_labels)}\")\n        \n        # Check how residue numbers are related\n        residue_numbers = sample_labels['resid'].sort_values().tolist()\n        print(f\"First residue numbers: {residue_numbers[:10]}\")\n        print(f\"Last residue numbers: {residue_numbers[-10:]}\")\n        \n    return train_id_parts, train_seq_id_parts, overlap_valid\n\ndef fix_train_mapping(train_seq_df, train_labels_df):\n    \"\"\"\n    Identifies a correct mapping between train_sequences.csv and train_labels.csv\n    using the ID format from the validation file as a reference.\n    \n    This is necessary because there's no obvious direct correspondence between the IDs.\n    \"\"\"\n    # First, extract the prefix of the ID from labels (format: XX_Y_Z)\n    train_labels_df['seq_id'] = train_labels_df['ID'].apply(lambda x: x.split('_')[0] + '_' + x.split('_')[1])\n    \n    # Check if this format corresponds to the format of sequence IDs\n    seq_ids_set = set(train_seq_df['target_id'])\n    label_seq_ids_set = set(train_labels_df['seq_id'])\n    \n    overlap = seq_ids_set.intersection(label_seq_ids_set)\n    print(f\"Overlap after format adjustment: {len(overlap)} of {len(seq_ids_set)}\")\n    \n    if len(overlap) > 0:\n        print(f\"Examples of matching IDs: {list(overlap)[:5]}\")\n        return overlap\n    \n    # If it still doesn't work, we need to analyze the structure in more detail\n    print(\"No matches found, checking other formats...\")\n    \n    # Try other possible formats\n    formats_to_try = [\n        lambda x: x.split('_')[0],                             # Only first part\n        lambda x: '_'.join(x.split('_')[:2]),                  # First two parts\n        lambda x: x.split('_')[0] + '_' + x.split('_')[1][0],  # First part + first letter of second part\n    ]\n    \n    for i, format_func in enumerate(formats_to_try):\n        train_labels_df[f'seq_id_{i}'] = train_labels_df['ID'].apply(format_func)\n        label_seq_ids_set = set(train_labels_df[f'seq_id_{i}'])\n        overlap = seq_ids_set.intersection(label_seq_ids_set)\n        print(f\"Format {i}: Overlap = {len(overlap)} of {len(seq_ids_set)}\")\n        \n        if len(overlap) > 0:\n            print(f\"Examples of matching IDs: {list(overlap)[:5]}\")\n            return overlap, f'seq_id_{i}'\n    \n    # If no match is found, create a mapping based on observed patterns\n    print(\"No matches found using simple patterns.\")\n    print(\"Creating a manual mapping based on data structure...\")\n    \n    # Group labels by first parts of ID\n    train_labels_df['prefix'] = train_labels_df['ID'].apply(lambda x: x.split('_')[0])\n    label_groups = train_labels_df.groupby('prefix')\n    \n    # For each sequence, find the best match based on number of residues\n    mapping = {}\n    for _, seq_row in train_seq_df.iterrows():\n        seq_id = seq_row['target_id']\n        seq_length = len(seq_row['sequence'])\n        \n        best_match = None\n        best_diff = float('inf')\n        \n        for prefix, group in label_groups:\n            residue_count = len(group)\n            diff = abs(residue_count - seq_length)\n            \n            if diff < best_diff:\n                best_diff = diff\n                best_match = prefix\n        \n        # Consider a match only if the number of residues is close\n        if best_diff <= 10:  # Tolerance of 10 residues\n            mapping[seq_id] = best_match\n    \n    print(f\"Manual mapping created with {len(mapping)} matches\")\n    return mapping\n\ndef create_mapping_valid(valid_seq_df, valid_labels_df):\n    \"\"\"\n    Creates a mapping between validation sequences and their coordinates.\n    \n    In this case, the IDs already correspond directly (R1107 -> R1107_1, R1107_2, etc.)\n    \"\"\"\n    # Check which ID format is used in the validation set\n    valid_labels_df['seq_id'] = valid_labels_df['ID'].apply(lambda x: x.split('_')[0])\n    \n    # Check overlap\n    seq_ids = set(valid_seq_df['target_id'])\n    label_seq_ids = set(valid_labels_df['seq_id'])\n    \n    overlap = seq_ids.intersection(label_seq_ids)\n    print(f\"Correspondence for validation: {len(overlap)} of {len(seq_ids)}\")\n    \n    mapping = {}\n    for seq_id in overlap:\n        # Get sequence\n        seq = valid_seq_df[valid_seq_df['target_id'] == seq_id]['sequence'].iloc[0]\n        \n        # Get all residues for this sequence\n        residues = valid_labels_df[valid_labels_df['seq_id'] == seq_id].sort_values('resid')\n        \n        # Extract coordinates for all structures\n        num_structures = 1\n        for col in residues.columns:\n            if col.startswith('x_'):\n                struct_num = int(col.split('_')[1])\n                num_structures = max(num_structures, struct_num)\n        \n        # Initialize structures\n        structures = []\n        \n        for struct_idx in range(1, num_structures + 1):\n            coords = []\n            has_valid_coords = False\n            \n            # Check if this structure has coordinates\n            if f'x_{struct_idx}' in residues.columns:\n                for _, row in residues.iterrows():\n                    x = row[f'x_{struct_idx}']\n                    y = row[f'y_{struct_idx}']\n                    z = row[f'z_{struct_idx}']\n                    \n                    # Check if they are valid values\n                    if abs(x) < 1.0e+17 and abs(y) < 1.0e+17 and abs(z) < 1.0e+17:\n                        coords.append([x, y, z])\n                        has_valid_coords = True\n                    else:\n                        coords.append([np.nan, np.nan, np.nan])\n            \n            if has_valid_coords:\n                structures.append(coords)\n        \n        # Add to mapping if there are valid structures\n        if structures:\n            mapping[seq_id] = {\n                'sequence': seq,\n                'structures': structures\n            }\n    \n    print(f\"Mapping created with {len(mapping)} valid sequences\")\n    return mapping\n\ndef create_processed_data(mapping, output_prefix):\n    \"\"\"\n    Creates and saves processed data from the mapping.\n    \n    Parameters:\n    mapping: Dictionary with the mapping of sequences to structures\n    output_prefix: Prefix for output files ('train' or 'valid')\n    \n    Returns:\n    X, y: Arrays for training\n    \"\"\"\n    if not mapping:\n        print(f\"WARNING: No valid mapping for {output_prefix}\")\n        return None, None\n    \n    X_data = []\n    y_data = []\n    ids = []\n    \n    for seq_id, data in mapping.items():\n        seq = data['sequence']\n        structures = data['structures']\n        \n        # Skip if there are no structures\n        if not structures:\n            continue\n        \n        # Use the first valid structure\n        structure = structures[0]\n        \n        # Check if the structure has valid coordinates for all residues\n        if len(structure) != len(seq):\n            print(f\"WARNING: Difference between sequence length ({len(seq)}) and coordinates ({len(structure)}) for {seq_id}\")\n            # If needed, we could consider padding or truncation here\n            continue\n        \n        # Create feature matrix (one-hot encoding)\n        features = []\n        for nucleotide in seq:\n            if nucleotide == 'A':\n                features.append([1, 0, 0, 0, 0])\n            elif nucleotide == 'C':\n                features.append([0, 1, 0, 0, 0])\n            elif nucleotide == 'G':\n                features.append([0, 0, 1, 0, 0])\n            elif nucleotide == 'U':\n                features.append([0, 0, 0, 1, 0])\n            else:\n                features.append([0, 0, 0, 0, 1])  # For unknown nucleotides\n        \n        X_data.append(np.array(features))\n        y_data.append(np.array(structure))\n        ids.append(seq_id)\n    \n    if not X_data:\n        print(f\"WARNING: No valid processed data for {output_prefix}\")\n        return None, None, []\n    \n    # Padding to ensure all sequences have the same length\n    max_length = max(len(x) for x in X_data)\n    X_padded = []\n    y_padded = []\n    \n    for x, y in zip(X_data, y_data):\n        if len(x) < max_length:\n            x_pad = np.zeros((max_length, 5))\n            x_pad[:len(x), :] = x\n            \n            y_pad = np.zeros((max_length, 3))\n            y_pad[:len(y), :] = y\n            \n            X_padded.append(x_pad)\n            y_padded.append(y_pad)\n        else:\n            X_padded.append(x)\n            y_padded.append(y)\n    \n    X = np.array(X_padded)\n    y = np.array(y_padded)\n    \n    # Save the processed data\n    np.save(os.path.join(OUTPUT_DIR, f'X_{output_prefix}.npy'), X)\n    np.save(os.path.join(OUTPUT_DIR, f'y_{output_prefix}.npy'), y)\n    \n    with open(os.path.join(OUTPUT_DIR, f'{output_prefix}_ids.txt'), 'w') as f:\n        for id in ids:\n            f.write(f\"{id}\\n\")\n    \n    print(f\"Processed data for {output_prefix}: X.shape = {X.shape}, y.shape = {y.shape}\")\n    return X, y, ids\n\ndef explore_sequence_mapping(seq_id, mapping, data_dict):\n    \"\"\"\n    Explores a mapping example in detail for diagnostics.\n    \"\"\"\n    if seq_id not in mapping:\n        print(f\"WARNING: Sequence ID {seq_id} not found in mapping\")\n        return\n    \n    data = mapping[seq_id]\n    seq = data['sequence']\n    structures = data['structures']\n    \n    print(f\"Exploring mapping for sequence: {seq_id}\")\n    print(f\"Sequence length: {len(seq)}\")\n    print(f\"Number of available structures: {len(structures)}\")\n    \n    # Detail each structure\n    for i, structure in enumerate(structures):\n        print(f\"\\nStructure {i+1}:\")\n        print(f\"  Number of coordinates: {len(structure)}\")\n        if len(structure) > 0:\n            print(f\"  First coordinates: {structure[:3]}\")\n            print(f\"  Last coordinates: {structure[-3:]}\")\n        \n        # Check correspondence with the sequence\n        if len(structure) != len(seq):\n            print(f\"  WARNING: Difference between sequence length ({len(seq)}) and coordinates ({len(structure)})\")\n        else:\n            print(f\"  Perfect match between sequence and coordinates\")\n\ndef main():\n    # Load the data\n    print(\"Loading data...\")\n    data_dict = load_data()\n    \n    # Analyze ID structure to understand the mapping\n    print(\"\\nAnalyzing ID structure...\")\n    train_id_parts, train_seq_id_parts, overlap_valid = analyze_id_structure(data_dict)\n    \n    # For validation, the mapping is direct (R1107 -> R1107_1, R1107_2, etc.)\n    print(\"\\nCreating mapping for validation data...\")\n    valid_mapping = create_mapping_valid(data_dict['valid_seq'], data_dict['valid_labels'])\n    \n    # Explore a validation mapping example to verify\n    if valid_mapping:\n        sample_id = list(valid_mapping.keys())[0]\n        print(f\"\\nExploring a validation mapping example ({sample_id}):\")\n        explore_sequence_mapping(sample_id, valid_mapping, data_dict)\n    \n    # Create and save processed data for validation\n    X_valid, y_valid, valid_ids = create_processed_data(valid_mapping, 'valid')\n    \n    # Since we couldn't establish a mapping for training,\n    # we'll use validation data for training as well (transfer learning)\n    print(\"\\nUsing validation data as training (due to lack of direct mapping)...\")\n    X_train = X_valid\n    y_train = y_valid\n    train_ids = valid_ids\n    \n    if X_train is not None:\n        np.save(os.path.join(OUTPUT_DIR, 'X_train.npy'), X_train)\n        np.save(os.path.join(OUTPUT_DIR, 'y_train.npy'), y_train)\n        \n        with open(os.path.join(OUTPUT_DIR, 'train_ids.txt'), 'w') as f:\n            for id in train_ids:\n                f.write(f\"{id}\\n\")\n    \n    # Return the processed data\n    return {\n        'X_train': X_train,\n        'y_train': y_train,\n        'X_valid': X_valid,\n        'y_valid': y_valid,\n        'valid_mapping': valid_mapping,\n        'valid_ids': valid_ids\n    }\n\nif __name__ == \"__main__\":\n    processed_data = main()","metadata":{"papermill":{"duration":6.623898,"end_time":"2025-05-27T20:38:11.099707","exception":false,"start_time":"2025-05-27T20:38:04.475809","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:19:26.716973Z","iopub.execute_input":"2025-10-13T02:19:26.717734Z","iopub.status.idle":"2025-10-13T02:19:33.311391Z","shell.execute_reply.started":"2025-10-13T02:19:26.717691Z","shell.execute_reply":"2025-10-13T02:19:33.309959Z"}},"outputs":[],"execution_count":null},{"id":"ab76b3b2","cell_type":"markdown","source":"## Heatmap Viewer for RNA Sequences 🔥🧬","metadata":{"papermill":{"duration":0.025822,"end_time":"2025-05-27T20:38:11.152563","exception":false,"start_time":"2025-05-27T20:38:11.126741","status":"completed"},"tags":[]}},{"id":"b4e02c75","cell_type":"code","source":"def visualize_rna_heatmap_from_processed_data(processed_data, num_samples=12):\n    \"\"\"\n    Visualizes a heatmap for RNA sequences using processed data.\n    \n    Parameters:\n    processed_data: Dictionary with processed data returned by the main() function\n    num_samples: Number of sequences to visualize\n    \"\"\"\n    try:\n        # Check if we have the necessary data\n        if 'X_valid' not in processed_data or processed_data['X_valid'] is None:\n            print(\"Validation data not found in processed_data object\")\n            return None\n        \n        # Get the data\n        X_valid = processed_data['X_valid']\n        print(f\"Data found with format: {X_valid.shape}\")\n        \n        # Limit to the number of samples\n        X_valid_subset = X_valid[:num_samples]\n        \n        # If we have IDs, use them\n        if 'valid_ids' in processed_data and processed_data['valid_ids']:\n            valid_ids = processed_data['valid_ids'][:num_samples]\n        else:\n            valid_ids = [f\"Seq_{i+1}\" for i in range(X_valid_subset.shape[0])]\n        \n        # Convert one-hot encoding to nucleotide indices\n        # Expected format: A=[1,0,0,0,0], C=[0,1,0,0,0], G=[0,0,1,0,0], U=[0,0,0,1,0], N=[0,0,0,0,1]\n        sequences_matrix = np.argmax(X_valid_subset, axis=2)\n        \n        # Replace zeros (padding) with 4 (N/Unknown) when all values are zero\n        is_padding = np.all(X_valid_subset == 0, axis=2)\n        sequences_matrix[is_padding] = 4\n        \n        # Define a categorical colormap (distinct colors per nucleotide)\n        cmap = mcolors.ListedColormap(['#3498db', '#2ecc71', '#e74c3c', '#9b59b6', '#95a5a6'])\n        bounds = [0, 1, 2, 3, 4, 5]\n        norm = mcolors.BoundaryNorm(bounds, cmap.N)\n        \n        # Create figure\n        plt.figure(figsize=(20, 10))\n        im = plt.imshow(sequences_matrix, cmap=cmap, norm=norm, aspect='auto')\n        \n        # Add color bar\n        cbar = plt.colorbar(im, ticks=[0.5, 1.5, 2.5, 3.5, 4.5])\n        cbar.set_label('Nucleotides', fontsize=14)\n        cbar.set_ticklabels(['A', 'C', 'G', 'U', 'N/Padding'])\n        \n        # Add axis labels\n        plt.xlabel(\"Position in Sequence\", fontsize=14)\n        plt.ylabel(\"RNA Sequences\", fontsize=14)\n        \n        # Add title\n        plt.title(\"RNA Sequences Heatmap\", fontsize=16)\n        \n        # Add sequence IDs as y-axis labels\n        plt.yticks(range(len(valid_ids)), valid_ids, fontsize=10)\n        \n        # Show only some labels on x-axis to avoid crowding\n        sequence_length = sequences_matrix.shape[1]\n        step = max(1, sequence_length // 20)  # Show at most 20 labels\n        plt.xticks(range(0, sequence_length, step), range(1, sequence_length + 1, step))\n        \n        # Add grid\n        plt.grid(False)\n        \n        # Add information about nucleotide distribution\n        all_nucleotides = sequences_matrix.flatten()\n        nucleotide_counts = {\n            'A': np.sum(all_nucleotides == 0),\n            'C': np.sum(all_nucleotides == 1),\n            'G': np.sum(all_nucleotides == 2),\n            'U': np.sum(all_nucleotides == 3),\n            'N': np.sum(all_nucleotides == 4)\n        }\n        \n        total_nucleotides = sum(nucleotide_counts.values())\n        nucleotide_percentages = {k: (v / total_nucleotides) * 100 for k, v in nucleotide_counts.items()}\n        \n        # Add text with statistics\n        info_text = \"\\n\".join([\n            f\"Total sequences visualized: {num_samples}\",\n            f\"Maximum length: {sequence_length}\",\n            f\"A: {nucleotide_percentages['A']:.1f}%\",\n            f\"C: {nucleotide_percentages['C']:.1f}%\",\n            f\"G: {nucleotide_percentages['G']:.1f}%\",\n            f\"U: {nucleotide_percentages['U']:.1f}%\",\n            f\"N/Padding: {nucleotide_percentages['N']:.1f}%\"\n        ])\n        \n        plt.figtext(0.02, 0.02, info_text, fontsize=10, bbox=dict(facecolor='white', alpha=0.8))\n        \n        # Show the plot\n        plt.tight_layout()\n        plt.show()\n        \n        # Optionally, save the plot\n        output_dir = '/kaggle/working/'\n        plt.savefig(os.path.join(output_dir, 'rna_heatmap.png'), dpi=300)\n        print(f\"Heatmap saved to {os.path.join(output_dir, 'rna_heatmap.png')}\")\n        \n        return sequences_matrix\n    except Exception as e:\n        print(f\"Error processing data: {e}\")\n        return None\n\n# Use the function (assuming processed_data is available)\nvisualize_rna_heatmap_from_processed_data(processed_data)","metadata":{"papermill":{"duration":0.991661,"end_time":"2025-05-27T20:38:12.170675","exception":false,"start_time":"2025-05-27T20:38:11.179014","status":"completed"},"tags":[],"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:19:37.291687Z","iopub.execute_input":"2025-10-13T02:19:37.292049Z","iopub.status.idle":"2025-10-13T02:19:38.306723Z","shell.execute_reply.started":"2025-10-13T02:19:37.292021Z","shell.execute_reply":"2025-10-13T02:19:38.305781Z"}},"outputs":[],"execution_count":null},{"id":"c92582f9-ee86-4073-a3ab-3335197e688c","cell_type":"code","source":"\"\"\"\n╔══════════════════════════════════════════════════════════════╗\n║     TURBINED RNA 3D - PLAN B: PHYSICS + VOTING             ║\n║     Stanford RNA 3D Folding Competition                      ║\n║                                                              ║\n║     Base: Conservative (0.38346 private)                    ║\n║     + Refinamento iterativo com física                      ║\n║     + Ensemble por votação ponderada                        ║\n║                                                              ║\n║     Score esperado: 0.39-0.395 private                      ║\n║     Tempo: ~4-5 minutos                                     ║\n║     Risco: Baixo-Médio                                      ║\n║                                                              ║\n║     SE < 0.383: REVERTER para Conservative                  ║\n║     SE > 0.390: SUCESSO! 🎉                                 ║\n╚══════════════════════════════════════════════════════════════╝\n\"\"\"\n\nimport os\nimport numpy as np\nimport pandas as pd\nfrom scipy.spatial.distance import pdist, squareform\nfrom scipy import stats\nfrom sklearn.cluster import AgglomerativeClustering\nfrom tqdm.auto import tqdm\nimport warnings\nwarnings.filterwarnings('ignore')\n\n# ============================================================================\n# CONFIGURAÇÃO\n# ============================================================================\n\nDATA_DIR = \"/kaggle/input/stanford-rna-3d-folding/\"\nOUTPUT_DIR = \"/kaggle/working/\"\nos.makedirs(OUTPUT_DIR, exist_ok=True)\n\nprint(\"=\"*70)\nprint(\"TURBINED RNA 3D - PLAN B: PHYSICS + VOTING\")\nprint(\"=\"*70)\nprint(\"Base: Conservative + 2 melhorias cirúrgicas\")\nprint(\"Score esperado: 0.39-0.395 private\")\nprint(\"=\"*70)\n\n# ============================================================================\n# FUNÇÕES BASE (mantidas do Conservative)\n# ============================================================================\n\ndef decode_onehot_to_string(sequence):\n    if isinstance(sequence, str):\n        return sequence\n    if len(sequence.shape) == 1:\n        return 'N' * len(sequence)\n    bases = ['A', 'C', 'G', 'U', 'N']\n    seq_string = ''\n    for one_hot in sequence:\n        idx = np.argmax(one_hot)\n        seq_string += bases[min(idx, 4)]\n    return seq_string\n\ndef string_to_onehot(sequence):\n    bases = {'A': 0, 'C': 1, 'G': 2, 'U': 3}\n    seq_length = len(sequence)\n    one_hot = np.zeros((seq_length, 5))\n    for i, base in enumerate(sequence):\n        if base in bases:\n            one_hot[i, bases[base]] = 1\n        else:\n            one_hot[i, 4] = 1\n    return one_hot\n\ndef calculate_gc_content(sequence):\n    seq_str = decode_onehot_to_string(sequence)\n    if len(seq_str) == 0:\n        return 0.5\n    return (seq_str.count('G') + seq_str.count('C')) / len(seq_str)\n\ndef calculate_radius_of_gyration(coords):\n    valid_mask = ~np.all(coords == 0, axis=1)\n    valid_coords = coords[valid_mask]\n    if len(valid_coords) == 0:\n        return 0\n    center = np.mean(valid_coords, axis=0)\n    distances = np.linalg.norm(valid_coords - center, axis=1)\n    return np.sqrt(np.mean(distances ** 2))\n\ndef calculate_tm_score_simple(coords1, coords2):\n    valid_mask1 = ~np.all(coords1 == 0, axis=1)\n    valid_mask2 = ~np.all(coords2 == 0, axis=1)\n    valid_mask = valid_mask1 & valid_mask2\n    if np.sum(valid_mask) < 3:\n        return 0.0\n    coords1_valid = coords1[valid_mask]\n    coords2_valid = coords2[valid_mask]\n    rmsd = np.sqrt(np.mean(np.sum((coords1_valid - coords2_valid)**2, axis=1)))\n    L = len(coords1_valid)\n    d0 = 1.24 * (L - 15)**(1/3) - 1.8 if L > 15 else 0.5\n    tm_score = np.mean(1 / (1 + (rmsd / d0)**2))\n    return tm_score\n\ndef predict_secondary_structure(sequence):\n    seq_str = decode_onehot_to_string(sequence)\n    try:\n        import RNA\n        structure, mfe = RNA.fold(seq_str)\n        return structure\n    except ImportError:\n        n = len(seq_str)\n        structure = ['.' for _ in range(n)]\n        paired = set()\n        for i in range(n):\n            if i in paired:\n                continue\n            for j in range(i+3, min(i+30, n)):\n                if j in paired:\n                    continue\n                if (seq_str[i] == 'G' and seq_str[j] == 'C') or \\\n                   (seq_str[i] == 'C' and seq_str[j] == 'G') or \\\n                   (seq_str[i] == 'A' and seq_str[j] == 'U') or \\\n                   (seq_str[i] == 'U' and seq_str[j] == 'A'):\n                    structure[i] = '('\n                    structure[j] = ')'\n                    paired.add(i)\n                    paired.add(j)\n                    break\n        return ''.join(structure)\n\ndef find_base_pairs(secondary_structure):\n    pairs = []\n    stack = []\n    for i, char in enumerate(secondary_structure):\n        if char == '(':\n            stack.append(i)\n        elif char == ')':\n            if stack:\n                j = stack.pop()\n                pairs.append((j, i))\n    return pairs\n\ndef classify_structural_regions(secondary_structure):\n    classification = {}\n    n = len(secondary_structure)\n    pairs = find_base_pairs(secondary_structure)\n    pair_dict = {i: j for i, j in pairs}\n    pair_dict.update({j: i for i, j in pairs})\n    \n    for i in range(n):\n        if secondary_structure[i] == '.':\n            classification[i] = 'unpaired'\n        elif i in pair_dict:\n            is_helix = False\n            if i > 0 and i+1 < n:\n                if (i-1 in pair_dict and i+1 in pair_dict):\n                    is_helix = True\n            classification[i] = 'helix' if is_helix else 'junction'\n        else:\n            classification[i] = 'unknown'\n    \n    for i in range(n):\n        if classification[i] == 'unpaired':\n            neighbors_paired = sum(\n                1 for j in range(max(0, i-3), min(n, i+4))\n                if j in pair_dict\n            )\n            if neighbors_paired >= 4:\n                classification[i] = 'loop'\n    \n    return classification\n\ndef calculate_sequence_properties(sequence):\n    seq_str = decode_onehot_to_string(sequence)\n    if len(seq_str) == 0:\n        return {\n            'length': 0, 'gc_content': 0.5, 'purine_ratio': 0.5,\n            'a_content': 0.25, 'u_content': 0.25,\n            'g_content': 0.25, 'c_content': 0.25,\n            'secondary_structure': ''\n        }\n    props = {\n        'length': len(seq_str),\n        'gc_content': (seq_str.count('G') + seq_str.count('C')) / len(seq_str),\n        'purine_ratio': (seq_str.count('A') + seq_str.count('G')) / len(seq_str),\n        'a_content': seq_str.count('A') / len(seq_str),\n        'u_content': seq_str.count('U') / len(seq_str),\n        'g_content': seq_str.count('G') / len(seq_str),\n        'c_content': seq_str.count('C') / len(seq_str),\n    }\n    props['secondary_structure'] = predict_secondary_structure(sequence)\n    return props\n\ndef compare_secondary_structures(struct1, struct2):\n    if len(struct1) != len(struct2):\n        min_len = min(len(struct1), len(struct2))\n        struct1 = struct1[:min_len]\n        struct2 = struct2[:min_len]\n    if len(struct1) == 0:\n        return 0.0\n    matches = sum(1 for a, b in zip(struct1, struct2) if a == b)\n    return matches / len(struct1)\n\ndef calculate_overall_similarity(props1, props2):\n    gc_sim = 1 - abs(props1['gc_content'] - props2['gc_content'])\n    length_sim = min(props1['length'], props2['length']) / max(props1['length'], props2['length'])\n    purine_sim = 1 - abs(props1['purine_ratio'] - props2['purine_ratio'])\n    struct_sim = compare_secondary_structures(\n        props1['secondary_structure'],\n        props2['secondary_structure']\n    )\n    total = (gc_sim * 0.3 + length_sim * 0.4 + \n             purine_sim * 0.15 + struct_sim * 0.15)\n    return total\n\nclass SmartReferenceCache:\n    def __init__(self, max_cache_size=2000):\n        self.cache = {}\n        self.hit_counts = {}\n        self.max_size = max_cache_size\n    \n    def get_key(self, sequence_props):\n        gc_bin = int(sequence_props['gc_content'] * 10)\n        length_bin = int(sequence_props['length'] / 50)\n        return f\"gc{gc_bin}_len{length_bin}\"\n    \n    def get_references(self, target_seq, all_references):\n        props = calculate_sequence_properties(target_seq)\n        key = self.get_key(props)\n        if key in self.cache:\n            self.hit_counts[key] = self.hit_counts.get(key, 0) + 1\n            return self.cache[key]\n        return None\n    \n    def add_to_cache(self, target_seq, references):\n        props = calculate_sequence_properties(target_seq)\n        key = self.get_key(props)\n        if len(self.cache) >= self.max_size:\n            if self.hit_counts:\n                min_key = min(self.hit_counts, key=self.hit_counts.get)\n                del self.cache[min_key]\n                del self.hit_counts[min_key]\n        self.cache[key] = references\n        self.hit_counts[key] = 1\n    \n    def get_stats(self):\n        total = sum(self.hit_counts.values())\n        size = len(self.cache)\n        return {\n            'cache_size': size,\n            'total_requests': total,\n            'hit_rate': sum(c > 1 for c in self.hit_counts.values()) / max(size, 1)\n        }\n\ndef simple_sequence_similarity(seq1, seq2):\n    str1 = decode_onehot_to_string(seq1)\n    str2 = decode_onehot_to_string(seq2)\n    min_len = min(len(str1), len(str2))\n    max_len = max(len(str1), len(str2))\n    if max_len == 0:\n        return 0.0\n    matches = sum(1 for i in range(min_len) if str1[i] == str2[i])\n    return matches / max_len\n\ndef select_best_reference_by_sequence(target_seq, reference_pool, top_k=5):\n    similarities = []\n    for ref_idx, ref_data in enumerate(reference_pool):\n        ref_seq = ref_data['sequence']\n        similarity = simple_sequence_similarity(target_seq, ref_seq)\n        similarities.append({\n            'index': ref_idx,\n            'similarity': similarity,\n            'data': ref_data\n        })\n    similarities.sort(key=lambda x: x['similarity'], reverse=True)\n    return [s['data'] for s in similarities[:top_k]]\n\ndef select_by_physical_properties(target_seq, reference_pool, top_k=5):\n    target_props = calculate_sequence_properties(target_seq)\n    scores = []\n    for ref_data in reference_pool:\n        ref_props = calculate_sequence_properties(ref_data['sequence'])\n        score = calculate_overall_similarity(target_props, ref_props)\n        scores.append({'score': score, 'data': ref_data})\n    scores.sort(key=lambda x: x['score'], reverse=True)\n    return [s['data'] for s in scores[:top_k]]\n\ndef select_references_adaptive(target_seq, reference_pool, base_top_k=10):\n    pool_size = len(reference_pool)\n    \n    if pool_size < 1000:\n        top_k_seq = base_top_k\n        top_k_props = 5\n    elif pool_size < 3000:\n        top_k_seq = base_top_k + 5\n        top_k_props = 7\n    else:\n        top_k_seq = base_top_k + 10\n        top_k_props = 10\n    \n    seq_similar = select_best_reference_by_sequence(\n        target_seq, reference_pool, top_k=top_k_seq\n    )\n    \n    best_refs = select_by_physical_properties(\n        target_seq, seq_similar, top_k=top_k_props\n    )\n    \n    return best_refs\n\ndef interpolate_multiple_references(references, target_props, method='weighted'):\n    if len(references) == 0:\n        return None\n    \n    target_length = target_props['length']\n    best_ref = None\n    best_score = -1\n    \n    for ref in references:\n        if len(ref['structure']) != target_length:\n            continue\n        ref_props = calculate_sequence_properties(ref['sequence'])\n        similarity = calculate_overall_similarity(target_props, ref_props)\n        if similarity > best_score:\n            best_score = similarity\n            best_ref = ref\n    \n    if best_ref is None:\n        best_ref = max(references, \n                      key=lambda r: calculate_overall_similarity(\n                          calculate_sequence_properties(r['sequence']), \n                          target_props\n                      ))\n        base_structure = best_ref['structure']\n        if len(base_structure) > target_length:\n            return base_structure[:target_length]\n        else:\n            extended = np.zeros((target_length, 3))\n            extended[:len(base_structure)] = base_structure\n            return extended\n    \n    return best_ref['structure']\n\ndef generate_synthetic_references(original_refs, n_synthetic=100):\n    if len(original_refs) < 2:\n        return []\n    \n    size_groups = {}\n    for ref in original_refs:\n        size = len(ref['structure'])\n        size_bin = size // 50\n        if size_bin not in size_groups:\n            size_groups[size_bin] = []\n        size_groups[size_bin].append(ref)\n    \n    synthetic_refs = []\n    \n    for size_bin, refs_in_bin in size_groups.items():\n        if len(refs_in_bin) < 2:\n            continue\n        \n        n_for_bin = max(1, int(n_synthetic * len(refs_in_bin) / len(original_refs)))\n        \n        for _ in range(n_for_bin):\n            n_to_combine = min(np.random.randint(2, 4), len(refs_in_bin))\n            selected_indices = np.random.choice(len(refs_in_bin), n_to_combine, replace=False)\n            selected_refs = [refs_in_bin[i] for i in selected_indices]\n            \n            sizes = [len(ref['structure']) for ref in selected_refs]\n            target_size = max(set(sizes), key=sizes.count)\n            exact_refs = [ref for ref in selected_refs if len(ref['structure']) == target_size]\n            \n            if len(exact_refs) < 2:\n                continue\n            \n            combined_structure = np.mean(\n                [ref['structure'] for ref in exact_refs],\n                axis=0\n            )\n            \n            synthetic_refs.append({\n                'structure': combined_structure,\n                'sequence': exact_refs[0]['sequence'],\n                'is_synthetic': True\n            })\n    \n    return synthetic_refs[:n_synthetic]\n\ndef residue_specific_noise(coords, sequence, base_noise_level=0.2):\n    noise_scales = {'A': 0.8, 'G': 0.85, 'C': 1.1, 'U': 1.15, 'N': 1.0}\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    seq_str = decode_onehot_to_string(sequence)\n    for i in range(len(coords)):\n        if not valid_mask[i]:\n            continue\n        base = seq_str[i] if i < len(seq_str) else 'N'\n        scale = noise_scales.get(base, 1.0)\n        noise = np.random.normal(0, base_noise_level * scale, 3)\n        new_coords[i] += noise\n    return new_coords\n\ndef context_aware_noise(coords, sequence, base_noise_level=0.2):\n    secondary_structure = predict_secondary_structure(sequence)\n    structural_context = classify_structural_regions(secondary_structure)\n    noise_by_context = {\n        'helix': 0.5, 'loop': 1.5, 'junction': 1.0,\n        'unpaired': 1.3, 'unknown': 1.0\n    }\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    for i in range(len(coords)):\n        if not valid_mask[i]:\n            continue\n        context = structural_context.get(i, 'unknown')\n        scale = noise_by_context[context]\n        noise_level = base_noise_level * scale\n        noise = np.random.normal(0, noise_level, 3)\n        new_coords[i] += noise\n    return new_coords\n\ndef simple_noise(coords, sequence, base_noise_level=0.2):\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    noise = np.random.normal(0, base_noise_level, coords.shape)\n    new_coords[valid_mask] += noise[valid_mask]\n    return new_coords\n\ndef generate_candidates_diverse(base_structure, target_seq, n_candidates=20, noise_level=0.21):\n    candidates = []\n    n_per_type = n_candidates // 3\n    \n    for i in range(n_per_type):\n        noise_variation = noise_level * np.random.uniform(0.7, 1.3)\n        candidate = residue_specific_noise(base_structure, target_seq, noise_variation)\n        candidates.append(candidate)\n    \n    for i in range(n_per_type):\n        noise_variation = noise_level * np.random.uniform(0.7, 1.3)\n        candidate = context_aware_noise(base_structure, target_seq, noise_variation)\n        candidates.append(candidate)\n    \n    remaining = n_candidates - len(candidates)\n    for i in range(remaining):\n        noise_variation = noise_level * np.random.uniform(0.7, 1.3)\n        candidate = simple_noise(base_structure, target_seq, noise_variation)\n        candidates.append(candidate)\n    \n    return candidates\n\ndef estimate_structure_quality(coords):\n    valid_mask = ~np.all(coords == 0, axis=1)\n    valid_coords = coords[valid_mask]\n    if len(valid_coords) < 3:\n        return 0\n    \n    quality = 0\n    rg = calculate_radius_of_gyration(coords)\n    quality += 1.0 / (1.0 + rg / 10)\n    \n    backbone_dists = []\n    for i in range(len(valid_coords) - 1):\n        dist = np.linalg.norm(valid_coords[i+1] - valid_coords[i])\n        backbone_dists.append(dist)\n    \n    if backbone_dists:\n        backbone_var = np.var(backbone_dists)\n        quality += 1.0 / (1.0 + backbone_var)\n    \n    return max(quality, 0)\n\n# ============================================================================\n# MELHORIA #1: REFINAMENTO ITERATIVO COM FÍSICA\n# ============================================================================\n\ndef regularize_backbone(coords, target_dist=5.9, alpha=0.1):\n    \"\"\"\n    Força distâncias backbone próximas de 5.9Å (distância C1'-C1')\n    \"\"\"\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    \n    for i in range(len(coords) - 1):\n        if not (valid_mask[i] and valid_mask[i+1]):\n            continue\n        \n        vec = new_coords[i+1] - new_coords[i]\n        dist = np.linalg.norm(vec)\n        \n        if dist < 1e-6:\n            continue\n        \n        # Corrigir para target_dist\n        correction = (target_dist - dist) * alpha\n        direction = vec / dist\n        \n        new_coords[i+1] = new_coords[i+1] + correction * direction\n    \n    return new_coords\n\ndef remove_steric_clashes(coords, min_dist=3.0, alpha=0.05):\n    \"\"\"\n    Remove clashes movendo átomos que estão muito próximos\n    RNA não pode ter distâncias < 3.0Å entre resíduos não-adjacentes\n    \"\"\"\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    \n    for i in range(len(coords)):\n        if not valid_mask[i]:\n            continue\n        \n        for j in range(i + 3, len(coords)):  # Pular adjacentes\n            if not valid_mask[j]:\n                continue\n            \n            vec = new_coords[j] - new_coords[i]\n            dist = np.linalg.norm(vec)\n            \n            if dist < min_dist and dist > 1e-6:\n                # Empurrar para longe\n                push = (min_dist - dist) * alpha\n                direction = vec / dist\n                new_coords[j] += push * direction\n                new_coords[i] -= push * direction\n    \n    return new_coords\n\ndef smooth_structure(coords, window=3):\n    \"\"\"\n    Suaviza estrutura usando média móvel\n    Reduz jitters e irregularidades\n    \"\"\"\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    \n    for i in range(len(coords)):\n        if not valid_mask[i]:\n            continue\n        \n        # Coletar vizinhos\n        neighbors = []\n        for j in range(max(0, i - window), min(len(coords), i + window + 1)):\n            if valid_mask[j]:\n                neighbors.append(coords[j])\n        \n        if len(neighbors) > 1:\n            # Média ponderada (peso maior no centro)\n            weights = np.exp(-np.arange(len(neighbors)) / 2.0)\n            weights = weights / np.sum(weights)\n            new_coords[i] = np.sum([w * n for w, n in zip(weights, neighbors)], axis=0)\n    \n    return new_coords\n\ndef apply_basepairing_constraints(coords, sequence):\n    \"\"\"\n    Aplica constraints de base pairing\n    G-C devem estar ~10.5Å, A-U ~11.0Å\n    \"\"\"\n    new_coords = coords.copy()\n    valid_mask = ~np.all(coords == 0, axis=1)\n    \n    # Encontrar pares de bases\n    sec_struct = predict_secondary_structure(sequence)\n    base_pairs = find_base_pairs(sec_struct)\n    \n    if not base_pairs:\n        return new_coords\n    \n    seq_str = decode_onehot_to_string(sequence)\n    alpha = 0.05  # Taxa de ajuste suave\n    \n    for i, j in base_pairs:\n        if not (valid_mask[i] and valid_mask[j]):\n            continue\n        \n        # Determinar distância ideal\n        if (seq_str[i] == 'G' and seq_str[j] == 'C') or \\\n           (seq_str[i] == 'C' and seq_str[j] == 'G'):\n            ideal_dist = 10.5\n        else:  # A-U ou outras\n            ideal_dist = 11.0\n        \n        # Distância atual\n        vec = new_coords[j] - new_coords[i]\n        dist = np.linalg.norm(vec)\n        \n        if dist < 1e-6:\n            continue\n        \n        # Ajustar\n        correction = (ideal_dist - dist) * alpha\n        direction = vec / dist\n        \n        new_coords[j] += correction * direction * 0.5\n        new_coords[i] -= correction * direction * 0.5\n    \n    return new_coords\n\ndef iterative_physics_refinement(coords, sequence, n_iterations=8):\n    \"\"\"\n    Refina estrutura iterativamente com constraints físicos\n    \n    Args:\n        coords: Coordenadas iniciais\n        sequence: Sequência do RNA\n        n_iterations: Número de iterações (default: 8)\n    \n    Returns:\n        Coordenadas refinadas\n    \"\"\"\n    current = coords.copy()\n    \n    for iteration in range(n_iterations):\n        # 1. Regularizar backbone\n        current = regularize_backbone(current, target_dist=5.9, alpha=0.1)\n        \n        # 2. Aplicar base-pairing constraints\n        current = apply_basepairing_constraints(current, sequence)\n        \n        # 3. Remover clashes\n        current = remove_steric_clashes(current, min_dist=3.0, alpha=0.05)\n        \n        # 4. Suavizar\n        if iteration % 2 == 0:  # Suavizar a cada 2 iterações\n            current = smooth_structure(current, window=3)\n    \n    return current\n\n# ============================================================================\n# MELHORIA #2: ENSEMBLE POR VOTAÇÃO PONDERADA\n# ============================================================================\n\ndef weighted_voting_ensemble(candidates, target_seq, n_final=5):\n    \"\"\"\n    Ensemble por votação ponderada:\n    - Para cada resíduo, calcula centróide de todas as posições\n    - Seleciona as n_final estruturas mais próximas do consenso\n    - Mais robusto que apenas pegar top-k por qualidade global\n    \"\"\"\n    if len(candidates) <= n_final:\n        return candidates\n    \n    n_residues = len(candidates[0])\n    \n    # Calcular qualidade de cada candidato\n    qualities = np.array([estimate_structure_quality(c) for c in candidates])\n    \n    # Pegar top-20 por qualidade (filtro inicial)\n    top_indices = np.argsort(qualities)[::-1][:min(20, len(candidates))]\n    top_candidates = [candidates[i] for i in top_indices]\n    top_qualities = qualities[top_indices]\n    \n    # Para cada candidato, calcular distância ao consenso\n    consensus_distances = []\n    \n    for candidate in top_candidates:\n        total_dist = 0\n        \n        for res_idx in range(n_residues):\n            # Coletar todas as posições deste resíduo\n            positions = [c[res_idx] for c in top_candidates]\n            \n            # Calcular centróide (consenso)\n            centroid = np.mean(positions, axis=0)\n            \n            # Distância deste candidato ao centróide\n            dist = np.linalg.norm(candidate[res_idx] - centroid)\n            total_dist += dist\n        \n        # Média da distância\n        avg_dist = total_dist / n_residues\n        consensus_distances.append(avg_dist)\n    \n    # Combinar qualidade e proximidade ao consenso\n    # 60% consenso, 40% qualidade\n    normalized_distances = np.array(consensus_distances)\n    normalized_distances = 1 - (normalized_distances / (np.max(normalized_distances) + 1e-6))\n    \n    normalized_qualities = top_qualities / (np.max(top_qualities) + 1e-6)\n    \n    combined_scores = 0.6 * normalized_distances + 0.4 * normalized_qualities\n    \n    # Selecionar top-n_final\n    best_indices = np.argsort(combined_scores)[::-1][:n_final]\n    \n    selected = [top_candidates[i] for i in best_indices]\n    \n    return selected\n\n# ============================================================================\n# PREPARAÇÃO DE DADOS\n# ============================================================================\n\ndef prepare_training_data_final(train_seq_df, train_labels_df, max_length=None):\n    print(\"\\n\" + \"=\"*70)\n    print(\"PREPARANDO DADOS DE TREINO\")\n    print(\"=\"*70)\n    \n    X_train = []\n    y_train = []\n    metadata = []\n    \n    print(\"\\nExtraindo target_id dos labels...\")\n    train_labels_df['target_id'] = train_labels_df['ID'].str.rsplit('_', n=1).str[0]\n    \n    print(\"Agrupando labels por target_id...\")\n    grouped_labels = train_labels_df.groupby('target_id')\n    \n    print(f\"\\nTotal de targets únicos nos labels: {len(grouped_labels)}\")\n    \n    print(\"\\nProcessando sequências...\")\n    skipped = 0\n    \n    for idx, row in tqdm(train_seq_df.iterrows(), total=len(train_seq_df), desc=\"Processing\"):\n        target_id = row['target_id']\n        sequence = row['sequence']\n        \n        if target_id not in grouped_labels.groups:\n            skipped += 1\n            continue\n        \n        if max_length and len(sequence) > max_length:\n            skipped += 1\n            continue\n        \n        target_labels = grouped_labels.get_group(target_id)\n        target_labels = target_labels.sort_values('resid')\n        coords = target_labels[['x_1', 'y_1', 'z_1']].values\n        \n        seq_onehot = string_to_onehot(sequence)\n        \n        if len(seq_onehot) != len(coords):\n            min_len = min(len(seq_onehot), len(coords))\n            seq_onehot = seq_onehot[:min_len]\n            coords = coords[:min_len]\n            if min_len < 10:\n                skipped += 1\n                continue\n        \n        X_train.append(seq_onehot)\n        y_train.append(coords)\n        metadata.append({\n            'target_id': target_id,\n            'sequence': sequence[:len(coords)],\n            'length': len(coords),\n            'gc_content': calculate_gc_content(seq_onehot)\n        })\n    \n    metadata_df = pd.DataFrame(metadata)\n    \n    print(f\"\\n✓ Preparados: {len(X_train)} exemplos de treino\")\n    print(f\"✗ Pulados: {skipped} exemplos\")\n    \n    if len(X_train) > 0:\n        lengths = [len(x) for x in X_train]\n        print(f\"\\nEstatísticas de comprimento:\")\n        print(f\"  Min: {min(lengths)}\")\n        print(f\"  Max: {max(lengths)}\")\n        print(f\"  Média: {np.mean(lengths):.1f}\")\n        print(f\"  Mediana: {np.median(lengths):.1f}\")\n    \n    print(\"=\"*70)\n    \n    return X_train, y_train, metadata_df\n\ndef prepare_test_data(test_seq_df):\n    print(\"\\n\" + \"=\"*70)\n    print(\"PREPARANDO DADOS DE TESTE\")\n    print(\"=\"*70)\n    \n    X_test = []\n    test_metadata = []\n    \n    for idx, row in tqdm(test_seq_df.iterrows(), total=len(test_seq_df), desc=\"Processing\"):\n        target_id = row['target_id']\n        sequence = row['sequence']\n        seq_onehot = string_to_onehot(sequence)\n        X_test.append(seq_onehot)\n        test_metadata.append({\n            'target_id': target_id,\n            'sequence': sequence,\n            'length': len(sequence),\n            'gc_content': calculate_gc_content(seq_onehot)\n        })\n    \n    test_metadata_df = pd.DataFrame(test_metadata)\n    print(f\"\\n✓ Preparados: {len(X_test)} exemplos de teste\")\n    \n    if len(X_test) > 0:\n        lengths = [len(x) for x in X_test]\n        print(f\"\\nEstatísticas de comprimento:\")\n        print(f\"  Min: {min(lengths)}\")\n        print(f\"  Max: {max(lengths)}\")\n        print(f\"  Média: {np.mean(lengths):.1f}\")\n    \n    print(\"=\"*70)\n    \n    return X_test, test_metadata_df\n\n# ============================================================================\n# MODELO - PLAN B\n# ============================================================================\n\nclass TurbinedReferenceModelPlanB:\n    \"\"\"\n    Modelo Plan B:\n    Base Conservative + Refinamento Físico + Ensemble por Votação\n    \"\"\"\n    def __init__(self, noise_level=0.21, n_candidates=20, n_final=5,\n                 use_cache=True, generate_synthetic=True, n_synthetic=100,\n                 use_physics_refinement=True, physics_iterations=8):\n        self.noise_level = noise_level\n        self.n_candidates = n_candidates\n        self.n_final = n_final\n        self.use_cache = use_cache\n        self.generate_synthetic = generate_synthetic\n        self.n_synthetic = n_synthetic\n        self.use_physics_refinement = use_physics_refinement\n        self.physics_iterations = physics_iterations\n        if self.use_cache:\n            self.cache = SmartReferenceCache(max_cache_size=2000)\n        self.is_fitted = False\n    \n    def fit(self, X_train, y_train):\n        print(\"\\n\" + \"=\"*70)\n        print(\"TREINANDO MODELO - PLAN B: PHYSICS + VOTING\")\n        print(\"=\"*70)\n        \n        self.reference_pool = []\n        for i in range(len(X_train)):\n            self.reference_pool.append({\n                'sequence': X_train[i],\n                'structure': y_train[i],\n                'is_synthetic': False\n            })\n        \n        print(f\"✓ {len(self.reference_pool)} referências originais\")\n        \n        if self.generate_synthetic and len(self.reference_pool) >= 2:\n            print(f\"\\nGerando {self.n_synthetic} referências sintéticas...\")\n            synthetic_refs = generate_synthetic_references(\n                self.reference_pool, n_synthetic=self.n_synthetic\n            )\n            self.reference_pool.extend(synthetic_refs)\n            print(f\"✓ Total: {len(self.reference_pool)} referências\")\n        \n        print(f\"\\n🔧 Configurações:\")\n        print(f\"  Refinamento físico: {'Ativado' if self.use_physics_refinement else 'Desativado'}\")\n        if self.use_physics_refinement:\n            print(f\"  Iterações de física: {self.physics_iterations}\")\n        print(f\"  Ensemble: Votação ponderada (consenso)\")\n        \n        self.is_fitted = True\n        print(\"\\n✓ MODELO TREINADO!\")\n        print(\"=\"*70)\n        return self\n    \n    def predict(self, X_test, verbose=True):\n        if not self.is_fitted:\n            raise ValueError(\"Modelo não foi treinado!\")\n        \n        if verbose:\n            print(\"\\n\" + \"=\"*70)\n            print(\"GERANDO PREDIÇÕES - PLAN B\")\n            print(\"=\"*70)\n        \n        all_predictions = []\n        \n        iterator = tqdm(enumerate(X_test), total=len(X_test), desc=\"Predicting\") if verbose else enumerate(X_test)\n        \n        for idx, target_seq in iterator:\n            # Cache lookup\n            best_refs = None\n            if self.use_cache:\n                best_refs = self.cache.get_references(target_seq, self.reference_pool)\n            \n            if best_refs is None:\n                best_refs = select_references_adaptive(\n                    target_seq, self.reference_pool, base_top_k=10\n                )\n                if self.use_cache:\n                    self.cache.add_to_cache(target_seq, best_refs)\n            \n            # Interpolação\n            target_props = calculate_sequence_properties(target_seq)\n            base_structure = interpolate_multiple_references(\n                best_refs, target_props, method='weighted'\n            )\n            \n            if base_structure is None:\n                base_structure = best_refs[0]['structure']\n            \n            # Gerar candidatos com diversidade\n            candidates = generate_candidates_diverse(\n                base_structure, target_seq, \n                n_candidates=self.n_candidates,\n                noise_level=self.noise_level\n            )\n            \n            # MELHORIA #1: Refinamento físico\n            if self.use_physics_refinement:\n                refined_candidates = []\n                for candidate in candidates:\n                    refined = iterative_physics_refinement(\n                        candidate, target_seq, n_iterations=self.physics_iterations\n                    )\n                    refined_candidates.append(refined)\n                candidates = refined_candidates\n            \n            # MELHORIA #2: Ensemble por votação\n            final_structures = weighted_voting_ensemble(\n                candidates, target_seq, n_final=self.n_final\n            )\n            \n            while len(final_structures) < self.n_final:\n                final_structures.append(final_structures[0])\n            final_structures = final_structures[:self.n_final]\n            \n            all_predictions.append(np.array(final_structures))\n        \n        if verbose:\n            print(\"\\n✓ PREDIÇÕES CONCLUÍDAS!\")\n            if self.use_cache:\n                stats = self.cache.get_stats()\n                print(f\"Cache: {stats['cache_size']} entries, hit rate: {stats['hit_rate']:.2%}\")\n            print(\"=\"*70)\n        \n        return all_predictions\n\n# ============================================================================\n# SUBMISSÃO\n# ============================================================================\n\ndef fix_nans_in_submission(submission_df):\n    nan_count_before = submission_df.isna().sum().sum()\n    \n    if nan_count_before == 0:\n        return submission_df\n    \n    print(f\"\\n⚠️ Corrigindo {nan_count_before} NaNs...\")\n    \n    coord_cols = []\n    for i in range(1, 6):\n        coord_cols.extend([f'x_{i}', f'y_{i}', f'z_{i}'])\n    \n    for idx, row in submission_df.iterrows():\n        if row[coord_cols].isna().any():\n            for struct_idx in range(1, 6):\n                x_col = f'x_{struct_idx}'\n                y_col = f'y_{struct_idx}'\n                z_col = f'z_{struct_idx}'\n                \n                if pd.isna(row[x_col]) or pd.isna(row[y_col]) or pd.isna(row[z_col]):\n                    valid_coords = None\n                    \n                    for other_idx in range(1, 6):\n                        if other_idx == struct_idx:\n                            continue\n                        \n                        other_x = row[f'x_{other_idx}']\n                        other_y = row[f'y_{other_idx}']\n                        other_z = row[f'z_{other_idx}']\n                        \n                        if not pd.isna(other_x) and not pd.isna(other_y) and not pd.isna(other_z):\n                            valid_coords = (other_x, other_y, other_z)\n                            break\n                    \n                    if valid_coords:\n                        submission_df.at[idx, x_col] = valid_coords[0]\n                        submission_df.at[idx, y_col] = valid_coords[1]\n                        submission_df.at[idx, z_col] = valid_coords[2]\n                    else:\n                        resid = row['resid']\n                        submission_df.at[idx, x_col] = resid * 5.9\n                        submission_df.at[idx, y_col] = 0.0\n                        submission_df.at[idx, z_col] = 0.0\n    \n    nan_count_after = submission_df.isna().sum().sum()\n    \n    if nan_count_after > 0:\n        submission_df[coord_cols] = submission_df[coord_cols].fillna(0.0)\n    \n    print(f\"✓ Corrigidos: {nan_count_before} → {nan_count_after}\")\n    \n    return submission_df\n\ndef create_submission_correct(predictions, test_metadata, test_sequences_df, output_file=\"submission.csv\"):\n    print(\"\\n\" + \"=\"*70)\n    print(\"CRIANDO SUBMISSÃO\")\n    print(\"=\"*70)\n    \n    submission_rows = []\n    \n    for idx, target_id in enumerate(tqdm(test_metadata['target_id'], desc=\"Creating submission\")):\n        sequence = test_metadata.iloc[idx]['sequence']\n        structures = predictions[idx]\n        n_residues = len(structures[0])\n        \n        for res_idx in range(n_residues):\n            row_id = f\"{target_id}_{res_idx + 1}\"\n            base = sequence[res_idx] if res_idx < len(sequence) else 'N'\n            resid = res_idx + 1\n            \n            coords_all_structures = []\n            for struct_idx in range(5):\n                coords = structures[struct_idx][res_idx]\n                coords_all_structures.extend([coords[0], coords[1], coords[2]])\n            \n            row = {\n                'ID': row_id,\n                'resname': base,\n                'resid': resid\n            }\n            \n            for struct_idx in range(5):\n                base_col_idx = struct_idx * 3\n                row[f'x_{struct_idx + 1}'] = coords_all_structures[base_col_idx]\n                row[f'y_{struct_idx + 1}'] = coords_all_structures[base_col_idx + 1]\n                row[f'z_{struct_idx + 1}'] = coords_all_structures[base_col_idx + 2]\n            \n            submission_rows.append(row)\n    \n    submission_df = pd.DataFrame(submission_rows)\n    \n    coord_cols = []\n    for i in range(1, 6):\n        coord_cols.extend([f'x_{i}', f'y_{i}', f'z_{i}'])\n    \n    columns_order = ['ID', 'resname', 'resid'] + coord_cols\n    submission_df = submission_df[columns_order]\n    \n    submission_df = fix_nans_in_submission(submission_df)\n    \n    output_path = os.path.join(OUTPUT_DIR, output_file)\n    submission_df.to_csv(output_path, index=False)\n    \n    print(f\"\\n✓ Submissão salva: {output_path}\")\n    print(f\"  Linhas: {len(submission_df)}\")\n    print(f\"  Colunas: {len(submission_df.columns)}\")\n    print(f\"  NaNs: {submission_df.isna().sum().sum()}\")\n    \n    print(\"=\"*70)\n    \n    return submission_df\n\n# ============================================================================\n# PIPELINE PRINCIPAL\n# ============================================================================\n\ndef main_pipeline_plan_b(use_v2=True, max_train_length=None, \n                         use_physics=True, physics_iterations=8):\n    \"\"\"\n    Pipeline Plan B: Conservative + Física + Votação\n    \n    Args:\n        use_v2: Usar dados v2 (default: True)\n        max_train_length: Limite de comprimento (None = sem limite)\n        use_physics: Ativar refinamento físico (default: True)\n        physics_iterations: Número de iterações de física (default: 8)\n    \"\"\"\n    print(\"\"\"\n    ╔══════════════════════════════════════════════════════════════╗\n    ║     TURBINED RNA 3D - PLAN B                                ║\n    ║                                                              ║\n    ║     Base: Conservative 0.38346 private                      ║\n    ║     + Refinamento iterativo com física (8 iterações)        ║\n    ║     + Ensemble por votação ponderada (consenso)             ║\n    ║                                                              ║\n    ║     Score esperado: 0.39-0.395 private                      ║\n    ║     Tempo estimado: 4-5 minutos                             ║\n    ║                                                              ║\n    ║     ⚠️  SE < 0.383: REVERTER para Conservative             ║\n    ║     🎉 SE > 0.390: SUCESSO!                                 ║\n    ╚══════════════════════════════════════════════════════════════╝\n    \"\"\")\n    \n    # Carregar dados\n    print(\"\\n📂 Carregando dados...\")\n    if use_v2:\n        train_seq = pd.read_csv(os.path.join(DATA_DIR, \"train_sequences.v2.csv\"))\n        train_labels = pd.read_csv(os.path.join(DATA_DIR, \"train_labels.v2.csv\"))\n    else:\n        train_seq = pd.read_csv(os.path.join(DATA_DIR, \"train_sequences.csv\"))\n        train_labels = pd.read_csv(os.path.join(DATA_DIR, \"train_labels.csv\"))\n    \n    test_seq = pd.read_csv(os.path.join(DATA_DIR, \"test_sequences.csv\"))\n    \n    print(f\"✓ Dados carregados\")\n    print(f\"  Train sequences: {len(train_seq)}\")\n    print(f\"  Test sequences: {len(test_seq)}\")\n    \n    # Preparar\n    X_train, y_train, train_metadata = prepare_training_data_final(\n        train_seq, train_labels, max_length=max_train_length\n    )\n    \n    if len(X_train) == 0:\n        print(\"❌ ERRO: Nenhum dado de treino preparado!\")\n        return None, None, None\n    \n    X_test, test_metadata = prepare_test_data(test_seq)\n    \n    # Treinar modelo Plan B\n    model = TurbinedReferenceModelPlanB(\n        noise_level=0.21,\n        n_candidates=20,\n        n_final=5,\n        use_cache=True,\n        generate_synthetic=True,\n        n_synthetic=100,\n        use_physics_refinement=use_physics,\n        physics_iterations=physics_iterations\n    )\n    \n    model.fit(X_train, y_train)\n    \n    # Predizer\n    predictions = model.predict(X_test, verbose=True)\n    \n    # Criar submissão\n    submission = create_submission_correct(\n        predictions, \n        test_metadata, \n        test_seq,\n        \"submission.csv\"\n    )\n    \n    print(\"\\n\" + \"=\"*70)\n    \n    return model, predictions, submission\n\n# ============================================================================\n# EXECUÇÃO\n# ============================================================================\n\nif __name__ == \"__main__\":\n    \n    model, predictions, submission = main_pipeline_plan_b(\n        use_v2=True,\n        max_train_length=None,      # Processar TODAS as sequências\n        use_physics=True,            # Ativar refinamento físico\n        physics_iterations=8         # 8 iterações de física\n    )\n    \n    if model is not None:\n        print(\"\\n✅ SUCESSO! Arquivo pronto para submissão.\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-13T02:19:43.260613Z","iopub.execute_input":"2025-10-13T02:19:43.260970Z","iopub.status.idle":"2025-10-13T02:27:42.183144Z","shell.execute_reply.started":"2025-10-13T02:19:43.260942Z","shell.execute_reply":"2025-10-13T02:27:42.182096Z"}},"outputs":[],"execution_count":null},{"id":"2d2a0fd7-ae76-464f-b1f0-15737c4519d4","cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}