Cartographic Interpretation of the Seafloor
Geomorphology Using GMT: a Case Study of the
Manila Trench, South China Sea

Polina Lemenkova

► To cite this version:

Polina Lemenkova. Cartographic Interpretation of the Seafloor Geomorphology Using GMT: a Case Study of the Manila Trench, South China Sea. Aksaray University Journal of Science and Engineering, 2020, 4 (1), pp.1 - 18. 10.29002/asujse.604761. hal-02873799

HAL Id: hal-02873799

https://hal.science/hal-02873799

Submitted on 18 Jun 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Cartographic Interpretation of the Seafloor Geomorphology Using GMT: a Case Study of the Manila Trench, South China Sea

Polina Lemenkova*

College of Marine Geo-sciences, Ocean University of China, Qingdao, Shandong Province 266100, People's Republic of China

*Received Date: Aug 09, 2019

•Revised Date: Jan 25, 2020

•Accepted Date: Feb 28, 2020

•Published Online: Jun 18, 2020

## Abstract

The study is geographically focused on the Manila Trench, located in the west Pacific Ocean, South China Sea, west Philippines. The research aims at the geological mapping, analysis and visualizing variations in the submarine geomorphology of the Manila Trench. Technically, the work was done using Generic Mapping Tools scripting toolset (GMT). A combination of various GMT modules was applied for geospatial modelling. Methodology includes cartographic data integration and interpretation through approaches of data analysis: topographic plotting, geophysical modelling, geological mapping and statistical analysis. The data included SRTM, ETOPO1, geoid and gravity grids (CryoSat-2, Jason-1). Two sets of the cross-section profiles of the trench were automatically digitized. The profile transects were compared and differences in the geomorphic shape in southern and northern parts revealed. Southern part has steeper slope on the western part. Northern part is steeper on the continental slope part. The submarine terraces are located on the northern segment at depths -2,000 m. The depth and geomorphology of the slope vary for the range -3,500 to -4,500 m: minimal for the northern part with 526 samples (18.2%) for the depths -4,000 to -4,200 m. The histogram for the northern part has bimodal distribution with two peaks. The southern part shows 142 values for the minimal -3,500 to -3400 m. The statistical analyses revealed that northern part of the trench is deeper. The GMT functionality shown in this paper enabled integration and interpretation of the multi-source data: automatically digitized profiles, geological mapping, 2D and 3D bathymetric modelling, statistical analysis, mathematical approximation of the trend modelling. The GMT proved to be capable of visualizing geodata that can significantly improve Earth studies and interpretation of submarine geomorphology of the oceanic trenches through the advanced data analysis.

## Keywords

Mapping, Cartography, GMT, Manila Trench, Pacific Ocean, South China Sea

<figure><img src="image_2.png" /></figure>

<footnote>*Corresponding Author: Polina Lemenkova, 238 Songling Road, Laoshan, 266100, Qingdao, Shandong Province, People's Republic of China. Tel.: +86-1768-554-1605. E-mail: pauline.lemenkova@gmail.com</footnote>
# 1. INTRODUCTION

Understanding a variety of the submarine landforms hidden from the human eyes due to their remote location, and factors affecting seafloor geomorphology, is a very complex task that requires a multi-disciplinary approach: geological data analysis, data processing and modelling by advanced algorithms, geostatistical analysis and visualization. The asymmetry of the deep-sea trenches reflects a phenomenon of tectonic plates subduction. Thus, as one plate bends down to the Earth's mantle, another plate is being deformed filling the growing empty space. The depths of trenches are influenced by many processes and factors controlling their actual shape and structure. Many attempts were undertaken to map and visualize general form of the submarine geomorphology and model oceanic deep-sea trenches. Rapid development of the IT technologies in XXI century facilitated geodata processing which contributes to our better understanding of the marine geomorphology.

The significance of the research presented in this paper is associated with the presented approach of GMT as an advanced toolset among geoinformation GIS technologies applied for geological mapping and geomorphic analysis. Having access to the machine learning technologies associated with rapid recent IT progress (scripting coding, geospatial datasets and DEMs), data analysis offered by GMT can be employed for marine geological research to provide new, more detailed insight into the seafloor bathymetry. Presented GMT codes enable to perform a speed yet quality mapmaking and data processing. In turn, this allows to study in details submarine geomorphic landforms of the oceanic seafloor. In turn, this allows an in-depth study of the Manila Trench.

The contribution of this research towards development of the methodologies of seafloor mapping and modelling geodata by GMT consists in the presented and explained GMT scripting codes and snippet. Rather than traditional GIS interface, a GMT is fully based on the scripting approach where cartographic methodology consists in executing the programming code which results in mapping. Geoinformation is crucial for understanding ocean seafloor and precise bathymetric mapping. Combining datasets of the raster maps with geomorphic modelling increases our knowledge of the factors affecting seafloor landforms.

