diff --git a/README.md b/README.md index 44a0eed..40b085a 100644 --- a/README.md +++ b/README.md @@ -1,30 +1,36 @@ # TNT Script used by Guoyi -These scripts follow MIT except for setk.run which belongs to Salvador Arias (Instituto Miguel Lillo, San Miguel de Tucumán, Argentina). +These scripts follow MIT, part of script is from setk.run belonging to Salvador Arias (Instituto Miguel Lillo, San Miguel de Tucumán, Argentina). ## Usage -### Method1 +- Place this script and your tnt matrix under the folder that you call `tnt` or `exe` file is placed -- OS: Arch Linux; +- Enter `tnt` -- `tnt run guoyi.run file.tnt`, `file.tnt` is your data; +- Enter command `guoyi filename;` -- Type `;` and enter; +## Functions -### Method2 +- Estimated implied weighting K value. -- Replace `%1` with your filename; +- Search trees via TBR Mult and Xmult. -- Run `tnt p guoyi.run` directly; +- Perform Strict consensus. + +- Calculate Relative Bremer support, jackknifing, and bootstrap. + +- Calculate TL, CI, and RI. ## Options -- `trees.tre`, `bt.tre` are trees with taxaname; +- Results instructions are at the end of `tnt.log`. -- `trees_no.tre`, `bt_no.tre` are trees without taxaname which can be put into Winclada with the `file.tnt`; +- `con.tre`, `con.tre` are trees with taxaname. -- `majority` can be replaced by `nelsen`; +- `trees_no.tre`, `bt_no.tre` are trees without taxaname which can be put into Winclada with the `file.tnt`. -- `boot` can be replaced by `jack`; +- `nelsen` can be replaced by `majority`. + +- `xmult` and `mult` replications and hold trees number can be adjusted. diff --git a/guoyi.run b/guoyi.run index 9f2fe68..8a744ba 100644 --- a/guoyi.run +++ b/guoyi.run @@ -5,29 +5,174 @@ mxram 10240; nstates 32; nstates NOGAPS; -procedure %1; -setk.run; +/*Arguments*/ +if ( argnumber != 1 ) + silent -console; + quote - /-----------------------------------------------\; + quote - | GUOYI TNT SCRIPT 2022-2023 MIT |; + quote - | You need to give your filename |; + quote - | tnt*> guoyi filename; |; + quote - \-----------------------------------------------/; + proc/; +end +/*Report what will be done*/ +quote - /-----------------------------------------------\; +quote - | Implied weighting will be estimated. |; +quote - | TBR Mult and Xmult will be performed. |; +quote - | Strict consensus will be used. |; +quote - | Relative bremer support, jackknifing and |; +quote - | bootstrap will be shown on the svg. |; +quote - | TL, CI and RI will be calculated finally. |; +quote - \------------------------------------------------/; + + +/*Set K*/ procedure %1; -hold 1000; +var: + actK + minK + maxK + gVal + a + b + Nval + maxIts +; +if (ntax<0) + quote NO HAY DATOS!; + proc/; +end +set gVal 1; +/** busca el valor de g **/ +loop 0 nchar + if (!isinfo[#1]) continue; end + if (!isact[#1]) continue; end + set a=maxsteps[#1]; + set b=minsteps[#1]; + if (('a'-'b')>'gVal') + set gVal 'a'-'b'; + end +stop +quote valor de g = 'gVal'; +set actK 10; +set minK 0; +set maxK 500; /** Asume 500 como el maximo posible valor de K **/ +set maxIts 0; +/** busca el mejor valor de K **/ +loop 0 1 + quote actual valor de k = 'actK'; + set a 1-('actK'/('actK'+1)); + quote 'a'; + set b ('actK'/('actK'+'gVal'-1))-('actK'/('actK'+'gVal')); + quote 'b'; + set Nval 'a'/'b'; + quote valor de N = 'Nval'; + if (('Nval'>14.8)&&('Nval'<15.2)) endloop; end /** N esta en el rango **/ + if ('Nval'>15.2) + set minK 'actK'; + else + set maxK 'actK'; + end + set a ('maxK'-'minK')/2; + set actK 'minK' + 'a'; + set maxIts ++; + if ('maxIts'==100) endloop; end /** salida de emergencia **/ + setloop 0; +stop +/** Si salio de emergencia **/ +if ('maxIts'==100) + quote NO SE TERMINARON LAS ITERACIONES; + quote mejor K encontrado: 'actK'; + proc /; +end +quote Valor de K: 'actK' (N='Nval'); +piwe='actK'; + +/*Reopen tnt*/ +procedure %1; +hold 5000; + +/*Search trees*/ xpiwe(*; -mult=replic 1000 tbr hold 1; +mult=replic 1000 tbr hold 10; +sect: slack 1000; +xmult=rss hit 1000 drift 10 ratchet 10 fuse 10; piwe&; + +/*Export trees*/ export= trees.tre; taxname-; export= trees_no.tre; taxname=; +/*Consensus tree*/ bbreak=tbr; -majority *; +nelsen *; +/*Get Bremer/jak/boot support*/ ttags=; +sub 1; hold +1000; bbreak=fill; +sub 2; hold +2000; bbreak=fill; +sub 3; hold +3000; bbreak=fill; +sub 4; hold +4000; bbreak=fill; +sub 5; hold +5000; bbreak=fill; +sub 6; hold +6000; bbreak=fill; +sub 7; hold +7000; bbreak=fill; +sub 8; hold +8000; bbreak=fill; +bsupport[; +resample jak replications 1000; resample boot replications 1000; -ttags & bt.svg thickness 7 italics fontsize 15; -export < bt.tre; -taxname-; -export - bt_no.tre; +/*Export consensus tree with supports*/ +ttags & con.svg thickness 7 italics fontsize 15; +export < con.tre; +taxname-; +export - con_no.tre; + +/*Caulculate TL/CI/RI score*/ length *; -stats.run; -quit ; +report-; +var-; +var = + 0 themin + 1 themax + + CIs [ ( ntrees + 1 ) ] + + RIs [ ( ntrees + 1 ) ] + + this +; +set themin minsteps; +set themax maxsteps; +loop 0 ntrees + progress #1 (ntrees+1) Calculating indices...; + set this length[#1]; + set CIs[#1] 'themin'/'this'; /*CI=1 means no homoplasy*/ + set RIs[#1] ('themax'-'this')/('themax'-'themin'); /*RI=1 character fits perfetcly*/ + stop +progress/; +report=; +macfloat 3; +maketable CIs Consistency index; +maketable RIs Retention index; + + +/*Report*/ + +quote - /----------------------------------------------\; +quote - | The analysis has been finished. |; +quote - | The file `tnt.log` contains |; +quote - | K, TL, CI and RI |; +quote - | The file `trees*.tre` contain |; +quote - | trees found by mult and xmult |; +quote - | The file `con*.tre` contain |; +quote - | strict consensus tree |; +quote - | The file `*_no.tre` contain |; +quote - | tree with `taxaname-` |; +quote - | The file `con.svg` contain |; +quote - | strict consensus tree with |; +quote - | relative bremer, jak and boot |; +quote - | support on the tree |; +quote - \----------------------------------------------/; + +/*Quit*/ +zzz; diff --git a/setk.run b/setk.run deleted file mode 100644 index 1cc962a..0000000 --- a/setk.run +++ /dev/null @@ -1,71 +0,0 @@ -macro= - -var: - actK - minK - maxK - gVal - a - b - Nval - maxIts -; - -if (ntax<0) - quote NO HAY DATOS!; - proc/; -end - -set gVal 1; - -/* busca el valor de g */ -loop 0 nchar - if (!isinfo[#1]) continue; end - if (!isact[#1]) continue; end - set a=maxsteps[#1]; - set b=minsteps[#1]; - if (('a'-'b')>'gVal') - set gVal 'a'-'b'; - end -stop - -quote valor de g = 'gVal'; - -set actK 10; -set minK 0; -set maxK 500; /* Asume 500 como el maximo posible valor de K */ - -set maxIts 0; - -/* busca el mejor valor de K */ -loop 0 1 - quote actual valor de k = 'actK'; - set a 1-('actK'/('actK'+1)); - quote 'a'; - set b ('actK'/('actK'+'gVal'-1))-('actK'/('actK'+'gVal')); - quote 'b'; - set Nval 'a'/'b'; - quote valor de N = 'Nval'; - if (('Nval'>14.8)&&('Nval'<15.2)) endloop; end /* N esta en el rango */ - if ('Nval'>15.2) - set minK 'actK'; - else - set maxK 'actK'; - end - set a ('maxK'-'minK')/2; - set actK 'minK' + 'a'; - set maxIts ++; - if ('maxIts'==100) endloop; end /* salida de emergencia */ - setloop 0; -stop - -/* Si salio de emergencia */ -if ('maxIts'==100) - quote NO SE TERMINARON LAS ITERACIONES; - quote mejor K encontrado: 'actK'; - proc /; -end - -quote Valor de K: 'actK' (N='Nval'); -piwe='actK'; -proc/; diff --git a/stats.run b/stats.run deleted file mode 100644 index dd28770..0000000 --- a/stats.run +++ /dev/null @@ -1,23 +0,0 @@ -macro= ; -report- ; -var = - 0 themin - 1 themax - + CIs [ ( ntrees + 1 ) ] - + RIs [ ( ntrees + 1 ) ] - + this -; -set themin minsteps ; -set themax maxsteps ; -loop 0 ntrees - progress #1 (ntrees+1) Calculating indices... ; - set this length[#1] ; - set CIs[#1] 'themin'/'this' ; - set RIs[#1] ('themax'-'this')/('themax'-'themin') ; - stop -progress/ ; -report= ; -macfloat 3 ; -maketable CIs Consistency index ; -maketable RIs Retention index ; -proc/ ;