Characteristics of secondary-ruptured faults in the Aso Caldera triggered by the 2016 Mw 7.0 Kumamoto earthquake

The 16 April 2016 Mw 7.0 Kumamoto earthquake caused prominent fault displacements and crustal deformation, not only around the main rupture faults but also around numerous secondary-ruptured faults. The physics and characteristics of such secondary faulting have not yet been studied in detail. We investigated a set of two secondary faults that appeared at the timing of the Mw 7.0 quake in the Aso Caldera by mainly using synthetic aperture radar interferometry and fault slip modeling. The two faults were found to be associated with surface displacement offsets of several centimeters or more, in the oblique sense of right-lateral and vertical motion. Fault slip inversions found that the slip was dominantly in normal sense with smaller contribution from the right-lateral component. The deeper limit of the slips was estimated to be around 1.3 km, which may coincide with the boundary between the superficial sediment layer and the basement rock. The shallowness of the slip and the difference in the dip angles of the main secondary fault and the Mw 7.0 seismogenic fault suggest separation of the two fault systems, although the fault strike and sense of motions were similar. The amount of slip on the two secondary faults was larger than that expected from the scaling law derived from seismogenic faults, which may indicate the difference in the physics of seismogenic and secondary faultings.


Introduction
The 16 April 2016 Mw 7.0 Kumamoto earthquake, the largest event of the 2016 Kumamoto earthquake sequence, caused prominent fault displacements and crustal deformation. Complex rupture patterns associated with the earthquake were revealed by InSAR (e.g., Fujiwara et al. 2016;Ozawa et al. 2016), SAR pixel offset analysis (e.g., Himematsu and Furuya 2016;Yue et al. 2017), differential Lidar (Scott et al. 2018), and field surveys (Goto et al. 2017;Kumahara and Research Group of Inter-University 2016;Shirahama et al. 2016). Most studies on fault slip estimation using geodetic and seismic data indicated that the maximum fault slip of the earthquake was around 5 m or larger (e.g., Himematsu and Furuya 2016;Kobayashi et al. 2017;Kubo et al. 2016;Tanaka et al. 2019;Yue et al. 2017). Multiple faults having different dip and rake angles were required to explain the geodetic and seismic data. Such complex faulting has been revealed by modern instruments and analysis techniques also for other crustal earthquakes of similar size (e.g., Elliott et al. 2012;Hamling et al. 2017;Oskin et al. 2012;Wei et al. 2011).
What makes the Kumamoto earthquake unique was numerous triggered secondary faultings captured by InSAR (Fujiwara et al. 2016;Ozawa et al. 2016). We define secondary-ruptured faults, or secondary faults, as those (1) that were triggered to slip in response to the main seismogenic faulting; (2) that are significantly smaller in their slip dimensions than those of the main faulting, and (3) whose continuity from the main seismogenic fault is not discernible. Fujiwara et al. (2016) found Page 2 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 approximately 230 secondary surface faultings, some along previously mapped fault traces, and many others along unmapped locations. Since the deformation around each of the secondary faults did not extend over a broad area, Ozawa et al. (2016) speculated that the fault slips must have been limited to shallow depths. In this study, we focus on a set of secondary-ruptured faults in the Miyaji district located in the northern floor of Aso caldera. Ishimura et al. (2017) reported that vertical and dextral displacements of less than 10 cm were found in the field along the faults in Miyaji. They also excavated a pit across a surface break and found a nearly vertical fault structure. As stated earlier, the secondary faults do not seem to be connected to the main source fault of the Kumamoto earthquake. If so, how do these secondary faults extend to deeper depths? What is the relation between the rupture depth and the structural properties of the rocks forming the shallow part of the crust? The purpose of this paper is to estimate the fault geometry and slip distribution of the secondary faults in the Aso caldera by inverting InSAR data observed from three different directions, in order to obtain insights on the nature of secondary faults. This paper pairs with Ishimura et al. (submitted to Earth, Planets and Space), who present and discuss the results of a trench survey conducted on the same secondary faults, along with some geomorphological discussion taking into account the detailed deformation field captured by InSAR. The set of two studies enables us to deduce insights not only on the ruptures that occurred in 2016, but also on past events including the long-term slip rates.

