diff --git a/README.md b/README.md index 2e1a1cf..2dba14e 100755 --- a/README.md +++ b/README.md @@ -2,7 +2,7 @@ This package provides a custom docker-free adaptation of AlphaFold 2.0 (Jumper et al, 2021) and ColabFold (Mirdita et al, 2021). -## Installation +## Installation of ComplexFold ### Get ComplexFold ```bash git clone https://github.com/muthoff/ComplexFold @@ -14,11 +14,11 @@ CFDIR=$PWD ```bash wget https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh bash Miniconda3-latest-Linux-x86_64.sh -``` -Let it initialize -```bash +# Let it initialize + echo "conda deactivate\n" >> ~/.bashrc source ~/.bashrc + conda update conda conda install pip rm -r Miniconda3-latest-Linux-x86_64.sh @@ -45,7 +45,7 @@ patch -p0 < $CFDIR/docker/openmm.patch conda deactivate ``` -## Genetic databases +## Download of genetic databases This step requires `aria2c` to be installed on your machine. @@ -59,7 +59,8 @@ AlphaFold needs multiple genetic (sequence) databases to run: * [PDB](https://www.rcsb.org/) (structures in the mmCIF format). There is a script `ComplexFold/scripts/download_all_data.sh` that can be used to download -and set up all of these databases: +and set up all of these databases. The scripts provided download the most recent database +vesions: * Default: @@ -118,7 +119,7 @@ $DOWNLOAD_DIR/ # Total: ~ 2.2 TB (download: 438 GB) `bfd/` is only downloaded if you download the full databasees, and `small_bfd/` is only downloaded if you download the reduced databases. -## Model parameters +## Downloaf of model parameters While the AlphaFold code is licensed under the Apache 2.0 License, the AlphaFold parameters are made available for non-commercial use only under the terms of the @@ -138,22 +139,24 @@ will download parameters for: Jumper et al. 2021, Suppl. Methods 1.9.7 for details). ## Final set up and simple run of ComplexFold -1. Go through the section "Defaults" `ComplexFold/run_complexfold.sh` and change appropriatly. +1. Go through the section "Defaults" in `ComplexFold/run_complexfold.sh` and change appropriatly 2. Run `run_complexfold.sh` pointing to a FASTA file containing the protein sequence for which you wish to predict the structure and a description line stating with '>'. You can provide multiple sequences per FASTA file in order to predict complexes. ```bash - python3 run_complexfold.sh -f path/to/fasta + bash run_complexfold.sh -f path/to/fasta ``` ## Features ### FASTA file Each FASTA file must contain a description line stating with '>' followed by the name of the protein. ComplexFold automatically cuts the name at the first space. The sequence -can encompass multiple line. +can encompass multiple lines. + For the predictions of complexes you have to provide the description and the sequence for each component, including homomers. The homomer sequence and description must be identical! + A trimer with two unique proteins would have a FASTA file like this for example: ```fasta @@ -172,41 +175,45 @@ PTWS ### Complex Mode ComplexFold auto detects heteromers and changes a few things for better prediction. This includes: -1. Usage of 'ptm' models (Equivalent to flag `-c`.) +1. Usage of 'ptm' models (Equivalent to flag `-c`) 2. Reducing the number of potential templates (10 per heteromer) and increaseing the number of actual templates used (4 per heteromer). At the moment, I cannot control - which templates are finally used. These changes increase the likelyhood of using a + which templates are finally used. These changes increase the likelihood of using a template specific for each heteromer though. -### MSA library +### MSA library and custom MSAs By default ComplexFold uses the same databases and tools like AlphaFold2.0. Three MSAs are computed for each heteromer individually and saved under the name of the tool and -heteromer in the 'msa' subdir of the output. This includes also the search for templates. -You can copy these files into the 'msa_library'. ComplexFold search this path for appropriate -files before computing new MSAs. This way you do not have to re-compute MSAs for new -predictions, like in complex ther one component in known. +the complex component description in the `msa/` subdir of the output. This includes also +the search for templates. In addition, custom MSAs can be provided in either `sto`or `a3m` +format. The file names must be `custom_{component description}.{sto/a3m}`. + +You can copy all these files into the `msa_library`. ComplexFold searches this path for +appropriate files before computing new MSAs. This way you do not have to re-compute MSAs +for new predictions, like in varying complexes where one component is always the same. Please delete any files in the msa_library after updating the databases! ### Re-run-in-place -You can re-run complexfold in the same ouput dir `-o` with differing parameters. All files -there will be moved into a new subfolder called `result_n`, where n is the n-th run. -Moreover, the `msa/` subfolder is searched for any matching MSA/Template file. Appropriate +You can re-run ComplexFold in the same ouput dir (`-o`) with differing parameters. All files +there will be moved into a new subfolder called `result_n`, where `n` is the `n`-th run. +Moreover, the `msa/` subfolder is searched for any matching MSA/template files. Appropriate files are read and missing ones are computed. ### Recycling AlphaFold uses its output and recycles the information into a new iteration of folding. By default 3 recycles are performed, however, this can be increased by the user by providing -the flag `-r `, where is an integer. Usually smaller -than 30 is sufficient. This can be combined with the flag `-l `, where - is a decimal number and declares when to stop recycling. This number -is based on a predicted Ca-RMS and thus and the smaller, the better the prediction. However, -this is not strictly true and pLDDT as well as pTM usually remain high in very difficult or -inheretenlty disorder proteins nonetheless. This happen even if the observed value drops below -the set threshold. +the flag `-r `, where `` is an integer. Usually `` smaller +than 30 is sufficient. + +This can be combined with the flag `-l `, where `` +is a decimal number and declares when to stop recycling. This number is based on a predicted +Ca-RMS and thus the smaller, the better the prediction. However, this is not strictly true +and pLDDT as well as pTM usually remain high in very difficult or inheritently disordered +proteins. This happen even if the observed value drops below the set threshold. -The smaller the more likely recycling happens exactly -times. is by default 0, usefull values are usually 0.25-0.33. +The smaller `` the more likely recycling happens exactly `` +times. `` is by default 0, usefull values are usually 0.25-0.33. ### Random seeds and sampling prediction multiple times AlphaFold utilises a random seed to initialize features for predictions and the predictins @@ -216,16 +223,18 @@ only a very shallow MSA (few sequences). However, GPU processing is sometimes (i case) not deterministic. Thus the entire AlphaFold pipeline is non-deterministic either and the same seed can give different results. -You can provide multiple seeds in a comma separated list `-s 123,567,...` or set `-x `. -The latter genrates as much random seed as requested. This lets ComplexFold run the prediction -multiple times with all given models. Only the best 5 models, based on plDDT (normal models) or -pTM (ptm models), are relaxed and given as output. Keep in mind that the scores are themselves -only prediction and a this way predicted model may just have an erroneously high score. +With ComplexFold you can provide multiple seeds in a comma separated list `-s 123,567,...` or +set `-x `. The latter genrates as much random seed as requested. This lets ComplexFold +run the prediction multiple times with all given models. Only the best 5 models, based on plDDT +(normal models) or pTM (ptm models), are relaxed and given as output. + +Keep in mind that the scores are themselves only prediction and this way predicted model may +just have an erroneously high score. ### Small bfd -You can either use the full bfd or the faster smaller than. `-p ` where the argument is -either full_dbs or reduced_dbs. This can especially help if MSA get too big (several GB) and -HHBlits crashes. +You can either use the full bfd or the faster smaller one. `-p ` where the argument is +either `full_dbs` or `reduced_dbs`. This can especially help if the MSA gets too big (several GB) and +HHBlits crashes. The MSAs often do not differ. ### Difficulty presets - thoroughness ComplexFold combines the listed features in difficulty presets: `-y `. By default @@ -237,35 +246,36 @@ ComplexFold uses "alphafold", which is equivalent to the original AlphaFold sett 5. `extreme`: `-r 20 -l 0.33 -x 30 -p full_dbs` ### Keep results with a good score in a certain region -I have got an a/b globular prediction which had always 10-15 amino acids unstructured while the -entire protein was around 70 amino acids long. I was afraid this bad region allows the other +Once, I have got an a/b globular prediction which had always 10-15 amino acids unstructured. The +entire protein was only 70 amino acids short. I was afraid this bad region allowed the other parts to be very good, while itself might even be physically impossible. Hence such a predictions might only be a local minimum and may prevent proper folding/prediction. `-i ` where the start and end position (x-y) is given, lets CompexFold only keep -the models with the best pLDDT in that region. This has to be combined with `-x `. +the models with the best pLDDT in that region. This has to be combined with larger `-x `. + +I was abel to get very different predictions while sampling 50+ seeds. This included a b-barrel (good +average pLDTT and good score at any position), further b-sheet-only-fold (bad average pLDDT) +and more a/b globular-folds which looked similar as the initial one but were differently arranged. -I was abel to get very different prediction sampling 50+ seed, including a b-barrel (good -average pLDTT and good score at any position), further b-sheet only fold (bad average pLDDT) -and further a/b globular fold which look similat as the initial one but were differentlya arranged. Eventually, the initial predction was confirmed to be true by crystallisation, though. ### AlphaFold output -The outputs will be in a subfolder of `output_dir`. They include the computed MSAs, -unrelaxed structures, relaxed structures, raw model outputs, some prediction metadata -and a comprehensive report of given argument, scores and timings. The `output_dir` -directory will have the following structure: +The outputs will be in a subfolder of `output_dir` named like the input FASTA file. It +includes the computed MSAs, unrelaxed structures, relaxed structures, raw model outputs, +some prediction metadata and a comprehensive report of given arguments, scores and timings. +The `output_dir` directory will have the following structure: ``` / msas/ - {protein description}_bfd_uniclust_hits.a3m - {protein description}_mgnify_hits.sto - {protein description}_uniref90_hits.sto - {protein description}_pdb70_hits.hhr + {heteromer description}_bfd_uniclust_hits.a3m + {heteromer description}_mgnify_hits.sto + {heteromer description}_uniref90_hits.sto + {heteromer description}_pdb70_hits.hhr heteromers/ - {protein description}.fa + {heteromer description}.fa parsed_results.pkl unrelaxed_model_{1,2,3,4,5}{,_ptm}-sx.pdb relaxed_model_{1,2,3,4,5}{,_ptm}-sx.pdb @@ -280,10 +290,10 @@ directory will have the following structure: The contents of each output file are as follows: -* `msas/` - A directory containing the files describing the various genetic +* `msas/` – A directory containing the files describing the various genetic tool hits that were used to construct the input MSA. Files are given for - each unique protein identified by it >protein description. -* `heteromers/` - A directory containing FASTA files for each heteromers. + each unique protein identified by the `description` given in the FAST file. +* `heteromers/` – A directory containing FASTA files for each heteromer. * `parsed_results.pkl` – A `pickle` file containing a class which collects raw_model outputs and auxiliary outputs. Please look up exact structure in the alphafold/complexfold.py `class Result` and `class Result_Handler`. @@ -295,17 +305,17 @@ The contents of each output file are as follows: structure prediction (see Jumper et al. 2021, Suppl. Methods 1.8.6 for details). `x` indicates the seed used, look it up in `report.json`. * `model_*-sx.png` – Graphical representation of the model by ColabFold. -* msa_coverage.png – Coverage of the used MSAs by ColabFold. -* PAE.png – Predicted aligned error by ColabFold. -* predicted_contacts.png – Predicted contacts by ColabFold. -* predicted_distogram.png – Predicted distogram by ColabFold. -* pLDDT.png – pLDDT of each residue by ColabFold. +* `msa_coverage.png` – Coverage of the used MSAs by ColabFold. +* `PAE.png` – Predicted aligned error by ColabFold. +* `predicted_contacts.png` – Predicted contacts by ColabFold. +* `predicted_distogram.png` – Predicted distogram by ColabFold. +* `pLDDT.png` – pLDDT of each residue by ColabFold. * `report.json` – A JSON format text file containing the scores, auxillary data, command arguments values. ## Citing this work -Please refer to the github but also acknowledge by citing: +Please refer to this github by Matthias Uthoff but also acknowledge by citing: ```bibtex @Article{AlphaFold2021, @@ -330,12 +340,11 @@ Please refer to the github but also acknowledge by citing: } ``` -5 - ## Acknowledgements -Great thanks goes to Deepmind and the ColabFold contributors. +Great thanks goes to [Deepmind](https://github.com/deepmind/alphafold) +and the [ColabFold](https://github.com/sokrypton/ColabFold) contributors. AlphaFold communicates with and/or references the following separate libraries and packages: @@ -369,6 +378,9 @@ We thank all their contributors and maintainers! This is not an officially supported Google product. +This is derivative of [AlphaFold](https://github.com/deepmind/alphafold) +and [ColabFold](https://github.com/sokrypton/ColabFold) + ### ComplexFold Code License Licensed under the Apache License, Version 2.0 (the "License"); you may not use diff --git a/alphafold/complexfold.py b/alphafold/complexfold.py index 0ff8728..7f2d5f7 100755 --- a/alphafold/complexfold.py +++ b/alphafold/complexfold.py @@ -17,6 +17,8 @@ from typing import Dict from absl import flags +from absl import logging +from alphafold.common import protein from alphafold.model import features from alphafold.data import colabfold diff --git a/imgs/casp14_predictions.gif b/imgs/casp14_predictions.gif deleted file mode 100755 index 0d9ce2c..0000000 Binary files a/imgs/casp14_predictions.gif and /dev/null differ diff --git a/imgs/header.jpg b/imgs/header.jpg deleted file mode 100755 index d4000a1..0000000 Binary files a/imgs/header.jpg and /dev/null differ diff --git a/notebooks/AlphaFold.ipynb b/notebooks/AlphaFold.ipynb deleted file mode 100755 index 6736838..0000000 --- a/notebooks/AlphaFold.ipynb +++ /dev/null @@ -1,694 +0,0 @@ -{ - "cells": [ - { - "cell_type": "markdown", - "metadata": { - "id": "pc5-mbsX9PZC" - }, - "source": [ - "# AlphaFold Colab\n", - "\n", - "This Colab notebook allows you to easily predict the structure of a protein using a slightly simplified version of [AlphaFold v2.0](https://doi.org/10.1038/s41586-021-03819-2). \n", - "\n", - "**Differences to AlphaFold v2.0**\n", - "\n", - "In comparison to AlphaFold v2.0, this Colab notebook uses **no templates (homologous structures)** and a selected portion of the [BFD database](https://bfd.mmseqs.com/). We have validated these changes on several thousand recent PDB structures. While accuracy will be near-identical to the full AlphaFold system on many targets, a small fraction have a large drop in accuracy due to the smaller MSA and lack of templates. For best reliability, we recommend instead using the [full open source AlphaFold](https://github.com/deepmind/alphafold/), or the [AlphaFold Protein Structure Database](https://alphafold.ebi.ac.uk/).\n", - "\n", - "Please note that this Colab notebook is provided as an early-access prototype and is not a finished product. It is provided for theoretical modelling only and caution should be exercised in its use. \n", - "\n", - "**Citing this work**\n", - "\n", - "Any publication that discloses findings arising from using this notebook should [cite](https://github.com/deepmind/alphafold/#citing-this-work) the [AlphaFold paper](https://doi.org/10.1038/s41586-021-03819-2).\n", - "\n", - "**Licenses**\n", - "\n", - "This Colab uses the [AlphaFold model parameters](https://github.com/deepmind/alphafold/#model-parameters-license) and its outputs are thus for non-commercial use only, under the Creative Commons Attribution-NonCommercial 4.0 International ([CC BY-NC 4.0](https://creativecommons.org/licenses/by-nc/4.0/legalcode)) license. The Colab itself is provided under the [Apache 2.0 license](https://www.apache.org/licenses/LICENSE-2.0). See the full license statement below.\n", - "\n", - "**More information**\n", - "\n", - "You can find more information about how AlphaFold works in our two Nature papers:\n", - "\n", - "* [AlphaFold methods paper](https://www.nature.com/articles/s41586-021-03819-2)\n", - "* [AlphaFold predictions of the human proteome paper](https://www.nature.com/articles/s41586-021-03828-1)\n", - "\n", - "FAQ on how to interpret AlphaFold predictions are [here](https://alphafold.ebi.ac.uk/faq)." - ] - }, - { - "cell_type": "code", - "execution_count": 0, - "metadata": { - "cellView": "form", - "id": "woIxeCPygt7K" - }, - "outputs": [], - "source": [ - "#@title Install third-party software\n", - "\n", - "#@markdown Please execute this cell by pressing the _Play_ button \n", - "#@markdown on the left to download and import third-party software \n", - "#@markdown in this Colab notebook. (See the [acknowledgements](https://github.com/deepmind/alphafold/#acknowledgements) in our readme.)\n", - "\n", - "#@markdown **Note**: This installs the software on the Colab \n", - "#@markdown notebook in the cloud and not on your computer.\n", - "\n", - "from IPython.utils import io\n", - "import os\n", - "import subprocess\n", - "import tqdm.notebook\n", - "\n", - "TQDM_BAR_FORMAT = '{l_bar}{bar}| {n_fmt}/{total_fmt} [elapsed: {elapsed} remaining: {remaining}]'\n", - "\n", - "try:\n", - " with tqdm.notebook.tqdm(total=100, bar_format=TQDM_BAR_FORMAT) as pbar:\n", - " with io.capture_output() as captured:\n", - " # Uninstall default Colab version of TF.\n", - " %shell pip uninstall -y tensorflow\n", - "\n", - " %shell sudo apt install --quiet --yes hmmer\n", - " pbar.update(6)\n", - "\n", - " # Install py3dmol.\n", - " %shell pip install py3dmol\n", - " pbar.update(2)\n", - "\n", - " # Install OpenMM and pdbfixer.\n", - " %shell rm -rf /opt/conda\n", - " %shell wget -q -P /tmp \\\n", - " https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh \\\n", - " && bash /tmp/Miniconda3-latest-Linux-x86_64.sh -b -p /opt/conda \\\n", - " && rm /tmp/Miniconda3-latest-Linux-x86_64.sh\n", - " pbar.update(9)\n", - "\n", - " PATH=%env PATH\n", - " %env PATH=/opt/conda/bin:{PATH}\n", - " %shell conda update -qy conda \\\n", - " && conda install -qy -c conda-forge \\\n", - " python=3.7 \\\n", - " openmm=7.5.1 \\\n", - " pdbfixer\n", - " pbar.update(80)\n", - "\n", - " # Create a ramdisk to store a database chunk to make Jackhmmer run fast.\n", - " %shell sudo mkdir -m 777 --parents /tmp/ramdisk\n", - " %shell sudo mount -t tmpfs -o size=9G ramdisk /tmp/ramdisk\n", - " pbar.update(2)\n", - "\n", - " %shell wget -q -P /content \\\n", - " https://git.scicore.unibas.ch/schwede/openstructure/-/raw/7102c63615b64735c4941278d92b554ec94415f8/modules/mol/alg/src/stereo_chemical_props.txt\n", - " pbar.update(1)\n", - "except subprocess.CalledProcessError:\n", - " print(captured)\n", - " raise" - ] - }, - { - "cell_type": "code", - "execution_count": 0, - "metadata": { - "cellView": "form", - "id": "VzJ5iMjTtoZw" - }, - "outputs": [], - "source": [ - "#@title Download AlphaFold\n", - "\n", - "#@markdown Please execute this cell by pressing the *Play* button on \n", - "#@markdown the left.\n", - "\n", - "GIT_REPO = 'https://github.com/deepmind/alphafold'\n", - "\n", - "SOURCE_URL = 'https://storage.googleapis.com/alphafold/alphafold_params_2021-07-14.tar'\n", - "PARAMS_DIR = './alphafold/data/params'\n", - "PARAMS_PATH = os.path.join(PARAMS_DIR, os.path.basename(SOURCE_URL))\n", - "\n", - "try:\n", - " with tqdm.notebook.tqdm(total=100, bar_format=TQDM_BAR_FORMAT) as pbar:\n", - " with io.capture_output() as captured:\n", - " %shell rm -rf alphafold\n", - " %shell git clone {GIT_REPO} alphafold\n", - " pbar.update(8)\n", - " # Install the required versions of all dependencies.\n", - " %shell pip3 install -r ./alphafold/requirements.txt\n", - " # Run setup.py to install only AlphaFold.\n", - " %shell pip3 install --no-dependencies ./alphafold\n", - " pbar.update(10)\n", - "\n", - " # Apply OpenMM patch.\n", - " %shell pushd /opt/conda/lib/python3.7/site-packages/ && \\\n", - " patch -p0 < /content/alphafold/docker/openmm.patch && \\\n", - " popd\n", - " \n", - " %shell mkdir -p /content/alphafold/common\n", - " %shell cp -f /content/stereo_chemical_props.txt /content/alphafold/common\n", - "\n", - " %shell mkdir --parents \"{PARAMS_DIR}\"\n", - " %shell wget -O \"{PARAMS_PATH}\" \"{SOURCE_URL}\"\n", - " pbar.update(27)\n", - "\n", - " %shell tar --extract --verbose --file=\"{PARAMS_PATH}\" \\\n", - " --directory=\"{PARAMS_DIR}\" --preserve-permissions\n", - " %shell rm \"{PARAMS_PATH}\"\n", - " pbar.update(55)\n", - "except subprocess.CalledProcessError:\n", - " print(captured)\n", - " raise\n", - "\n", - "import jax\n", - "if jax.local_devices()[0].platform == 'tpu':\n", - " raise RuntimeError('Colab TPU runtime not supported. Change it to GPU via Runtime -> Change Runtime Type -> Hardware accelerator -> GPU.')\n", - "elif jax.local_devices()[0].platform == 'cpu':\n", - " raise RuntimeError('Colab CPU runtime not supported. Change it to GPU via Runtime -> Change Runtime Type -> Hardware accelerator -> GPU.')" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "W4JpOs6oA-QS" - }, - "source": [ - "## Making a prediction\n", - "\n", - "Please paste the sequence of your protein in the text box below, then run the remaining cells via _Runtime_ > _Run after_. You can also run the cells individually by pressing the _Play_ button on the left.\n", - "\n", - "Note that the search against databases and the actual prediction can take some time, from minutes to hours, depending on the length of the protein and what type of GPU you are allocated by Colab (see FAQ below)." - ] - }, - { - "cell_type": "code", - "execution_count": 0, - "metadata": { - "cellView": "form", - "id": "rowN0bVYLe9n" - }, - "outputs": [], - "source": [ - "#@title Enter the amino acid sequence to fold ⬇️\n", - "sequence = 'MAAHKGAEHHHKAAEHHEQAAKHHHAAAEHHEKGEHEQAAHHADTAYAHHKHAEEHAAQAAKHDAEHHAPKPH' #@param {type:\"string\"}\n", - "\n", - "MIN_SEQUENCE_LENGTH = 16\n", - "MAX_SEQUENCE_LENGTH = 2500\n", - "\n", - "# Remove all whitespaces, tabs and end lines; upper-case\n", - "sequence = sequence.translate(str.maketrans('', '', ' \\n\\t')).upper()\n", - "aatypes = set('ACDEFGHIKLMNPQRSTVWY') # 20 standard aatypes\n", - "if not set(sequence).issubset(aatypes):\n", - " raise Exception(f'Input sequence contains non-amino acid letters: {set(sequence) - aatypes}. AlphaFold only supports 20 standard amino acids as inputs.')\n", - "if len(sequence) < MIN_SEQUENCE_LENGTH:\n", - " raise Exception(f'Input sequence is too short: {len(sequence)} amino acids, while the minimum is {MIN_SEQUENCE_LENGTH}')\n", - "if len(sequence) > MAX_SEQUENCE_LENGTH:\n", - " raise Exception(f'Input sequence is too long: {len(sequence)} amino acids, while the maximum is {MAX_SEQUENCE_LENGTH}. Please use the full AlphaFold system for long sequences.')" - ] - }, - { - "cell_type": "code", - "execution_count": 0, - "metadata": { - "cellView": "form", - "id": "2tTeTTsLKPjB" - }, - "outputs": [], - "source": [ - "#@title Search against genetic databases\n", - "\n", - "#@markdown Once this cell has been executed, you will see\n", - "#@markdown statistics about the multiple sequence alignment \n", - "#@markdown (MSA) that will be used by AlphaFold. In particular, \n", - "#@markdown you’ll see how well each residue is covered by similar \n", - "#@markdown sequences in the MSA.\n", - "\n", - "# --- Python imports ---\n", - "import sys\n", - "sys.path.append('/opt/conda/lib/python3.7/site-packages')\n", - "\n", - "import os\n", - "os.environ['TF_FORCE_UNIFIED_MEMORY'] = '1'\n", - "os.environ['XLA_PYTHON_CLIENT_MEM_FRACTION'] = '2.0'\n", - "\n", - "from urllib import request\n", - "from concurrent import futures\n", - "from google.colab import files\n", - "import json\n", - "from matplotlib import gridspec\n", - "import matplotlib.pyplot as plt\n", - "import numpy as np\n", - "import py3Dmol\n", - "\n", - "from alphafold.model import model\n", - "from alphafold.model import config\n", - "from alphafold.model import data\n", - "\n", - "from alphafold.data import parsers\n", - "from alphafold.data import pipeline\n", - "from alphafold.data.tools import jackhmmer\n", - "\n", - "from alphafold.common import protein\n", - "\n", - "from alphafold.relax import relax\n", - "from alphafold.relax import utils\n", - "\n", - "from IPython import display\n", - "from ipywidgets import GridspecLayout\n", - "from ipywidgets import Output\n", - "\n", - "# Color bands for visualizing plddt\n", - "PLDDT_BANDS = [(0, 50, '#FF7D45'),\n", - " (50, 70, '#FFDB13'),\n", - " (70, 90, '#65CBF3'),\n", - " (90, 100, '#0053D6')]\n", - "\n", - "# --- Find the closest source ---\n", - "test_url_pattern = 'https://storage.googleapis.com/alphafold-colab{:s}/latest/uniref90_2021_03.fasta.1'\n", - "ex = futures.ThreadPoolExecutor(3)\n", - "def fetch(source):\n", - " request.urlretrieve(test_url_pattern.format(source))\n", - " return source\n", - "fs = [ex.submit(fetch, source) for source in ['', '-europe', '-asia']]\n", - "source = None\n", - "for f in futures.as_completed(fs):\n", - " source = f.result()\n", - " ex.shutdown()\n", - " break\n", - "\n", - "# --- Search against genetic databases ---\n", - "with open('target.fasta', 'wt') as f:\n", - " f.write(f'>query\\n{sequence}')\n", - "\n", - "# Run the search against chunks of genetic databases (since the genetic\n", - "# databases don't fit in Colab ramdisk).\n", - "\n", - "jackhmmer_binary_path = '/usr/bin/jackhmmer'\n", - "dbs = []\n", - "\n", - "num_jackhmmer_chunks = {'uniref90': 59, 'smallbfd': 17, 'mgnify': 71}\n", - "total_jackhmmer_chunks = sum(num_jackhmmer_chunks.values())\n", - "with tqdm.notebook.tqdm(total=total_jackhmmer_chunks, bar_format=TQDM_BAR_FORMAT) as pbar:\n", - " def jackhmmer_chunk_callback(i):\n", - " pbar.update(n=1)\n", - "\n", - " pbar.set_description('Searching uniref90')\n", - " jackhmmer_uniref90_runner = jackhmmer.Jackhmmer(\n", - " binary_path=jackhmmer_binary_path,\n", - " database_path=f'https://storage.googleapis.com/alphafold-colab{source}/latest/uniref90_2021_03.fasta',\n", - " get_tblout=True,\n", - " num_streamed_chunks=num_jackhmmer_chunks['uniref90'],\n", - " streaming_callback=jackhmmer_chunk_callback,\n", - " z_value=135301051)\n", - " dbs.append(('uniref90', jackhmmer_uniref90_runner.query('target.fasta')))\n", - "\n", - " pbar.set_description('Searching smallbfd')\n", - " jackhmmer_smallbfd_runner = jackhmmer.Jackhmmer(\n", - " binary_path=jackhmmer_binary_path,\n", - " database_path=f'https://storage.googleapis.com/alphafold-colab{source}/latest/bfd-first_non_consensus_sequences.fasta',\n", - " get_tblout=True,\n", - " num_streamed_chunks=num_jackhmmer_chunks['smallbfd'],\n", - " streaming_callback=jackhmmer_chunk_callback,\n", - " z_value=65984053)\n", - " dbs.append(('smallbfd', jackhmmer_smallbfd_runner.query('target.fasta')))\n", - "\n", - " pbar.set_description('Searching mgnify')\n", - " jackhmmer_mgnify_runner = jackhmmer.Jackhmmer(\n", - " binary_path=jackhmmer_binary_path,\n", - " database_path=f'https://storage.googleapis.com/alphafold-colab{source}/latest/mgy_clusters_2019_05.fasta',\n", - " get_tblout=True,\n", - " num_streamed_chunks=num_jackhmmer_chunks['mgnify'],\n", - " streaming_callback=jackhmmer_chunk_callback,\n", - " z_value=304820129)\n", - " dbs.append(('mgnify', jackhmmer_mgnify_runner.query('target.fasta')))\n", - "\n", - "\n", - "# --- Extract the MSAs and visualize ---\n", - "# Extract the MSAs from the Stockholm files.\n", - "# NB: deduplication happens later in pipeline.make_msa_features.\n", - "\n", - "mgnify_max_hits = 501\n", - "\n", - "msas = []\n", - "deletion_matrices = []\n", - "full_msa = []\n", - "for db_name, db_results in dbs:\n", - " unsorted_results = []\n", - " for i, result in enumerate(db_results):\n", - " msa, deletion_matrix, target_names = parsers.parse_stockholm(result['sto'])\n", - " e_values_dict = parsers.parse_e_values_from_tblout(result['tbl'])\n", - " e_values = [e_values_dict[t.split('/')[0]] for t in target_names]\n", - " zipped_results = zip(msa, deletion_matrix, target_names, e_values)\n", - " if i != 0:\n", - " # Only take query from the first chunk\n", - " zipped_results = [x for x in zipped_results if x[2] != 'query']\n", - " unsorted_results.extend(zipped_results)\n", - " sorted_by_evalue = sorted(unsorted_results, key=lambda x: x[3])\n", - " db_msas, db_deletion_matrices, _, _ = zip(*sorted_by_evalue)\n", - " if db_msas:\n", - " if db_name == 'mgnify':\n", - " db_msas = db_msas[:mgnify_max_hits]\n", - " db_deletion_matrices = db_deletion_matrices[:mgnify_max_hits]\n", - " full_msa.extend(db_msas)\n", - " msas.append(db_msas)\n", - " deletion_matrices.append(db_deletion_matrices)\n", - " msa_size = len(set(db_msas))\n", - " print(f'{msa_size} Sequences Found in {db_name}')\n", - "\n", - "deduped_full_msa = list(dict.fromkeys(full_msa))\n", - "total_msa_size = len(deduped_full_msa)\n", - "print(f'\\n{total_msa_size} Sequences Found in Total\\n')\n", - "\n", - "aa_map = {restype: i for i, restype in enumerate('ABCDEFGHIJKLMNOPQRSTUVWXYZ-')}\n", - "msa_arr = np.array([[aa_map[aa] for aa in seq] for seq in deduped_full_msa])\n", - "num_alignments, num_res = msa_arr.shape\n", - "\n", - "fig = plt.figure(figsize=(12, 3))\n", - "plt.title('Per-Residue Count of Non-Gap Amino Acids in the MSA')\n", - "plt.plot(np.sum(msa_arr != aa_map['-'], axis=0), color='black')\n", - "plt.ylabel('Non-Gap Count')\n", - "plt.yticks(range(0, num_alignments + 1, max(1, int(num_alignments / 3))))\n", - "plt.show()" - ] - }, - { - "cell_type": "code", - "execution_count": 0, - "metadata": { - "cellView": "form", - "id": "XUo6foMQxwS2" - }, - "outputs": [], - "source": [ - "#@title Run AlphaFold and download prediction\n", - "\n", - "#@markdown Once this cell has been executed, a zip-archive with \n", - "#@markdown the obtained prediction will be automatically downloaded \n", - "#@markdown to your computer.\n", - "\n", - "# --- Run the model ---\n", - "model_names = ['model_1', 'model_2', 'model_3', 'model_4', 'model_5', 'model_2_ptm']\n", - "\n", - "def _placeholder_template_feats(num_templates_, num_res_):\n", - " return {\n", - " 'template_aatype': np.zeros([num_templates_, num_res_, 22], np.float32),\n", - " 'template_all_atom_masks': np.zeros([num_templates_, num_res_, 37, 3], np.float32),\n", - " 'template_all_atom_positions': np.zeros([num_templates_, num_res_, 37], np.float32),\n", - " 'template_domain_names': np.zeros([num_templates_], np.float32),\n", - " 'template_sum_probs': np.zeros([num_templates_], np.float32),\n", - " }\n", - "\n", - "output_dir = 'prediction'\n", - "os.makedirs(output_dir, exist_ok=True)\n", - "\n", - "plddts = {}\n", - "pae_outputs = {}\n", - "unrelaxed_proteins = {}\n", - "\n", - "with tqdm.notebook.tqdm(total=len(model_names) + 1, bar_format=TQDM_BAR_FORMAT) as pbar:\n", - " for model_name in model_names:\n", - " pbar.set_description(f'Running {model_name}')\n", - " num_templates = 0\n", - " num_res = len(sequence)\n", - "\n", - " feature_dict = {}\n", - " feature_dict.update(pipeline.make_sequence_features(sequence, 'test', num_res))\n", - " feature_dict.update(pipeline.make_msa_features(msas, deletion_matrices=deletion_matrices))\n", - " feature_dict.update(_placeholder_template_feats(num_templates, num_res))\n", - "\n", - " cfg = config.model_config(model_name)\n", - " params = data.get_model_haiku_params(model_name, './alphafold/data')\n", - " model_runner = model.RunModel(cfg, params)\n", - " processed_feature_dict = model_runner.process_features(feature_dict,\n", - " random_seed=0)\n", - " prediction_result = model_runner.predict(processed_feature_dict)\n", - "\n", - " mean_plddt = prediction_result['plddt'].mean()\n", - "\n", - " if 'predicted_aligned_error' in prediction_result:\n", - " pae_outputs[model_name] = (\n", - " prediction_result['predicted_aligned_error'],\n", - " prediction_result['max_predicted_aligned_error']\n", - " )\n", - " else:\n", - " # Get the pLDDT confidence metrics. Do not put pTM models here as they\n", - " # should never get selected.\n", - " plddts[model_name] = prediction_result['plddt']\n", - "\n", - " # Set the b-factors to the per-residue plddt.\n", - " final_atom_mask = prediction_result['structure_module']['final_atom_mask']\n", - " b_factors = prediction_result['plddt'][:, None] * final_atom_mask\n", - " unrelaxed_protein = protein.from_prediction(processed_feature_dict,\n", - " prediction_result,\n", - " b_factors=b_factors)\n", - " unrelaxed_proteins[model_name] = unrelaxed_protein\n", - "\n", - " # Delete unused outputs to save memory.\n", - " del model_runner\n", - " del params\n", - " del prediction_result\n", - " pbar.update(n=1)\n", - "\n", - " # --- AMBER relax the best model ---\n", - " pbar.set_description(f'AMBER relaxation')\n", - " amber_relaxer = relax.AmberRelaxation(\n", - " max_iterations=0,\n", - " tolerance=2.39,\n", - " stiffness=10.0,\n", - " exclude_residues=[],\n", - " max_outer_iterations=20)\n", - " # Find the best model according to the mean pLDDT.\n", - " best_model_name = max(plddts.keys(), key=lambda x: plddts[x].mean())\n", - " relaxed_pdb, _, _ = amber_relaxer.process(\n", - " prot=unrelaxed_proteins[best_model_name])\n", - " pbar.update(n=1) # Finished AMBER relax.\n", - "\n", - "# Construct multiclass b-factors to indicate confidence bands\n", - "# 0=very low, 1=low, 2=confident, 3=very high\n", - "banded_b_factors = []\n", - "for plddt in plddts[best_model_name]:\n", - " for idx, (min_val, max_val, _) in enumerate(PLDDT_BANDS):\n", - " if plddt >= min_val and plddt <= max_val:\n", - " banded_b_factors.append(idx)\n", - " break\n", - "banded_b_factors = np.array(banded_b_factors)[:, None] * final_atom_mask\n", - "to_visualize_pdb = utils.overwrite_b_factors(relaxed_pdb, banded_b_factors)\n", - "\n", - "\n", - "# Write out the prediction\n", - "pred_output_path = os.path.join(output_dir, 'selected_prediction.pdb')\n", - "with open(pred_output_path, 'w') as f:\n", - " f.write(relaxed_pdb)\n", - "\n", - "\n", - "# --- Visualise the prediction & confidence ---\n", - "show_sidechains = True\n", - "def plot_plddt_legend():\n", - " \"\"\"Plots the legend for pLDDT.\"\"\"\n", - " thresh = [\n", - " 'Very low (pLDDT < 50)',\n", - " 'Low (70 > pLDDT > 50)',\n", - " 'Confident (90 > pLDDT > 70)',\n", - " 'Very high (pLDDT > 90)']\n", - "\n", - " colors = [x[2] for x in PLDDT_BANDS]\n", - "\n", - " plt.figure(figsize=(2, 2))\n", - " for c in colors:\n", - " plt.bar(0, 0, color=c)\n", - " plt.legend(thresh, frameon=False, loc='center', fontsize=20)\n", - " plt.xticks([])\n", - " plt.yticks([])\n", - " ax = plt.gca()\n", - " ax.spines['right'].set_visible(False)\n", - " ax.spines['top'].set_visible(False)\n", - " ax.spines['left'].set_visible(False)\n", - " ax.spines['bottom'].set_visible(False)\n", - " plt.title('Model Confidence', fontsize=20, pad=20)\n", - " return plt\n", - "\n", - "# Color the structure by per-residue pLDDT\n", - "color_map = {i: bands[2] for i, bands in enumerate(PLDDT_BANDS)}\n", - "view = py3Dmol.view(width=800, height=600)\n", - "view.addModelsAsFrames(to_visualize_pdb)\n", - "style = {'cartoon': {\n", - " 'colorscheme': {\n", - " 'prop': 'b',\n", - " 'map': color_map}\n", - " }}\n", - "if show_sidechains:\n", - " style['stick'] = {}\n", - "view.setStyle({'model': -1}, style)\n", - "view.zoomTo()\n", - "\n", - "grid = GridspecLayout(1, 2)\n", - "out = Output()\n", - "with out:\n", - " view.show()\n", - "grid[0, 0] = out\n", - "\n", - "out = Output()\n", - "with out:\n", - " plot_plddt_legend().show()\n", - "grid[0, 1] = out\n", - "\n", - "display.display(grid)\n", - "\n", - "# Display pLDDT and predicted aligned error (if output by the model).\n", - "if pae_outputs:\n", - " num_plots = 2\n", - "else:\n", - " num_plots = 1\n", - "\n", - "plt.figure(figsize=[8 * num_plots, 6])\n", - "plt.subplot(1, num_plots, 1)\n", - "plt.plot(plddts[best_model_name])\n", - "plt.title('Predicted LDDT')\n", - "plt.xlabel('Residue')\n", - "plt.ylabel('pLDDT')\n", - "\n", - "if num_plots == 2:\n", - " plt.subplot(1, 2, 2)\n", - " pae, max_pae = list(pae_outputs.values())[0]\n", - " plt.imshow(pae, vmin=0., vmax=max_pae, cmap='Greens_r')\n", - " plt.colorbar(fraction=0.046, pad=0.04)\n", - " plt.title('Predicted Aligned Error')\n", - " plt.xlabel('Scored residue')\n", - " plt.ylabel('Aligned residue')\n", - "\n", - "# Save pLDDT and predicted aligned error (if it exists)\n", - "pae_output_path = os.path.join(output_dir, 'predicted_aligned_error.json')\n", - "if pae_outputs:\n", - " # Save predicted aligned error in the same format as the AF EMBL DB\n", - " rounded_errors = np.round(pae.astype(np.float64), decimals=1)\n", - " indices = np.indices((len(rounded_errors), len(rounded_errors))) + 1\n", - " indices_1 = indices[0].flatten().tolist()\n", - " indices_2 = indices[1].flatten().tolist()\n", - " pae_data = json.dumps([{\n", - " 'residue1': indices_1,\n", - " 'residue2': indices_2,\n", - " 'distance': rounded_errors.flatten().tolist(),\n", - " 'max_predicted_aligned_error': max_pae.item()\n", - " }],\n", - " indent=None,\n", - " separators=(',', ':'))\n", - " with open(pae_output_path, 'w') as f:\n", - " f.write(pae_data)\n", - "\n", - "\n", - "# --- Download the predictions ---\n", - "!zip -q -r {output_dir}.zip {output_dir}\n", - "files.download(f'{output_dir}.zip')" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "lUQAn5LYC5n4" - }, - "source": [ - "### Interpreting the prediction\n", - "\n", - "Please see the [AlphaFold methods paper](https://www.nature.com/articles/s41586-021-03819-2) and the [AlphaFold predictions of the human proteome paper](https://www.nature.com/articles/s41586-021-03828-1), as well as [our FAQ](https://alphafold.ebi.ac.uk/faq) on how to interpret AlphaFold predictions." - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "jeb2z8DIA4om" - }, - "source": [ - "## FAQ & Troubleshooting\n", - "\n", - "\n", - "* How do I get a predicted protein structure for my protein?\n", - " * Click on the _Connect_ button on the top right to get started.\n", - " * Paste the amino acid sequence of your protein (without any headers) into the “Enter the amino acid sequence to fold”.\n", - " * Run all cells in the Colab, either by running them individually (with the play button on the left side) or via _Runtime_ > _Run all._\n", - " * The predicted protein structure will be downloaded once all cells have been executed. Note: This can take minutes to hours - see below.\n", - "* How long will this take?\n", - " * Downloading the AlphaFold source code can take up to a few minutes.\n", - " * Downloading and installing the third-party software can take up to a few minutes.\n", - " * The search against genetic databases can take minutes to hours.\n", - " * Running AlphaFold and generating the prediction can take minutes to hours, depending on the length of your protein and on which GPU-type Colab has assigned you.\n", - "* My Colab no longer seems to be doing anything, what should I do?\n", - " * Some steps may take minutes to hours to complete.\n", - " * If nothing happens or if you receive an error message, try restarting your Colab runtime via _Runtime_ > _Restart runtime_.\n", - " * If this doesn’t help, try resetting your Colab runtime via _Runtime_ > _Factory reset runtime_.\n", - "* How does this compare to the open-source version of AlphaFold?\n", - " * This Colab version of AlphaFold searches a selected portion of the BFD dataset and currently doesn’t use templates, so its accuracy is reduced in comparison to the full version of AlphaFold that is described in the [AlphaFold paper](https://doi.org/10.1038/s41586-021-03819-2) and [Github repo](https://github.com/deepmind/alphafold/) (the full version is available via the inference script).\n", - "* What is a Colab?\n", - " * See the [Colab FAQ](https://research.google.com/colaboratory/faq.html).\n", - "* I received a warning “Notebook requires high RAM”, what do I do?\n", - " * The resources allocated to your Colab vary. See the [Colab FAQ](https://research.google.com/colaboratory/faq.html) for more details.\n", - " * You can execute the Colab nonetheless.\n", - "* I received an error “Colab CPU runtime not supported” or “No GPU/TPU found”, what do I do?\n", - " * Colab CPU runtime is not supported. Try changing your runtime via _Runtime_ > _Change runtime type_ > _Hardware accelerator_ > _GPU_.\n", - " * The type of GPU allocated to your Colab varies. See the [Colab FAQ](https://research.google.com/colaboratory/faq.html) for more details.\n", - " * If you receive “Cannot connect to GPU backend”, you can try again later to see if Colab allocates you a GPU.\n", - " * [Colab Pro](https://colab.research.google.com/signup) offers priority access to GPUs. \n", - "* Does this tool install anything on my computer?\n", - " * No, everything happens in the cloud on Google Colab.\n", - " * At the end of the Colab execution a zip-archive with the obtained prediction will be automatically downloaded to your computer.\n", - "* How should I share feedback and bug reports?\n", - " * Please share any feedback and bug reports as an [issue](https://github.com/deepmind/alphafold/issues) on Github.\n", - "\n", - "\n", - "## Related work\n", - "\n", - "Take a look at these Colab notebooks provided by the community (please note that these notebooks may vary from our validated AlphaFold system and we cannot guarantee their accuracy):\n", - "\n", - "* The [ColabFold AlphaFold2 notebook](https://colab.research.google.com/github/sokrypton/ColabFold/blob/main/AlphaFold2.ipynb) by Sergey Ovchinnikov, Milot Mirdita and Martin Steinegger, which uses an API hosted at the Södinglab based on the MMseqs2 server ([Mirdita et al. 2019, Bioinformatics](https://academic.oup.com/bioinformatics/article/35/16/2856/5280135)) for the multiple sequence alignment creation.\n" - ] - }, - { - "cell_type": "markdown", - "metadata": { - "id": "YfPhvYgKC81B" - }, - "source": [ - "# License and Disclaimer\n", - "\n", - "This is not an officially-supported Google product.\n", - "\n", - "This Colab notebook and other information provided is for theoretical modelling only, caution should be exercised in its use. It is provided ‘as-is’ without any warranty of any kind, whether expressed or implied. Information is not intended to be a substitute for professional medical advice, diagnosis, or treatment, and does not constitute medical or other professional advice.\n", - "\n", - "Copyright 2021 DeepMind Technologies Limited.\n", - "\n", - "\n", - "## AlphaFold Code License\n", - "\n", - "Licensed under the Apache License, Version 2.0 (the \"License\"); you may not use this file except in compliance with the License. You may obtain a copy of the License at https://www.apache.org/licenses/LICENSE-2.0.\n", - "\n", - "Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an \"AS IS\" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License.\n", - "\n", - "## Model Parameters License\n", - "\n", - "The AlphaFold parameters are made available for non-commercial use only, under the terms of the Creative Commons Attribution-NonCommercial 4.0 International (CC BY-NC 4.0) license. You can find details at: https://creativecommons.org/licenses/by-nc/4.0/legalcode\n", - "\n", - "\n", - "## Third-party software\n", - "\n", - "Use of the third-party software, libraries or code referred to in the [Acknowledgements section](https://github.com/deepmind/alphafold/#acknowledgements) in the AlphaFold README may be governed by separate terms and conditions or license provisions. Your use of the third-party software, libraries or code is subject to any such terms and you should check that you can comply with any applicable restrictions or terms and conditions before use.\n", - "\n", - "\n", - "## Mirrored Databases\n", - "\n", - "The following databases have been mirrored by DeepMind, and are available with reference to the following:\n", - "* UniRef90: v2021\\_03 (unmodified), by The UniProt Consortium, available under a [Creative Commons Attribution-NoDerivatives 4.0 International License](http://creativecommons.org/licenses/by-nd/4.0/).\n", - "* MGnify: v2019\\_05 (unmodified), by Mitchell AL et al., available free of all copyright restrictions and made fully and freely available for both non-commercial and commercial use under [CC0 1.0 Universal (CC0 1.0) Public Domain Dedication](https://creativecommons.org/publicdomain/zero/1.0/).\n", - "* BFD: (modified), by Steinegger M. and Söding J., modified by DeepMind, available under a [Creative Commons Attribution-ShareAlike 4.0 International License](https://creativecommons.org/licenses/by/4.0/). See the Methods section of the [AlphaFold proteome paper](https://www.nature.com/articles/s41586-021-03828-1) for details." - ] - } - ], - "metadata": { - "accelerator": "GPU", - "colab": { - "collapsed_sections": [], - "name": "AlphaFold.ipynb" - }, - "kernelspec": { - "display_name": "Python 3", - "name": "python3" - }, - "language_info": { - "name": "python" - } - }, - "nbformat": 4, - "nbformat_minor": 0 -} diff --git a/run_complexfold.py b/run_complexfold.py index c03b8cf..90607d8 100755 --- a/run_complexfold.py +++ b/run_complexfold.py @@ -39,7 +39,6 @@ # - Added colabfold.plot_protein """Full AlphaFold protein structure prediction script.""" -import json import os import pathlib import pickle @@ -61,16 +60,10 @@ from alphafold.model import model from alphafold.model import features from alphafold.relax import relax -from alphafold.data import colabfold from alphafold import complexfold import jax import numpy as np -import matplotlib.pyplot as plt -#import dill -#from pathos.multiprocessing import ProcessingPool as Pool -from multiprocessing import Pool -from multiprocessing import Process # Internal import (7716). flags.DEFINE_list('fasta_paths', None, 'Paths to FASTA files, each containing '