Answers.Body: <p>So, as the first reply said, there's sort of multiple questions here.</p>
  
  <p>For filtering things out from a table, I would absolutely use pandas. <a href="https://swcarpentry.github.io/python-novice-gapminder/08-data-frames/" rel="nofollow noreferrer">https://swcarpentry.github.io/python-novice-gapminder/08-data-frames/</a> Pandas is a way of managing dataframes in python. It looks like what you want from your table is to remove certain rows, if the content of those rows at a certain column matches certain IDs from a list (pathogen IDs). That's the sort of thing the pandas library is good at.</p>
  
  <p>For managing workflows, I use either bash or snakemake. Snakemake is a python thing. It has some advantages, such as being able to start from an intermediary file (for example, if a workflow is "sequences -> multiple alignment -> phylogenetic tree", it can easily start from the multiple alignment). Also deals with branched workflows quite well. Can be easily installed with conda (anaconda is absolutely a lifesaver for bioinformatics work imo). I still use bash for simple workflows, because it's slightly simpler/faster.</p>
  
  <p>For potential filtering - If you want to keep the top hits as they are, but filter the lower-scoring hits based on whether they have a PDB ID or not, I'd personally look into Entrez queries. They're a way of making your search more specific. Here they're shown for WebBlast, but they work quite well for command line blast too. <a href="https://www.youtube.com/watch?v=bxx5uaKjMa8" rel="nofollow noreferrer">https://www.youtube.com/watch?v=bxx5uaKjMa8</a></p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-04-27T07:52:33.677Z
