Impact of crustal deformation detection by the DSI (difference of split-band interferograms) method with PALSAR-2 data: a case study on the 2016 Kumamoto Earthquake

Interferometric synthetic aperture radar (InSAR) is a useful tool for detecting surface deformations at high spatial resolutions. When InSAR is applied to large surface deformations, clear fringes with complicated phase gaps often appear in the interferograms. Although the surface deformations in such areas provide valuable data for earthquake research and disaster investigation, it is difficult to convert the complicated interferometric phase to surface deformation information because of the difficulties associated with phase unwrapping. To resolve these difficulties, we created multiple SAR pairs with different frequencies using a bandpass filter and calculated the difference between the interferograms generated from these SAR pairs (referred to as the DSI analysis in this study). Generally, the obtained difference corresponds to SAR observations using long-wavelength radar. Therefore, phase wrap is less likely to occur, simplifying phase unwrapping. We applied the DSI analysis to PALSAR-2 data pairs for the 2016 Kumamoto Earthquake and successfully identified large crustal deformations with complicated phase gaps in the vicinity of the surface ruptures. Comparing these results with the crustal deformations observed with global navigation satellite system measurements, the root-mean-squares of the differences were found to be approximately 4 cm. Although this accuracy was lower than that of conventional InSAR, it was nearly equivalent to that of offset-tracking analysis. It should be noted that the spatial resolution of the DSI analysis was significantly improved compared to that of offset-tracking analysis. A disadvantage of this method is that its detection accuracy is significantly degraded in zones with low coherence owing to noise amplification. The standard deviation of the noise component was approximately 2 cm for pixels with coherence > 0.7. However, for pixels with a coherence < 0.2, the standard deviation was > 10 cm, and the noise component occasionally exceeded 1 m. Despite its disadvantages, this method is effective for detecting large crustal deformations with high spatial resolution in areas where conventional InSAR processing is inappropriate.


Introduction
Interferometric synthetic aperture radar (InSAR) is a useful tool for detecting surface deformations associated with earthquakes, volcanic activities, landslides, and so on. However, InSAR may not be able to detect deformations that are extremely large. Figure 1 shows an interferogram depicting crustal deformation around the surface rupture associated with the 2016 Kumamoto Earthquake (April 16, 2016, M JMA 7.3) derived from PALSAR-2 data. Clear fringes were obtained near the fault, but many phase gaps exceeding half a cycle were observed around the fault. Crustal deformations in such areas are important for earthquake research and disaster investigation; however, the phase information around such areas is not precisely converted to the corresponding amount of crustal deformation. This is due to the difficulty of solving modulus 2π ambiguities of the phases in the interferogram (phase unwrapping). Considering simple 1-D phase unwrapping, it is impossible to solve the ambiguities in cases where the phase gaps of neighboring pixels exceed half a cycle because of the Nyquist criterion. Many phase unwrapping algorithms for InSAR processing extended to 2-D data have been proposed. One of these is the Brunch-cut method, which solves for modulus 2π ambiguities at each pixel while searching for an adequate path that connects phases while considering residues and coherence (Goldstein et al. 1988). However, large crustal deformations often cause phase gaps exceeding half of the cycle; therefore, it is difficult to obtain unwrapped phases in areas surrounded by phase gaps. Although many other advanced phase unwrapping algorithms, such as minimum cost flow (Constantini 1998) and statistical cost flow (Chen and Zebker 2000), have been proposed, it remains difficult to unwrap phases in such areas. Himematsu and Furuya (2016) applied offset-tracking analysis to PALSAR-2 data and obtained the crustal deformation near the fault. This demonstrates that offset-tracking analysis is an effective tool, but a large correlation window is necessary to obtain high accuracy. Thus, there is a trade-off between detection accuracy and spatial resolution in this method. Madsen and Zebker (1992) proposed a method known as the split-spectrum method (SSM) to derive the unwrapped phase from the phase differences of interferograms with different frequencies. SAR observations use a chirp pulse that varies the radar frequency within the pulse width to enable a high image resolution in the range (satellite to pixel) direction. SSM creates SAR pairs with different frequencies by dividing the chirp band using a bandpass filter (Additional file 4: Fig. S1a) and generates interferograms from their sub-band pairs. Bamler and Eider (2005) and De Zan et al. (2015) discussed the theoretical performance of SSM. In the 1990s, there were limited datasets with sufficiently high coherence for SSM applications. Recently, an L-band SAR with a wide bandwidth and sufficiently high coherence has become available. SSM has recently been used as a general technique to estimate ionospheric delays (Brcic et al. 2010;Rosen et al. 2010;Gomba et al. 2016;Furuya et al. 2017;Wegmüller et al. 2018). This suggests that recent L-band SAR data may be useful for obtaining the unwrapped phase using SSM. Although this is an application of SSM, in this study, we referred to it as the difference of split-band interferograms (DSI) to distinguish it from the SSM used for estimation of the ionospheric delay.
We applied DSI analysis to PALSAR-2 data for the 2016 Kumamoto Earthquake (Fig. 1) and successfully detected large crustal deformations near the fault. In this paper, we describe the principles of crustal deformation detection using the DSI analysis with practical processing procedures and discuss its detection accuracy, as well as its advantages and disadvantages.

