Robust Vector MSSA for SNR Enhancement of Seismic Robust Vector MSSA for SNR Enhancement of Seismic RecordsRecords

—A robust vector MSSA algorithm for denoising seismic data is introduced, and initial results of applying this algorithm to denoise synthetic and real seismic data records are presented. In particular, the MSSA algorithm, originally applied to denoise scalar seismic waveﬁelds, is generalized to denoise vector seismic wave ﬁelds. We also introduce a robust rank reduction algorithm within the MSSA denoising algorithm which attenuates erratic signals in the input waveﬁeld without distorting the output.


I. INTRODUCTION
S IGNAL PROCESSING of seismic records which enhances the signal-to-noise ratio (SNR) is important in various applications of geophysics, including computing velocity models of the Earth's interior, rockburst monitoring for improving mine safety, and hydraulic fracturing monitoring for identifying fracturing hazards and optimizing extraction of oil and natural gas resources from producing reservoirs [1], [2], [3]. In each of these applications, if arrays of receivers are used to record noisy seismic signals, array signal processing techniques may be utilized to achieve better SNR improvement than if signals from individual receivers are processed in isolation.
For seismic monitoring applications which require locating seismic events from noisy receiver data, there are several factors affecting the accuracy of computed locations which must be understood to optimize results [4]. These factors include random and coherent noise contaminating received signals, receiver array geometry, accuracy of the subsurface velocity model, and the choice of mathematical method for computing the event hypocenter location from received signals [5]. Here we are principally concerned with signal processing for removing random noise from signals received during seismic monitoring, prior to event location computation.
More specifically, the objective of this article is to assess the utility of singular spectrum analysis (SSA) to removing random noise from seismic records. In the past, SSA has been applied in the f-x domain to denoise signals received by 1 dimensional borehole receiver arrays, and multichannel singular spectrum analysis (MSSA) has been applied in the f-x-y domain to denoise signals received by 2 dimensional surface receiver arrays [6]. Here we extend these results by examining the performance of robust vector f-x SSA and fx-y MSSA algorithms for denoising linear and surface array signals. II. ROBUST VECTOR SSA/MSSA ALGORITHMS We introduce the robust vector SSA/MSSA algorithms by recalling the application of f-x SSA and f-x-y MSSA to linear and surface array signal processing, and then explaining how these algorithms may be modified by vector processing and robust rank reduction of Hankel matrices [6]. Assuming the receiver array is an N x by N y rectangular grid in which each receiver (i, j) records a single component displacement signal S(i, j, n) at discrete times t n = n∆t, the discrete Fourier transform coefficientŜ(i, j, f ) of the trace recorded by each receiver in a particular window of time at a particular frequency f can be computed to determine the matrix: (1) In turn, the elements of this Fourier transform coefficient matrix can be used to write an L x L y by K x K y level 2 block Hankel matrix M f : in terms of the L x by K x Hankel matrices: where: Note that in the event N y = 1 and the receiver array is a linear array, the level 2 block Hankel matrix M f equates to the level 1 Hankel matrix M 1,f Given M f , the f-x SSA algorithm, applicable when N y = 1, and the f-x-y MSSA algorithm, applicable when N y > 1, proceed by computing a rank reduced approximation of M f via singular value decomposition and performing successive antidiagonal averaging operations to compute a new (i.e. denoised) Fourier coefficient for each receiver in the array. These denoised Fourier coefficients, defined across all frequencies f , determine denoised traces at each receiver via inverse Fourier transform. In passing, we note that for surface array signals determined by a noiseless superposition of k 2D traveling waves with different wave vectors, the rank of M f is k, a fact which can be used to determine appropriate rank reduction of the block Hankel matrix in the denoising algorithm if the number of seismic event signals in a noisy seismic record is known [7].
To remove random noise from 3 component (x, y, z) receiver array displacement signals, it is also possible, rather than applying an f-x-y MSSA algorithm to each component separately, to denoise the 3 components simultaneously by writing a single block Hankel matrix 3 M f using the 3 component Fourier coefficient data of the surface array in which each matrix elementŜ(i, j, f ) in equation (3) is replaced by the 3 component column vector: This 3 component variation of single component f-x-y MSSA, referred to as vector f-x-y MSSA, is expected to offer improved performance over single component f-x-y for removing random noise from surface arrays when the component signals are correlated in time [8].
When receiver array signals are contaminated by non-Gaussian erratic noise in isolated receivers, as may occur because of power line interference or transient electronic performance of receivers, the results of f-x-y MSSA denoising may be negatively affected because the result of singular value decomposition of the block Hankel matrix is sensitive to such noise [9]. To improve stability of the results of applying the denoising algorithm in such cases, rank reducing the block Hankel matrix by robust matrix product approximation instead of singular value decomposition may be performed [11]. More specifically, each block Hankel matrix M f of dimension m×n may be approximated by a matrix product U V H of rank at most p, where U is an m×p matrix, V is an n×p matrix, and an iteratively reweighted least squares algorithm is applied to minimize the expression: for a set of weights w t ij , thereby determining updated values of the U and V matrix elements u t+1 ik and v t+1 jk from their matrix elements at iteration t, and these updated values are used to repeatedly update the weights w t+1 ij using a weight function W via the relation: for σ equal to the M.A.D. of the matrix with elements: Using this robust rank reduction method allows for f-x-y MSSA denoising of seismic displacement signals containing erratic trace noise, because erratic traces are associated with outlier matrix elements in block Hankel matrices written in the f-x-y domain whose weights in the rank reduction process are iteratively decreased. Note that for the purpose of this article, we focus on an f-x-y MSSA denoising algorithm which implements robust rank reduction of 3 component Hankel matrices 3 M f in an algorithm we've called robust vector MSSA. Different weight functions W can be selected for computing weights from the matrix elements of U and V . Following previous research, we select the Tukey weight function: for updating the weights w t ij , and choose the number of weight iterations to be 9, recalling that convergence of the weight iterations is not mathematically proven [9]. For each particular set of weights w t ij , new values of the U and V matrix elements are determined by alternating solution of weighted least squares problems for the values u t+1 ik and v t+1 jk (i.e. criss-cross regression), where the alternating solution is computed for a single iteration to avoid problems with crisscross regression returning non-optimal values of the matrices U and V [10].