The research aim is comparative geomorphic analysis of deep-sea oceanic trench: Manila Trench located in the geologically complex region of the western part of the Pacific Ocean, a part of the 'Ring of Fire' where active movements of the tectonic plates and volcanism take place. Complex geophysical settings affect the formation of the trench, high seismicity and geodynamic instability visualized on the thematic maps of geological and tectonic settings as the most important causes of the ocean trench formation. Cartographic objective of the research was to visualize bathymetry, geomorphology, tectonic and geological settings of the trenches through 2D and 3D modelling and mapping. The statistical analysis aimed to compare datasets of the trench and show differences in geomorphic and bathymetric structure.
## 2. STUDY AREA AND MATERIALS

The study geographically focuses on the Manila Trench, an oceanic trench located in the South China Sea, west Pacific Ocean, west of the islands of Luzon and Mindoro in the Philippines (Fig. 1).

The study presents geomorphic mapping and the use of data models of various spatial resolution to visualize certain geophysical behaviors of the crust in the region. The area of the Manila Trench is located in the Manila subduction zone at the Philippine Sea plate boundary where it moves in the northwest direction toward the Eurasian plate with a high convergence rate [1].

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Bathymetric map of the study area based on the SRTM 15 sec grid. Image source: author, made using GMT. Applied color palette: 'geo.cpt' by GMT, with colors adjusted for global bathymetry/topography relief [R=-8000/8000, H=0, C=RGB]. SRTM15_PLUS is a unique dataset, introduced in 2019 from the SRTM mission. Publicly available SRTM15_PLUS dataset was used in the current research. Hence, the base map for this study is a high resolution (15 arc second) raster grid: a data fusion, SRTM15_PLUS combining SRTM and other models for land and other datasets for bathymetry. Image source: author, made using GMT.</figcaption></figure>

The Manila Trench with its maximal depth of 5,400 m [2] stretches in almost vertical North-South direction. It is created by the subduction of the Eurasian Plate (through its part Sunda Plate) under the Philippine Sea Plate. Specific geological structure of the Manila Trench causes potential repetitive
earthquakes [3]. Another frequent hazard caused by the subducting slab of the Eurasian plate beneath the Manila Trench (Fig. 2) is destructive tsunamis causing catastrophic damages in the source area and along the coastlines of the Luzon Island.

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Geologic map of the Manila Trench based on the ETOPO1 1 min grid. Geologic map of the PSB. Solid red circles indicate active volcanoes in the study region. Yellow fronts denote trench axes (Manila Trench). Tectonic slab contours are depicted by the orange lines. Color palette is visualized by 'globe.cpt' embedded in GMT with following specifications: Colors for global bathymetry/topography relief [R=-10000/10000, H=0, C=RGB]. Image source: author, made using GMT.</figcaption></figure>

The hypocenters of the tsunami of the Manila Trench are located at the depths <100 km [4]. The tsunami hazard from the Manila Trench source has been assessed in more details in several research papers [5-7]. The zone of the Eurasian Plate subduction explains the belt of volcanoes in the Manila Trench area (Fig. 2), on the west side of the Philippine island of Luzon. The area between the northernmost Manila subduction zone and southern Taiwan is considered a
transition from the subduction to initial collision, and a weakly coupled nature for the northern part of the Manila Trench [8].

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Model of the geoid along the Manila Trench. Density variations of the Earth's interior cause geoid undulations. Because the structure and submarine geomorphology of the ocean trench is affected by Earth's interior through tectonic plates subduction and movements, the visual shape of the geoid is distinctly repeating the isolines, as shown on image above. Image source: author, made using GMT.</figcaption></figure>

Submarine regions, in contrast with the terrestrial ones, still remain the least accessible areas on the Earth due to their remote location. Therefore, modelling oceanic trenches is one of the most complicated issues in the marine geology explained by their inaccessibility. Oceanic trenches can only be studied using computer based modelling, advanced mapping and algorithms of the data analysis.

Therefore, using special software and tools designed for geospatial data modelling enables to get an insight into the deepest regions of the ocean and to visualize areas with the most difficult access [9].
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Model of the gravimetry along the Manila Trench. The accuracy of the marine gravity measurements is important for reliable results and refines approximation of the models with Earth ellipsoid. Gradient changes in gravity anomalies have visual correlations with geoid as caused by the indirect effects of the geoid undulations. Hence, visualized image of the gravity is reflected by the structure and submarine geomorphology of the trench affected by plates subduction. Image source: author, made using GMT.</figcaption></figure>

The morphology of the seafloor is caused by the uneven distribution of the elevations on the Earth with distinctly uneven hypsography. Thus, the majority of the depths on the Earth is occupied by the deep basins (ca. 4–6.5 km), while relatively few areas are covered by the shallow zones [10].

