Storminess in North West Europe and volcanic activity during the Holocene

Evidence from observational records and model simulations suggest that volcanic eruptions can strengthen mid- to high-latitude atmospheric circulation and enhance westerly wind strength, with recent proxy data-model assimilations supporting this. However, assessments of Holocene variability in storminess rarely consider whether major volcanic eruptions could be a possible driver of reconstructed periods of enhanced storminess. This research presents a new reconstruction of past storminess from a coastal peatbog situated in western Ireland spanning the last ∼ 7 ka. The record is based on the measurement of the sand content along the core, with XRF core scanning analysis also applied to test whether variations in quartz sand, shell sand and sea spray can be detected by variations in silica, calcium and bromine respectively. Ca measurements were similar to the long-term changes in sand content along the core, however, peaks in sand content were not detected, while Si reflected increases in sand content only within the last millennium when the inorganic content was highest. Br concentrations appear to have been influenced primarily by humification. We also compared sand-based storminess records from northwest Europe. Six multi-decadal to centennial periods with enhanced storminess are common to records from Ireland and Wales during the last 2.5 ka BP, centred at c. 2.25, 2, 1.4, 1.1, 0.5 and 0.2 ka BP, with less agreement between records before this time. The storm periods at 2.8, 2.2–2, 1.1 and 0.5 ka BP are more widespread events and agree with records from Sweden and Scotland. Each of the episodes of increased storminess coincide roughly with major volcanic eruptions during the late Holocene, as well as with periods of enhanced North Atlantic ice-rafting. We hypothesise, therefore, that both enhanced storminess and ice-rafting may have resulted from the climate and environmental impacts of these eruptions, aligning with the findings of recent observational and modelling studies on the climate response to eruptions. Challenges remain, however, in testing this hypothesis, given chronological uncertainties in peatland records and uncertain interpretations of the factors influencing sand deposition. Therefore, to provide an independent assessment of the influence of explosive eruptions on storminess for Ireland's northeast Atlantic position, we draw upon the rich tradition of annalistic record keeping on the island, including many reports of major storms and windy seasons, to develop a windiness index running from the sixth to seventeenth centuries CE. A set of superposed epoch analyses shows that the ice-core-based dates of explosive volcanic eruptions are statistically significantly associated with the dates of documented storms and windy seasons in Ireland, suggesting avenues for future research.

1 Introduction

There have been advancements recently in our understanding of atmospheric circulation during the last millennium. A jet stream reconstruction has been created based on European tree ring records (Xu et al., 2024), while other jet stream reconstructions have been developed by the assimilation of annually resolved palaeoclimate proxies, documentary evidence, and observations with model simulations of the past climate (e.g. Sjolte et al., 2018; Osman et al., 2021; Brönnimann et al., 2025). Annually resolved North Atlantic Oscillation (NAO) reconstructions have also been developed using proxy records from the North Atlantic validated with past millennium climate simulations (Ortega et al., 2015). These studies of annual and seasonal atmospheric circulation of the last millennium include some discussion of the natural drivers of jet stream position, with volcanic activity often mentioned. While one study finds only weak evidence for volcanic influence on the North Atlantic jet stream (Osman et al., 2021), others have found evidence of jet stream strengthening following major (tropical and northern hemisphere) eruptions (Brönniman et al., 2025), extreme southward shifts following the tropical Mayon and Tambora eruptions in 1814/1815 CE (Xu et al., 2024) and positive NAO circulation after tropical eruptions (Ortega et al., 2015; Sjolte et al., 2021). A strengthening of atmospheric circulation is supported by evidence of high storminess in the NEEM-2011-S1 ice core following major volcanic eruptions (Andreasen et al., 2024). These support observational and modelling evidence that large tropical eruptions (which can inject substantial volumes of sulphur dioxide gas into the stratosphere in both hemispheres, where it oxidizes to form sulphate aerosol particles that absorb outgoing terrestrial long-wave radiation and incoming near-infrared radiation) can act to warm the tropical lower stratosphere, creating a meridional temperature gradient that promotes a stronger circumpolar vortex circulation, positive NAO anomalies, enhanced westerlies and wetter conditions in northern Europe in successive years (Robock and Mao, 1992; Shindell et al., 2004; Fischer et al., 2007; Stenchikov et al., 2006; Kodera, 1994; Zambri et al., 2017; Sjolte et al., 2021). A more recent modelling study has shown eruptions are associated with reduced storm frequency in the mid-latitudes, with increases instead in the subtropics and high latitudes (Andreasen et al., 2024).

Our understanding of longer-term (over multiple millennia) variations in storminess is primarily based on evidence from sedimentary and organic deposits. Storminess records of decadal to centennial temporal resolution can be derived from coastal depositional settings, such as lagoons, lakes and peatlands, where variations in the sediment characteristics (e.g. sand content, grain size) reflect the energy and/or frequency of past storms (e.g. Sabatier et al., 2012; Sorrel et al., 2012; Orme et al., 2015, 2016, 2017; Goslin et al., 2018; Vandel et al., 2019; Kylander et al., 2020; Sjöström et al., 2024). These records and others have advanced our knowledge of past storminess in Europe, although one challenge that remains is explaining the origins of observed differences between such reconstructions.

Several authors have compiled multiple reconstructions of storminess to identify common signals in northwest Europe. Based on multiple records from different proxies (i.e. from sand dunes, estuaries and peatlands), Sorrel et al. (2012) found centennial periods of enhanced storminess at 5.8–5.5, 4.5–4, 3.3–2.4, 1.9–1.1 and 0.6–0.3 ka BP. Kylander et al. (2023) created regional storm stacks by combining peatland-based storm records from southern Sweden and western Scotland. The southern Sweden storm stack shows high storminess centred at approximately 4.5–4.3, 3.9–3.8, 2.87, 2.3–2, 1.2–1.1 and 0.7–0.4 ka BP, while the west Scotland storm stack shows high storminess centred at 3.4, 2.7–2.6, 2.3–2 and 0.9–0.4 ka BP (Kylander et al., 2023). There are also other periods of high storminess identified from single sites in northwest Europe, including from Wales (Orme et al., 2015), western Denmark (Goslin et al., 2018), Estonia (Vaasma et al., 2025) and Ireland (Sjöström et al., 2024). While there are common periods of higher storminess identified in several of these records, including at 2.8–2.7, 1.1 and 0.5 ka BP, there are also differences. The differences may result from local anthropogenic disturbances, the influence of other climate and environmental variables (e.g. snow cover, vegetation cover), site geography (e.g. distance/position of the sand sources in relation to the bogs) and the proxies used to measure sand content. On the other hand, they may reflect real differences in storminess between regions caused by shifts in the storm track or pressure centres (Orme et al., 2016; Nielsen et al., 2024; Kylander et al., 2023; Sjöström et al., 2024).

Several authors have found higher storminess at times when ice-rafting in the North Atlantic was more extensive, which occurred at c. 6.3–5.1, 4.7–3.9, 3.4–2.7, 1.5–1.1 and 0.6–0.1 ka BP (Bond et al., 2001; Sorrel et al., 2012; Orme et al., 2015; Kylander et al., 2020; Goslin et al., 2018; Sjöström et al., 2024). However, the explanation for this varies between studies with uncertainty about cause and effect. Some have suggested that high ice-rafting reflects changes in the North Atlantic latitudinal temperature gradient, which may have enhanced past storm intensity (Lamb, 1995; Raible et al., 2007; Trouet et al., 2012; Sorrel et al., 2012). Sjöström et al. (2024) concluded that a shared driving mechanism could explain the similarity of oscillations in storminess and IRD, while Goslin et al. (2018) has suggested that enhanced storminess caused the changes in the ice rafting, by modulating the strength of the subpolar gyre and ice export from the Arctic. Solar variability has been discussed as a possible driver of the changes in storminess and ice rafting (e.g. Martin-Puertas et al., 2012; Goslin et al., 2018; Kylander et al., 2020), although other work has found no evidence of solar influence, instead positing that internal oceanic oscillations are a key driver of ocean and atmospheric circulation (e.g. Sorrel et al., 2012). While these studies have considered solar and oceanic drivers of past storminess, the role of volcanic eruptions has received little attention, most likely due to the decadal-centennial resolution of available storm records.

Figure 1 Vellum page (H 42r b) from TCD MS1282, the Annals of Ulster, covering the years 852 to 858 CE. The top and bottom insets present the original language text (Irish and Latin) plus English translations of entries reporting severe cold and storminess in consecutive years, shortly after the 852/3 CE eruption of Mount Churchill, Alaska (Mackay et al., 2022). Image reproduced with permission from the Board of Trinity College Dublin. Text and translation adapted from Mac Airt and Mac Niocaill (1983).

