{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.11","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":87793,"databundleVersionId":12024591,"sourceType":"competition"}],"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"This notebook based on trRosettaRNA : https://yanglab.qd.sdu.edu.cn/trRosettaRNA/\n","metadata":{}},{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\n# for dirname, _, filenames in os.walk('/kaggle/input'):\n#     for filename in filenames:\n#         print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-05-02T08:14:57.001776Z","iopub.execute_input":"2025-05-02T08:14:57.002027Z","iopub.status.idle":"2025-05-02T08:14:57.350778Z","shell.execute_reply.started":"2025-05-02T08:14:57.001970Z","shell.execute_reply":"2025-05-02T08:14:57.349711Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!wget -c https://yanglab.qd.sdu.edu.cn/trRosettaRNA/download/trRosettaRNA_v1.1.zip","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:39:11.647718Z","iopub.execute_input":"2025-05-02T06:39:11.648077Z","iopub.status.idle":"2025-05-02T06:42:02.349165Z","shell.execute_reply.started":"2025-05-02T06:39:11.648055Z","shell.execute_reply":"2025-05-02T06:42:02.347950Z"},"collapsed":true,"jupyter":{"outputs_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!unzip ./trRosettaRNA_v1.1.zip\nos.chdir('/kaggle/working/trRosettaRNA_v1.1/')\n!git clone https://github.com/jaswindersingh2/SPOT-RNA.git\n!wget 'https://www.dropbox.com/s/dsrcf460nbjqpxa/SPOT-RNA-models.tar.gz' \n# !wget -O SPOT-RNA-models.tar.gz 'https://app.nihaocloud.com/f/fbf3315a91d542c0bdc2/?dl=1'\n!tar -xvzf SPOT-RNA-models.tar.gz && rm SPOT-RNA-models.tar.gz\n!pip install pyrosetta-installer \n!mv '/kaggle/working/trRosettaRNA_v1.1/SPOT-RNA-models/' '/kaggle/working/trRosettaRNA_v1.1/SPOT-RNA'\n!pip install biopython\n!python -c 'import pyrosetta_installer; pyrosetta_installer.install_pyrosetta(serialization=True, silent=False)'\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:42:02.608843Z","iopub.execute_input":"2025-05-02T06:42:02.609148Z","iopub.status.idle":"2025-05-02T06:43:58.168080Z","shell.execute_reply.started":"2025-05-02T06:42:02.609120Z","shell.execute_reply":"2025-05-02T06:43:58.166523Z"},"collapsed":true,"jupyter":{"outputs_hidden":true}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:43:58.170539Z","iopub.execute_input":"2025-05-02T06:43:58.171006Z","iopub.status.idle":"2025-05-02T06:43:58.308975Z","shell.execute_reply.started":"2025-05-02T06:43:58.170964Z","shell.execute_reply":"2025-05-02T06:43:58.307707Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import os\nos.getcwd()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:43:58.310893Z","iopub.execute_input":"2025-05-02T06:43:58.311185Z","iopub.status.idle":"2025-05-02T06:44:03.573265Z","shell.execute_reply.started":"2025-05-02T06:43:58.311157Z","shell.execute_reply":"2025-05-02T06:44:03.572214Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%writefile fold.py  \nimport tempfile\nimport glob\nimport sys\nimport os, time\nfrom pathlib import Path\nfrom folding.arguments import get_args\nfrom folding.utils_cst import npz2cst\nfrom folding.utils_ros import fold_from_cst\nimport os\nimport argparse\nimport tempfile\ndef main():\n    args = get_args()\n    \n    # Ensure output directory exists\n    os.makedirs(os.path.dirname(os.path.abspath(args.OUT)), exist_ok=True)\n    \n    # Use -tmp directly as args.tmpdir\n    args.tmpdir = args.TMPDIR if hasattr(args, 'TMPDIR') else args.tmp\n    os.makedirs(args.tmpdir, exist_ok=True)\n    print('temp folder:     ', args.tmpdir)\n    \n    # Parse NPZ into Rosetta-format restraint files\n    npz2cst(args)\n    \n    # Perform energy minimization\n    fold_from_cst(args)\n\nif __name__ == '__main__':\n    main()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:44:19.440726Z","iopub.execute_input":"2025-05-02T06:44:19.441054Z","iopub.status.idle":"2025-05-02T06:44:19.450114Z","shell.execute_reply.started":"2025-05-02T06:44:19.441022Z","shell.execute_reply":"2025-05-02T06:44:19.448779Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import os\n\n# Physical cores (excludes hyperthreading)\nphysical_cores = os.cpu_count() // 2  # For hyperthreaded CPUs\nprint(f\"Physical cores: {physical_cores}\")\n\n# Total logical processors\nprint(f\"Logical cores: {os.cpu_count()}\")\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:44:19.824680Z","iopub.execute_input":"2025-05-02T06:44:19.825000Z","iopub.status.idle":"2025-05-02T06:44:19.831096Z","shell.execute_reply.started":"2025-05-02T06:44:19.824978Z","shell.execute_reply":"2025-05-02T06:44:19.830065Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# %% [markdown]\n# # trRosettaRNA for RNA 3D Structure Prediction Challenge\n# This notebook runs the trRosettaRNA pipeline on the competition's test dataset,\n# converting FASTA MSAs to A3M, predicting 3D structures, and generating a submission CSV.\nfrom tqdm import tqdm\n# %% [code]\nimport os\nimport pandas as pd\nimport numpy as np\nimport subprocess\nimport tempfile\nimport torch\nfrom Bio import SeqIO\nfrom Bio.Seq import Seq\nfrom Bio.SeqRecord import SeqRecord\nimport re\nfrom pathlib import Path\n\n# Paths to trRosettaRNA scripts and data\nTRROSETTA_DIR = \"/kaggle/working/trRosettaRNA_v1.1/\"\nPREDICT_PY = os.path.join(TRROSETTA_DIR, \"predict.py\")\nFOLD_PY = os.path.join(TRROSETTA_DIR, \"fold.py\")\nMODEL_DIR = os.path.join(TRROSETTA_DIR, \"params/model_1\")\nDATA_DIR = \"/kaggle/input/stanford-rna-3d-folding/\"#/app/Kaggle/StanfordRNA/\"\nMSA_DIR = os.path.join(DATA_DIR, \"MSA\")\nOUTPUT_DIR = \"/kaggle/working/output/\"\nPath(OUTPUT_DIR).mkdir(exist_ok=True)\n\n# Function to convert FASTA MSA to A3M\ndef fasta_to_a3m(fasta_file, a3m_file):\n    \"\"\"\n    Convert FASTA MSA to A3M format using HH-suite's reformat.pl or Python parsing.\n    \"\"\"\n    try:\n        # Option 1: Use HH-suite's reformat.pl (recommended)\n        subprocess.run(\n            [\"reformat.pl\", \"fas\", \"a3m\", fasta_file, a3m_file],\n            check=True, capture_output=True\n        )\n    except FileNotFoundError:\n        # Option 2: Python-based conversion (simplified)\n        records = list(SeqIO.parse(fasta_file, \"fasta\"))\n        with open(a3m_file, \"w\") as f:\n            for rec in records:\n                # Remove gaps (except '-') and convert to uppercase\n                seq = re.sub(r\"[^A-Z-]\", \"-\", str(rec.seq).upper())\n                f.write(f\">{rec.id}\\n{seq}\\n\")\n        print(f\"Converted {fasta_file} to {a3m_file} using Python\")\n\n# Function to create FASTA file from sequence\ndef create_fasta(target_id, sequence, fasta_file):\n    \"\"\"\n    Create a FASTA file for a single RNA sequence.\n    \"\"\"\n    record = SeqRecord(Seq(sequence), id=target_id, description=\"\")\n    with open(fasta_file, \"w\") as f:\n        SeqIO.write(record, f, \"fasta\")\n\n# Function to extract C1' coordinates from PDB\ndef extract_c1p_coordinates(pdb_file, sequence):\n    \"\"\"\n    Extract C1' atom coordinates from a PDB file for each residue.\n    Returns a list of [x, y, z] for C1' atoms, ordered by resid.\n    \"\"\"\n    coords = []\n    L = len(sequence)\n    with open(pdb_file, \"r\") as f:\n        for line in f:\n            if line.startswith(\"ATOM\"):\n                atom_name = line[12:16].strip()\n                resname = line[17:20].strip()\n                resid = int(line[22:26].strip())\n                if atom_name == \"C1'\" and resname in [\"A\", \"C\", \"G\", \"U\"] and resid <= L:\n                    x = float(line[30:38].strip())\n                    y = float(line[38:46].strip())\n                    z = float(line[46:54].strip())\n                    coords.append([resid, resname, x, y, z])\n    # Sort by resid to ensure correct order\n    coords.sort(key=lambda x: x[0])\n    return [[x[2], x[3], x[4]] for x in coords]  # Return [x, y, z] list\n\n# Load test sequences\ntest_df = pd.read_csv(os.path.join(DATA_DIR, \"test_sequences.csv\"))\n\n# Submission dataframe\nsubmission_rows = []\nimport glob\n\n# Process each test sequence\nfor idx, row in tqdm(test_df[:1].iterrows()):\n    target_id = row[\"target_id\"]\n    sequence = row[\"sequence\"]\n    L = len(sequence)\n    \n    print(f\"Processing {target_id} (length {L})...\")\n    \n    # Step 1: Convert MSA\n    fasta_msa = os.path.join(MSA_DIR, f\"{target_id}.MSA.fasta\")\n    a3m_msa = os.path.join(OUTPUT_DIR, f\"{target_id}.a3m\")\n    if not os.path.exists(fasta_msa):\n        print(f\"MSA file {fasta_msa} not found, skipping {target_id}\")\n        continue\n    try:\n        fasta_to_a3m(fasta_msa, a3m_msa)\n    except Exception as e:\n        print(f\"MSA conversion failed for {target_id}: {e}\")\n        continue\n    \n    # Step 2: Create FASTA\n    fasta_file = os.path.join(OUTPUT_DIR, f\"{target_id}.fasta\")\n    create_fasta(target_id, sequence, fasta_file)\n    \n    # Step 3: Run predict.py\n    npz_file = os.path.join(OUTPUT_DIR, f\"{target_id}.npz\")\n    print('predict.py')\n    predict_cmd = [\"python\", PREDICT_PY, \"-i\", a3m_msa, \"-o\", npz_file, \"-mdir\", MODEL_DIR, \"-cpu\", \"6\"]\n    try:\n        result = subprocess.run(predict_cmd, check=True, capture_output=True, text=True)\n        print(f\"predict.py output: {result.stdout}\")\n    except subprocess.CalledProcessError as e:\n        print(f\"predict.py failed for {target_id}: {e.stderr}\")\n        continue\n    \n    # Step 4: Run fold.py\n    pdb_out = os.path.join(OUTPUT_DIR, f\"{target_id}.pdb\")\n    temp_dir = os.path.join(OUTPUT_DIR, f\"temp_{target_id}\")\n    os.makedirs(temp_dir, exist_ok=True)\n    print('temp_dir',temp_dir)\n    fold_cmd = [\"python\", FOLD_PY, \"-npz\", npz_file, \"-fa\", fasta_file, \"-out\", pdb_out, \"-nm\", \"5\", \"-tmp\", temp_dir]\n    try:\n        result = subprocess.run(fold_cmd, check=True, capture_output=True, text=True)\n        print(f\"fold.py output: {result.stdout}\")\n    except subprocess.CalledProcessError as e:\n        print(f\"fold.py failed for {target_id}: {e.stderr}\")\n        continue\n    \n    # Step 5: Check for PDB files\n    pdb_files = glob.glob(f\"{temp_dir}/model_*.pdb\")\n    if len(pdb_files) < 5:\n        print(f\"Found only {len(pdb_files)} PDB files in {temp_dir}, expected 5\")\n        # Fallback: Search subdirectories\n        pdb_files = glob.glob(f\"{temp_dir}/**/model_*.pdb\", recursive=True)\n        if len(pdb_files) < 5 and os.path.exists(pdb_out):\n            print(f\"Falling back to primary PDB {pdb_out} for {target_id}\")\n            coords = extract_c1p_coordinates(pdb_out, sequence)\n            if len(coords) == L:\n                all_coords = [coords] * 5\n            else:\n                print(f\"Warning: Incomplete coordinates in {pdb_out} for {target_id} ({len(coords)}/{L})\")\n                continue\n        elif len(pdb_files) < 5:\n            print(f\"Still insufficient PDB files ({len(pdb_files)}) for {target_id}\")\n            continue\n            \n\n    # Step 6: Extract C1' coordinates\n    all_coords = []\n    for i in range(1, 6):\n        pdb_file = next((f for f in pdb_files if f.endswith(f\"model_{i}.pdb\")), None)\n        if pdb_file and os.path.exists(pdb_file):\n            coords = extract_c1p_coordinates(pdb_file, sequence)\n            if len(coords) == L:\n                all_coords.append(coords)\n            else:\n                print(f\"Warning: Incomplete coordinates in {pdb_file} for {target_id} ({len(coords)}/{L})\")\n        else:\n            print(f\"Warning: {pdb_file} not found for {target_id}\")\n    \n    # Step 7: Format submission rows\n    if len(all_coords) == 5:\n        for resid in range(1, L + 1):\n            ID = f\"{target_id}_{resid}\"\n            resname = sequence[resid - 1]\n            row_dict = {\"ID\": ID, \"resname\": resname, \"resid\": resid}\n            for i, coords in enumerate(all_coords, 1):\n                x, y, z = coords[resid - 1]\n                row_dict[f\"x_{i}\"] = x\n                row_dict[f\"y_{i}\"] = y\n                row_dict[f\"z_{i}\"] = z\n            submission_rows.append(row_dict)\n    else:\n        print(f\"Warning: Insufficient structures for {target_id} ({len(all_coords)}/5), skipping submission\")\n\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-02T06:44:23.227670Z","iopub.execute_input":"2025-05-02T06:44:23.228770Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"submission_df = pd.DataFrame(submission_rows)\nsubmission_cols = [\"ID\", \"resname\", \"resid\"] + \\\n                 [f\"{c}_{i}\" for i in range(1, 6) for c in [\"x\", \"y\", \"z\"]]\nsubmission_df = submission_df[submission_cols]  # Reorder columns\nsubmission_df.to_csv(\"submission.csv\", index=False)\nprint(\"Submission saved to submission.csv\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# !pip install jupytext  \n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-01T17:30:09.569954Z","iopub.execute_input":"2025-05-01T17:30:09.572283Z","iopub.status.idle":"2025-05-01T17:30:09.577494Z","shell.execute_reply.started":"2025-05-01T17:30:09.572237Z","shell.execute_reply":"2025-05-01T17:30:09.575834Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# jupytext --update-metadata '{\"jupytext\": {\"cell_markers\": \"\\\"\\\"\\\"\"}}' --to py notebook.ipynb  \n# %run fold.py --OUT '/kaggle/working/trRosettaRNA_v1.1' --TMPDIR /kaggle/working/tmp\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-01T17:36:40.403423Z","iopub.execute_input":"2025-05-01T17:36:40.403809Z","iopub.status.idle":"2025-05-01T17:36:40.409300Z","shell.execute_reply.started":"2025-05-01T17:36:40.403780Z","shell.execute_reply":"2025-05-01T17:36:40.408216Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# ls# %tb","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-01T19:10:48.254639Z","iopub.execute_input":"2025-05-01T19:10:48.255018Z","iopub.status.idle":"2025-05-01T19:10:48.260670Z","shell.execute_reply.started":"2025-05-01T19:10:48.254991Z","shell.execute_reply":"2025-05-01T19:10:48.259362Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}