The tectonic plate boundaries are the hotspot areas where the most of the largest earthquakes take place and oceanic slabs descend beneath the continental lithosphere causing trench migration [11]. The largest earthquakes mostly occur in the shallow part of the subduction zones [12]. The earthquakes are sometimes accompanied by strong tsunami waves. Slab dynamics is one of the most important drivers for the trench formation, dynamics and migration [13]. Many factors affect speed and direction of the trench migration. These include plate-mantle coupling, slab interactions with the mantle transition, plate geometry, kinematics, strain rate, temperature,
fluid pressure, deformation mechanisms, as well as mantle rheology [14]. The geological complexity of the Manila Trench is also expressed by the connection with the Philippine Trench where a horizontal mantle flow exists between the Manila Trench and the Philippine Trench. It is caused by the collision between the Palawan block and the Philippine Mobile Belt, and movement of the South China Sea slab [15]. Various factors cause the formation of the Manila Trench.

<figure><img src="image_7.png" /><figcaption>**Figure 5.** 3D map of the Manila Trench visualizing geomorphology. A detailed 3D mesh plot on top of the 2D geoid contour image as a block diagram was created as a combination of bathymetric and geoid image. 3D bathymetry model shows selected segment of the Manila Trench to the west of the trench. The margin is characterized by steep slope adjacent close to the Luzon Island. Image source: author, made using GMT.</figcaption></figure>

Being an oceanic trench, it presents a special area of the ocean seafloor with distinct geomorphological structures characterized by the notable depths and steep gradient angles, located in the zones of the continental margin tectonic plates bending [16]. The convergence between the two plates forming the Manila Trench is roughly northwestward. There is a high-
velocity zone present in the crust and upper mantle beneath the Luzon arc where the Manila Trench is stretching [17].

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Transect of the cross-section profiles digitized along the central and southern regions of the Manila Trench. Cross-sectional transects of the Manila Trench were taken in the central and southern regions of the trench, through key features such as deep valley, slopes, moraines, South China Sea basin, and adjacent shelf area. The selected cpt color palette include 'terra.cpt' colors for global bathymetry/topography relief [R=-7000/7000, H=0, C=RGB]. Image source: author, made using GMT.</figcaption></figure>

All these factors briefly mentioned to illustrate the tectonic settings of the study area indicate high frequency of earthquakes and special tectonic zone where the Manila Trench has been formed. The seafloor of the deep-sea trench presents a complex structure combined of various landforms: mid-ocean ridges, transform faults, ocean plains complicated by chains of
seamounts, minor ridges, trenches and plateaus [18]. Factors influencing geomorphic
development, structure and bathymetric patterns of the trench are diverse. To mention some of
them: geological, hydro-chemical, geothermal, climatic, tectonic and bathymetric determinants
[19]. The steepness of the Manila geomorphology varies: the dip angles along the Manila
Trench increasing gradually southwards from ca. 25° at the latitudes between 18° and 21° N,
32° at the latitude 17° N, and nearly vertical at the latitude 14° N [20]. The tectonic front of the
Manila Trench continuing northward into the frontal thrust faults in the western Taiwan shows
active plate boundary between the Eurasian and Philippine Sea plates in the Luzon-Taiwan
region [21].

<figure><img src="image_9.png" /><figcaption>Histograms of the bathymetry, Manila Trench: northern (A) and southern (B) parts</figcaption></figure>

**Figure 7.** Histograms of the bathymetry of the two selected segments of the Manila Trench. The statistical histograms showing frequency of data distribution for the transecting profiles of the Manila Trench bordering the South China Sea basin are visualized on the figure above. The shape of the histogram varies by the two segments of the trench reflecting their geomorphology. Image source: author, made using GMT.

The structure of the northern Manila Trench has been studied in various papers focused on the
problems of crustal structure and deformation in the north of the Manila Trench [22-24]. Among
other findings, the increasing dip angle of the Manila Trench from north to south has been discovered. The Manila trench loses its surface roughness in the early collision zone and gradually becomes a less well-defined deformation front [25]. However, the comparative analysis of its southern and northern parts is still missing. Therefore, current paper contributes to better understanding of the Manila Trench by comparative analysis of its geomorphology in the northern and southern segments. The data used in the current study included SRTM, ETOPO1, geoid and gravity grids (CryoSat-2, Jason-1). Topographic data used in this research include ETOPO1 Global Relief Model grid containing elevation data (topography/bathymetry) for the Earth.

## 3. METHODOLOGY

Data analysis applied for marine geology may be supported by programming languages, e.g. R [26], Octave, AWK [27], statistical libraries, e.g. Gretl [28] and other advanced tools, such as statistical software, e.g. SPSS [29]. For this study, a Generic Mapping Tool (GMT) scripting toolset [30] was selected as a main tool. It enabled to perform all the steps of the research: data analysis, modelling, mapping and statistical analysis and visualization. The functionality of the GMT and its powerful cartographic possibilities explain the choice for the GMT. The use of the GMT specifically for the Manila Trench was presented as mapping in various works, e.g. [31]. Examples of other advanced tools for geological data analysis include PHASER diffractometer, direct observations received from the Research Vessel (R/V) cruises, observations from Digital Elevation Models (DEM) [32]. The mapping was based on the available data sets embedded in the GMT [33] as well as Shuttle Radar Topographic Mission (SRTM) data [34-35]. The geologic characteristics was performed by modelling seafloor geomorphology, mapping tectonic slabs and submarine volcanoes.

