PhilippJunk / homelette

A unified and modular interface to homology modelling software
MIT License
10 stars 2 forks source link

Sequence name of fasta file can not exceed 14 character #5

Closed wt12318 closed 1 year ago

wt12318 commented 1 year ago

Hi,

It seems that the sequence name of fasta file can not exceed 15 character.

The name of first fasta file is 14 character long:

>AAAAAAAAAAAAAA
DIELTQSPTTLSVTPGDRVSLSCRASQSISDYLHWYQQQSHESPRLLIKYASQSISGIPSRFSGSGSGSDFTLSINSVEPEDVGVYYCQNGHSFPLTFGAGTKLELKGGGGSGGGGSGGGGSQVQLQQSGAELVQPGASVKLSCKASGYTFTNYDINWVRQRPEQGLEWIGWIFPGDGSTQYNEKFQGKATLTTDTSSSTAYMQLSRLTSEDSAVYFCARQTTATWFAYWGQGTTVTVSS

Code can run normal:

print("####read fasta file#####")
gen = hm.alignment.AlignmentGenerator_hhblits.from_fasta("fasta/B7-H3_8H9_10_kq.fa",template_location="templates/B7-H3_8H9_10_kq/")
print("#####search template######")
gen.get_suggestion(database_dir="hhsuite_dbs/", n_threads=3)

gen.select_templates(list(temp.template)[0:10])
gen.alignment.sequences["AAAAAAAAAAAAAA"].sequence

'DIELTQSPTTLSVTPGDRVSLSCRASQ-----S-----I------S----DYLHWYQQQSHESPRLLIKYASQSISGIPSRFSGSGSGSDFTLSINSVEPEDVGVYYCQNGHSFPLTFGAGTKLELKGGGGSGGGG---------SGGGGSQVQLQQSGAELVQPGASVKLSCKASGYTFTNYDINWVRQRPEQGLEWIGWIFPGDGSTQYNEKFQGKATLTTDTSSSTAYMQLSRLTSEDSAVYFCARQTTAT-----WFAYWGQGTTVTVSS'

gen.alignment.print_clustal(70)

AAAAAAAAAA  DIELTQSPTTLSVTPGDRVSLSCRASQ-----S-----I------S----DYLHWYQQQSHESPRLLIKY
6EJG_C      -IVLTQSPASLSVSLGQRATISCRASK-----S-----V------STSIYSYMHWYQQKPGQPPKLLIKY
3WBD_B      -VVMTQTPLSLPVSLGDQASISCRSSQ-----SLVHSNG------N----TYLYWYLQKPGQSPKPLIYR
5F3J_C      ---ITQDELSNPVTSGESVSISCRSSK-----SLLYQDG------K----TYLNWFLQRPGQSPQLLIYL
1H8O_B      --VLTQSHKFMSTSVGDRVSITCKASQ-----D-----V------G----TAVAWYQQKPGQSPKLLIYW
1H8N_A      --VLTQSHKFMSTSVGDRVSITCKASQ-----D-----V------G----TAVAWYQQKPGQSPKLLIYW
5KVE_L      -IVMSQSPSSLAVSVGEKITMSCKSSQ-----S-----LLYSNNEK----NYLAWYQQKPGQSPKLLIYW
3ESU_F      --VLIQSTSSLSASLGDRVTISCRASQ-----D-----I------R----NYLNWYQQKPDGTVKLLIYY
1KTR_L      -ILMTQTPSSLPVSLGDQASISCRSSQSIVHSN-----G------N----TYLEWYLQKPGQSPKLLIYK
3ESV_G      --QMTQTTSSLSASLGDRVTVSCRASQ-----D-----I------R----NYLNWYQQKPDGTVKFLIYY
4F9L_C      -IVLTQSPSSLSASLGDTITITCHASQ-----N-----I------N----LWLSWYQQRPGNIPKLLIYR