The primary aim of this study is to add to our knowledge about Holocene storminess in northwest Europe, through a new storm record from western Ireland and a synthesis of records from the region. This will enable further identification of particularly stormy periods in recent millennia in northwest Europe and discussion of past drivers, including the possible role of volcanic eruptions. As an independent source to assess the potential eruption-storminess linkage in Ireland, we also exploit Ireland's rich written tradition of medieval annalistic record keeping that stretches from the early medieval to the early modern periods (c. 1.6–0.4 ka BP). These records comprise annually arranged lists of events deemed of importance by their authors (Fig. 1). Up to the thirteenth century, most recording took place in major monastic centres across Ireland, such as Clonmacnoise on the bank of the River Shannon in the Irish midlands. Thereafter, the disruption to Irish monasticism initiated by the Anglo-Norman invasion and colonisation in the twelfth and thirteenth centuries promoted a move to recording further west and north in Gaelic dominated areas where the learned families of historians in the employment of the Gaelic aristocracy continued the annalistic tradition (McCarthy, 2008). Of all natural meteorological extremes experienced in Ireland, major windstorms are among the most feared, given the island's exposure to deep low-pressure zones and associated cyclones swept across the island by the (particularly winter season) westerlies (Sweeney, 1997, 2000). Medieval Irish society was notably vulnerable (agriculturally and infrastructurally) to such hazards (Ó Corráin, 2005; Campbell and Ludlow, 2020), as the following report from the Annals of Ulster for 804 CE attests: “Violent thunder, accompanied by wind and fire, on the night before [i.e. 16 March, Julian Calendar] St. Patrick's Day, which destroyed many persons, i.e. one thousand and ten in Corcu Baiscinn [Thomond, Co. Clare]; and the sea divided the island of Fita [Mutton Island, Co. Clare] into three parts, and covered the land of Fita with sand ….” (Mac Airt and Mac Niocaill, 1983). Such events could also clearly mobilize considerable volumes of sediment. This vulnerability, alongside an interest in the incidence of spectacular and damaging weather extremes as potential divine portents or instances of divine retribution for social misbehaviours (McCarthy and Breen, 1997; Williams, 2010; Baker et al., 2017), has ensured their abundant reporting in the medieval Irish records.

A secondary aim of this paper is to assess whether the results from X-Ray Fluorescence core scanning (XRF-CS) analysis provides useful information about past storminess, either through sand deposition, with silica as a proxy for quartz sand and calcium carbonate for shell sand, or indicators of sea spray such as bromine. XRF core scanners rapidly measure relative elemental concentrations along sediment cores at up to sub-millimetre resolution (Croudace et al., 2006). An XRF beam is focused at the sediment surface and then the fluorescent X-rays generated are measured as “counts” for each element (Croudace et al., 2006). There are limitations, however, to the method that need to be accounted for. The counts of individual elements can be influenced by changes in the organic and water contents and variations in other elements (the dilution effect and closed sum effect, respectively; Löwemark et al., 2011). Furthermore, the XRF-CS struggles to detect elements with low atomic numbers (Gregory et al., 2019), and counts are also influenced by uneven core surfaces and changes in water content and grain size (Croudace et al., 2006; Weltje and Tjallingii, 2008). These issues can be mitigated by following recommended data processing techniques, such as removing data that has high Mean Squared Error (MSE) and large deviations in the total counts (Löwemark et al., 2019), normalising against a conservative element using a log-ratio (Löwemark et al., 2011; Weltje and Tjallingii, 2008; Croudace et al., 2019) and a centred log-ratio approach (Bertrand et al., 2024).

Comparisons between XRF-CS and ICP-OES/ICP-MS measurements along peat cores have shown that both Ca and Si can be measured effectively by core scanners (Poto et al., 2015; Kern et al., 2019), although this is likely to be sample specific, where sand concentrations are abundant enough to exceed the detection limits of the core scanner. While XRF-CS measurements of Si and Ca along peat cores have not been widely used to reconstruct past storminess, Si has been used as a proxy for dust deposition to peatlands (Lim et al., 2015; Longman et al., 2017). Bromine has been measured on cores from coastal peatlands to provide a proxy for sea spray deposition and therefore past storminess (Turner et al., 2014; Orme et al., 2015; Stewart et al., 2017), with instrumental measurements demonstrating that there is an exponential relationship between wind speed and the deposition of sea spray aerosols across coastal regions (Franzén, 1990; Gustafsson and Franzén, 1996; Meira et al., 2008). However, there remain questions about this storminess proxy, as it is thought that <5 % of deposited bromine is retained within peat (Shotyk, 1997). Furthermore, studies have found that over 90 % of bromine in peat is organically bound, meaning that the degree of peat decomposition may alter the concentration of bromine along the core, with higher (lower) peat humification leading to enhanced (reduced) concentrations of bromine (Biester et al., 2004, 2012; Martínez Cortizas et al., 2007; Kylander et al., 2020). In this research, comparison of the XRF-CS results with measured sand concentrations and humification will enable an assessment of their potential as proxies for storminess research.

Figure 2 Map showing the location of Bullsmouth Bog on Achill Island on the west coast of Ireland and locations for other regional storm reconstructions. Map (A) shows the British Isles including the location of Achill Island (orange square) and Cors Fochno bog in Wales (CF; Orme et al., 2015; red square). Map (B) shows Achill Island, including the location of Bullsmouth Bog (orange square), nearby Roycarter Bog (RC; Sjöström et al., 2024; red square), the islands elevation (thin and thick contours at 50 and 200 m intervals respectively and peak elevations in m above sea level), water bodies (blue shading), beaches (yellow shading) and forests (green shading). Map (C) shows the main features of the area around Bullsmouth Bog and the core site location (red square).

2 Study Area

Achill Island is located on the west coast of County Mayo, Ireland (Fig. 2). The nearest meteorological station record from Belmullet, 25 km north of Achill Island, shows that between 1981 and 2010 the mean annual temperature was 10.3 °C (annual range from 6.3 °C in January to 15 °C in August), the mean annual precipitation was 1245 mm (mean monthly totals ranging from 70 to 146 mm) and the annual mean wind speed was 6.6 m s⁻¹ (ranging from 5.7 m s⁻¹ in July to 7.9 m s⁻¹ in January), with almost 30 days with gales each year (Met Éireann, 2020a). The prevailing winds are from a west and southwesterly direction (Walsh, 2012). The west coast of Ireland typically experiences multiple severe storms each year, for example there are on average eight named storms each year (Met Éireann, 2025b). Storm Éowyn, in January 2025, was the strongest storm on record for Ireland; at Belmullet station there were sustained (10 min) wind speeds of 91 km h⁻¹ (25 m s⁻¹; Storm Force 10) and a maximum gust of 137 km h⁻¹ (38 m s⁻¹) (Met Éireann, 2025c).

Achill Island has mountains reaching ∼ 670 m above sea level, although much of the island is low lying and relatively flat at <100 m elevation (Fig. 2B). The geology of Achill island is part of the precambrian Dalradian Supergroup, a siliciclastic metasedimentary succession, which includes quartzites, gneisses and schists (Chew, 2003; Prave et al., 2009). Along the south western side of Achill island there are beaches including the 3 km long Keel Beach, which is a large sandy beach with pebble storm beach, sand dunes and a rare machair environment (DAHG, 2013). The lower slopes and valleys on the island are widely covered by Atlantic blanket bog, some of which have been previously cut for fuel, as evident from satellite images. Bullsmouth Bog is a small area of blanket bog (∼ 0.2 km²) on the eastern side of the island close to the village of Bullsmouth (Fig. 2). While there are houses, roads and drainage ditches around the edges of the bog (Fig. 2C), satellite images show no signs of disturbance across the central parts, across an area approximately 250 m in diameter. The peat in the centre of the bog is ∼ 4 m thick. The closest beach to Bullsmouth Bog is 400 m to the northeast of the core site, where a narrow sandy beach and ridge of pebbles and stones extend along the coastline for approximately 700 m. The beach and machair at Keel is 7.9 km to the southwest, while to the northwest there are sandy beaches at a distance of 3.8 km.

3 Materials and Methods

3.1 Fieldwork

Two cores (BM-1 and BM-2) were extracted 1 m apart from the centre of Bullsmouth Bog, Achill Island, in 2019, from a longitude and latitude of 53°59.791^′° N and 9°56.352^′° W and at an elevation of 5.5 m above sea level. The cores were taken using a Russian corer with a barrel length of 1 m and width of 5 cm. The cores were sampled from the surface to the bottom of the bog at 4 m depth, although the bottom 25 cm of BM-2 was not retained successfully by the corer. The cores were sampled in 1 m contiguous sections, which were transferred into guttering and wrapped in plastic film. The cores were stored in a cool room at Maynooth University prior to analysis.

3.2 Chronology

The core chronology was constrained by accelerator mass spectrometry radiocarbon dates from six bulk peat samples taken from BM-1, analysed by Beta Analytic Inc., Florida. Bulk peat was sampled as previous research has shown no significant difference between dated bulk peat and macrofossils (Blaauw et al., 2004; Holmquist et al., 2016). The ¹⁴C ages were calibrated and the age-depth model constructed in Bacon version 2.5.7. (Blaauw and Christen, 2011) using the IntCal20 calibration dataset (Reimer et al., 2020). Given the close sampling proximity of the two cores it is likely that the accumulation rates of the peat will have been similar, although differences may have occurred due to microtopography (e.g. Foster et al., 1988; Harris et al., 2020). Following this consideration, a single chronology was developed for core BM-1 and tentatively applied to both cores.

3.3 Loss-on-ignition

Loss-on-ignition (LOI) analysis was applied to core BM-1 (Dean, 1974; Heiri et al., 2001) to establish the inorganic and organic content of the peat. Sampling was conducted at contigous 1 cm intervals along core BM-1, at each depth using 5 cm³ of wet peat. The dry weight of the samples was measured by drying in an oven overnight at 102 °C, before cooling in a dessicator and weighing. The weight of organic material in each sample was measured by ignition in a furnace at 550 °C for 4 h, before again being cooled and weighed to determine the mass loss relative to the dry sample, which reflects the organic content.

