Reconstructing Sound Speed Profile From Remote Sensing Data: Nonlinear Inversion Based on Self-Organizing Map

By establishing a linear regression relationship between the projection coefficient of the empirical orthogonal function (EOF) of the sound speed profile (SSP) and remote sensing parameters of the sea surface, the single empirical orthogonal function regression (sEOF-r) method was used to reconstruct the underwater SSP from satellite remote sensing data. However, because the ocean is a complex dynamical system, the parameters of the surface and the subsurface did not conform to the linear regression model in strict sense. This paper proposes a self-organizing map (SOM)-based nonlinear inversion method that used satellite observations to obtain anomalies in data on the sea surface temperature and height, and combined them with the EOF coefficient from an Argo buoy to train and generate a map. The SSP was then reconstructed by obtaining the best matching neuron. The results of SSP reconstruction in the northern part of the South China Sea showed that the relationship between the parameters of the sea surface and the subsurface could be adequately expressed by the nonlinear neuronal topology. The SOM algorithm generated a smaller inversion error than linear inversion and had better robustness. It improved the average accuracy of reconstruction by 0.88 m/s and reduced the mean-squared reconstruction error to less than 1.19 m/s. It thus offered significant promise for acoustic applications.


I. INTRODUCTION
As the basic physical parameter used to describe the acoustic characteristics of water columns, the sound speed profile (SSP) is crucial to applications of marine acoustics, such as underwater target recognition [1], marine environmental monitoring [2], and underwater communication [3]. The most direct way to obtain the SSP is by directly measuring the sound speed on site, which is time consuming and laborious. Owing to the significant temporal and spatial variations in underwater sound speed, it is nearly impossible to obtain a large-scale real-time SSP through direct measurement. Marine satellite remote sensing can satisfy the demands of a large scale and immediacy, but the parameters measured by The associate editor coordinating the review of this manuscript and approving it for publication was Haiyong Zheng .
using it are limited to the surface area of water. With recent developments in observation technologies, such as Argo and remote sensing, data samples of the underwater profile and ocean surface are being continually accumulated, and a strong relevance has been noted between the underwater profiles and remote sensing data. This has been used to obtain empirical relations between the surface and subsurface parameters at sea, and for the real-time remote sensing-based monitoring of underwater profiles [4]- [11].
As the sound speed is a function of temperature, salinity, and pressure (depth), the change in it reflects the temporal and spatial variations in the characteristics of sea water. Acoustic tomography is an important means of ocean exploration, and the SSP is an important parameter to this end. Munk and Wunsch proposed a method for the large-scale inversion of the SSP by using acoustic rays time differences and used this to observe global water temperatures perturbation [12]. To limit the dimensionality of the parameters in SSP-related Inversion, the empirical orthogonal function (EOF) has been introduced to provide a compact representation of the SSP. LeBlanc proved that the EOF as shape function yields the smallest error in describing the SSP [13], and its first few orders can be used to reconstruct it [14]. To reconstruct the underwater profile from surface observation, empirical relationships between surface quantities and subsurface perturbation shape were studied. Stammer and Wunsch found that variations in the height of the sea surface and its surface velocity are generally correlated with barotropic, and the surface kinetic energy of water is controlled by the first baroclinic model [15], [16]. Carnes obtained the statistical regression relationship between the EOF of the vertical temperature structure and the dynamic surface height [17], and then used the EOF with a linear relationship to invert the parameters of the water columns, this is called single empirical orthogonal function regression (sEOF-r) [4]. The US Navy has adopted this method as part of modular data prediction for the marine environment [5]. In the framework of the sEOF-r, MODAS provides a dynamic climatology that can obtain data on the height and temperature of the sea surface to predict underwater structures by constructing synthetic profiles generated by regression analysis [6]. Based on the assumption that the empirical relationships are static in time, the vertical structure of different parameters can be reconstructed [7]- [10]. The SSP were also reconstructed directly in a global scale via the sEOF-r method [11]. Although the above methods have been proven effective, they are all based on a linear framework, whereas the ocean system is a complex and nonlinear system. This makes it difficult to describe it accurately by using a linear system.
With advances in computation power, machine learning has been gradually proposed for similar problems. Hjelmervik used k-means clustering and gradient search to estimated real-time profiles from the surface parameters [18], [19]. Analyses based on the self-organizing map (SOM) have been used to extract patterns of the ocean current featuring dynamic and unique spatial and temporal structures from data obtained by HF radars [20]. Charantonis applied the SOM to data on underwater autonomous gliders to reconstruct thermohaline profiles [21]. In recent studies, a generalized regression neural network that uses the fruit fly optimization algorithm has shown a better capability to reconstruct the salinity profile than the linear method [22]. Su proposed extreme gradient boosting (XGBoost), a new ensemble learning algorithm, to retrieve anomalies in the temperature and salinity at altitudes above 2000 m in the ocean [23]. Compared with traditional inversion algorithms, machine learning can be used to obtain optimal solutions without presetting the linear model to improve accuracy. This paper proposes a SOM-based inversion scheme and evaluates it through comparative tests on SSP reconstruction in the northern part of the South China Sea. The South China Sea is a semi-enclosed marginal sea with a typical East Asian monsoon climate. Owing to the alternating rotation of southwesterly winds in summer and northeasterly winds in winter, the marine and atmospheric environments of the South China Sea have clear seasonal differences. Due to the complexity of the ocean-atmosphere system and the unique topographic effect of the South China Sea, dynamic marine activities in this region are complex [24], which makes it difficult to reconstruct the SSP. There is also a lack of data samples owing to political tensions and suboptimal experimental conditions. Because of the complex perturbation mechanism of the SSP and the lack of a large number of samples, a robust reconstruction of the SSP from remote sensing data in this region could be of immense benefit to the technical progress in SSP acquirement. By introducing the artificial neural network, the linear constraint between the parameters of the sea surface and those of the subsurface can be eliminated, and the relationship between them can be established through neuronal topology. Compared with the sEOR-r method, the results of SSP reconstruction thus obtained have a smaller error and better robustness, which improves the accuracy and applicability of the method of SSP inversion using parameters of remote sensing.