AAAAAAAAAA  ASQSISGIPSRFSGSGSGSDFTLSINSVEPEDVGVYYCQNGHSFPLTFGAGTKLELKGGGGSGGGG----
6EJG_C      ASYLESGVPARFSGSGSGTDFTLNIHPVEEEDAATYYCEHSREFPFTFGTGTKLEIKGGGGSGGGGSGGG
3WBD_B      VSNRFSGVPDRFSGSGSGTDFTLKISRVEAEDLGVYFCFQGTHVPYTFGGGTRLEIKGGGGSGGGG----
5F3J_C      MSTRASGVSDRFSGSGSGTDFTLEISRVEAEDVGVYYCQQLVEYPLTFGAGTKLELKRADGGGGSGGGG-
1H8O_B      ASTRHTGVPDRFTGSGSGTDFTLTISNVQSEDLADYFCQQYSSYPLTFGAGTKLELKRGGGGSGGGGSGG
1H8N_A      ASTRHTGVPDRFTGSGSGTDFTLTISNVQSEDLADYFCQQYSSYPLTFGAGTKLELKRGGGGSGGGGSGG
5KVE_L      ASARDSGVPDRFTGSGSGTDFTLTISSVKAEDLAVFYCQQYYSYPYTFGGGTKLEIKGGGGSGGGG----
3ESU_F      TSRLQSGVPSRFSGSGSGTDYSLTISNLEQEDIGTYFCQQGNTLPWTFGGGTKLEIRRGGGGSGGGGSGG
1KTR_L      VSNRFSGVPDRFSGSGSGTDFTLKISRVEAEDLGVYYCFQGSHVPFTFGSGTKLEIKRGGGGSGGGGSGG
3ESV_G      TSRLQPGVPSRFSGSGSGTDYSLTINNLEQEDIGTYFCQQGNTPPWTFGGGTKLEIKRGGGGSGGGGSGG
4F9L_C      ASNLHTGVPSRFSGSGSATGFTLTISSLQPEDIATYYCQQGHSYPYTFGGGTKLDIKRADAGGGGSGGGG

AAAAAAAAAA  -----SGGGGSQVQLQQSGAELVQPGASVKLSCKASGYTFTNYDINWVRQRPEQGLEWIGWIFPGDGSTQ
6EJG_C      G----SGGGGSQVQLQQSGPELVKPGASVKISCKASGYTFSSSWMNWVKQRPGKGLEWIGRIYSGDGDAI
3WBD_B      -----SGGGGSQIQLQQSGPELVRPGASVKISCKASGYTFTDYYIHWVKQRPGEGLEWIGWIYPGSGNTK
5F3J_C      -----SGGGGSQVQLQQSGPELVKPGASVKISCKASGYAFISSWMNWVKQRPGKGLEWIGRIYPGDGDTH
1H8O_B      GG---SGGGGSQVQLQEPGGELVRPGASVKLSCKASGYTFTSYWINWVKQRPGQGLEWIGNIYPSDSYTN
1H8N_A      GG---SGGGGSQVQLQESGGELVRPGASVKLSCKASGYTFTSYWINWVKQRPGQGLEWIGNIYPSDSYTN
5KVE_L      -----SGGGGSQVQLQQPGAELLKPGASVKLSCKASGYSFSNYWMHWVKQRPGQGPEWIGMIHPNSGNTK
3ESU_F      GG---SGGGGSEVQLQQSGPELVKPGASVKISCKDSGYAFSSSWMNWVKQRPGQGPEWIGRIYPGDGDTN
1KTR_L      GG---SGGGGSQVQLQQSGPEDVKPGASVKISCKASGYTFTDYYMNWVKQSPGKGLEWIGDINPNNGGTS
3ESV_G      GG---SGGGGSEVQLQQSGPELVKPGASVKISCKDSGYAFNSSWMNWVKQRPGQGLEWIGRIYPGDGDSN
4F9L_C      SGGGGSGGGGSQVQLQESGAELVKPGASVKLSCKASGYTFTRYYLYWVKQRPGQGLEWIGEINPNNGGTK

AAAAAAAAAA  YNEKFQGKATLTTDTSSSTAYMQLSRLTSEDSAVYFCARQTTAT-----WFAYWGQGTTVTVSS
6EJG_C      YNGKFKGKATLTADKSSSTAYMQLSSLTSEDSAVYFCAREG-KTGDL--LLRSWGQGSALTVSS
3WBD_B      YNEKFKGKATLTVDTSSSTAYMQLSSLTSEDSAVYFCARGG-KF-----AMDYWGQGTSVTVSS
5F3J_C      YNGKFKGKATLTADKSSSTAYMQLSSLTSEDSAVYFCAREE-TAQTG--GFDYWGQGTTLTVS-
1H8O_B      YNQKFKDKATLTVDKSSSTAYMQLSSLTSEDSAVYFCARWG-Y----------WGQGTLVTVS-
1H8N_A      YNQKFKDKATLTVDKSSSTAYMQLSSLTSEDSAVYFCARWG-Y----------WGQGTLVTV--
5KVE_L      YNEKFKNKATLTVDKSSSMVYMQLSSLTSEDSAVFYCARLG-ND------MDYWGQGTSVTVSS
3ESU_F      YNGKFKGKATLTADKSSSTAYMQLSSLTSVDSAVYFCARSG-LLRY---AMDYWGQGTSVTVSS
1KTR_L      YNQKFKGRATLTVDKSSSTAYMELRSLTSEDSSVYYCESQS-G--------AYWGQGTTVTVSA
3ESV_G      YNGKFEGKAILTADKSSSTAYMQLSSLTSVDSAVYFCARSG-LLRY---AMDYWGQGTSVTVSS
4F9L_C      FNEKFKSKATLTVDKSSRTTYIQLSSLTSEDSAVYYCSRED-DYDGTPDAMDYWGQGTAVTVS-