Variability in sand content in coastal peatbogs can be used as a proxy for storminess (e.g. Björck and Clemmensen, 2004; De Jong et al., 2006; Orme et al., 2015, 2016), as supported by comparisons between peatland-based Aeolian Sediment Influx (ASI) reconstructions and instrumental storm records (Vandel et al., 2019). The ignition residue (IR) reflects the inorganic content of the peat sample, and was calculated using the weight of the ashed sample as a percentage of the dry weight. Previous research on proxies for sand content in coastal peatbogs supports the conclusion that the IR predominantly reflects changes in sand content, with strong similarities between the measured IR and sand grain counts (Orme et al., 2015) and sand fraction weights (Orme et al., 2016).

3.4 XRF Core Scanner Analysis

Core BM-2 was scanned with an ITRAX X-ray fluorescence (XRF) core scanner (Croudace et al., 2006) at University College Dublin. Non-destructive element analysis was conducted using a molybdenum (Mo) anode X-ray tube (settings: 30 kV, 45 mA, dwell time 10 s) at 0.5 mm increments along the core.

XRF core scanning results are semi-quantitative and can be influenced by the dilution and closed-sum effect, as detailed above (Löwemark et al., 2011) as well as uneven core surfaces and changes in water content and grainsize (Croudace et al., 2006; Weltje and Tjallingii, 2008). To address these limitations, sections of the XRF-CS record with anomalous peaks in the Mean Squared Error (>2.5 and <1.5) and anomalous changes in the total “counts per second” (kcps) were removed, as these are indicators of inaccurate measurements that can introduce artefacts to the data (Löwemark et al., 2019). The results were normalised using a centred log-ratio approach, which is a multivariate log-ratio transformation of the data that removes issues resulting from sediment properties and non-linear matrix effects (Bertrand et al., 2024). The individual element counts were divided by the mean count of all well detected elements (those without zeros) and expressed as a centred log ratio (Bertrand et al., 2024). 1 cm moving average records were calculated from the centred log ratio (clr) element records.

3.5 Humification

The degree of humification was measured on 1 cm layers at 5 cm intervals on core BM-2, following the method of Chambers et al. (2011). This involved simmering dried and ground peat samples in 8 % NaOH for an hour before filtering them through Whatmann number 1 grade filter paper. Three liquid subsamples taken from each sample were then analysed on a Hitachi U-1100 spectrophotometer and the percentage light transmission and absorbance recorded, with averages calculated for each sample.

3.6 Statistical Analysis

Pearson correlation coefficients were calculated between the centred log-ratio elements, the IR and humification (average percentage light transmission) results to test their similarity. The records are at different resolutions and therefore the element results were downsampled to the resolutions of the IR and humification records prior to calculation of the correlations.

Principle Component Analysis (PCA) was conducted on the centred log-ratio element results. The transformation of the element data removes issues associated with the constant-sum constraint in raw element data, whereby changes in one element cause changes in others (Aitchison, 1983; Kucera and Malmgren, 1998). PCA was performed in R using the PCA function from the “FactoMineR” package (Lê et al., 2008) and visualised using the fviz_pca_var function of the “factoextra” package (Kassambara and Mundt, 2020).

3.7 Historical Windiness Index

A survey was conducted of the surviving corpus of medieval Irish annalistic texts (McCarthy, 2008), collectively known as the “Irish Annals”, with the purpose of identifying all references to storms that likely arose from large-scale Atlantic cyclones (excluding more transient thunder and hailstorms arising solely from mesoscale weather systems), as well as references to years described as notably windy or unsettled. The survey included 22 surviving texts amounting to approximately 1.15 million words, spanning the fifth to seventeenth centuries CE (c. 1.6–0.4 ka BP). This identified 99 years credibly experiencing either one or more large individual storms, or prolonged windy and unsettled conditions. These were combined into a simple storminess index in binary 0.1 format (i.e. a year either has credible documentary evidence of such conditions, or it does not). A date range of 600 to 1616 CE was chosen for further analysis, with the loss of one storm event, dated 563 CE. This range was chosen in recognition of the changing density of reporting in the surviving texts, with at least one event (of any type) being reported effectively each year from 600 to 1616 CE (Ludlow, 2012). The former date represents a convenient cut-off, with the preceding years having increasingly thinner (i.e. lower density) coverage, including outright multi-year gaps, whilst the latter date represents the cessation of coverage of the last major annalistic text, the Annals of the Four Masters (McCarthy, 2008). Within the 600–1616 range of effectively continuous reporting, any years experiencing relevant weather have, in principle, the possibility of being reported, and hence there is likely to be a meaningful climatic signal in both the years assigned a 0 (not windy) and a 1 (windy), important for later statistical testing (Sect. 4.5).

The above does not mean that the record of windiness is complete, even within the 600–1616 range. The likelihood of any given storm being included in the surviving record will vary according to the subjective experience and interests of the original authors (and later copyists), perhaps relative to the magnitude of the event in question (Ludlow, 2010). It has also been shown that a positive correlation exists between years with greater coverage (i.e., more texts available and more density of coverage per text) and the number of weather reports available (Ludlow, 2010, 2012). In other words, some unknown number of years assigned a value of 0 within the 600–1616 range will in fact have experienced notable storminess but will be absent from the annalistic record because of a (temporarily) lower interest in recording storms below a certain perceived severity or because of a lower density of coverage for the year in question (either because of genuinely less active recording or because of loss in transmission through time, as detailed in McCarthy, 2008). However, we posit that because (1) the relative rarity of notably stormy years and (2) the great interest in the reporting of storms and adverse weather in the surviving texts (Ludlow and Travis, 2019 count 313 reported meteorological and related phenomena and conditions from 431 to 1649 CE, averaging one every 3.9 years), a year assigned 0 for windiness within the 600–1616 range will more often reflect such conditions than not. By contrast, years assigned a 1 can be understood to more confidently reflect notably stormy conditions, although again imperfectly, given subjectivity on the part of the authors and some variability in the spatial extent of storms experienced across the island. For example, a major storm impacting one region of the island may not have been experienced as sufficiently severe to merit reporting in another region, although the texts often do report events beyond their specific centres and regions of recording, likely drawing information from networks of parent and daughter monastic foundations across the island and from the hosting of visiting ecclesiastics, dignitaries and pilgrims; Ludlow, 2010).

Similarly, some years assigned a 1 in our index may be included based on misleading, fabricated or exaggerated reporting, but this can be considered rare. The field of historical climatology, which specialises in the use of written evidence for palaeoclimatological reconstruction and human impact studies, emphasizes the importance of applying historical source criticism (e.g. MacNeil, 2000) to assess the reliability of historical sources, a process that begins by developing an understanding of the context in which the evidence was created (Brázdil et al., 2005). This includes a consideration of the worldview of the authors, the intended audience and purpose of its creation, necessary to properly interpret the evidence and assess potential biases. The abundant meteorological content of the Irish Annals has previously been assessed in this way and most reports shown to be historically credible (on a three-point scale of apparent reliability, 74.5 % of all relevant reports were scored by Ludlow (2010) at the highest (most reliable) point of the scale). Years assigned a value of 1 in the windiness index constructed here exclude reports deemed apparently unreliable (lowest point on the scale) by Ludlow (2010). Comparison to independent evidence such as tree-ring and ice-core records has also pointed to the general reliability of the annalistic weather evidence (e.g. McCarthy, 2008; Ludlow et al., 2013; Gao et al., 2016; Ludlow and Travis, 2019; Kostick and Ludlow, 2022). Moreover, when applying chronological corrections summarized by McCarthy (2008, 2017, 2018) and verifiable by the comparison of eclipses recorded in the Irish Annals to their independently retro-calculated dates (for which see also McCarthy and Breen, 1997), the chronology of the Irish Annals can be shown to be accurate to within a very small uncertainty (at most ±1 year) for most of the span covered by these sources.

To independently assess the credibility of a possible association between explosive volcanism and storminess in Ireland's northeast Atlantic location, a selection of superposed epoch analyses (SEAs) was conducted. This is a compositing approach that examines whether there is a coherent response (e.g. increase or decrease) in the average frequency of a potential response variable (here the frequency of years experiencing storms and windiness as captured in our windiness index) relative in time to a set of dated potential trigger or causal events (here the polar ice-core-based dates of major explosive eruptions; Sigl et al., 2015). The approach is widely used in palaeoclimatic studies and allows the testing of the likely randomness (statistical significance or confidence), scale and temporal signature of any potential response (i.e. lagged or multi-year responses) (e.g. Sigl et al., 2015; Manning et al., 2017; Rao et al., 2019; Gao et al., 2021). To test whether the value in any observed year is elevated beyond that which might be expected by chance alone (i.e., statistically significant), we adopted a Monte Carlo randomisation and resampling approach, following that of Gao et al. (2021).

4 Results

4.1 Chronology

