A review of Resolving the Geochemical Provenance of Stonehenge Bluestones: A Sample-Size-Independent Multivariate Framework (G. W. Taylor, June 2026)
Taylor's paper, DOI:10.13140/RG.2.2.17039.96164 , sets out to demonstrate that the Stonehenge bluestones are geochemically identical to outcrops in the Mynydd Preseli, using rare earth element ratio data and three complementary techniques: PERMANOVA, SIMPER, and Euclidean nearest-neighbour distance mapping. The paper is generous with its data. Table 1 gives the full REE ratio matrix for all forty-two samples, Table 2 the complete pairwise PERMANOVA output, and Table 3 the entire 20 × 22 distance matrix with nearest-neighbour assignments. That openness is unusual and welcome, and it is what makes a substantive review possible. Everything below is derived from those three tables. Nothing in this review depends on data the author has not published.
The conclusion the paper reaches is, in outline, correct. The problem is that the method used to reach it cannot distinguish that conclusion from its opposite — and the paper's own data contain the proof.
1. What the paper argues
The argument runs in three stages.
A global PERMANOVA on the pooled data fails to reject the null hypothesis of no difference between Stonehenge and Welsh samples (p = 0.85, pseudo-F = 0.074). The author treats this as positive support for common origin. Because the pairwise PERMANOVA matrix contains a handful of rejections, and because archaeological sample sizes are small, the author judges the pairwise results unreliable and sets them aside.
To bypass the sample-size problem, all Stonehenge samples are pooled into one group and all Welsh samples into another, and SIMPER is used to compare group means. The means agree closely across all twelve ratios — La/Lu at 15.2 against 15.3, La/Yb at 2.12 against 2.12 — which the paper presents as demonstrating homogeneity.
Finally, Euclidean distances are computed between every Stonehenge sample and every Welsh sample. Several minima are very small, and these are offered as sample-level confirmation, with the two smallest presented as the headline results.
Stated at its strongest, the argument is: three independent methods at three different scales all fail to find a difference, so there is no difference to find.
2. The conclusion is almost certainly right
It is worth being clear about this before going further. That the great majority of the Stonehenge bluestones derive from the Mynydd Preseli has been the settled position since Thomas in 1923, and the last fifteen years of work by Bevins, Ixer, Pearce and colleagues has narrowed several of the lithologies to individual outcrops — Carn Goedog for the spotted dolerites, Craig Rhos-y-felin for the foliated rhyolites. Nobody reading this review needs persuading of the Preseli connection.
That is precisely why the paper is worth examining carefully. We already know the answer. A method that is applied to a case with a known answer, and which cannot recover that answer reliably, has been shown not to work — and the demonstration is much cleaner than it would be on an open question. What follows is not a defence of some rival provenance. It is an argument that the framework proposed here would have produced the same confident result had the stones come from anywhere.
3. What a provenance method has to demonstrate
A geochemical fingerprint is only useful if it discriminates. Showing that a Stonehenge sample resembles a Welsh outcrop establishes nothing on its own; the question is always whether it resembles that outcrop more than it resembles the alternatives, and by a margin larger than the measurement error.
This gives two requirements that any provenance study must meet:
- A comparison set. At least one candidate source outside the favoured region, so that a match can be shown to be selective rather than universal.
- A resolution limit. Some estimate of how much two measurements of the same rock differ, so that "close" can be distinguished from "indistinguishable given the noise."
The paper meets neither. No non-Welsh source appears anywhere in the analysis, so the specificity of the match is never tested. And no error estimate is offered — indeed the framework is presented as a way of avoiding variance estimates rather than as a way of quantifying them.
The second gap turns out to be recoverable from the author's own data, and doing so is the substance of this review.
4. The replicates
The dataset contains two analytical replicate pairs. On the Stonehenge side, OU10 and OU10 rpt are the same sample analysed twice. On the Welsh side, PCM7 and PCM7 rpt are the same sample analysed twice. The author includes both, and treats each pair as two independent observations.
They are better used as a control. Two measurements of one rock should, if the method has any resolving power, be closer to each other than either is to a genuinely different rock. The distance between the replicates is therefore an estimate of the method's noise floor.
Table 3 does not report Welsh-to-Welsh or Stonehenge-to-Stonehenge distances, so the replicate distances are not given directly. They can nonetheless be bounded from below without any additional data. For any three samples A, B and C, the triangle inequality requires
d(A, B) ≥ | d(A, C) − d(B, C) |
so any column of Table 3 in which the two replicates differ places a lower bound on the distance between them. Taking the largest such difference in each case:
| Replicate pair | Reference sample | Distances | Lower bound on d |
|---|---|---|---|
| OU10 / OU10 rpt | PCG20 | 1.092765, 0.917049 | ≥ 0.1757 |
| PCM7 / PCM7 rpt | OU15 | 0.750, 1.022 | ≥ 0.2720 |
Now set those figures against the paper's results.
The headline match is smaller than the noise floor. The abstract, the methodology and the conclusion all cite OU11 to CGD2 at d = 0.143 as the exemplary near-zero pairing. A sample in this dataset sits at least 0.176 from itself. The flagship match is closer than a rock is to its own second analysis.
Half the reported matches fall inside the noise floor. Of the twenty-two nearest-neighbour distances in the "Closest Euclidean" row, eleven are below 0.272: OU11, OU12, OU14, OU19A, SH33, SH37, SH65, OU10, OU10 rpt, SH43 and SH61. For these samples the assignment carries no information at all.
The replicates disagree with each other about provenance. This is the clearest demonstration available. OU10 and OU10 rpt are one stone:
| Sample | Nearest Welsh outcrop | Group | d | Runner-up | Group | d |
|---|---|---|---|---|---|---|
| OU10 | CGD2 | PO1 | 0.203463 | PCM30 | PO3 | 0.317994 |
| OU10 rpt | PCM30 | PO3 | 0.222013 | CGD2 | PO1 | 0.345259 |
The same physical stone is assigned to Carn Goedog on one analytical run and to a PO3 locality on the other. On the map in the paper's Figure 1 these are different red boxes. The ranking is reversed by the difference between two measurements of one rock.
The Welsh-side replicate behaves the same way. Against SH33, PCM7 rpt ranks first at 0.196369 while PCM7 — the same rock — ranks fifth at 0.422311, behind PCM31, PCGF29 and PCM32.
The margins between competing sources are far smaller than the noise. SH61 sits at 0.195945 from PCC11 (group PO1) and 0.218942 from PCAW49 (group PO3): a separation of 0.023, roughly a tenth of the replicate bound. SH37 has four candidates within 0.08 of each other — PCM32 at 0.1886, PCM7 rpt at 0.1969, PCM7 at 0.2009, PCM31 at 0.2677 — spanning two different Welsh groups.
None of this requires any judgement about geology. It is arithmetic on Table 3, and any reader can repeat it.
5. The assignments contradict the paper's own groups
Reading the bottom two rows of Table 3 together produces the paper's actual provenance result: each Stonehenge sample matched to a Welsh source group.
| Stonehenge group | n | Welsh groups assigned |
|---|---|---|
| SLF1 | 7 | PO1, PO2ii, PO1, PO2ii, PO2ii, PO1, PO2iv |
| SODC1 | 4 | PO2iv, PO3, PO1, PO3 |
| SLF3 | 7 | PO1, PO3, PO3, PO3, PO1, PO1, PO2iii |
Each Stonehenge group is scattered across three Welsh source groups. If SLF1 is a genuine petrological grouping — and the paper treats it as one throughout — its members came from one place. A method that distributes a single group across three outcrop clusters is not recovering provenance; it is sorting noise. The map in Figure 1, which presents the assignments as tight and localised, does not reflect what Table 3 says.
6. What the distances are actually measuring
The methodology states that the distances are computed on normalised REE ratios. They are not. The calculation can be reproduced exactly from the untransformed values in Table 1: taking OU11 and CGD2 and summing the squared differences across all twelve ratios gives 0.020460, whose square root is 0.14304 — the reported 0.143. The full working is in the appendix below.
This matters because the twelve variables are on wildly different scales. Across the dataset La/Lu ranges from about 13.8 to 16.8, a spread of 3.1, while La/Ce ranges from 0.355 to 0.407, a spread of 0.05. In an unstandardised Euclidean distance the contribution of each variable goes as the square of its difference, so La/Ce can contribute at most about 0.0025 to a squared distance while La/Lu can contribute over 9.
The consequence is visible in the largest distance in the matrix, SH49 to PCDL26 at 3.7325. Decomposing it: La/Lu supplies 75.5% of the squared distance, La/Tb a further 13.0%, and La/Ho 6.1%. Three of the twelve variables account for 94.5% of the result. The remaining nine are, for practical purposes, absent.
So the twelve-variable multivariate distance is largely a single-variable comparison of La/Lu, dressed in twelve dimensions. Two further points follow.
First, all twelve variables share La as numerator. Dividing a common quantity by twelve different denominators induces strong correlation among the resulting ratios by arithmetic alone — the spurious correlation of ratios described by Chayes. The PCA reporting 87.47% of variance on the first component is not evidence that the projection is dependable; it is largely a measurement of that induced correlation. (There is also an unresolved inconsistency here: an unstandardised PCA on these ratios would put considerably more than 87% on PC1, so the PCA appears to have been run on the correlation matrix while the distances were not standardised at all. Whichever is intended, the two analyses are not on the same footing and cannot corroborate one another.)
Second, ratios of compositional parts are not amenable to ordinary Euclidean geometry. Centred log-ratio transformation exists for exactly this case, and the paper's appendix presents the absence of transformation as a methodological virtue.
7. The statistical framework
Three problems, stated compactly.
Failing to reject H₀ is not evidence for H₀. This is the load-bearing assumption of the entire paper, and the reasoning is inverted at the point where it matters most. The author's stated justification for trusting the non-rejection is that sample sizes are small and variances tight — but those are precisely the conditions under which a non-significant result is least informative. Low power makes p > 0.05 the expected outcome whether or not a difference exists.
The appendix result should have made this visible. A pooled pseudo-F of 0.074 means that variation between the two groups is about 7% of the variation within them. That is not a signature of common origin. It is a statement that the grouping explains essentially nothing — which is what happens when heterogeneous populations are pooled.
The pairwise matrix is uninformative in both directions, for reasons of design. In a permutation test the smallest attainable p-value is fixed by the group sizes: with n₁ and n₂ samples there are C(n₁+n₂, n₁) distinct partitions, so no p-value below 1/C(n₁+n₂, n₁) is achievable, however large the true difference. Reading the group sizes off Table 1 and comparing with Table 2:
| Comparison | Group sizes | Partitions | Minimum possible p | Reported p | Reported F |
|---|---|---|---|---|---|
| SOF1 vs SODC1 | 1, 4 | 5 | 0.200 | 0.197 | 5.317 |
| SOF1 vs PO2iii | 1, 3 | 4 | 0.250 | 0.245 | 37.24 |
| SOF1 vs PO2iv | 1, 3 | 4 | 0.250 | 0.251 | 92.41 |
| SODC1 vs PO2ii | 4, 3 | 35 | 0.029 | 0.029 | 16.3 |
| SLF3 vs PO2ii | 7, 3 | 120 | 0.008 | 0.009 | 9.532 |
| SLF3 vs PO2iv | 7, 3 | 120 | 0.008 | 0.009 | 7.611 |
Every one of these results sits at the boundary of what the test could return. The comparison the paper singles out as its showcase anomaly — SOF1 against SODC1, where a pseudo-F of 5.317 accompanies a non-significant p of 0.197 — is a test in which no value below 0.2 was arithmetically possible. It is read as evidence of a geochemical match. Meanwhile the three rejections at 0.029 and 0.009 are the most extreme outcomes the design permits, and they are the results discarded as artefacts.
Two of the groups, SOF1 and SLF2vi, contain a single sample. PERMANOVA on a group of one has no within-group variance to estimate, and every entry in those rows should be removed rather than interpreted.
The paper's diagnosis of this problem is also worth correcting. It argues that a single critical F value of 2.1532 cannot apply to pairs of differing sizes because the degrees of freedom differ. The deeper point is that the pseudo-F in PERMANOVA does not follow an F distribution at all — that is the reason for permuting in the first place. There is no critical F for any pair, at any sample size. The recommendation to work from permutation p-values is right; the reasoning offered for it is not.
Pooling for SIMPER is contradicted by the paper's own results. Table 2 records SODC1 differing from SLF3 at p = 0.035 — both Stonehenge groups — and PO2ii differing from PO3 at p = 0.020 — both Welsh groups. By the author's own analysis, neither pool is internally homogeneous. Averaging each into a single mean profile produces two composite figures that correspond to no actual rock, and their agreement is guaranteed by the mixing rather than by any shared origin. It should also be said that SIMPER is built on Bray–Curtis dissimilarity, which is defined for abundance data; applied to element ratios it has no clear interpretation, and the paper's own observation that the high-magnitude ratios dominate the dissimilarity profile is a symptom of this rather than a finding.
8. Two errors of fact
The text states that SH65 matches CGD1 at d = 0.1585. In Table 3 the distance from SH65 to CGD1 is 1.3004. The figure 0.1585 belongs to PCM30 — a different outcrop in a different dolerite group, and more than eight times nearer than CGD1. One of the paper's two headline pairings names the wrong source.
The note beneath Table 3 states that only two or three pairings reject the null hypothesis. Table 2 contains five p-values below 0.05, two of which are within-side comparisons.
Separately, the submitted document contains material that reads as unedited drafting: passages addressed in the second person within a first-person paper, unrendered LaTeX in the body and appendix, two alternative titles both retained, and one sentence in the appendix that describes an intention to draw the reader's attention away from a discrepancy in the author's own figure. That sentence cannot be what the author meant to publish, and it should be removed.
9. What a working version would look like
The instinct behind the sample-level analysis is sound. Nearest-neighbour matching in geochemical space is a legitimate provenance technique, and preferring it to group-level tests when groups are small is a reasonable judgement. The execution is what fails. A version that would carry weight would need:
- Concentrations, centred-log-ratio transformed, rather than twelve ratios sharing a numerator.
- Standardisation before any distance is computed, so that all variables contribute and the result is not a proxy for La/Lu.
- Groups kept separate. The pairwise structure is the provenance signal; pooling destroys it. The rejections in Table 2 are the informative entries, because exclusion is what geochemistry can actually establish.
- Candidate sources outside Preseli, so that the specificity of a match can be demonstrated rather than assumed.
- A permutation null for the nearest-neighbour distances, answering whether the observed minima are smaller than would arise by chance from twenty candidate outcrops.
- The replicate distance reported as the resolution limit, with any assignment whose margin over the runner-up falls below it declared undetermined.
Applied honestly, that last step alone would leave most of the assignments in Table 3 unresolved. That is not a failure of the study; it is the correct result for this dataset, and stating it would be a genuine contribution.
10. Why this matters beyond one paper
The pattern here is not unusual. "No statistically significant difference" is quietly doing the work of a positive finding across a good deal of provenance literature, archaeological and geological alike, and the smaller the sample the more confident the claim tends to become. The framework in this paper is an unusually explicit version of a common move: replacing tests that can fail with descriptive statistics that cannot, and describing the absence of an error estimate as independence from sample size.
The replicate check offers a cheap and general guard against it. Most analytical programmes run duplicates already. Computing the distance between two analyses of the same sample, and refusing to report any assignment whose margin is smaller than that distance, costs nothing and would prevent a great deal of overclaiming. It is the one thing this paper's data do establish, and the author supplied the means to establish it.
Appendix: worked arithmetic
Reproducing d(OU11, CGD2) from untransformed Table 1 values.
| Ratio | OU11 | CGD2 | Difference | Squared |
|---|---|---|---|---|
| La/Ce | 0.382173 | 0.373181 | 0.008992 | 0.0000809 |
| La/Pr | 2.464567 | 2.409396 | 0.055171 | 0.0030438 |
| La/Nd | 0.473525 | 0.460256 | 0.013269 | 0.0001761 |
| La/Sm | 1.462617 | 1.436000 | 0.026617 | 0.0007085 |
| La/Eu | 3.556818 | 3.626263 | −0.069445 | 0.0048226 |
| La/Gd | 1.298755 | 1.291367 | 0.007388 | 0.0000546 |
| La/Tb | 7.113636 | 7.039216 | 0.074420 | 0.0055383 |
| La/Dy | 1.046823 | 1.052786 | −0.005963 | 0.0000356 |
| La/Ho | 5.396552 | 5.439394 | −0.042842 | 0.0018354 |
| La/Er | 1.956250 | 1.983425 | −0.027175 | 0.0007385 |
| La/Yb | 2.100671 | 2.124260 | −0.023589 | 0.0005564 |
| La/Lu | 14.904760 | 14.958330 | −0.053570 | 0.0028698 |
| Σ | 0.0204605 |
√0.0204605 = 0.14304, matching the reported 0.143. The distances are therefore computed on untransformed ratios, not on normalised values as the methodology states.
Decomposition of the largest distance, d(SH49, PCDL26) = 3.7325.
| Ratio | Difference | Squared | % of total |
|---|---|---|---|
| La/Lu | 3.24276 | 10.5155 | 75.5 |
| La/Tb | 1.34784 | 1.8167 | 13.0 |
| La/Ho | 0.92424 | 0.8542 | 6.1 |
| La/Eu | 0.51740 | 0.2677 | 1.9 |
| La/Yb | 0.38592 | 0.1489 | 1.1 |
| La/Er | 0.37596 | 0.1413 | 1.0 |
| remaining six | 0.1870 | 1.3 | |
| Σ | 13.9313 |
√13.9313 = 3.7325. Three variables account for 94.6% of the result.
Lower bounds on the replicate distances.
By the triangle inequality, d(A, B) ≥ |d(A, C) − d(B, C)| for any C. Scanning all columns of Table 3 for the largest such difference:
- OU10 / OU10 rpt, against PCG20: |1.092765 − 0.917049| = 0.175716, so d ≥ 0.1757
- PCM7 / PCM7 rpt, against OU15: |0.750 − 1.022| = 0.272, so d ≥ 0.2720
These are lower bounds; the true replicate distances may be larger.
No comments:
Post a Comment
Comments welcome on fresh posts - you just need a Google account to do so.