heche-psb / wgd

wgd v2: a suite of tools to uncover and date ancient polyploidy and whole-genome duplication
https://wgdv2.readthedocs.io/en/latest/
GNU General Public License v3.0
21 stars 0 forks source link

`wgd viz` throws error when including self-comparison as --speciespair #12

Closed taprs closed 10 months ago

taprs commented 10 months ago

Hi again! I found a few bugs which are not absolutely prohibitive but would be good to fix at some point. I think I will post them here one by one.

I run substitution rate correction with wgd viz in wgd 2.0.22 as follows:

wgd viz -d wgd_ksd/global_MRBH.tsv.ks.tsv \
    -epk species1_ksd/species1.fa.tsv.ks.tsv -sp speciestree.txt -rw -ap species1_syn/iadhore-out/anchorpoints.txt  \
    -sr "species1.fa;species2.fa" -sr "species1.fa;species3.fa" -sr "species1.fa;species4.fa" -sr "species1.fa;species1.fa" \
    -gs wgd_ksd/gene_species.map --plotkde --plotelmm

and get the following StopIteration output:


14:15:25 INFO     This is wgd v2.0.22                                                                                                               cli.py:32
14:15:28 INFO     Implementing node-averaged Ks analysis                                                                                           viz.py:511
Traceback (most recent call last):
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd_2.0.22/bin/wgd", line 11, in 
    load_entry_point('wgd', 'console_scripts', 'wgd')()
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd_2.0.22/lib/python3.8/site-packages/click/core.py", line 829, in __call__
    return self.main(*args, **kwargs)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd_2.0.22/lib/python3.8/site-packages/click/core.py", line 782, in main
    rv = self.invoke(ctx)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd_2.0.22/lib/python3.8/site-packages/click/core.py", line 1259, in invoke
    return _process_result(sub_ctx.command.invoke(sub_ctx))
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd_2.0.22/lib/python3.8/site-packages/click/core.py", line 1066, in invoke
    return ctx.invoke(self.callback, **ctx.params)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd_2.0.22/lib/python3.8/site-packages/click/core.py", line 610, in invoke
    return callback(*args, **kwargs)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/cli.py", line 557, in viz
    _viz(**kwargs)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/cli.py", line 597, in _viz
    multi_sp_plot(df,spair,gsmap,outdir,onlyrootout,title=prefix,ylabel=ylabel,viz=True,plotkde=plotkde,reweight=False,sptree=speciestree,ap = anchorpoints, extraparanomeks=extraparanomeks,plotapgmm=plotapgmm,plotelmm=plotelmm,components=components,max_EM_iterations=em_iterations,num_EM_initializations=em_initializations,peak_threshold=prominence_cutoff,rel_height=rel_height, na=True,user_xlim=xlim,user_ylim=ylim)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd/viz.py", line 527, in multi_sp_plot
    df_perspair,allspair,paralog_pair,corrected_ks_spair,Outgroup_spnames = getspair_ks(spair,df,reweight,onlyrootout,sptree=sptree,na=na,spgenemap=spgenemap)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd/viz.py", line 72, in getspair_ks
    if sptree != None and len(paralog_pair) !=0 : corrected_ks_spair,Outgroup_spnames = correctks(df,sptree,paralog_pair[0],reweight,onlyrootout,na=na)
  File "/netscratch/dep_mercier/grp_novikova/software/wgd/wgd/viz.py", line 180, in correctks
    focussp_clade = next(tree.find_clades({"name": focusp}))
StopIteration

Removing the -sr "species1.fa;species1.fa" portion of the command results in a successful run. Can you reproduce this?

heche-psb commented 10 months ago

Hi, could you show me the content of your speciestree.txt?

taprs commented 10 months ago

Here it is:

(zostera:0.261188,((polygonifolius:0.0158094,coloratus:0.0182822)0.806137:0.0266174,(pusillus:0.0230147,(obtusifolius:0.0181003,(trichoides:0.0109874,acutifolius:0.013523)0.386264:0.0102274)0.345251:0.0108846)0.679006:0.0221993)1:0.261188);

I should have thought about matching the tree labels with the input data! Maybe it would be good to have a check for that for people as absent-minded as I. Also I noticed the file did not have a newline in the end.

heche-psb commented 10 months ago

Yes, the label in your speciestree.txt should match the parameters given by -sr option.

taprs commented 10 months ago

Inputting the correct species tree solved the issue for me.