The geoid gravitational regional modeling was done through the GMT modules 'psbasemap', 'grdcontour' and 'grdimage'. Plotting contours of the terrestrial and water areas, geoid and bathymetry, net grids and basic cartographic elements such as titles, scales, annotations, was performed for the mapping (Fig. 3). Each GMT module consists of the set of code lines, which, taken together as a script, produce and visualize maps. The main code line to generate geoid image with shading was as follows (Code 1):

```
gmt grdimage geoid.egm96.grd -l+a45+nt1 -R105/123/8/24 -JD114/13/16/19/6i -Cgeoid.cpt -P -K > $ps
```

The -R command shows the WESN coordinates of the main area: 105/123/8/24.
The gravity mapping shows free-air gravity anomalies covering study area including terrestrial and marine regions. Modelling was based on the data sets based on the CryoSat-2 and Jason-1 satellite-derived data [36], for the Manila Trench and surrounding regions (South China Sea, the Philippines). The gravity modelling revealed ground information for the geologic and geophysical settings of the study area, as free-air gravity anomalies reflect subsurface density variations [37]. The gravity dataset covering research area is based on the map gridded with a resolution of one arc-minute (Fig. 4). The gravity map was plotted using following GMT code (Code 2):

```
gmt grdimage gravRT.grd -I+a45+nt1 -R105/123/8/24 -JY114/12/6.5i -CgravRT.cpt -P -K > $ps
```

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Transect of the cross-section profiles digitized along the Manila Trench. Calculated statistical trends were compared for both segments of the trench to show mathematic approximations of the curvature trends of their bathymetric shapes of the cross-sections. The statistical values of the geomorphic elevations show variations in the both segments. It enables to reveal general trends and helps to highlight correlations in the data distribution. Image source: author, made using GMT.</figcaption></figure>

The GMT modules used for plotting 3D model (Fig. 5) of the mesh grid showing geomorphological shape of the study area. Specifically, various GMT modules have been used: grdcut, grd2cpt, grdcontour, pscoast, grdview, logo, psconvert. The mesh gridding was selected
to visualize and highlight topographic and bathymetric landforms of the region, through
dividing the input grid model into smaller shapes formed after discretization of the geometric
model of the area. The main GMT code for modelling include the following code (Code 3):

```
grdview mnt_relief.nc -J -R -JZ3.0c -CrainbowRT.cpt -p215/30 -Qsm -N-7500+glightgray -Wm0.07p -Wf0.1p,red -Wc0.1p,magenta -B2/2/4000:"Bathymetry and topography (m)":eSWZ -S5 -Y6.5c -O -K >> $ps
```

Using this code, the 2D initial meshing was divided into small simple polygons with curved contour lines and color palette highlighting the topography.

The geomorphological modelling consisted in plotting two sets of the cross-section profiles and then comparing their parameters (Fig. 6). The raster image as a background (Fig. 6C) was mapped using following GMT code (Code 4):

```
gmt grdimage mnt_relief.nc -CmyoceanMnT.cpt -R105/123/8/24 -JM6i -P -I+a15+ne0.75 -Xc -K >> $ps
```

The cross-track profiles were then generated with specific parameters. These are 400 km length, 20 km distance, 2km sampling. The plotting and stacking these was done using the statistical mean, through the following code (Code 5):

```
gmt grdtrack trenchMTn.txt -Gmnt_relief.nc -C400k/2k/20k+v -Sa+sstackMTn.txt > tableMTn.txt
```

Two segments were plotted with the following coordinates (Code 6):

```
cat << EOF > trenchMTn.txt 119.2 16.4 119.2 14.0 EOF
cat << EOF > trenchMTs.txt 119.2 14.0 120.0 13.2 EOF
```

The graphs (Fig. 6A and Fig. 6B) were plotted using following GMT code (Code 7):

```
gmt psxy -R -J -W1.0p -Ey+p0.2p stackMTn.txt -O -K >> $ps
gmt psxy -R -J -W1.0p,red stackMTn.txt -O -K >> $ps
```

Minor GMT modules used for plotting cartographic elements included 'psxy', 'logo' 'psconvert', as well as Unix 'echo' utility.

The statistical modelling of the digitized cross-section profiles aimed at comparison of two sets in southern and northern segments, respectively, were done using module 'pshistogram' by following GMT code snippet (Code 8):

```
gmt pshistogram tableMTn.txt -i4 -R-5500/1500/0/20 -JX4.8i/2.4i -Y7.5c -Bpxg1000a1000f100+l"Bathymetry (m)" -Bpyg5a5f2.5+l"Frequency"+u" %" -Bsyg2.5 -BWSne+t"Histograms of the bathymetry, Manila Trench: northern (A) and southern (B) parts"+gsnow1 -Glightsteelblue1 -D+f7p,Times-Roman,black -L0.1p,dimgray -Z1 -W250 -N0+pred -N1+pblue -N2+pgreen -O -K >> $ps
```
## 4. RESULTS AND DISCUSSION

