Optimal interpolation method to predict the bathymetry of Saldanha Bay

Accurate interpolation when compiling bathymetric maps is essential in any water depth study. In the case of Saldanha Bay, continuous dredging operations are constantly altering the ocean floor, which has a detrimental effect on sedimentation and coastal hydrodynamics. If the integrity of the coastline is to be secured, accurate bathymetry predictions would be invaluable in determining the effect of dredging operations on coastal erosion. Inverse distance weighting (IDW) and ordinary kriging (OK) are two well‐known and commonly used interpolation methods to produce surfaces through spatial autocorrelation for numerous applications, inter alia, to estimate bathymetry. This study aims to analyse and compare the efficiency of the IDW and OK interpolation methods to predict the bathymetry of Saldanha Bay. Three comparative interpolation tests, which vary according to the decrease in the quantity of sounding points, are conducted. SPSS statistical software was used to assess the performance of the interpolation methods. Firstly, 2D scatterplots were used to show the correlation between predicted and measured sounding values for each interpolation method. Secondly, analysis of variance was employed to investigate whether the difference between the IDW and OK interpolation methods was statistically significant, and to determine which method was best suited for determining the bathymetry of Saldanha Bay. Findings revealed a strong linear relationship between predicted and measured sounding values for both IDW and OK when 100% of the sounding points are used. Conversely, for medium and small quantities of sounding points, a weak correlation exists. Clear similarities exist in the way that IDW and OK estimate and generate the continuous surface of bathymetry. However, IDW consistently performed better than OK across all interpolation tests. The findings of this study will assist in selecting the most suitable interpolation method for future bathymetry surveys of Saldanha Bay.

coastline is to be secured, accurate bathymetry predictions would be invaluable in determining the effect of dredging operations on coastal erosion. Inverse distance weighting (IDW) and ordinary kriging (OK) are two well-known and commonly used interpolation methods to produce surfaces through spatial autocorrelation for numerous applications, inter alia, to estimate bathymetry. This study aims to analyse and compare the efficiency of the IDW and OK interpolation methods to predict the bathymetry of Saldanha Bay.
Three comparative interpolation tests, which vary according to the decrease in the quantity of sounding points, are conducted. SPSS statistical software was used to assess the performance of the interpolation methods. Firstly, 2D scatterplots were used to show the correlation between predicted and measured sounding values for each interpolation method. Secondly, analysis of variance was employed to investigate whether the difference between the IDW and OK interpolation methods was statistically significant, and to determine which method was best suited for determining the bathymetry of Saldanha Bay. Findings revealed a strong

| INTRODUC TI ON
Since the construction of its breakwater and causeway in the early 1970s, Saldanha Bay has continuously expanded its iron ore, gas and oil export capabilities (Saldanha Bay Municipality, 2021). During the initial port construction in the 1970s, substantial dredging operations took place which altered the ocean floor and increased the average depth by approximately 1.5 m (Henrico & Bezuidenhout, 2020). Existing anthropogenic changes, such as the construction of the Saldanha Bay harbour and consequent continuous maintenance dredging operations, already impacted the structure of the ocean floor. Anticipated future development projects are likely to contribute to structural changes. These changes will affect oceanic forces (e.g., refraction and energy of currents and waves) which will again contribute to coastal erosion within Saldanha Bay. The aim of this study is to analyse and compare the efficiency of the inverse distance weighting (IDW) and ordinary kriging (OK) interpolation methods to predict the bathymetry of Saldanha Bay when applying different quantities of sampling points. Three comparative interpolation tests, with varying quantities of sounding points, were conducted to achieve the aim of this study. The original set of sounding points was divided randomly to create an additional two sets of sounding points. The first data set contained 100% of the original sounding points (n = 1,653), the second only 66% (n = 1,091) and the third 33% (n = 546). This study applied geographic information system (GIS) techniques to interpolate the bathymetry of Saldanha Bay. GIS techniques were assessed using SPSS statistical software. The respective performance and accuracy of IDW and OK interpolation methods were analysed by means of bathymetry map construction, graphic illustration and statistical calculation.
Interpolation methods consider the spatial correlation of known measurement points to estimate locations that lack measurements. Interpolation methods are therefore based on the principle that points in close proximity are more closely related than those that are further apart (Ajvazi & Czimber, 2019). Spatial interpolation is an important GIS technique for creating and analysing spatial surfaces. This study used ArcGIS 10.5.1 to conduct all interpolation tests (GIStandards.eu, 2020). Interpolation tools, which provide various options to conduct surface modelling, are available in the Geostatistical Analyst extension. linear relationship between predicted and measured sounding values for both IDW and OK when 100% of the sounding points are used. Conversely, for medium and small quantities of sounding points, a weak correlation exists. Clear similarities exist in the way that IDW and OK estimate and generate the continuous surface of bathymetry. However, IDW consistently performed better than OK across all interpolation tests. The findings of this study will assist in selecting the most suitable interpolation method for future bathymetry surveys of Saldanha Bay.