The results of the AMS ¹⁴C dating of core BM-1 are shown in Table 1. The developed age-depth model and 2σ age uncertainty for the Bullsmouth Bog core (Fig. 3) indicates that the record spans the interval from 6.41 to 0.03 ka BP, with an average temporal resolution (at 1 cm sampling resolution) of 17 years. The BM-2 core was not dated directly and the bottom 25 cm of this core was not retained during the coring process. Assuming that the sediment accumulation history of BM-1 and BM-2 has been broadly consistent, the humification and XRF-CS records developed on BM-2 span the interval from 5.9 to 0.03 ka BP.

Table 1 Radiocarbon dates and calibrated ages for samples from core BM-1.

Figure 3 Stratigraphy and age-depth model measured on core BM-1, including the calibrated radiocarbon dates and 2σ errors (orange points and error bars), modelled median age (black line) and modelled two-sigma age uncertainty along the core (blue shading).

Figure 4 Ignition residue (BM-1 core), humification (BM-2 core) and clr-normalised element (BM-2 core) results from Bullsmouth Bog. The element results show the high resolution results measured at 0.5 mm (light blue line) and the 1 cm moving average results (black line).

4.2 Loss-on-ignition

The reconstructed mineral influx to Bullsmouth Bog is shown by variability in the IR (Fig. 4; Orme, 2026). The record has low values (2 %–3 %) and a slight increasing trend between ∼ 6.4 and 2.1 ka BP with some peaks in the IR at 3.4, 3.2 and 2.4 ka BP. After 2.1 ka BP the IR steadily increases until 1.1 ka BP (from 2.6 % to 4.1 %) with a reduction in the IR between 1.1 and 0.8 ka BP. Between 0.8 to 0.65 ka BP there is a large peak in IR (reaching a maximum of 9.5 %) after which it remains high (5 %–6 %).

4.3 Humification

The average light transmission record is a proxy for peat humification, with higher values indicating low peat humification, and lower values indicating well-humified peat. Values range from 21 %–47  % through most of the record, with higher values at c. 5.5, 5.1–4.8, 4.6, 4.4–4.1, 3.9, 3.5–3.1, 2.2, 1.4 and 1.1 ka BP and lower values at c. 6–5.6, 3.9–3.5, 3.1–1.8, 0.8–0.2 ka BP (Fig. 4; Orme, 2026). After 2.4 ka BP there is an increasing trend in the results. The peat in the upper 35 cm (representing the period since 0.17 ka BP) has higher values of 50 %–60 %, reflecting the partially decomposed peat of the acrotelm.

4.4 XRF Core Scanner Analysis

XRF-CS data is semi-quantitative as the measure of different element counts can be influenced by the organic and water content of the peat, surface roughness, the abundance of other elements and the sediment grainsize (Croudace et al., 2006; Weltje and Tjallingii, 2008; Löwemark et al., 2011). To address some of these issues the raw XRF-CS results (Supplement Fig. S1) were processed to remove sections of the core where the scan results were compromised, as indicated by anomalous changes in the MSE and total kcps (Fig. S2; Löwemark et al., 2019). The sections that contained significantly erroneous data were the ends of the core sections, therefore all data from these sections was removed from further analysis as the resulting data would likely contain artefacts. The centred-log ratio element results are shown in Fig. 4 (see also Orme, 2026), including the moving average calculated with a 1 cm window.

The clr-Si record (Fig. 4) has small variations along much of the core, with small peaks at approximately 2.6, 1, 0.6 and 0.2 ka BP. The clr-Ca record shows higher variability, including low values at c.4.7 ka BP, an increasing trend to 3.1 ka BP, a decreasing trend to 0.9 ka BP and increasing values from 0.9 to 0.2 ka BP. The clr-Br record shows an increasing trend through the record. There are higher values at approximately 5.9, 5.5, 4.9, 4.4–4.2, 4.1, 3.7, 3.4, 3.2, 2.3, 1.5, 0.8 and 0.2 ka BP and lower values at approximately 5.6, 5.2, 4.5, 3.9, 2.9–2.3 and 0.7–0.4 ka BP.

4.5 Statistical Analysis

The correlations between the elements, IR and humification results are shown in Table 2. Considering the element to element correlations, the strongest correlations are between clr-Fe and clr-Cu ($R=-\mathrm{0.37}$), clr-Ca and clr-Br (R=0.31) and clr-S and clr-Ba ($R=-\mathrm{0.3}$). The IR results correlated most strongly with the clr-Br (R=0.48) and clr-Ca (R=0.39) results. Humification correlated most strongly with clr-Br (R=0.46) and clr-Zr ($R=-\mathrm{0.42}$). It should be noted that while the element and humification analysis were conducted on BM-2, the IR analysis was conducted on BM-1 which, given potential variations in accumulation rates between the two cores, may influence some of these results.

The PCA results are shown in Fig. 5. The first and second eigenvectors each explained 15 %–16 % of the variability. PC1 is associated with positive loadings for clr-Zr, clr-S, clr-Al and clr-Si and negative loadings for clr-Cu, clr-Ba, clr-Ar and clr-Fe. PC2 has positive loadings for clr-Br, clr-Ca and clr-Cu and negative loadings for clr-Fe.

Table 2 Correlation matrix of normalised elements (clr-element) and the IR and humification results. The sample number (n) for the element-element correlations is 6520. For correlations with the IR n=340 and with the humification results n=66.

Figure 5 PCA results for the centred-log ratio element data.

4.6 Documentary evidence

To begin the superposed epoch analyses, we selected all 25 eruptions with a Greenland sulphate deposition signal of at least 20 kg km⁻² dated to between 600 and 1616 CE (Sigl et al., 2015). This choice was informed by the Sigl et al. (2015) Greenland sulphate deposition estimate of 21.4 kg km⁻² for the 1991 (June) eruption of Mount Pinatubo in the Philippines, after which evidence was observed for an enhanced westerly flow over Northern Hemisphere mid-latitudes (e.g. Robock and Mao, 1992; Robock, 2000; Paik et al., 2023), a general effect reproduced in several modelling efforts (e.g., Zambri and Robock, 2016; Zambri et al., 2017). Figure 6a thus shows the mean windiness index value in each of the five years before these 25 eruption dates (i.e. superposed years −5 to −1, respectively, on the horizontal axis), the mean value during these eruption years (superposed year 0), and the mean in each of the following 5 years (i.e. superposed years 1 to 5). The two highest observed values are seen at superposed years −1 and 1, representing the set of first years preceding and the set of first years following our volcanic dates, respectively, with twenty percent of years in each set reported by contemporaries as being notably stormy (i.e. registering in the windiness index with a value of 1). Superposed year 0 also has suggestively, if not uniquely, elevated reporting, with 12 % of such years registering.

To test whether the value in any observed year is statistically significant a Monte Carlo randomisation and resampling approach was applied (Gao et al., 2021). In this, our 25 eruption dates were randomised 10 000 times, as were our windiness index values, with the mean frequency of windiness index values being re-calculated for each superposed year after each randomisation. This allowed us to create a randomised reference distribution against which to compare our actually observed values. This exercise revealed the mean frequencies in superposed years −1 and 1 to be statistically significantly elevated at 94.3 % (i.e. p=0.057) and 94.5 % (i.e. p=0.055) confidence, respectively. Superposed year 0 did not meet our 90 % confidence threshold. The elevated frequencies at superposed years −1 and 1 are consistent with an association between substantial explosive eruptions and increased storminess for Ireland when considering the small remaining dating uncertainties in the ice-core timescale over our period of interest (approximately ±1 year; Sigl et al., 2015), plus a similarly small uncertainty within the corrected Irish annalistic chronologies (McCarthy et al., 2008, 2017, 2018), as well as potentially variable lag-times in the deposition of sulphate onto the Greenland ice (such that eruptions may occasionally have occurred in the calendar year before their deposition-based dates) and multi-year or lagged impacts (such that a climatic response may be expressed, or continue to be expressed, in one or more post-eruption years).

We note that doubt has recently been cast on the causality (and even reality) of an enhanced extratropical Northern Hemispheric westerly flow following the 1991 Pinatubo eruption, and large eruptions more generally, including on the basis that internal climate variability may overwhelm a volcanically forced climate response that might otherwise promote such an enhancement (e.g. Polvani et al., 2019; Tejedor et al., 2024). We therefore also examined whether the association between windiness index values and eruptions changed notably depending upon their climate forcing potential. To do so, we repeated our analyses, varying the number of chosen eruption dates according to their Greenland sulphate deposition values, ranging from the inclusion of all identifiable volcanic dates irrespective of deposition volume (n=85) up to a highly restrictive selection of only those eruptions with at least 60 kg km⁻² deposition (n=4; Sigl et al., 2015). Figure 6b shows the SEA output for all variants (11 total). While some variability is evident, there is a coherence in response captured by the mean and median (of all variants combined), in which superposed years −1, 0 and 1 again exhibited the most elevated windiness values, on average, from all superposed years examined (−5 through 0 to 5). Of the individual variant eruption lists, superposed year −1 exhibited the highest mean windiness index values, with one variant (comprising eruptions with at least 40 kg km⁻² deposition) achieving a mean of 0.375 (i.e. 37.5 % of these eruptions registered an association with windy years) at 98.12 % confidence (i.e. p=0.0188).

