The manuscript titled ’ Holocene deglaciation drove rapid genetic diversification of Atlantic walrus’ by Ruiz-Puerta et al. sequence ancient DNA from walrus bones to analyze full mitochondrial genome (mitoegenome) sequences. The mitogenome sequences are combined with already published mitogenomes and analyzed to explore phylogenetic structure, as well as demographic history of Atlantic walrus. The results show strong signatures of genetic population structure that matches well with the known contemporary observations. Divergence time estimation (Time of Most Recent Common Ancestor [TMRCA]) points to that divergence of all populations coincides with the time following the last Glacial Maximum (LGM), which provides compelling evidence for rapid colonization after the melting of the Atlantic ice sheets.

The manuscript is very well written and the research and analysis overall well conducted. As such, I find that the manuscript should be of general interest for the readership of the journal and people interested in Arctic phyleography.

I am very positive about the manuscript, however, I still have a few comments and suggestions that I believe should be addressed before the manuscript can be published. These are mostly of technical character and something I believe the authors should be able to include without any major revisions.

General.
I think it is somewhat unclear how many sequences are being analyzed. As far as I understand the authors generated 75 new sequences, but they also mention that they included some from earlier published work. However, from the text it is not clear how many of these sequences are being analyzed in the phylogeny. I suggest that the authors somewhere in the manuscript states exactly how many of the analyzed sequences are generated by the authors themselves and how many are from earlier published work.

Abstract
In the abstract the authors state that ‘Our results indicate that the significant levels of genetic diversity in contemporary Atlantic walrus populations was directly shaped by late Pleistocene-Holocene ice loss’. I think the sentence is somewhat unclear. I guess that authors refers to that the inferred population structure of the ancient samples matches well with what we see in comtemporary population studies. Maybe this can be emphasized more clearly.

Introduction
Page 4. ‘Clearly, the ways and rapidity with which…’ Suggestion: replace ‘rapidity’ with ‘speed’ or ‘pace’
Page 4. ‘ …abundance and distribution of Arctic marine biota is highly nuanced, possibly involving both “tabula rasa” and/or “polynya” scenarios’… Suggestion: replace ‘nuanced’ with ‘complex’.

Material and methods.
Can the authors elaborate on why they chose to use ANGSD to analyze the data. Since they use a sequence threshold of 3 reads to call bases, I wonder if ANDSD has any advantages compared to other base calling software?

Phylogenetic analysis.
I would like the authors to elaborate a bit on the choice of priors used in their phylogenetic analysis in BEAST 2. I think it is especially important to understand the reasoning behind using the Coalescent exponential population model and a relaxed clock exponential model. Is this something that is based on some prior knowledge or analyses? For example, did the authors test other models to understand the potential sensitivity of the data to different population and/or clock models? Such exploration of the data is normally suggested when applying Bayesian statistic like BEAST, and is something the authors developing the software in general suggest.

I further think it is important the authors present the mutation rate estimate somewhere in the paper. It may be useful for other studies. In such a case, it would further be interesting to compare the estimated mutation rate to the general literature. That should allow a validation of the mutation rate estimate and hence the overall conclusion of the paper of a post LGM divergence. Additionally, such a rate would further be useful in future genetic research on walrus.

Last, the authors use a sequence from the pacific walrus as outgroup in the analysis. The divergence between this sequence and the Atlantic walrus sequences is estimated to only roughly 35,000 years. This date is very far from other previous estimates based on shorter fragments of the mitogenome and would correspond to divergence within the last ice age, which is very unlikely. I expect that this could be explained by potential time-dependency (as suggested e.g. by Ho et al. 2005; 2007), and is something the authors need to discuss in the text to avoid any misunderstandings.

Bayesian skyline plots.
I think it is important that the authors consider the potential bias when analyzing the clades instead of the more confined populations. Heller et al. (2013) showed that BSPs are sensitive to population structure, which can lead to appearance of population decrease when not accounted for. As their data in general show strong evidence for tight population structure, even within the two clades, this possibility needs to be addressed somewhere in the text.

The authors should also write how time was estimated for the BSPs. Did the authors use the sequence ages as done for the phylogeny, or did they use the mutation rate estimate from this analysis as a strict prior for the BSPs?

Last, I wonder why the BSPs show the same TMRCA for both the Eastern and Western clade? From the phylogeny the TMRCA of the Eastern clade is almost twice as old as for the Western clade, which should be reflected by the BSPs. Can the authors elaborate on this? I wonder if it could be a mistake in the datasets, although it might also account for differences when applying a different model on the data.

The y-axis should be female effective population size. If the authors scaled the estimates to include full effective size, they need to mention this somewhere in the text.

References:
Heller R, Chikhi L, Siegismund HR (2013) The Confounding Effect of Population Structure on Bayesian Skyline Plot Inferences of Demographic History. PLoS ONE 8(5): e62992. https://doi.org/10.1371/journal.pone.0062992
Ho, S.Y.W., Phillips, M.J., Cooper, A. & Drummond, A.J. 2005. Time dependency of molecular rate estimates and systematic overestimation of recent divergence times. Mol Biol Evol 22: 1561-1568.
Simon Y W Ho, Beth Shapiro, Matthew J Phillips, Alan Cooper, Alexei J Drummond, Evidence for Time Dependency of Molecular Rate Estimates, Systematic Biology, Volume 56, Issue 3, June 2007, Pages 515–522, https://doi.org/10.1080/10635150701435401

To the editor and authors,

I have now read through the revised manuscript titled ’Holocene deglaciation drove rapid genetic diversification of Atlantic walrus’ by Ruiz-Puerta et al. I think that the authors have done a brilliant job in revising the manuscript and answered well to the comments by all the referees. As such, I feel that the manuscript is in very good shape and recommend publication after minor revision.

Comment #1. Line 232-2
The authors discuss the potential for an early wave of colonisers to Iceland. This is based on the observation of two Icelandic sequences, which show an older origin. I like the author to elaborate on whether the positions of these two sequences might have been impacted be sequencing error or sequence length of the two sequences.

Comment #2. Line 265. Remove duplicated ‘using’.

Comment #3. Line 309-311. The authors test for potential demographic changes using the bayesian skyline plots in BEAST 2 and discuss their findings. Specifically they mention the decline seen in the western clade and discuss potential reasons for it. While I very much acknowledge the observation then I wonder whether the signal in the data is actually strong enough to warrant this discussion, or whether it (perhabs) is too speculative. I suggest that the authors reanalyze for demographic change within the two clades using the extended bayesion skyline plot methods (EBSP) implemented in BEAST 2 (or potentially a similar method). This analysis estimates the number of demographic changes directly from the data and thereby tests for deviations from constant size. The outcome of such analysis would allow a statistical framework and more objective test for whether potential changes are indeed supported by the data.

Comment #4. Line 455-457. In this revision, the authors have tested both a strict and relaxed clock model and assessed the best model by comparing the likelihood. While I am not an expert in this then I think that a more appropriate way to do this is to estimate the marginal likelihood and then taking the difference to get the BF, which then can be used to assess the confidence of support. In the ‘old’ version of the ‘BEAST package’ there was a way to calculate this in TRACER, however, as far as I can see the estimator was not reliable and there are other more accurate ways to do this now (https://www.beast2.org/model-comparison/).