Delineation of shallow volcanic structures from audio-frequency magnetotelluric data beneath Ulleung Island, East Sea (Sea of Japan)

As volcanoes are closely related to the living environment of humans, in particular via natural hazards and geothermal energy, scientific studies on volcanic edifice structures are required. Ulleung Island is a Quaternary volcanic island located in the mid-western East Sea (Sea of Japan). In this study, we conducted an audio-frequency magnetotelluric (AMT) survey to image the substructure of the Ulleung volcanic edifice. In addition, the electrical structure was interpreted from the geothermal system perspective, geochemical compositions of volcanic rocks, and the possibility of the presence of a magma reservoir. AMT data were obtained from 25 stations and processed using the remote reference technique. Then, the three-dimensional (3-D) approach was reasonably adopted according to dimensionality analysis. Before conducting the 3-D inversion, the effects on topography and ocean were analyzed using a simplified 3-D synthetic model, because Ulleung Island is surrounded by sea and the topography is undulating. Most AMT stations on Ulleung Island are distorted by topographical and oceanic effects; in particular, oceanic effects are significant at frequencies lower than 10 Hz. The 3-D inversion was conducted with full impedance components and vertical magnetic transfer functions in a frequency range of 10,000–1 Hz. The results show that the Ulleung volcanic edifice is characterized by two layers of electrical structures as follows: near-surface resistive anomalies and an underlying conductive layer. Considering the high geothermal gradient on Ulleung Island, we suggest that the conductive layer of the volcanic edifice is due to hydrothermal alteration of basaltic rocks, with a potential heat source underneath. Based on the geochemical characteristics of the Ulleung volcanic rocks, the possibility of heat transfer from a trachytic magma reservoir within the shallow crustal is suggested. In summary, this study presents a geological interpretation of shallow volcanic structures beneath Ulleung Island and the possibility of an active magma reservoir as a potential heat source for the volcanic edifice.