Principles of DSI analysis
In this section, we describe the principles of the DSI analysis for obtaining an unwrapped phase. Generally, the interferometric phase is expressed as where Δρ non-disp is the non-dispersive component of the slant-range change due to crustal deformation Δρ deform , orbit difference Δρ orb , topography Δρ topo , and atmospheric delay Δρ atm , as Eq.
(2) (e.g., Rosen et al. 2010). The second term on the right-hand side of Eq. (1) is called the dispersive component, which is mainly due to ionospheric delay, where f is the radar frequency, c is the speed of light, ΔTEC is the total electron content along the radar propagation path, K is a constant (40.28 m 3 / s 2 ), and C 1 is the constant for all phases to be treated as relative values in the image. The crustal deformation component can be extracted by subtracting the orbital, topographic, atmospheric delay, and ionospheric delay components from the initial interferogram. However, it �ρ non -disp = �ρ deform +�ρ orb + �ρ topo + �ρ atm is difficult to obtain an appropriate unwrapped phase for large crustal deformations with complex phase gaps, as previously described. For example, in InSAR applications using PALSAR-2 data (radar wavelength: approximately 24 cm), a slant-range change of approximately 6 cm causes a phase change of half a cycle, making it difficult to apply phase unwrapping. To overcome this, DSI analysis was employed to obtain the unwrapped phase. DSI analysis creates SAR pairs with different frequencies by dividing the chirp band using a bandpass filter and generating sub-band interferograms from their pairs. After subtracting the orbital, topographic, atmospheric, and ionospheric delay components from their sub-band interferograms, the phase difference of the interferograms for sub-bands with the highest and lowest frequencies (Δφ diff ) was estimated. This process can be expressed as follows: where f high and f low are the center frequencies of the highest and the lowest sub-bands, respectively; Δφ high and Δφ low are the phase differences of their respective subband interferograms; and C 2 is the constant for all phases to be treated as relative values in the image. Δφ diff corresponds to the interferometric phase obtained by the radar frequency of f high -f low . The wavelength of this interferogram is c/(f high -f low ). Phase wrapping does not occur when the slant-range change is less than a quarter of the wavelength; therefore, 2-D phase unwrapping is not necessary. When an SAR image with a bandwidth of 80 MHz is divided into two sub-bands, f high -f low is 40 MHz and the phase difference of Eq. (3) corresponds to the interferometric phase observed at a wavelength of approximately 7.5 m. The size of slant-range change required to cause phase change exceeding half a cycle is approximately 1.9 m. Even if the maximum slant-range change exceeds its threshold, application of 2-D phase unwrapping is simple due to the low fringe rate. However, in this study, we attempt to obtain an unwrap phase without 2-D phase unwrapping. When the chirp band is divided into N subbands (center frequency for each sub-band is shifted to f 1~N ), the interferometric phases for the sub-bands are φ 1~N . The bandwidth of one sub-band B sub is B/N, where B is the bandwidth of the full-band single-look complex (SLC) image. The phase difference between the interferometric phases for neighboring sub-bands Δφ i -Δφ i-1 corresponds to the interferometric phase in which the frequency is equal to the bandwidth of sub-band B sub and its wavelength is c/B sub . When the slant-range change does not exceed a quarter of its wavelength, phase wrap does not appear. As Δφ 1 ~ Δφ N are wrapped phases, Δφ N relative to Δφ 1 can be derived by 1-D phase unwrapping Page 4 of 11 Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 along f 1 ~ f N (a schematic diagram of 1-D phase unwrapping along the frequency is shown in Additional file 4: Fig. S1b) and Δφ diff in Eq. (3) was obtained from Δφ N -Δφ 1 . For this analysis, the sub-band bandwidth must be limited such that the maximum slant-range change Δρ max does not exceed a quarter of the corresponding wavelength, and the number of band divisions N must be larger than 4BΔρ max /c.

