Issues Estimating Theta

[Sean Canfield]

Aloha Peter et al.,

I am analyzing 4-5 populations using a single mitochondrial locus. No matter how much I reduce the upper bound for theta (uniform prior), I get the triangular posterior distribution (see attached). Does anyone have any tips or advice for this scenario?

Thanks,

Sean

[Sean Canfield]

I should note that while for this run I used a random starting value for theta, in other runs I use the 10% prior value and get the same result (I was just experimenting a bit here).

[Peter Beerli]

Sean,

could you add a positive number random seed for the subsampling, I do not believe that is the problem, but it should a positive number.

And could you show me a run with an upper bound for theta of 0.1, I assume you know how variable the dataset is [and it is >0 ;-)]

mtDNA dataset of behave difficult without rate heterogeneity taken into account, do you know the alpha parameter for the gamma distribution?

if so then you can use the menu to set up, say, 3 mutation classes [see menu Data and there  6   One region of substitution rates?

change that to 3 and then use your alpha].

Peter

--
You received this message because you are subscribed to the Google Groups "migrate-support" group.
To unsubscribe from this group and stop receiving emails from it, send an email to migrate-suppo...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/migrate-support/999a1778-da75-4639-89e1-20ff036ba9c6n%40googlegroups.com.
<migrate-n_het_fullmodel_noGI_out3.pdf>

[Sean Canfield]

Hi Peter,

I'll re-run the program with a positive seed number and an upper bound of 0.1 for theta (I got rid of the last result with a 0.1 upper bound for some reason). I'll post that when it's completed (shouldn't take more than a couple of hours on my laptop).

The dataset is variable, but not highly...there are 12 variable sites among the 724 basepairs, with a lot of singleton haplotypes. I'm not opposed to using all of the individuals for the analysis simultaneously, but there are unequal sample sizes ranging from 12 to 131 (if I exclude the 12 sample population, the next smallest is ~ 30). I'm not sure how much this will affect the analysis. But I'm worried that excluding the singletons, which comprise > 50% of the total haplotype diversity, might also be misleading. 

Do you have any advice for selecting alpha? I'm a bit in the weeds here, but I love tinkering with this stuff, so I'm eager to learn.

Thanks as always for your help,

Sean

[Peter Beerli]

Sean,

To estimate alpha I usually use paup*, Essentially

#Load the data

#create a Neighborjoining tree (nj)

set crit=like

# set the HKY model, + basefreq= estimate and the allow to estimate the rate parameter alpha

lscore

then use that value in migrate, but commonly with mtdna simply using an alpha=0.2 will do because all the mtdna dataset I have seen

use use an alpha similar to that. If you need more instruction within migrate about alpha let me know.

Concerning samples, I think 30 should be a good sample size per population, no problem with that.

When you send the output send me also the parmfile

best,

Peter

To view this discussion on the web visit https://groups.google.com/d/msgid/migrate-support/deabd9ac-8463-43c1-b493-c4d7ad6f73ebn%40googlegroups.com.

[Sean Canfield]

Thanks Peter! I'll get on that analysis right now. I think I can intuitively see why ~ 0.2 would work in this case. Will post again when the analyses are completed!

- Sean

[Sean Canfield]

Hi again,

Here are the results and the parmfile. Alpha was set to 0.2, priors for theta were uniform [0, 0.1].

The program is still struggling with theta, with a triangular posterior distribution and equal estimates for each population. I'm at a bit of a loss!

Thanks,

Sean

[Peter Beerli]

Dear Sean,

after some probing I found the problem!

your parmfile specifies the inheritance scalar option as

inheritance-scalars=YES:{1.00}

this is correct except that I messed up at one point in the code, originally in version 4.x the option was 

inheritance-scalars={1.00}

but while working on version 5.x I moved all these type of options to 

inheritance-scalars=YES:{1.00} or

inheritance-scalars=NO

the confusion arose when I started to edit the comments on the parmfile to the new version but did not update the corresponding code that reads the option, but unfortunately I updated the rewriting parmfile in the menu to the newer version.

TL;DR inheritance-scalars=YES:{1.00} leads to a problem and the inheritance scalar becomes ~0.0 which leads to a the result in the posterior and wrong results. The fix is easy: use

inheritance-scalars={1.00}

and it will work,

hope this helps all who see weird triangular posteriors.

Peter

To view this discussion on the web visit https://groups.google.com/d/msgid/migrate-support/2f913779-63cf-4db7-aa6f-7b1d1153852cn%40googlegroups.com.
<migrate-n_gamma_out.pdf><parmfile_prelim_noGI_gamma>

[Sean Canfield]

Wonderful, thanks Peter! I've been running it myself over here and that has definitely taken care of the problem.

- Sean