Lamarc is a tool for studying populations. It can estimate:

• Effective population size
• Migration rates between populations
• Changes to population size with time (growth rates)
• Recombination rates
• Variations of mutation rates across genomic regions (gamma)
• Fine-scale mapping of trait data

First off, you need genetic data. This might be standard genetic data like DNA or SNP sequences, or microsatellite counts, or it could (with version 2.0) be observations with genetic underpinnings, like phenotypic observations or electrophoretic data. More data will result in more accurate results, but it's important to note that the addition of more individuals, after a certain point (10-15 per population), is much less helpful than the addition of more independent genetic regions. The history of a known region in an additional individual is probably similar to that of its cohorts. The history of a new unlinked region from the same (or new) set of individuals is almost certainly very different from that of the old region, and therefore more likely to hold new information.

Secondly, you need a computer. In particular, you probably need a computer with a lot of memory. A good LAMARC work-out can use hundreds of megabytes of RAM, and for particularly ambitious analyses with lots of data, gigabytes are probably called for. For learning purposes, LAMARC can run on more limited computers, if you're careful to give it a limited data set and not allow it to record too many observations.

Finally, you need some time. The total amount of time can vary from data set to data set, but it's not unusual for a complete, solid run of LAMARC to take a week or two. Larger data sets requiring full error analysis (profiling) can take a month or longer. However, it would be a mistake to throw all your data at LAMARC in January, then write to us in May asking when you could expect results. The answer would probably be a tragic, "LAMARC ran out of memory, and inexplicably failed to crash." Much better, instead, to test LAMARC on subsets of your data and get preliminary results in a day, which will then help you determine how best to analyze your full data set. This tutorial is designed to help you do just that.

Hang on, let's take this slowly. If this is your first time using LAMARC, I strongly recommend using only a subset of your data, to get a feel for how the program works, first.

If you have multiple populations, let's select two of them to start with. If you know something about the migration patterns of the populace already, pick two populations you know share migrants. Then select around 10 individuals from each population. If you're only studying one population, pick 10-15 individuals. If you want accurate results, it's vitally important that you pick these individuals randomly from your data set, and not on the basis of them being the 'most different' or 'most interesting'. You can grossly distort your results should you rely on such a scheme.

Next, pick one region to analyze. If you have DNA or SNP data, that means one continuous stretch of sequence. If you have microsat data, that probably means just one microsatellite, though if you happen to have linked microsat data, it could be those linked microsats.