A. Synthetic Data
In this section we characterize the difference between single component f-x-y MSSA denoising and robust vector f-x-y denoising of a synthetic microseismic event signal received by a 2D surface array during seismic monitoring. To do this, we begin by comparing the SNR improvement achieved with vector f-x-y MSSA denoising to the SNR improvement achieved with single component f-x-y MSSA when the input noise is Gaussian, and then compare the erratic noise attenuation of robust vector f-x-y MSSA to that of vector f-x-y MSSA when the input signal is contaminated by erratic noise in one of the receiver channels.
Here we consider a synthetic seismic event with magnitude −2.0 on the Richter scale, and a hypocenter location whose distance from the surface array is on the order of hundreds of meters, as is typical of microseismic events induced by hydraulic fracturing [13]. More specifically, we suppose the surface array is a 28 by 28 square grid of receivers with distance 1.5 between each pair of receivers in (x, y, z) locations (0, 0, 0) − (40.5, 40.5, 0), the seismic event occurs at (200, 300, −200), where distances are measured in meters, and that the source is a left lateral strike slip double couple with moment tensor: where the seismic moment M 0 = 10 6 N m is determined by the Richter scale relation: and the source time function s(t) is a Ricker wavelet centered on f 0 = 60Hz [14]: (15) Figure 1 shows plots of vector f-x-y MSSA denoising and single component f-x-y MSSA denoising of a noisy microseismic event signal (P and S waves) received by a column of receivers in the surface array. Specifically, having defined a microseismic signal contaminated by Gaussian white noise with SNR 0.8, plots A and C show an SNR improvement of 13.38 when vector f-x-y MSSA denoising is applied, and plots B and D show an SNR improvement of 9.99 when single component f-x-y MSSA denoising is applied. Note that this 3.39dB difference in SNR improvement between vector and single component processing decreases is similar to results of denoising multicomponent seismic data obtained previously, and decreases as the dimensions of the rectangular surface array are decreased.
In this example, each component of the P and S wave signal at any receiver in the surface array has a maximum magnitude in the Fourier domain at 60Hz, and the magnitude squared coherence between any 2 of the displacement signals x,y, and z components, in the Fourier domain, can be computed across any row or column of the surface array. Figure 2 shows a plot of the magnitude squared coherence between x and y components of the signals received by surface array column (0, 0, 0) − (0, 0, 40.5), in the f-x-y domain, at 60Hz. Note that since the distance between receivers in the column is 1.5m, the expected phase differences between adjacent received 60Hz P and S wave signals in the column are: e ik y,P ∆y ≈ e i0.73·2π(1.5m·60Hz/2500m/s)∆y (16) ≈ e 2πi·0. 026 (17) e ik y,S ∆y ≈ e i0.73·2π(1.5m·60Hz/1800m/s)∆y (18) ≈ e 2πi·0.037 , implying peaks in the magnitude squared coherence should occur at normalized frequencies 0.026 and 0.037. When erratic noise is present in one or more of the surface array receiver channels, application of the vector f-x-y MSSA denoising algorithm with regular singular value decomposition may yield poor results, as shown in Figure 3.

