FAQ¶

My species is not human — what do I need to change?¶

The prior on the population size. It is the one default that is genuinely tuned to humans, and it is easy to miss because it is never printed as a warning.

SweepLink samples log10N under a Gamma prior in shape–rate form, whose hyperparameters default to:

Argument	Default
`--log10N.alpha`	`16.0`	shape
`--log10N.lambda`	`4.0`	rate

That puts the mean of log10N at 16/4 = 4, i.e. N = 10,000, with a standard deviation of 1 — a full order of magnitude. 95% of the prior mass sits between roughly N = 190 and N = 1,500,000.

If your species sits comfortably inside that range you can leave it alone. If it does not — an inbred crop line, an experimental population, a domesticated animal, or anything with N in the hundreds or in the tens of millions — the prior will pull the estimate towards 10,000 and you should re-centre it.

To centre the prior on N = 10^k with a standard deviation of σ on the log10 scale:

alpha  = k² / σ²
lambda = k  / σ²

Some worked values:

Target N	σ (log10)	`--log10N.alpha`	`--log10N.lambda`	95% of prior mass
500	0.5	`29.14`	`10.8`	N = 64 – 5,800
1,000	0.5	`36`	`12`	N = 126 – 11,400
100,000	0.5	`100`	`20`	N = 11,700 – 1,060,000
1,000,000	0.75	`64`	`10.67`	N = 42,000 – 36,000,000

$ ./sweeplink infer ... --log10N.alpha 36 --log10N.lambda 12

Warning

--N does not change the prior. It only sets the value the chain starts from. To move the prior you must set --log10N.alpha and --log10N.lambda.

Two things worth knowing before you pick numbers:

N is the diploid effective size. The drift term is x(1-x)/2N, so there are 2N chromosomes. Do not halve or double your literature value to compensate.
You can look before you choose. The calculateLLForN task evaluates the likelihood over a grid of population sizes without running an MCMC, which shows where your data alone put N:
```
$ ./sweeplink calculateLLForN --counts counts.txt --meta meta.txt \
      --mu_a_A 1e-8 --minN 100 --maxN 1e7 --numN 100 --logN --sumLL --out probe
```
The result is a single row in probe_LL.txt; the largest value marks the best-supported N. Use it to sanity-check the prior, not to replace it.

Everything else in the model is species-agnostic. You will of course also want a mutation rate appropriate to your organism (--mu_a_A), and generation times that match how you filled the meta file.

I do not know which allele is ancestral¶

This is common — without an outgroup there is no way to polarise, so you have only REF and ALT.

Use whichever assignment you have, and change two defaults, both of which quietly assume that you do know the direction:

$ ./sweeplink infer ... \
      --mu_a_A 1e-8 --mu_A_a 1e-8 \
      --forwardPrior uniform --newMutPrior uniform

--mu_A_a 1e-8: Make the mutation rates symmetric. The default --mu_A_a is 0, which hard-codes “the derived allele can arise but never revert” — a statement about a direction you have not established. Setting both rates equal removes it.
--newMutPrior uniform: The default is neutral, a 1/x density describing a derived variant entering at low frequency. That is also a polarisation assumption, so make both priors uniform.

With those set, any REF/ALT assignment gives a valid run. What changes is the interpretation:

Important

The sign of s becomes arbitrary — it tells you which of your two arbitrary labels is rising, not whether the derived allele is favoured. The magnitude |s| is unaffected, and so is the evidence that the locus is under selection at all. Both panels of sweeplink-plot remain meaningful; the bottom one already plots |s|.

Does my data need to be phased?¶

No. SweepLink never resolves haplotypes. Its only input is a count c/n per locus, time point and population, so unphased genotypes are fine — you just need to count copies of an allele at each site.

How long does a run take, and how to make it faster?¶

Runtime scales with loci × time points × iterations. The tutorial’s 413 loci and 100,000 iterations take about 13 minutes on 4 threads.

--numThreads is the main lever — the likelihood parallelises across loci:

$ ./sweeplink infer ... --numThreads 8      # or 'all', or 'allButOne'

Warning

The default is one thread, not all cores. That is deliberate: the parallel sums do not accumulate in a fixed order, so a multi-threaded run is not reproducible bit-for-bit even with --fixedSeed. Use threads for speed, and one thread when you need an exactly repeatable run.

For a whole genome, run one job per chromosome rather than one enormous job. The companion scripts are built for that layout — put each run in its own directory with a shared name and use --dirs:

$ sweeplink-plot --dirs chr1 chr2 chr3 --run-name sweeplink

How many iterations do I need?¶

The defaults (100,000 iterations, 10 burn-in rounds of 1,000, thinning 10) are a reasonable starting point. To judge whether they were enough:

Watch the acceptance rates printed after each burn-in round. SweepLink retunes proposals between rounds, so you do not need to touch anything — but rates near 0 or near 1 after tuning mean the chain is struggling.
Look at the trace, not just the summary. <out>_trace.txt holds the retained samples; plot log10N with sweeplink-plot-posteriors and check the posterior is a sensible unimodal shape rather than drifting.
Compare two seeds. Two runs with different --fixedSeed that disagree have not converged.

Raise --iterations first; raise --numBurnin if the chain is still moving towards its final region when the main chain starts.

My run was killed. Do I have to start over?¶

Not from scratch. Every run writes <out>_state.txt, holding the current value and jump size of every parameter, so you can carry on from there — and skip burn-in, because the tuned jump sizes come back with it:

$ ./sweeplink infer ... --initVals dead_state.txt --numBurnin 0 --out rerun

Two things to know. If the state file is 0 bytes the run died during burn-in and there is nothing to resume from. And the dead run’s _statePosteriors.txt is empty regardless — those are only written at the end — so the resumed run still needs the full number of iterations you want the posterior built from. It saves the converging, not the sampling.

Full details in resuming an interrupted run.

What counts as missing data?¶

Write it as 0/0 — a locus with n = 0 at a time point simply contributes nothing there.

sweeplink-vcf2input writes 0/0 for you. It also offers filters for the cases where too much is missing to be useful — --max-missing-recent, --max-missing-total, --min-segregating-timepoints. Ancient DNA time series are usually thinnest at the oldest time points, which is worth checking in the table the script prints before you rely on those samples.

Common errors¶

Argument 'mu_a_A' must be defined: --mu_a_A is mandatory; there is no sensible default mutation rate.
String 'NA' is not a number!: Missing data in the counts file must be 0/0, not NA.
Size of grid (...) does not match the number of states for abs(s): --grid_abs_s must have exactly --numStatesAbsS entries.
Neither or both 'betaMean' and 'betaESS' should be specified!: Give both, or neither and let them be estimated from the data.
the state posteriors are not multiples of 1/<n>: A companion script was told the wrong number of MCMC iterations. Normally it reads this from <prefix>.parameters; pass --iterations only for runs made by older versions that did not record it.
no locus of the second run is present in the first: --round2-prefix points at a run that is not a refinement of --prefix.