Introduction
The East Sea (Sea of Japan) is a semi-enclosed back-arc basin that comprises three major sub-basins (Ulleung, Japan, and Yamato) with intervening topographic highs of crustal and volcanic origins (Chough et al. 2000).During the late Miocene-Quaternary, volcanic islands and seamounts formed along the northern Ulleung Basin.Ulleung Island is a Quaternary volcanic island located in the mid-western East Sea (Sea of Japan); it erupted from the Pliocene until the late Holocene.The most recent volcanic activity on the island was approximately 5 ka, according to K-Ar dating (Kim et al. 1999), and the island can be classified as an active volcano.As such, Ulleung Island has potential geothermal resources; however, volcanic gases and hot springs have not been directly observed on the surface.Therefore, borehole investigations were conducted to examine potential geothermal resources.In the temperature logging results of four boreholes drilled from 2011 to 2016 (MKE 2012;KEPRI 2015KEPRI , 2016)), the highest geothermal gradient was approximately 95 °C/km, which is higher than the general geothermal gradient of 25.1 °C/ km in Korea (Kim and Lee 2007).In addition, in a three-dimensional (3-D) full-waveform seismic tomography model, a low-velocity anomaly was observed at a depth of approximately 100 km beneath Ulleung Island (Simutė et al. 2016).These findings suggest that a hot shallow magma reservoir or intrusive rocks may still exist beneath Ulleung Island.
Based on various physical properties, geophysical techniques can be used to detect and identify the substructures of volcanic systems.The magnetotelluric (MT) method, a passive electromagnetic technique that uses natural time variations of the Earth's magnetic and electric fields (Vozoff 1990), was used to obtain the electrical structures beneath the Ulleung volcanic edifice.As electrical resistivity is sensitive to temperature, permeability, the presence of fluids, and clay alteration minerals (Manzella et al. 2004;Nurhasan et al. 2006;Gasperikova et al. 2015), it is commonly used to investigate the structural setting, geothermal setting, and hydrothermal circulation of volcanic areas (Uchida and Sasaki 2006;Aizawa et al. 2009;Ingham et al. 2009;Ogawa et al. 2014;Gresse et al. 2021).Thus, MT surveys can be used as a powerful tool to obtain substructure images beneath volcanic edifices.
Previous marine geophysical studies have been conducted extensively in the vicinity of Ulleung Island and Ulleung basin (Kim et al. 2005(Kim et al. , 2011(Kim et al. , 2022;;Park et al. 2006Park et al. , 2009;;Kim and Park 2010;Kim and Yoon 2017).These studies have mainly focused on the tectonic structure of the Ulleung basin and the eruptive history of Ulleung Island.Ground geophysical studies have been conducted on the Ulleung volcanic edifice using gravity surveys, magnetic surveys, electrical resistivity tomography (ERT), and MT surveys (Kwon et al. 1995;Kim andKim 2018, 2019;Lee et al. 2019).Kwon et al. (1995) interpreted the geological structure of the volcanic edifice using two-dimensional (2-D) forward modeling by adding shipborne gravity data to ground gravity and magnetic data.Kim and Kim (2018) interpreted the presence of faults and intrusive rocks in the volcanic edifice from a magnetic anomaly in a ground magnetic survey.In addition, craters and volcanic vents were estimated from the differences in magnetic anomalies according to the types of volcanic rock.Kim and Kim (2019) identified the caldera sedimentary layers and geological structure using ERT, focusing on the Nari caldera formed by volcanic eruptions.Meanwhile, Lee et al. (2019) conducted an MT survey to explain the substructure and fractures for geothermal development, and interpreted the electrical structure of the volcanic edifice through 3-D modeling and inversion.The results suggested that a shallow conductive anomaly in the east-west direction is due to the circulation of fluid through connected fractures.In addition, a conductive anomaly at a depth of 2-4 km was interpreted as a seawater intrusion or a deep-seated remnant heat source.
In this study, we used a high-frequency MT technique, audio-frequency magnetotelluric (AMT) surveying, to analyze on the shallow volcanic structures that could explain the geological structure of the volcanic edifice and the high geothermal gradient on Ulleung Island.The results of this study reveal the integrated interpretation of electrical structures using AMT data, the electrical structure of the geothermal system, and geochemical characteristics.First, a reasonable 3-D electrical structure was identified through the acquisition, analysis, and processing of AMT data, based on which the shallow volcanic structures were inferred.Second, the 3-D electrical structure was interpreted from the conceptual model of the geothermal system, and the possibility of the presence of a magma reservoir was determined.Finally, the geochemical characteristics of the volcanic rocks were evaluated to determine the nature of the magmatic heat source at shallow crustal depths.

Geological setting
Ulleung Island is approximately 130 km from the eastern coast of the Korean Peninsula (Fig. 1a); it has a pentagonal form with a part of the volcanic edifice exposed above sea level (a.s.l.).Mt.Seonginbong (984 m) sits at the center of the island (Fig. 1b), Nari caldera (with a diameter of 3.4 km across the topographical rim) is on the north side of Mt.Seonginbong, and the Albong lava dome/tephra cone complex is on the northwest side of the Nari caldera (Fig. 1c).The area is approximately 73 km 2 , and the lengths of the major and minor axes are approximately 12 and 8 km, respectively.However, the volcanic edifice extends to 2200 m below sea level (b.s.l.), and at this depth, the major axis is approximately 30 km.The Ulleung volcanic edifice forms an Aspite up to 500 m b.s.l.from the seafloor, but above it forms a Tholoide; it is estimated that approximately 30 parasitic cones are scattered between sea level and 1500 m b.s.l.(Won and Lee 1984).
The volcanic stratigraphy of Ulleung Island was first divided into five periods by Harumoto (1970): stratovolcano (Periods 1, 2, and 3), caldera (Period 4), and erupting pumices (Period 5).Subsequently, several researchers (Kim and Lee 1983;Won and Lee 1984;Kim 1985a;Min et al. 1988;Kim and Lee 2008) have divided the volcanic stratigraphy into five stages or five periods according to the occurrences and lithofacies of volcanic rock based on Harumoto (1970).Although the chemical composition of each eruption period varies, the Ulleung volcanic rocks are considered products of extensive fractional crystallizations (Won and Lee 1984;Kim 1985bKim , 1986;;Song et al. 1999;Choi 2021;Choi et al. 2022).They have compositions from alkali olivine basalt or trachybasalt-basaltic trachyandesite to trachyandesite-trachyte-phonolite, belonging to the alkaline series (Choi 2021).In this study, the volcanic stratigraphy was divided into four groups based on volcanic cycles and unconformities: basaltic rocks, Ulleung Group, Seonginbong Group, and Nari Group (Fig. 1c; Hwang et al. 2012).The Ulleung Group consists of volcaniclastic and lower-trachytic rocks.The Seonginbong Group comprises upper-trachytic and phonolitic rocks.The Nari Group consists of calderaforming volcanic rocks including the so-called Albong trachyandesites.
Although the original form of the Ulleung volcanic edifice has not been preserved owing to weathering, the caldera formed at the end of volcanic activity and Albong formed upon lava eruption and retain their original form (Kwon et al. 1995).The volcanic rocks below sea level are not well known; however, basaltic rocks are presumed to be widely distributed (Hwang et al. 2012).Faults are observed on a small scale and cause small displacements of the geological boundaries; these faults developed owing to tensile stress accompanied by later volcanism (Hwang et al. 2012).
Four borehole investigations were carried out in the eastern, western, southern, and northern parts of Ulleung Island (Fig. 1c), and the results are summarized in Table 1.Based on temperature logging, geothermal gradients of the four boreholes are higher than the typical regional geothermal gradient (25-30 °C/km), especially in the eastern and western parts.

Data acquisition and processing
AMT surveys were conducted in 2019 and 2020 at regular intervals throughout Ulleung Island.AMT data were acquired using an MTU-5C system (Phoenix Geophysics Ltd.).The system consisted of an ultra-wideband magnetotelluric receiver that provided simultaneous data measurements in a frequency range of 10,000-0.0002Hz.We measured the two horizontal components ( E x and E y ) of the electric field using Pb-PbCl 2 non-polarizable electrodes and three orthogonal components ( H x , H y , and H z ) of the magnetic field using MTU-150 induction coils.The time series data for the five measured components Fig. 2 Apparent resistivity of the sum of the squared impedances (SSQ impedance; grey dot) at each station and the regional average of 1-D impedance (black dot).At a frequency of 1.02 Hz, skin depth was estimated to be approximately 5 km were processed using the EMPower software (Phoenix Geophysics Ltd.).The recording times for AMT surveys were approximately 15 h (from 17:00 to 08:00 KST) to ensure high-quality data (including night time data) at each station.This is because the influence of local noise at night time is minimized and the magnetic field signal of the AMT dead band (5-1 kHz) is sufficiently strong (Garcia and Jones 2002).Impedance tensors were estimated in an electromagnetic field between 10,000 and 1 Hz using the robust processing of the remote reference technique to minimize local noise effects (Gamble et al. 1979).To obtain the data of horizontal magnetic field components for remote reference processing, Sawauchi station in 2019 (operated by Nittetsu Mining Consultants Co., Ltd. in Japan) and Jeju Island station in 2020 (obtained by Kangwon National University in Korea) were used (Fig. 1a).However, remote reference stations of Korea and Japan, which are far ≥ 4 times the skin depth from Ulleung Island, may have low correlation at high frequencies.To solve this problem, AMT data were obtained from two stations every day.Since the simultaneously measured AMT data showed a high correlation above approximately 3000 Hz, AMT data obtained at two stations on Ulleung Island were used simultaneously as remote reference data for each other.
To obtain the regional one-dimensional (1-D) impedance and skin depth, we used the sum of the squared elements of the impedance tensor (SSQ impedance; Rung-Arunwan et al. 2016).This was because SSQ impedance is less affected by downward bias due to distortion than the determinant of impedance proposed by Berdichevsky et al. (1980).Regional 1-D impedance can be obtained by geometric average of the SSQ impedance.Figure 2 shows the apparent resistivity of the SSQ impedances at each station and regional average of 1-D impedance.The apparent resistivity of individual stations showed a consistent trend; the apparent resistivity decreased toward lower frequencies and increased below 5 Hz.Skin depth was calculated from regional 1-D impedance.At a frequency of 1.02 Hz, apparent resistivity was 98 Ωm and the skin depth was approximately 5 km.

Dimensionality analysis
A dimensionality analysis of AMT data was conducted before AMT modeling and inversion to determine whether a 1-D, 2-D, or 3-D approach was required.In general, 2-D interpretation of 3-D MT data can lead to reasonable geological interpretation; however, proper dimensional interpretation is important because it may not be acceptable in some cases (Ledo 2005).Therefore, the WALDIM code (Martí et al. 2009) was applied to evaluate the dimensionality of the AMT data.
WALDIM, a FORTRAN code that conducts dimensionality analysis of MT data, is a numerical approach based on WAL invariant criteria (Weaver et al. 2000); it analyzes seven independent ( I 1 ,I 2 ,I 3 , I 4 , I 5 , I 6 , I 7 ) invari- ants and one dependent ( Q ) invariant.A set of eight invariants was defined such that the invariants could be represented in a Mohr circle diagram, except for two ( I 1 and I 2 ).Invariants I 1 and I 2 are non-dimensional and normalized to unity, and their vanishing has a physical interpretation, specifically related to geoelectric dimensionality.They also serve to normalize the remaining invariants and provide information on the 1-D magnitude and phase of geoelectrical resistivity.Invariants I 3 to I 7 and Q can be used to determine the dimensionality type and identify galvanic distortion.
Figure 3 shows the results of the WALDIM application for the four bands at each AMT station.In this study, a threshold value used to distinguish between types of dimensionality was tested on real data, where a value of 0.15 was used.The error of the MT tensor components was set to 5%.The results of the dimensionality analysis varied from 1-D to 3-D for the four bands.The high-frequency (10,000-1000 Hz) and low-frequency (10-1 Hz) bands were predominantly 3-D, except for a few 1-D bands and anisotropy (Fig. 3a and d).Anisotropy occurs in inferred 2-D cases when the strike direction is mismatched or the relationships between the tensors do not imply isotropic media (Martí et al. 2009).The intermediate bands (1000-100 Hz and 100-10 Hz) were generally 3-D, except for complex patterns near the coastline and caldera rim (Fig. 3b and c).In short, these results suggested that a 3-D approach was necessary for a consistent interpretation of the entire dataset.

Topographical and oceanic effects
Ulleung Island is a small island surrounded by sea and with undulating topography; therefore, the influences of topographical and oceanic effects cannot be excluded.Accordingly, the AMT response was analyzed using a simplified 3-D synthetic model of the island to identify topographical and oceanic effects on AMT data.To obtain the 3-D computation domain, topography was obtained from the Shuttle Radar Topography Mission (SRTM; 30 m resolution) and bathymetry was obtained from a standard bathymetry map (ETOPO1; 450 m resolution).Simplified land and ocean models were used as synthetic models for the AMT response analysis.The ocean model included topography and bathymetry data, with seawater resistivity of 0.33 Ωm and homogeneous underground resistivity of 100 Ωm.For the land model, the resistivity value corresponding to seawater was replaced by 100 Ωm from the ocean model.The AMT forward response was obtained using the Geotools software (CGG) based on the finite integration technique applied to orthogonal Cartesian grids.
The entire computation domain covered an area 98 × 95 km (east and north directions, respectively), corresponding to Ulleung Island and the surrounding ocean.Figure 4a shows the core region and part of the external cell in the 3-D computation domain, including topography and bathymetry data.The core region covered an area of 14 × 11 km (east and north directions, respectively).Along the horizontal direction, the entire mesh was discretized into 85 × 70 cells.The core mesh was constant at 200 × 200 m for 55 × 40 cells and gradually increased for the external cells.The vertical direction was discretized into 103 layers, and layers of ~ 1 km from the 1st to the 37th layer were constant at 25 m; subsequent layers were increased by a factor of 1.05 to a depth of 2 km and 1.1 to a depth of 8 km.Later layers were increased as the factor gradually increased to a maximum depth of 544 km.
Figure 4b shows a plan view of the AMT stations, topography, and bathymetry on the mesh area.The locations of the stations used in the synthetic model were the same as the field survey stations, and the two stations (sites A and B) were analyzed according to the distance from the coastline and topographical slope.Site A was located close to the coastline to the north, with the upslope to the south.The distance from the coastline was approximately 0.3 km.Site B was located on the Nari caldera with flat topography surrounded by the caldera rim; the distance from the coastline was approximately 2.7 km. Figure 5 shows the apparent resistivity, phase, and vertical magnetic transfer function (or tipper) for 17 frequencies between 10,000 and 1 Hz, considering the topographical and oceanic effects in the synthetic model.The impedance tensor represented the off-diagonal components ( Z xy and Z yx ), and tipper represented real and imaginary parts for T x and T y components.For Site A close to the coastline, the impedance tensors between the land and ocean model differed below about 30 Hz.The tipper was distorted by the oceanic effect from about 1000 Hz to lower frequencies.Here, T x showed an oce- anic effect from a relatively higher frequency than T y , because Site A is close to the northward conductive seawater.If the station was located east or west of Ulleung Island, the ocean effect would be greater at T y than at T x .For Site B located in the southern part of the Nari caldera, the impedance tensors between the land and ocean models were almost identical in all frequencies.However, the tipper was distorted by oceanic effects below 10 Hz.Considering both sites, the oceanic effect on the AMT data of Ulleung Island affected the tipper more than the impedance tensor.Taken together, in the 3-D inversion of Ulleung Island, the 3-D computation domain including seawater should be considered for both full impedance tensor and tipper.

3-D inversion
In a previous MT study by Lee et al. (2019), 3-D inversion was conducted using the off-diagonal components of the impedance tensor, and the ocean effect was considered, but topography and bathymetry data were not utilized.We used topography and bathymetry data to consider topography and oceanic effects in the AMT data and conducted 3-D inversion using full components of the impedance tensor and tipper to delineate shallow volcanic structures.
The 3-D resistivity model was obtained using the ModEM inverse algorithm (Egbert and Kelbert 2012;Kelbert et al. 2014) with the nonlinear conjugate gradient method (NLCG) at the 10,000-1 Hz frequency band.The computational domain for inversion was the same as the mesh of the synthetic model (Fig. 4).The 3-D inversion at 54 frequencies was conducted with full components of impedance ( Z xx , Z xy , Z yx and Z yy ) and tipper ( T x and T y ).The starting model used for inversion was a homogenous model of 100 Ωm for land area, and seawater was fixed at 0.3 Ωm.This model referred to the previous MT study (Lee et al. 2019) and was set through the regional average of SSQ impedances.Model regularization employed a smoothing parameter of 0.3 for all directions.The initial normalized root mean square (nRMS) was 5.3 while the final one was 1.1 after 54 iterations using the error floor.The error floor was set to For tipper, the error floor was set to 3% for T x and T y .
Figure 6 shows the final nRMS values at each AMT station, which were around 1 except for a few stations that exhibited poor data fit due to local noise.
Figure 7 shows the observed data used for inversion, and calculated AMT responses of the land and ocean models at Sites A and B. For the impedance tensors, the AMT responses calculated from the land and ocean models show differences below 10 Hz, especially in Site A adjacent to the seawater.In the case of tipper, there was a difference between the land and ocean models below 10 Hz.At Site A, in particular, the difference below 10 Hz was clear.In summary, the 3-D inversion of the ocean model including the seawater provided a satisfactory geoelectric structure for delineating shallow volcanic structures.

3-D resistivity model
The final 3-D electrical resistivity model is shown in Figs. 8 and 9. Figure 8 shows cross-sections in the eastwest (E-W) and north-south (N-S) directions, while Fig. 9 shows four horizontal plan views between sea level and 1 km b.s.l.The main electrical structures observed in the cross-section of this model were classified into two layers (Fig. 8): (1) near-surface resistive anomalies above sea level (R1 and R2); (2) a conductive layer between sea level and 1 km b.s.l.(C1 and C2).
The near-surface formation consisted of mostly the Ulleung, Seonginbong, and Nari groups, mainly composed of trachytic rocks throughout Ulleung Island and trachyandesitic rocks in the Nari caldera.For the N-S cross-section (Fig. 8b), the two areas of high resistivity (R1 and R2) were divided by the caldera rim, and these features represented trachyandesitic rocks (R1) in the Nari caldera and trachytic rocks (R2) in the south.The E-W cross-section (Fig. 8c and d) showed that the Ulleung volcanic edifice above sea level was mostly composed of high-resistivity trachyandesitic (R1) and trachytic rocks (R2).The high resistivity at the Nari caldera appeared continuously up to approximately 500 m b.s.l.(Fig. 9a-c).This may be associated with the volcanic activity of the Nari caldera and Albong.
Considering that the collapse depth of the caldera is at least 600 m (Hwang et al. 2012), it was interpreted as a resistive zone caused by the presence of caldera sediments and trachyandesitic rocks (R1).In addition, considering that the high resistivity (R2) in the southern Nari caldera extended to the deep part, it is likely associated with the volcanic activity that produced the trachytic rocks.In the horizontal plan view, the high resistivity observed north of the Nari caldera is not prominently evident near the surface, but exhibits a continuous presence to an elevation of 800 m b.s.l.The present asymmetric geometry of the Nari caldera was formed through repeated eruptive episodes occurring in the northern region during the last 19,000 years (Kim et al. 2014).Therefore, this high resistivity is likely related to volcanic activity associated with the caldera collapse.
The conductive layer between sea level and 1 km b.s.l. was mainly located outside the Nari caldera and was distributed throughout Ulleung Island.These features appeared continuously down to 500 m b.s.l. in a ring shape similar to the coastline (Fig. 9c), and almost disappeared at 800 m b.s.l.(Fig. 9d).The thickness of this anomaly is approximately 500 m and it is most prominent at 200 m b.s.l.(Fig. 9b).In addition, along the E-W cross-section (Fig. 8c), the low-resistivity anomalies (C1 and C2) likely exhibit connectivity southward at a depth slice of 500 m b.s.l.(Fig. 9c).Considering that these anomalies are volcanic islands with high geothermal gradients, they are likely associated with seawater intrusion or convection of hot water.
We conducted a sensitivity test to evaluate the robustness of the low-resistivity anomalies (C1 and C2) that were clearly visible in the inverted 3-D resistivity model.We replaced C1 and C2 with a resistivity value of 100 Ωm and performed forward modeling for C1 and C2, individually.Figure 10 shows the modified model used for the sensitivity test and the resulting AMT responses of the inversion and modified model.There were differences in the AMT responses for both anomalies.The C1 anomaly, which is relatively close to the surface, shows a difference at a frequency of 3000 Hz and lower, and the C2 anomaly, which is located deeper, shows a difference at a frequency of 30 Hz and lower.Therefore, the sensitivity test proved that the inverted 3-D resistivity model was sensitive to C1 and C2.

Geoelectric structure of the geothermal system
Geothermal systems around the world are classified into four groups based on geological and tectonic relationships: volcanism and tectonics, continental collision zones, continental rift systems, and active volcanism (Chandrasekharam and Bundschuh 2008).Although volcanic gases and hot springs are not directly observed, the existence of a geothermal system is likely owing to the new magma injection that formed Albong, the high geothermal gradient, and the eruption cycle, which is an active volcano with a subvolcanic magma plumbing system (Kim et al. 2022).Therefore, the electrical structure from the AMT data on Ulleung Island can be interpreted from the perspective of a geothermal system.A standard approach for the exploration of geothermal systems is to describe them in conceptual models (Leeuwen 2016), which consist of a qualitative description of at least three elements: (1) heat source, (2) reservoir, and (3) cap rock.These structures can be identified through AMT data, where the resistivity is directly related to the parameters characterizing geothermal reservoirs, such as temperature, alteration, salinity, and porosity (Hersir and Björnsson 1991).There is generally a conductive cap rock over the resistive region where the reservoir is located (Spichak and Manzella 2009).This suggests that AMT surveys are an effective tool for identifying conductive cap rocks and detecting potential reservoirs.
Figure 11 shows the lithology distribution by depth of the four boreholes.Most of the rock types consist of basaltic rocks in GH-1 (east) and GH-2 (west), as well as trachytic rocks in GH-3 (south) and GH-4 (north) (MKE 2012; Lee and Lee 2016).The 3-D resistivity model shows low-resistivity anomalies at boreholes GH-1 and GH-2, and relatively high resistivity at GH-3 and GH-4.These conductive anomalies might be associated with seawater intrusion, hydrothermal solution, and hydrothermal alteration, considering that Ulleung Island is surrounded by the sea; a high geothermal gradient and basaltic rocks are to be expected below sea level.Basaltic rocks in GH-1 and GH-2 show argillic alteration due to hydrothermal alteration, without distinguishing lava flow and volcanic breccia (MKE 2012).The temperature logging results also show a high geothermal gradient (Fig. 12), providing further evidence that the alteration zone formed in a high-temperature environment.If the electrical structure is interpreted based on these facts, the conductive layer below the sea level may be associated with the cap rock.As cap rocks are generally clay, the presence of clay was determined by physical parameters such as relatively low resistivity values in the range of 1-100 Ωm (Musset and Khan 2000).This suggests that a high-temperature heat source can form a hydrothermal solution beneath basaltic rocks.This heat source may have been a magma reservoir or intrusion beneath the Ulleung volcanic edifice.
Reservoirs and heat sources were difficult to confirm owing to the limited resolution of the AMT data for deeper targets; however, the presence of cap rocks and a high geothermal gradient support the presence of a heat source.Furthermore, in a numerical heatconduction model, assuming a magma reservoir under Ulleung Island, the presence of a hot magma reservoir in a molten or solidified state under the volcanic edifice was reported by Lee et al. (2016).This is explained by the heat transfer from the heat source, which is consistent with the electrical structure results.In addition, although there was no resolution below 1 km b.s.l. in our resistivity model and it differs from previous MT studies (Lee et al. 2019), the presence of a deep-seated remnant heat

Geochemical characteristics of Ulleung volcanic rocks
Subaerial volcanic activity on Ulleung Island has been divided into five lithostratigraphic stages (Kim et al. 1999(Kim et al. , 2014;;Hwang et al. 2012;Brenna et al. 2014;Choi 2021).Basal Stage 1 occurred between 1.37 and 0.97 Ma, with eruptions of mainly basaltic lava and agglomerate.Stages 2 and 3 consisted of trachytic rocks and occurred between 0.83 and 0.77 Ma, and between 0.73 and 0.24 Ma, respectively.Stage 4 occurred between ~ 19 and 5.6 ka, and involved phonolitic-trachytic pyroclasts deposited during explosive caldera-forming eruptions.Stage 5 is the post-caldera activity within the Nari caldera, which commenced at ~ 5 ka and has composed of leucite-bearing trachyandesite lava and spatter (the socalled Albong lava).Caldera-forming eruptions brought to the surface abundant coarse-grained or intrusive lithic fragments of a cogenetic suite of extrusive rocks (Brenna et al. 2014).
Representative major and trace element concentrations for the five stages of the Ulleung volcanic rocks are shown in Fig. 13 against SiO 2 as an index of differentiation.MgO, CaO, FeO*, TiO 2 , Ni, Co, and Cr contents decrease with increasing SiO 2 content, indicating fractionation of ferromagnesian phases such as olivine, clinopyroxene, kaersutite, and titanomagnetite.Al 2 O 3 , Ba, and Eu contents increase up to ~ 55 wt.% SiO 2 , followed by a decrease, indicating feldspar fraction at ~ 55 wt.% SiO 2 .Immobile incompatible elements, such as Nb, increase with increasing SiO 2 content.Meanwhile, for Sm, a rare earth element compatible with amphibole, there is an offset between Stage 1 and differentiated stages 2-5, indicating the presence of at least two discrete magmatic batches during the volcanic evolution of Ulleung Island.Furthermore, Stage 4 phonolites have much lower MgO, CaO, TiO 2 , Ba, and Eu contents, and higher Sm and Nb contents than stage 2 and 3 trachytes at a given SiO 2 content, implying magma differentiation at different depths (Brenna et al. 2014;Choi 2021).It is noteworthy that the Stage 5 volcanic rocks are less silicic than those of Stages 3 and 4, suggesting that a new magma batch might have formed after the explosive volcanic activity that formed the Nari caldera in Stage 4.
In addition, the application of clinopyroxene-melt thermobarometry for submarine Ulleung volcanic rocks shows two magma reservoirs at different depths: sub-Moho (~ 17-32 km depth) for basaltic rocks and shallow crustal (~ 2-9 km depth) for trachytes (Choi et al. 2022).These may have been independent reservoirs associated with two discrete magmatic batches, as shown above.Although these results cannot be directly compared with the 3-D AMT inversion model because of the low resolution below 1 km b.s.l., there is likely a heat source for these magma reservoirs within shallow crustal depths below Ulleung Island.

Conclusions
An AMT survey was conducted on Ulleung Island to reveal the substructure of the volcanic edifice from its 3-D electrical properties.AMT data were obtained from 25 stations, and impedance tensors were estimated in an electromagnetic field between 10,000 and 1 Hz.For the analysis of AMT data, the 3-D approach was reasonably adopted from dimensionality analysis.The 3-D synthetic forward modeling was conducted to analyze the effects of topography and the ocean.Topographical and oceanic effects mostly affected all stations; in particular, oceanic effects were common at frequencies lower than 10 Hz.As a result, the inverted 3-D resistivity model obtained from the 3-D inversion was interpreted from the electrical structure of the Ulleung volcanic edifice, a conceptual model of the geothermal system, and the geochemical characteristics of Ulleung volcanic rocks.
From the cross-sections and depth slices of the inverted 3-D resistivity model, the volcanic edifice can be interpreted as follows: (1) near-surface resistive anomalies above sea level, and (2) a conductive layer between sea level and 1 km b.s.l.The near-surface formation, composed mostly of trachytic and trachyandesitic rocks, is a relatively high resistivity layer that extends continuously from the Nari caldera and south of the caldera to 500 m b.s.l.This may be associated with the volcanic activity of the Nari caldera and Albong.The conductive layer under the near-surface resistive layer is associated with the cap rock from the perspective of the geothermal system, demonstrating the possible presence of a heat source under this layer.However, it is difficult to interpret the electrical structure at greater depths ( ≥ 1 km b.s.l.) because of the limitation of AMT resolution.Nevertheless, the AMT results are consistent with geochemical evidence showing at least two discrete magmatic batches, possibly at different depths: sub-Moho for basaltic rocks and shallow crustal for trachytes.The trachytic magma reservoir may play a role as a heat source for hydrothermal alteration, resulting in a cap rock.In addition, the heat source below the conductive layer inferred in this study is consistent with the conductive anomaly of a deep-seated remnant heat source in a previous MT study (Lee et al. 2019).
In conclusion, the 3-D electrical structure of Ulleung Island reveals a conductive layer below sea level that is consistent with borehole and temperature logging data.This structure can be explained based on the high geothermal gradient and geochemical characteristics of volcanic rocks, and it is concluded that a hot magma reservoir still exists beneath the volcanic edifice.However, the existence of this magma reservoir is only inferred, and the specific deep structure was not revealed.Therefore, long-term measurements for obtaining high-quality data below 1 Hz will be required to image deep structures on Ulleung Island.

Fig. 1
Fig. 1 Study area maps.a Location of Ulleung Island in the East Sea (Sea of Japan) and remote reference stations on Jeju Island, Korea, and in Sawauchi, Japan.b Bathymetry and topography of the study area.c Simplified geological map of Ulleung Island (modified from Hwang et al. 2012)

Fig. 5
Fig. 5 Comparison of apparent resistivity and phase of land (dot) and ocean (solid line) models for Site A and Site B. The impedance tensor represents the off-diagonal components ( Z xy and Z yx ), and tipper represents real and imaginary parts of the T x and T y components.The AMT data are distorted by the influence of topography and the ocean at all stations, especially close to the coastline 20% of Z xy * Z yx for the diagonal components Z xx and Z yy , 10% of Z xy for the off-diagonal component Z xy , and 10% of Z yx on the off-diagonal component Z yx .

Fig. 6
Fig. 6 Distribution of final normalized root mean square (nRMS) value at each AMT station.For most stations, the value is around 1

Fig. 7
Fig. 7 Apparent resistivity, phase, and tipper curves of 3-D inversion at Site A and B. Dots indicate observed data used for inversion, the green line indicates calculated data from the land model, and the blue line indicates calculated data from the ocean model.The AMT response of the ocean model provided satisfactory geoelectric structure

Fig. 8
Fig. 8 Resistivity model cross-sections.a Map view of vertical cross-sections: b north-south section across Albong, c east-west section across Albong, and d east-west section across GH-1 and GH2 in the 3-D resistivity model.Red circles and black lines indicate the locations of borehole investigations.Yellow triangles indicate Mt.Seonginbong and Albong

Fig. 9
Fig. 9 Horizontal plan view of the inverted 3-D resistivity model below sea level.a 0 m, b 200 m, c 500 m, and d 800 m.Red circles indicate the locations of borehole investigations.Yellow triangles indicate Mt.Seonginbong and Albong, and the white dotted line indicates the caldera rim

Fig. 10
Fig. 10 Sensitivity tests of low-resistivity anomalies: C1 and C2. a Modified model in which low-resistivity anomalies were replaced with a block of 100 Ω•m.AMT response before and after sensitivity test for b C1 and c C2

Fig. 13
Fig. 13 Representative major oxide and trace element variations in Ulleung volcanic rocks.Data from Choi (2021) and references therein.Total Fe as FeO*

Table 1
Results of borehole investigations and temperature logging conducted on Ulleung Island