Based on the patterns of the free-air gravity map (Fig. 4), the high density materials are distributed continuously from the northern to the southern Taiwan vertically southwards to the Luzon Island, the Philippines from -50 mGal to -200 mGal. The majority of the values in the South China Sea areas lie between the 10 to 30 mGal, while values reaching >60 mGal can be noted in the island areas and southern Philippines. The significance on the free-air gravity differences for the two profiles consists in the principle of the free-air gravity map, which contains signals from the seafloor topography, sediments, and crust and mantle density anomalies. More specifically, the free-air anomaly is dominated by short wavelength variations which reflect the density contrast at the seafloor. Therefore, the differences highlight the unevenness in the geophysical properties in the study areas. Thus, the area crossed by the southern profile (Fig. 6, green transect) has dominant values of -50 to -220 mGal along the profile while northern profile (Fig. 6, yellow transect) has more values of -20 to 20 mGal in the surrounding coastal area. Asymmetric profile shapes and more gentle slope in the northern profile (Fig. 6) can point at local geophysical variations, e.g. pulsing, dehydrating, and radially flowing mantle plume. Furthermore, bathymetry and gravity data demonstrate that V-shaped profiles in the deepest parts of the cross-section (-4,000 m) may be the result of sublithospheric flow. It is formed because of the displacement of large volumes of magma at the trench axis which is reflected in its geomorphology as east-west asymmetry (Fig. 6).

The statistical results (Fig. 7) show differences in the depth ranges of the Manila Trench along the Luzon island. The depth and morphological properties of the slope vary significantly for the range -3,500 to -4,500 m. Maximal values are notable for the northern part of the trench where 526 observation points (18.2%) are recorded for the depths -4,000 to -4,200 m, as shown on Fig. 7A. The histogram for the northern part has a clear bimodal distribution with two peaks. The 2nd peak shows areas crossing the Luzon Island. On the contrary, southern part demonstrates 142 values for the minimal bathymetry with values -3,500 to -3400 m. The 1st peak also shows areas of the Luzon Island. The comparison of the statistical data analysis clearly shows that the northern part of the trench is deeper.

The southern part of the trench has steeper slope from the oceanward part and on the contrary, the northern part is steeper in the continental slope part (Figs. 8A and 8F). The submarine terraces can be noted on the northern segment of the trench (Fig. 8F) at depths -2,000 m. Generally, mapping northern part shows steeper degree of the oceanward slope (pairwise comparing on Fig. 8: subplots C and H, D and I, E and J).
However, a gradually slanting terrace on the southern part of the trench can be seen on the fragment 50-100 of the cross-section with depths from -1,100 to -1,300 m. Hence, the comparison between two segments of the trench located on the northern and south-eastern part of the trench show that northern part of the trench has shallower values of depths comparing to the southern yet more steep slope on the Luzon Island side.

The approaches to the data analysis, visualization and geological modelling are diverse. To mention some of them: tomographic, seismic and bathymetric 2D and 3D modelling, geologic cross-sectioning, seismic cross-section profiling, R programming language for analysis of correlation between variables [38, 39], ArcGIS based assessment and calculation of measured data [40-42], ILWIS GIS [43], Python statistical libraries, such as Matplotlib, NumPy, SciPy, Seaborn, Pandas, StatsModels [44, 45]. Using special software for automatization of the geodata can be illustrated by Autotrace [46], R packages [47], [48] or geological software [49, 50]. In view of a variety of the approaches and methods illustrated above, the advantage of the GMT consists in its embedded statistical module that enable to visualize data and to perform descriptive statistical visualization and modelling as histograms showing data distribution.

## 5. CONCLUSION

The presented research sets out to demonstrate the utility of the Generic Mapping Tool (GMT) for cartographic visualization in support of interpreting geophysical data of the seafloor geomorphology of the Manila Trench. This investigation of the tectonic dynamics of the Manila Trench is an important problem in the field of geodynamics and relevant for understanding the potential origins of hazards related to earthquakes and tsunamis in the region. The study produces useful visual cues to interpret the tectonics of the crustal structure based on the derived data models such as the geoid, free-air gravity anomalies, bathymetry, and geology aimed at modelling and visualizing tectonic lithospheric settings and geological situation in the study area.

The GMT cross-section stacking methodology [51-53] proved to be a successful means of visualizing and plotting geomorphological models applied for the submarine bathymetry by effectively minimizing the hand-made cartographic routine. Using GMT enables to automatically digitize profiles, perform statistical analysis by plotting histograms and model trends for the general shape of the profiles using mathematical models of the lines approximation. Several GMT modules were tested and applied for the cartographic visualization of the marine free-air gravity, geoid, bathymetric mapping, geomorphic modelling
and statistical data analysis. Two different trench segments were analyzed aimed at detecting their spatial changes for an area covering Manila Trench and adjusting Luzon Island (northern Philippines). The results provide useful recipes on GMT scripts that produce visualization outputs for interpretation and explanation of the tectonic dynamics of the Manila Trench. The application of GMT scripts, produced in this paper, can easily be adapted to other investigations using similar datasets.