Figure 6 Superposed epoch analyses (SEAs) assessing the temporal association between windiness index values (derived from reporting in the Irish annals) and ice-core-based volcanic eruption dates (Sigl et al., 2015). (a) Initial SEA showing the mean windiness index value in each of the 5 years before (−5 to −1 on x-axis), the years of (0 on x-axis), and each of the five years following (1 to 5 on x-axis) 25 eruption dates between 600 and 1616 CE with a minimum Greenland sulphate deposition of 20 kg km⁻² (Sigl et al., 2015). (b) Set of 11 SEAs varying according to the number of eruption dates included as based on differing minimum Greenland sulphate deposition volumes (grey lines). Black (solid) and red (dashed) lines represent the mean and median, respectively, of all variants. (c)–(e) Mean windiness index values and corresponding statistical (Monte Carlo-based) confidence level in superposed years −1 (c), 0 (d) and 1 (f) relative to each variant set of eruption dates, as based on Greenland deposition volumes (Sigl et al., 2015).

Figure 6c–e focus on the three central superposed years (−1, 0, 1) and show the mean windiness frequency for each variant eruption list, ranked by increasing minimum deposition volume, along with the associated statistical significance. It was observed that the mean frequency of stormy years was low and statistically insignificant when including all identifiable volcanic events (i.e. imposing no minimum sulphate deposition requirement). Thereafter, for superposed years −1 and 1, the mean was seen to broadly increase as eruptions with smaller Greenland deposition values were iteratively excluded. For superposed year −1 (Fig. 6c), the association was deemed statistically significant at better than 90 % confidence (p<0.10) once any eruptions with <10 kg km⁻² deposition were excluded (leaving n=49). Statistical significance then persisted for seven further eruption lists (i.e. progressing through eruptions with increasingly greater Northern Hemispheric climate forcing potential as reflected in increasing Greenland sulphate deposition). A similar pattern was observed for superposed year 1 (Fig. 6e), although the 90 % confidence threshold was not breached until all eruptions with <20 kg km⁻² deposition were excluded (leaving n=25), and statistical significance was less persistent. Superposed year 0 never formally breached the 90 % threshold, coming closest when eruptions with less than 25 kg km⁻² deposition were excluded (80.6 % confidence, p=0.194). Ultimately, as the minimum deposition threshold was further increased and the sample size of eruptions in each list became very small, the mean windiness index value for each superposed year was observed to fall (including to below the 90 % confidence threshold for superposed years −1 and 1).

Figure 7 Comparison of potential element proxies for storminess with the ignition residue and humification results. Clr-Si (A, D, G), clr-Ca (B, E, H) and clr-Br (C, F, I) results are shown normalised (light blue) and smoothed with a moving average window of 1 cm (black line). The purple lines show the normalised IR results and the orange lines in plots (C), (F) and (I) show the normalised humification for comparison with the bromine proxy. Plots (A)–(C) present the results for 6–0 ka BP, plots (D)–(F) show the results for 1–0 ka BP and (G)–(I) show the results for 4–2 ka BP.

5 Discussion

5.1 Assessment of storminess proxies

In this section the potential storm proxies from XRF-CS analysis will be discussed, including an assessment of whether Si and Ca are potential indicators of sand, and if Br is a potential indicator of sea spray. These elements are of interest as the local sand is rich in quartz and shell fragments, with bromine a major component of sea spray. The correlation matrix shows that Br and Ca correlate more strongly with the IR storminess proxy than other elements (Table 2), warranting closer examination. The normalised element, IR and humification results are compared in Fig. 7, with a focus on the periods 1–0 and 4–2 ka BP to aid comparison.

The clr-Si record is rather noisy without large variations along the core. Through much of the record changes in clr-Si do not mirror the peaks observed in the IR record (Fig. 7), with which it has a weak negative correlation (Table 2). On the other hand, within the last millennium there is a greater similarity between the smoothed clr-Si reconstruction and the IR record, with moderate increases in both the clr-Si and IR records at c.0.6 and 0.2 ka BP (Fig. 7D), highlighting the fact that sand was detected in this section of the core where the mineral content of the peat exceeded 5 %–6 %. Elements with low atomic numbers and/or low concentrations are less well measured by the XRF core scanners due to “noise” from the matrix (Gregory et al., 2019). Therefore, the low concentrations of sand within the peat prior to 1–2 ka BP, as well as somewhat poor detection of Si (atomic number 14), can explain why in this study Si was not a suitable storm proxy.

The clr-Ca variations are similar to the larger magnitude variations in the IR record particularly during the last 2 ka, with the normalised comparison in Fig. 7 showing that both the clr-Ca and IR results show lower values at ∼ 1 ka BP followed by higher values from 0.7–0.2 ka BP (Fig. 7B and E). A significant correlation of R=0.39 between the clr-Ca and IR records supports the contention that Ca is detecting variations in the sand content to some extent and thus is a potential storm proxy. On the other hand, clr-Ca does not show increased values coinciding with most of the IR peaks, for example between 4–2 ka BP there are a number of low magnitude peaks in the IR record that are not clearly aligned to peaks in the clr-Ca record, which instead shows elevated but less variable values (Fig. 7H). These differences may in part arise because the IR and element analyses were measured on BM-1 and BM-2, respectively; slight differences in depths between the two cores may have occurred if there were variations in deposition rates across the site (i.e. caused by microtopography and the evolution of pools and hummocks). The large magnitude variations shown by the clr-Ca records may be reflecting environmental changes rather than large storms, for example fine Ca-rich silt may be transported by lighter winds following deforestation or sand dune mobilisation. Bullsmouth Bog is located 8 km from an area of machair at Keel, which is a highly calcerous ecosystem containing high concentrations of shell sand (NPWS, 2018, “Keel Machair/Manaun Cliff SAC”), therefore we speculate that the variations in Ca may reflect changes in this source environment.

Br has been interpreted as a proxy for storminess previously (e.g. Turner et al., 2014; Orme et al., 2015; Stewart et al., 2017). However, there remain questions about this storminess proxy, particularly because changes in the rates of peat decomposition may affect the concentration of Br along the core (Biester et al., 2004, 2012; Martínez Cortizas et al., 2007). The comparison of clr-Br and humification along the Bullsmouth Bog core supports the contention that Br concentrations are altered by humification (see Fig. 7C), as several of the most pronounced changes in the humification record are reflected in the clr-Br record. For example, the humification record shows low values at 3.1–2.4 ka BP (indicating well humified peat) followed by an increasing trend towards the surface (increasingly poorly humified peat), which is mirrored by the changes in the clr-Br record. The statistical analysis showed a significant positive correlation between the clr-Br and humification records (R=0.46). While these results support the interpretation that peat decomposition influences the Br concentration in the peat, this is the opposite relationship to previous studies (Biester et al., 2004; Martnez Cortizas et al., 2007). Furthermore, the clr-Br record also positively correlates with the IR record (R=0.48) and the PCA analysis shows clr-Br and clr-Ca have similar variability (both with positive PCA 2 loadings), indicating that some of the bromine preserved within the peat could reflect past periods of high storminess when both sand and sea spray were deposited. Despite these results, the peaks in IR do not appear to coincide with many of the increases in clr-Br aside from during recent centuries, where both the IR and clr-Br peak at 0.2 ka BP (Fig. 7C and F). We speculate that in this youngest section of core the peat has not been as decomposed as deeper sections, preserving the Br storm signature within the sediments.

Si and Ca show promise for future storminess reconstructions, particularly in sites where the sand content of the peat is higher than the present study. Given the apparent complex interplay between Br deposition and humification, the former may not be a useful proxy for either past storminess or peat humification, although future research could potentially use a high-resolution humification record to normalise the Br record, with residual Br showing variability in deposition and thus storminess. The storminess reconstruction in the following sections is therefore based on the IR results, rather than the element results.

Figure 8 Summary of selected reconstructions of northwest European storminess (b–e). (a) Ignition residue results from Roycarter Bog (County Mayo), Ireland (Sjöström et al., 2024). (b) Ignition residue results from Bullsmouth Bog (County Mayo), Ireland (this study). (c) Western Scotland storm stack developed by Kylander et al. (2023) based on the results of Orme et al. (2016) and Kylander et al. (2020). (d) Southern Sweden storm stack based on the results of Kylander et al. (2023). (e) Ignition residue results from Cors Fochno Bog, Wales (Orme et al., 2015). Frequently used climate intervals are marked at the top; including the Roman Warm Period (RWP; Hu et al., 2022), Dark Ages Cold Period (DACP; Helama et al., 2017), Medieval Climate Anomaly (MCA; Graham et al., 2011) and Little Ice Age (LIA; Wanner et al., 2022).

5.2 Storminess in northwest Europe