HENRICO
Over 40 different types of interpolation methods exist, many of which are integrated into GIS software to create, among others, digital elevation models (Ferrando, De Rosa, Federici, & Sguerso, 2016;Ferreira et al., 2017;Meng, Liu, & Borders, 2013). Only a few of these interpolation methods are widely used (e.g., IDW, kriging, splines, nearest neighbours, local and global polynomials, radial basis functions, topo-to-raster), of which IDW and kriging are the most attractive (Curtarelli, Leão, Ogashawara, Lorenzzetti, & Stech, 2015;Ferreira et al., 2017;Kalivas, Kollias, & Apostolidis, 2013). This study considers IDW and kriging, specifically OK, to identify the optimal spatial interpolation method to use for mapping the bathymetry of Saldanha Bay from sounding measurements. It should be noted that no single interpolation method is ideal for all types of data sources. The accuracy of such estimations is largely dependent on sounding density, distribution and accuracy. Low density and heterogeneous sample points that are inaccurate in terms of spatial location and referencing will negatively influence the accuracy of surface interpolation (Aguilar, Agüera, Aguilar, & Carvajal, 2005;Diaconu et al., 2019;Erdoğan, 2010;Gao, 2009;Setianto & Triandini, 2013).
Numerous studies have compared the efficiency of the IDW and kriging interpolation methods to predict bathymetric surfaces accurately (Azpurua & Ramos, 2010;Bello-Pineda & Hernández-Stefanoni, 2007;Diaconu et al., 2019;Ferreira et al., 2017;Meng et al., 2013;Murphy, Curriero, & Ball, 2010;Zarco-Perello & Simões, 2017;Zimmerman, Pavlik, Ruggles, & Armstrong, 1999). These two methods differ a great deal in their approach towards calculating and predicting the correlation in geospatial data. Furthermore, the choice of interpolation method is dependent on the specific research application and the characteristics of the spatial data. Most studies show that kriging outperforms IDW, especially when the number of sampling points varies and when measurements are sparsely distributed (Diaconu et al., 2019;Ferreira et al., 2017;Zimmerman et al., 1999). However, a few studies have found the IDW method to produce better results (Azpurua & Ramos, 2010;Gong et al., 2014;Meng et al., 2013;Zarco-Perello & Simões, 2017). However, Murphy et al. (2010) found that results are site-specific, that there is "no consensus as to a superior or preferred method". Given these contradictory findings, and the possibility that results may be site-specific, it is imperative to find the most suitable method for interpolation of the bathymetry of Saldanha Bay. Regarding Saldanha Bay, various studies had been conducted on the effects of dredging on changes in sedimentation (Flemming, 1977;Henrico & Bezuidenhout, 2020;Luger, Schoonees, Mocke, & Smit, 1999;Wiese, 2013). Some studies focused on ocean movements and the impact of wind on the Saldanha Bay environment (Flemming, 2015;Weeks, Boyd, Monteiro, & Brundrit, 1991), while other studies concentrated on the characteristics of the water and marine life within the Bay (Probyn, Pitcher, Pienaar, & Nuzzi, 2001;Shannon & Stander, 1977;Van der Merwe & Lohrentz, 2001). Some unpublished reports (dated 2000, 2004, 2010 and 2013) include selective bathymetry mapping, mainly of Inner (Saldanha) Bay. These reports were acquired from the Saldanha Bay port captain (A. Miya, personal communication, 9 September 2019). They consist mainly of bathymetry surveys of selected areas within the Bay, including the harbour approach channels and turning basin.
The Council for Scientific and Industrial Research (CSIR) created these reports, which carry restricted access, for Transnet Capital Projects. No other studies could be found that compared IDW and kriging interpolation methods to determine the most suitable method to predict the bathymetry of Saldanha Bay.
Finding the most suitable interpolation method to map the ocean floor of Saldanha Bay will not only inform local and provincial government of the changes occurring within the Bay, but will also inform other researchers of the preferred interpolation method for future analysis of the bathymetry of Saldanha Bay.

