Magmatic–hydrothermal system of Aso Volcano, Japan, inferred from electrical resistivity structures

We investigated the magmatic–hydrothermal system of Aso Volcano, Japan, using broadband magnetotelluric (MT) data. To establish the nature of the shallow crust, a previous resistivity model based on data from 100 measurement sites in and around Aso volcano was revised using data from 9 additional sites near Naka-dake crater, which is located in the central part of the volcano. The components of MT impedance and the tipper vector were used to obtain the resistivity structure by three-dimensional inversion. The resistivity structure shows a subvertical low-resistivity (< 1 Ωm) column-shaped body beneath Naka-dake crater that extends from − 600 m to 10 km below sea level (BSL) and dips steeply to the north-northeast. The position of the upper part of the column is displaced eastward compared with the previous model and does not overlap the position of the presumed magma reservoir inferred previously from seismic and geodetic observations underneath the western side of Naka-dake crater at a depth of 5 km. We interpret this low-resistivity column to be a magmatic–hydrothermal system composed of brine and magma that were transported to Naka-dake crater from the main deep-seated magma reservoir. A horizontal low-resistivity (< 10 Ωm) layer occurs beneath post-caldera cones at the depths of 0–2 km BSL, and this layer extends laterally from the upper part of the low-resistivity column. We interpret this low-resistivity layer as representing a shallow hydrothermal system that has developed around the central column-shaped magmatic–hydrothermal system.