But when I change the FASTA sequence name to 15 character length:

>AAAAAAAAAAAAAAA
DIELTQSPTTLSVTPGDRVSLSCRASQSISDYLHWYQQQSHESPRLLIKYASQSISGIPSRFSGSGSGSDFTLSINSVEPEDVGVYYCQNGHSFPLTFGAGTKLELKGGGGSGGGGSGGGGSQVQLQQSGAELVQPGASVKLSCKASGYTFTNYDINWVRQRPEQGLEWIGWIFPGDGSTQYNEKFQGKATLTTDTSSSTAYMQLSRLTSEDSAVYFCARQTTATWFAYWGQGTTVTVSS

Run the same code, this time the cluster is empty and all templates identity is 0:

print("####read fasta file#####")
gen = hm.alignment.AlignmentGenerator_hhblits.from_fasta("fasta/B7-H3_8H9_10_kq.fa",template_location="templates/B7-H3_8H9_10_kq/")
print("#####search template######")
gen.get_suggestion(database_dir="hhsuite_dbs/", n_threads=3)
gen.select_templates(list(temp.template)[0:10])
gen.alignment.sequences["AAAAAAAAAAAAAAA"].sequence

'DIELTQSPTTLSVTPGDRVSLSCRASQSISDYLHWYQQQSHESPRLLIKYASQSISGIPSRFSGSGSGSDFTLSINSVEPEDVGVYYCQNGHSFPLTFGAGTKLELKGGGGSGGGGSGGGGSQVQLQQSGAELVQPGASVKLSCKASGYTFTNYDINWVRQRPEQGLEWIGWIFPGDGSTQYNEKFQGKATLTTDTSSSTAYMQLSRLTSEDSAVYFCARQTTATWFAYWGQGTTVTVS'

gen.alignment.print_clustal(70)

AAAAAAAAAA  DIELTQSPTTLSVTPGDRVSLSCRASQSISDYLHWYQQQSHESPRLLIKYASQSISGIPSRFSGSGSGSD
6EJG_C      ----------------------------------------------------------------------
3WBD_B      ----------------------------------------------------------------------
5F3J_C      ----------------------------------------------------------------------
1H8O_B      ----------------------------------------------------------------------
1H8N_A      ----------------------------------------------------------------------
5KVE_L      ----------------------------------------------------------------------
3ESU_F      ----------------------------------------------------------------------
1KTR_L      ----------------------------------------------------------------------
3ESV_G      ----------------------------------------------------------------------
4F9L_C      ----------------------------------------------------------------------

AAAAAAAAAA  FTLSINSVEPEDVGVYYCQNGHSFPLTFGAGTKLELKGGGGSGGGGSGGGGSQVQLQQSGAELVQPGASV
6EJG_C      ----------------------------------------------------------------------
3WBD_B      ----------------------------------------------------------------------
5F3J_C      ----------------------------------------------------------------------
1H8O_B      ----------------------------------------------------------------------
1H8N_A      ----------------------------------------------------------------------
5KVE_L      ----------------------------------------------------------------------
3ESU_F      ----------------------------------------------------------------------
1KTR_L      ----------------------------------------------------------------------
3ESV_G      ----------------------------------------------------------------------
4F9L_C      ----------------------------------------------------------------------

AAAAAAAAAA  KLSCKASGYTFTNYDINWVRQRPEQGLEWIGWIFPGDGSTQYNEKFQGKATLTTDTSSSTAYMQLSRLTS
6EJG_C      ----------------------------------------------------------------------
3WBD_B      ----------------------------------------------------------------------
5F3J_C      ----------------------------------------------------------------------
1H8O_B      ----------------------------------------------------------------------
1H8N_A      ----------------------------------------------------------------------
5KVE_L      ----------------------------------------------------------------------
3ESU_F      ----------------------------------------------------------------------
1KTR_L      ----------------------------------------------------------------------
3ESV_G      ----------------------------------------------------------------------
4F9L_C      ----------------------------------------------------------------------

