Estimating Bathymetry from Multi-Mission Satellite Altimetry Data Using Gravity-Geologic Method over Malaysian Seas

Altimetry Data Using Gravity-Geologic Method over Malaysian Seas Astina Tugi, Ami Hassan Md Din 1, , Nornajihah Mohammad Yazid, Abdullah Hisam Omar, Amalina Izzati Abdul Hamid and Noor Anim Zanariah Yahaya, Geospatial Imaging and Information Research Group (GI2RG), Geoscience and Digital Earth Centre (INSTeG), Faculty of Built Environment and Surveying, Universiti Teknologi Malaysia, 81310 Johor Bahru, Johor, Malaysia. Geomatics Innovation Research Group (GnG), Faculty of Built Environment and Surveying, Universiti Teknologi Malaysia, 81310 Johor Bahru, Johor, Malaysia. This manuscript is submitted to "Terrestrial, Atmospheric and Oceanic Sciences"


INTRODUCTION
Topography is essential in understanding the process of the Earth. It varies from small mountain valleys to large continental landmasses, which causes weather and climate variations. Due to tectonic activity, erosion, and sedimentation transfer, changes in the land have stimulated the need for detailed topography to investigate geological occasions. While, in the ocean, marine administrations can be organised, and marine geology, biology, and physical oceanography can be discovered with detailed bathymetry information given (Sandwell et al., 2001;Rosmorduc et al., 2006;Hell, 2011;Annan and Wan, 2020). Hence, it is vital to have knowledge of ocean bathymetries.
Traditionally, accurate mapping of the bathymetric ocean floor uses echo-sounding approaches. In recent days, modern multi-beam echo sounder (MBES) techniques have made the conventional single-beam echo sounder (SBES) obsolete whereby, the use of MBES has greatly improved the accuracy, efficiency, and spatial resolution of coastal and ocean mapping (Hell, 2011). Nevertheless, this technique is challenging to use in mapping large ocean floor areas, as it is time-consuming (Carron et al., 2001;Sandwell and Smith, 2001;Smith et al., 2005;Jena et al., 2012;Minzhang et al., 2014). According to Jena et al. (2012), collecting MBES bathymetry data for unexplored offshore areas is a challenging task. As airborne, space-borne optical, or hyperspectral sensors can reveal bottom topography in shallow areas, and these systems are applicable only at a depth of < A c c e p t e d M a n u s c r i p t > Page 5 in 66 pages water less than tens of metres (Smith et al., 2005;Kim et al., 2019). Optical images must be analysed using the attenuation of sunlight in water, the reflectance of the bottom of the ocean, and water properties in order to predict depths using optical images (Hsiao et al., 2016;Hugue et al., 2016;Tourian et al., 2017;Kasvi et al., 2019). All in all, optical images can obtain around 20 m of maximum bathymetry depth, according to the previous studies.
Ocean surface height anomalies for uniform reconnaissance of deep-sea floor topography and bathymetry modelling on a global scale can be retrieved from space-borne radar altimetry (Smith et al., 2005;Minzhang et al., 2014;Tugi et al., 2016;Breda et al., 2019;Annan and Wan, 2020).
Currently, the sounding of ships and satellite altimetry gravity anomalies has constructed all highresolution global bathymetry models, which depend on the 20-200 km waveband gravity anomalies. Therefore, researchers must be cautious when analysing the isostatic mechanisms of the seafloor (Minzhang et al., 2014). According to Xu et al. (2009), knowledge of the global ocean before the employment of satellite altimeter missions is spatially and temporally separated by scattered observations. Subsequently, this reflected insufficient information in the global components of ocean observation. However, the measurement of sea surface height from global ocean circulation can be obtained reliably and consistently with the implementation of satellite altimeter measurements. Satellite gravity missions have provided information on the gravity of the Earth, allowing the derivation of marine gravity anomalies to explore the ocean basin (Yildiz, 2012;Sandwell et al., 2014;Yazid et al., 2016).
The measurement of the sea surface height is the concept of satellite altimeter, which can be determined from the time back and forth of the microwave pulse emitted by the satellite altimeter to the sea surface, and the power of the microwave pulse radiation is acknowledged (Chelton et al., 2001;Din et al., 2014;Zulkifle et al., 2019;Hamden et al., 2021). The primary measurement taken by the satellite altimeter is the measurement range, satellite height above the sea surface height as well as the altitude of the satellite with relation to the ellipsoid (orbital height). Current progress in satellite altimetry has led to improvements in high-resolution marine gravity fields (Andersen et al., 2010;Yazid et al., 2016) and global bathymetric models, which provide refined depth resolutions for the South China Sea (SCS) (Sandwell et al., 2014). This paper focuses on ocean floor bathymetry estimation for the Malaysian Seas using the space-borne technique, namely, satellite altimeters, to derive the marine gravity anomalies. Subsequently, the ocean floor bathymetry is estimated from the derived marine gravity anomalies. To the best of the authors' knowledge, no such study has been done in the Malaysian region, particularly in estimating bathymetry from a space-based approach using Gravity-Geological Method and comparing it to Global Bathymetry Models.