Answers.Text: @Mike I'm afraid that I haven't used EMBL-EBI much, so I'm not sure how to help. Entrez blast queries are mostly an NCBI thing. Still, the 'sequences associated with solved protein structures' should be the same in both cases, because they represent links to the PDB protein database. So I guess there are different options, from re-doing the search in NCBI, to downloading the sequences of interest from EMBL-EBI and then blasting in the NCBI with Entrez (either web-based or command line), to trying to get the NCBI code for EMBL sequences...
Answers.UserDisplayName: Laura
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-05-26T21:18:37.503Z
Answers.Text: Thank you! I swear I replied to this, so sorry Laura
Answers.UserDisplayName: Mike 
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-04-25T12:16:19.680Z
Answers.Text: Thanks so much for this; my issue is, or my question is, will entrez have the CDS-based proteome that I'm using? I have a proteome from deposited into EMBL-EBI
Answers.UserDisplayName: Mike 
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-04-24T07:17:35.613Z
Answers.LastActivityDate: 2020-04-24T07:17:35.613Z
Answers.UserDisplayName: Laura
Answers.Body: <p>The question is still in numerous parts,</p>
  
  <ul>
  <li>parsing a blast output</li>
  <li>obtaining sequences </li>
  <li>feeding this as input into clustalo</li>
  </ul>
  
  <p>If I was automating this I would parse the Blast genbank codes [regex either Perl one-liner (scripts aren't trendy ;-( ), or Pythons "re"] and feed them into NCBI's <code>efetch</code> to hook out the sequence data. Then I would create a subprocess that feed the files into clustalo. In Perl you simply issue the linux command wrapped in "`" either side, in Python use </p>
  
  <pre><code>import os
  os.system('clustalo -in file.in -out file.out')
  </code></pre>
  
  <p>Ensuing your 'bin' is correctly configured to execute clustalo</p>
  
  <p>or </p>
  
  <pre><code>import subprocess # use popen
  </code></pre>
  
  <p>To my knowledge there is no wrapper for <code>clustalo</code> unlike several other phylogenetics programs in BioPython, which enormously simplifies pipelines. Others might recommend using shell rather than e.g. Python.</p>
  
  <p>BTW I personally use <code>Muscle</code> having been a long term user of clustalo, thats a long discussion. </p>
  
  <hr>
  
  <p>To be honest you should look at the power of NBCI's online blast utility. The download facilities are certainly powerful particularly if you are pulling down e.g. 1500 sequences from blast hits, because it automatically can splice the hits of interest from their genes/genomes. However, its context specific.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-03-26T08:30:01.600Z
Answers.Text: To be clear, I already have the sequence data for these uniprot IDs on my desktop in another file and they are linked to their sequences. I wish to do pairwise alignments on sequences detected by jackhmmer
Answers.UserDisplayName: Mike 
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-03-26T09:38:10.407Z
Answers.Text: Hmmmmm  ... anyway subprocessing is the specific answer you requested on how to automate clustalo ... unless you are using the bash shell
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-03-26T10:34:09.517Z
Answers.Text: Currently doing some bipython tutorials at the moment, but yeah I do a lot through the command line
Answers.UserDisplayName: Mike 
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-03-26T13:15:26.203Z
Answers.Text: Biopython does include a wrapper for ClustalO (and MAFFT, MUSCLE etc.)
Answers.UserDisplayName: Chris_Rands
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-03-25T23:59:47.137Z
Answers.LastActivityDate: 2020-03-26T00:05:44.343Z
Answers.LastEditDate: 2020-03-26T00:05:44.343Z
Answers.UserDisplayName: M__
Body: <p>I've detected homology between targets of ligands in drugbank and proteins in the proteome of a pathogen. I've parsed the output very rudimentary and calculated my query coverage. This exists in an excel file (see attached)</p>
  
  <p>I'd like to further filter these hits. For query coverages > 70% alone I have 1500+ hits. Some redundancy exists and some hits are for drugbank targets for which no data has been deposited in the PDB. </p>
  
  <p>So there are ways of further filtering these hits, but I'd like to do pairwise alignments first before I do anything else, as I think this is the most rational filter I can think of. </p>
  
  <p>I do not want to filter based on lack of PDB data necessarily, as when I check my high query coverage hits, I'm getting predictions that are actually backed up by in vitro data. This includes high coverage hits lacking x-ray data</p>
  
  <p>My question:</p>
  
  <p>I have CLUSTALO installed on my desktop: How do programmatically send it jobs if I have the Drugbank sequences and the pathogen sequences in the same file? I would also like to parse out the pathogen IDs from the excel file. There’s 1500 and they are a bit scattered throughout the columns. </p>
  
  <pre><code>0.980634         -                    572         sp|P06672|PDC_ZYMMO                         -                     568                               1.2e-263             874.5          0.0            1        1      6.6e-267  1.4e-263   874.2      0.0       5         562    15     570    12     572    0.98   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_005678  PE=3             SV=1
  
  0.950704         -                    580         sp|P06672|PDC_ZYMMO                         -                     568                               1.2e-258             858.0          0.0            1        1      6.8e-262  1.5e-258   857.7      0.0       6         546    10     560    7      578    0.98   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_005014  PE=3             SV=1
  
  0.950704         -                    691         sp|P06672|PDC_ZYMMO                         -                     568                               8.2e-162             538.3          0.1            1        1      2.4e-164  5.1e-161   535.6      0.1       1         541    89     650    89     658    0.89   Acetolactate                            synthase                  OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_006830  PE=3             SV=1
  
  0.492958         -                    585         sp|P06672|PDC_ZYMMO                         -                     568                               8e-87                290.7          0.0            1        2      3.3e-83   7e-80      267.8      0.0       3         283    2      287    1      292    0.94   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_001897  PE=3             SV=1
  
  0.107394         -                    585         sp|P06672|PDC_ZYMMO                         -                     568                               8e-87                290.7          0.0            2        2      1.7e-08   3.6e-05    21.1       0.0       413       474    436    505    376    583    0.72   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_001897  PE=3             SV=1
  </code></pre>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T10:04:58.227Z
Comments.Text: you will need to structure your questions so it is one question per post by breaking-down your goals into components. At present your question is not answerable I'm afraid, the site rules are designed to enable search engines to easily identify the post.
Comments.UserDisplayName: M__
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T10:24:33.677Z
Comments.Text: Is this an automated comment? I shall parse this out into multiple questions (pun intended)
Comments.UserDisplayName: Mike 
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T11:31:23.500Z
Comments.Text: Thanks @Mike (thats not automated "thanks")
Comments.UserDisplayName: M__
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T13:04:08.733Z
Comments.Text: Also, make sure to tell us what operating system you are using, and why excel (!?) is involved. Finally, explain if you need every possible combination of pairwise alignments from your sequence file or only some.
Comments.UserDisplayName: terdon
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T13:17:35.090Z
Comments.Text: Thank you, Michael! Terdon, I just wanted to view the output in a consolidated table so I could send to others for review, though of course the awk one-liners and GREP commands were done on the command line (gnome emulator); I'm using the latest version of ubuntu
Comments.UserDisplayName: Mike 
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T13:19:15.677Z
Comments.Text: I'm using CLUSTALO mainly for single alignments; I have a few MSAs I need to run, but not many. Sorry for the confusion. I'm using CLUSTALO for the single alignments we've agreed it's gotten us optimal alignments
Comments.UserDisplayName: Mike 
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T15:39:46.803Z
Comments.Text: If you can give us a few entries from your input files and explain what you want to do with them, we can try and answer. But I don't understand why you would be using an excel file. Can you export to a simple text-based format? What do you have in the excel? We can't really help you process data if you don't show us the data you need to process.
Comments.UserDisplayName: terdon
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T19:17:01.533Z
Comments.Text: Sure thing, see an edit to the original post. The original output is in a .txt, I just copied and pasted it into the excel sheet to send to a colleague.
Comments.UserDisplayName: Mike 
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T20:02:38.553Z
Comments.Text: Peronsally although I automate alot and automate clustalo I don't automate downloads because Blast online is pretty slick and has very neat tools, but its context specific.
Comments.UserDisplayName: M__
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-03-25T20:19:50.430Z
Comments.Text: Sure understood, automating clustalo, which I have locally, will do me well right now. I would eventually like to look at automating what looks to be a simple pipeline
Comments.UserDisplayName: Mike 
Tags: <ngs><python><protein-structure><pdb><protein-protein-interaction>
Title: Can I automate CLUSTALO and output alignment sequence identity?
ThreadText: # Can I automate CLUSTALO and output alignment sequence identity?
  
  (Asked by: username_0 on 2020-03-25)
  
  I've detected homology between targets of ligands in drugbank and proteins in the proteome of a pathogen. I've parsed the output very rudimentary and calculated my query coverage. This exists in an excel file (see attached)
  
  I'd like to further filter these hits. For query coverages > 70% alone I have 1500+ hits. Some redundancy exists and some hits are for drugbank targets for which no data has been deposited in the PDB.
  
  So there are ways of further filtering these hits, but I'd like to do pairwise alignments first before I do anything else, as I think this is the most rational filter I can think of.
  
  I do not want to filter based on lack of PDB data necessarily, as when I check my high query coverage hits, I'm getting predictions that are actually backed up by in vitro data. This includes high coverage hits lacking x-ray data
  
  My question:
  
  I have CLUSTALO installed on my desktop: How do programmatically send it jobs if I have the Drugbank sequences and the pathogen sequences in the same file? I would also like to parse out the pathogen IDs from the excel file. There’s 1500 and they are a bit scattered throughout the columns.
  
  ```
  0.980634         -                    572         sp|P06672|PDC_ZYMMO                         -                     568                               1.2e-263             874.5          0.0            1        1      6.6e-267  1.4e-263   874.2      0.0       5         562    15     570    12     572    0.98   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_005678  PE=3             SV=1
  
  0.950704         -                    580         sp|P06672|PDC_ZYMMO                         -                     568                               1.2e-258             858.0          0.0            1        1      6.8e-262  1.5e-258   857.7      0.0       6         546    10     560    7      578    0.98   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_005014  PE=3             SV=1
  
  0.950704         -                    691         sp|P06672|PDC_ZYMMO                         -                     568                               8.2e-162             538.3          0.1            1        1      2.4e-164  5.1e-161   535.6      0.1       1         541    89     650    89     658    0.89   Acetolactate                            synthase                  OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_006830  PE=3             SV=1
  
  0.492958         -                    585         sp|P06672|PDC_ZYMMO                         -                     568                               8e-87                290.7          0.0            1        2      3.3e-83   7e-80      267.8      0.0       3         283    2      287    1      292    0.94   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_001897  PE=3             SV=1
  
  0.107394         -                    585         sp|P06672|PDC_ZYMMO                         -                     568                               8e-87                290.7          0.0            2        2      1.7e-08   3.6e-05    21.1       0.0       413       474    436    505    376    583    0.72   Uncharacterized                         protein                   OS=Lomentospora            prolificans        OX=41688             GN=jhhlp_001897  PE=3             SV=1
  
  ```
  
  **Tags:** ngs, python, protein-structure, pdb, protein-protein-interaction
  
  > **Comments:**
  > **username_1**: you will need to structure your questions so it is one question per post by breaking-down your goals into components. At present your question is not answerable I'm afraid, the site rules are designed to enable search engines to easily identify the post.  
  > **username_0**: Is this an automated comment? I shall parse this out into multiple questions (pun intended)  
  > **username_1**: Thanks @username_0 (thats not automated "thanks")  
  > **username_2**: Also, make sure to tell us what operating system you are using, and why excel (!?) is involved. Finally, explain if you need every possible combination of pairwise alignments from your sequence file or only some.  
  > **username_0**: Thank you, Michael! username_2, I just wanted to view the output in a consolidated table so I could send to others for review, though of course the awk one-liners and GREP commands were done on the command line (gnome emulator); I'm using the latest version of ubuntu  
  > **username_0**: I'm using CLUSTALO mainly for single alignments; I have a few MSAs I need to run, but not many. Sorry for the confusion. I'm using CLUSTALO for the single alignments we've agreed it's gotten us optimal alignments  
  > **username_2**: If you can give us a few entries from your input files and explain what you want to do with them, we can try and answer. But I don't understand why you would be using an excel file. Can you export to a simple text-based format? What do you have in the excel? We can't really help you process data if you don't show us the data you need to process.  
  > **username_0**: Sure thing, see an edit to the original post. The original output is in a .txt, I just copied and pasted it into the excel sheet to send to a colleague.  
  > **username_1**: Peronsally although I automate alot and automate clustalo I don't automate downloads because Blast online is pretty slick and has very neat tools, but its context specific.  
  > **username_0**: Sure understood, automating clustalo, which I have locally, will do me well right now. I would eventually like to look at automating what looks to be a simple pipeline  
  
  ---
  
  ## 2 answers
  
  ### Answer by username_3 on 2020-04-24. Score: 1
  
  So, as the first reply said, there's sort of multiple questions here.
  
  For filtering things out from a table, I would absolutely use pandas. <https://swcarpentry.github.io/python-novice-gapminder/08-data-frames/> Pandas is a way of managing dataframes in python. It looks like what you want from your table is to remove certain rows, if the content of those rows at a certain column matches certain IDs from a list (pathogen IDs). That's the sort of thing the pandas library is good at.
  
  For managing workflows, I use either bash or snakemake. Snakemake is a python thing. It has some advantages, such as being able to start from an intermediary file (for example, if a workflow is "sequences -> multiple alignment -> phylogenetic tree", it can easily start from the multiple alignment). Also deals with branched workflows quite well. Can be easily installed with conda (anaconda is absolutely a lifesaver for bioinformatics work imo). I still use bash for simple workflows, because it's slightly simpler/faster.
  
  For potential filtering - If you want to keep the top hits as they are, but filter the lower-scoring hits based on whether they have a PDB ID or not, I'd personally look into Entrez queries. They're a way of making your search more specific. Here they're shown for WebBlast, but they work quite well for command line blast too. <https://www.youtube.com/watch?v=bxx5uaKjMa8>
  
  > **Comments:**
  > **username_0**: Thanks so much for this; my issue is, or my question is, will entrez have the CDS-based proteome that I'm using? I have a proteome from deposited into EMBL-EBI  
  > **username_3**: @username_0I'm afraid that I haven't used EMBL-EBI much, so I'm not sure how to help. Entrez blast queries are mostly an NCBI thing. Still, the 'sequences associated with solved protein structures' should be the same in both cases, because they represent links to the PDB protein database. So I guess there are different options, from re-doing the search in NCBI, to downloading the sequences of interest from EMBL-EBI and then blasting in the NCBI with Entrez (either web-based or command line), to trying to get the NCBI code for EMBL sequences...  
  > **username_0**: Thank you! I swear I replied to this, so sorry username_3  
  
  ---
  
  ### Answer by username_1 on 2020-03-25. Score: 0
  
  The question is still in numerous parts,
  
  * parsing a blast output
  * obtaining sequences
  * feeding this as input into clustalo
  
  If I was automating this I would parse the Blast genbank codes [regex either Perl one-liner (scripts aren't trendy ;-( ), or Pythons "re"] and feed them into NCBI's `efetch` to hook out the sequence data. Then I would create a subprocess that feed the files into clustalo. In Perl you simply issue the linux command wrapped in "`" either side, in Python use
  
  ```
  import os
  os.system('clustalo -in file.in -out file.out')
  
  ```
  
  Ensuing your 'bin' is correctly configured to execute clustalo
  
  or
  
  ```
  import subprocess # use popen
  
  ```
  
  To my knowledge there is no wrapper for `clustalo` unlike several other phylogenetics programs in BioPython, which enormously simplifies pipelines. Others might recommend using shell rather than e.g. Python.
  
  BTW I personally use `Muscle` having been a long term user of clustalo, thats a long discussion.
  
  ---
  
  To be honest you should look at the power of NBCI's online blast utility. The download facilities are certainly powerful particularly if you are pulling down e.g. 1500 sequences from blast hits, because it automatically can splice the hits of interest from their genes/genomes. However, its context specific.
  
  > **Comments:**
  > **username_0**: To be clear, I already have the sequence data for these uniprot IDs on my desktop in another file and they are linked to their sequences. I wish to do pairwise alignments on sequences detected by jackhmmer  
  > **username_1**: Hmmmmm  ... anyway subprocessing is the specific answer you requested on how to automate clustalo ... unless you are using the bash shell  
  > **username_0**: Currently doing some bipython tutorials at the moment, but yeah I do a lot through the command line  
  > **username_4**: Biopython does include a wrapper for ClustalO (and MAFFT, MUSCLE etc.)  
  
  ---
  
  