II. METHODOLOGY
The SSP c(z, t) at a specific depth z and time t can be represented by a 1 × n column vector. Each element of the matrix represents a sampling point at depth z. To provide constraints on the number of parameters in applications, the c(z, t) is usually expressed as: where c 0 (z) is the constant part in the SSP representing the background profile of the sea that is stable in the long term.
∞ s=1 a s (t)K s (z) represents the abnormal value of the SSP that results from ocean energy and matter exchange. It has timevarying characteristics but does not affect the background profile. K (z) is the shape function of the perturbation that is determined by the characteristics of the sea, and a s is the projection coefficient representing changes in the weight of the shape function with time. In the SSP inversion, the result can be reconstructed by the shape functions and their corresponding amplitudes.
The EOF is the most commonly used perturbation shape function to solve the problem of SSP inversion. It can identify the main modes of ocean perturbation while reducing noise [25]. EOFs are obtained by extracting the principal components of the SSP sample matrix. The SSP anomaly matrix X is calculated by subtracting the background mean SSP from the m × n SSP sample matrix, where m is the depth of the sample profile and n is the number of sample profiles. The EOF vectors can be derived by principal component analysis [26], [27]: where R is the covariance matrix, K is the eigenvector, and E is the diagonal matrix of the eigenvalues. The non-zero element of E represents the variance in perturbation described by the corresponding column vector of K , which is the EOF vector. The EOF vector corresponds to the largest eigenvalue that explains the greatest variance R. In general application, the first five leading modes with the highest variance are used to reconstruct the SSP with reasonable accuracy.
In applying remote sensing parameters to reconstruct the SSP, the sea area is divided into grids (2 • × 3 • grids are considered in this paper) and SSP samples are calculated by introducing temperature, salinity profiles of Argo to the empirical formula for sound speed. The perturbation of the SSP in a grid is considered to follow a law, that is, it is composed of a unique EOF mode of the sea area. The key to solving the inversion problem is to establish the relationship between the remote sensing parameters and the amplitude of K through historical samples. Based on this, the real-time remote sensing values can be used to obtain a s and then reconstruct the SSP through (1).