B. Real Data
Real seismic data from the IRIS community wavefield experiment in Oklahoma was also used to test robust vector MSSA denoising [15]. Figure 4 shows the 2D surface receiver array, and 60 seismic traces (z displacements) recorded by a linear subarray for a 4 second time interval, at sampling interval 4ms. The distance between stations in this subarray is approximately 24.5m. Figure 5 shows the result of applying robust vector f-xy MSSA to the seismic record obtained by receivers 6-24. More specifically, seismic records DP1 Out (x displacement), DP2 Out (y displacement), and DPZ Out (z displacement) show the result of iteratively applying robust vector f-x-y MSSA denoising to the noisy input seismic records 26 times, replacing the missing input (i.e. 0 valued) traces with their corresponding denoised traces in the noisy input after each iteration. Because these receivers constitute a linear subarray of the 2D receiver array, the algorithm is equivalent to robust vector f-x SSA denoising in this case. Improvement of SNR between input and output signals is not clearly noticeable in this example, since Gaussian random noise is not clearly present in the input signals, but erratic noise present in isolated traces in the input signals is successfully removed by denoising.

IV. CONCLUSION
A robust vector MSSA algorithm for denoising seismic data has been introduced, and initial results of applying this algorithm to synthetic and real seismic data records have been presented. In the case of synthetic data denoising, for a noisy seismic signal of SNR 0.8 received by a 28 by 28 surface array, application of vector f-x-y MSSA improved the signalto-noise ratio by 13.38dB, a better result than applying single component f-x-y MSSA, which improved the signal-to-noise ratio by 9.99dB. In addition, for a 16 by 16 subarray of the 2D receiver array, robust vector MSSA was demonstrated to remove erratic noise from a different synthetic signal contaminated with impulsive noise. In the case of real data denoising..
Results may be improved in several respects. Firstly, the comparison between vector and single component f-x-y MSSA denoising depends on severable adjustable parameters which have been set to fixed values to obtain initial results, including: dimensions of surface array, spacing of surface array, source location relative to surface array, source event magnitude and moment tensor, noisy signal SNR, and magnitude squared coherence of the 3C components. Variation of these parameters could be more thoroughly tested, with results compared against previous results to verify their validity. Similarly, the performance of the robust vector f-x-y MSSA algorithm (including algorithm processing time) could be tested for different specifications of the erratic noise, algorithm settings (e.g. convergence criteria), and surface array dimensions.

ACKNOWLEDGEMENT
Thanks to the Incorporated Research Institutions for Seismology for publicizing seismic data obtained by the community wavefield experiment in Oklahoma, and Dr. Mirko van der Baan for helpful direction in using this data.