Study Area
The research area in this study encompasses the Malaysian Seas, i.e., the Malacca Straits, South China Sea, Sulu Sea, and Celebes Sea (refer to Figure 1). The limitations of the study area are ranged from the latitude and longitude of 0º 0' 0" N to 14º 0' 0" N and 95º 0' 0" E to 126º 0' 0" E, respectively.

Satellite Altimeter Data
Several satellite altimeter missions are found globally, but only six satellite altimeter missions, namely ERS-2, Jason-1, Envisat-1, Jason-2, CrySat-2, and SARAL, are used in this paper. For satellite altimetry data processing, the Radar Altimeter Database System (RADS) is utilised to cover altimetry data from 2005 to 2015. RADS also provides the measurements of satellite orbit altitude and altimeter range from the time taken by the return radar pulse from the satellite antenna returns to the satellite receiver to achieve the sea surface height (SSH) in RADS (Yahaya et al., 2016). The MSSH is then calculated from the SSH processing. A multi-mission satellite altimeter approach was applied in this study, where altimeter missions were used, and its information is listed in Table 1. The altimeter has conducted the measurement by emitting a short radiation pulse to measure the required time to reflect the pulse from the sea surface. This measurement refers to the altimeter range, R, and this altimeter range provides the height of the instrument above the surface of the sea.
From the round-trip travel time of the pulse, an estimated value of R from satellite to mean sea surface (MSS) can be acquired. Figure 3 illustrates the principle of satellite altimeter in deriving sea level data.
The data extraction for SSH can be completed after the required correction has been made. Two corrections are found in the study, namely range and geophysical correction, to obtain the SSH (Din et al., 2014;Abazu et al., 2017). The range correction for obtaining the SSH was automatically applied in this system by using Equation 1. In geophysical correction, the effect of tidal height variations, hT, and the ocean surface response to atmospheric loading, ha, play a significant part in obtaining the SSH. There are four components to the tidal height variations effect: ocean tide, solid Earth tide, pole tide, and load tide. As stated by Andersen and Scharroo (2011), the ocean tide is regarded to be the most critical factor in generating fluctuations in the ocean surface compared to other effects, and it is caused by the attraction forces of the astronomical bodies (moon and sun) towards the Earth, as well as the Earth's gravitational pull and rotation. Like the ocean tides, the solid Earth tide compensates for temporal variation in land and ocean surface elevations caused by tidal deformation of the underlying non-rigid Earth, including the ocean floor, due to the moon and sun attraction (Fok, 2012).
On the other hand, the pole tides a tidal effect of both the solid Earth and the oceans to the centrifugal potential caused by minor disturbances in the Earth's rotating axis, the polar motion which occurring durations mainly of 365 and 433 days (Desai, 2002 andFok, 2012). Furthermore, the load tide is known as a correction for variations in height caused by changes in the tide-induced forces operating on the Earth's surface and can be computed using mathematical models such as the GOT4.10 model in RADS. The effects of tidal height variations, hT (solid earth tide, pole tide, ocean tide, and load tide effects), and the ocean surface response to atmospheric loading, ha, must be considered in order to achieve the dynamic sea surface height, SSHd (see Equation 3). Equation 2 is the simplification of Equation 3.

