Contrary to the success of 3D imaging, Rosas‐Carbajal et al. (2017) point out a non‐negligible offset between the density estimated from muography and that from gravity. In the case of La Soufrière de Guadeloupe lava dome, such offset amounts to 0.47 g/cm3. Specifically, muography tends to yield lower density compared to the apparent density estimated from gravity data. In other words, muography detectors placed on the mountain slope perceive flux higher than expected. As an immediate solution, Rosas‐Carbajal et al. (2017) propose to separate the density sensed by gravity data ( ρ g ) and the density sensed by muography data ( ρ μ ) and relate the two vectors with an unknown constant bias (Δρ), written as ρ g = ρ μ + (Δρ, Δρ, ⋯, Δρ)T. This constant density bias is solved together with ρ μ . Lelièvre et al. (2019) derived a comprehensive formulation for estimating this sort of constant density bias.
3.6 REGULARIZATION PARAMETERS
Smoothing constraints introduced as in equation 3.13 and 3.14 successfully reduce the number of unknown parameters. However, one practical difficulty is to give proper correlation length (λ) and variance of prior density (σ ρ ), because these parameters should be given “before” the inversion is performed. Barnoud et al. (2019) use the cross‐validation method to avoid overfitting in parameter estimation. The idea of the cross‐validation is to randomly pick one data out of n muon + n grav data (leave‐one‐out cross‐validation) and apply the inversion analysis to the remaining dataset for various sets of (σ ρ ,λ); the optimal set of (σ ρ ,λ) that reproduces the originally extracted data is chosen in the end. Instead of picking one data, it is possible to divide the dataset into small k subsets and extract the one subset out of k (k‐fold cross‐validation). They examined this idea for an actual observation setup at Puy de Dôme (France), where muography is performed from three directions and 650 points of gravity data are available (Ambrosino et al., 2015; Portal et al., 2016). From their numerical test, the above two cross‐validation methods are proved to be feasible when the density profile in the mountain is randomly distributed in a Gaussian way.
3.7 DISCUSSION
The addition of surface gravity data helps to constrain the solution when the number of muography stations is small. In the case of Fig. 3.5, where only one muography station is given, the gravity data contributes greatly to locating the high‐density lava block. However, the impact of the gravity data would be not so impressive when the target is covered by several muography detectors. According to the numerical simulation by Jourde et al. (2015), the addition of gravity data is only effective when the target is surveyed by one or two muography detectors; gravity data are almost useless when the target is sampled by more than two muography sites. However, significant improvement of the resolution is obtained in deeper regions which are not covered by muography paths, because a fraction of information brought by muography data is transferred to the deep regions coupled to gravity data (Jourde et al., 2015). As technology for particle detection improves, it will be easier to install multiple muography detectors surrounding the target volcano. If the density distribution of the mountain above the muography detectors is accurately determined from only muon tomography, the gravitational contribution can be straightforwardly calculated and removed for density distribution analysis of deeper parts of volcanoes. How does the volcanic conduit extend to the deep reservoir? Is it filled with solidified magma or vacant, thus complying with deformation? These questions could be addressed by the joint analysis. One such measurement is ongoing at Mount Omuro scoria cone, Izu Peninsula, Japan, where more than ten muography detectors made of emulsion films are installed around the volcano (Nagahara & Miyamoto, 2018).
The major technical problem of joint inversion so far could be the existence of density bias, reported, for instance, by Rosas‐Carbajal et al. (2017). In most cases, muography tends to yield lower density compared to the density estimated from gravity data. This indicates that muography detectors perceive higher flux than expectation. This sort of flux excess has been observed in many research groups (Ambrosino et al., 2015; Carbone et al., 2014) and it originates from the following three possible causes: combinatorial background (Saracino et al., 2017); low‐energy muon, electromagnetic, and hadronic components of cosmic rays (Gómez et al., 2017; Nishiyama et al., 2016); and upward‐going particles entering from the rear side of the detector (Jourde et al., 2013). Muography detectors must be designed so that they are capable of rejecting these noise signals before the data is combined with gravity data. Rosas‐Carbajal et al. (2017), Barnoud et al. (2019), and Lelièvre et al. (2019) introduced a constant density bias between density values felt by muography and gravimetry. However, as the authors of these works acknowledge, the density bias may not be constant and it may vary depending on the direction of lines of sight and the topography around the detector. Figs. 3.6a and 3.6b represent the maps of average density inside Showa‐Shinzan lava dome along muons’ trajectories from two types of detectors with different energy thresholds; the top panel rejects the contribution of low‐energy (< 1 GeV) particles from analysis, whereas the bottom panel includes all the low‐energy particles and the dataset can be regarded as contaminated by density bias. The density offset between the two datasets is not constant and correlates with the true density. The large offset near the horizon may originate from the upward‐going particles. The author believes a physical treatment on the background noise reduction would be essentially important for solving the discrepancy between the two methods.
Another serious limit is that the inversion of gravity data is inherently non‐unique. In other words, an infinite number of possible density profiles can produce an identical external gravity field. Especially, locating the depth of mass anomaly is inherently difficult. This affects the accuracy when extracting the gravitation contribution from the volume of interest as described in equation 3.5. Nishiyama et al. (2014) assumes the regional trend Δg trend is constant over the entire survey region; however, this is not true in general. One possible strategy is to remove a long‐wavelength regional trend obtained by fitting a planar surface to the regional Bouguer anomaly map (e.g., Roy et al., 2018). Such detrending would cause additional uncertainty that would probably spatially correlate, which should be properly introduced in non‐diagonal elements of the data covariance matrix ( C d ).
3.8 SUMMARY
Joint inversion of gravity and muography data takes advantage of different sensitivities for the reconstruction of the three‐dimensional density distribution of volcanoes. It has been successfully applied in volcanological studies, such as the shape of the solidified lava dome, the extension of fracture zone and its hydrothermal activity, the identification and location of volcanic deposits. However, some technical difficulties are recognized, i.e., optimization of a priori parameters; density bias between muography and gravity methods; and detrending of gravity signals. Overcoming these problems would lead to improve the resolution and reliability of three‐dimensional imaging and provide a better understanding of volcanic plumbing systems.