A. INVERSION BASED ON sEOF-R
Based on the parameters of remote sensing and SSP samples obtained the same time and in the same grid, a linear relation between the parameters of the surface and subsurface of the sea can be calculated. The regression relationship between the projection coefficients of each order of the EOF and the parameters of remote sensing can be described by [11] where H is the sea surface height anomaly (SSHA), T is the sea surface temperature anomaly (SSTA), a s is the projection coefficient of EOF when s is greater than zero, a 0 is a constant term coefficient, and A is an unknown coefficient to be determined. By using the regression expression (4), the projection coefficient can be inverted using the H and T as input parameters, and the profile can be reconstructed with the help of (1). It is easy to see that the sEOF-r is based on a linear regression between the remote sensing parameters of the sea surface and the projection coefficient. This linear relationship is based on the statistical results of a large number of samples in a specific sea area [15]- [17]. The differences between individual and statistical characteristics may lead to errors.

B. SOM-BASED INVERSION
Given that individuals do not fully conform to the linear relationship, and a large number of samples are not even strictly linear, the accuracy of inversion can be improved by violating the constraints of linear inversion. To improve inversion accuracy, we introduce the self-organizing neural network, which is a side inhibition phenomenon in a simulated biological nervous system. That is, nerve cells that win in the nerve cell competition are excited and the failed nerve cells are inhibited. The SOM is a vector projection with two layers, an input layer and a output layer (competition layer), as shown Fig. 1. It defines the nonlinear projection from the input layer to the low-dimensional competition layer. In the training process, the prototype vector moves to follow the probability density of the input data and keeps the topological structure of the data unchanged. The processing flow of SOM-based inversion is shown in Fig. 2. In the input layer, EOF coefficients, anomalies in the sea surface height and temperature, and latitude and longitude are input to train the SOM, which maps the training data to a discreet set of classes. In the inversion, real-time remote sensing data are imput to calculate the distance between pairs of neurons, and the missing projection coefficients are extracted from neurons from the best matching class with the minimum distance from input data [28]. The training of the SOM algorithm can be described as follows: 1) Use linear initialization, where the i-dimensional weight vector (prototype vector) [m 1 , . . . , m i ] is initialized in order along the linear subspace formed by the two main feature vectors of the input dataset. The feature vector can be calculated using the Gram-Schmidt method. 2) Randomly select a sample vector x from the input training dataset, and calculate the degree of similarity between it and all weight vectors on the map. The best matching unit (BMU) is denoted by m c , and its weight vector has the highest similarity with the input sample x. This similarity is defined by the Euclidean distance metric, and the BMU is defined as a neuron: 3) After finding the BMU, use the batch algorithm for the prototype vector of the SOM to simply replace it with the weighted average of the samples, where the weight factor is the value of the neighborhood function: where t denotes time, c(j) is the BMU of the sample vector x j , h c(j) i is the neighborhood function (weighting factor), and n is the number of sample vectors. 4) Select another learning data input to the input layer of the network, and return to step 2 until all input data have been provided to the network. 5) Let t = t + 1, return to step 2, and stop when training times T are reached. Once the SOM has been trained, The missing paremeters can be then obtained by extracting the corresponding dimensions from the BMU. Chapman proposed a formula that uses these data to find the BMU with the closest Euclidean distance in the SOM, and then employs this unit to fill in the missing profile parameter [28]: where X is the input data item, p is the index of each class, Y is the reference vector (BMU), d E is the Euclidean distance between the input vector and the SOM unit, C i,j is the correlation matrix between the input data and the output data, a is the set of known information, and b is the set of unknown information. Finally, the SSP is reconstructed by obtaining the set of projection coefficients.

