dauparas / ProteinMPNN

Code for the ProteinMPNN paper
MIT License
910 stars 278 forks source link

Indexing with Missing Residues #77

Open brucejwittmann opened 9 months ago

brucejwittmann commented 9 months ago

Hey! Another quick question: How does indexing work for ProteinMPNN when there are missing residues present in the PDB file? It looks like it will add "X" to the output sequences where it identifies gaps in the residue indexing, but how is that reflected in indexing? For instance, if I had a PDB file with 2 residues present but 5 missing in between them, I would get sequences of the below format out:

NXXXXXN

The first "N" will have index 1, but what about the second? Is it "2" (because it's the second residue in the PDB file) or "7" (because it's the seventh index in the PDB file if we account for the missing residues). I'm getting some weird outputs when I tie residues together assuming the first indexing scheme (where the Ns are indices 1 and 2 rather than 1 and 7).

Extending the above, what if I have missing residues at the beginning or end of the PDB sequence? For instance, say in the same example above I have 2 residues missing from the beginning and 3 missing from the end in addition to the 3 in the middle. I believe I would get exactly the same result as the first case (NXXXXXN), but would the indices of the first "N" change? Would it be "3" now to account for the added residues at the front?

tony-res commented 8 months ago

From what I've seen the second "N" would be at position 2. What I do is I just replace the "X" with "" and then count.

sequence_no_gaps = sequence_with_gaps.replace("X","")

That said, I forget if the index starts at 0 or at 1. I think it starts at 0. So

--position_list "0 1"

rather than

--position_list "1 2"

As a sanity check I always use a python script to go through the output .fa file and calculate the changes. A good Python library for that is from collections import Counter.

For example,

>>> from collections import Counter
>>> word = "mississippi"
>>> counter = {}

>>> for letter in word:
            counter[letter] = counter.get(letter, 0) + 1

>>> counter
{'m': 1, 'i': 4, 's': 4, 'p': 2}

Some cool code to handle everything: https://huggingface.co/spaces/simonduerr/ProteinMPNN/blob/main/app.py#L382

jadolfbr commented 8 months ago

protein_mpnn_utils.py, line 125; I commented out a few of the places where 'gaps' are being added for sequence numbering. Becomes this:

try:
    for resn in range(min_resn,max_resn+1):
      if resn in seq:
        for k in sorted(seq[resn]): seq_.append(aa_3_N.get(seq[resn][k],20))
      #else: seq_.append(20)
      if resn in xyz:
        for k in sorted(xyz[resn]):
          for atom in atoms:
            if atom in xyz[resn][k]: xyz_.append(xyz[resn][k][atom])
            else: xyz_.append(np.full(3,np.nan))
      #else:
      #  for atom in atoms: xyz_.append(np.full(3,np.nan))
    print("Seq_", seq_)
brucejwittmann commented 8 months ago

Thanks everyone. @jadolfbr, based on the code you referenced, it looks like the logic is as follows for indexing:

  1. Find the smallest and largest residue indices in the PDB file, extracting from the coordinates section.
  2. Scan over the range of integers between the smallest and largest integers, recording masks/pads for missing residues.

Based on that logic, in my example with "NXXXXXN", the first "N" would be index "1" and the second index "7" assuming that there were 5 indices between them. If I were to add missing residues at the beginning or end, I would get the same indexing scheme as the first case as those residues will not show up in the coordinates section of the PDB file.

@tony-res, I believe indexing should start at 1 for ProteinMPNN. If you look at example 4 you can see that noted in a comment. Also you can see when masking fixed positions that "1" is subtracted from the input array.

I think it's important to note that the index in the code @jadolfbr copied is not the same as the index for freezing positions. The index in the copied code is direct from the PDB file, while the indices in the example I linked are 1-indexed to the relative position of a residue in the file.

jadolfbr commented 8 months ago

Yes, nice.

Also, it’s 1-indexed for each chain - not necessarily for the whole PDB.

On Mon, Dec 4, 2023 at 8:47 PM Bruce Wittmann @.***> wrote:

Thanks everyone. @jadolfbr https://github.com/jadolfbr, based on the code you referenced, it looks like the logic is as follows for indexing:

  1. Find the smallest and largest residue indices in the PDB file, extracting from the coordinates section.
  2. Scan over the range of integers between the smallest and largest integers, recording masks/pads for missing residues.

Based on that logic, in my example with "NXXXXXN", the first "N" would be index "1" and the second index "7" assuming that there were 5 indices between them. If I were to add missing residues at the beginning or end, I would get the same indexing scheme as the first case as those residues will not show up in the coordinates section of the PDB file.

@tony-res https://github.com/tony-res, I believe indexing should start at 1 for ProteinMPNN. If you look at example 4 https://github.com/dauparas/ProteinMPNN/blob/8907e6671bfbfc92303b5f79c4b5e6ce47cdef57/examples/submit_example_4.sh#L23C1-L23C88 you can see that noted in a comment. Also you can see when masking fixed positions https://github.com/dauparas/ProteinMPNN/blob/8907e6671bfbfc92303b5f79c4b5e6ce47cdef57/protein_mpnn_utils.py#L310 that "1" is subtracted from the input array.

I think it's important to note that the index in the code @jadolfbr https://github.com/jadolfbr copied is not the same as the index for freezing positions. The index in the copied code is direct from the PDB file, while the indices in the example I linked are 1-indexed to the relative position of a residue in the file.

— Reply to this email directly, view it on GitHub https://github.com/dauparas/ProteinMPNN/issues/77#issuecomment-1839860158, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAZDHREF2YATJSRL43GQZK3YHZ4JLAVCNFSM6AAAAAA7QLOLHSVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTQMZZHA3DAMJVHA . You are receiving this because you were mentioned.Message ID: @.***>