A reconstruction of storminess has already been developed from Roycarter Bog (Sjöström et al., 2024), located just 10 km to the north of Bullsmouth Bog, allowing regional patterns of storminess to be explored. Storm periods at Roycarter Bog, inferred from grain size distributions, mostly corresponded to IR increases in this record, except for two relatively long storm periods at 6.15–5.5 and 4.97–4.13 ka BP. From c.4.3 ka BP there is a similarity between the IR from Bullsmouth Bog and the IR record from Roycarter Bog (Fig. 8), with peaks in both records at approximately 4.35–4.1, 3.45 and 3.2 ka BP. Between 3 and 1.5 ka BP there are some differences between the two records, including a pronounced peak in the Bullsmouth Bog reconstruction at 2.4 ka BP whereas Roycarter Bog exhibits a peak at 2.1–1.8 ka BP. From 1.5 to 1 ka BP the records both show enhanced storminess that reaches a maximum at c. 1.1 ka BP, followed by a decrease in storminess at ∼ 1 ka BP lasting two to three hundred years. During the last 800 years the IR results from Bullsmouth Bog are significantly higher than the older part of the record, an increase which is not observed in the Roycarter record. Bullsmouth Bog may have received more soil from local agriculture or land disturbance after 0.8 ka BP, as the site is positioned on the eastern side of Achill Island where the prevailing westerly winds would first cross the island. This is supported by increased detection of Fe and K during the last millennium (Fig. S1), as Fe₂O₃ and K₂O are constituents of the podzol soils that occur on the areas of Achill Island that are not covered by peatland (Gardiner, 1980; ISRIC, 2025). Nevertheless, during this time the two records show the same peaks in IR at c. 0.6–0.5 and 0.2 ka BP that are likely indicative of enhanced storminess. An earlier potential interval of storminess has been identified in a peat core situated on the western side of Achill Island spanning the mid-Holocene, as there is a mineral (silt) layer dated to c. 5.2–5.1 ka BP, which is thought to have resulted from a series of storms (Caseldine et al., 2005). This coincides with rather small peaks in the Bullsmouth Bog and Roycarter records (<1 % increase in IR) supporting the interpretation that this was a storm event, with the closer proximity to beaches on Achill Island a likely explanation for the larger amount of sand and silt deposition in the Caseldine record (>3 % increase in the IR).

We next review the storm events identified in similar research that has been conducted across northwest Europe, to assess whether the episodes of apparently enhanced storminess in Ireland had impacts more widely (Fig. 8). While storminess records have been developed from different depositional settings in northwest Europe (i.e. estuaries, sand dunes, cliff top boulder deposits; Sorrell et al., 2012; Hansom and Hall, 2009; Goslin et al., 2018; Wilson and Braley, 1997) using a range of methods and proxies, we here focus on those records that span the mid-late Holocene interval and have measured the aeolian sand influx to coastal peatlands. This is to avoid differences between records caused by the choice of proxies and depositional settings, as they can reflect different processes (i.e. wind vs. wave energy), and likely have different thresholds for the transport and deposition of material, and thus sensitivity to storm events.

In Wales a reconstruction of past storminess, likewise derived from the IR of a coastal peat sequence (Orme et al., 2015), shows similar timing in storminess peaks to those in western Ireland, particularly during the last 2500 years, with high storminess at 2.3–2, 1.4, 1.1, 0.55 and 0.2 cal ka BP (Fig. 8). The storm stacks from western Scotland and Southern Sweden (Kylander et al., 2023; Fig. 8) show that storminess was high at c.2.2–2 ka BP and 0.5 ka BP, highlighting the occurrence of high storminess more widely across Sweden, Scotland, Wales and Ireland at these times. The South Sweden stacked record furthermore shows peaks in storminess at c. 1.1 and 2.8 ka BP that may reflect the storm events in the records from Ireland and Wales.

This comparison highlights the benefit of developing multiple reconstructions of storminess across a region using the same method, as well as some challenges. Shared storm events emerge through this comparison particularly during the last 2500 years, including at c. 2.25, 2, 1.4, 1.1, 0.5 and 0.2 ka BP. There are intervals of lower storminess across the records at 1–0.6 and 1.8–1.5 ka BP, corresponding to the Medieval Climate Anomaly (1.05–0.7 ka BP; Graham et al., 2011) and during the late Roman Warm Period (2.15–1.55 ka BP; Hu et al., 2022), suggesting that broadly warmer intervals had reduced storm intensity and/or frequency. Prior to 2.5 ka BP, while some storm events are seen across many of the records others appear non-synchronous or absent in one record but present in others, and sand peaks in the Bullsmouth record are muted prior to 4 ka BP. Differences between records could result from dating uncertainties and local factors, for example, changes in forest cover as a barrier to sand transport. A mid-Holocene pollen record from Achill Island shows that around 4.6 ka BP pine pollen reduced and grass pollen increased (Caseldine et al., 2005), suggesting that forest cover before 4.6 ka BP may have blocked sand transport across the island, and during the late Holocene a more open landscape allowed this transport process to occur. Sea level changes and associated shifts in coastal beaches and dunes are another factor that can complicate comparisons between records, for example on the coast of County Mayo there has been a gradual transgression, from sea levels that were approximately 7 m below present at 7 ka BP (Shennan et al., 2018). This may have caused the lower magnitude of the IR peaks in the earlier parts of the Bullsmouth and Roycarter Bogs, as a lower sea level during the mid-Holocene likely caused a greater distance between the core site and coastal beaches and sand dunes than seen at present. It is also not known when the sand dunes formed at Keel Beach to the west of Bullsmouth Bog, or if there were transgressive phases when sand may have been more easily entrained. Nevertheless, this record comparison and others (Sorrell et al., 2012; Vandel et al., 2019; Kylander et al., 2023; Leszczyńska et al., 2025) demonstrate that greater certainty in past storm events can be achieved where simultaneous storm signals are identified across regional sites, which helps to remove some local factors.

Figure 9 Summary of the sedimentary based storminess records from Ireland and possible drivers of annual to decadal storm variability. (a) Reconstruction of ice-rafted debris (IRD) from the North Atlantic (percentage of haematite stained grains from 4-stacked records from the North Atlantic; Bond et al., 2001). (b) Total Solar Irradiance reconstruction (Steinhilber et al., 2009). (c) Reconstructions of storminess from the Bullsmouth Bog record (this study). (d) Reconstruction of volcanic activity from the West Antarctic Ice Sheet Divide deep ice core based on annual non-sea salt sulphate flux (Cole-Dai et al., 2021), with red circles showing the 6 major tropical eruptions since 2.5 ka BP (426 BCE, 44 BCE (or 43 BCE; McConnell et al., 2020), 540 CE, 1257 CE – Samalas, 1458 CE, 1815 CE – Tambora) identified using ice core records from Greenland and Antarctica (Sigl et al., 2015). The purple circle marks the 852/3 CE Churchill Eruption considered a high intensity eruption but without a strong signature in the ice core records (Mackay et al., 2022). The periods of Northern European storminess identified in Fig. 8 are shown by the blue shading. The lower plot shows the median age and 2σ uncertainties for the peaks in the IR record for Bullsmouth Bog and Roycarter Bog, plotted with the age of the 7 major eruptions identified by Sigl et al. (2015) and Mackay et al. (2022).

5.3 Storminess, sea ice and major volcanic eruptions

The drivers of past storminess over Europe have been discussed widely as new sedimentary based storminess records have been developed, with the influence of solar (Martin-Puertas et al., 2012; Kylander et al., 2020; Goslin et al., 2018) and internal oceanic variations (Sorrel et al., 2012) often suggested. A comparison of the storminess records from northwest Europe with records of potential drivers is presented in Fig. 9.

Many of the identified peaks in storminess, particularly during the last 2500 years, coincide with major volcanic eruptions (Sigl et al., 2015; Cole-Dai et al., 2021; Mackay et al., 2022). There is an apparent correspondence between the timing of major volcanic eruptions at 0.14, 0.49, 0.69, 1.10, 1.41, 1.99, 2.38 ka BP (426 BCE, 44 BCE, 540 CE, 852/3 CE, 1257 CE, 1458 CE, 1815 CE) and many of the peaks in storminess observed in the Bullsmouth Bog record and the selected records shown in Fig. 9, as highlighted by the blue shading. While the chronologies of the storm reconstructions have relatively large age uncertainties, there is a good alignment between the median age estimates of the maximum IR peaks and the eruption years (Fig. 9). Eruptions can potentially cause a strengthening of the European jet stream, a positive NAO equating to a more vigorous atmospheric circulation and associated enhanced storminess, as identified by studies of eruptions during the last millennium (Andreasen et al., 2024; Brönnimann et al., 2025) and by modelling simulations (Shindell et al., 2004; Fischer et al., 2007; Zambri and Robock, 2016; Zambri et al., 2017). Documentary evidence shows that in Ireland following the 852/3 CE Churchill eruption there were several stormy years (Ludlow et al., 2013; Mackay et al., 2022; Fig. 1), which may be associated with the peak in IR at 1.1 ka BP observed in the Irish and Welsh reconstructions. Although some studies have questioned the reality (or likelihood) of an invigorated Northern Hemispheric westerly circulation following volcanic eruptions (e.g. Polvani et al., 2019; Tejedor et al., 2024), our SEA analysis (Fig. 6) suggests that enhanced storminess has repeatedly occurred post-eruption for Ireland's northeast Atlantic location. This increase in storminess does not seem to result from any systematic uni-directional shift in the prevailing NAO pattern, at least as based upon a comparison of the timing of major eruptions with the NAO reconstructions of Trouet et al. (2009) and Olsen et al. (2012), which show that eruptions and the identified periods of high storminess were associated with positive, negative and neutral NAO states (Fig. S3). However, these NAO reconstructions may not capture shorter-term interannual excursions in the NAO as they are intended to reconstruct decadal to multi-decadal and longer variations.