III. EXPERIMENT OF SSP RECONSTRUCTION A. DATA
The SSP samples were obtained from Argo buoy data. The international Argo plan was implemented in 2000, VOLUME 9, 2021  whereby more than 30 countries and international organizations currently deployed roughly 4000 automatic profile buoys throughout the world's oceans to measure temperature and salinity profiles. The Argo buoy data used in this paper were from the Global Ocean Argo Scattering Data Collection (V3.0) (1997/07-2020/06) (ftp://ftp.argo.org.cn/pub/ARGO)/global/) [29]. Del Grosso's empirical formula for sound speed was used to convert the temperature, salinity, and pressure data into a SSP [30]. For matrix processing, the inversion methods require the same depth sampling in profiles. As the number of samples decreases very quickly due to exclusion of shallower profiles when the maximum depth exceeded 900 m, a maximum depth of 900 m was chosen in order to include more variability of sound speed and SSP samples. All the profiles were   The background profile employed the climatic profile data of the WOA13, which contained climatic hydrological data on the ocean (https://www.nodc.noaa.gov/OC5/woa13/) provided by the US National Oceanic and Atmospheric Administration's Ocean Data Research Center (NODC). It consisted of averaged data from a variety of datasets on the global temperature, salinity, and density. The data are divided into annual average data, seasonal average data, and monthly average data. The spatial resolutions provided are 5 • , 1 • , and 0.25 • . We used the annual average data from 1995 to 2012 at a spatial resolution of 1 • × 1 • . The SSP anomalies were obtained by subtracting the climate profile from the SSP samples, that is, the buoy profile minus the WOA13 profile.   Fig. 3. This area is located in the northern basin of the South China Sea. As shown in Fig. 4, a total of 490 profiles were finally obtained. The training dataset consisted of 439 profiles from 2008 to 2016, and 51 profiles from 2017 were used as the test set. Once the map had been generated, it was used to iteratively train the algorithm by gradually extending the duration of training and reducing the radius of influence.

B. EOF ANALYSIS
The EOF provided a large number of SSP samples with a compact description by finding spatial patterns that explained much of the variance in the SSP perturbation. The first five shape functions are shown in Fig. 5. From the perspective of modal shape, the perturbation in sound speed mainly occurred in the range of depth of 600 m, and the amplitude was close to zero at around a depth of 900 m. We focused on reconstruction within 900 m. Table 1 shows the effects of reconstruction of the first five EOFs. The rates of contribution of the first five eigenvectors to the variance were 81.35%, 8.35%, 4.10%, 1.81% and 1.21%, respectively, and their cumulative contributions to the rate of variance wsa 96.82%. The samples were reconstructed using the first five modes, and reconstruction values were compared with the measured values. The root mean-squared error was 0.42 m/s, indicating that the EOF shape function can adequately represent most of the variance in perturbation in the region, thus ensuring a reasonably accurate reconstruction.

IV. SSP RECONSTRUCTION
To evaluate the reconstruction effect of SOM-based inversion and sEOF-r, the root mean-squared errorσ is was used: where c m is the measured SSP, c r is the profile reconstructed by inversion that varies with the depth z and time t, M corresponds to the number of depth sampling points, and N is the total number of SSP samples of the test set. A comparison of the reconstruction results of the two methods is shown in Fig. 6. Based on the same EOF mode and the same number of samples used to construct the relationship between the parameters of the sea surface and the projection coefficients, the accuracy of reconstruction of the SOM method was significantly higher than that of the sEOF-r. The SSPs reconstructed by using the SOM had a mean error of 1.19 m/s with a coefficient of determination R 2 of 0.99, which can satisfy the needs of conventional applications. The error in the linear regression algorithm fluctuated significantly, and most of the errors were larger than 2 m/s. The average error of linear regression was 2.07 m/s with a R 2 of 0.98. The errors of reconstruction at different depths are shown in Fig. 7. Because the SSP near the surface is directly controlled by parameters of the sea surface, the linear relationship which is based on the physical mechanism was more effective near the sea surface and the sEOF-r errors were smaller near the sea surface. The optimization processing of the SOM resulted in low errors in the thermocline, but increased the error near surface. The large errors for the linear regression algorithm were obtained at depths close to 500 m, which may be caused by the water exchange between South China Sea and Pacific Ocean on the East of Luzhou Strait [31]. As the linear regression was based on a large number of samples, differences were obtained between individual and statistical characteristics that led to uncertainty in inversion using sEOF-r, especially in the experimental sea area that featured a strong water transport. Fig. 8 shows the reconstructed samples of the SSP. The results are consistent with those in Fig. 6 and Fig. 7. The result of the SOM method was thus in good agreement with the actual profile, whereas the linear method yielded larger errors in many samples. Under the constraints of linear inversion, even if the statistical law had an exact statistical relationship with (4), the differences between individual characteristics and the statistical law of a large number of samples inevitably introduced errors. In the SOM method, each sample had a unique transformation relationship and no linear constraint. The results show that the error introduced by the inversion relationship was relatively small. Fig. 9 shows the relationship between the inverted projection coefficients and input remote sensing parameters. The distributions of both show that they did not have a strict linear relationship, and the actual relationship between the EOF coefficients and remote sensing parameters was rather nonlinear. The results of inversion of the SOM method were significantly better than those of linear regression. The remote sensing coefficients enabled the accurate inversion of the first two orders of coefficients, but that of coefficients of the third order and higher incurred a certain error. The sea surface temperature and height affected changes in the first five orders of the coefficients, that is, the first five orders accounted for a large part of the variance in the results. The analysis of the neural map showed that the SSP can be retrieved using the remote sensing parameters. In practicality, inversion of five shape modes is a proper compromise between accuracy and noise avoidance.

V. CONCLUSION
In this paper, a nonlinear inversion method based on self-organizing map is proposed. The remote sensing parameter was used as input to the BMU to obtain the projection coefficients for SSP reconstruction. In tests on data from the northern part of the South China Sea, the SOM-based yielded higher accuracy than the sEOF-r method.
The analysis of the neuronal topology of the SOM also preliminarily explained the perturbation transfer-based relationship between parameters of the sea surface and the subsurface profile. The linear relationship embodied in a large amount of data was not strictly applicable to some sea areas. At the same time, the first five orders of the projection coefficients of the EOF had a higher correlation with the parameters of remote sensing. Choosing the first five order projection coefficients can satisfy the accuracy-related requirements of the method of inversion without introducing excessive noise in the highorder coefficients.
The proposed method also has certain limitations. The method used to acquire the EOF to establish the inversion relationship is highly dependent on the number of samples, which limits its general application. The machine learning method can reduce the limitations of simple linear fitting so that the relationship between the parameters can be explored more accurately. Moreover, as the ocean features the complex effects of multiple parameters, the neuronal structure is conducive to introduce other parameters, such as heat flux and wind speed. Currently available methods, including the one proposed here and the sEOF-r, are statistical, and require the introduction of physical mechanisms to further improve their performance.
HAIPENG LI is currently pursuing the B.S. degree with the College of Electronics and Information Engineering, Guangdong Ocean University, China. His current research interests include neural networks and ocean color remote sensing.
KE QU received the M.S. and Ph.D. degrees in signal and information processing from Shanghai Acoustic Laboratory, Chinese Academy of Sciences, Shanghai, China, in 2011 and 2014, respectively. He is currently an Associate Professor with the School of Electronics and Information Engineering, Guangdong Ocean University, Zhanjiang, China. His current research interests include acoustic inversion and joint ocean-acoustic modeling.
JIANBO ZHOU received the B.S. and Ph.D. degrees in underwater acoustic engineering from the School of Underwater Acoustic Engineering, Harbin Engineering University, Harbin, China, in 2013 and 2018, respectively. He is currently holding a postdoctoral position with the Department of Acoustic and Information Engineering, School of Marine Science and Technology, Northwestern Polytechnical University. His current research interests include modeling and application of ocean ambient noise and related fields.