Overview of the 2016 Kumamoto Earthquake
The 2016 Kumamoto Earthquake began with a foreshock (M W 6.1) at 21:26 on April 14, 2016 (JST) in the Kumamoto Prefecture of the Kyushu district, western Japan. F-net, the broadband seismometer network of the National Research Institute for Earth Science and Disaster Resilience (NIED) (Okada et al. 2004), derived a lateral-slip mechanism with a north-south tension axis at a depth of 17 km. The mainshock (M W 7.1) occurred 2 days later on April 16 at 1:25 JST. F-net determined a lateral-slip mechanism similar to that of the foreshock at a depth of 11 km (Fig. 1a). In our previous study, we detected the crustal deformation associated with an earthquake by applying conventional InSAR analysis to PALSAR-2 data and constructed a fault model that explained the observed crustal deformation (Ozawa et al. 2016). The fault model consisted of four rectangular fault planes. Three of them were right-lateral faults along the Hinagu and Futagawa faults, and the remaining plane was a normal fault with a low dip angle parallel to the Futagawa Fault. Dense fringes with complicated phase gaps appeared in the vicinity of the Futagawa Fault ( Fig. 1d), where surface ruptures with offsets exceeding 2 m were found (e.g., Shirahama et al. 2016, Toda et al. 2016. Complicated phase gaps suggesting shallow fault ruptures were also observed in the western section of the Aso Caldera (Fig. 1c); however, crustal deformation was not quantitatively obtained because of difficulties in phase unwrapping in our previous study (Ozawa et al. 2016). In this case study, we attempted to detect this crustal deformation using DSI analysis.

Conventional InSAR analysis
We analyzed the PALSAR-2 data obtained from Path 28 (descending orbit, left looking) and Path 126 (ascending orbit, left looking) (Additional file 1: Table S1). Most of the scenes were mountainous and covered by vegetation.
The observation mode was SM1 for both the paths. The center frequency was 1.2575 GHz, and the bandwidth of the SLC image was 80 MHz. The incidence angles were 32° and 24° at the center of the scene. The observation dates were April 15, 2016, and April 29, 2016, for both paths, indicating that the interferograms constructed from these SAR pairs did not include the crustal deformations associated with the foreshock on April 14.
Here, we describe the procedure for conventional InSAR analysis. We used RINC software for the analysis, which is an InSAR analysis tool developed by NIED (Ozawa et al. 2016). The number of looks in the interferogram generation was set to 11 pixels in the range direction and 15 pixels in the azimuth (along-track) direction for Path 28, and to 8 pixels in the range direction and 12 pixels in the azimuth direction for Path 126. This ensured that the ground resolution of imagery from each path was almost equal, ranging from 25 to 30 m. Although processing with a smaller number of looks is possible, we adopted a larger number of looks to reduce the decorrelation noise. The spectral shift filter described by Gatelli et al. (1994) was used for the interferogram generation. For the simulation of the topographic phase, a numerical ellipsoidal height model created from a 10 m mesh digital elevation model published by the Geospatial Information Authority of Japan (GSI) and the EGM96 geoid model (Lemoine et al. 1997) was used. Ionospheric delay was estimated using the SSM method described by Wegmüller et al. (2018). Details of the procedure are described in the following subsections. We simulated atmospheric delay based on the results of a numerical weather model (Ozawa and Shimizu 2010;Ozawa et al. 2019). The data used in this simulation were the results of the Meso-Scale Model of the Japan Meteorological Agency. The simulated atmospheric delays and interferograms before and after correction are shown in Additional file 4: Fig. S2. After subtraction of the orbital phase, topographic phase, atmospheric delay, and ionospheric delay components from the initial interferograms, the spectral enhancement filter described by Goldstein and Werner (1998) was applied with a filter coefficient of 0.8 and a window size of 32 pixels. By applying phase unwrapping to the obtained interferograms using SNAPHU software (Chen and Zebker 2000, 2001, 2002, we derived crustal deformation maps using procedure for conventional InSAR analysis. However, unwrapping errors are expected to occur in the areas where many phase gaps appear.

Estimation of the ionospheric delay component by SSM
The estimation of ionospheric delay using SSM described by Brcic et al. (2010), Rosen et al. (2010), and Gomba et al. (2016) is necessary to unwrap the phase that includes the crustal deformation component. It is difficult to estimate the ionospheric delay in areas where unwrapping errors are likely to occur. Therefore, we used the method described by Wegmüller et al. (2018) to estimate the ionospheric delay component. This method creates multiple SLC image pairs using SSM. In this study, we Page 5 of 11 Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 divided the data into four sub-bands for use in the subsequent DSI analysis. Interferograms for the highest and the lowest sub-bands were generated by applying InSAR analysis to each SAR pair. Orbital and topographic contributions were removed from the initial interferograms.
To reduce noise, we applied a spectral enhancement filter (Goldstein and Werner 1998) to the interferograms and completed further multilook processing (called the "2nd multilook" in this study). The window size adopted for the 2nd multilook was 5 × 5 pixels. Additional file 4: Figure S3a, b shows the resulting interferograms. The dispersive component (Δφ disp ) was estimated as follows: where f 0 is the center frequency of the full-band SLC image, Δφ obs denote the phase of the full-band interferogram, and Δφ diff is the phase difference of the interferograms for the highest and lowest sub-bands (Additional file 4: Fig. S3c). Outliers of Δφ diff were removed using the method described by Gomba et al. (2016) (Additional file 4: Fig. S3d) and 2Δφ disp was calculated by substituting them into Eq. (4) (Additional file 4: Fig. S3f ). As Δφ disp is a wrapped phase, it was necessary to apply 2-D phase unwrapping. However, its spatial distribution was smooth except under extreme conditions, and an appropriate unwrapped phase was easily obtained (Additional file 4: Fig. S3g). A smoothed ionospheric delay map was obtained by applying a Gaussian filter to its unwrapped phase (Additional file 4: Fig. S3h). Finally, the image was resampled to the size of the original interferogram. It should be noted that this is double the value of the ionospheric delay. Additional file 4: Figure S3i indicates that the interferogram with ionospheric component subtracted.

DSI analysis and results
For the DSI analysis, we used four sub-band SLC image pairs created in the estimation of the ionospheric delay and derived sub-band interferograms by applying the same procedure for the conventional InSAR analysis.
The frequency of f i -f i-1 was 20 MHz, corresponding to a wavelength of 15 m. A phase wrap did not appear when the slant-range change was less than a quarter of the wavelength (approximately 3.75 m). From the differential interferograms produced for all sub-bands, we obtained the phase difference of the DSI, Δφ diff , in Eq. (4) by 1-D phase unwrapping along f 1 ~ f 4 , as described in the previous section. The procedure above is shown in Additional file 4: Fig. S4.
The slant-range change distributions obtained from DSI analysis are shown in Fig. 2a, d. The observation times for Paths 28 and 126 differed by approximately half a day, but it can be assumed that the crustal deformation that occurred during the time difference in the observations was negligible because the maximum magnitude of the aftershocks that occurred during this period was M JMA 4.2. We then derived the quasieast-west (QEW) and quasi-up-down (QUD) components from their results using the method proposed by Fujiwara et al. (2000) (Fig. 3). The QEW direction was rotated approximately 3° clockwise from the east and the QUD direction inclined approximately 8° northward from the vertical direction. In general, the results showed an eastward displacement in the northern area of the Futagawa Fault and a westward displacement in the southern area, which is consistent with the focal mechanism of the main shock. The boundary separating the area of the eastward and westward displacements clearly coincided with the Futagawa Fault. The contrast between the positive and negative QEW displacements became obscure around the western rim of the Aso Caldera. In the western area, the boundary extends to the Hinagu Fault. Conversely, a clear boundary of the QUD was observed in the western area of the Futagawa Fault, whereas in the eastern area, the boundary was consistent with the Idenokuchi Fault, which is parallel to the Futagawa Fault (Fig. 3g, h). Thus, the QEW and QUD boundaries were different in this area, corresponding to the partitioning of a strike-slip and normal-slip, as described by Toda et al. (2016) and Himematsu and Furuya (2016). In addition, this study clarified that subsidence gradually increased from the west and exceeded 2 m in the area between the Futagawa and Idenokuchi faults. A maximum westward displacement of 0.3-0.4 m was also detected in this area.
Crustal deformation in west of the Aso Caldera, where a complicated fringe pattern appeared (Fig. 1c), was also clearly observed in the DSI analysis (Fig. 3d, e). As pointed out by Fujiwara et al. (2016), the surface displacement obtained by DSI analysis showed that northdown dislocations dominated in the southern area and south-down dislocations dominated in the northern area, similar to a graben (Fig. 3f ). In addition, this result showed that the deformation gap was large in the southern area, as reported by Fujiwara et al. (2016). Clear deformation gaps can be observed at points where surface ruptures have been found in field surveys (Awata et al. 2019; National Institute of Advanced Industrial Science and Technology 2019), and their extensions could also be identified in these results (thin lines in Page 6 of 11 Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 Fig. 3). Although the locations of the gaps could also be seen from the results of the conventional InSAR analysis, the DSI analysis was more resistant to phase unwrapping errors, allowing a precise amount of surface displacement to be determined. Artificial patterns in the differences between the DSI and conventional InSAR results appeared in such locations, indicating unwrapping errors in the conventional InSAR analysis (Additional file 4: Fig. S5).

Discussion
As described in the previous section, the application of DSI analysis to the 2016 Kumamoto Earthquake successfully detected detailed crustal deformation, which is consistent with previous researches. Although new seismological knowledge was not obtained in this case study, it is important to demonstrate accurate results using DSI analysis. In this section, we focus on the detection accuracy and limitations of the DSI analysis.  Page 7 of 11 Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 First, crustal deformations obtained from the DSI analysis were compared with those obtained using offset-tracking method with image amplitude. We used the GAMMA SAR processor for offset-tracking analysis. In the co-registration of SAR images for the offset-tracking analysis, we corrected for the parallax effect due to the satellite orbit difference and estimated the affine transformation coefficient by masking areas with large deformations. In the offset-tracking analysis, we set a correlation window of 64 × 64 pixels for every 24 pixels and calculated the offset with the maximum correlation. As crustal deformation can be calculated with higher spatial resolution using a smaller correlation window, we also attempted to perform an analysis using a correlation window of 32 × 32 pixels. However, the spike noise significantly affected these results, and the results of the 64 × 64 correlation window were adopted for this comparison. Figure 2b, e shows the slant-range change maps obtained from the offset-tracking analysis, and Fig. 2c, f shows the differences between these maps and the results of the DSI analysis. Generally, the slant-range change maps obtained by the offset-tracking analysis produced   (g, h) Distribution of the QEW and QUD components, respectively, enlarged around the Futagawa Fault, corresponding to Fig. 1d Page 8 of 11 Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 results similar to those obtained from the DSI analysis. Although a slight bias and ramp of phase can be observed in Fig. 2f, we suspect that this was due to errors in the SLC image matching for the offset-tracking analysis. The spatial distributions of the complicated deformations west of the Aso Caldera, derived from the DSI analysis, are sharper than those obtained from the offset-tracking analysis (Fig. 2g, h). This can be observed more clearly by comparing the profiles of the slant-range changes (Fig. 2i), indicating that the DSI analysis detected crustal deformations at a higher spatial resolution than that in the offset-tracking analysis. While DSI analysis can obtain a high spatial resolution with a pixel size, the spatial resolution in offset-tracking analysis depends on the size of the correlation window. Next, we compared the crustal deformations obtained by DSI analysis with those from GEONET, which is the global navigation satellite system (GNSS) observation network operated by the GSI (Sagiya et al. 2000). Figure 4a shows the locations of the GNSS sites used in this comparison. Crustal deformations from the GNSS were derived from the difference in the average coordinates of the F3 solutions (Nakagawa et al. 2009) for April 8 and 29, 2016. We adopted a window size of 11 days to reduce noise. As crustal deformation from GNSS included the effects of the foreshock on April 14, we subtracted the crustal deformation calculated from the two-fault model with uniform slip estimated by Kobayashi (2017). Crustal deformation from the DSI analysis was picked up from a single pixel, which included the location of the GEONET station. Figure 4b shows a comparison between the crustal deformations obtained from DSI analysis and GEONET (these values are shown in Tables S2 and  S3). The reference site for their deformation was the GEONET site "960699" (open diamond in Fig. 4a). Except for GEONET site "960701, " the difference in slant-range changes calculated by DSI and GNSS does not exceed 10 cm, and the root-mean-squares (RMS) of the differences between them were 4.6 and 3.8 cm for Paths 28 and 126, respectively. Using the same process, we compared the crustal deformations derived from the offset-tracking analysis and those from GEONET and obtained RMS values of 4.4 and 4.2 cm for Paths 28 and 126, respectively (Additional file 2: Table S2 and Additional file 3: S3). These values were very similar to the results of the DSI analysis.
We also investigated the results of the conventional InSAR analysis using the same method. The phaseunwrapping error did not appear around the GEONET stations (Additional file 4: Fig. S5). The RMS values of the conventional InSAR results were 2.5 and 2.2 cm for Paths 28 and 126, respectively (Additional file 2: Table S2 and Additional file 3: S3), indicating that the precision of the DSI analysis was lower than that of the conventional InSAR. The SM1 mode of PALSAR-2 uses a radar wave with a frequency of 1.2575 GHz, whereas the DSI analysis in this study corresponds to a radar wave with a frequency of 60 MHz, indicating that it is much less sensitive. Furthermore, phase noise due to decorrelation at almost the same level as that of conventional InSAR  Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 appeared on each sub-band interferogram, and it seems that this was not compensated for by the phase difference between them. The decorrelation noise in DSI results is amplified by f 0 /(f high -f low ) times the amount of decorrelation noise present in conventional InSAR data. The decorrelation noise in the DSI results of this case study was amplified to 21 times that in the conventional InSAR. Therefore, the accuracy of the DSI analysis is lower than that of conventional InSAR. In particular, for GEONET site "960701, " the slant-range change in the DSI analysis for Path 126 differed by more than 1 m from that in the conventional InSAR (Additional file 3: Table S3). This is much larger than the usual accuracy of GNSS observations, meaning that the difference of 1 m was due to an error in the DSI analysis. The coherence for this pixel was 0.04, which is extremely low. In addition, large amounts of random noise in the DSI analysis results were obtained near the vicinity of the fault, and coherence in this area was very low. The phase difference between sub-band interferograms for such pixels was not appropriately estimated owing to noise amplification, and the errors may have been enlarged by a misestimation of the phase difference Δφ diff in Eq. (4).
To investigate the relationship between the coherence and accuracy of the DSI analysis, we analyzed the postseismic pair of Path 126 (April 29, 2016, andOctober 14, 2016). For the conventional InSAR analysis of this pair, phase unwrapping could be applied precisely because of the low fringe rate. The observation mode and procedures of the DSI and conventional InSAR analyses were the same as those used for the co-seismic pair. Figure 5a, b shows the slant-range change maps obtained using DSI and conventional InSAR analyses. Figure 5c, d shows the interferogram and coherence values obtained from conventional InSAR analysis. In the slant-range change map produced by DSI analysis, the area where high-frequency noise is dominant corresponded to the area with low coherence. The differences between the slant-range change obtained from the DSI and conventional InSAR analyses increased with decreasing coherence (Fig. 5e, f ). The standard deviation of the pixels increased as coherence decreased. Pixels with a coherence larger than 0.7 had a standard deviation of 2 cm or less, whereas pixels with a coherence lower than 0.6 had a standard deviation of 3.5 cm. This corresponds to the RMS of the comparison with GEONET, suggesting a moderate level of precision   . 5 a, b Slant-range change maps of the post-seismic pair of Path 126 obtained using the difference of split-band interferograms (DSI) analysis and conventional Interferometric Synthetic Aperture Radar (InSAR), respectively. c, d Differential interferogram and coherence obtained from conventional InSAR. e Scatter plot showing the differences in slant-range changes obtained using DSI (y-axis) versus coherences (x-axis). f Relationship between the coherence and standard deviation of the differences in the slant-range changes Page 10 of 11 Ozawa and Himematsu Earth, Planets and Space (2022) 74:72 for these coherence measurements. An RMS of approximately 10 cm was calculated for a pixel coherence of 0.2, but the difference occasionally exceeded 1 m (Fig. 5e). The insufficient reduction of atmospheric and ionospheric delays also affected the DSI analysis results. However, it should be emphasized that the atmospheric delay was not amplified because it was included as a nondispersive component. Noise levels due to ionospheric delays were similar to those in the conventional InSAR analysis. By converting the interferometric phase into slant-range change, the ionospheric delay component was K/f high /f low ΔTEC for DSI analysis and K/f 0 2 ΔTEC for conventional InSAR from Eq. (1).

Summary
In this study, we described the use of DSI analysis to detect large crustal deformation. As a case study, we applied DSI analysis to PALSAR-2 data for the 2016 Kumamoto Earthquake and detected crustal deformation in areas where the unwrapped phase could not be determined by conventional InSAR analysis because of complicated phase gaps. Compared with conventional InSAR, the detection accuracy of DSI analysis is lower; however, DSI can detect crustal deformation in areas with large crustal deformation. Furthermore, DSI analysis is useful for detecting crustal deformation with comparable accuracy and a higher spatial resolution than those of the offset-tracking method. The detection accuracy of DSI analysis depends on the coherence of the analyzed pixels. An accuracy of approximately 2 cm was estimated for pixels with a coherence of more than 0.7-0.8, but an estimation error exceeding 1 m was occasionally obtained for pixels with a coherence of less than 0.2. The disadvantage of DSI analysis is the dependence of the detection accuracy on coherence, as low coherence produces low accuracy. However, this disadvantage will soon become insignificant, as several new L-band SAR missions, PALSAR-3/ALOS-4 and NISAR, will commence in the near future, facilitating significant improvements in the observation frequency over that of PALSAR-2/ ALOS-2. As these missions will reliably provide InSAR results with high coherence, it is expected that DSI analysis will become an effective tool for detecting large crustal deformations.