This paper also aims to close the knowledge gap on the geomorphological differences between the northern and southern part of the Manila Trench. Therefore, the goal of the GMT-based mapping was to relate geomorphological mapping to underlying tectonics conditions. These include the collision of two tectonic plates which causes instability in the region: repetitive earthquakes, submarine volcanism and high frequency of earthquakes in this region. Spatial variation in the south-east directed second segment (plotted as green line on Fig. 6C) shown steeper slopes on the oceanward side of the trench and deeper bathymetric records. On the contrary, northern segment (plotted as almost vertical yellow-colored line on Fig. 6C) demonstrated softer geomorphic shapes yet steeper gradient from the Luzon Island side. The fine-resolution data such as SRTM grid was used to map bathymetry and topography of the study area and 1-min resolution from the gravity data (Jason-1, CryoSat).

Future work should consider applying other GMT modules and their combination. Developing a complex approaches of the cartographic-geophysical investigations include earthquake intensity modelling, mapping detailed volcanic arcs. The combination and overlapping of various geospatial data as layers aimed to examining and quantifying geomorphological submarine landforms in connection with tectonic settings would enable to get more information of these data using advanced geospatial data analysis. Other recommendation for the future research could be integration of the GMT with programming and advanced statistical packages that absent in GMT: clustering, dendrogram plotting, factor analysis, Principal Component Analysis, etc.

## Acknowledgement

This research was funded by the China Scholarship Council (CSC), State Oceanic Administration (SOA), Marine Scholarship of China, Grant Nr. 2016SOA002, Beijing, People's Republic of China. I cordially thank two anonymous reviewers for their critics and comments that helped to improve the original manuscript.
References

[1] C.-L. Lo, W.-B. Doo, K.-C. Hao, S.-K. Hsu, *Plate coupling across the northern Manila subduction zone deduced from mantle lithosphere buoyancy*. Physics of the Earth and Planetary Interiors, 273 (2017) 50-54.

[2] Y. Liu, A. Santos, S.M. Wang, Y. Shi, H. Liu, D.A. Yuen, *Tsunami hazards along Chinese coast from potential earthquakes in South China Sea*. Physics of the Earth and Planetary Interiors 163 (1-4) (2007) 233-244.

[3] T.-R. Wu, H.-C. Huang, *Modeling tsunami hazards from Manila trench to Taiwan*, Journal of Asian Earth Sciences 36 (2009) 21-28.

[4] A. Ruangrassamee, N. Saelem, *Effect of Tsunamis generated in the Manila Trench on the Gulf of Thailand*, Journal of Asian Earth Sciences 36 (2009) 56-66.

[5] P.H. Nguyen, Q.C. Bui, P.H. Vu, T.T. Pham, *Scenario-based tsunami hazard assessment for the coast of Vietnam from the Manila Trench source*, Physics of the Earth and Planetary Interiors, 236 (2014) 95-108.

[6] C.B. Bautista, M.L.P. Bautista, K. Oike, F.T. Wu, R.S. Punongbayan, *A new insight on the geometry of subducting slabs in northern Luzon, Philippines*. Tectonophysics 339 (2001) 279-310.

[7] S.L. Soloviev, C.N. Go, *A Catalogue of Tsunamis on the Western Shore of the Pacific Ocean*. Academy of Sciences of the USSR. Nauka Publishing, Moscow, (1974) 439.

[8] S.-J. Chin, J.-Y. Lin, Y.-C. Yeh, K.-C. Hao, C.-W. Liang, *Seismotectonic characteristics of the Taiwan collision-Manila subduction T transition: The effect of pre-existing structures*. Journal of Asian Earth Sciences 173 (2019) 113-120.

[9] B.R. Calder, L.A. Mayer, *Automatic processing of high-rate, high-density multibeam echosounder data*. Geochemistry Geophysics Geosystems, 4(6) (2003) 1-24.

[10] P. Lemenkova, *Testing Linear Regressions by StatsModel Library of Python for Oceanological Data Interpretation*. Aquatic Sciences and Engineering 34, (2019) 51-60.

[11] L. Zhang, C. He, Y. Liu, J. Lin, *Frictional properties of the South China Sea oceanic basalt and implications for strength of the Manila subduction seismogenic zone*. Marine Geology 394, 16-29 (2017).

[12] H. Yu, Y. Liu, H. Yang, J. Ning, *Modeling earthquake sequences along the Manila subduction zone: Effects of three-dimensional fault geometry*. Tectonophysics 733 (2018) 73-84.

[13] P. Lemenkova, *R scripting libraries for comparative analysis of the correlation methods to identify factors affecting Mariana Trench formation*. Journal of Marine Technology and Environment 2 (2018) 35-42.