AAAAAAAAAA  EDSAVYFCARQTTATWFAYWGQGTTVTVS
6EJG_C      -----------------------------
3WBD_B      -----------------------------
5F3J_C      -----------------------------
1H8O_B      -----------------------------
1H8N_A      -----------------------------
5KVE_L      -----------------------------
3ESU_F      -----------------------------
1KTR_L      -----------------------------
3ESV_G      -----------------------------
4F9L_C      -----------------------------

gen.show_suggestion()
template | coverage | identity -- | -- | -- 6EJG_C | 0.0 | 0.0 3WBD_B | 0.0 | 0.0 5F3J_C | 0.0 | 0.0 1H8O_B | 0.0 | 0.0 1H8N_A | 0.0 | 0.0 5KVE_L | 0.0 | 0.0 3ESU_F | 0.0 | 0.0 1KTR_L | 0.0 | 0.0 3ESV_G | 0.0 | 0.0 4F9L_C | 0.0 | 0.0
PhilippJunk commented 1 year ago

Thanks for reporting this bug, this is definitely not intended behavior. I should have some time on the weekend to fix this.

wt12318 commented 1 year ago

Hi,

Thanks for updating. I used the up-to-date docker, but the code not change:

import homelette as hm
hm.__version__
##'1.4'

import inspect
with open('hm.txt', 'w', encoding='utf-8') as f:
    print(inspect.getsource(hm.alignment), file=f)
f.close()
cat hm.txt  | grep "limit"
## self.target = target[:14]  # limit length because of hhblits

But the commit history said this changed to self.target = target :

image

And the error still exist:

###fasta file:
>tttttttttttttttttttttttttttttttttttttttttt
EVQLVESGGGLVQPGGSLRLSCAASGYTFTDYYIHWVRQAPGKGLEWMAWISPHTGGTIYADSVKGRFTISADTSKNTAYLQMNSLRAEDTAVYYCARGPDDWNDGDAFDIWGQGTLVTVSSGGGSGGSGGGSDIQMTQSPSSLSASVGDRVTITCKSSQSVLYSSNNENFLAWYQQKPGKAPKLLIYWASTRESGVPSRFSGSRSGTDFTLTISSLQPEDFATYYCQQYYSTPITFGQGTKVEIK

##error
Traceback (most recent call last):
  File "./get_patch_V2.py", line 729, in <module>
    t.execute_routine(
  File "/usr/local/src/homelette-1.4/homelette/organization.py", line 256, in execute_routine
    r.generate_models()
  File "/usr/local/src/homelette-1.4/homelette/routines.py", line 944, in generate_models
    self.alignment.rename_sequence(self.target, "Target")
  File "/usr/local/src/homelette-1.4/homelette/alignment.py", line 575, in rename_sequence
    self.sequences[new_name] = self.sequences.pop(old_name)
KeyError: 'tttttttttttttttttttttttttttttttttttttttttt'
PhilippJunk commented 1 year ago

It should have changed from

self.target = target

to

self.target = target[:14]

so the code in the docker container seems to be correct.

See here for the current version.


My documentation of the change has been not great though, I will try to improve that. Basically, the issue was related to how hhblits under the hood limits the name of the sequence to search to 14 characters because of the output format.

Before the fix, when parsing the output file of hhblits, since the name got changed, the alignment ended up empty because the regex of the longer name would be no-where to be found in the file.

Now, the target name is automatically (and silently) limited to 14 characters on the generation of a AlignmentGenerator object. This should guarantee that the parsing of the hhblits output file succeeds every time.


For me, the following code runs flawlessly in the docker container:

import homelette as hm
gen = hm.alignment.AlignmentGenerator_hhblits.from_fasta(
    "fasta_file.fa")
gen.get_suggestion(database_dir="hhsuite_dbs/")
gen.show_suggestion()

# select only a few to make output more comprehensive
gen.select_templates(['6VUN_A', '7CU5_A'])
gen.show_suggestion()
gen.alignment.print_clustal(70)

# pull PDB files
gen.get_pdbs()

# initialize task and perform modelling
task = gen.initialize_task()
task.execute_routine(
    tag = 'example_modeller',
    routine = hm.routines.Routine_automodel_default,
    templates = ['6VUN_A'],
    template_location = './templates/')
task.models

Regarding your actual error, it seems to be generated from inside Task.execute_routine. Without the actual code to recreate this, it is kind of difficult to talk about the error, but my guess would be that you initialized the Task object while still expecting the long name, not the automatically shortened one?

wt12318 commented 1 year ago

Sorry, it's my fault. I still used the original target name, not the shortened one. Thanks a lot.

PhilippJunk commented 1 year ago

No worries, glad it is working now.