A factor that may have amplified, and perhaps prolonged, the deposition of sand to the Irish and Welsh peatlands following eruptions may have been vegetation responses to the eruptions. Colder temperatures and reduced sunlight caused by stratospheric aerosols potentially reduced coastal vegetation cover and may have promoted sand dune instability, leading to enhanced sand erosion and transport inland (Jackson et al., 2019). Following major tropical eruptions in recent centuries the summer temperatures in northern Europe were lower (e.g. Fischer et al., 2007; Sigl et al., 2015) and severe winter temperatures with reports of snow and ice cover in Ireland also followed major eruptions (Ludlow et al., 2013). Tree ring records show summer temperatures decreased by an average of ∼ 0.6 °C and took 5–10 years to fully recover, with tropical eruptions of similar or larger magnitude to the Tambora eruption in 1815 CE causing cooling of 0.89 °C (Sigl et al., 2015). Tree-ring maximum latewood density measurements, rather than ring widths (as used by Sigl et al., 2015), are increasingly preferred for palaeoclimatic reconstructions in better capturing interannual temperature variance and such records often imply a somewhat shorter post-volcanic recovery time, but also a more severe temperature impact (e.g. Esper et al., 2015). Following the 1991 CE Pinatubo eruption, the amount of direct (diffuse) sunlight decreased (increased), with an overall 2.5 % decrease in total sunlight (Proctor et al., 2018). Episodes of explosive volcanism have also been associated with cooling of the North Atlantic during recent centuries (resulting in part from enhanced evaporative cooling from stronger westerly winds) (Knudson et al., 2014; Birkel et al., 2018), potentially enhancing the cooling over surrounding regions such as western Ireland. The impacts of these changes in temperature and sunlight on vegetation is a topic of debate; evidence has shown changes in sunlight reduces crop yields globally following eruptions (Proctor et al., 2018), and following past eruptions in recent millennia there has widely been reduced tree growth (Krakauer and Randerson, 2003; Sigl et al., 2015). However, enhanced diffuse light following eruptions (caused by increased scattering of incoming direct solar radiation by volcanic aerosols) is also thought to promote photosynthesis (e.g. Gu et al., 2003). Assessing whether vegetation changes occurred on Achill Island at times of eruptions is difficult, as there are no late Holocene pollen records, and records from elsewhere in Ireland are complicated by human alterations to vegetation (e.g. Brown et al., 2005).

The cooler temperatures following major eruptions may have also enhanced the transport of sand to the peatlands through niveo-aeolian processes. Studies have suggested that snow- and ice-covered landscapes can enhance the transport of sand (Björck and Clemmensen, 2004 and references therein; Vaarsma et al., 2025). The reports in the Irish Annals of severe snow and ice cover in years following eruptions (Ludlow et al., 2013) supports that this may be an additional factor amplifying sand transport and deposition from storms following eruptions.

The patterns of storminess identified in the Welsh and Irish sedimentary storm records resemble the ice-rafting record from the North Atlantic (Bond et al., 2001; Fig. 9), with peaks in IRD also coinciding with the major eruptions particularly during the last 2.5 ka BP. We hypothesise that the peaks in IRD may also have been triggered by these eruptions, with cooler Atlantic SST's, an expansion of Arctic sea ice and greater export of sea ice southward in the North Atlantic. Model simulations, reconstructions and observations have shown cooler Atlantic SST's and enhanced sea ice extent following eruptions (e.g. Knudson et al., 2014; Gagné et al., 2017; Birkel et al., 2018; Slawinska and Robock, 2018; Sicre et al., 2013; Zhong et al., 2011; Miller et al., 2012). Many of these studies indicate sea ice enhancement for up to a decade. Moreover, Lehner et al. (2013) proposed a feedback where greater sea ice export could cause freshening in the North Atlantic, reduced convection and weakening of the AMOC and northward heat transport, resulting in further cooling and prolonged sea ice expansion.

The general cooling of the LIA is thought to have been caused by a combination of volcanic eruptions and grand solar minima, which were amplified by ocean and atmospheric feedbacks including weakening of the thermohaline circulation (see references in Wanner et al., 2022). Visual comparison of the Total Solar Irradiance record (Steinhilber et al., 2009) with the storminess records (Fig. 9) shows that some solar minima coincide with peaks in storminess while others do not. Therefore, while solar minima likely contributed to lower LIA temperatures, it is not clear from these results whether there is an influence on storm intensity or frequency.

A challenge of assessing the impact of volcanic eruptions on Holocene storminess for this region is that the sedimentary based storm records considered here have decadal resolutions, whereas eruptions are short-lived events with climate impacts from aerosols generally lasting <10 years, although ocean temperature and circulation impacts might play out over longer timescales (e.g. Pausata et al., 2015b) with uncertain outcomes on storminess. While the peatland-based storm records appear to show decadal to centennial intervals of increased storminess, it is not presently possible to determine if the associated sand was deposited during one major event or a series of more moderate storms. While the age resolution of the records is decadal to multi-decadal (averaging 11, 16 and 27 yr cm⁻¹ for the Bullsmouth Bog, Cors Fochno and Roycarter records, respectively, since 2.5 ka BP), it is possible that a single storm during these intervals may have deposited the measured sand. The peaks in sand span between 1 and 10 cm in the three records, meaning that either periods of high storminess lasted 10's–100's of years, or there was some bioturbation or infiltration of the sand grains within the sediment column that spread sand from short-lived stormy periods vertically within the peat column. This latter outcome is possible, as experiments with tephra deposits have shown shards move up to 15 cm down through the peat, with upward movement occurring with plant growth and water table variability (Payne and Gehrels, 2010).

The hypothesis that volcanic eruptions during the Holocene have been a driver of episodes of enhanced storminess in northwest Europe ultimately requires testing further. Payne and Egan (2019) highlight some challenges of investigations into volcanic impacts in palaeoecological studies, which are also applicable to this study on the impacts on storminess. They highlight the difficulty in proving cause-effect relationships, uncertain interpretations and chronological uncertainties. While this study follows the advice of Payne and Egan (2019) to replicate between studies and base interpretations on better understood modern volcanic impacts, the sedimentary based storminess records have large age uncertainties (e.g. the Bullsmouth Bog storm events since 2.5 ka BP have 2σ age ranges of between 280 and 430 years; Fig. 9) and several possible interpretations (i.e. enhanced storminess or vegetation responses). Nonetheless, Ireland's long written medieval record of storms and stormy years has in principle shown an association with ice-core-based dates of explosive volcanic eruptions. This association becomes clearer and more statistically significant as smaller (in terms of sulphate output) and potentially less climatically impactful eruptions are iteratively excluded from the analysis. This finding supports the credibility of a volcanic influence on storminess for Ireland's northeast Atlantic position and shows the potential of employing environmental data from written archives alongside palaeoenvironmental proxies.

Nonetheless, some aspects of this evidence bear further consideration. Whilst the observed association is a recurrent one, it is worth noting explicitly that storms are reported in the Irish annals without associated ice-core evidence of explosive volcanism. This is expected as the westerly winds prevail here and characteristically expose Ireland to accompanying Atlantic depressions, including periodically severe and damaging events as reported in the annals (Sweeney, 1997, 2000). Conversely, some large eruptions registering in the polar ice do not have any associated reporting of storms or stormy conditions. This is also expected, not least because of the character of the annalistic record. Despite providing a remarkable continuity of coverage spanning more than a millennium, multiple texts that once existed, and which doubtless included further weather reporting, are now lost, whilst the geography and temporal density of the surviving coverage varies (McCarthy, 2008). Such factors have been shown to influence the density of weather reporting (Ludlow, 2010, 2012). While the surviving record offers a sample sufficient to test hypothesised associations between phenomena (here, volcanic eruptions and storminess, or elsewhere, cold winters and eruptions; Ludlow et al., 2013; Kostick and Ludlow, 2022), it cannot provide a complete accounting of whether notable storminess followed any given eruption. The same considerations make the use of the annalistic records to infer longer-term trends in storminess similarly challenging, without normalizing reported storm frequencies to account for changing record densities. Addressing this in future work (requiring a count of total annalistic entries through time or other measures of coverage density) would allow a credible comparison of trends between the annalistic and sedimentary based records, as well as a windiness index weighted by annual record density for application in SEA approaches.

A further potential explanation of the absence of reported storminess following some (even large) eruptions is the inherent variability in post-volcanic climatic responses, such that not all eruptions can be expected to induce increased storminess for Ireland's high mid-latitude location. Many variables beyond sulphate emission volume govern the climatic response. These include the timing (seasonality) of the eruptions (e.g. Kravitz and Robock, 2011; Zhuo et al., 2021), which relates to the importance of initial conditions including the prevailing state of the ocean-atmospheric system and associated major modes of variability such as the NAO (plus polar vortex) and ENSO (e.g. Pausata et al., 2015a; Fuglestvedt et al., 2024; Zhuo et al., 2024), with internal climatic variability feasibly overwhelming the forcing from some eruptions at least partly depending on such factors (e.g. Polvani et al., 2019; Tejedor et al., 2024). The altitude of sulphate injected into the atmosphere is also critical (Burke et al., 2019; Toohey et al., 2019), with some apparently large Greenland sulphate deposition signals associated with a mainly tropospheric (not stratospheric) aerosol presence, and hence a shorter and less severe climatic impact. Guillet et al. (2023) have, for example, identified a substantial Greenland sulphate deposition signal in 1182 CE that may have been a tropospheric-only event (i.e. the eruption failed to inject a substantial volume of its sulphate output into the stratosphere where it would have had a longer residence time and ultimately larger spatial coverage). Perhaps significantly, the Irish annals provide no obvious evidence of any corresponding storm or windiness.