[14] C. Faccenna, A.F. Holt, T.W. Becker, S. Lallemand, L.H. Royden, *Dynamics of the Ryukyu/Izu-Bonin-Marianas double subduction system*. Tectonophysics 746 (2018) 229-238.

[15] J. Fan, D. Zhao, P-wave anisotropic tomography of the central and southern Philippines, Physics of the Earth and Planetary Interiors 286 (2019) 154-164.

[16] P. Lemenkova, *Factor Analysis by R Programming to Assess Variability Among Environmental Determinants of the Mariana Trench*. Turkish Journal of Maritime and Marine Sciences, 4 (2018) 146-155.

[17] W.-B. Doo, K.-C. Hao, D. Brown, C.-L. Lo, S.-K. Hsu, Y.-S. Huang, *Serpentinization of the fore-arc mantle along the Taiwan arc-continent collision of the northern Manila subduction zone inferred from gravity modeling*. Tectonophysics 691 (2016) 282-289.

[18] P. Lemenkova, *Hierarchical Cluster Analysis by R language for Pattern Recognition in the Bathymetric Data Frame: a Case Study of the Mariana Trench, Pacific Ocean*. Virtual Simulation, Prototyping and Industrial Design. 2(5) (2018) 147-152.
[19] P. Lemenkova, *Scatterplot Matrices of the Geomorphic Structure of the Mariana Trench at Four Tectonic Plates (Pacific, Philippine, Mariana and Caroline): a Geostatistical Analysis by R.* Problems of Tectonics of Continents and Oceans 1 (2019) 347-352.

[20] Z. Cheng, W. Ding, M. Faccenda, J. Li, X. Lin, L. Ma, P. Fang, H. Ding, Geodynamic effects of subducted seamount at the Manila Trench: Insights from numerical modeling, Tectonophysics 764 (2019) 46-61.

[21] H.-S. Yu, Nature and distribution of the deformation front in the Luzon Arc-Chinese continental margin collision zone at Taiwan, Marine Geophysical Researches 25 (2004) 109-122.

[22] F. Li, Z. Sun, D. Hu, Z. Wang, *Crustal structure and deformation associated with seamount subduction at the north Manila Trench represented by analog and gravity modeling.* Marine Geophysical Research 34 (2013) 393.

[23] E. He, M. Zhao, X. Qiu, J.-C. Sibuet, J. Wang, J. Zhang, *Crustal structure across the post-spreading magmatic ridge of the East Sub-basin in the South China Sea: Tectonic significance.* Journal of Asian Earth Sciences 121 (2016) 139-152.

[24] C.-Y. Ku, S.-K. Hsu, *Crustal structure and deformation at the northern Manila Trench between Taiwan and Luzon islands.* Tectonophysics 466 (3-4) (2009) 229-240.

[25] R. Lester, K. McIntosh, H.J.A. Van Avendonk, L. Lavier, C.-S. Liu, T.K. Wang, *Crustal accretion in the Manila trench accretionary wedge at the transition from subduction to mountain-building in Taiwan.* Earth and Planetary Science Letters, 375 (2013) 430-440.

[26] P. Lemenkova, *An Empirical Study of R Applications for Data Analysis in Marine Geology.* Marine Science and Technology Bulletin, 8(1) (2019) 1-9.

[27] P. Lemenkova, *AWK and GNU Octave Programming Languages Integrated with Generic Mapping Tools for Geomorphological Analysis.* GeoScience Engineering 65(4) (2019) 1-22.

[28] P. Lemenkova, *Regression Models by Gretl and R Statistical Packages for Data Analysis in Marine Geology.* International Journal of Environmental Trends 3(1) (2019) 39-59.

[29] P. Lemenkova, *Numerical Data Modelling and Classification in Marine Geology by the SPSS Statistics.* International Journal of Engineering Technologies 5(2) (2019) 90-99.

[30] P. Wessel, W.H.F. Smith, *New, improved version of the generic mapping tools released.* Eos Transactions American Geophysical Union 79 (1998) 579.

[31] Y.-J. Hsu, S.-B. Yu, A.T.R. Song, T. Bacolcol, *Plate coupling along the Manila subduction zone between Taiwan and northern Luzon.* Journal of Asian Earth Sciences 51 (2012) 98-108.

[32] S.J. Giletycz, A.T.S. Lin, C.-P. Chang, J. Shyu, *Relicts of mud diapirism of the emerged wedge-top as an indicator of gas hydrates destabilization in the Manila accretionary prism in southern Taiwan (Hengchun Peninsula).* Geomorphology 336 (2019) 1-17.

[33] W.H.F. Smith, D.T. Sandwell, *Global seafloor topography from satellite altimetry and ship depth soundings.* Science 277 (1997) 1957-1962.

[34] T.G. Farr, P.A. Rosen, E. Caro, R. Crippen, R. Duren, S. Hensley, M. Kobrick, M. Paller, E. Rodriguez, L. Roth, D. Seal, S. Shaffer, J. Shimada, J. Umland, M. Werner, M. Oskin, D. Burbank, D. Alsdorf, *The Shuttle Radar Topography Mission.* AGU Review of Geophysics, 45(2) (2007).