InSAR analysis and features of the displacement field
We investigate two faults located in the Miyaji district, Aso City, Kumamoto Prefecture, within the Aso caldera (Fig. 1). The targeted faults run in the direction of N50 • E, roughly parallel to the main rupture fault of the Mw 7.0 Kumamoto earthquake. Hereafter, we call these faults as Miyaji faults. No aftershocks were located around the Miyaji faults ( Fig. 1), which is a piece of supporting evidence that the faults are not connected to the main rupture fault of the Kumamoto earthquake.
In order to capture displacements around the faults from different directions, three sets of ALOS-2 SAR images obtained along three different line-of-sight (LOS) directions were used to form interferograms (Table 1). According to the local residents, the faulting in Miyaji occurred simultaneously with the Mw 7.0 rupture. In principle, measurements from three independent directions enable us to obtain the three-dimensional displacement field, but the north-south component has relatively large uncertainty (Wright et al. 2004;Morishita et al. 2016). This limitation is inherent to InSAR, originating from the fact that the LOS horizontal directions are close to east or west (Table 1). Standard deviation of 1 cm in each interferograms leads to standard deviation of 1.3 cm, 6.5 cm, and 0.8 cm for east-west, north-south, and vertical components, respectively, under the configuration of the data used in this study.
We conducted a 2-pass interferometry with RINC software (Ozawa et al. 2016). We multi-looked the images by 4 and 7 in the range and azimuth directions,   Page 3 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 respectively, which resulted in the spatial resolution of approximately 15-20 m in the two directions. Digital ellipsoidal height model developed by M. Tobita (Geospatial Information Authority of Japan, GSI) and T. Ozawa (National Research Institute for Earth Science and Disaster Resilience of Japan) based on the digital elevation model of GSI was used for geocoding and computation of topographic phase. We used an adaptive spectral filter (Baran et al. 2003), followed by median filter ( 5 × 5 ) and mean filter ( 3 × 3 ) to reduce noise. The resulting interferograms show high coherence in the target area (Fig. 2). While similar fringe patterns were obtained for Paths 126 and 23 (obliquely from ENE and ESE, respectively), the interferogram of Path 130 (obliquely from WSW) exhibits considerably different fringe pattern.
For phase unwrapping of the interferograms (conversion from phase to displacements), SNAPHU (Chen and Zebker 2000) was used. We subsequently converted the relative displacements mapped by the interferograms to absolute (referenced) displacements by using the GNSS displacements as a reference. For that purpose, we used the coseismic GNSS displacements in the SAR line-of-sight (LOS) direction, observed at two stations Aso (960703) and Takamori (960704) of the GNSS Earth Observation Network System (GEONET) operated by GSI, which were the nearest from our target area ( Fig. 2a-c). The offset for each interferogram was estimated in such a way that the misfit between the corrected InSAR displacements and GNSS displacements was minimized at the GNSS sites. Namely, we estimated the offset l for each interferogram by where u InSAR i and u GNSS i denote the InSAR and GNSS displacements in the LOS direction at the ith GNSS station location, respectively. For each unwrapped interferograms, the estimated offset was subtracted to produce referenced LOS displacements (Fig. 3).
We decomposed the three referenced LOS displacements along different directions into eastward, northward and upward components by (Wright et al. 2004) where u is the column vector containing the vector components of displacement (u east , u north , u up ) T , P is formed by a stack LOS vectors (unit row vectors pointing from the ground to the satellite in a local east, north, and up directions) for the Path 126, Path 130 and Path 23 acquisitions, and the vector r = (r 126 , r 130 , r 23 ) T is composed of the referenced InSAR displacements (Fig. 3). Although we used different LOS vectors for each pixel taking into account the small differences in the pointing directions to the satellite, the LOS vectors at the center of the Miyaji fault location are shown in Table 2. To better extract the displacement features related to the faultings, the obtained horizontal displacements were further converted to a horizontal component along N50 • E, approximately parallel to the strike of the Miyaji faults, by The displacement field along the N50 • E horizontal direction (Fig. 4a) shows localized displacements around the Miyaji faults, superimposed by broad regional right-lateral movement that mainly resulted from the rupture of the main fault of the Mw 7.0 Kumamoto earthquake. The vertical displacements (Fig. 4b), on the other hand, exhibits north-side-down vertical displacements across two strands of faults. More closely looking, each of two strands of faults is composed of multiple smaller faults located en echelon (see Ishimura et al. (submitted to Earth, Planets and Space) for detailed interpretation). The N50 • E horizontal component is noisier than the vertical one because of the relatively large uncertainty in the NS displacements.
The profiles of displacements across the faults show clear displacement offsets at the location of the faults (Fig. 5). The offsets approximately amount to 5 cm and 8 cm, for the northern and southern faults, respectively, in both vertical and right-lateral profiles. On the northern fault, the maximum vertical and dextral surface slip measured from the InSAR displacements were 8 cm and 19 cm, respectively, whereas those for the southern fault were 12 cm and 19 cm, respectively (Ishimura et al., submitted to Earth, Planets and Space).

