search

CartwrightLab / dawg

Simulating Sequence Evolution
GNU General Public License v2.0
11 stars 3 forks source link

Sequence with no Adenine(A) using 'GTR' model to simulate #61

Closed yukimayuli-gmz closed 2 years ago

yukimayuli-gmz commented 4 years ago

Hi, I meet a problem during I used the dawg to simulate some sequence using the GTR model by some random trees. The input file for the dawg is(one of the random trees):

Tree.Tree = "( ( ( ( ( ( sp.19 : 0.0274309776261303 , sp.20 : 0.0363371779192595 ) : 0.0231000164056752 , ( ( ( sp.17 : 0.0154588571141289 , sp.18 : 0.00146999226656243 ) : 0.008786572370401 , ( sp.27 : 0.00179051684256983 , sp.28 : 0.00755231767160561 ) : 0.0635450661666089 ) : 0.00618192060060325 , sp.12 : 0.0337847168172262 ) : 0.0555383411224361 ) : 0.002542308222857 , sp.10 : 0.0382953974492483 ) : 0.00719380989169898 , ( sp.5 : 0.0756897850378422 , sp.6 : 0.0197184969210785 ) : 0.0228941653868867 ) : 0.00301641855849631 , ( sp.1 : 0.0951581028441154 , ( ( sp.21 : 0.00049993230701374 , sp.22 : 0.000456900812670007 ) : 0.202166829497331 , ( sp.7 : 0.0263569567773822 , sp.8 : 0.0346980961911223 ) : 0.00514439857638284 ) : 0.0570595001942826 ) : 0.00836196371217841 ) : 0.164987526960894 , ( ( sp.3 : 0.0686698084012003 , ( ( sp.29 : 0.0106409125382184 , sp.30 : 0.0049558033984565 ) : 0.0568440949960711 , sp.9 : 0.0148929885248412 ) : 0.0381298668406275 ) : 0.0279709721036675 , ( sp.2 : 0.0401266752960971 , ( ( sp.4 : 0.00874337946377212 , ( sp.15 : 0.017538496785908 , sp.16 : 0.00967390252271784 ) : 0.103838353162518 ) : 0.00111384934694459 , ( ( sp.11 : 0.0799675125536995 , ( ( sp.13 : 0.00190281104678998 , sp.14 : 0.126332344299632 ) : 0.00545928077950286 , ( sp.23 : 0.0170462282898584 , sp.24 : 0.00524201356675814 ) : 0.284812760193555 ) : 0.00368770510917451 ) : 0.00230321546575498 , ( sp.25 : 0.0128043584869294 , sp.26 : 0.0229539442894333 ) : 0.193254767154063 ) : 0.0204965286993815 ) : 0.0163225158532997 ) : 0.144957139198762 ) : 0.0309704561882311 ) ;" Subst.Model = "GTR" Subst.Params = {1, 1, 1, 1, 1, 1} Subst.Freqs = {0.25, 0.25, 0.25, 0.25} Root.Length = 10000

But in the output sequence file of dawg. There no Adenine(A) in the sequence and C, G, T are about 1/3 in the sequence respectively. Even I change the Subst.Freqs and Subst.Params by referring some other paper, there also no A in the output sequence.
And my computer system is CentOS Linux 7.

Best Wishes, Guo Mengzhen

reedacartwright commented 4 years ago

@jgarciamesa can you check if this affects dawg 2.0?

reedacartwright commented 4 years ago

@jgarciamesa Can you check if this bug affects on dawg 2.0?

yukimayuli-gmz commented 4 years ago

@jgarciamesa Can you check if this bug affects on dawg 2.0?

yes, the "dawg --version" will get "dawg 2-current-rUnknown".

jgarciamesa commented 4 years ago

@jgarciamesa Can you check if this bug affects on dawg 2.0?

As @yukimayuli pointed out it does affect the current version of dawg. I was able to reproduce this with the latest version of the develop branch. It seems the issue happens here

https://github.com/reedacartwright/dawg/blob/e471fdb1c47d443779c4d1c2a289b79ab339c020/src/lib/ma.cpp#L48

The values returned by read_section(*secit) replace the first values of substitution rates (1) and frequencies (0.25) both with a 0. Not sure why yet, I'll keep debugging.

jgarciamesa commented 4 years ago

Given they way the arguments are parsed, having a curly bracket { character in front of a value disrupts parsing it as in {0.25,0.25,0.25,0.25} or {1,1,1,1,1,1}. For now, please remove the curly brackets when specifying parameter values. @yukimayuli, let me know if after this dawg still misbehaves.

yukimayuli-gmz commented 4 years ago

Given they way the arguments are parsed, having a curly bracket { character in front of a value disrupts parsing it as in {0.25,0.25,0.25,0.25} or {1,1,1,1,1,1}. For now, please remove the curly brackets when specifying parameter values. @yukimayuli, let me know if after this dawg still misbehaves.

When I remove the curly brackets, the dawg works good. No problem I find.