hmorlon / PANDA

Phylogenetic ANalyses of DiversificAtion
24 stars 15 forks source link

Problematic lambda values "RPANDA" #20

Closed majonavarrete closed 5 years ago

majonavarrete commented 5 years ago

Dear all,

I am trying to fit an environmental-dependent model using paleoelevation data of the tropical Andes. I am also testing different sampling strategies since we know that there are some cryptic species within my group. The clade is composed of 90 species, so I am trying with f=90/90 (100%), 81/90 (90%), and 72/90 (80%). I run 10 models per sampling strategy (total=30 analyses) in RPANDA in which rates were modified by both an exponential and linear function of paleo-environment: (i) speciation rate varies with the environmental variable and no extinction; (ii) speciation rate varies with the environmental variable and extinction rate is constant; (iii) speciation rate is constant and extinction rate varies with the environmental variable; (iv) speciation and extinction rates both vary with the environmental variable. These models were compared to simpler diversification models, (v) namely a pure speciation model (Yule) and a (vi) constant rate birth-death models.

Everything was ok until I saw the results for the linear function of the model where speciation and extinction rates both vary with the environmental variable (which is also the best model according to delta AIC), specifically for the f=81/90 sampling strategy (90%). All the lambda values in the other 29 runs vary between 0.14343217 and 1.432733, but for some reason that I can't understand, the "problematic" model is giving me a lambda value of 10.830502. I would like to know what am I doing wrong. Thank you in advance, Majo

linLambdadeponaltandlinmudeponalt

f.lamb<-function(t,x,y){y[1] + (y[2] x)} f.mu<-function(t,x,y){y[1] + (y[2] x)} lamb_par<-c(0.05,0.01) mu_par<-c(0.005, 0.009) tot_time<-max(node.age(prist)$ages) result_exp20<- fit_env(prist, elev2, tot_time, f.lamb, f.mu, lamb_par, mu_par, df= dof, f = 81/90, meth = "Nelder-Mead", cst.lamb = FALSE, cst.mu = FALSE, expo.lamb = FALSE, expo.mu = FALSE, fix.mu = FALSE, dt=1e-3, cond = "crown")

I also have tried with 79/90, 80/90, 82/90, and 83/90. These are all of the outputs that I got:

72/90 model (80%): "environmental birth death"

LH: -228.9846 aicc: 466.4397 lamb_par: 1.432733 -0.321531 mu_par: 0.6416556 -0.1663267

79/90 model : "environmental birth death"

LH: -223.9706 aicc: 456.4119 lamb_par: 12.837443 -3.275568 mu_par: 10.750152 -2.734816

80/90 model: "environmental birth death"

LH: -230.036 aicc: 468.5427 lamb_par: 1.1869399 -0.2533093 mu_par: 2.5183912 -0.6622045

81/90 (90%) model: "environmental birth death"

LH: -222.3847 aicc: 453.24 lamb_par: 10.830502 -2.768822 mu_par: -8.996084 2.307149

82/90 model: "environmental birth death"

LH: -223.7253 aicc: 455.9212 lamb_par: 4.0083097 -0.9913247 mu_par: -4.765867 1.224885

83/90 model: "environmental birth death"

LH: -232.0119 aicc: 472.4944 lamb_par: 0.42181768 -0.05442131 mu_par: -1.1430388 0.3022878

90/90 model (100%): "environmental birth death"

LH: -228.9846 aicc: 466.4397 lamb_par: 1.432733 -0.321531 mu_par: 0.6416556 -0.1663267

hmorlon commented 5 years ago

Hello Majo, there is not necessarily something wrong. Note that the lambda values you are referring to are the lambda values when the environmental variable has a value of 0, which does not necessarily happen in your data (i guess you never have an elevation of 0 for the Andes in the timscale of your phylogenies?). Maybe you can plot the estimated lambda values as a function of time and see if the estimates are reasonable or not. I hope this helps. Best, Hélène

majonavarrete commented 5 years ago

Thank you for your time Helene!