Method
Our main interests on the secondary faulting can be summarized by two questions: (1) what are the depths of the slipped portion of the faults? (2) How are the two secondary faults related to each other? To answer these questions, we estimated the fault slip models of the two secondary faults using the inversion procedure developed and conducted in previous studies (Fukushima et al. 2013(Fukushima et al. , 2018Ghayournajarkar and Fukushima 2020). For more detailed explanation on the inversion procedure, see Fukushima et al. (2013).
For preparation of the inversion dataset, we first removed the regional long-wavelength displacements caused mainly by the fault rupture of the Kumamoto earthquake from each of the unwrapped interferograms (Fig. 6a). The regional displacements were modeled by a two-dimensional cubic polynomial function expressed by: where a 1 , . . . a 10 are constants estimated by the leastsquares method, and x and y are longitude and latitude. This estimation was performed in the areas shown in Fig. 6. The resulting two-dimensional polynomial function ( Fig. 6b) was subtracted from the original (Fig. 6a) to extract the displacements due to the rupture of Miyaji faults (Fig. 6c). An alternative approach for removing the long-wavelength displacements would be to model the coseismic displacements of the main rupture of the Kumamoto earthquake and subtract them, but this approach requires caution because complete removal of the long-wavelength displacements is not guaranteed. Our approach has an advantage that all the long-wavelength variations can be eliminated including long-wavelength noise due to orbital inaccuracies and atmospheric delay.
The detrended interferograms show localized displacements around the Miyaji faults and little displacements far from them except for localized depressions seemingly unrelated to the faultings, indicating the adequacy of the detrending processing. After removal of the long-wavelength displacements, the data were subsampled with concentric grids (Fig. 6d).
We used triangular dislocation elements to express curved fault surfaces, and faults were discretized with an average interval of 200 m. We solved for the geometry (5) d trend = a 1 + a 2 x + a 3 y + a 4 x 2 + a 5 xy + a 6 y 2 + a 7 x 3 + a 8 x 2 y + a 9 xy 2 + a 10 y 3 ,   Fig. 6 a Original, b estimated trend, c detrended, and d subsampled data for the three data sets. For a-c, the data are shown rewrapped for better visibility. Black curves correspond to the top of the assumed fault planes. Local variations in areas isolated from the faults in c were assumed to be unrelated to the faulting and were not used in the subsampled data sets (d) Page 7 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 and slip distribution on the northern and southern secondary faults simultaneously. We assumed that, for each of the northern or southern fault, the small en echelon structures mentioned above were connected to a single fault at depth, and approximated the top of the fault to be a smooth curve close to a straight line. The bottom side of the fault was assumed to be a straight line parallel to the top side. The lengths of the faults were taken large enough to cover the extent where deformation was observed in the interferograms (black curves in Fig. 6). The displacement fields (e.g., Fig. 4) seem to indicate smaller amount of faulting extending beyond the assumed fault location, to the northeast from the southern Miyaji fault en echelon, and perhaps to the southwest of it as well, but we did not attempt to model these displacements because they are small and the displacement gradients across the fault that are important for constraining the depth of the faulting were hardly discernible. The widths of the faults were taken large enough so that the slipped area does not reach the bottom of the assumed faults. We solved for the amount of slip on the triangular elements, dip and rake angles of the two faults, and the weight of the smoothing of the slip distribution. Slips on the fault elements have linear relationship to the observed data, whereas the other parameters have nonlinear relationship to the data. This is hence a mixed linear-nonlinear inverse problem.
The likelihood function p was defined as follows (Fukuda and Johnson 2008): where C is an arbitrary constant, d is the data vector (length N), s is a model vector composed of linear parameters (slip on the triangular elements, length M), σ 2 is the variance of the data which was set to 1.0 × 10 −4 m 2 based on the calculation in the far-field, α 2 is the weight of the smoothing matrix L(m) , m is a model vector composed of the dip and rake angles of the two faults, G(m) is the design matrix relating the slip on the fault elements to displacements on the ground surface, and −1 d is the inverse of the data covariance matrix.
Since L(m) and G(m) are functions of the geometrical model parameters, these matrices are fixed once the values of m are fixed. For the design matrix G(m) , we used solutions of angular dislocations in a homogeneous elastic half-space (Thomas 1993) with an assumption of Poisson's ratio of 0.25. We adopted the smoothing matrix defined by Maerten et al. (2005) p(s, m, σ 2 , α 2 |d) = C(σ 2 ) N /2 (α 2 ) −M/2 , where the operator concerning the ith element and the three elements in contact of it is expressed by: Here, h ij is the distance between the center of the ith and jth elements, s i and s j are the ith and jth elements in s , and l i = M j=1 h ij , where M is the number of elements, is the sum of the distance between element centers.
We assumed that the inverse of the data covariance matrix −1 d was an identity matrix, meaning that we assumed no spatial correlation of noise. Trials with data correlation taken into account led to models that had slip deeper than a few kilometers and large residuals in the vicinity of the faults. This problem occurred probably because taking the data correlation into account is equivalent to placing less weight into short-scale variations in data; that is, the inversion algorithm attempted to search for models that fit with relatively small amplitude of large wavelength noise variations rather than the short-wavelength signals around the faults.
The nonlinear parameters were solved using particle swarm optimisation (PSO) (Eberhart and Kennedy 1995). PSO is a Monte Carlo search algorithm and it evaluates the likelihood function for many sets of nonlinear parameter values. For each set, the slip distribution (values of slip on the fault elements) on each fault was solved with a damped least-squares method with non-negativity constraints. The misfit of the leastsquares inversion was used to evaluate the likelihood function.
As for the model uncertainties, neighborhood approximation and Monte Carlo integration of the likelihood function (Eq. 6) was performed to form the probability density function (Sambridge 1999), and the standard deviation of the slip distribution was calculated from the diagonal elements of the model covariance matrix. The 95% confidence intervals of the nonlinear parameters were calculated by discarding the top and bottom 2.5% of the one-dimensional marginal probability density functions (Fukushima et al. 2013).

Results
The convergence of the nonlinear parameters (Fig. 7 top) and the one-dimensional marginal probability distribution of the parameters (Fig. 7 bottom) show that wide ranges of model parameters were searched and the search converged toward one single set of model parameters, leaving a clear peak in each marginal probability distribution. The dip angle of the northern fault was less constrained compared with the southern one (Table 3), but the difference in the confidence intervals indicate Page 8 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 that the northern fault is dipping steeper than the southern one. The dip angles of the northern and southern faults were 63 • and 42 • , respectively, for the maximum likelihood model. As for the rake angles, similar values of around 240 • were estimated for the two faults, indicating dominantly normal faulting with smaller amount of right-lateral slip. The geometry and slip distributions of the faults of the maximum likelihood model and the uncertainty of slip are shown in Figs. 8 and 9. Larger dip angle for the northern fault means that the two faults are converged at depth, forming splay fault branching. It is worth pointing out that the convergence depth coincides with the deeper limit of the slipped area of the southern fault. One should note, however, that our results do not warrant the connection of the two faults at depth because the shallow slip on the northern fault does not extend down to the depth where two fault surfaces are converged. This set of fault slip models well explains the data (Fig. 10).
The estimated slip extends from the surface to ∼0.5 km with a peak slip of 11 cm, and to ∼1.3 km with a peak slip of 19 cm, for the northern and southern faults, respectively. For the northern fault, the maximum standard deviation was 3.7 cm, whereas the standard deviation for most of the southern fault elements shallower than 1.5 km was less than 5 cm (Fig. 9). The estimated slip at most area along the faults can be considered real, given the small amount of standard deviation of slip.
The moment magnitudes for the northern and southern faults amount to 3.52 and 4.36, respectively. Here, we assumed the rigidity of 6.2 GPa calculated from the density of 1900 kg/m 3 (Komazawa 1995) and the shear velocity of 1.8 km/s (Abe et al. 2010). As far as the authors are aware, the smallest moment magnitude of surface-rupturing earthquakes found so far in the world is Mw 4.7 of the 2007 Katanning earthquake in Australia (Dawson et al. 2008;King et al. 2019), highlighting the peculiarity of the secondary ruptures.
We calculated the Coulomb stress changes caused by the Kumamoto mainshock on the Miyaji faults using the   Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 Coulomb software (version 3.3.01) (Lin and Stein 2004;Toda et al. 2005) and the mainshock slip model of Ozawa et al. (2016). In this calculation, we assumed a homogeneous half-space with Poisson's ratio 0.25 and rigidity of 6.2 GPa, effective fault friction coefficient of 0.4, and the rake directions same as the estimated rake angles (227 • and 233 • for the northern and southern Miyaji faults, respectively). The Coulomb stress changes on the two faults were positive and amounted up to 32 kPa and 85 kPa for the northern and southern faults, respectively. This result indicates that the secondary faultings could have been statically triggered by the mainshock.

Relation to the subsurface structure
The superficial structure beneath the northern part of the Aso caldera floor, also called the Aso Valley, consists of sediments deposited in the caldera lake underlain by the pumice flow from the Kuju volcano (Matsumoto and Fujimoto 1969). Further down, basement rock of granite has been identified deeper than ∼500 m near the northern rim of the caldera (Matsumoto and Fujimoto 1969).
Beneath the location of the Miyaji faults, the data of gravity (Komazawa 1995;Miyakawa et al. 2016), seismic velocity structure (Sudo and Kong 2001;Abe et al. 2010), and electrical resistivity (Handa et al. 1998;Hase et al. 2005;Hata et al. 2016Hata et al. , 2018Matsushima et al. 2020)   Page 10 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 altogether indicate a low-density (porous having high water content) superficial layer of 1 to 2 km, which can be attributed to the lake deposits and pumice flow. The consistency between the thickness of the superficial sediment layer and the depth extent of the estimated slip ( ∼ 1.3 km) indicates that the fault ruptured the sediment layer. Although our results cannot exclude the possibility that the fault extends further down, it is more likely that Miyaji faults are "root-less" and cause minor ruptures within the superficial layer, considering the fact that the slip of the penultimate event was almost the same as the 2016 one (Ishimura et al., submitted to Earth, Planets and Space).
The optimal dip angle of the northern and southern faults were 63 • and 42 • , respectively. Since the dip angles found in the trench were nearly vertical (Ishimura et al., submitted to Earth, Planets and Space), it is implied that the Miyaji faults are listric. The dip angle of the main southern fault is smaller than the Futagawa fault responsible for the Mw 7.0 Kumamoto main rupture, estimated in a range of 65 • -79 • (Kubo et al. 2016;Himematsu and Furuya 2016;Ozawa et al. 2016;Yue et al. 2017). This supports the lack of continuation of Miyaji faults to the main Futagawa fault structure.
Considering the fact that numerous secondary faultings of the Kumamoto earthquake appeared not only in the Aso caldera, but also outside of it, we speculate that the volcanic products (lava and pyroclastic material), commonly constituting the superficial layer in the region, facilitate secondary faultings. Around the Unzen volcano and the Beppu-Haneyama fault zone (10-20 km from the Kuju volcano), both located in the Beppu-Shimabara graben in Kyushu, there are areas with densely clustered active faults of a few to several kilometers (Fig. 11a). Similarly, the Tanna fault connected to the Hakone volcano in central Japan is surrounded by a cluster of minor faults of similar lengths (Fig. 11b). The latter resembles the tectonic setting of the Kumamoto area, having a major fault (Futagawa fault) surrounded by minor faults, and hence it is possible that the minor faults in the Tanna fault area are activated by secondary faulting.

Effect of Earth's layered structure
The deformation field in the crust due to faulting is dependent on the properties of the surrounding medium. It has been known that ignoring the effect of layering biases the estimated fault depth to be shallower than the actual depth, when the fault is located beneath a compliant superficial layer (Cattin et al. 1999;Savage 1998;Segall 2010). In our case, the faults are inferred to be located in a compliant superficial layer overlain by stiffer bedrock. We performed a synthetic test to evaluate the degree of bias in the depth and dip angle estimated under the simplified assumption of homogeneous half-space.
As the actual fault slip, we assumed a uniform slip of 1 m on a rectangular area reaching the ground surface  Fig. 10 a Modeled data corresponding to the maximum likelihood fault slip model and b residual data, calculated by subtracting (a) from Fig. 6c, for the three data sets. The data are shown rewrapped with the same colorscale as Fig. 6 for better visibility Page 11 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 with length of 2 km, dip angle of 45 • , and rake angle of 240 • . We tested with different widths, hence with different bottom depths of the faulting. We assumed that the Earth is formed by a compliant layer of 2 km thickness having rigidity of 6.2 GPa and a half-space beneath it having rigidity of 33 GPa; the latter was given based on the shear velocity of 3.5 km/s (Abe et al. 2010) and the density of 2700 kg/m 3 (Komazawa 1995) (Fig. 12a).
The observation points were assumed to be distributed in a concentric grid, similar to the data we used for the inversion. The eastern, northern, and vertical displacements at the observation points were calculated by the EDGRN/EDCMP package (Wang et al. 2003) (Fig. 12b-d).
The three-dimensional displacements were inverted to estimate the fault width and dip angle using PSO for each case of actual fault width. The results revealed that  the effect of layering was minor for our studied case (Fig. 13). When the actual width was smaller than or equal to 1.0 km, the estimated width was nearly equal to or slightly smaller than the actual width. When the actual width was larger than 1.0 km, the estimated width was slightly larger than the actual width. The maximum degree of overestimation was 7.5% (for the case where the actual width was 1.5 km). As for the dip angle, nonsystematic deviation of less than 2.5° from the actual dip angle was observed when the actual width was smaller than 0.5 km, which probably resulted from a relatively small number of points on which the deformation was observed. Systematic underestimation of approximately 1° or less was observed when the actual width was larger than or equal to 0.5 km.
From the synthetic test, we conclude, for the purpose of this research, that the bias was limited and that the assumption of homogeneous half-space and hence the discussion on the rupture depth made in the previous section were valid. Wells and Coppersmith (1994) derived empirical relationships among magnitude, rupture length and width, rupture area, and average displacement, from 244 historical events in the world. Although the relationships are basically for seismogenic faults, here we compared them with our findings on the Miyaji faults to shed light on the differences in the seismogenic and secondary faultings. Wells and Coppersmith (1994) gave, for strike-slip and normal faults, and where MD x and SRL x (x = s for strike-slip and x = n for normal slip) correspond to the maximum surface slip (m) and surface rupture length (km), respectively. We assumed here that the relationships were applicable to the Miyaji faults even though their parameter values lie around the smaller limit of the data range used to derive the empirical relationships of Eqs. (8) and (9) (minimum of MD x and SRL x for both strike-slip and normal faults used by Wells and Coppersmith (1994) were 2 cm and 3.8 km, respectively).

Fault scaling
As the surface fault lengths, we used 1.5 km and 3.2 km for the northern and southern faults, each of which correspond to the extent where slip was estimated on the (8) log(MD s ) = −1.69 + 1.16 log(SRL s ), (9) log(MD n ) = −1.98 + 1.51 log(SRL n ), surface elements by the inversion (Fig. 8). Substituting these values to the equations, the expected amount of strike-slip would be 3.2 cm and 7.9 cm, respectively, and that of normal slip would be 1.9 cm and 6.1 cm, respectively. The amount of maximum dextral surface slip measured from the InSAR displacements on the northern and southern faults was, on the other hand, 8 cm and 19 cm, respectively, and that of normal slip was 12 cm and 19 cm, respectively (Ishimura et al., submitted Page 13 of 15 Fukushima and Ishimura Earth, Planets and Space (2020) 72:175 to Earth, Planets and Space). The large deviations of the observed values from the relationship of Eqs. (8) and (9) may indicate that secondary faults tend to slip more than the average of seismogenic faulting. Whether this tendency is a ubiquitous characteristic is an intriguing question to be investigated in future studies.

Conclusions
The Mw 7.0 Kumamoto earthquake triggered numerous secondary faultings. We obtained detailed fault-parallel and vertical displacement fields around a set of two secondary-ruptured faults in Miyaji located in the Aso caldera by conducting InSAR analysis using the ALOS-2 SAR data. The results revealed clear dextral and vertical surface displacement offsets across the two faults. By conducting fault slip inversions, we further found that the northern and southern faults ruptured down to ∼0.5 km and ∼1.3 km, with the optimal dip angles of 65 • and 42 • , respectively. The two faults seemed to form a splay structure whose connection depth lay at ∼1.3 km. The deeper limit of the rupture may correspond to the base of the superficial sediment layer. The shallowness of the ruptures as well as the difference in the dip angles of the main southern fault and the steeper dip angle of the Futagawa fault, the main seismogenic fault responsible for the Mw 7.0 quake, suggest that the secondary faults are separated from the seismogenic fault system.
Our analysis suggested that the slip on the Miyaji secondary faults was larger than expected from the scaling law proposed for seismogenic faults (Wells and Coppersmith 1994). Systematic studies of secondary faulting, being detected for the Kumamoto earthquake case (Fujiwara et al. 2016;Ozawa et al. 2016) and other earthquakes by InSAR and differential Lidar techniques, will be able to elucidate the characteristics of secondary faults and provide insights into the physics of faulting.