You are ready to run LAMARC. The very first thing LAMARC asks you for is the name of a file, and the converter-created file is what it wants. You can edit the contents of that file if you wish (it's all XML), but it's probably much easier to change things in the menu. If you like, you can even enter the name of your file at the command line ('lamarc yourfile.xml'), and if you don't want to do anything in the menu, you can run LAMARC in 'batch mode' by adding a '-b' flag ('lamarc yourfile.xml -b').

For a first crack at your data, yes. It'll use a lot of default values that are perhaps not the best choice for your data set, but they'll work for now. After you see how it runs, we'll get around to changing the defaults.

OK, hang on, I'm going to go try it.

I'll be here.

After creating a tree, LAMARC is rearranging it to see if it can find better ones. Each time it does so, it increments the counter. Using the defaults, it first counts to zero from -1000 without sampling anything (this is the phase known as 'burn-in'), then it starts counting up to 10,000, sampling every 20th tree. Those trees will be used to estimate parameters, which you'll see in a moment.

OK, it finished counting and now it's stopped.

It's busy calculating the max--

No, never mind, there we go. OK, it gave me a lot of information and started counting again. What does all this mean?

OK, your screen probably looks something like this:

15:39:34  Initial chain   1:  [====================]        10000 steps
 
15:40:09  Predicted end of chains for this region:  Tue Feb  8 16:06:43
          2005
 
15:40:09  Accepted  7.39% | Posterior lnL 27.0287215 | Data lnL -627.243228
Trees discarded due to too many events:        2
Trees discarded due to too small population sizes:        0
Trees discarded due to an infinitesimal data likelihood:        0
Trees discarded due to extremely long branch lengths:        0
Tree-Arranger accepted            284/8268 proposals
Tree-Size-Arranger accepted       455/1732 proposals
  Class                  Theta
  Population 1          0.039795
  Population 2          0.065072
   
  Population                     Mig
  Population 1            --------  164.6673
  Population 2            142.1566  --------

15:40:09  Initial chain   2:  [=|                  ]           35

Let's take a look at each piece in turn.

Predicted end of chains

LAMARC's first estimate of when it will be done. This estimate gets revised at the end of every chain, and hopefully gets more accurate as it goes.

Accepted 7.39%

The number of proposed changes to the tree that were accepted by LAMARC. This number is usually in the 10% range--if it's very low, it means that it's latched on to some form of a tree, and is very reluctant to vary from that, which means that it probably is not exploring tree-space very effectively. Too high, and it's accepting almost everything that comes its way, which probably means it's not proposing radical-enough changes, or else that your data give it no reason to prefer one tree over another (such data are likely mostly noise, so they aren't very informative about population parameters - if you run into this, go review your input data for obvious problems like missing data).

Posterior lnL 27.0287215

The posterior log likelihood. This number (which should always be positive) indicates how much better the maximized parameters are than the initial, driving set of parameters. In early chains (like this one), it may be quite high, but it will drop down to the single digits in later chains. Please note that this is a relative likelihood, valid only within this specific run; it cannot be used in log-likelihood tests or other statistical applications.

Data lnL -627.243228

The data log likelihood. This number will probably be very negative, should increase for the first few chains, then level off for the last few chains. It represents the log of the probability of the data given the currently accepted tree. Since a large data set is highly improbable on any tree, it's not a problem that this number is so low. However, if it is still rising rapidly by the end of your run, you haven't yet found the best trees and your run is too short.

One other thing you should note is that if you're comparing this run to a new run with more of your data, the data log likelihood will go down. This is due to the fact that a tree with more tips will generally fit the data at those tips more poorly than a tree with fewer tips.

Trees discarded due to...

Sometimes, LAMARC will discard trees because the trees themselves are inherently too tricky to deal with. Almost always, these trees would also be rejected from having too low a likelihood, so you shouldn't worry about this too much unless one of these numbers gets very high (say, larger than 5% of the total number of proposed trees). If that happens, your starting parameters might be too extreme, or you might be calling two populations different when they are actually genetically identical (rejections due to 'too many events' can have this cause). Generally, though, it's just LAMARC being efficient.

If absolutely no trees are rejected, you'll see the message "No trees discarded due to limit violations." which means you're fine.

Arranger accepted

This is a more detailed breakdown of the 'Accepted 7.39%', above. The two arrangers on by default in LAMARC are the Tree-Arranger, which breaks a branch of a tree and then re-attaches it, and the Tree-Size-Arranger, which preserves the topology of the tree but picks new sizes for some or all of the branches. The Tree-Arranger is absolutely required, since it's the only one to change tree topologies. The Tree-Size-Arranger is more of a helper function, which is why (by default) it is only attempted 1/5 as often as the Tree-Arranger.

These numbers can vary fairly widely, but pay attention if they get too large or too small, just as you would to the overall 'Accepted X%' line, above.

Theta and Mig

These are the parameters LAMARC is trying to estimate. You will be estimating one theta value for every population in your data, and two migration rates for every pair of populations in your data. In this case, with two populations, that means a theta for each (0.039795 for population 1 and 0.065072 for population 2), and two migration rates (164.6673 for the rate from pop2 to pop1, and 142.1566 for the rate from pop1 to pop2). These are LAMARC's best estimates of the parameters, as judged from one 'chain'-worth of trees. They will be used as driving values for the next chain, as will the last tree sampled in the chain. These reported values for the last chain will be reported as LAMARC's estimates of your parameters.

Initial Chain 2

LAMARC is starting over, essentially. Except that this time, it hopefully has a halfway-decent tree to start with, and better driving values than it did the first time. It'll repeat this process (by default) until it's done it 10 times.

By default, LAMARC runs 10 'Initial' chains and 2 'Final' chains. The initial chains take 10,000 steps, and the final chains take 200,000 steps. This is LAMARC's attempt to solve a basic problem in likelihood analysis: LAMARC does a great search in the area right around the driving values, but not so great a search elsewhere. By using multiple chains, LAMARC moves to driving values as close to the true maximum as possible, so that it can find that maximum efficiently. The first 10 chains are an attempt to find the right neighborhood, and the first final chain is an attempt to narrow it down even further. The last final chain is used to actually estimate your parameters, and to provide support intervals for those estimates.

If you can see the message "Beginning profiling, please be patient", it means that LAMARC has entered a computationally-intense phase where it calculates the support intervals for its point estimates. We call this phase 'Profiling', since it provides profiles of the likelihood surface as seen from the perspective of each parameter in turn. After a while, your screen will look something like this:

03:44:15  Beginning profiling, please be patient
04:06:09  Finished profile 1 of 7.  Predicted end of this set of profiles:
          Sat Apr  9 06:17:33 2005
 
04:36:23  Finished profile 2 of 7.  Predicted end of this set of profiles:
          Sat Apr  9 06:46:43 2005
 
04:37:26  Finished profile 3 of 7.  Predicted end of this set of profiles:
          Sat Apr  9 05:48:20 2005

It can take an extremely variable amount of time to finish profiling, but a good rule of thumb is 'about half again as long as it's taken so far'. LAMARC provides estimates of when it expects to be done as it finishes each profile in turn, predicting that each profile it has to go will take about as long as the the average profile so far. For an exploratory run like we've assembled here, a few hours is not unreasonable.

If you have more than one genetic region in your data, it first goes through the chains for the first region, then calculates the profiles for that region, then goes on to the next region. So what you're seeing is the beginnings of LAMARC calculating the data from the second region. Back when we were getting started, I suggested only using data from one region, and this is why--each region goes through the same process, which can be rather lengthy. But, no harm done--it just means a longer wait.

OK, it's finally done. So, what can I tell from this 'outfile' it produced?

Sorry. OK, here's the basic questions these sections answer:

What are the values for my parameters?

The values for your parameters are the first row of the 'Maximum Likelihood Estimates' tables; the ones labeled 'MLE'.

What are the confidence intervals for my parameters?

Technically for this kind of analysis they're called 'support intervals', so the standard 95% support intervals can be found by looking for the rows that start with '95%'. So, in our Theta1 example above, the lower 95% point is at 0.004314, the upper is at 0.007735, and the MLE is at 0.005697. One way to write this is:

0.005697 +.002083/-.001383

Since this is clearly too many digits of precision, you might round off to:

0.0057 +.0021/-.0013

If you thought this was confusing and were willing to sacrifice precision, you could average this error to:

0.0057 +/-.0017

even though it's, you know, wrong. Also keep in mind that results from likelihood LAMARC are believed to estimate too-narrow confidence intervals due to its reliance on the driving values. Runs with replication get away from this problem.

What else can I say about my estimated parameters?

You can comment on parameter correlations by looking for trends in the 'Profile Likelihoods' section. For example, you could look at M21 as compared to M12, and if a high M21 correlated with a low M12 and vice versa, you could talk about the total amount of migration between those two populations.

OK, I think I get it. Let's go to the next sections.

Gladly! Actually the next bits don't need a lot of explanation--the next section ("User Specified Options") is a summary of your starting values, search strategy, and filenames. The next section ("Data summary") is an overview of what your data looked like, what data models you used, and an actual copy of the raw data. And the final section ("Run Reports by Region") is a copy of most of what was sent to the screen during the LAMARC run (sans the progress bar, of course).

OK, then! I'm ready to publish!

Ha ha.

Er, right. So, what do I do next?