Eruptive history and 40 Ar/ 39 Ar geochronology of the Milos volcanic 1 field, Greece 2

. High-resolution geochronology is essential for determining the growth-rate of volcanoes, which is one of the key 7 factors for establishing the periodicity of volcanic eruptions. However, there are less high-resolution eruptive histories (>10 6 8 years) determined for long-lived submarine arc volcanic complexes than for subaerial complexes, since submarine volcanoes 9 are far more difficult to observe than subaerial ones. In this study, high-resolution geochronology and major element data are 10 presented for Milos Volcanic Field (VF) in the South Aegean Volcanic Arc, Greece. The Milos VF has been active for over 3 11 Ma, and the first two million years of its eruptive history occurred in a submarine setting that has been emerged above sea 12 level. The long submarine volcanic history of the Milos VF makes it an excellent natural laboratory to study the growth-rate 13 of a long-lived submarine arc volcanic complex. This study reports twenty-one new high-precision 40 Ar/ 39 Ar ages and major 14 element compositions for eleven volcanic units of the Milos VF. This allows us to divide the Milos volcanic history into at 15 least three periods of different long-term volumetric volcanic output rate (Q e ). Period I (submarine, ~3.3-2.13 Ma) and III 16 (subaerial, 1.48 Ma-present) have low Q e of 0.9 ± 0.5 × 10 -5 km 3 ·yr -1 and 0.25 ± 0.05 × 10 -5 km 3 ·yr -1 , respectively. Period II 17 (submarine, 2.13 - 1.48 Ma) has a 3-12 times higher Q e of 3.0 ± 1.7 × 10 -5 km 3 ·yr -1 . The Q e of the Milos VF is 2-3


Introduction
Short-term eruptive histories and compositional variations of lavas and pyroclastic deposits of many arc volcanic fields are well established.However, high-resolution eruptive histories that extend back > 10 5 -10 6 years have been determined only for a handful of long-lived subaerial arc volcanic complexes.Some examples are: Mount Adams (Hildreth and Lanphere, 1994), Tatara-San Pedro (Singer et al., 1997), Santorini (Druitt et al., 1999), Montserrat (Cole et al., 2002), Mount Baker (Hildreth et al., 2003a), Katmai (Hildreth et al., 2003b), and Ceboruco-San Pedro (Frey et al., 2004).To establish the growth-rate of volcanic complexes and disentangle the processes responsible for the eruption, fractionation, storage, and transport of magmas over time, comprehensive geological studies are required.These include detailed field mapping, sampling, high-resolution geochronology and geochemical analysis.Based on these integrated studies, the growth-rate of volcanoes can be determined to establish the periodicity of effusive and explosive volcanism.
The Milos Volcanic Field (VF) is a long-lived volcanic complex that has been active for over 3 Ma.The Milos VF erupted for a significant part of its life below sea level, similar to the other well studied volcanic structures in the eastern Mediterranean (Fytikas et al., 1986;Stewart and McPhie, 2006).The eruptive history of the Milos VF has been examined with a broad range of chronostratigraphic techniques such as K-Ar, U-Pb, fission track, 14 C and biostratigraphy (e.g.Angelier et al., 1977, Fytikas et al., 1976, 1986, Traineau and Dalabakis, 1989, Matsuda et al., 1999, Stewart and McPhie, 2006, Van Hinsbergen et al., 2004and Calvo et al., 2012).However, most of the published ages have been measured using the less precise K-Ar or fission track methods, and modern, high precision 40 Ar/ 39 Ar ages for the Milos VF have not been published so far.In this study, (1) we provide high-precision 40 Ar/ 39 Ar geochronology of key volcanic units of the Milos VF and (2) refine the stratigraphic framework of the Milos VF with the new high-precision 40 Ar/ 39 Ar ages and major element composition.(3) We also quantify and constrain the compositional and volumetric temporal evolution of volcanic products of the Milos VF.