SSHAlt = h -(R' + ΔRDry + ΔRWet + ΔRIon + ΔRSSB)
(1) where R = Altimeter Range  Table 2. Crossover adjustment is carried out to limit errors in the satellite track and other long wavelength errors by minimising height differences between the ascending and < A c c e p t e d M a n u s c r i p t > Page 11 in 66 pages descending tracks at a crossover location (Hamid et al., 2018;Hamden et al., 2021). Subsequently, data filtering and data gridding are also implemented. The data provided by RADS is daily solution data is then processed to obtain climatology data in the form of average data. Afterward, the data obtained is filtered to ensure that the obtained data is free from errors of un-modelled tides, orbit error, or residual sea surface height variability (Andersen, 2008;Din et al., 2019) and then gridded.
This paper employed a filter value of 2.0 and a block size or cell size of 0.25 degrees for the data grid.

Marine Gravity Anomaly from Satellite Altimeter
Three essential steps are employed to obtain the final output of the derived Free Air Gravity Anomaly (FAGA) from satellite altimetry data, which are the selection of the appropriate global After the geoid computation is completed, the gravity anomaly is derived using the Fast Fourier Transformation (FFT) technique.
According to Salam (2005), the gravity anomaly can be derived from the MSSH using Fast Fourier Transformation (FFT) technique and to implement this technique, GRAVSOFT software is used to derive the gravity anomaly of satellite altimeter under the GEOFOUR program. The FFT formula used in deriving the MSSH to the gravity anomaly is denoted in Equations 4 and 5 (Schwarz et al., 1990;Forsberg and Tscherning, 2003;Sandwell, 2004;Salam, 2005) as follows: From Equations 4 and 5, the value of u, v, and w represent the frequency domain data of the spatial coordinate of x, y, and z, while i is the imaging unit (i= 1  ). The symbol F represents the FFT function, while the F -1 is the inverse FFT function. The FFT technique is implemented in the GEOFOUR program to derive the gravity anomaly to the geoid or vice versa. The information of MSSH can be assumed as the geoid undulation, while N can be estimated using the Stokes Integral (Heiskanen and Moritz, 1967;Salam, 2005). The FFT implemented in the GEOFOUR program is the planar Stokes Integral, which can be clarified using Equation 6 as follows:

< A c c e p t e d M a n u s c r i p t >
where Δg = the gravity anomaly go = the average of the gravity or the normal gravity The cross-validation process is applied for the chosen satellite-derived free-air gravity anomaly (FAGA) for further computation. Hence, in this study, Surfer Software has been used to execute the cross-validation process in calculating the gridded error of the observation. The original data set is then given the known values for the geoid (N) observation location, and this information has been used to evaluate the grid error for the relative quality of the grid data (Golden Software, 2002). The gridded error is calculated by eliminating the first data set observation, and the remaining data are used to interpolate the value of the eliminated observation by using a specific interpolation algorithm at the same location of the deleted point. The new satellite-derived FAGA (filtered satellite-derived FAGA) dataset will be used to estimate the bathymetry after the interpolation process. Figure 5 shows the flowchart of the gravity anomaly derivation from the satellite altimeter.
NGGM represents the geoid from the global geopotential model (GGM), while Δg, Δgr, ΔgGGM, and MDTDTU15 represent the altimeter-derived gravity anomaly, the residual gravity, the gravity anomaly from the GGM model, and the mean dynamic topography model (MDT) 15 from DTU, respectively.
In this study, three filtering processes are used with seven different interpolation methods to acquire better results. Residual filtering is classified in Table 3,  This interpolation is required in the grid creation of the long-wavelength gravity field, represented by the generation of regional gravity anomalies gREG (i) using the Gravity-Geologic Method in the bathymetry estimation. Therefore, the interpolation method has to be selected, and in doing so, seven interpolation methods are tested to narrow down the selection of interpolation.

Computation of Bathymetry using Gravity-Geologic Method
A total of 12,362 points of shipborne bathymetry obtained from the National Geophysical Data Centre (NGDC) are employed to estimate the bathymetry with 6584 points used in the computation, while the remaining 5778 points are used in the validation process. A total of 6584 points of shipborne bathymetry are used to compute the residual gravity by obtaining the regional gravity for the study area. Figure 6 shows the shipborne track of bathymetry, which was divided into two different colours. The green point data is utilised in the computation process, while the yellow points represent the data used in the validation process. The computation and validation data classifications were chosen randomly according to the coverage of the shipborne area. In this paper, the Gravity-Geologic Method is used to predict the bathymetry. The input data for the bathymetry prediction is the derived-FAGA from multi-mission satellite altimetry data.
According to Kim et al. (2010), the residual of the gravity effects is essential information for the acquisition of short-wavelength signals, and these signals are used to predict the depth of the bedrock. From this computation, the depth undulation can be predicted. Figure 7 shows the Gravity-Geologic Method geometry and the computed parameter related to the computation (Hsiao et al., 2010).
Two types of gravity field signals are produced, namely the shorter wavelength and the longer wavelength. The information from the shorter wavelength gravity field, residual gravity (gRES), is generated based on local depth variations from the shipborne measurement (control points). The deepest depth of the shipborne measurement (D) represents the reference datum used in the Bouguer plate measurement. Meanwhile, a longer wavelength gravity field, regional gravity (gREG) is produced from the deeper mass variations. Then, the final depth was computed from the residual gravity (Hsiao et al., 2010;Wei et al., 2021).
Besides, by combining residual gravity (gRES) with the regional gravity (gREG), the observed Bouguer gravity anomaly (gOBS) can be obtained. Equation 7 shows the gOBS linear computation and the connection between gREG and gRES. The subscript j indicates the control points with the measured depth; while subscript i represent the unmeasured depth of control points.
According to Kim et al. (2010), the measured depth at the jn control point is assumed to be used in the estimation of the residual gravity field of the bedrock (gRES), which produces the effect of the shorter wavelength from a simple Bouguer slab formula listed in Equation 8. Figure 7 shows the geometry of the Gravity-Geologic Method and the downward continuation method.
where; G = Gravitational constant, 6.672 × 10-8 cm 3 /g sec 2 Δρ = Density contrast between seawater and bedrock (g/cm 3 ) E (j) = Depth at the j th control point (in meter) D = The deepest depth of the control points as a reference datum (in meter) The importance of the control points is to determine the shorter wavelength effect in order to generate the regional gravity at the location of the control point. The residual of gravity represents the effect of the bedrock surface at the control point. This residual has been eliminated from the observed gravity Bouguer to obtain the estimated regional gravity (gREG) that will represent the effect of the longer wavelength. Equation 9 shows the computation of the regional gravity field at the control point (gREG (j)).
After gREG (j) has been acquired, the regional gravity of the unmeasured depth i, gREG (i) can be determined by interpolating the regional gravity gREG (j) generated at the measured depth, j (control points), creating a gridded regional surface for the longer wavelength effect. Subsequently, the regional gravity that has been estimated at site i (gREG (i)) is deducted from the observed gravity, gOBS (i) in order to achieve the estimated residual of the gravity (gRES (i))(see Equation 10).

gRES (i) = gOBS (i) -gREG (i)
A mass change below the ocean floor has caused the density contrast between seawater and the ocean floor to vary differently. It is known that the density (ρ) is 1.03 g/cm 3 and 2.70 g/cm 3 for the seawater and ocean floor bedrock. In the application of the Gravity-Geologic Method, the density contrast between the seawater and the ocean floor bedrock (Δρ) with the assumed value of 1.67 g/cm 3 is applied. According to Kim et al. (2010), the ocean floor bedrock density in the crust can be assumed to be between the mean range of 2.67 g/cm 3 and 2.73 g/cm 3 . Essentially, the theoretical value of the density contrast (1.67 g/cm 3 ) may be higher than the specified value, as the derived-FAGA from satellite altimeter that has been compiled on the sea surface in Gravity-Geologic Method is affected by the ocean bottom sources (Roman, 1999;Kim et al., 2010).
An assessment of the global bathymetry model has been executed before the bathymetry computation using the Gravity-Geologic method is done where it is executed using the validation points obtained from NGDC (5778 points) of the shipborne bathymetry. In addition to that, this assessment has been conducted to determine the ocean depth differences between each model and apart from the validation points from the shipborne bathymetry. In this assessment, four bathymetry models, namely DTU10, ETOPO1, GEBCO, and Sandwell V18.1, are extracted from the global bathymetry models. The validation points recorded in this study are around 5778 points, which subsequently followed the validation points of the shipborne bathymetry location.

Statistical Analysis of the Global Bathymetry Models and the Shipborne Bathymetry
The data analysis of the global bathymetry models and the shipborne bathymetry from NGDC (refer to Figure 6) are listed in Table 4, including the minimum and maximum depth, the average of the data, and the standard deviation for both types of data used. Table 4 shows that the minimum depth value for the shipborne bathymetry data is -9.00 m, and the standard deviation of the data for the

Statistical Analysis of the Type of Filter and the Interpolation Method
Earlier in Section 2.3, seven interpolation methods were used to grid the satellite altimeter FAGA. Figure 8 shows the illustration of the gridded satellite altimeter FAGA. As shown in Figure 8, the Local Polynomial method (Figure 8c) is not appropriate to be employed for deriving the satellite-derived FAGA as it only gives values to grid nodes using a weighted least-squares fit on data contained inside the grid node's search area (Golden Software, 2002). The gridded factor used was the First Order polynomial with a power of 2. As for the computed bathymetry using Gravity-Geologic Method, the estimation is made according to the type of filter used, with a different interpolation method effect, with a cell size of 0.1º x 0.1º (~11 km). The interpolation size is chosen to be smaller than the original data size of SSH from the satellite altimeter, which is 0.25°, in order to get a better resolution of the estimated bathymetry. This assessment is executed to visualise the depths difference and the distinction between the bathymetry models and the ground truth data used (NGDC shipborne data).
All the statistical analyses explained in this section are done to examine the accuracy of the estimated bathymetry obtained using the Gravity-Geologic Method. Figure 9,  Nevertheless, the minimum curvature is chosen since it gives a better RMSE value.
Surfer software executes the interpolation process. Several explanations concerning the minimum curvature interpolation are suggested in the literature, where according to Yang et al. (2004), the use of the minimum curvature is common in the earth sciences. From this interpolation, the resulting surface is analogous to a thin, linearly elastic plate that passes through each data value with a minimum amount of bending. Apart from that, previous studies conducted by Yang et al. (2004) and Vohat et al. (2013) suggest that interpolation produces the smoothest possible surface while respecting the data as closely as possible.
The grid from this interpolation is generated by repeatedly applying an equation over the grid in order to smooth the grid (Vohat et al., 2013). Hence, this method is believed to be the most suitable interpolation method here. However, the interpolation method yields the poorest RSME values in the radial basis function. All filters that are applied to radial basis function interpolation gave RMSE values of more than 200 meters. Conversely, the radial basis function fits with the source data, producing a smooth surface. Despite that, minimum curvature interpolation has been chosen as it practically yields the best-estimated bathymetry output in this study, and RMSE value is used as the indicator in choosing the best-estimated bathymetry. This minimum curvature method is also known as exact interpolators, as they seek to respect the data used in the study (Yang et al., 2004;Vohat et al., 2013).
From previous studies, minimum curvature interpolation has been mainly used to estimate the bathymetry using the Gravity-Geologic method (Sandwell and Smith, 2001;Roman, 1999;Kim et al., 2011). Therefore, towards selecting the final estimated bathymetry, the minimum curvature with Filter 1 is chosen as the best output. Concerning the final estimated bathymetry, the estimated bathymetry Filter 1 from minimum curvature interpolation, an alternative name has been made, namely the UTM18 Bathymetry Model. The UTM18 Bathymetry Model has been used to simplify the final estimated bathymetry in this study and to match with other global bathymetry models.
Since the best-estimated bathymetry goes to the minimum curvature with the implementation of Filter 1, the gravity anomaly information used in the computation is charted. Figure 10 illustrates the three components of gravity used in this study, namely the regional, residual and observed gravity anomaly. The comparison of all these components are plotted to visualise the pattern of gravity anomaly in order to calculate the bathymetry specifically at the latitude of 6.5ºN and along the longitude 108ºE to 112ºE. Analysis along the longitude was made in accordance with the variability of the seabed from shallow water to deep water. Figure 10 shows the g_obs(i) and the computed gravity anomaly, namely g_reg(i) and g_res(i), in which these three components of gravity anomaly are used in estimating the bathymetry. Based on the line graph plotted, the g_res(i) represents the short-wavelength information since it is computed

< A c c e p t e d M a n u s c r i p t >
Page 23 in 66 pages based on the observation (g_obs(i)) and the regional gravity (g_reg(i)) data. Nevertheless, g_reg(i) is smoother than the residual gravity (g_res(i)) since the regional gravity represents the longwavelength effects (Hsiao et al., 2010). With the readiness of these three gravity components, the production of the estimated bathymetry can be done smoothly.
The density contrast value used in the computation is approximately 1.67 g/cm 3 . With the implementation of the Gravity-Geologic method computation, a stronger and homogeneous density contrast between the seawater and the bedrock is applied. The depth to the marine density contrast also tends to be more profound. Therefore, the marine Bouguer slab calculations are not significantly affected by the effects of subsurface density errors (Kim et al., 2011). Figure 11 (a) illustrates the best-estimated bathymetry using the Gravity-Geologic Method in which it is executed in a small area at the South China Sea, since this area has recorded several shipborne bathymetry data in various depth levels, from shallow to the deep ocean. Regarding Figure 11 (b), the bathymetry computed shows the maximum depth of approximately -3400 metres, while the minimum depth is less than -100 metres. The bathymetry is deeper as the longitude increases, from the west to the east part of the South China Sea, and is considered the deep ocean. Figure 11 (b)

< A c c e p t e d M a n u s c r i p t >
Page 24 in 66 pages reveals that the deep basin of this study area lies between 110ºE to 114ºE. The statistical analysis of the estimated bathymetry afterward is then evaluated.
The assessment of the selected estimated bathymetry is executed by comparing the UTM18 Bathymetry Model with the global bathymetry models to examine their depth differences with the ground truth data (validation data). Figure 12 shows the bathymetric profile along 6.5ºN latitude and from 108ºE to 112ºE of longitude. In comparison, 4 global bathymetry models found in this study are used: ETOPO, GEBCO, Sandwell V18.1, and DTU10. From the results collected, the estimated bathymetry from Gravity-Geologic Method seemed to fit better with the NGDC shipborne bathymetry compared to other bathymetry models. This is because the data used for the computation is taken directly from the shipborne and considered the high-frequency bathymetry (short-wavelength) (Roman, 1999).
Moreover, the bathymetry predicted from Gravity-Geologic Method was influenced by the variation of the seabed terrain and the relative depth of the predicted area. Based on Figure 12, predicted bathymetry at longitude 111˚E to 112˚E shows the largest deviation compared to another longitude.
Meanwhile, the bathymetry predicted using Gravity-Geologic Method (UTM18) is closer to NGDC data than the other global models. As the Gravity-Geologic method used the deepest depth from the computation data (NGDC depth) as a reference datum, this information was relatively used in predicting bathymetry. The increase of the depth may improve the predicted bathymetry. Another

< A c c e p t e d M a n u s c r i p t >
factor for enhancing the predicted bathymetry is the estimation of the density contrast between the seawater and seabed topography.
Correspondingly, from the results shown in Figures 10 and 12, the pattern of the g_res(i) and the estimated bathymetry from the Gravity-Geologic Method are reasonably similar since the gravity starts to decrease as the bathymetry goes deeper. The depth variation for the UTM18 Bathymetry Model and the global models is listed in Table 4, presenting the minimum, maximum, mean, and standard deviation of the bathymetries. This information is described to capture the variation of the Gravity-Geologic Method bathymetry with the NGDC shipborne and the global models used.
Importantly, the statistical analysis of the differences between the estimated bathymetry and the global bathymetry models with the ground truth data (NGDC data) is presented in Table 5 Based on the results shown in Table 5, it can be concluded that the estimated bathymetry from the UTM18 Bathymetry Model yields a better RMSE value compared to other bathymetry models. The  Table   6. The analysis included the minimum, maximum, and mean differences, standard deviation, RMSE, and correlation. Additionally, these analyses were executed to attest to the consistency of the UTM18 Bathymetry Model with other global bathymetry models for bathymetry mapping.
Based on Table 6, the largest depth difference is recorded from the ETOPO1 bathymetry model, The histogram (in Figure 13) and scatter plot (in Figure 14) on the depth differences are also plotted, which subsequently illustrates the separation between the estimated bathymetry with the existing global bathymetry models, as well as the validation points (checkpoints). Figure 13 consists Generally, most of the depth differences lie between ±500 m. Nevertheless, according to the final result of the computation, the estimated bathymetry of UTM18 shows a good correlation with the NGDC shipborne bathymetry. It can be seen that the depth differences for the UTM18 Bathymetry Model lie between ±200 m, but most of the differences are in the range of ±50 m. There are more than 2500 points, which are recorded to have a depth difference between mentioned ranges. When compared to another global bathymetry, it appears that the depth differences from Sandwell global model lie between -200 m to 500 m. However, more than 2500 points have small depth differences compared to other global models. As for ETOPO1 bathymetry, some of the points have depth differences that exceeded ±500 m.
From the depth differences recorded, their distributions corresponding to the depth are plotted on the scatter plot in Figure 14. The distributions of the differences are also combined into the exact figure for the sake of understanding the distinction between each bathymetry model. As shown in Figure 14, most of the point differences lie between ±500 m. However, the difference recorded is significant between the depths of -500 m to -3000 m, especially for the Sandwell and ETOPO1 bathymetry models. Overall, the separation of the depths lies between ±500 m.
As mentioned earlier in this section, the estimated bathymetry in this study has a relatively large RMSE value, which is ±96.949 m. This may be due to the dynamic surface of the ocean may have influenced the accuracy of the estimated bathymetry due to the fluctuation of the ocean.
Subsequently, the comparison between the global bathymetry models is conducted to evaluate the estimated bathymetry computed based on Gravity-Geologic Method (UTM18 Bathymetry Model), not to contend with the existing global bathymetry model. From the scatter plot, it can be summarised that the UTM18 Bathymetry Model has a good relationship with the NGDC shipborne data. Moreover, since 6584 points of shipborne measurement are included in the Gravity-Geologic Method computation, the accuracy of the estimated bathymetry can be enhanced and improvised. As the final output of the study, the map of the bathymetry for the Malaysian Seas has been plotted successfully. Figure 15 illustrates the UTM18 Bathymetry Model, covering the Malaysian Seas. Apart from that, the depth assimilation from Sandwell V18.1 Bathymetry Model is exploited for the area with no shipborne data and terrain area coverage due to sparse shipborne data. Figure 16 portrays the topography of the ocean over the Malaysian Seas, including the terrain area. From this figure, the undulation of the ocean topography can clearly be seen. The maximum depth for each sea is also gathered from the bathymetry computation. Table 7

CONCLUSION
Over the years, the conventional technique in acquiring bathymetry has provided the most accurate bathymetry. However, this technique has its limitation, which can be overcome by using satellite altimeters. Satellite altimeter data can be used to obtain the gravity anomaly and can be used to predict the ocean bathymetry in extensive coverage without consuming the time it takes to obtain the data. This study proves that the gravity anomaly derived from satellite altimeter and satellite gravity mission is reliable in estimating the bathymetry of the ocean. With the combination of gravity anomaly from the satellite gravity missions, the acquired gravity anomaly data are dense, thus improving the result of this study, aiming to get a good estimation of bathymetry.

< A c c e p t e d M a n u s c r i p t >
Page 31 in 66 pages The use of Filter 1 (removal of the > ±100 m residual) also supports the computation of the bathymetry towards a better accuracy. Based on this computation, the RMSE value for the estimated bathymetry is ±96.949 m. Using the minimum interpolation method in the Gravity-Geologic method has improved the computed bathymetry compared to other interpolation methods applied. UTM18 Bathymetry Model has also revealed that its accuracy has improved by 69% compared to the ETOPO1 global bathymetry model. This condition demonstrates the reliability of the Gravity-Geologic method and bathymetry from the space-based technique.
In general, the computation of bathymetry using the Gravity-Geologic method required shipborne bathymetry as the ground reference or as a datum. Nevertheless, the estimated bathymetry produced by using the Gravity-Geologic method is not recommended for the navigational purpose since this computation includes the interpolation procedures, and also, the bathymetry is not referred to as the lowest astronomical tide (LAT). The consistency of the Gravity-Geologic method predictions with bathymetry estimated from other models attests to the effectiveness of the Gravity-Geologic method for bathymetry mapping. Undeniably, the predictions are not unique, and thus, the estimates must compute with care. All in all, the analytical simplicity of the Gravity-Geologic method makes it easy to update the bathymetry estimation once the new and improved data on the FAGA and bathymetric observations have emerged. Individual satellite track and the combination track from multi-mission Satellite Altimeter (ERS-2, Jason-1, Envisat-1, Jason-2, CryoSat-2 and SARAL). 3

Figure legends
Satellite Altimeter measurement in obtaining the sea surface height 4 Flowchart of the RADS Processing. 5 Flowchart of the derivation of gravity anomaly from MSS 6 Shipborne points used to estimate the bathymetry. 7 Geometry of the Gravity-Geologic Method (Hsiao et al., 2010). 8 The gridding method available in Surfer 8. a) IDTP, b) Kriging c) Local Polynomial, d) Minimum Curvature, e) Nearest Neighbour, f) Radial Basis Function g) Triangulation with Linear Interpolation (Units are in mGal).

9
The RMSE value for the estimated bathymetry based on the type of filter and the interpolation method used (green box shows the interpolation method used in obtaining the lowest RMSE value). 10 The regional (g_reg(i)), the observed (g_obs(i)), and the residual (g_reg(i)) gravity anomaly at the latitude of 6.5ºN. 11 Estimated bathymetry with the minimum curvature interpolation from Filter 1; a) Map of the estimated bathymetry, b) The topography surface of the estimated bathymetry. 12 Comparison for bathymetry along the latitude of 6.5ºN. 13 Bathymetry differences between UTM18 Bathymetry Model (using Gravity-Geologic Method) and a) NGDC shipborne, b) DTU10 bathymetry model, c) ETOPO1 bathymetry model, d) GEBCO bathymetry model, e) Sandwell V18.1 bathymetry model. 14 Scatter plot for the distribution of the depth differences. 15 UTM18 Bathymetry Model. 16 Topograsurface of UTM18 estimated bathymetry.

No.
Filtering of Cross-Validation Filter 1 Removing the Δg residual that is more than 100 mGal (> 100 mGal) Filter 2 Removing the Δg residual that is more than 50 mGal (> 50 mGal) Filter 3 Removing the Δg residual that is more than 20 mGal (> 20 mGal)