pytorch chemistry biology translation license:apache-2.0 region:us
---
license: apache-2.0
pipeline_tag: translation
tags:
- chemistry
- biology
inference: false
---

# **Contributors**

- Sebastian Lindner (GitHub [@Bienenwolf655](https://github.com/Bienenwolf655); Twitter [@lindner_seb](https://twitter.com/lindner_seb)) 
- Michael Heinzinger (GitHub [@mheinzinger](https://github.com/mheinzinger); Twitter [@HeinzingerM](https://twitter.com/HeinzingerM))
- Noelia Ferruz (GitHub [@noeliaferruz](https://github.com/noeliaferruz); Twitter [@ferruz_noelia](https://twitter.com/ferruz_noelia); Webpage: [www.aiproteindesign.com](https://www.aiproteindesign.com) )


# **REXzyme: A Translation Machine for the Generation of New-to-Nature Enzymes**
**Work in Progress**

REXzyme (Reaction to Enzyme) (manuscript in preparation) is a translation machine -similar to Google Translator-
for the generation of enzymes that catalize user-defined reactions.

![Inference of REXzyme](./rexzyme3.png)

It is possible to provide fine-grained input at the substrate level. 
Akin to how translation machines have learned to translate between complex language pairs with great success, 
often diverging in their representation at the character level (Japanese - English), we posit that an advanced architecture will
be able to translate between the chemical and sequence spaces. REXzyme was trained on a set of 2480 reactions 
and ~32M enzyme pairs and it produces sequences that are predicted to perform their intended reactions. A second 
version of the model with 14k more reactions will be uploaded to this repository shortly.

 you will need to provide a reaction in the SMILES format
(Simplified molecular-input line-entry system). A useful online server to convert from molecules to SMILES
can be found here: https://cactus.nci.nih.gov/chemical/structure.

After converting each of the reaction components you should combine them in the following scheme: ```ReactantA.ReactantB>AgentA>ProductA.ProductB```<br/>
Additionally, one should prepend the task suffix ```r2s``` and append the eos token ```</s>```
e.g. for the carbonic anhydrase reaction: ```r2sO.COO>>HCOOO.[H+]</s>```

We provide this python script to convert reactants to the required reaction format, but 
we always recommend to draw and double-check the structures in a server like [cactus](https://cactus.nci.nih.gov/chemical/structure)

```python
#  left reactants (seperated by '.') seperated by a equal sign from the products (also seperated by '.')
reactions =  "CO2 . H2O =  carbonic acid . H+"
# agents (seperated by .) 
agent = ""

# https://stackoverflow.com/questions/54930121/converting-molecule-name-to-smiles
from urllib.request import urlopen
from urllib.parse import quote

def CIRconvert(ids):
    try:
        url = 'http://cactus.nci.nih.gov/chemical/structure/' + quote(ids) + '/smiles'
        ans = urlopen(url).read().decode('utf8')
        return ans
    except:
        return 'Did not work'

reagent = [CIRconvert(i) for i in reactions.replace(' ','').split('=')[0].split('.') if i != ""]
agent = [CIRconvert(i) for i in agent.replace(' ','').split('.') if i != ""]
product = [CIRconvert(i) for i in reactions.replace(' ','').split('=')[1].split('.') if i != ""]
f"r2s{'.'.join(reagent)}>{'.'.join(agent)}>{'.'.join(product)}</s>"
```

We are still working in the analysis of the model for different tasks, including experimental testing. 
See below in this documentation information about the models' performance in different in-silico tasks and how to generate your own enzymes. 

## **Model description**

REXzyme is based on the [Efficient T5 Large Transformer](https://huggingface.co/google/t5-efficient-large) architecture (which in turn is very similar to the current version of Google Translator)
and contains 48 (24 encoder/ 24 decoder) layers with a model dimensionality of 1024, totaling 770 million parameters. 

REXzyme is a translation machine trained on portion the [RHEA database](https://www.rhea-db.org/) containing 31,970,152 reaction-enzyme pairs. 
A second dataset with >14k reactions is being trained and will be uploaded soon.
The pre-training was done on pairs of SMILES and amino acid sequences, tokenized with a char-level 
Sentencepiece tokenizer. Note that two seperate tokenizers were used for input (./tokenizer_smiles) and labels (./tokenizer_aa).

REXzyme was pre-trained with a supervised translation objective  i.e., the model learned to process the continous
representation of the reaction from the encoder to autoregressively (causual language modeling) produce the output. 
The output tokens (amino acids) are generated one at a time, from left to right, and the model learns to match the original enzyme sequence.
Hence, the model learns the dependencies among protein sequence features that enable a specific enzymatic reaction.

There are stark differences in the number of members among reaction classes.
However, since we are tokenizing the reaction SMILES on a character level, 
the model has learnt dependencies among molecules and enzyme sequence features, and it can transfer learning from more to less populated
reaction classes.

## **Model Performance**

- **Dataset curation**
We converted the reactions from rxn format to smile string including only left-to-right reactions.
The enzyme sequences were truncated to 1024.
Enzymes catalyzing more than one reaction appear in multiple enzyme-reaction pairs.
<br/><br/>
- **General descriptors**

    | Method                | Natural     | Generated <sup>[1]</sup> |
    | :---                  |    :----:   |          ---:            |
    | **IUPRED3 (ordered)** | 99.9%       | 99.9%                    |
    | **ESMFold (avg. plddt)**            | 85.03       | 79.82  |
    | **FlDPnn**            | 0.0878      | 0.0929                  |
<sup>[1]|</sup> We excluded sequences with %identities ≥ 70% and pLDDTs < 60%.
<br/><br/>


- **Functional classification**
<br/><br/>
<table>
  <tr>
    <td><b>Method </b></td>
    <td colspan="2"> <a href="https://google-research.github.io/proteinfer/">ProteInfer</a></td>
    <td colspan="2"> <a href="https://www.science.org/doi/10.1126/science.adf2465">CLEAN</a></td>
  </tr>
  <tr>
    <td><b>Dataset</b></td>
    <td >Natural (%) </td>
    <td >Generated (%) </td>
    <td >Natural (%) </td>
    <td >Generated (%) </td>
  </tr>
    <tr>
    <td><b> EC: Level 1 </b></td>
    <td >81</td>
    <td >80</td>
    <td >80</td>
    <td >79</td>
  </tr>
    <tr>
    <td><b> EC: Level 2 </b></td>
    <td >78</td>
    <td >77</td>
    <td >79</td>
    <td >78</td>
  </tr>
    <tr>
    <td><b> EC: Level 3 </b></td>
    <td >76</td>
    <td >75</td>
    <td >78</td>
    <td >77</td>
  </tr>
    <tr>
    <td><b> EC: Level 4 </b></td>
    <td >62</td>
    <td >58</td>
    <td >70</td>
    <td >65</td>
  </tr>
    <tr>
    <td><b> No EC predicted </b></td>
    <td >10</td>
    <td >7</td>
    <td >0</td>
    <td >0</td>
  </tr>
    <tr>
    <td><b> GO-Terms </b></td>
    <td >41</td>
    <td >39</td>
    <td >-</td>
    <td >-</td>
  </tr>
    <tr>
    <td><b> No GO predicted </b></td>
    <td >1</td>
    <td >1</td>
    <td >-</td>
    <td >-</td>
  </tr>
</table>
<br/><br/>
- **PGP pipeline** [(see GitHub)](https://github.com/hefeda/PGP)

    | Method      | Natural | Generated |
    | :---        | :----   |      :--- |
    | **Disorder**      | 11.473       | 11.467   |
    | **DSSP3**   | L: 42%, H: 41%, E:18% | L: 45%, H: 39%, E: 16%|
    | **DSSP8**   | C:25%, H:38% T:10%, S:5%, I:0%, E:19%, G:2%, B:0% | C:29%, H:38% T:10%, S:4%, I:0%, E:17%, G:3%, B:0%|
    | **CATH Classes**   | Mainly Beta: 6%, Alpha Beta: 78%, Mainly Alpha: 16%, Special: 0%, Few Secondary Structures: 0% |  Mainly Beta: 4%, Alpha Beta: 87%, Mainly Alpha: 9%, Special: 0%, Few Secondary Structures: 0%|
    | **Transmembrane Prediction**  | Membrane: 9%, Soluble: 91% | Membrane: 9%, Soluble: 91%|
    | **Conservation**      | High: 37%, Low: 33%      | High: 38%, Low: 33%  |
    | **Localization**   | Cytop.: 66%, Nucleus: 4%, Extracellular: 6%, PM: 4%, ER: 11%, Lysosome/Vacuole: 1%, Mito.: 6%, Plastid: 1%, Golgi: 1%, Perox.: 1% | Cytop.: 85%, Nucleus: 2%, Extracellular: 6%, PM: 1%, ER: 6%, Lysosome/Vacuole: 0%, Mito.: 4%, Plastid: 0%, Golgi: 0%, Perox.: 0%|
<br/><br/>


## **How to generate from REXzyme**
REXzyme can be used with the HuggingFace transformer python package.
Detailed installation instructions can be found [here](https://huggingface.co/docs/transformers/installation).

Since REXzyme has been trained on the objective of machine translation, users have to specify a chemical reaction, specified in the format of SMILES.

Disclaimer: Although the perplexity gets computed here it is not the best selection criteria.
Usually the BLEU score is deployed for translation evaluation, 
but this score would enforce a high sequence similarity (thus not *de novo* design, which is what we tend to go for). 
We recommend generating many sequences and selecting them by plDDT, as well as other metrics.

```python
from datasets import load_from_disk
from transformers import AutoTokenizer
from transformers import T5Tokenizer, T5ForConditionalGeneration
import math
import torch
from tqdm import tqdm
import pickle
tokenizer_aa = AutoTokenizer.from_pretrained('/path/to//tokenizer_aa')
tokenizer_smiles = AutoTokenizer.from_pretrained('/path/to//tokenizer_smiles')

model = T5ForConditionalGeneration.from_pretrained("/path/to/REXzyme").cuda()
print(model.generation_config)
reactions = ["NC1=NC=NC2=C1N=CN2[C@@H]1O[C@H](COP(=O)([O-])OP(=O)([O-])OP(=O)([O-])[O-])[C@@H](O)[C@H]1O.*N[C@@H](CO)C(*)=O>>NC1=NC=NC2=C1N=CN2[C@@H]1O[C@H](COP(=O)([O-])OP(=O)([O-])[O-])[C@@H](O)[C@H]1O.[H+].*N[C@@H](COP(=O)([O-])[O-])C(*)=O"]

def calculatePerplexity(inputs,model):
    '''Function to compute perplexity'''
    a=tokenizer_aa.decode(inputs)
    b=tokenizer_aa(a, return_tensors="pt").input_ids.to(device='cuda')
    b = torch.stack([[b[b!=tokenizer_aa.pad_token_id]] for label in b][0])
    with torch.no_grad():
        outputs = model(b, labels=b)
    loss, logits = outputs[:2]
    return math.exp(loss)


for idx,i in tqdm(enumerate(reactions)):
    input_ids = tokenizer_smiles(f"r2s{i}</s>", return_tensors="pt").input_ids.to(device='cuda')
    print(f'Generating for {i}')
    ppls_total = []
    for _ in range(4):
        outputs = model.generate(input_ids,
                top_k=15,
                top_p = 0.92,
                repetition_penalty=1.2,
                max_length=1024,
                do_sample=True,
                num_return_sequences=25)
        ppls = [(tokenizer_aa.decode(output,skip_special_tokens=True), calculatePerplexity(output, model),len(tokenizer_aa.decode(output))) for output in tqdm(outputs)]
        ppls_total.extend(ppls)
```

## **A word of caution** 

- We have not yet fully tested the ability of the model for the generation of new-to-nature enzymes, i.e.,
  with chemical reactions that do not appear in Nature (and hence neither in the training set). While this is the intended objective of our work, it is very much work in progress. We'll uptadate the model and documentation shortly.