| S TUDY ARE A
Saldanha Bay (hereafter "the Bay") was selected as the area of study (see Figure 1). The Bay is approximately 105 km north-north-west of Cape Town, in the Western Cape province of South Africa and is bordered by the town of the same name, Saldanha Bay, to the west and Langebaan to the south-east. The Bay is the deepest natural harbour in Southern Africa. The town of Saldanha Bay is located within a very unique ecosystem, characterised by West Coast flora and fauna that host various endangered and endemic plant and animal species (CSIR, 2014).
The Bay covers an area of approximately 85 km 2 . It has experienced significant changes to its bathymetry over the years, caused not only by anthropogenic influences, but also by the severe weather conditions experienced from time to time along the West Coast region of South Africa (Flemming, 2015). The West Coast region has a Mediterranean type climate characterised by long, dry summers, with an average rainfall of 270 mm per year that mostly occurs during the months of winter (CSIR, 2014). This region is characterised by a seasonal reversal of wind direction and speed. The on-shore wind blows south-westerly during the summer months (relatively towards the entrance of the Bay) and off-shore north-easterly during the winter months. Wind speed, especially during summer, often exceeds 11 m/s (Transnet, 2018).
Construction of the Saldanha Bay harbour, which included the building of a causeway with an iron ore and oil loading jetty and a breakwater, started in May 1973. The breakwater, which connected Hoedjiespunt with Marcus Island, was mostly constructed from soil collected from the massive dredging operations that took place to build the port.
The construction of the breakwater effectively divided the Bay into Inner Bay and Outer Bay (Smith & Pitcher, 2015).
The causeway further divided Inner Bay into Small Bay to the west of the causeway and Big Bay to the east. It is estimated that approximately 30 million cubic metres of sediments were removed by dredging operations (Zwemmer & Van't Hof, 1979). The port started its operations in September 1976. Recent plans to expand port facilities and infrastructure by developing the marine and ocean services industry should now take effect with the recent procurement of multi-million-dollar investments for this purpose (Saldanha Bay Industrial Development Zone, 2020).
The consequent changes to the bathymetry of Saldanha Bay had a foreseeable effect on wave refraction and direction and tidal movement (predominantly from the south-west and south-south-west). These changes altered F I G U R E 1 Saldanha Bay was selected as the area of interest for this study. The spatial location of Saldanha Bay on the West Coast of South Africa is indicated by the red dot ( • ) on the inset map. Source: South African National Hydrographer | 1995 HENRICO sediment movement and beach sand deposit patterns within the Bay. Because of this, several beaches within the Bay have been and continue to be eroded. These processes are described and illustrated by wave refraction patterns that were reconstructed from aerial photography 1 as before and after images of harbour construction ( Figure 2).
The general wave direction and tidal movement before harbour construction were naturally spread across the entire Inner Bay. Before construction, the strongest refraction was observed in the north-eastern corner of the Bay where coastal and beach formation areas are most exposed. This picture changed considerably after harbour construction, when the direction of wave refraction shifted more towards the north-north-eastern and eastern parts of the Bay. The altered direction of wave refraction shifted the erosion potential towards Langebaan beaches. Here the most exposed coastal and beach formation areas correspond to the erosion area indicated by the red hatch fill in Figure 2. These changes provide the rationale for this study, namely, to identify the most efficient interpolation method to accurately predict the bathymetry of Saldanha Bay, an imperative to satisfy the requirement of continuously analysing the hydrodynamics of the Bay.
Direction changes in wave refraction and energy dispersion, as well as occasional storm surges experienced along the West Coast of South Africa, caused the Saldanha Bay Municipality (SBM) to construct a bolder beach front and a set of groynes to protect the affected beaches and stop coastal erosion (Legg, 2013). However, according to Mr Jaco Kotze, the chairperson of the Langebaan Ratepayers' Association, the 43 million rand spent by the SBM "was money wasted", because the Langebaan beaches are still experiencing notable erosion and are slowly disappearing (Legg, 2013). The unique natural characteristics of Saldanha Bay make it ideally situated to be a world-class freeport service location. This is evident in the recent multi-million-dollar investment to expand Saldanha Bay Port facilities and infrastructure (Saldanha Bay Industrial Development Zone, 2020). The challenge, however, lies in achieving a balance between protecting the Saldanha Bay ecosystem and welcoming the potential economic growth associated with development projects. Some of these development projects will influence the bathymetry significantly, and will ultimately put more strain on the Saldanha Bay coastline and beaches. In the case of Saldanha Bay, an accurate and current representation of the ocean floor is particularly important for the continuous assessment of the rate of anthropogenic erosion along the Saldanha Bay coastline, and its impact on coastal morphology. This will allow relevant authorities to take timely corrective and preventive actions to protect and sustain the beaches of the Bay.
F I G U R E 2 Wave refraction patterns and exposure zones before and after harbour construction, at a scale of 1:50,000. Source: Adapted from Flemming (2015, p. 63)

| ME THODOLOGY
The entire bathymetric data set used during this study (n = 1,653) was received from the South African National Hydrographer. These soundings were collected during the last known survey of the entire Saldanha Bay in 1977.
The lack of updated bathymetric data on Saldanha Bay was communicated to the Hydrographer during a visit to the South African Navy Hydrographic Office in Tokai, Cape Town, on 2 July 2020. The Hydrographer recognised the urgency of new and updated bathymetric data on Saldanha Bay. Subsequent to the meeting, it was planned to conduct a complete survey of the entire Saldanha Bay by the end of 2020.
Through both the IDW and OK methods three interpolation comparison tests consisting of different quantity soundings were conducted to create bathymetric profiles of Saldanha Bay. Each interpolation test aimed to determine the efficiency of the IDW and OK methods in predicting the bathymetry of Saldanha Bay. The IDW and OK interpolation methods were selected because they are commonly used and are considered the most commonly favoured interpolation methods to estimate continuous surfaces (Curtarelli et al., 2015;Ferreira et al., 2017;Kalivas et al., 2013).

| Inverse distance weighting
IDW applies a deterministic and rather simple algorithm. It is described as a non-statistical interpolation method that ignores the spatial distribution of the data, and predicts unknown values by considering the proximity of known values (Li, 2013;Wang et al., 2014). IDW follows a nearest neighbour approach that gives influential weights to data points. Known data points carry more weight when in close proximity to unknown points and estimates are obtained as a weighted average of known neighbours. The general expression for IDW, as described by Wang et al. (2014, p. 3746), is: where Z is the estimated value of the interpolation point, Z i is the value of sample point i, n is the number of sample points, d i is the distance between the sample and interpolated points, and p is the power parameter (positive real number).
According to De Souza, Krueger, and Sluter (2003), cited in Ferreira et al. (2017), this algorithm has the advantage of creating smooth interpolated surfaces, and accounts for dimension parameters, number of sampling points and the power parameter, which controls the weighted neighbouring points on the interpolated points.
However, this method does not consider data trends. It also has the disadvantage that it prefers the closest neighbouring points where weighted averages are essentially similar for all points in close proximity (Ferreira et al., 2017). Another disadvantage relates to the number of sample points used. An increase in sample points will smooth out the interpolated surface and a lack of sample points will affect the accuracy of the results unfavourably and create what is known as the "bull's-eye" effect (Setianto & Triandini, 2013). The bull's-eye effect is an artefact of IDW. It causes concentric circles to appear around known measured points. This happens when there is a lack of sample points around known values, which results in the occurrence of peaks or channels (Li, Wang, Ma, & Wu, 2018). The value of the power parameter is another factor which adversely influences the accuracy of the results created through the IDW method. A lower power parameter will smooth out the interpolated surface. A high power value will place more weight on the nearest known points. The surface will therefore be less smooth (Esri, 2020a;GIS Geography, 2021). According to Pham, Van Huynh, Tran, and Chau (2016), the best choice of power parameter is 2.

| Kriging
In contrast to deterministic models, kriging is a geostatistical interpolation method that accounts for directional trends in data. Kriging is similar to deterministic methods in that it also applies linear regression to predict unknown values by considering weights assigned to known points. It accounts for the spatial autocorrelation between known values (Luo, Taylor, & Parker, 2008). Optimal interpolation weights are determined by a best-fit semivariogram model that considers both distance and direction to determine the spatial relationship of the data.
The semivariogram applies the following unbiased equation to compute both distance and direction of the data (Aziz, Yusof, Daud, Yusop, & Kasno, 2019, p. 161): where (h) denotes the semivariogram, h is the lag/distance, n (h) is the number of pairs of data values, Z is an intrinsic random function, n is the number of sample points, and the x i are data values at point i.
The ArcGIS Geostatistical Analyst offers various types of kriging methods to choose from, namely ordinary, simple, universal, indicator, probability, disjunctive, empirical Bayesian, and areal interpolation, each with its own intrinsic capacity to handle different data types (Esri, 2019). Ordinary kriging is the most widely used interpolation method and is also the default selection in ArcGIS software (Esri, 2020a;Kis, 2016). OK was therefore applied during this study. ArcGIS allows the user to choose among different functions (e.g., circular, spherical, exponential, Gaussian and linear) to model the empirical semivariogram when applying the Kriging tool. Each function influences the modelling of data differently, and the choice of the optimal model is dependent on the spatial autocorrelation and prior knowledge of the data to achieve a best-fit representation of the data. This normally entails a trial-and-error approach to find the parameters that optimise the model (GIS Geography, 2017). The typical semivariogram and its properties (Azavea, 2016;Esri, 2020a;Ferreira et al., 2017) are graphically illustrated in Figure 3. The range is the distance at which the semivariogram meets the variance of the data set. Points located within the area of the range are spatially correlated. The sill is the variance of the entire data set. Sample points that fall outside this point are considered not to be spatially dependent. The nugget is the error that exists between zero separation distance and the position where the semivariogram intercepts the y-axis. This random effect can be attributed to measurement errors in the data. Ideally, the nugget value should be y(0) = 0, but this seldom happens because of random sampling errors and short-scale reliability. Large errors in location measurements and sparsely distributed data contribute to a larger nugget effect. However, Camana and Deutsch (2020, p. 4) state that "a higher nugget effect in kriging will lead to smoother estimates" which are "realistic and correct".
The partial sill is the extent of variation of the sample points; it is derived by subtracting the nugget from the sill.
Naturally, a strong spatial correlation exists between the sample points. When compared to the nugget, the partial sill is higher (Aspexit, 2020). The semivariance is the spatial correlation between data points plotted on the y-axis.
Finally, the lag distance is the distance between pairs of sampling points plotted on the x-axis. Two rules of thumb for selecting the lag distance are highlighted by Esri (2018): first, to have at least 30-50 pairs of sampling points for every semivariogram point measured; and second, half the largest distance between all sampling points should be equivalent to the number of lags multiplied by the lag size.
For this study, the OK interpolation method was applied, with a Gaussian function to model the semivariogram. The mathematical formula normally used to express OK is as follows (Setianto & Triandini, 2013, p. 23): (2) where Z S 0 is an unsampled location, n is the number of sample points, the i are weights assigned to each observed sample point, Z is a random function, and the S i are the sample locations.
OK is considered the best model when data samples are similar at short distances (GIS Geography, 2016), a characteristic of the Saldanha Bay bathymetry data used in this study. OK is a simple geostatistical method that uses a linear-weighted technique to predict values from a stationary random field based on the assumption that the constant mean of the data set is unknown (Adhikary & Dash, 2017;Setianto & Triandini, 2013).

| Bathymetric data and interpolation tests
Measurements of water depth (bathymetry) are called soundings and are commonly used to map the bottom of oceans, rivers or lakes through interpolation to create high-resolution bathymetric maps. Various methods are available in mapping bathymetry, such as remote sensing platforms (aircraft, satellites and drones), ships and underwater vehicles, using techniques such as SoNAR and LiDAR (Giordano, Mattei, Parente, Peluso, & Santamaria, 2015).
A multi-beam echosounder is currently most commonly used to capture soundings. This SoNAR instrument is typically attached to the bottom or side of ships to record reflected sound pulses from the ocean floor to calculate bathymetry (Hasan, Ierodiaconou, Laurenson, & Schimel, 2014).
In this study, pre-processing included the conversion of soundings to metres above mean sea level with the transverse Mercator projection and WGS 1984 UTM Zone 34S coordinate system (Henrico & Bezuidenhout, 2020).
To simulate the use of a varied quantity of soundings to conduct the three comparison tests, the ArcGIS 10.5.1 Subset Features tool was used to extract sampling points randomly from the entire sounding data set. First 66% and then 33% of the entire data set was extracted, so that, besides the original sounding data set, two additional data sets were created to constitute the three different quantity data sets used to conduct the three interpolation comparison tests. Test 1 used the entire sounding data set (Figure 4a), which consisted of 1,653 points; test 2 used F I G U R E 3 Illustration of a typical semivariogram. Source: Adapted from Scheeres (2016) and Ferreira et al. (2017) 66% of the sounding data set (Figure 4b), which consisted of 1.091 points; and test 3 used 33% of the sounding data set (Figure 4c), which consisted of 546 points.
In this study the sounding data set was divided into three sets with varying quantity of sounding points to conduct the interpolation (comparison) tests to create bathymetric maps with a ground sampling distance of 2 m. For each test, both the IDW and OK methods were applied to determine its efficiency when the number of sampling points differed. The ArcGIS 10.5.1 Geostatistical Analyst wizard was used to conduct all interpolation tests. The parameters selected for both the IDW and OK methods were the same for all tests. They are described below.
The IDW interpolation method is based on the principle that points closer to each other are more alike than those further apart and therefore assigns greater weights to points closer to the interpolated location. For all IDW tests, the Optimize function available in the wizard was selected to optimise the output for the model (Table 1).
However, for the IDW test only the power value is optimised (i.e., between 2.19 and 2.9). This indicates that lower weights were assigned to points further from the interpolated location. Therefore, immediately surrounding points were more inclined to influence the predictions. The maximum and minimum neighbouring points to consider during these tests were 15 and 10 (default selection), respectively. The one-sector parameter (a single ellipse) was selected for the neighbourhood sector type.
OK bases the estimation of interpolated points on the assumption that the mean of the local neighbourhood is constant and unknown. The ArcGIS Geostatistical Analyst wizard follows a two-phase approach in conducting kriging tests. Firstly, a semivariogram is applied to analyse and optimally model the data. Secondly, the semivariogram model is used to create the interpolation map. For all OK tests the Optimize function was selected to optimise the output for the model. The properties of the semivariogram (nugget, lag size and number, partial sill and range) were optimised, and their values are shown in Table 2. The maximum and minimum neighbours were 5 and 2 (default selection) respectively, and the neighbourhood sector type was a four-sector one with 45° offset.
The Geostatistical Analyst wizard allows one to examine the spatial distribution and relationship of measured points to make changes to the parameters before creating the interpolation map. This happens after the OK method is selected and the semivariogram is created and displayed. Figure 5 illustrates the Gaussian semivariance models created for tests 1, 2 and 3.

| RE SULTS AND PERFORMAN CE OF THE INTERP OL ATI ON ME THODS
The results and performance of the IDW and OK comparison tests are highlighted in this section via the interpolation maps created by each method. The performance of each method is analysed by assessing the original F I G U R E 4 The different quantities of sounding points used during this study to conduct the three interpolation comparison tests. (a) Test 1 was conducted using the entire sounding data set, which consisted of 1,653 points. (b) Test 2 consisted of 1,091 points. (c) Test 3 consisted of 546 points measured and predicted values of the three interpolation tests. Consequently, test 1 used the entire data set of 1,653 inputs points for performance assessment, test 2 used 1,091 input points and test 3 used 546 input points.
Results are statistically quantified and comparisons are made between the difference and performance of IDW and OK. The interpolation maps ( Figure 6) are displayed with a red to blue colour map ramp which indicates deep to shallow waters. The performance of the interpolation methods was assessed through SPSS statistical software (version 27) (https://www.ibm.com/za-en/analy tics/spss-stati stics -software). Firstly, 2D scatterplots were created to show the relationship between predicted and measured sounding values for each of the interpolation methods. Secondly, an analysis of variance (ANOVA) test was conducted to determine whether the difference between the IDW and OK interpolation methods was statistically significant, and to deduce which method is most suitable to determine the bathymetry of Saldanha Bay. The following null and alternative hypotheses are defined:   Through visual comparison of these results, it is evident that there are explicit similarities in the way IDW and OK estimate and generate the continuous surface of bathymetry. It is also visually evident that the OK method produced a smoother representation of the interpolated surface compared to the rougher depth contours created by IDW, which may be attributed to the estimation approach that each interpolation method follows. The bull'seye effect (described in Section 3.1) is clearly visible in the IDW result of Figure 6e (indicated by black circles).
Even though this effect is also found on some of the OK results, it is more evident in the IDW results.

| Relationship between predicted and measured sounding values
Two-dimensional scatterplots (Figure 7) were used to evaluate the relationship between the predicted and measured sounding values for each of the interpolation methods through the use of a different quantity of sounding points.
The results (see Table 3) show that a strong linear relationship exists between predicted and measured sounding values for both IDW (r = 0.9560, p < 0.001) and OK (r = 0.9349, p < 0.001) when 100% (n = 1,653) of the sounding points are used. In Table 3, Spearman's correlation coefficient (r s ) is 0.95 (IDW) and 0.93 (OK). This is statistically significant (p < 0.001). However, this picture changes dramatically when 66% (n = 1,091) and 33% (n = 546) of the sounding values are used. Table 3 shows that during IDW interpolation r s is respectively 0. The results indicate that the relationship between the predicted and measured sounding values for each of the interpolation methods is highly influenced by sample size. This means that for a large quantity of sounding samples (n = 1,653) to determine the bathymetry of the Saldanha Bay, the strong correlation can be classified as significant.
Conversely, for medium (n = 1,091) and small (n = 546) quantities of sounding samples the weak correlation can be classified as not significant. It is also evident by evaluating the scatterplots that the OK prediction showed an underprediction for all three interpolation tests. Underprediction is a property of OK, as confirmed by other studies (Bradaï, Douaoui, Bettahar, & Yahiaoui, 2016;Esri, 2020b;Thiart, Ngwenya & Haines, 2014).

| Statistical analysis of most suitable interpolation method between IDW and OK
Next, an ANOVA test was conducted to determine if the means of the IDW and OK interpolation methods are significantly different (Liu & Wang, 2020;Singh, 2018). Bland and Altman plots were used to plot the measured versus the predicted values for IDW and OK and to show the relationship between these values to determine the variability of differences between the two interpolation methods. The Bland and Altman method is good for visually checking that "the approach was reasonable and that the data were 'well-balanced'" (Flegal, Graubard, & Ioannidis, 2020, p. 1312).
An analysis of all sounding data indicated that between the two interpolation methods, IDW performed better than OK in predicting Saldanha Bay depth values. is accepted.
The Bland and Altman plots for measured and predicted analysis of IDW and OK are shown in Figure 9. The majority of the values fall within the "limits of agreements", but the effectiveness of IDW is confirmed with a mean of 0.01 compared to a mean of 0.32 for OK. The Bland-Altman Limits of Agreement (BA LoA) are statistical limits, which are "calculated by using the mean and the standard deviation/s of the differences between two measurements" (Giavarina, 2015, p. 143).
TA B L E 3 Two-dimensional scatterplot statistics for the IDW and OK interpolation methods

| D ISCUSS I ON AND CON CLUS I ON S
The accuracy prediction of any interpolation method, including IDW and OK, to estimate bathymetry (of Saldanha Bay) is dependent on various factors of influence, including data variations (e.g., accuracy, density and distribution), type of terrain and spatial resolution of the output grids (Li & Heap, 2011 (2011), 53 comparative studies were assessed and a total of 72 interpolation methods-with IDW and OK among the most commonly used-were compared and analysed. Li and Heap (2011) stated that "there are no consistent findings about" the effects that the influencing factors have on the performance of the interpolation methods. It was also stated that sampling density was found to be insignificant and that it does not affect the performance of the methods.
In this study, the influencing factors include: highly accurate sounding points, which were collected using single-beam echosounders; varied density soundings (evident in the three interpolation tests that were conducted) and a good distribution of sample points, covering the entire study area. The terrain researched was the ocean floor of Saldanha Bay. Its terrain type is characterised by a slope profile ranging from 0 m (coastline) to approximately 50 m (where the Bay meets the open Atlantic Ocean). The bathymetric maps created had a very high spatial resolution of 2 m which generated smooth bathymetric maps. In terms of factors that could affect the accuracy prediction of the IDW and OK interpolation methods to estimate the bathymetry of Saldanha Bay, all factors had an equal effect on both IDW and OK to produce optimal prediction results, because the same software tools, processes, settings and data sets were utilised to conduct each comparative test.
The aim of this study was to analyse and compare the efficiency of the IDW and OK interpolation methods to predict the bathymetry of Saldanha Bay through a varying quantity of sampling points. Visual comparison reveals definite similarities in the way IDW and OK estimate and generate the continuous surface of bathymetry.
The results show that using more sample points produced better interpolation results, which is confirmed by the prediction errors that increased as the quantity of soundings decreased. IDW consistently performed better than OK across all interpolation tests. However, only a large quantity of sounding samples (n = 1,653) to determine the bathymetry of Saldanha Bay showed a strong correlation that was classified as statistically significant. Based on the results obtained from this study, IDW is identified as the most suitable interpolation method to predict the bathymetry of Saldanha Bay. Compared to OK, it produced better accuracy predictions for all tests conducted.
These findings were confirmed by the ANOVA statistical method.
Although this study evaluated two well-known and common interpolation methods using ArcGIS software, some limitations do exist. It might be feasible for future studies to determine the effectiveness of other interpolation methods to predict the bathymetry of Saldanha Bay. Scientists agree that no superior or preferred method exists and that results are largely site-specific and data-specific (Murphy et al., 2010). It might also be interesting to perform similar interpolation tests with the use of different GIS software, such as QGIS. Such a study, to compare the software functionalities and interpolation tools between open-source and proprietary software, would be invaluable given the current global GIS market trend to move towards and adopt open-source solutions (Maida, 2020).
The effectiveness of IDW as opposed to OK in predicting the bathymetry of Saldanha Bay was illustrated by this study. Although the Optimize function available in the ArcGIS 10.5.1 Geostatistical Analyst wizard was applied to conduct the interpolation tests, and default parameters were applied, it is the opinion of the author that altering the parameters and settings of either of the models would not have had a significant influence on the outcome of these tests. The approach followed during this study established that IDW consistently outperformed OK when sampling quantities and patterns varied. This identifies IDW as the preferred interpolation method to use for future analysis of the bathymetry of Saldanha Bay.

ACK N OWLED G EM ENTS
The author would like to thank the South African National Hydrographer for supplying the bathymetric data of Saldanha Bay, Professor Martin Kidd for the statistical consultation provided, and Dr Gerhard Van Zyl (Stellenbosch University) for providing language help and proofreading the manuscript. The insightful comments of reviewers helped to significantly improve the manuscript.