The differing climatic impacts of low-latitude (tropical) versus extratropical (particularly high-latitude) Northern Hemispheric eruptions has been a major focus of study. Relevant here is the contention that low-latitude eruptions can induce an invigorated winter-season westerly flow at mid and high Northern Hemispheric latitudes. It is thought that this occurs when the meridional (latitudinal) thermal gradient is steepened through the preferential warming of the low-latitude stratosphere where the (at least initial) concentration of sulphate aerosols acts to absorb outgoing longwave radiation and incoming near-infrared radiation, producing a relative “winter warming” effect over Northern Hemispheric landmasses that derives from an enhanced westerly flow (e.g. Robock and Mao, 1992; Robock, 2000; Fischer et al., 2007; Schneider et al., 2009; Zambri et al., 2017). If correct (and see Paik et al., 2023 for a review of opposing findings), this offers one explanatory hypothesis for the apparent association between volcanic eruptions and storminess for Ireland, and we note that periods of Northern European storminess identified from our sedimentary based records mainly appear associated with some of the largest low-latitude eruptions (Fig. 9). However, the earlier finding by Ludlow et al. (2013) of an association between eruptions and severe winter cold events (also reported in the Irish annals) suggests considerable potential variability in responses because some of these events, particularly those involving severe and prolonged frosts, are likely to have arisen from large-scale anticyclonic (negative NAO) conditions less associated with storminess (though some, particularly those involving heavy snowfall, are not incompatible with strong cyclonical development). Overall, this variability is consistent with the interplay of multiple factors determining the outcome of volcanic climatic forcing, as noted above. It is also potentially consistent with modelling by Zambri et al. (2017) that suggests any enhanced westerly flow can be confined mainly to the first post-eruption winter, after which sulphate aerosols may become distributed more broadly across the Northern Hemisphere, thereby diminishing the enhanced thermal gradient thought to drive this dynamical response, allowing the aerosol-induced tropospheric radiative cooling response to dominate. Supportive evidence is offered by Galvin et al. (2011) who examined the winter-season response to five explosive tropical eruptions in long instrumental data (post-1800) for Ireland and found minor warming in the first and second post-eruption winters, switching to cooling thereafter.

The GISP2 Greenland ice-core sulphate record used by Ludlow et al. (2013) did not, however, distinguish low-latitude (tropical) from extratropical eruptions, whereas the Sigl et al. (2015) bipolar volcanic forcing history employed here identifies probable low-latitude eruptions (when sulphate is deposited coterminously in Antarctica and Greenland) from extratropical eruptions (when sulphate is deposited in only one polar region). It may be important to note, therefore, that in conducting our SEAs (Fig. 6), the largest set of eruption dates, in which no minimum Greenland sulphate deposition limit was imposed, and for which no statistically significant association with storminess was observed, is dominated by extratropical (including high latitude) eruptions. These comprise 54 (63.5 %) of the 85 eruptions between 600 and 1616 CE, including many apparently small events, and are unlikely to have led to a steepening of the meridional thermal gradient. While Robock and Mao (1992) did see evidence of (and posited potential mechanisms for) enhanced westerly flow even after some mid- and high-latitude Northern Hemispheric eruptions, recent work has found that such eruptions can instead promote a negative NAO (Sjolte et al., 2021), associated with decreased windiness (or storm intensity; Matthews et al., 2016) at Ireland's high mid-latitude location (Sweeney, 1997). Greenland ice-cores also preferentially capture small Icelandic eruptions due to their proximity (Abbott and Davies, 2012), many of which many have been less climatically impactful with tropospheric-only sulphate release (although multi-year periods of prolonged release may be an exception; Gabriel et al., 2024).

When imposing a minimum Greenland deposition volume of 10 kg km⁻² (the first variant eruption list for which a statistically significant association with windiness index values is observed, at 93.7 % confidence, for superposed year −1 (Fig. 6)), the percentage share of extratropical eruptions drops to 51 % (25 of 49 eruptions), with low-latitude eruptions almost at parity, and when imposing a minimum of 20 kg km⁻² (for which both superposed years −1 and 1 first become simultaneously significant at >90 % confidence), the share drops marginally further to 48 % (12 of 25 eruptions), indicating the importance of low-latitude eruptions in the observed association. Nonetheless, extratropical eruptions remain in all variant eruption lists, including the list (n=8) that achieves the highest statistical significance (at 98.1 % confidence in year −1, excluding all events with <40 kg km⁻² deposition), which comprises 5 extratropical events (two with positive windiness index values) and 3 low-latitude events (one with a positive value). Similarly, the variant eruption list (n=19) that achieves the second highest statistical significance (at 98.02 % confidence for year −1, excluding all eruptions with <25 kg km⁻² deposition) comprises 9 low-latitude eruptions (two with positive windiness index values) and 10 extratropical events (three with positive values). Our results thus also imply a role for extratropical eruptions in storminess for Ireland, suggesting a need to further explore relevant mechanisms (e.g. a “local stratospheric top-down mechanism” proposed by Guðlaugsdóttir et al., 2025, p. 3961).

Further work can benefit from the inclusion of documentary storm records from a broader region of Europe (e.g. Lamb and Frydendahl, 1991; De Kraker, 1999, 2005; Garnier et al., 2018) to provide a fuller chronology of storm events against which to assess the respective influence of eruptions not only by deposition volume, but also by eruption location and the identification of sulfur mass-independent fractionation signatures to discriminate between eruptions that did and did not achieve a stratospheric sulphate injection (e.g. Burke et al., 2019, 2023; Gautier et al., 2019; Hutchison et al., 2024). For our sedimentary proxies, future work in northwest Europe may benefit from using cryptotephras for Icelandic eruptions to better constrain the timing of the storm events and allow correlation between the available sedimentary-based records. Furthermore, pollen analysis on the same cores may allow synchronous local vegetation responses to eruptions to be detected. While the focus here has been on the possibility of volcanic eruptions as drivers of past storminess for Ireland's northeast Atlantic (and proximate) locations, we cannot exclude other possible drivers, nor the important role of internal variability. There are similar challenges when investigating whether solar minima caused high storminess during the Holocene, as chronological uncertainties make it challenging to align the TSI and storminess records.

6 Conclusions

The Bullsmouth Bog storminess reconstruction from the west coast of Ireland spans the last 7 ka. A comparison with peatland-derived storminess records from nearby in Ireland and the west coast of Wales have enabled a regional storm event sequence to be developed, showing synchronous periods with high storminess at 2.8, 2.25, 2, 1.4, 1.1, 0.5 and 0.2 ka BP. The storm events at c. 2.8. 2.2–2. 1.1 and 0.5 ka BP are observed more widely in records from Ireland, Wales, Sweden and Scotland. There are also periods with widely reduced storminess coinciding with warmer climate intervals, notably at 1–0.6 and 1.8–1.5 ka BP during the Medieval Climate Anomaly and late Roman Warm Period. The records from western Ireland show earlier shared storm peaks at 4.35–4.1, 3.45 and 3.2 ka BP, although these are not observed in the other reconstructions from northwest Europe. Prior to ∼ 4.5 ka BP the Bullsmouth Bog storminess reconstruction had low magnitude variability, which may have been caused by more distal or reduced sand sources as a result of lower sea levels at this time, or greater forest cover inhibiting sand transport.

Consideration of possible drivers reveal that many of the stormy periods coincide with major volcanic eruptions and periods of enhanced ice-rafting in the North Atlantic, particularly when looking at the last 2.5 ka BP. We hypothesise that peaks in storminess may have often resulted from an enhancement in the strength of the jet stream following eruptions, in agreement with the atmospheric response shown by observations and models. An additional factor contributing to greater sand transport may have been a reduction in coastal vegetation resulting from reduced sunlight and temperatures following eruptions. The coincidence with high ice-rafting to the North Atlantic may also reflect a shared driver, as volcanic eruptions contribute to colder LIA-like events where ice-rafting is more prevalent. While the evidence points towards a volcanic driver for past storm events, confirming this is challenging given the significant chronological uncertainty of the storminess reconstructions and a lack of clarity on the environmental conditions and processes responsible for the sand deposition (i.e. the role of landscape change or snow cover following eruptions and the time-frame for sand deposition). This requires improved chronologies and investigation of whether the peatland-based storm records reflect storminess or other environmental responses related to the eruptions. Nevertheless, the volcanic driver hypothesis is supported by an apparent statistical association between storms and stormy years reported in Irish medieval annalistic texts and explosive volcanic eruptions dated between 600 and 1616 in ice-core-based volcanic forcing histories. These results suggest that written records have potential to contribute further here, including in a more direct future comparison of documented trends in storminess and the evidence of the available sedimentary proxies.

An additional aim of this research was to test whether XRF core scanning can be usefully utilised as a fast and high resolution method of reconstructing past storminess, with Ca, Si and Br targeted as potential indicators of carbonate shells, quartz sand grains and sea spray deposition, respectively. The Si and Ca records were compared with the ignition residue record; Ca showed the low magnitude variations but not the peaks in sand content, while Si showed peaks in sand content only within the last millennium where the mineral content was highest. The findings support previous research concluding that retained bromine concentrations along the core are primarily influenced by variations in peat decomposition, although some influence from storm variability cannot be discounted.