Sequence order in the input and parallelism affect reproducibility

[apcamargo]

I conducted an experiment to evaluate the variance in phylogenetic trees reconstructed from a gene. I observed that trees reconstructed from two files with identical MSAs but different sequence orders produced different trees. This was assessed using tree distance metrics and visual inspection. Despite setting a constant seed (via --seed), the two files still resulted in distinct trees. Additionally, when using the same input file (therefore, the same sequence order) and the same seed, the resulting trees were still different.

The same substitution model was used across all tests:

iqtree2 -T 8 -s test_msa.afa -m WAG+G4 --prefix test_tree --nstop 50 --seed 8787

Using a single thread (-nt 1) resulted in identical trees when the input was identical. However, shuffling the order of the sequences in the input FASTA led to distinct trees.

I could not find any documentation indicating that using multiple threads would make results irreproducible or that changing the sequence order in the input file would lead to different results, regardless of the seed and number of threads. Am I missing something?

Version:

IQ-TREE multicore version 2.3.4 COVID-edition for Linux x86 64-bit built Apr 26 2024

[bqminh]

There are two questions:

• Can sequence order input change the tree topology? YES, even with the same seed number. IQ-TREE at the beginning uses random stepwise addition technique to construct initial parsimony trees, where different orders may result in different parsimony trees, and subsequent tree search may diverge in different trees.

• Can different thread number result in different trees? YES, even with the same seed number. The explanation is a bit involved but revolved around this inequality:

(a+b)+c != a +(b+c)

While they are mathematically equivalent, they are not always the same in a computer due to rounding errors (there are many articles about this on the web). As an example, let's say we have an alignment with 10 sites. The tree log-likelihood is the sum over site log-likelihoods: L1+L2+...+L10. However, the sum might be (very slightly) different depending on which order you do the summation. If there is one thread, the sum is taken from left to right: ((((L1+L2)+L3)+...+L9)+L10. However, if there are two threads, thread no. 1 will compute L1+L2+L3+L4+L5 from left to right, whereas thread no. 2 will do L6+L7+L8+L9+L10 (left to right), and the results of two threads will be added together. That's where the likelihood calculation can be different between runs of different thread numbers, but the difference can be expected to be very tiny. However, for a dataset with little phylogenetic information (e.g. when you have few sites compared with number of sequences - how many sites do you have?), there can be two trees T1 and T2 with very similar log-likelihoods, and a small difference might actually make L(T1)<L(T2) under 1 thread but L(T1)>L(T2) under 2 threads. And the tree search might subsequently diverge into different trees.

BONUS:

• Can two single-threaded runs with the same seed number result in different trees in different computers? YES, they can. The explanation is again a bit involved: IQ-TREE makes use of the vector instruction set to speed up the likelihood calculation. There are mainly two vector instruction sets: SSE which allows CPU to perform two arithmetic calculations at once, whereas AVX allows for four concurrent calculation. Under SSE, the above log-likelihood calculation will become (L1+L3+L5+L7+L9)+(L2+L4+L6+L8+L10). Under AVX it becomes (L1+L5+L9)+(L2+L6+L10)+(L3+L7)+(L4+L8) (the four sub-sums within parenthesis are calculated concurrently). And as explained above, these two results can be slightly different. IQ-TREE automatically detects which instruction sets does your computer have. Thus, three computers, one without vector instructions, one with SSE, and one with AVX, may result in different trees.

I believe this is the explanation for a paper (https://doi.org/10.1038/s41467-020-20005-6) about irreproducibility in IQ-TREE and RAxML.

However, you don't really have to worry about this, because it happens only for dataset with little phylogenetic information. But for such dataset, there is much bigger issue in the reliability of the trees (e.g. low support values) and users need to do multiple runs with different seed numbers anyway and other analyses to rectify the results.

Hope that helps.

[bqminh]

PS: At the beginning of the log file, IQ-TREE print several lines like this:

Host:    xxxxxx (AVX2, FMA3, 16 GB RAM)
Command: iqtree2 -s example.phy 
Seed:    241511 (Using SPRNG - Scalable Parallel Random Number Generator)
Time:    Mon May 20 17:08:01 2024
Kernel:  AVX+FMA - 1 threads (8 CPU cores detected)

The results are only reproducible if Command, Seed and Kernel with 1 thread are identical (and of course identical input alignment).

[apcamargo]

Beautiful explanation! Thank you!

However, you don't really have to worry about this, because it happens only for dataset with little phylogenetic information. But for such dataset, there is much bigger issue in the reliability of the trees (e.g. low support values) and users need to do multiple runs with different seed numbers anyway and other analyses to rectify the results.

My case was precisely this. I was using alternative MSAs to estimate the confidence interval of a metric that depends on the tree topology. When I decided to double check if the variance was due only to the MSA or also to the phylogeny reconstruction, I ended up finding that there's variation in the tree topology even with a constant seed.

This isn't an issue to me, since I was precisely measuring variance, but I decided to open up a bug anyway because (1) I feel that it could be beneficial to users if this was mentioned in the documentation somewhere (as there's little explanation about the seed parameter); (2) wanted to have this discussed somewhere in case someone googles this in the future.

[roblanf]

Since this isn't an issue, but expected behaviour (noting that it's not expected by everyone, particularly if you haven't spent a long time thinking about how computers really work), I suggest we move this to the discussion forum.

[apcamargo]

Sure! Thanks again for the answer.

[ktmeaton]

Just wanted to say thanks @apcamargo and @bqminh for documenting this so thoroughly! I was hoping to get some exact reproducibility of trees for unit testing and continuous integration tests. The combination of 1) Identical Input Files, 2) random seed --seed, and 3) single threaded (-T 1) was exactly what I needed to get the same output across multiple runs on the same computer. Thanks so much!

[apcamargo]

I appreciate that :)