Geological setting
The Milos VF is part of the South Aegean Volcanic Arc (SAVA), an arc which was formed in the eastern Mediterranean by subduction of the African plate beneath the Aegean microplate (Figure 1, Nicholls, 1971;Spakman et al., 1988;Duermeijer et al., 2000;Pe-Piper and Piper, 2007;Rontogianni et al., 2011).The present-day Benioff zone is located approximately 90 km underneath Milos (Hayes et al., 2018).The upper plate is influenced by extensional tectonics (e.g.McKenzie, 1978;Pe-Piper and Piper, 2013), which is evident on the island of Milos as horst and graben structures (Figure 2).
The Milos VF is exposed on the islands of the Milos archipelago: Milos, Antimilos, Kimolos and Polyegos.The focus of this study is Milos which has a surface area of 151 km 2 .The geology and volcanology of Milos have been extensively studied in the last 100 years.The first geological map was produced by Sonder (1924).This work was extended by Fytikas et al. (1976) and Angelier et al. (1977) and the subsequent publications of Fytikas et al. (1986) and Fytikas (1989).Interpretations based on volcanic facies of the complete stratigraphy were made by Stewart andMcPhie (2003, 2006).More detailed studies of single volcanic centres (e.g.Bombarda volcano and Fyriplaka complex) were published by Campos Venuti and Rossi (1996) and Rinaldi et al. (2003).Milos has also been extensively studied for its epithermal gold mineralization, summarized by Alfieris et al. (2013).Milos was known during the Neolithic period for its export of high-quality obsidian.Today the main export product is kaolinite mined from hydrothermally altered felsic volcanic units in the centre of the island (e.g.Alfieris et al., 2013).
The geology of Milos can be divided into four main units: (1) metamorphic basement, (2) Neogene sedimentary rocks, (3) volcanic sequences and ( 4) the alluvial cover.The metamorphic basement crops out at the southwest, south and southeast of Milos (Figure 3) and is also found as clasts in many volcanic units.The metamorphic rocks include lawsonite-free jadeite eclogite, lawsonite eclogite, glaucophane schist, quartz-muscovite-chlorite and chlorite-amphibole schist (Fytikas et al., 1976(Fytikas et al., , 1986;;Grasemann et al., 2018;Kornprobst et al., 1979).The exposed units belong to the Cycladic Blueschist Unit (Lower Cycladic nappe), whereas eclogite pebbles in the phreatic eruption products called "green lahar" by Fytikas (1977) are derived from the Upper Cycladic Nappe (Grasemann et al., 2018).On top of this metamorphic basement, Neogene fossiliferous marine sedimentary rocks were deposited (e.g.Van Hinsbergen et al. 2004).This sedimentary sequence can be divided into a lower unit A and upper unit B that is unconformably overlain by volcaniclastic sediments (Van Hinsbergen et al., 2004).Unit A is 80 m thick and consists of fluviatile-lacustrine, brackish and shallow marine conglomerate, sandstone, dolomite and limestone.Unit B is 25-60 m thick and consists of sandstone overlain by a succession of alternating marls and sapropels, suggesting a deeper marine setting (Van Hinsbergen et al., 2004).Five volcanic ash layers that contain biotite are found in this Neogene sedimentary sequence, either suggesting that volcanic eruptions in small volume already occurred in the Milos area or that these ash layers are derived from larger eruptions of volcanic centres further away from Milos (van Hinsbergen et al., 2004).Age determinations by bio-magneto-and cyclostratigraphy suggested that deposition of Unit A started at approximately 5 Ma, and that Milos subsided 900 m in 0.6 Ma (Van Hinsbergen et al. 2004) due to extension.This subsidence happened ca 1.0-1.5 Ma before the onset of the main phase of Pliocene-recent volcanism on Milos.
The Pliocene-recent volcanic sequence of Milos has been subdivided into different units by Angelier et al. (1977) and Fytikas et al. (1986).In addition, Stewart and McPhie (2006) provided a detailed facies analysis of the different volcanic units.The subdivision by Angelier et al. (1977) is not constrained well due to their limited amount of age data.The subdivision of volcanic units by Fytikas et al. (1986) and facies descriptions of Stewart and McPhie (2006) are summarized below.It is important to note that according to Stewart and McPhie (2006), the five volcanic cycles described by Fytikas et al. (1986) are difficult to match with existing age data and the continuous progression in volcanic construction (Fig. 4).For example, the first phase of Fytikas et al. (1986), the Basal Pyroclastic Series, contains the large pumice cone-crypto dome volcanoes according to Stewart and McPhie (2006).Two of these pumice-cone crypto dome volcanoes are much younger and intercalated between the Complex of Domes and Lava Flows (CDLF) of Fytikas et al. (1986).
Graded sandstone and bioturbated and fossil rich (in-situ bivalve shells) mudstone are intercalated, indicating a marine environment and a water depth of several hundreds of meters (e.g.Stewart, 2003;Stewart and McPhie, 2006), whereas later degassed magmas with a similar composition intruded as sills and cryptodomes.The BPS has been strongly affected by hydrothermal fluids, especially the proximal deposits (e.g.Kilias et al., 2001).The second volcanic unit was named the Complex of Domes and Lava Flows (CDLF, Fytikas et al., 1986) and the volcanic facies of this unit are described as the submarine dacitic and andesitic domes by Stewart and McPhie (2006).This phase of effusive submarine volcanism was predominantly andesitic/dacitic in composition and produced microcrystalline rocks with phenocrysts of pyroxene, amphibole, biotite and plagioclase.The eruption centres were mainly located along NNE faults and formed up to 300 m thick deposits extending over areas of 2.5 to 10 km 2 around the eruption centres.In the north-eastern part of Milos, an andesitic scoria cone provided scoria lapilli and bombs to deeper water settings.Sandstone intercalated in the CDLF contains both igneous and metamorphic minerals suggesting input from the basement.Rounded pebbles of rhyolite and dacite indicate that some of the volcanic deposits were above sea level, or in very shallow, near shore environments (e.g.Stewart and McPhie, 2006).
The third volcanic unit is called the Pyroclastic Series and Lava Domes (PSLD) by Fytikas et al. (1986) and belongs to submarine-to-subaerial dacitic and andesitic lava domes of Stewart and McPhie (2006).This highly variable group is dominated by rhyolitic, dacitic and andesitic lavas, domes, pyroclastic deposits and felsic pumiceous sediments (Stewart and McPhie, 2006).Thickness varies between 50-200 m, and the deposits are located in the eastern and northern parts of Milos (Figure 2 and 3).The initial pyroclastic layers were subaqueously deposited and the extrusion of a dome resulted in the deposition of talus around the margins by mass flow.On top of the dome sand-and siltstone with fossils (Ostrea fossil assemblage) and traction-current structures suggest that the top of the dome was above wave base.The youngest deposits of this unit are dacitic and andesitic lavas and domes.These domes generated subaerial block-and-ash flow and surge deposits.
Paleosols within these deposits are a clear indicator that some areas were above sea level.The last unit of the PSLD is represented by large subaerial rhyolitic lava that contains quartz and biotite phenocrysts and is found near Halepa in the southcentral part of Milos.Angelier et al. (1977) do not provide sample names, only numbers for the sample locations.Here the location is given after "Angelier_" (Angelier et al. 1977, their Fig. 3).Abbreviations: BPS=Basal pyroclastic series; CDLF=Complex of domes and lava flows; PSLD=Pyroclastic series and lava domes; CTF=Complexes of Trachilas and Fyriplaka.See more details in Figure .4.
The fourth unit consists of the subaerially constructed rhyolitic Complexes of Trachilas and Fyriplaka (CTF) (Fytikas et al., 1986), which Stewart and McPhie (2006) interpreted as subaerial rhyolitic lava-pumice cones.These two volcanic complexes are built from rhyolitic pumice deposits and lavas that contain quartz and biotite phenocrysts (10-20 modal %).The deposits have a maximum thickness of 120 m and decrease to several meters thickness in the distal parts.Basement-derived schist is found as lithic clasts (Fytikas et al., 1986).In addition, the Kalamos rhyolitic lava dome, which outcrops on the southern coast of Milos, produced lava that spread westwards to the Fyriplaka beach (Figure 2).This lava belongs to this fourth phase and is probably derived from an older volcano and not the Fyriplaka complex (Campos Venuti and Rossi, 1996).
The fifth volcanic unit comprises deposits from phreatic activity, especially in the northern part of the Zefiria Graben and near Agia Kiriaki (Figure 2 of Stewart and McPhie, 2006).Many overlapping craters are surrounded by lithic breccias that are composed of variably altered metamorphic basement clasts and volcanic clasts.This phreatic activity has continued into historic times (Trainau and Dalabakis, 1989).Fytikas et al. (1986) referred to this unit as "green lahar", although it indicated that this deposit is not a lahar but the product of phreatic eruptions in the last 0.2 Ma.

Mineral separation and sample preparation
Samples were collected from all major volcanic units on Milos island based on the studies of Fytikas et al. (1986), Stewart and McPhie (2006) and our own observations in the field.Photos of the sample locations and thin sections can be found in supplementary material I. Approximately 2 kg of fresh pumice clasts or lava was sampled from each unit.Samples were cut into ~5 cm 3 cubes using a diamond saw to remove potentially altered surfaces and obtain the fresh interior parts.These cubes were ultra-sonicated for 30 minutes in demi-water to remove dust and seawater and dried in an oven overnight at 50 °C.Dry sample cubes were crushed in a steel jaw crusher, and this fraction was split into two portions of roughly equal size.One of them was powdered in an agate shatter box and agate ball mill to a grain size of less than 2 µm for the major-element analysis.
The second fraction was sieved to obtain a grain size of 250-500 µm for 40 Ar/ 39 Ar dating.

40 Ar/ 39 Ar dating
The mineral and groundmass samples were wrapped in either 6-or 9-mm aluminium foil and packed in 20 mm aluminium cups, that were vertically stacked.Based on stratigraphy and previous geochronological constraints >1 Ma samples and the <1 Ma samples were irradiated for 7 and 1 hours respectively in irradiation batches VU108 and VU110 in the Cadmium-Lined in-Core Irradiation Tube (CLICIT) facility of the Oregon State University Training Research, Isotopes, General Atomics (TRIGA) reactor.The neutron flux for all irradiations was monitored by standard bracketing using the Drachenfels sanidine (DRA; 25.52 ± 0.08 Ma, modified from Wijbrans et al., 1995 andcalibrated relative to Kuiper et al., 2008) and Fish Canyon Tuff sanidine (FCs; 28.201 ± 0.023 Ma, Kuiper et al., 2008) with Min et al. (2000) decay constants.
In total, 24 samples (8 groundmasses, 15 biotites and 2 amphiboles, for sample G15M0026 both biotite and amphibole were analysed) were measured by either 40 Ar/ 39 Ar fusion and/or incremental heating techniques.For incremental heating experiments, 80-100 grains per sample were loaded into a 25-hole (surface per hole ~36 mm 2 ) copper tray together with single grain standards in ~12 mm 2 holes.The tray was prebaked in vacuum (10 -5 -10 -6 mbar) at 250 °C overnight to remove atmospheric argon and subsequently baked overnight at 120 °C in the ultra-high vacuum sample chamber (<5*10 -9 mbar) and purification system connected to a Thermo Scientific Helix MC mass spectrometer.
Samples and standards were heated with a focused laser beam at 8 % power using a 50W CW CO2 laser.The released gas was cleaned by exposure to a cold trap cooled by a Lauda cooler at -70 °C, a SAES NP10 at 400 °C, Ti sponge at 500 °C and cold SAES ST172 Fe-V-Zr sintered metal.The five isotopes of argon were measured simultaneously on five different collectors: 40 Ar on the H2-Faraday, 39 Ar on the H1-Faraday or the H1-CDD, 38 Ar on the AX-CDD, 37 Ar on the L1-CDD and 36 Ar on the L2-CDD for 15 cycles with 33 seconds integration time (CDD: compact discrete dynodes).The Faraday cups on H2 and H1 were equipped with 10 13 Ω amplifiers.Procedural blanks were measured every 2 or 3 analyses in different sequences, and airshots were measured every 8-12 hours to correct the instrumental mass discrimination.Gain between different collectors was monitored by measuring CO2 on mass 44 in dynamic mode on all collectors.Gain was generally stable over periods of weeks.
Note that because samples, standards and air calibration runs are measured during the same period, gain correction does not substantially change the final age results.The raw mass spectrometer data output was converted by an in-house designed Excel macro script to be compatible with the ArArCalc 2.5 data reduction software (Koppers, 2002).The 40 Ar/ 36 Ar atmospheric air value of 298.56 from Lee et al. (2006) is used in the calculations.The correction factors for neutron interference reactions are (2.64 ± 0.02) x10 -4 for ( 36 Ar/ 37 Ar)Ca , (6.73 ± 0.04) x10 -4 for ( 39 Ar/ 37 Ar)Ca, (1.21 ± 0.003) x10 -2 for ( 38 Ar/ 39 Ar)K and (8.6 ± 0.7) x10 -4 for ( 40Ar/ 39 Ar)K.All uncertainties are quoted at the 1σ level and include all analytical errors (i.e.blank, mass discrimination and neutron interference correction and analytical error in J-factor, the parameter associated with the irradiation process).
A reliable plateau age is defined as experiments with at least 3 consecutive steps overlapping at 2-sigma, containing >50% of the 39 ArK, a Mean Square Weighted Deviate (MSWD) value<2.5, and with a 40 Ar/ 36 Ar inverse isochron intercept that does not deviate from atmospheric argon at 2-sigma.All the inverse isochron ages used the same steps as used in the weighted mean ages, and all relevant analytical data for the age calculations following standard practices (Schaen et al., 2020) can be found in supplementary material II.

Whole-rock major element analysis by XRF
Major-element concentrations were measured by X-ray fluorescence spectroscopy (XRF) on a Panalytical AxiosMax.A Panalytical Eagon2 was used to create 40mm fused glass beads of Li2B4O7/LiBO2 (65.5:33.5%,Johnson & Johnson Spectroflux 110) with a 1:6 dilution sample-flux ratio that were molten at 1150 °C.Sample powders were ignited at 1000 °C for 2 hours to determine loss on ignition (LOI) before being mixed with the Li2B4O7/LiBO2 flux.Interference corrected spectra intensities were converted to oxide-concentrations against a calibration curve consisting of 30 international standards.The precision, expressed as the coefficient of variation (CV), is better than 0.5%.The accuracy, as measured on the international standards AGV-2, BHVO-2, BCR-2 and GSP-2 was better than 0.7% (1 RSD) (supplementary material III).

Eruption volume calculation
The minimum and/or maximum eruption volume of each volcano during each eruption period is derived from the ranges of thickness and surface areas that are reported in Campos and Rossi (1996) and Stewart and McPhie (2006).We converted these volumes to Dense Rock Equivalent (DRE) based on the magma type of different deposits.This analysis only includes the onshore deposits and results in a smaller estimate for larger pyroclastic volumes.The DRE volume is calculated using the equation of Crosweller et al. (2012): Tephra density is assumed to be 1000 kg/m 3 (Crosweller et al., 2012).Magma density varies depending on the magma type.
Here we used 2300 kg/m 3 for rocks with a SiO2 range of 65-77 wt.% and 2500 kg/m 3 for all samples with SiO2 < 65 wt.%.

40 Ar/ 39 Ar age results
In this section, we present our groundmass, biotite and amphibole 40 Ar/ 39 Ar results for eleven volcanic units of Milos.The 40 Ar/ 39 Ar ages range from 0.06 to 4.10 Ma and cover most of the major volcanic units of Milos.Table 2 and 3 show the 40 Ar/ 39 Ar results of incremental heating steps and single grain fusion analyses, respectively.Note that the Irr-ID column in these two Tables represents the irradiation ID of the analytical experiment (e.g.VU108-, VU110-) and the top right superscripts (G, B, A, O) in the sample IDs (e.g., G15M0029 G , G15M0021 B ) refer to groundmass, biotite, amphibole and obsidian.

Groundmass 40 Ar/ 39 Ar plateau and/or isochron ages
All groundmass samples yielding 40 Ar/ 39 Ar plateau and isochron ages with more than 50% 39 ArK and less than 2.5 MSWD included in their age spectrum are shown in Figure 4 and reported in Table 2.The 40 Ar/ 36 Ar isochron intercepts do not deviate from atmospheric argon at the 2-sigma level, unless stated otherwise (Table 3).Sample G15M0016 was collected from a dyke at Kleftiko in the southwest of Milos (Figure 2).Three incremental heating experiments were performed on the groundmass of this sample (Figure 5A).Two lava samples, G15M0019 and G15M0020, were collected from Kontaro in north-eastern Milos (Figure 2).Three replicate incremental heating step experiments on groundmass from sample G15M0019 (VU108-Z6a_4; VU108-Z6a_5 and VU108-Z6b_1, Figure 5B) were performed that are not reproducible.Their plateau ages range from 1.55 Ma to 1.62 Ma with relatively high MSWD (3.8-4.5),56-95% of the total 39 ArK, 34-53% of radiogenic 40 Ar, 0.88-1.02 of K/Ca and an atmospheric isochron intercept of 297-315.We consider the isochron age from the last experiment (VU108-Z6b_1) as the reliable age (1.48 ± 0.02 Ma, MSWD 0.44) because its MSWD value is the only one smaller than 2.5 in this experiment, and therefore the best estimate for the eruption age.Three replicate incremental heating step experiments on groundmass from sample G15M0020 (VU108-Z5a_5; VU108-Z5b_1 and VU108-Z5b_2, Figure 5C) were analysed.These experiments are similar at the lower temperature heating steps.They produced statistically meaningful plateau ages ranging from 1.52-1.56Ma with 41-62% of the total 39 ArK, 18-48% of radiogenic 40 Ar, 1.51-1.73 of K/Ca and an atmospheric isochron intercept of 295-300.Their combined weighted mean age is 1.54 ± 0.01 Ma (MSWD 3.06; 39 ArK 57.32%) with 25.31% of 40 Ar*.
The age of 1.825 ± 0.002 Ma is considered the best estimate for the eruption age of the Demeneghaki obsidian.The results shown in Figure 5 did not yield weighted mean plateau ages according to standard criteria including 39 ArK > 50%, but still provide some useful age information.Sample G15M0017 was collected from a cryptodome of the Profitis Illias volcano of southwestern Milos (Figure 2).Three replicate incremental heating experiments, VU108-Z7a, VU108-Z7a_4 and VU108-Z7b_1, have been performed on this sample which resulted in disturbed age spectra (Figure 6A).The consecutive lower temperature steps of all experiments define ages of <2.5 Ma, which is much younger than the ages of the submarine pyroclastic products of the lower series at Kleftiko and/or Profitis Illias (3.0-3.5 Ma, Fytikas et al., 1986 andStewart andMcPhie, 2006).At the consecutive higher temperature heating steps, these experiments yielded 3.64 ± 0.08 Ma ( 40 The amount of radiogenic 40 Ar of both the 40 Ar/ 39 Ar result from our sample and the K-Ar age data from previous studies (Fytikas et al., 1986) is rather low (<15%) for a sample of this age based on our laboratory experience.Therefore, the estimated age range for the oldest volcanic products of the Milos VF should be confirmed by other dating techniques.Sample G15M0015 is also a cryptodome breccia from Profitis Illias (Figure 2).Two replicate incremental step heating experiments were performed on the groundmass of this sample (VU108-Z9a and VU108-Z9b_1, Figure 6B).Experiment VU108-Z9a groundmass shows a disturbing age spectrum and ages increase from ~3 Ma in the initial heating steps to ~3.2

Groundmass
Ma, followed by a decrease to ~3 Ma in the high temperature heating steps.The consecutive heating steps only exist at the lower temperature steps yielding a "plateau" of 3.12 ± 0.02 Ma (MSWD 9.07).Due to the excess argon (  (Fytikas et al. 1986).
Sample G15M0029 is an andesite collected from Korakia in the northeast of Milos (Figure 2).Two incremental heating experiments (VU108-Z16a and VU108-Z16b_1, Figure 6C) were performed on this sample.Sample G15M0023 and G15M0024 are from the Triades lava dome northeast of Milos (Figure 2).A mafic enclave G15M0022

Single biotite grain
(host rock G15M0021) was collected from a lava near Cape Vani (Figure 2).The total fusion experiments of the biotites show that their initial 40 Ar/ 36 Ar estimates overlap with air (296-300).The total fusion ages gave the best estimates for their eruption ages of 2.10-2.13Ma using 22 out of 31 fusions with a range of radiogenic 40 Ar between 30-36% (Figure 7B).
Sample G15M0013 is from the rhyolitic Halepa lava dome in the south of Milos (Figure 2).The total fusion experiment (VU108-Z13, Figure 7C) on biotite of this sample produced a weighted mean age of 1.04 ± 0.01 Ma (MSWD 1.62; 9/10, 40  Sample G15M0034 and G15M0035 were collected from a lava dome located southeast of the Trachilas cone (Figure 2).Nine total fusion experiments (VU108-Z21, Figure 7C) were performed on biotite of sample G15M0035 and yielded the age of 0.63 ± 0.02 Ma (MSWD 1.26; 6/9; 40 Ar* 4.9%; inverse isochron age 0.77 ± 0.13 Ma).The atmospheric isochron intercept overlaps with air at 2-sigma (296.4 ± 1.7).The 4.9% of radiogenic 40 Ar is so low that we should consider the age of 0.63 ± 0.02 Ma with caution.For biotite of sample G15M0034 (VU108-Z20, Figure 7C) one total fusion experiment produced a weighted mean age of 0.51 ± 0.02 Ma (MSWD 0.95; 6/10; 40 Ar* 3.5%; inverse isochron age 0.61 ± 0.08 Ma) with an atmospheric isochron intercept.The age of 0.51 ± 0.02 Ma also needs to be considered as possibly suspect due to the low amount of radiogenic 40 Ar.Sample G15M0033 was collected from the Kalamos lava along the coast of the southwest of the Fyriplaka rhyolitic complex (Figure 2).Biotite of this sample (VU108-Z19, Figure 7C) yielded 0.412 ± 0.004 Ma (MSWD 1.10; 8/10; inverse isochron age 0.39 ± 0.02 Ma) with ~22.2% of radiogenic 40 Ar which is considered as the eruption age for the Kalamos lava.Sample G15M0007 was collected from the rhyolitic Trachilas complex in the north of Milos (Figure 2).Twenty-two total fusion (VU110-Z12, Table 3) and two incremental heating experiments (VU110-Z12a and 12b, Figure 8B) were performed on biotite of this sample.The total fusion experiments did not result in a reliable age due to the large errors of single steps (± 0.19 Ma on average) and the rather low amount of radiogenic Three pumice clasts (G15M0008-9 and G15M0012) were sampled from different layers of the Fyriplaka complex (Figure 2).

Multiple biotite grain
The first incremental step heating experiment on biotite from sample G15M0009 (VU110-Z23a, Figure 8C) gave negative ages at the lower temperature heating steps.Four consecutive higher temperature heating steps seem to define a "plateau" of 0.11 ± 0.02 Ma (MSWD 1.37) only using 18.33% of the total 39 ArK with 1.65% of radiogenic 40 Ar.The second experiment (VU110-Z23b) also yielded a "plateau" of 0.11 ± 0.03 Ma (MSWD 6.77) at higher temperature heating steps including 41.05% of the total 39 ArK and 3.13% of radiogenic 40 Ar.The significantly larger error of the isochron age may be due to the clustering of data close to zero on the y-axis.8C).The clustering of data points of experiment VU110-Z24a could result in the lower initial estimate of 40 Ar/ 36 Ar (285.98 ± 4.76).However, the combined age of 0.07 ± 0.01 Ma, using 43.53% of the total 39 ArK with an atmospheric isochron intercept (295.67 ± 7.39), could be the representative age of eruption.
Biotite of sample G15M0008 did not result in a reliable plateau in the first incremental step heating experiment (VU110-Z22a, Figure 8C) but shows a very disturbed age spectrum.The second experiment (VU110-Z22b) yielded 0.062 ± 0.003 Ma (MSWD 0.91) using 71.81% of the total 39 ArK with 2.69% of radiogenic 40 Ar as the best estimate of the eruption age. 40Ar/ 39 Ar multi-grain incremental heating plateau and/or isochron ages

Multiple amphibole grain
There are only two amphibole samples that yielded 40 Ar/ 36 Ar plateau and/or isochron ages (Figure 9A and B).Sample G15M0004 was collected from the pyroclastic series of Adamas from the PSLD (Fytikas et al., 1986), to the north of Bombarda (Figure 2).Two replicate heating experiments of G15M0004 amphibole (VU108-Z10_1 and VU108-Z10_2) were performed  The age in bold is considered as the best estimate of the eruptive age.
The 40 Ar* (%) is the average radiogenic 40 Ar of the analyses included in the weighted mean.
*The K/Ca ratio is calibrated by removing the total fusion with excess 37 Ar (Ca) (fA>1).

B
The experiment was analyzed on biotite of the sample.
The same steps were used for the calculation of isochron ages as used in the weighted mean ages.yielding 2.99 ± 0.11 Ma (MSWD 1.00; 39 ArK 87.31%,

Major element results
Major-element results are given in Table 4.The major element compositions range from 54 to 78 wt.% SiO2 (basaltic-andesiterhyolite to dacite-rhyolite, see Figure 10A).The most felsic samples (SiO2>75 wt.%) belong to the Fyriplaka and Trachilas complexes.Our data overlap with those of previous studies and display a similar range in SiO2-K2O (Francalanci and Zellmer, 2019 and reference therein).The samples of Polyegos are similar to the Fyriplaka and Trachilas complexes, whereas the older Milos samples overlap with Kimolos and Antimilos (Fytikas et al., 1986, Francalanci et al., 2007).
Although some samples of Antimilos are tholeiitic, all of the Milos volcanic units belong to the calc-alkaline and medium to high-K series (Figure 10B).A mafic inclusion, sample G15M0022, has high K2O (6%), similar to sample G15M0021 (7.2 wt.%).Both of them were collected from the Vani Cape area (Fig. 2).The SiO2 wt.% versus our 40 Ar/ 39 Ar ages diagram (Figure 11A) shows that there is a tendency of the volcanic units to become more felsic over time.In the diagram with K2O/SiO2 versus age there is no significant change (Figure 11C).

Variations of eruption volume with ages
Figure 11a shows the cumulative volcanic output volume of the Milos VF over time.This diagram shows that the Milos VF can be separated into three periods: Periods I (~3.3-2.13Ma) and III (1.48-0.00Ma) are characterised by low volcanic output volumes, whereas Period II (2.13-1.48Ma) shows a rapid increase in volcanic output volume.Period I and II are build up in submarine settings, whereas Period III is in a subaerial setting.The Milos VF was largely (~85% by volume) constructed in submarine before ~1.48 Ma (Period I and II) (Figure 11A).During Period III (1.48 Ma-present), only a small volume (~15%) of rhyolitic magma was added from different eruption vents.See the details of Period I-III in section 4.3.2.

Discussion
In this section, our 40 Ar/ 39 Ar results are compared with previously published geochronological data, and subsequently used to refine the stratigraphy of the Milos VF.In the last part, we will discuss the temporal variations in major elements and the volumetric volcanic output rate of the Milos VF.

Comparison with the previous geochronological studies on the Milos VF
K-Ar ages may show undesirable and unresolvable scatter due to various problems: (1) inaccurate determination of radiogenic argon due to either incorporation of excess argon or incomplete degassing of argon during the experiments; (2) inclusion of cumulate or wall rock phenocrysts in bulk analyses; (3) disturbance of a variety of geological processes such as slow cooling, thermal reheating; (4) unrecognized heterogeneities due to separate measurements of potassium and argon content by different methods; (5) requirement of relatively large quantities (milligrams) of pure sample (e.g.Lee, 2015).In addition to these methodological issues, in the case of Milos we observe that hydrothermal alteration caused substantial kaolinitisation, in particular the felsic volcanic samples, that most likely has affected the K-Ar systematics.Some of these issues are also valid for the 40 Ar/ 39 Ar method.However, the K-Ar method does not allow testing if ages are compromised. 40Ar/ 39 Ar ages only need isotopes of argon to be measured from a single aliquot of sample with the same equipment that can eliminate some of the problems with sample inhomogeneity.Furthermore, step heating and multiple single fusion experiments can shed light on sample inhomogeneity due to partial alteration effects.The high sensitivity of modern noble gas mass spectrometers for 40 Ar/ 39 Ar measurements results in very small sample amounts needed for analysis, that can yield more information on the thermal or alteration histories than larger samples.Moreover, other argon isotopes ( 36 Ar, 37 Ar and 38 Ar) can be used to infer some information about the chemical compositions (i.e.Ca and Cl) of samples.A high-resolution laser incremental heating method of 40 Ar/ 39 Ar dating allows us to resolve the admixture of phenocryst-hosted inherited 40 Ar in the final temperature steps of the incremental step heating experiments.More than half of our 40 Ar/ 39 Ar ages derived for this study are based on this method.All incremental step heating experiments are reproducible, except for the sample G15M0017 which gave the oldest age.The total fusion experiments of this study gave at least five times smaller analytical uncertainty (1SE on average ≤0.01 Ma) than the previous studies using conventional K-Ar (Angelier et al., 1977;Fytikas et al., 1976Fytikas et al., , 1986;;Matsuda et al., 1999) and SHRIMP U/Pb zircon methods (Stewart and McPhie, 2006).Fission track dating on obsidians of the Milos VF produced two ages (Bigazzi and Radi, 1981;Arias et al., 2006) which seems to overlap with the K-Ar and 40 Ar/ 39 Ar ages, but with larger uncertainty.U/Pb zircon ages could indicate the timing of zircon formation at high temperature (>1000 °C) in magma chambers significantly prior to volcanic eruption (e.g.Flowers et al., 2005).On the other hand, the lower closure temperature of K-rich minerals (<700 °C) makes the K-Ar and 40 Ar/ 39 Ar ages better suited to determine the timing of extrusion of volcanic products (e.g.Grove and Harrison, 1996;Cassata and Renne, 2013).
The MSWD value, as a measure of the scatter of the individual step ages, is based on the error enveloping around the data point.The decrease in error will automatically cause an increase in MSWD (e.g.York, 1968;Wendt and Carl, 1991).The MSWD values reported in this study are relatively high.In part this is caused by the fact that modern multi-collector mass spectrometers used for 40 Ar/ 39 Ar dating can measure the isotope ratios very precisely, which in turn would increase the MSWD.
It will be more valuable and challenging to find a plateau or isochron age which meets the MSWD criteria (<2.5) by modern multi-collector 40 Ar/ 39 Ar dating than by K-Ar or 40 Ar/ 39 Ar dating using a single detector instrument (e.g.Mark et al., 2009).
Potential drawbacks of the 40 Ar/ 39 Ar method are its dependence on neutron irradiation causing the production of interfering argon isotopes that need to be corrected for.The uncertainty in the ages of standards that are required to quantify the neutron flux also needs to be incorporated in the final ages as are uncertainties related to decay constants (supplementary material II).
Finally, recoil can occur during irradiation.Minerals such as biotite can be prone to recoil, yielding slightly older ages (e.g.Hora et al., 2010).
Figure 12 compares previous published K-Ar, U/Pb zircon and fission track ages from the same volcanic units with the new 40 Ar/ 39 Ar data of this study.In general, there is a good agreement, however, six ages out of twenty-three differ significantly from previous studies and will be discussed below.
The obsidian fission track ages (Bigazzi and Radi, 1981;Arias et al., 2006) for the Dhemeneghaki volcano are 0.25 My younger than the K-Ar ages (1.84 Ma, Angelier et al., 1977) and the 40 Ar/ 39 Ar age of this study (1.825 Ma,G15M0032B).The good agreement between the K-Ar and 40 Ar/ 39 Ar ages suggests that the fission track ages record another, lower temperature event, than the K-Ar and 40 Ar/ 39 Ar ages.In addition, the larger uncertainty of fission track ages (>0.05 Ma) also overlaps with the 40 Ar/ 39 Ar age at 2-sigma.We assume that the 40 Ar/ 39 Ar age is the correct extrusion age for the obsidian of the Dhemeneghaki volcano.Angelier et al. (1977) reported one dacite sample in the northwest of Milos with an age of 1.71 Ma (Angelier_3, location 3 on Figure 3 of Angelier et al., 1977).Argon loss could result in these ages (Angelier_3-5 in Figure 12) being younger than our 40 Ar/ 39 Ar groundmass ages of 1.97 ± 0.01 .
The amphibole of sample G15M0004 of the Adamas dacitic lava dome, located ~1 km north of rhyolitic Bombarda volcano, gave an inverse isochron age of 1.95 Ma ± 0.45 Ma.This age overlaps with the K-Ar age for the Adamas lava dome of 2.03 ± 0.06 Ma (dacite M 66) of Fytikas et al. (1986).The large analytical uncertainty of our sample G15M0004 is caused by a combination of low 40 Ar* yields and clustering of data points that define the inverse isochron showing excess argon was identified by the 40 Ar/ 39 Ar method ( 40Ar/ 36 Ar 319.51 ± 14.70), whereas the presence of excess argon cannot be tested by the K-Ar technique, implying that the Fytikas et al. (1986) might be slightly old.The Korakia andesite has an age of 1.59 ± 0.25 Ma (M 103, Fytikas et al., 1986) and was deposited in a submarine-subaerial environment on top of the Sarakiniko Formation that was dated based on paleomagnetic polarity in combination with a K-Ar age (1.80-1.85 Ma, Stewart and McPhie, 2003 and reference therein).The much older 40 Ar/ 39 Ar groundmass age (2.68 ± 0.01 Ma) of Korakia andesite sample G15M0029 is unreliable and it could indicate the emplacement age of the Kalogeros cryptodome (2.70 ± 0.04 Ma, Stewart and McPhie, 2006) or represents a geological meaningless age with only 23-27% of the total 39 Ar released in the "plateau".In this case, the K-Ar age of 1.59 ± 0.25 Ma is considered as the likely eruption age for the Korakia andesite although its argon loss or excess Ar component is unknown.
We obtained 40 Ar/ 39 Ar ages of 3.41-4.10Ma and 3.06 ± 0.02 Ma, respectively, from the groundmasses of dacite samples G15M0017 and G15M0015 in the southwest of Milos (Figure 2 and 13B).Both of these samples are derived from the coherent dacite facies of the rhyolitic Profitis Illias volcano based on the Figure 11 of Stewart and McPhie (2006).Sample G15M0015 yielded much higher radiogenic 40 Ar (41.77%) than that of sample G15M0017 (<10% of 40 Ar*), and the rhyolite sample M 164 from Fytikas et al. (1986) (23.5% of 40 Ar*) gave an estimate the eruptive age of 3.08 ± 0.08 Ma to the Profitis Illias volcano which is much younger than that given by our sample G15M0017 (Figure 12).Therefore, we consider our 40 Ar/ 39 Ar ages of 3.06 ± 0.02 Ma as the best estimate of the emplacement age of the coherent dacite facies of Profitis Illias volcano.
A basaltic andesite dyke near Kleftiko on the south-western coast of Milos has a K-Ar age of 3.50 ± 0.14 Ma which only gave 13.9% of 40 Ar* (Fytikas et al. 1986).This age is significantly older than the eruptive ages of Profitis Illias volcano which the dyke intruded (Stewart, 2003).Although containing relatively low 40 Ar* (16.87%), our 40 Ar/ 39 Ar age of 2.66 ± 0.01 Ma with 67.27% of 40 Ar* from the groundmass of basaltic andesitic sample G15M0016 of the dyke near Kleftiko is probably an accurate intrusion age.
Figures 13 and 14 summarize our new 40 Ar/ 39 Ar ages in combination with previously published stratigraphic, biostratigraphic, fission track, 14 C, K-Ar and U-Pb age data.We did not consider the Matsuda et al. (1999) data as the fission-track ages seem to be offset to other dating techniques ages obtained from the same deposits (see section 4.1 above).The exact start of volcanism in the Milos VF is still unclear since these older deposits are strongly hydrothermally altered.Van Hinsbergen et al. (2004) reported five ash layers in the Pliocene sedimentary rocks of southern Milos, ranging between 4.5-3.7 Ma in age, based on biostratigraphy, magnetostratigraphy and astronomical dating.In a slightly wider circle around Milos island, the 6.943 ± 0.005 Ma a1-tephra event recorded in several locations on nearby Crete (Rivera et al., 2011) shows that explosive volcanism along the Aegean arc, possibly on Milos, already occurred during the Messinian.These ash beds cannot be traced to currently exposed centres in the Milos VF and could conceivably be related to volcanic centres further north (Antiparos and Patmos), which were active during this time interval (Vougioukalakis et al., 2019).
Biostratigraphy shows that the youngest layer with dateable fossils (bio-event, the last common occurrence of Sphenolithus spp., Van Hinsbergen et al., 2004) in the Neogene sedimentary rocks is 3.61 Ma old (GTS2020, Raffi et al., 2020).The diatomite Unit II from Calvo et al. (2012) on top of the oldest volcaniclastic deposit from the north-eastern coast of Milos is constrained within 2.83-3.19Ma.These data suggest that the oldest products must be older than 2.83 Ma and younger than 3.61 Ma.Our oldest 40 Ar/ 39 Ar ages of this study displayed a wide range of 3.41-4.10Ma that is probably not correct due to alteration of the samples.Alteration might induce Ar loss and that would imply that the age is even older than 3.4-4.1 Ma.The age of 3.50 ± 0.14 Ma given by Fytikas et al. (1986) for an andesitic pillow lava or dyke has been discussed above and probably belongs to a series of basaltic andesite intrusions in the younger dacitic-rhyolitic deposits of Profitis Illias (~ 3.08 Ma, Fytikas et al., 1986), and therefore the 3.5 Ma age is probably not correct (e.g.Stewart, 2003).Fytikas et al. (1986) measured one sample from Kimolos (Figure 2 and 3) with an age of 3.34 Ma. Furthermore, Ferrara et al. (1980)

Periods with different volumetric output
The volume estimates of the Milos VF are hampered by limited exposure of several volcanic units and unknown age relationships.Therefore, not all units can be attributed to a certain volcano.Furthermore, we also do not know how much the volcanic products were lost through transport by air, sea currents and erosion.Therefore, the discussion here only provides a first order estimate of the onshore extruded magma volume.Taken into account all these limitations, our age data and the volume estimates by Stewart and McPhie (2006) indicate at least three periods of different long-term volumetric volcanic output rates (Qe) from ~3.3 to 0.0 Ma.We define a "Period" as a time interval were the Qe is significantly different from the average output rate (Qe average=1.0×10 - km 3 •yr -1 ) of the Milos VF over the last 3.3 Ma. Figure 11 shows that the Qe can be subdivided into two slow-growth periods (I and III) and one period (II) during which the Qe was significantly larger.
The lower boundary of Period I is based on our estimate of the oldest volcanic units of Milos at ~3.3 Ma.These oldest units were deposited in the southwest of Milos between ~3.3 and 3.08 Ma and include the BPS of Fytikas et al. (1986) and the felsic pumice cone/crypto dome facies of Stewart and McPhie (2006).These deposits have a minimum thickness of 120 m.The estimates of the DRE volume and Qe of these earliest volcanic deposits are hampered by the lack of precise age information, the high degree of alteration and structural complexities.Therefore, we only calculated the Qe of Period I from 3.08 Ma for which the eruption products are mainly dacitic-rhyolitic in composition (Table 5,Fig 11), and the first products that can be reliably dated are cryptodomes (3.06 Ma, sample G15M0015) and dykes (2.66 Ma, sample G15M0016) into the BPS of Fytikas et al. (1986) or the units of Profitis Illias volcano of Stewart and McPhie (2006, 3.08 Ma) in the southwest of Milos.These cryptodomes and dykes were followed by the formation of the submarine Fylakopi pumice cone volcano at 2.66 Ma (Stewart and McPhie, 2006) and Kalogeros cryptodome at 2.62 Ma (sample G15M0006) in the north-eastern part of Milos.These two pumice cone volcanoes contributed 3-11 km 3 DRE in volume to the Milos VF.The last two volcanic activities of Period I occurred in the southwest (Mavro Vauni, 2.50 Ma, Angelier et al., 1977) and west of Milos (Mavros Kavos, 2.36 Ma, this study), respectively, which produced two high-aspect-ratio andesitic-dacitic lava domes with a total volume of 1-3 km 3 DRE (Stewart and McPhie, 2006).During the submarine Period I, which lasted ~ 1.2 Ma, the estimated Qe is 0.9 ± 0.5×10 -5 km 3 •yr - 1 .The change from Period I to II is based on the sharp increase in Qe at 2.13 Ma (Fig. 11).During this period the Qe (3.0 ± 1.7×10 -5 km 3 •yr -1 ) increased by a factor of ~3 compared to Period I and III.Period II began with the submarine extrusions of the dacitic-rhyolitic Triades lava dome in the north-west and dacitic Adamas lava dome in the north-east of Milos and was followed by the rhyolitic Dhemeneghaki pumice cone/cryptodome and the Bombardo volcano in the north-east of Milos.For the Bombarda centre a large age range is reported in the literature (1.71-2.15Ma, Fig. 13B).We did not successfully date samples from the Bombarda centre, but Rinaldi and Campos Venuti (2003) reported that an age of 1.71 Ma is the best approximation based on other stratigraphic information.For the Dhemeneghaki centre, we obtained a 40 Ar/ 39 Ar age of 1.825 ± 0.002 Ma from obsidian.The Triades, Adamas, Dhemeneghaki and Bombarda centres all developed in submarine settings, as the intercalated sediments from the northern coast of Milos show (Calvo et al., 2012;Fig. 14).The last two volcanic expressions in Period II consist of two submarine-to-subaerial lava dome extrusions, Kantaro (1.59 Ma, Fytikas et al., 1987) and Korakia (1.48 Ma,this study) in the north-west and north-east of Milos, respectively.The products of these two centres are andesitic-dacitic in composition.All volcanic centres of Period II produced 8-30 km 3 DRE in volume for the Milos VF.
Period III began with a time interval of 0.4 Ma with no eruptions and has a very low Qe of 0.25 ± 0.05×10 -5 km 3 •yr -1 .The boundary between Period II and III can be placed at the last eruption of Period II, at the start of the first eruption in the low output interval, or halfway in between.The difference between those options is not significant, given the large uncertainties of the volume estimates (Fig. 12), and therefore we have decided to start Period III directly after the last eruption of the high Qe of Period II.The composition of nearly all Period III volcanic products is rhyolitic, an exception is the dacitic Plakes lava dome (Fig. 12).The Plakes lava dome is probably the last volcano erupting at ~0.97 Ma (Fytikas et al., 1987) in a submarine environment in the north of Milos, whereas the other lava dome in Period III, Halepa, produced rhyolitic lavas in a subaerial setting in the south (Stewart and McPhie, 2006).The Halepa and Plakes domes contributed 1-3 km 3 DRE in volume to the Milos VF and were followed by a 0.3 Ma interval with no or limited volcanic eruptions.Two subaerial pumice cone volcanoes with biotite bearing rhyolites were constructed during the last 0.6 Ma, the Trachilias and Fyriplaka complexes.The Trachilas complex was active for approximately 300 kyr (0.63-0.32 Ma) in the northern part of Milos.The evolution of this complex began with phreatic eruptions which became less explosive over time (Fytikas et al., 1986).During the last eruption (0.317 ± 0.004 Ma) of the Trachilas complex rhyolitic pumices filled up the crater area and did breach the northern tuff cone walls.The Trachilas complex only added a small volume (1-2 km 3 DRE) to the Milos VF.The Kalamos lava dome was also extruded in the south of Milos (Fig. 2) contemporaneously with the Trachilias complex.
The youngest volcanic activity of Milos (0.11 Ma-present) is characterized by subaerial eruptions of biotite phyric rhyolite from the Fyriplaka complex in the south of Milos, and was studied in detail by Campos Venuti and Rossi (1996).This complex is constructed on a paleosol that developed in a phreatic deposit ("Green Lahar", Fytikas et al., 1986) or lies directly on the metamorphic basement.Campos Venuti and Rossi (1996) indicated that the stratigraphic order is: Fyriplaka and Gheraki tuff rings, Fyriplaka lava flow, tuff cone of Tsigrado-Provatas.The total estimated volume of volcanic material is 0.18 km 3 DRE.
The boundary between the Fyriplaka and Tsigrado tuff cones is characterized by a marked erosive unconformity.The composition of these young volcanic products is very constant (Fig. 10-11), as noted by Fytikas et al. (1986) andCampos Venuti andRossi (1996).The products from Fyriplaka and Tsigrado cones are covered by a paleosol rich in archaeological remains and a phreatic deposit consisting largely of greenschist metamorphic fragments.According to Campos Venuti and Rossi (1996), the Fyriplaka cone was quickly built by phreatic and phreatomagmatic eruptions, as there are no paleosols observed between the different units.However, our data do suggest a large range in ages between 0.11 and 0.06 Ma. Fytikas et al. (1986) also reported a range between 0.14 and 0.09 Ma.These ages are inconsistent with the "Green Lahar" age of 27 kyrs (Principe et al., 2002), suggesting that the "Green Lahar" deposit consists of many different phreatic eruption layers that were formed during a time interval of more than 0.4 Ma, as the Kalamos lava is underlain by a green phreatic eruption breccia (Campos Venuti and Rossi 1996).We, therefore, conclude that phreatic eruptions occurred for more than 400 kyr, predominantly in the eastern part of Milos until historical times (200 BC -200 AD, Traineau and Dalabakis, 1989).

Temporal evolution of the magma flux and composition
Figure 11 shows temporal major-element variations during the evolution of the Milos VF.The volcanic units of Period III are dominantly rhyolitic in composition, whereas during Period I and II the compositions of volcanic units range between basalticandesite to rhyolite.However, the K2O/SiO2 ratio is constant (0.05 ±0.02) over the 3.3 Ma evolution of the Milos VF, with one exception, sample G15M0021 collected near Cape Vani which is altered by hydrothermal processes (e.g.Alfieris et al.

2013
). Period I and III contain large explosive pumice cone volcanoes, whereas Period II is dominated by effusive dome extrusions.The difference in volcanic structures is not observed in the SiO2 content and the K2O/SiO2 ratio of the volcanic products.
It is noteworthy that the value of the Qe (0.2-4.7×10 -5 km 3 •yr -1 ) for the Milos VF is at least 2-3 orders lower than the average for rhyolitic systems (4.0Í10 -3 km 3 •yr -1 ) and the mean for continental arcs (~70Í10 -3 km 3 •yr -1 ) (White et al., 2006).Milos overlaps with the lowest Qe values of the study of White et al. (2006).No data are available for the ratio between intruded magma in the crust below Milos and extruded volcanic units (I:E).White et al. (2006) argued that a ratio of 5:1 (I:E) is probably a realistic estimate for most volcanic centres and that this ratio can be higher in volcanic centres constructed on continental crust.A magma supply rate from the mantle beneath the Milos VF could be estimated in the order of 0.1-3.3Í10 - km 3 •yr -1 .
Druitt et al. ( 2019) reported a long-term average magma supply rate of approximately 1Í10 -3 km 3 •yr -1 beneath the Kameni islands of Santorini, which is comparable to that of the Milos.Besides the case of Santorini VF, no other information on the long-term average magma supply rate of other volcanic centres of the SAVA is available to our knowledge.
Milos is approximately 15 km long (W-E), a magma production rate of approximately 0.7-22 km 3 •km -1 •Ma -1 can be estimated over the last ~3.34Ma.Although this magma production rate per km arc length is the onshore estimate for the Milos VF, it is still significantly lower than for oceanic arcs: 157-220 km 3 •Ma -1 •km -1 (Jicha and Jagoutz, 2015).For continental arcs, the longterm magma production rate is more difficult to establish because magmatism is cyclic, and short periods (5-20 Ma) of intense magmatism ("flare ups") with 85 km 3 •km -1 •Ma -1 being alternated with periods of 25-50 Ma of low magma production rate of 20 km 3 •km -1 •Ma -1 (e.g.Jicha and Jagoutz, 2015).The periods of low magma production overlap with the magma production rates beneath the Milos VF over the past ~3.34Ma.

Conclusions
This study reports twenty-one new 40 Ar/ 39 Ar ages and major element data for 10 volcanic units of the Milos Volcanic Field.
In combination with previously published age data, geochemistry and facies analysis the following points can be made.
(1) The exact age of the start of volcanism in the Milos VF is still unclear due to the high degree of alteration of the oldest deposits.The best estimate based on our new 40 Ar/ 39 Ar ages, published K-Ar data and nannofossil biozones is between 3.5 and 3.15 Ma.
(2) Based on the long-term volumetric volcanic output rate, the volcanic history of the Milos VF can be divided into two slow growth periods, Period I (~3.3-2.13Ma) and III (1.48 Ma-present), and one relatively fast growth period, Period II (2.13-1.48Ma).
(3) Period I and II are characterised by andesitic to rhyolitic lavas and pyroclastic units, whereas those of Period III are dominantly rhyolitic.The K2O/SiO2 ratio is constant over the 3.3 Ma history of the Milos VF.
(4) The long-term volumetric volcanic output rate of Milos is 0.2-4.7Í10 - km 3 •yr -1 , two-three orders of magnitude lower than the average for rhyolitic systems and continental arcs.

Figure 1 .
Figure 1.Map of the South Aegean Volcanic Arc (SAVA).Red triangles indicate Volcanic Fields (VF): Susaki, Methana and Milos VFs in the western SAVA, Santorini VF in the centre and Nisyros VF in the eastern SAVA.Red contour lines show the depth to the Benioff zone (Hayes et al., 2018).The white arrow represents the GPS-determined plate velocity of the Aegean microplate relative to the African plate from Doglioni et al. (2002).

Figure 2 .
Figure 2. Distribution of the proximal and medial facies of the submarine pumice cone/crypto dome volcanoes, submarine, submarine-subaerial and subaerial domes and rhyolitic complexes (tuff cone and associated lava) of Milos, modified after Fytikas et al. (1986) and Stewart and McPhie (2006).The distal facies of Stewart and McPhie (2006) is not shown.

Figure 3 .
Figure 3. Simplified geological map of Milos with our 40 Ar/ 39 Ar ages and sample locations of key volcanic deposits, modified after Stewart and McPhie (2006) and Grasemann et al. (2018).The stratigraphic units of Milos are from Fytikas et al. (1986).Age data from this study are in black, published ages are shown in red (Angelier et al., 1977, Fytikas et al., 1986, Traineau and Dalabakis, 1989, and Stewart and McPhie, 2006).The "Green Lahar" (Fytikas, 1977) consists of deposits from multiple phreatic explosions and contains fragments of metamorphic, sedimentary and volcanic rocks.

Figure 6 .
Figure 6.Groundmass 40 Ar/ 39 Ar plateau or inverse isochron ages for samples G15M0017 (A), G15M0015 (B) and G15M0029 (C).Individual steps and final age calculation are reported with 1σ errors.The Profitis Illias volcano (A, B) and dacitic Korakia dome (C) are located in the south-western and north-eastern parts of Milos VF, respectively (Fig. 2).See the individual steps of sample G15M0015 and G15M0029 in supplementary material II.

Figure 7 .
Figure 7. Biotite 40 Ar/ 39 Ar total fusion ages for samples G15M0006 (A) and G15M0025-26(B, C), G15M0022-24 (D-F), G15M0013 (G) and G15M0033-35 (H-J).Data outside the shaded area are not included in the weighted mean.Individual steps and final age calculation are reported with 1σ errors.The Kalogeros cryptodome and Mavros Kavos lava dome are located in the north-eastern and south-western parts of Milos VF, respectively, and Triades lava dome, Halepa lava dome, Trachilias complex and the Kalamos lava are situated in the southern, northern and south-eastern parts of Milos VF, respectively (see Fig. 2).

Figure 9 .
Figure 9. Amphibole 40 Ar/ 39 Ar plateau or inverse isochron ages for samples G15M0004 (A) and G15M0026 (B).Final age calculation is reported with 1σ errors.The Adamas and Mavros Kavos lava domes are located in the northern and south-western parts of Milos VF, respectively (see Fig. 2).See the individual steps of sample G15M0004 and G15M0026 in supplementary material II.

Figure 10 .
Figure 10.SiO2 versus K2O (A) and AFM (B) diagrams for the Milos volcanic field with data of this study as solid circles.Published data are represented by shaded fields (Francalanci and Zelmer, 2019 and reference therein).Fields for the tholeiite, calc-alkaline, high-K calc-alkaline and shoshonitic series are from Peccerillo and Taylor (1976).Vertical lines defining fields for basalt, basaltic-andesite, andesite, dacite and rhyolite are from Le Bas et al. (1986).The solid line dividing tholeiitic and calcalkaline fields is from Irvine and Baragar (1971).

Figure 11 .
Figure 11.Eruption age versus (A) cumulative eruption volume for the volcanic deposits of Milos, (B) SiO2 wt.%, (C) K2O%/SiO2%, of Milos volcanic units of this study and previous studies.The maximum (Max; red line) and minimum (Min; dashed red line) cumulative eruption volume curves were estimated from Campos et al. (1996) and Stewart and McPhie (2006).Qe is the long-term volumetric volcanic output rate (see discussion).The exact volume of volcanic products between 4.1 and 3.08 Ma is not well constraint and indicated with a question mark.The major element data of the old pumices of Filakopi volcanoes (2.66 Ma) are from Stewart (2003).The major element data of the Plakes lava dome is from Fytikas et al. (1986).Geochemical data of the old pumices of the Profitis Illias (~3.08 Ma) is lacking due to the severe alteration.
reported an age of 3.15 Ma for a lithic clast derived from the Petalia intrusion in the Kastro volcaniclastics of Polyegos.If we assume that this reported age is a cooling age, volcanism in the Milos VF must have started before 3.15 Ma.Although age constraints for the start of volcanism on Milos both from the Neogene sedimentary rocks and the dated volcanic samples are poor, the evidence at this stage would suggest that volcanism in the Milos VF started ~3.3 Ma ago.

Figure 14 .
Figure 14.Diagram presenting three periods of different long-term volumetric volcanic output rate on Milos volcanic field based on the new 40 Ar/ 39 Ar data of this study and published data.The location of the different volcanoes is given in Fig 2 and indicated in the left panel (from left to right: SW, W, NW, N, NE, E, SE and S of Milos).The right panel corresponds to published age data: [A]=Fytikas et al., 1976, [B]=Angelier et al., 1977, [C]=Fytikas et al., 1986, [D]= Bigazzi & Radi, 1981, [E]=Matsuda, 1999, [F]=Stewart and McPhie (2006), [G]= Trainau and Dalabakis, 1989, and Biostratigraphic data of the Neogene sediments (NG) is from [H]=Calvo et al. (2012) and [I]=Van Hinsbergen et al. (2004) calibrated to Raffi et al. (2020) (LCO of Sphenolithus spp.and FO of D. tamalis).The number in the left panel represents the volcanic centres of Milos (see details in Table 5).The start of volcanism (3.08-3.61Ma) on Milos and the basement of the other Islands (Antimilos, Kimolos and Polyegos) are not well constrained and indicated with question marks (see text for discussion).The simplified basement cross-section (NS: Neogene sedimentary rock; MB: Metamorphic basement) under Milos volcanic units is based on Fytikas et al. (1989).We used the filled symbols as the best estimate for the eruption ages at the different volcanic centres, and the open symbols are not used as the best estimate due to their relatively large uncertainties.

Ar/ 39 Ar fusion and/or isochron ages
MSWD 0.99).Note that excess argon (40Ar/ 36 Ar 310.2 ± 4.0) is present.Hence the inverse isochron age is younger compared to the weighted mean age.The isochron age of 2.62 ± 0.04 Ma is considered as the best estimate for the emplacement age.Sample G15M0025 was collected from the Mavros Kavos lava dome located in the west of Milos (Figure2).The biotite of this sample (VU108-Z2, Figure7B) shows a weighted mean age of 2.36 ± 0.01 Ma (MSWD 0.70; 9/10;40Ar* 37.60%, inverse isochron age 2.34 ± 0.04 Ma) with an 40 Ar/ 36 Ar intercept of 300.6 ± 3.5.The age of 2.36 ± 0.01 Ma is considered the best eruption age estimate for this sample.
Results of nine single fusion experiments are given in Figure7.Nine or ten replicate single fusion experiments were conducted on 5-10 grains biotite per fusion.Sample G15M0006 is from dacite with columnar joints from the Kalogeros cryptodome in the northeast of Milos (VU108-Z11, Figure7A).The sample shows a weighted mean age of 2.72 ± 0.01 Ma for 9 out of 10 total fusion experiments (MSWD 1.95; 9/10) with an average 47.9% of radiogenic40Ar.The inverse isochron age is 2.62 ± 0.04Ma ( Ar* 26.3%; inverse isochron age 1.02 ± 0.04 Ma) with an initial 40 Ar/ 36 Ar estimate of 299.8 ± 4.1.The best estimate for the eruption age of the Halepa rhyolite is 1.04 ± 0.01 Ma.
40Ar (9.1%).On the other hand, the first incremental heating experiment produced a plateau age of 0.30 ± 0.01 Ma (MSWD 4.61; 39 ArK 56.60%; inverse isochron age 0.28 ± 0.05 Ma) including 14.51% of radiogenic40Ar.The second incremental heating experiment yielded a plateau of 0.317 ± 0.004 Ma (MSWD 1.29; 39 ArK 74.05%; inverse isochron age 0.31 ± 0.03 Ma) with a higher amount of radiogenic40 Ar (18.30%).The isochron intercepts of both incremental heating experiments are atmospheric.The second experiment is the best estimate for the eruption age, since it contained the largest amount of radiogenic40Ar and has a better reproducibility of single heating steps.

Table 3 . 40 Ar/ 39 Ar results of single grain fusion analyses on the Milos volcanic field.
.95 ± 0.45 Ma (MSWD 1.17;40Ar/ 36 Ar 319.51 ± 14.70) is considered the best estimate, but ideally this age should be checked by other techniques.Sample G15M0026 is from the same location as sample G15M0025, which gives us the opportunity to compare the biotite age with the amphibole age.One total fusion experiment on biotite (VU108-Z1b) yielded a weighted mean age of 2.35 ± 0.01 Ma (MSWD 1.36;40Ar* 38.6%).The atmospheric isochron intercept is low ( 40 Ar/ 36 Ar 292.01 ± 2.92), the inverse isochron age of 2.42 ± 0.04 Ma (MSWD 0.93) is considered the best result from the biotite.Two incremental heating experiments for amphibole (VU108-Z1b_1 and VU108-Z1b_2) gave plateau ages of 2.67-2.70Mawhicharemuchhighervaluesthan the biotite inverse isochron ages (2..This result could be caused by the high40Ar/ 36 Ar isochron intercepts (>320) with large uncertainties of ~29.Therefore, on the basis of the remarkable similarity of the two experiments, the combined inverse isochron age of 2.31 ± 0.28 Ma (MSWD 0.93, 39 ArK 71.36%,40Ar* 34.97%) is considered as the best estimate from amphibole which overlaps with the biotite age of 2.42 ± 0.03 Ma.This biotite age of 2.42 ± 0.03 Ma is considered to the best approximation of the eruption age.

Table 4 . Major-element composition of volcanic samples from the Milos Volcanic Field.
The classification of rock type for each sample is on the basis of field observation and SiO2 versus K2O plot of LeBas et al. (1986).All iron expressed as Fe2O3T(otal).