[35] C.J. Olson, J.J. Becker, D.T. Sandwell, *A new global bathymetry map at 15 arcsecond resolution for resolving seafloor fabric: SRTM15_PLUS.* AGU Fall Meeting Abstracts (2014).

[36] D.T. Sandwell, R.D. Müller, W.H.F. Smith, E. Garcia, R. Francis, *New global marine gravity model from CryoSat-2 and Jason-1 reveals buried tectonic structure.* Science, 346(6205) (2014) 65-67.

[37] W.-B. Doo, C.-L. Lo, S.-K. Hsu, C.-H. Tsai, Y.-S. Huang, H.-F. Wang, S.-D. Chiu, Y.-F. Ma, C.-W. Liang, *New gravity anomaly map of Taiwan and its surrounding regions with some tectonic interpretations.* Journal of Asian Earth Sciences 154 (2018) 93-100.
[38] M.-H. Chang, S. Jan, V. Mensah, M. Andres, L. Rainville, Y.J. Yang, Y.H. Cheng, *Zonal migration and transport variations of the Kuroshio east of Taiwan induced by eddy impingements*. Deep-Sea Research, I 131 (2018) 1-15.

[39] P. Lemenkova, *Statistical Analysis of the Mariana Trench Geomorphology Using R Programming Language*. Geodesy and Cartography 45(2) (2019) 57-84.

[40] I.A. Suetova, L.A. Ushakova, P. Lemenkova, *Geoinformation mapping of the Barents and Pechora Seas*. Geography and Natural Resources 4 (2005) 138-142.

[41] M. Klaučo, B. Gregorová, U. Stankov, V. Marković, P. Lemenkova, *Determination of ecological significance based on geostatistical assessment: a case study from the Slovak Natura 2000 protected area*. Central European Journal of Geosciences 5(1) (2013) 28-42.

[42] M. Klaučo, B. Gregorová, U. Stankov, V. Marković, P. Lemenkova, *Land planning as a support for sustainable development based on tourism: A case study of Slovak Rural Region*. Environmental Engineering and Management Journal 2(16) (2017) 449-458.

[43] P. Lemenkova, I. Elek, *Clustering Algorithm in ILWIS GIS for Classification of Landsat TM Scenes: a Case Study of Mecsek Hills Region, Hungary*. Geosciences and Environment, Section ‘Near-Surface Geophysics’ (2012).

[44] P. Lemenkova, *Processing oceanographic data by Python libraries NumPy, SciPy and Pandas*. Aquatic Research 2 (2019) 73-91.

[45] P. Lemenkova, *Geospatial Analysis by Python and R: Geomorphology of the Philippine Trench, Pacific Ocean*. Electronic Letters on Science and Engineering 15(3) (2019) 81-94.

[46] H.W. Schenke, P. Lemenkova, *Zur Frage der Meeresboden-Kartographie: Die Nutzung von AutoTrace Digitizer für die Vektorisierung der Bathymetrischen Daten in der Petschora-See*. Hydrographische Nachrichten 25(81) (2008) 16-21.

[47] P. Lemenkova, *Plotting Ternary Diagrams by R Library ggtern for Geological Modelling*. Eastern Anatolian Journal of Science 5(2) (2019) 16-25.

[48] P. Lemenkova, *K-means Clustering in R Libraries {cluster} and {factoextra} for Grouping Oceanographic Data*. International Journal of Informatics and Applied Mathematics 2(1) (2019) 1-26.

[49] S. Gauger, G. Kuhn, K. Gohl, T. Feigl, P. Lemenkova, C. Hillenbrand, *Swath-bathymetric mapping, The expedition ANTARKTIS-XXIII/4 of the Research Vessel 'Polarstern' in 2006*. Reports on Polar and Marine Research, 557 (2007) 38-45.

[50] G. Kuhn, C. Hass, M. Kober, M. Petitat, T. Feigl, C.D. Hillenbrand, S. Kruger, M. Forwick, S. Gauger, P. Lemenkova, *The response of quaternary climatic cycles in the South-East Pacific: development of the opal belt and dynamics behavior of the West Antarctic ice sheet*. Expeditionsprogramm Nr. 75 ANT XXIII/4, AWI for Polar and Marine Research (2006).

[51] P. Lemenkova, *Geomorphological modelling and mapping of the Peru-Chile Trench by GMT*. Polish Cartographical Review 51(4) (2019) 181-194.

[52] P. Lemenkova, *GMT Based Comparative Analysis and Geomorphological Mapping of the Kermadec and Tonga Trenches, Southwest Pacific Ocean*. Geographia Technica 14(2) (2019) 39-48.

[53] P. Lemenkova, *Topographic surface modelling using raster grid datasets by GMT: example of the Kuril-Kamchatka Trench, Pacific Ocean*. Reports on Geodesy and Geoinformatics, 108 (2019) 9-22.