Introduction
Determining the source and location of a magmatichydrothermal system is important for understanding volcanic activity. Geophysical methods such as seismic data analysis, crustal deformation monitoring, and magnetotelluric (MT) surveys are effective tools for investigating the existence of magma within a volcano and associated hydrothermal systems (e.g., Nakamichi et al. 2007;Heise et al. 2010;Siena et al. 2010;Seccia et al. 2011;Lin et al. 2013;Luttrell et al. 2013;Farrell et al. 2014;Aizawa et al. 2014;Ogawa et al. 2014).
Various geophysical studies of Aso Volcano in Japan have suggested the presence of magma reservoirs, including on the basis of seismic data (Sudo and Kong 2001;Sudo 2004, Abe et al. 2010;Abe et al. 2016Huang et al. 2018, gravity data (Komazawa 1995;Miyakawa et al. 2016), and leveling and global navigation satellite system (GNSS) data (Sudo et al. 2006;Ohkura and Oikawa 2012). Hata et al. (2016Hata et al. ( , 2018 investigated the spatial relationships between the resistivity structure obtained by broadband MT data and the regions beneath Aso Volcano thought to be occupied by magma as inferred from existing seismic and leveling/GNSS data.

Open Access
*Correspondence: matsushima-n@aist.go.jp 1 Geological Survey of Japan, National Institute of Advanced Industrial Science and Technology, Central 7, Higashi 1-1-1, Tsukuba, Ibaraki 305-8567, Japan Full list of author information is available at the end of the article Matsushima et al. Earth, Planets and Space (2020) 72:57 Although MT surveys had previously been conducted at Aso Volcano (Handa et al. 1998;Takakura et al. 2000;Asaue et al. 2005;Utsugi et al. 2009), Hata et al. (2016Hata et al. ( , 2018 were the first to analyze the entire area of Aso volcano at depths of 0-20 km. Hata et al. (2016Hata et al. ( , 2018 identified an extremely low-resistivity (< 1 Ωm) columnar body at depths of 0-15 km below sea level (BSL). The low-resistivity column overlaps with the magma reservoirs. The sill-like inflation source of crustal deformation detected using GNSS data (Geospatial Information Authority of Japan (GSI 2004)) and the deep, low-frequency earthquakes observed by the Japan Meteorological Agency (JMA) are located just below this low-resistivity column. On the basis of these observations, Hata et al. (2016Hata et al. ( , 2018 interpreted the low-resistivity column as a magma supply system connecting the magma reservoirs. Some geophysical phenomena associated with the magma hydrothermal system have been detected at shallow depths beneath the post-caldera cones by analyzing seismic data (Yamamoto et al. 1999;Takagi et al. 2009;Ichihara et al. 2018) and by conducting audio-frequency magnetotelluric (AMT) ) and MT (Handa et al. 1998;Takakura et al. 2000;Asaue et al. 2005;Utsugi et al. 2009) surveys.
To clarify the structure and nature of the shallow part of the low-resistivity column detected by Hata et al. (2016Hata et al. ( , 2018 and its relationship to hydrothermal activity beneath the post-caldera cones, we revised the previous model of the resistivity structure, focusing on the shallower parts beneath the post-caldera cones by using new data from nine additional sites near Naka-dake crater, which is a currently active crater located in the central part of Aso volcano. Our aims were to identify the spatial relationships between the low-resistivity column, Nakadake crater, and geothermal manifestations at the surface, and to determine the origin of the low-resistivity column by additionally considering other geophysical and geochemical observations.

Geological background
Aso Volcano is located in central Kyushu, southwestern Japan (Fig. 1a). The Beppu-Shimabara Graben in central Kyushu is bounded by major faults (Matsumoto 1979). The southern boundary of the graben, which passes near Aso Volcano, divides the surface geology in this area, whereby Mesozoic and Paleozoic accretionary complexes crop out south of the boundary, and plutonic rocks, metamorphic rocks, Mesozoic sedimentary rocks, and Paleogene volcanic products are exposed mainly to the north of the boundary (GSJ 2002). Tectonic stress in the graben is represented by normal faults that indicate N-S extension, and by dextral strike-slip faults along the southern boundary of the graben (Matsumoto et al. 2015).
The volcanic activity of Aso Volcano started at ~ 300 ka, with four caldera-forming eruptions occurring between that time and 90 ka (Ono and Watanabe 1985). The last caldera-forming eruption at 90 ka, named Aso4, was the most active, erupting more than 600 km 3 of material (Machida and Arai 2003). The caldera is approximately oval shaped and measures 25 km × 18 km in area. More than 17 post-caldera cones have formed since the caldera-forming eruptions.
Eruptions occur intermittently from the Naka-dake post-caldera cones. There have been many eruptions from Naka-dake crater in historic times, including more than 70 eruptions since the beginning of the twentieth century. Naka-dake crater is characterized by cyclic activity involving: (1) crater lake formation (with a lake area of 46,000 m 2 ; Terada et al. 2012) with acidic, hot water (temperature > 50 °C, pH < 1.0; Ohsawa et al. 2010) as well as high-temperature volcanic gas emissions from a fumarole outside the lake; (2) drying up of the crater lake following the emergence of another high-temperature fumarole; (3) eruptions of mud as well as phreatic and phreato-magmatic eruptions; and (4) Strombolian eruptions (Ono et al. 1995;Ikebe et al. 2008). During this cyclic activity, the crater lake is maintained by a balance between the inflow of hydrothermal waters, seepage of lake water, precipitation, and evaporation, and the system indicates the presence of a persistent hydrothermal system just below the crater (Terada et al. 2012).
Surface manifestations of geothermal activity appear inside Aso caldera (Fig. 1). Active fumaroles (Yunotani; Fig. 1), hot springs (Jigoku and Tarutama; Fig. 1), and hydrothermal alteration zones are distributed along the western foot of the central cones. Wells drilled in this area (W1 and W2; Fig. 1) emit superheated steam, showing that a vapor-dominated system is developed in the hydrothermal field. The hot spring waters in this area have very low Cl, but relatively high SO 4 and HCO 3 concentrations, which indicates steam-heated surface waters (Yamasaki et al. 1978). In contrast to this activity, abundant spring waters that are distributed widely on the caldera floor around the post-caldera cones show high concentrations of Cl of magmatic origin (Takahashi et al. 2018). The Cl flux from these spring waters is estimated to be 30 ton/day, which is two times greater than that emitted by volcanic gases from the summit crater. The Cl and SO 4 fluxes are highest at the post-caldera cone on the northern caldera floor (Honzuka; Fig. 1).

Magnetotelluric measurements and modeling
We acquired broadband MT data for 100 sites in and around Aso caldera in 2015 and 2016 (Hata et al. 2016Matsushima et al. 2018). Five-component electromagnetic data, including the horizontal components of the electrical field and the horizontal and vertical components of the magnetic field, were processed and remotely referenced to produce the MT response functions, which included the four components of the MT impedance tensor and the two components of the vertical magnetic field transfer function (tipper vector). We obtained the previous 3D electrical resistivity model that was based on inversion of the MT response functions (Hata et al. 2016. However, the dataset of the model did not cover the summit area of Aso Volcano owing to enhanced volcanic activity at that time of data collection, and we therefore obtained additional Monthly volcanic activity report of Aso Volcano; https ://www.data.jma.go.jp/svd/vois/data/tokyo /STOCK / month ly_v-act_doc/month ly_vact_vol.php?id=503).
The 3D electrical resistivity model was then revised using these data to clarify the structure beneath the post-caldera cones at depths of < 10 km BSL. For the model revision, we used the code of Siripunvaraporn and Egbert (2009) which uses complex MT response functions for 16 total periods ranging from 0.005 to 2380 s for the respective sites, following Hata et al. (2018). The differences between our study and that of Hata et al. (2018) are as follows: 1. We excluded sites at which observations were discontinued halfway through the entire acquisition period as a result of high artificial noise. This resulted in the revised model using data from 99 sites (90 sites from previous observations in 2015 and 2016, and 9 sites used in 2017 and 2018). 2. The spatial dimensions of the model were 1400 km × 1400 km × 1001.7 km, constituting 114 cells × 108 cells × 70 cells in the north (x), east (y), and down (z) directions of the model, respectively. The vertical grid spacing was set to be uniform (50 m) above sea level, but increasing from 0.05 to 300 km with depth below sea level. The thickness of the model domain above sea level was 1.7 km, whereas that below sea level was 1000 km. The minimum horizontal cell size was set at 0.1 km × 0.1 km in the area of the central post-caldera cones and was increased to a maximum of 150 km × 150 km with increasing distance from the center of the model domain. The finer cell size was used to reproduce the steep topography in the summit area of the volcano, including the depression of Naka-dake crater. Figure 1b shows the central part of the model domain (maximum cell size, 2 km × 2 km). 3. We reconstructed the topography of the model domain by setting fixed resistivities for the air (100 MΩm) and sea (0.33 Ωm) in the corresponding grids. The positions of the boundary between land and air or land and sea were estimated using 50 m and 10 m mesh GSI digital elevation models for land areas in Japan, and gridded bathymetry data (15 arc second) from The General Bathymetric Chart of the Oceans (https ://www.gebco .net/) for sea and land areas in other countries.
The model fit to the observed data was validated using the root mean square (RMS) residual (Siripunvaraporn and Egbert 2009). The best-fit model was selected when the RMS residual was minimized within 10 iterations. A fixed prior model, corresponding to the initial model, was used in the inversion. The inversion comprised two steps. In the first step, the initial resistivity for land was uniform at 80 Ωm, and error floors defined by Siripunvaraporn and Egbert (2009) were set to 20% for the MT impedance tensor and 30% for the tipper vector to obtain an approximate model. For the initial model, we checked the response of MT impedance. The apparent resistivity and phase of MT impedance were close to the observations at shorter periods, especially the diagonal components of MT impedance, which closely matched the observations when topography was included in the model. In the second step, the best-fit model obtained in the first step was used as the initial model, and error floors were set to 5% for the MT impedance tensor and 10% for the tipper vector to obtain a well-constrained final model. The best-fit model was obtained when the RMS residual in the first step was 0.99, and in the second step when it was 2.19. The latter model was the final model used in this study.
The fit between observed and synthetic data for the final model is shown in Figs. 2 and 3. Figure 2 shows the sounding curves of the apparent resistivity and phase for the four components (Z xx , Z xy , Z yx , and Z yy ) and the real and imaginary parts of tipper vector at three representative sites (sites 401, 406, and 409) around Naka-dake crater. The sounding curves of the six other sites around Naka-dake crater are presented in Additional file 1: Figs. S1 and S2. A comparison of the observed and synthetic data shows that their off-diagonal components (Z xy and Z yx ) agree well, whereas the diagonal components (Z xx and Z yy ) show some disagreement. Sounding curves of the response calculated from the model of Hata et al. (2018) are also shown in Figs. 2, Additional file 1: Figs. S1, and S2, which confirm that those authors' model does not reproduce the observed data. Figure 3 shows phase tensor ellipses (Caldwell et al. 2004) and induction vectors (Perkinson convention) for periods of 1, 31, and 182 s for the observed data (left panel) and the synthetic data from the final model (right panel). The high absolute skew angle of the observed data indicates that it is preferable to describe the subsurface in terms of 3D structure rather than 2D. The induction vectors of the observed data converge on the post-caldera cones at periods of 1 s and 31 s, whereas they tend to be oriented southwest at a period of 182 s. The induction vectors trending southwest at longer periods correspond to the regional trend  in central Kyushu (Handa et al. 1992) and are not relevant to the particular structure of Aso Volcano. Overall, the synthetic data adequately reproduce the observed data. From these comparisons and the RMS residual, we concluded that the final model was reasonable.

Results
The 3D electrical resistivity model obtained in this study is shown in Fig. 4. Figure 4a and b shows crosssections of resistivity along the lines in Fig. 4c oriented ESE-WNW and NNE-SSW, respectively. The predominantly low-resistivity (< 1 Ωm) body extends subvertically beneath Naka-dake crater at depths of 0-10 km BSL and dips steeply to the north-northeast. A low-resistivity (< 10 Ωm) layer is identified extending from the lowresistivity column at depths shallower than 2 km BSL and extends laterally beneath the post-caldera cones. High-resistivity (> 500 Ωm) blocks are distributed under the northern caldera boundary and east of the model domain. No structure linked to caldera formation, such as a depression of the caldera floor, can be identified in the resistivity model. The differences in low-resistivity areas between our new model and the previous model of Hata et al. (2018) are shown in Fig. 5. The low-resistivity areas are more scattered at shallow depths in the previous model compared with our new model. In addition, the location of the low-resistivity column is different, in that the head of the column is positioned just beneath Naka-dake crater in the new model, whereas it is located on the west side of the crater in the previous model.

Validity of the resistivity model
We evaluated the sensitivity of the observed data to the low-resistivity column in the final model. After changing the resistivities of cells within the column, MT responses were calculated through 3D forward modeling. For test models 1 to 3, cells with resistivities of < 100 Ωm within the low-resistivity column were replaced with a constant value of 100 Ωm from the surface to a depth of 10.0 km BSL in test model 1, from the surface to a depth of 2.5 km BSL in test model 2, and through a depth range of 5.0 to 10.0 km BSL in test model 3. The test models and sounding curves of the tests are shown in Additional file 1: Figs S3, S4, and S5. RMS residuals were 3.12 for test model 1, 3.06 for test model 2, and 2.45 for test model 3. The differences between the final model (RMS residual = 2.19) and the test models are significant, meaning that the sounding curves of the test models are distinct from those of the final model. The disparity is obvious at all selected sites for test model 1 (Additional file 1: Fig. S3), whereas it is obvious at summit sites for test model 2 (sites 409 and 410 in Additional file 1: Fig. S4), and a mountainside site for test model 3 (site 022 in Additional file 1: Fig. S5). The significantly higher RMS residuals of the test models relative to the final model indicate that the resistivityvalue-replacement models are inferior compared with the final model. We conclude, therefore, that the lowresistivity column is sensitive to the observed data at all depths. These sensitivity tests indicate that the low-resistivity column (< 100 Ωm) is highly certain to be located beneath Naka-dake crater and that the column extends to a depth of at least 10 km BSL.
Because Naka-dake crater is steep and sometimes filled with acid lake water, which has extremely low resistivity (Kanda et al. 2019), we reexamined the existence of the low-resistivity column using the magnetotelluric phase tensor, which is insensitive to galvanic distortion (Caldwell et al. 2004). Figure 6 shows enlarged representations of the central area of phase tensor ellipses for periods of 1, 31, and 182 s for the observed data (left panel), the synthetic data from the final model (central panel), and the synthetic data from test model 1 (right panel). The geometric mean of the maximum and minimum phase tensor (Ф 2 ) of the observations exceeds 45° on the northern side of Naka-dake crater. The value of Ф 2 indicates the vertical conductivity gradient, with values greater than 45° indicating decreasing resistivity with increasing depth (Hill et al. 2009). Although the distribution of Ф 2 of the final model is somewhat different from observations at periods of 31 s and 182 s, the model results show the same distribution at a period of 1 s. However, the distribution of Ф 2 of test model 1 does not match the observations at any period. In conclusion, the similarity of the observed pattern of Ф 2 to that of the final model suggests that the low-resistivity column is located under the northern part of Naka-dake crater.

Comparisons with other observations 2D resistivity model
We assessed the consistency of our 3D model with other observations. Figure 7a, b shows cross-sections from the 2D model of Takakura et al. (2000) and our final model along lines A-A' and B-B' in Fig. 7c. The cross-sections of our model retain the characteristic features of the 2D model. The shallow conductive layer in the cross-sections was interpreted by Takakura et al. (2000) as a rock layer in which pores had been filled with hydrothermal fluids or as hydrothermally altered rock.

Well log data
In Fig. 8, the distribution of resistivity with depth for cells corresponding to wells W1 and W2 (Fig. 1) drilled by the New Energy and Industrial Technology Development Organization (NEDO 1994a) is superimposed on the electrical log data. The resistivity calculated from our model follows the general trend of the well log data except for the short-wavelength variation. As shown in the well log data, there is little difference between the resistivities of the products of the caldera-forming eruption (Aso4) and those of the products associated with post-caldera activity, which indicates that the distribution of Aso4 products is difficult to detect from the resistivity model alone. Where the alteration product montmorillonite is abundant in the rock, the resistivity is low in our model as well as in the well log data. Temperature increases with depth owing to the vapor-dominated geothermal reservoir in this area (NEDO 1994a), but resistivity shows little change with depth.

S-wave velocity structure
The S-wave velocity structure in and around Aso volcano was obtained by Huang et al. (2018). Figure 9 compares our resistivity structure with the velocity structure of Huang et al. (2018). The grid spacing of Huang et al. (2018) is 0.05°, which corresponds to ~ 5.6 km and ~ 4.7 km in the north-south and east-west directions, respectively. Although the depths are displaced upward by ~ 1 km, the distributions of velocity anomalies L1, L2, L4, and H2 identified by Huang et al. (2018) coincide well with the areas of low resistivity (low-velocity anomalies L1, L2, and L4) and high resistivity (high-velocity anomaly H2) in our model. In particular, anomaly L2, which was interpreted by Huang et al. (2018) as the top of the magma reservoir, corresponds to the low-resistivity column of our model. However, the high-velocity H3 and low-velocity L3 anomalies do not correlate clearly with our resistivity data.

Hypocenters
A comparison between the hypocenters of earthquakes and our resistivity model for several horizontal slices shows that seismicity is absent in areas around those parts of the low-resistivity column with resistivities of < 30 Ωm (Fig. 10). If we assume that earthquakes do not occur in the plastic region, this spatial correlation suggests that the plastic region has a resistivity of < 30 Ωm. Generally, the transition from brittle to plastic deformation occurs across a narrow zone of rock, typically at temperatures in the range of 370 to 400 °C, although the transition temperature depends on the strain rate and the properties  Fig. 1 of the rock (Fournier 1999). If the low-resistivity column were located within the plastic region, it would have an overall temperature of > 370 °C. Such a high temperature would exclude the possibility that the low resistivity arose from acidic thermal waters produced by condensation of volcanic gases at temperatures below the critical temperature, and it would also exclude the possibility that the low resistivity was produced by alteration products, such as montmorillonite, because such products are destroyed at temperatures in excess of 200 °C.

Source of crustal deformation
The source of crustal deformation of the Mogi model reported by Ohkura and Oikawa (2012) is located on the western side of Naka-dake crater at a depth of 5 km (Fig. 1), which corresponds to the P-and S-wave velocity perturbation detected by Sudo and Kong (2001). Ohkura and Oikawa (2012) and Sudo and Kong (2001) interpreted this anomaly as a magma reservoir. Tsutsui and Sudo (2004) identified a region without seismic reflectors and interpreted it as a magma conduit connecting the magma reservoir to Naka-dake crater. From those observations, an active magma plumbing system is inferred beneath the west side of Naka-dake crater.
The resistivity anomaly is absent from the region of the source of crustal deformation in our final model (Fig. 10). The ability to detect the anomaly was checked using sensitivity test model 4 (Additional file 1: Fig. S6). Assuming an andesite melt composition and ~ 10% melt fraction, resistivity values were replaced by a constant value of 5 Ωm at cells within the supposed magma reservoir, which is a cylinder of 4 km diameter located at a depth range of 4.25 to 10.0 km.  observed data, and we cannot detect the anomaly from the observed data. This insensitivity may be caused by an electrically isolated low-resistivity anomaly that does not connect to the subsurface low-resistivity layer. A similar circumstance has been identified for Usu volcano (Matsushima et al. 2001), where MT observations are unable to detect the magma reservoir located beneath the subsurface low-resistivity layer as predicted by petrological investigations.
The location of the source of crustal deformation is clearly different from that of the low-resistivity column in our final model. The discrepancy in the locations makes it unlikely that the low-resistivity column is a part of magma supply system associated with the magma reservoir inferred from geodetic and seismic studies.

Resistivity of magma
The low-resistivity column identified in the present study is subvertical, dips steeply to the north-northeast, and trends toward the sill-like source of crustal deformation of GSI (2014) and the distribution of low-frequency earthquakes that occur below 15 km BSL (JMA earthquake catalog). Hata et al. (2018) interpreted the low-resistivity column as a pathway for magma that originated in a deep magma reservoir inferred from receiver function analysis (Abe et al. 2016).
Persistent volcanic gas emissions from Naka-dake crater (Shinohara et al. 2015) and frequent eruptions of magma require a system that supplies magma and volcanic gas from a deep-seated magma reservoir to the surface. On the basis of melt inclusions estimated to have been entrapped at depths shallower than 3 km BSL, Saito et al. (2018), inferred the presence of a magma conduit beneath Naka-dake crater, considering the differing densities of the volatile-bearing magma and the host crust, with the head of the magma being at sea level. Their estimated depth of the magma head (i.e., sea level) agrees well with the position of the top of the low-resistivity column, which is located at a depth of − 600 m BSL just beneath Naka-dake crater. If we consider a resistivity value of 0.50-0.52 Ωm for silicate melt having the  (Saito et al. 2018), the low resistivity of the column (~ 1 Ωm) can be explained in terms of a mixture of melt and crystals with a high melt fraction of > 60%. The estimated melt fraction coincides with the melt fraction of the same scoria samples (61-65%) determined petrologically by Saito et al. (2018). Details of the relationship between resistivity and melt fraction have been presented by Hata et al. (2018). The correspondence between the melt fraction estimated from the resistivity and that from petrological study suggests that the low-resistivity column represents another magma pathway originating from a deep-seated magma reservoir, which is different from the magma plumbing system deduced from geodetic and seismic studies. However, if we consider that the whole of the low-resistivity column represents the magma conduit, the diameter of the column would be too large for a conduit.  (Huang et al. 2018; contours in each section) and resistivity of our final model (colors in each section) along lines a-d on the topographic map. L and H indicate low-and high-velocity anomalies marked by Huang et al. (2018), respectively. The depth is relative to BSL. The color scale for topography is as in Fig. 1  a  b c d e f 30 9.00−10.00 km 8.00−9.00 km 7.00−8.00 km 6.00−7.00 km 5.00−6.00 km 4.25−5.00 km 10 (Ωm) Fig. 10 Horizontal sections of the final resistivity model obtained from 3D inversion at a 4.25-5.00 km, b 5.00-6.00 km, c 6.00-7.00 km, d 7.00-8.00 km, e 8.00-9.00 km, and f 9.00-10.00 km below sea level. Black circles show the epicenters of earthquakes at the corresponding depths detected by JMA from January 1983 to March 2016, and red circles show the epicenters of aftershocks of the 2016 Kumamoto Earthquake from April 2016 to September 2016. Black triangles indicate the five main post-caldera cones. The red star in (a) indicates the position of the source of deformation at 5 km BSL detected by Ohkura and Oikawa (2012). In the color scale, the value of 30 Ωm is indicated by an arrow because the plastic region has a resistivity of < 30 Ωm (see text)

Radius of the magma conduit
Persistent volcanic gas emissions from Naka-dake crater indicate convection of magma in a conduit, as proposed by, Kazahaya et al. (1994) allowing volcanic gas to be supplied continuously without accumulation of magma near the surface. To satisfy the condition that the degassing magma in the conduit circulates without accumulating at the degassing depth, the convection model allows prediction of the radius of the conduit from the degassing rate of the magma (Kazahaya et al. 1994;Stevenson and Blake 1998). Estimates of conduit radius for various volcanoes have typically ranged from 3 to 60 m (Kazahaya et al. 1994;Stevenson and Blake 1998;Kazahaya et al. 2002;Kazahaya et al. 2004;Shinohara 2008;Ohwada et al. 2013). These estimates are much smaller than the radius of the low-resistivity column in our final model (about 300 m). This disparity indicates the need for another source of the low-resistivity column.

Resistivity of brine
Another possible source of low resistivity is a brine within an H 2 O-NaCl system, which has a critical temperature higher than that of pure water. For instance, at 125 MPa, which is roughly equal to a depth of 7 km assuming a lithostatic pressure gradient and a temperature of 800 °C, an H 2 O-NaCl fluid can be separated into a vapor with 4 wt % NaCl and a brine with 50 wt % NaCl (Shinohara 1994), although the mass ratio of vapor to brine exsolved from the magma depends on their unknown initial concentrations (Shinohara and Kazahaya 1995). Such brine and gas exsolved from a magma moves into the plastic region owing to the pressure difference between the magma and the surrounding fluid, which subjects it to a hydrothermal pressure gradient (Fournier 1999).
The bulk resistivity of a mixture of low-salinity vapor and high-salinity brine has a low resistivity (0.1 to 1.0 Ωm; Afanasyev et al. 2018) that is consistent with the resistivity of the column in our final model. The envelope of saline brine-vapor surrounding the columnar magma was modeled on the basis of a geochemical study of White Island, New Zealand (Giggenbach 1987). Because high-salinity brine has a low resistivity corresponding to that of the low-resistivity column, the presence of a brine envelope can explain the low-resistivity column. If it is assumed that brine exsolution is caused by deepseated solidifying magma (Shinohara 1994), then such brine would need to ascend and accumulate in shallower depths to explain the low-resistivity column.

Source of the shallow hydrothermal system below Naka-dake crater
The volcanic tremors recorded before the 2014 ash-gas emissions were located between the surface and a depth of ~ 400 m below Naka-dake crater (Ichihara et al. 2018), which corresponds to depths of − 1150 to − 750 m BSL. Thin low-resistivity layers are observed at the depth of volcanic tremors, and these layers have been interpreted to reflect infiltration of acid hydrothermal water (Minami et al. 2018;Kanda et al. 2019). Tanaka (1993) observed geomagnetic changes that show heat accumulation or consumption associated with hydrothermal conditions at a depth of about 200 m below the crater rim of Nakadake crater. These observations indicate the occurrence of hydrothermal activity related to interactions between volcanic gas and groundwater beneath Naka-dake crater. The cyclic activity in Naka-dake crater, including the appearance of crater lakes, is controlled by hydrothermal activity. Our resistivity model shows that the top of the low-resistivity column is located at a depth of -600 m BSL. This spatial relationship suggests that the region corresponding to the low-resistivity column supplies volcanic gas to the hydrothermal system beneath Naka-dake crater. On the basis of the magma degassing, we conclude that the low-resistivity column represents a brine envelope and magma that constitutes the origin of the degassing activity.

Shallow hydrothermal system within Aso volcano
The low-resistivity layer (< 10 Ωm) extending laterally at depths shallower than 2 km BSL is associated with hydrothermal manifestations at the surface. The Yunotani fumarole and Tarutama hot spring are located along the western foot of the post-caldera cones (Fig. 11a). The resistivity section along the line through the Yunotani fumarole (line A-A' in Fig. 11a) shows a thin, horizontal low-resistivity layer at a depth of < 0 km BSL (Fig. 11b). This low-resistivity layer is thought to be composed of rocks containing the alteration product montmorillonite, as suggested from well log data near the fumarole (W2; Fig. 8). The log data also show that a vapor-dominated reservoir is located beneath the low-resistivity layer, suggesting that an altered layer serves as the caprock of the geothermal reservoir. Vertical, elongate geothermal reservoirs in a depth range of 1 to 3 km BSL associated with the Yunotani fumarole were identified by Asaue et al. (2005) but these reservoirs were not detected in our resistivity model owing to insufficient spatial resolution.
Honzuka spring and Uchinomaki hot spring are found on the northern side of the post-caldera cones (Fig. 11a). The resistivity section along the line through these springs (line B-B' in Fig. 11a) shows that the low-resistivity layer extends laterally from Naka-dake crater at a depth range of 0 to 2 km BSL (Fig. 11b). The thickness of the low-resistivity layer beneath the cones increases with horizontal distance from Naka-dake crater to the Honzuka area, which is located 7 km north-northwest of the crater. Hot spring water that is drawn from well W3 (Fig. 11a) has a temperature of 60 °C at a depth of 1000 m BSL (Kohei Kazahaya, pers. comm.). The base of the well is located within the horizontal low-resistivity layer, suggesting that the layer reflects the presence of hydrothermal water. Kagiyama et al. (2012) found that spring water from Honzuka has an anomalously high electrical conductivity compared with other springs along the caldera floor, possibly a result of concentrated infiltration of volcanic fluids at Honzuka. As Honzuka is located above the termination of the horizontal low-resistivity layer, we conclude that the low-resistivity layer indicates subsurface lateral flow of acidic hydrothermal fluids that were formed by condensates of volcanic gas that originated in the region corresponding to the low-resistivity column. Similar subsurface lateral flow of acidic hydrothermal fluids from a central magma conduit has been reported from Iou-dake volcano at Satsuma-Ioujima (Matsushima 2011).
High-resistivity bodies (> 500 Ωm) are distributed below the northern and eastern boundaries of the caldera (Fig. 11b). Three wells located within the northern caldera boundary penetrate granite at depths of 150 m, 482 m, and 420 m below the surface (W4, W5, and W6, respectively; Fig. 11a; NEDO 1994b). The resistivity section along line B-B' in Fig. 11b shows that the granitic rocks observed at W4 correspond to the northern highresistivity body. Uchinomaki hot spring on the northern caldera floor (Fig. 11a) has a chemical composition that suggests interaction with granitic rocks (Kohei Kazahaya, pers. comm.). From these observations, we interpret the high-resistivity bodies (> 500 Ωm) to be granitic rocks.

Conclusions
We constructed a revised 3D resistivity model for the region in and around the caldera of Aso Volcano using existing data from 100 measurement sites and data from 9 new sites near the Naka-dake crater. Our new model clearly reveals a shallow crustal structure located beneath the post-caldera cones in the center of Aso caldera, as inferred from a low-resistivity (1-10 Ωm) layer widely distributed under the cones at depths shallower than 2 km below sea level (BSL). This low-resistivity layer correlates with hydrothermal surface manifestations such as fumaroles, hot springs, and cold springs. Comparisons between well log data and the resistivity model suggest that the low-resistivity layer comprises rocks containing high-temperature hydrothermal waters and montmorillonite that formed during hydrothermal alteration. The condensates of volcanic gas and admixing with meteoric water have caused thermal water to migrate from the summit of Aso Volcano, generating the low-resistivity layer below the post-caldera cones.
Our model also identifies a subvertical low-resistivity columnar-shaped body located under Naka-dake crater at a depth range of −0.6 to 10 km BSL. The earthquakefree region around the low-resistivity column has a resistivity of < 30 Ωm, which suggests that it represents a plastic region with temperatures above 370 °C. We suggest that salty brine, exsolved from the solidifying magma, exists in the plastic region and contributes to the low resistivity of the column. Furthermore, the head of the low-resistivity column is located just beneath Nakadake crater, where volcanic gas degassing from magma generates the subsurface hydrothermal activity at Nakadake crater and at the foot of the volcano. On the basis of these considerations, we conclude that the low-resistivity column represents a magmatic-hydrothermal system composed of a brine envelope and magma that constitute the origin of the degassing activity. However, the relationship between this magma and recent eruptions is unknown because the source of crustal deformation associated with the eruptions is not observed within the low-resistivity column. In addition, the position of the low-resistivity column does not overlap with the location of the magma reservoir inferred from previous geodetic and seismic observations. The relationship between the low-resistivity column and the magma reservoir is poorly constrained and requires further, more detailed work. Recording changes in the condition of the low-resistivity column associated with eruptions will help constrain its relationship to the magma reservoir and lead to a better understanding of volcanic activity at Aso Volcano.
Additional file 1: Fig. S1 Sounding curves of the apparent resistivity and phase for four components (Z xx , Z xy , Z yx , and Z yy ), and the real (upper panel) and imaginary (lower panel) parts of tipper vector, respectively, at three sites (402, 403, and 405) around Naka-dake crater. Black dots with error bars and open red circles indicate observed and synthetic data, respectively. Blue crosses indicate responses calculated from the model of Hata et al. 2018. Fig. S2 Sounding curves of the apparent resistivity and phase for four components (Z xx , Z xy , Z yx , and Z yy ), and the real (upper panel) and imaginary (lower panel) parts of tipper vector, respectively, at three sites (407, 408, and 410) around Naka-dake crater. Black dots with error bars and open red circles are observed and synthetic data, respectively. Blue crosses are responses calculated from the model of Hata et al. 2018. Fig. S3 Results of sensitivity test model 1. The uppermost panel shows a NNE-SSW-oriented cross-section passing through Naka-dake crater for the test model. The three lower panels show a comparison of the sounding curves for observed data (black dots), final model (red circles), and test model (blue triangles) at sites 022, 409, and 410. Fig. S4 Results of sensitivity test model 2. The uppermost panel shows a NNE-SSW-oriented cross-section passing through Naka-dake crater for the test model. The three lower panels show a comparison of the sounding curves for observed data (black dots), final model (red circles), and test model (blue triangles) at sites 022, 409, and 410. Fig. S5 Results of sensitivity test model 3. The uppermost panel shows a NNE-SSW-oriented cross-section passing through Naka-dake crater for the test model. The three lower panels show a comparison of the sounding curves for observed data (black dots), final model (red circles), and test model (blue triangles) at sites 022, 409, and 410. Fig. S6 Results of sensitivity test 4. The upper panels show a horizontal section for the test model at 4.25-5.00 km BSL and the locations of the sites (028, 032, and 035) selected for this test (left figure), and a E-W-oriented cross-section passing through Naka-dake crater for the test model (right figure). The three lower panels show a comparison of the sounding curves for observed data (black dots), final model (red circles), and test model (blue triangles) at sites 028, 032 and 035.