Quantifying deformation in North Borneo with GPS

The existence of intra-plate deformation of the Sundaland platelet along its eastern edge in North Borneo, South-East Asia, makes it an interesting area that still is relatively understudied. In addition, the motion of the coastal area of North-West Borneo is directed toward a frontal fold-and-thrust belt and has been fueling a long debate on the possible geophysical sources behind it. At present this fold-and-thrust belt is not generating significant seismic activity and may also not be entirely active due to a decreasing shelfal extension from south to north. Two sets of Global Positioning System (GPS) data have been used in this study; the first covering a time period from 1999 until 2004 (ending just before the Giant Sumatra–Andaman earthquake) to determine the continuous Sundaland tectonic plate motion, and the second from 2009 until 2011 to investigate the current deformations of North Borneo. Both absolute and relative positioning methods were carried out to investigate horizontal and vertical displacements. Analysis of the GPS results indicates a clear trend of extension along coastal regions of Sarawak and Brunei in North Borneo. On the contrary strain rate tensors in Sabah reveal that only insignificant and inconsistent extension and compression occurs throughout North-West Borneo. Moreover, station velocities and rotation rate tensors on the northern part of North Borneo suggest a clockwise (micro-block) rotation. The first analysis of vertical displacements recorded by GPS in North-West Borneo points to low subsidence rates along the western coastal regions of Sabah and inconsistent trends between the Crocker and Trusmadi mountain ranges. These results have not been able to either confirm or reject the hypothesis that gravity sliding is the main driving force behind the local motions in North Borneo. The ongoing Sundaland–Philippine Sea plate convergence may also still play an active role in the present-day deformation (crustal shortening) in North Borneo and the possible clockwise rotation of the northern part of North Borneo as a micro-block. However, more observations need to be collected to determine if the northern part of North Borneo indeed is (slowly) moving independently.

methods were carried out to investigate horizontal and vertical displacements. Analysis of the GPS results indicates a clear trend of extension along coastal regions of Sarawak and Brunei in North Borneo. On the contrary strain rate tensors in Sabah reveal that only insignificant and inconsistent extension and compression occurs throughout North-West Borneo. Moreover, station velocities and rotation rate tensors on the northern part of North Borneo suggest a clockwise (microblock) rotation. The first analysis of vertical displacements recorded by GPS in North-West Borneo points to low subsidence rates along the western coastal regions of Sabah and inconsistent trends between the Crocker and Trusmadi mountain ranges. These results have not been able to either confirm or reject the hypothesis that gravity sliding is the main driving force behind the local motions in North Borneo. The ongoing Sundaland-Philippine Sea plate convergence may also still play an active role in the present-day deformation (crustal shortening) in North Borneo and the possible clockwise rotation of the northern part of North Borneo as a micro-block. However, more observations need to be collected to determine if the northern part of North Borneo indeed is (slowly) moving independently.

Introduction
The Sundaland block is surrounded by three main plates known as the Philippine Sea, Australian and Indian plates (Simons et al. 2007). The Eurasian plate and Yangtze or South China block is attached to the north of Sundaland block. Generally speaking North Borneo is considered to be part of Sundaland near its eastern margins, close to highly active deformation zones bordering the south Philippine and Sulawesi Islands (Fig. 1). North Borneo is located at the junction between the North-West Borneo Trough and the Sulu Sea fold-and-thrust belt, to the north of the Celebes Sea. The way these two structures interact or connect unfortunately is still poorly understood (Hall 2013). At the northern end of the Crocker Range in North Borneo, Mount Kinabalu is the highest mountain peak in SE Asia and was formed by melting rocks deep in the crust about eight million years ago. It is composed of granite and has a ∼4100 m summit height above local mean sea level (Hall et al. 2009). Mount Kinabalu is assumed to be still rising at a long term rate of about 0.5 mm/year ). The region is characterized by low to moderate seismicity and most earthquakes occurred at depths of less than 50 km. There are also no significant major active faults that cut throughout Borneo. Nonetheless, North Borneo nowadays seems to be still actively deforming in possibly different ways. The region thereby appears to be slowly moving and/or deforming differently from the rigid Sundaland block with a clockwise rotation (Rangin et al. 1999;Simons et al. 2007;Sapin et al. 2013). In addition, a longlived debate about the cause of the present motions of North-West Borneo as driven by gravity sliding (King et al. 2010;Sapin et al. 2013;Hall 2013) or crustal shortening (Simons et al. 2007;Franke et al. 2008;Hesse et al. 2009) is still ongoing. In previous studies, the lack of a dense Global Navigation Satellite System (GNSS) network was a critical constraint to further analyze the behavior of North-West Borneo.
The recently proven ongoing (surface) deformations in North Borneo make it therefore a very appealing area to study in more detail. Moreover, one expected source of deformation in this area is a deep water fold-and-thrust belt. Franke et al. (2008) and Hesse et al. (2009) suggested that the unbalanced compression-extension observed offshore in the northern part of North Borneo may accommodate remaining convergence between Sundaland and Philippine Sea plates in this region. Meanwhile, Sapin et al. (2013) proposed the collapsing of the Crocker Range as the main mechanism to balance this offshore extension deformation deficit. The status of deformation is an important input for hydrocarbon accumulation studies and useful for future explorations (Ingram et al. 2004;Morley et al. 2011).
GNSS has become a great space-geodetic tool for geodynamics studies and has been used extensively in the past three decades. Hereby the Global Positioning System (GPS) is the pioneer GNSS satellite positioning system. GPS campaigns in Eastern Malaysia with repeated observations gradually were replaced with continuous observations. Nowadays, many continuous stations have been established, providing a better opportunity for understanding the geophysical sources behind the ongoing tectonic motions.
We therefore examine the recent deformation of North Borneo using the latest continuous GPS data set. Special attention is given to the representation of the computed GPS station positions that are used to define the Sundaland block motion in the International GNSS Service (IGS) realization (IGS08, Rebischung et al. 2012) of the International Terrestrial Reference Frame 2008 (ITRF2008, Altamimi et al. 2011) and use it as a local reference to investigate the motions in North and North-West Borneo. Individual baseline analyses are carried out to identify both spreading and shortening in North-West Borneo. Strain and rotation rates were also produced to provide a clearer view on the ongoing deformations in this region. First estimates of the recorded vertical motions in North-West Borneo are addressed after which conclusions and further recommendations are presented.

GPS data
The recent collaborations between European (Netherlands and France) and Asian (Malaysia, Thailand and Indonesia) universities and government agencies (e.g., in the EU-ASEAN SEAMERGES and GEO2TECDI projects) provided the main input to gather the data in this study. The continuous GPS stations were used to determine the absolute inter-seismic motion of the Sundaland plate before the tremendous impact of the 2004 Mw 9.2 Sumatra-Andaman earthquake. Available GPS data between 1999 and Boxing Day 2004 from 31 local continuous stations in Malaysia, Thailand and Indonesia (Fig. 2) were used. These three countries along with Laos, Cambodia, Singapore and Vietnam make up the main part of the Sundaland block. In Malaysia, 18 continuous stations cover the Peninsular Malaysia (ten stations) and Borneo (eight stations). The other stations are located in Thailand and Indonesia with six and seven stations, respectively. Final GPS solutions consist of 3-dimensional (3D) station coordinates. These coordinates are time dependent due to plate motion and any earthquake (co-and post-seismic) displacements. Representing these coordinates in the ITRF2008 (version IGS08) will produce useful solutions in a common global reference frame. To facilitate this and to also obtain velocities (compatible with ITRF2008) in this global reference frame, 74 IGS stations were added to the SE Asia GPS data set (blue, red and yellow diamonds in the inset in Fig. 2; listing per station name is given in Table M of the auxiliary material). This selection was based on data availability of the 91 IGS08 core stations in 1999 (59 stations) plus 11 additional IGS08 stations (ALGO, FAIR, HYDE, KIT3, KOSG, ONSA, PERT, SHAO, VILL, WUHN, YELL) and also 4 non-usable IGS stations in SE Asia: BAKO, CUSV, NTUS and SAMP (shown as yellow diamonds in Fig. 2). The latter 4 are not included in the IGS08 as they experienced strong post-seismic motions due to the 2004 Mw 9.2 and following earthquake events west of Sumatra. The included IGS sub-network contains stations on Sundaland  Bird (2003). Small blocks are indicated by twoletter abbreviations; Burma (BU), Molucca Sea (MS), Banda Sea (BS), Timor (TI) and Birds Head (BH). Relative velocities on selected boundaries are with respect to Sundaland and calculated from Morvel56-NNR (Argus et al. 2011). Absolute velocities are in IGS08 (ITRF2008) based on Sundaland motion from this study and Morvel56-NNR plate vectors relative to Sundaland. Major sources of geophysical processes on Sundaland are the converging motions of the Australian, Indian and Philippine Sea plates and its interaction with smaller blocks such as the BH, MS and BS. The biggest relative motions occur along the Philippine trench with values up to 10 cm/year. Black lines on land area indicate country boundaries. Malaysia comprises Peninsular Malaysia, Sabah and Sarawak. Meanwhile, Indonesia consists of Sumatra, Java, Kalimantan, Sulawesi, Papua and Maluku islands. Borneo Island is shared by Malaysia (Sabah and Sarawak), Indonesia (Kalimantan) and Brunei and surrounding regions. As such the IGS sub-network has a twofold purpose: to provide sufficient and globally welldistributed reference stations and to densify the regional GPS network. Historically, the initial few IGS stations in SE Asia together with newer additions have been repeatedly affected by numerous earthquakes (and their post-seismic aftermath) since end-2004. Therefore, their actual positions may differ (again) beyond the latest time epoch (2009.5) of the GNSS solutions included in the ITRF2008 (IGS08 version) (further on referred to as IGS08). Hence, 70 IGS stations outside this region with large data time span (1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011) and stable linear motion were selected as candidate reference stations (illustrated by the 32 blue and 38 red diamonds in the inset in Fig. 2; the latter symbol indicates those that finally were not used) for mapping in IGS08 as is discussed in Sect. 3.2. A well-distributed global IGS reference net-work is required to avoid unfavorable network effects during transformation of the GPS coordinates into global reference frames such as IGS08, especially when also vertical velocities will be derived (e.g., Bergeot et al. 2009;Collilieux et al. 2012;Legrand and Bruyninx 2009;Legrand et al. 2012). But likewise consistent use (both in size and geometry) of such reference networks throughout different time spans of GPS data from the same region is important. Our IGS sub-network selection for GPS data analysis in SE Asia considered both requirements (as in Simons et al. (2007) but expanded with ∼40 additional IGS08 core stations) and was used alongside the newer capability of the scientific GNSS software to also directly compute precise GPS solutions in IGS08.
The Department of Surveying and Mapping Malaysia (DSMM) is a government agency, which from 2003 on has been responsible for establishing and maintaining the  (Nordin et al. 2009; some examples are given in the auxiliary material) and therefore well suited for geophysics (e.g., Vigny et al. 2005;Simons et al. 2007). The latest expansion also includes stations in Borneo especially along the coastal area of North Borneo. These new GNSS stations have been providing reliable data to further study the geodynamics of North Borneo. Unfortunately, some older MASS stations such as KINA and LABU were shut down here due to urban development and monument instability. The possibil-ity that some (new) MyRTKnet stations are affected likewise is considered in this first detailed scientific analysis of the expanded network in this region. In 2008, many MyRTKnet stations also have a large data gap and exhibit noisy GPS results due to instability of date storage systems and communication failures. In order to have significant and consistent results, a continuous data span of 3 years was used from January 2009 to December 2011. The selected time period also minimizes the impact of any ongoing post-seismic deformations following the 2004 Sumatra-Andaman earthquake that even co-seismically displaced North Borneo up to 1 centimeter toward the west (Vigny et al. 2005). The data span consists of 20 continuous stations in North Borneo as shown in Table A (available in the paper auxiliary materials) that also provides the GPS antenna and receiver configurations. These data have been processed along with the other 58 MyRTKnet and 74 IGS stations (Fig. A in the auxiliary materials).

GPS processing
Dual frequency GPS data have been processed using the scientific GIPSY-OASIS II software version 6.1.2 (Jet Propulsion Laboratory 2013). In this study, the Precise Point Positioning (PPP) was implemented to derive precise coordinate results (Zumberge et al. 1997). Precise ephemeris of satellites along with Earth rotation parameters were obtained from the Jet Propulsion Laboratory (JPL) enabling consistent derivation of highly accurate daily absolute GPS positioning results. Furthermore (relative positioning) baseline analyses were performed to study local motions along specific coastal regions of North-West Borneo.

Daily solutions
The GPS data were decimated to 5-min intervals. Then the data were processed zero-differenced into daily coordinates making use of the ionospheric free combination of the observables. An elevation mask angle of 7 • was used. Absolute IGS antenna phase center corrections were applied to take into account the systematic errors of different elevation angle and azimuth of the observed satellites. The corresponding IGS08 absolute antenna phase center table (igs08_gpsweek.atx) was downloaded along with the weekly JPL orbits and clock products. Elevation mask angles can be below 10 • if absolute phase center variation tables are used (Schmid et al. 2005). The tropospheric mapping function used is known as the Vienna Mapping Functions (VMF1) and based on numerical data obtained from actual weather conditions. The database of the VMF1 mapping coefficient grid was downloaded from the Global Geodetic Observing System website maintained by the Vienna University of Technology (Boehm et al. 2006). Lyard et al. (2006) assimilated satellite altimetry (TOPEX/POSEIDON) data into their hydrodynamic tide model to produce FES2004. This model is one of the standard global ocean tide models besides GOT4.7 (Ray 2008). Hence, FES2004 was applied in the processing and ocean loading parameters were retrieved from the Onsala Space Observatory website (Bos and Scherneck 2014).
To enhance the coordinate solutions in the east-west direction in this region, ambiguity fixing has been applied by using the Ambizap algorithm (Blewitt 2008). Ambizap thereby provides unique ambiguity-fixed single network solutions. Moreover, Ambizap provides a very straightforward and fast technique for adding extra consistent station data to an already existing and processed network, due to the characteristic of PPP solutions that are not correlated with any other stations due their zero-differenced processing properties.

Weekly combinations in IGS08
For our own IGS subset mapping approach, seven-parameter Helmert transformations (3 rotations, 3 translations and one scale factor) were used to combine the daily (fiducial free) solutions into weekly averaged positions. The weekly averaging was performed to screen for any outliers and thereby improve the reliability of the coordinate solutions. Overall repeatability (precision) of the final daily station coordinates is 1.6, 1.9, 6.0 mm for the 1999-2004 (6 years) and 1.5, 1.6, 5.6 mm for the 2009-2011 (3 years) data periods in, respectively, the north, east and vertical components.
Next step is to map the weekly averaged coordinate solutions into the global reference frame IGS08. We compared transformation solutions using both the JPL X-files and an own well-defined IGS sub-network. For the X-file mapping approach, the corresponding X-files were applied, after which weekly averaged positions are computed. The X-file contains daily transformations parameters to IGS08 that are provided by JPL along with the ephemeris package (Jet Propulsion Laboratory 2013). The seven transformation parameters of each X-file are based on a global fit of IGS stations as daily computed by JPL and facilitate the direct transformation of single station position solutions into IGS08. Besides a reference station quality assessment, the IGS08 uses an updated set of absolute antenna calibrations (igs08.atx) from the one used for ITRF2008 (igs05.atx) (Rebischung et al. 2012). A further update released in October 2012 (IGb08, but not applicable to the data periods analyzed here) also includes coordinate updates for IGS08 stations that underwent discontinuities since the last (2009.5) data inputs into ITRF2008. Our data processing setup is compatible with the one used at JPL and hence the X-files provide us with daily coordinates in the IGS08. In general, X-files are very well suited to obtain both single point and network solutions globally but the user is fully dependent on both the quality and availability of these X-files.
For our independent mapping realization weekly reference station coordinate sets were extracted (using the valid reference position and velocity estimate at each epoch) from the IGS08 Solution INdependent EXchange Format (SINEX) provided by IGS since the changes in absolute antenna calibrations resulted in position corrections for 65 out of the 232 IGS stations in ITRF2008 between −1.3 ± 2.8, 0.1 ± 1.1 and 2.2 ± 3.3 mm in, respectively, the north, east and vertical components (Rebischung et al. 2012). Prior to the mapping, the IGS sub-network stations were screened to avoid network effects due to including IGS stations with position and/or velocity discontinuities in the 2009.5-2012.0 period by making use of the JPL GPS time series (in IGS08) and break estimates available at http://sideshow.jpl.nasa. gov/post/series.html as well as the ITRF2014 discontinuities file available at http://itrf.ign.fr/ITRF_solutions/2014/ doc/ITRF2014-soln-gnss.snx. That discontinuity screening already removed 24 candidate IGS08 reference stations. Also we did not use IGS stations for which no valid IGS08 solution (AUCK, FAIR, GOUG) or GPS data (ASC1, KWJ1, MALI, PETP) was available in both data periods. After a first IGS08 mapping round, we checked the weekly IGS08 coordinate residuals and further removed any stations that exhibited gradual/instant position offset (MAC1, COCO, IISC, PETP, GUAM) or big seasonal variations (MCM4, TIXI). As a result and also to conserve the IGS sub-network geometry in each data period (which was also an important criteria in its selection) 38 IGS stations (marked with red diamonds in the inset in Fig. 2) were not used in both solutions. Stations BAKO, CUSV, NTUS and SAMP (yellow diamonds in Fig. 2) have unstable time series due to co-and post-seismic events and are also not included the IGS08. The remaining 32 IGS stations (blue diamonds in Fig. 2) all have well determined positions and velocities in IGS08 to serve as connecting links with this global reference frame and also have a good distribution in the northern (19) and southern (13) hemisphere.
The pre-verified subset of IGS stations was used to independently estimate weekly seven-parameter Helmert transformations to the IGS08. Hereby an additional outlier detection criterion with upper margins of ±6, ±8 and ±20 mm in, respectively, the latitudinal, longitudinal and vertical direction was applied to remove any remaining (temporary) weekly outliers. As such a well-defined global IGS subset mapping network that weekly included 27 ± 5 and 30 ± 2 stations in, respectively, the 1999-2004 and 2009-2011 solutions were finally used (details available in Tables D and E in the auxiliary material). At the start of 1999, fewer IGS08 core stations initially were available, but later in the year all 32 became operational. Satirapod et al. (2013) using a similar IGS sub-network for ITRF mapping GPS stations in Thailand tested multiple geometry variations and concluded that 15 IGS stations in general suffice to get an accurate representation in the global reference frame. In general, IGS08 station position misfits with root-mean-square (RMS) values of 1.6, 1.8, 6.6, and 1.5, 2.0, 6.5 mm for, respectively, the north, east and vertical components were observed in mapping the 1999-2004 and 2009-2011 solutions (close to the repeatabilities of the daily coordinate solutions). A more extensive summary of the RMS results per year is available in the auxiliary materials (Tables F and G). The rejected IGS stations were mapped along with the local network in SE Asia. Figure 3 shows the IGS08 mapped coordinate differences of these two mapping strategies at station TAWX whereby the own IGS sub-network station mapping in IGS08 serves as reference. Other local stations showed a similar pattern of variation. In the 1999-2004 solution, variation of the latitude component also seems to indicate that a yearly cycle might be present. Contrarily, latitude components of the 2009-2011 solution appear to have more likely random variations also in the longitude and vertical position components. In latitude component, the differences fluctuate within −2.7 to 3.0 and −1.7 to 3.5 mm for, respectively, the 1999-  Table J of the auxiliary materials). This might indicate that although the automated X-files are well suited to apply anywhere globally, the manual use of a consistent (both in size and geometry) global IGS08 sub-network may result in a more consistent representation of coordinates in IGS08, especially if two different time periods are analyzed and compared. Always some aliasing of (small) position offsets in the Helmert transformations will take place during the (numerical) matching of the predicted (IGS08) and computed GPS position of each included IGS08 reference station. This is somewhat similar to the more extreme case of exclusively using a regional IGS sub-network where seasonal signals are reduced as these are aliased in the reference frame transformation parameters as has been observed by several other studies (Freymueller 2009;Legrand et al. 2012;Zou et al. 2014). In our case this aliasing seems to be more of a (seasonal) position misfit balancing nature that depends on the size and geometry of the chosen global IGS08 sub-network. Analysis of the 2009-2011 X-files also reveals that the JPL X-files (IGS08) not always use the same IGS08 network size and geometry (e.g., variations in size of ±20 reference stations on subsequent days in the same month) and also used some IGS08 stations at the end of 2011 that actually already had position discontinuities (e.g., KARR, IRKJ, BOGT, HARB, QAQ1 were still used on day 365 of 2011 for the daily mapping into IGS08). On average the number of used IGS08 stations to derive the daily X-files decreases from 47, 39 to 33 in, respectively, 2009, 2010 and 2011. This is due to the degradation of the IGS08 reference frame due to station position discontinuities and the main reason the IGb08 was released (Rebischung et al. 2012). In contrast, we use 32 (weekly 30±2) IGS08 stations throughout each time period (as well as in 1999-2004) and we do not make use any IGS08 stations with discontinu-  tively, ±8, ±9 and ±25 mm for the north, east and vertical component to identify any outliers in the velocity estimation. These (absolute position) outlier criteria are based on the error propagation of the daily coordinate repeatabilities and the IGS station position misfits in IGS08. They include extra tolerances of 2 and 4 mm for observed seasonal variations in, respectively, the horizontal and (most influenced) vertical position time series. The standard errors (99% confidence level) of each linear estimation are obtained from a Weighted Root Mean Square (WRMS) procedure defined as (2*WRMS)/T whereby T is the total observation period in years (Simons et al. 2007). This procedure has proven to be both a simple and safe method of ensuring the velocity errors are realistic, since formal errors coming from least-squares adjustment methods are typically too optimistic, and was validated in Simons et al. (2007) against procedures based on empirical noise models. For the 1999-2004 solution, the residuals with respect to the linear trends in Fig. 4 show that the IGS sub-network mapping approach marginally gives a better horizontal fit compared to the X-file mapping ( all components for the IGS sub-network mapping, whereas the X-file approach yields a slightly better fit of the vertical positions in IGS08. However, the estimated station velocities obtained using either one of these methods display only small differences at the sub-millimeter level per year as can be seen for the six year long time series at TAWX (Fig. 4). The differences between both methods for the other analyzed stations in Sundaland (1999Sundaland ( -2004  sub-network appears more in agreement with those of the BNKK and KMIT stations in Bangkok, whereas the X-file approach results in a 1.5 mm/year faster eastward motion estimate. The velocity differences for the 32 IGS08 stations used for mapping both data periods (Tables I and K in the auxiliary materials) have similar 3D values that vary in each direction due to their geometrical location in the IGS08 sub-network.
The estimated velocities for these IGS stations computed using the IGS sub-network approach were also compared with their nominal values in IGS08 (Tables I and K in the auxiliary materials). With values well below ±0.5 mm/year and ±2 mm/year for, respectively, the horizontal and vertical directions, the velocity differences remain within their error estimates (3-sigma) at all IGS08 stations in the 1999-2004 period. A few stations exhibit higher differences due to the position discontinuity jumps estimated in their 6-year time series (whereas in ITRF2008/IGS08 up to 12.5 years of position solutions are used to estimate both the velocity trend and any discontinuity jumps at once). The largest differences are typically observed in the vertical direction and most likely result from significantly less weekly solutions and temporary and/or gradual variations of velocity trends. The shorter 2009-2011 data period is more susceptible to this especially as it continues beyond the latest data inputs into ITRF08/IGS08. However, mainly a larger spread/variation in the velocity differences in each direction is observed here (Table K). Thereby the overall RMS values of the position misfits from the 32 IGS08 stations with respect to their linear trend estimates are 1.5, 1.7, 5.3 mm for the 1999-2004 and Hence, the velocity estimates derived with the own IGS08 global sub-network mapping approach exhibit no aliasing of seasonal signals in the Helmert transformation parameters as reported by studies that used a purely regional IGS sub-network approach (e.g., Freymueller 2009;Legrand et al. 2012;Zou et al. 2014).
Another quality indicator for the mapping procedure is the Root Mean Square (RMS) of the linear trend fitting for the stations in SE Asia (Fig. 5). The lower RMS values in the horizontal directions are most significant in the longitudinal direction for the 2009-2011 solution [panel (b) of Fig. 5] indicate that the IGS sub-network mapping method here performs better for almost all stations. The vertical direction shows a marginally better fit for the X-file solutions in both periods. As such the X-file solution appears slightly less consistent especially during the 2009-2011 period. We earlier have identified some quality control issues with them in this period and also the IGS08 sub-network used to compute them at JPL sometimes significantly varies both in size and geometry in short time periods, not to be confused with the gradual degradation of the IGS08 reference network due to occurrence of new undefined position discontinuities after epoch 2009.5 discussed by Rebischung et al. (2012). By the end of 2011, both the X-file and our own IGS08 sub-network approach make use of the same amount (∼ 32) of IGS08 stations but with a different geometrical global distribution. As explained before, we used the same consistent set of 32 IGS08 stations throughout mapping all coordinate solutions in both data periods and avoided using any IGS08 stations with new discontinuities after epoch 2009.5. The X-files results for the 2009-2011 period appear less reliable in the longitudinal direction for estimating velocities in IGS08 for stations in SE Asia. The X-files thereby do not make use of a consistent (size and geometry) IGS08 sub-network and (automated) IGS08 station discontinuity detection might not always have been optimal. The own IGS08 sub-network method uses the same 32 IGS08 stations in both data periods and hence is not plagued by degradation of the IGS08 reference network. That doesn't guarantee that it delivers more accurate velocity estimates, but any biases due to (small) misfits of IGS08 positions will be similar in both data periods, and cancelled out in the mostly relative motion analysis presented here. In the vertical velocity estimates, (small) constant differences were observed for stations in SE Asia that most probably result from different seasonal averaging/balancing during the Helmert parameter estimations as different IGS08 reference stations are used in both mapping approaches. The own IGS08 sub-network mapping approach also enabled independent verification and quality control of the IGS08 mapping approach that may also benefit from the multi-day averaging of IGS position solutions (including outlier removal) before reference frame transformations are estimated to the middle IGS08 epoch of each week.
Therefore, velocity estimates obtained with the IGS subnetwork mapping approach were finally used for the discussion of the GPS results in this paper (available from Tables B and C in the auxiliary materials for each data period). However, interpretation of the GPS velocities from the Xfile mapping approach will not provide different insights, as the velocity estimates are primarily biased by sub-millimeter per year values in both the horizontal and vertical directions. They also overlap with the IGS sub-network approach well within their formal error estimates. Hence, the impact of such small horizontal velocity differences on plate motion estimates is minimal and virtually absent in strain and rotation rate tensors that are based on relative velocity differences. Because the velocity estimates with the X-file method were concurrently computed, they have been included in the auxiliary materials (Tables H, I, J and K on SE Asia and IGS stations velocity comparisons for each data period).

Baseline processing
Baseline (relative positioning) analyses were also performed to study the interaction between two stations, because typically common atmospheric effects and mapping errors decrease in relative coordinate solutions especially for short baselines. In the baseline processing, the same software configuration was used with the exception that fixed coordinates were employed for the base station and satellite orbit clocks were estimated. For ambiguity fixing, double differences of wide and narrow lane phase biases were computed. The

Sundaland motion
The boundaries of the Sundaland platelet were defined from strain-rate tensor analysis except to the northwest, where the geological boundary of Sundaland and the South China blocks was used (Simons et al. 2007). A plate(let) motion can be described by a vector and a pole of rotation using Euler's rotation theorem. By using velocity vectors of GPS stations on Sundaland, the block motion can be estimated. In this paper, position time series of six years were used to estimate the velocity of stations from 1999 to 2004 in IGS08. No data on and after 26 December 2004 were used to avoid inclusion of station positions affected by co-and post-seismic motions due to the Sumatra-Andaman earthquake event (Vigny et al. 2005). These (inter-seismic) velocities were subsequently used as input to determine an IGS08 (ITRF2008) compatible rotation rate and pole of Sundaland. The Sundaland motion can then be removed from the IGS08 velocity field to highlight the relative motions of North and North-West Borneo with respect to the center core of Sundaland.

Previous estimates
There are many estimates on the Sundaland block available of which the majorities were published (Table 1) after 2002 when the ITRF2000 was released. All these estimates make use of GPS data on Sundaland during the inter-seismic phase (1994)(1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004). The solutions before 2002 were imprecise due to the less well determined ITRF1994 and ITRF1997 reference frames. Kreemer et al. (2003) developed a new no-net-rotation (NNR) model by combining geodetic site velocities and geological strain rates. Their NNR-GSRM1 reference frame implies an alignment with ITRF2000 in this region (Simons et al. 2007). As can be seen, the solutions by Bock et al. (2003), Simons et al. (2007), DeMets et al. (2010) and Argus et al. (2011) indicate quite similar results although using different selection of stations and reference frame. The estimate from Prawirodirdjo and Bock 2004 gives the lowest pole latitude and fastest rotation rate but made use of only two stations BAKO (Java) and NTUS (Singapore) resulting in a less accurate solution. Station BAKO also exhibits a noisier position time series and is located relatively close to the Australian-Sundaland plate boundary. Simons et al. (2007) have defined both the motion and boundaries of the Sundaland block. They used 28 stations to estimate the Sundaland block motion in ITRF2000; the highest number of stations so far. Their combination of campaign and continuous stations resulted in an exemplary outcome. Their motion was also updated (using the same GPS data set and sites on Sundaland) in ITRF2005 and illustrates the different velocity estimates obtained in ITRF2005 versus ITRF2000. Later, DeMets et al. (2010) used 18 station velocities by Simons et al. (2007) in their Sundaland estimation. This is 10 stations less than used by Simons et al. due to the location of these stations close to the Sumatra trench and one of them being located on the island of Java. However, their Sundaland pole location and rotation rate estimates remain almost identical. Then, Argus et al. (2011) used the outcome (MORVEL model) from DeMets et al. (2010) to determine their NNR-MORVEL56 model. The NNR-MORVEL56 is considered to be the best plate model at present, covering most (99.8%) of the surface on earth; 25 major and 31 small plates (Bird 2003;DeMets et al. 2010).
We would like to emphasize that any attempts to reestimate the Sundaland motion using or including GPS data after 26 December 2004 from Cambodia, Laos, Thailand, Vietnam, Peninsular Malaysia, western parts of Indonesia and West Borneo will need to take into account (by geophysical modelling) the post-seismic motions as a result of the Mw 9.2 Sumatra-Andaman earthquake. Stations in Peninsular Malaysia and Singapore still exhibit significant horizontal displacement rate offsets ranging from 10 to 25 mm/year in WSW direction when their inter-seismic (1999)(2000)(2001)(2002)(2003)(2004) and post-seismic (2009-2011) velocity estimates are compared (e.g., comparison of velocity estimates of ARAU and NTUS in Tables B and C in the auxiliary materials). Although their magnitudes diminish over time, post-seismic motions may remain noticeable here for decades.
The most recent solution was determined by Altamimi et al. (2012) as part of their study on plate motions in ITRF2008 to investigate the suitability of NNR-MORVEL56 in defining future ITRF solution releases. Their motion estimates for Sundaland could only make use of the two stations on Peninsular Malaysia (GETI and NTUS) included in the ITRF2008. Therefore, this solution is not constrained by the northern and eastern parts of Sundaland. Prior to the 2004 Mw 9.2 Sumatra-Andaman earthquake the GETI station also had a faster steady-state velocity due to inter-seismic loading of the locked subduction interface in front of Sumatra that accommodates convergence of Indian and Australian plates with Sundaland (Simons et al. 2007) which is also visible in Fig. 7. This became evident as the co-seismic displacement at GETI due to the Mw 9.2 earthquake was significantly larger than the one observed at NTUS (Vigny et al. 2005 Their overall findings did illustrate that ITRF2008 plate motion estimates have better agreement to NNR-MORVEL56 rather than NNR-NUVEL1A. NNR-NUVEL1A is a wellknown no-net rotation frame consisting of a bundle of geologically angular velocities from 14 plates (Demets and Gordon 1991). It is less significant for our work as it does not contain a definition for Sundaland (considered there to be still part of Eurasia). However, the large difference of rotation rates between NNR-MORVEL56 and ITRF2008 indicates that alignment accuracy of ITRF2008 to NNR-MORVEL56 is not better than 2 mm/year . Finally ITRF2000 velocities are not fully compatible with either ITRF2005 and ITRF2008 as velocities in latitude component (at least within SE Asia) were found to differ up to a few mm/year (Satirapod et al. 2011) and hence also affect the Euler pole and rotation rate estimations (e.g., Sundaland ITRF2000 and ITRF2005 results by Simons et al. (2007) in Table 1).

Pole and rate estimation
An updated absolute rotation pole and rate for the Sundaland plate needs to be computed, as the latest station velocities are in IGS08 (ITRF2008). It is obvious that the resulting rotation solution depends on the selection of stations assumed to be part of undeformed Sundaland. Therefore, selection of these stations also requires determination of the Sundaland block boundaries. However, not all of the available stations in Sundaland are located on the center core of this platelet but also located inside its wide boundaries with surrounding plates that are affected by deformation. Therefore, stations close to active faults, plate boundary zones and/or affected by local ground motions need to be excluded from the computations. A threshold limit of 3 mm/year for the residual velocities was selected as a rejection criterion in subsequent rotation vector estimations. This criterion was chosen as the majority of the estimated (1999)(2000)(2001)(2002)(2003)(2004) horizontal velocities in SE Asia have a 99% confidence level of ∼1.0 mm/year, while intra-plate deformation of ±1 mm/year may occur (Simons et al. 2007). Therefore, two stations can differ up to 3 mm/year in absolute velocity while still moving as part of the same quasi undeformed tectonic entity. Several iterations were carried out to detect and remove stations affected by local motions or exhibiting corrupted velocity estimates by any other errors. The iteration process continued until all the remaining stations passed the criterion. Finally, 14 stations were used in the final estimate of the updated Sundaland rotation pole and rate (Fig. 6). Four of the included stations have residuals above 2 mm/year with more than 70% of used stations exhibiting residuals below 2 mm/year. Four stations located on Thailand (SIS2 and KMIT) and Borneo (KUCH and BINT) exhibit a residual component in south direction. The remaining stations have the residual components toward the north direction and are located on Peninsular Malaysia and Singapore. The similar directions of the residuals in each of these regions indicate that intra-plate deformation of Sundaland is occurring far from its plate boundaries (e.g., due to elastic loading) up to 750 km away from its boundary with the Fig. 6 East (in blue) and north (in green) components of residual velocities for the 14 stations used in estimating the Sundaland rotation pole and rate in IGS08 (ITRF2008). Stations are listed per their longitude (east to west) whereby the purple bar represents the total horizontal magnitude of the residuals subducting Indian and Australian plates to the west (Simons et al. 2007). Therefore, it would be better to also have spatially well-distributed GNSS stations available on especially Natuna Island in the middle of the South China Sea, Vietnam, Laos and southwest Borneo. These areas have less deformation as shown in previous strain-rate analyses (Simons et al. 2007).
The final solution in IGS08 (ITRF2008) gives a Sundaland pole located at 48.05 ± 1.1 • N, 88.51 ± 0.4 • W with a rotation rate of 0.341 ± 0.015 • /Myr. Residual velocities of all available stations with respect to this solution are shown in Fig. 7. Intra-plate deformations have clearly contributed to the residual velocities of stations SAMP (Sumatra) and BAKO (Java) near the plate boundary zones of Sundaland. The velocities on Peninsular Malaysia show similar characteristics of deformation at 1-2 mm/year due to the impact of the Indian/Australian plate subduction interface (strain accumulation prior to the 2004 Sumatra-Andaman earthquake) beneath Sumatra. At least five stations (MIRI, LABU, KINA, SAND and TAWX) on North Borneo also illustrate clearly that intra-plate deformations are occurring here (Rangin et al. 1999;Kreemer et al. 2003;Simons et al. 2007;Sapin et al. 2013).

North Borneo study
SE Asia exhibits a complex crustal plate convergence configuration. Sundaland is wedged in between (inbound) plate convergence regions to the west, east and south, including a complex Sundaland/Australia/Philippine Sea plate triple junction in the southeast. The overall motion of the Sunda- land region is eastwards, at a rate of up to 1 cm faster than the Eurasian plate to the north (Fig. 1) and its likely separation from the Eurasian plate was an eventual result from the collision of the Indian and Eurasian plates followed by redistribution of crustal flow from the Tibetan plateau around the Eastern Himalayan Syntaxis (Simons et al. 2007) accompanied by a predominantly diffuse plate boundary to the north of Sundaland. Junctions between the converging plates create additional complexity of tectonic behavior mainly in the southeast of Sundaland (Sulawesi, Indonesia) at the triple junction with the Australian and Philippine Sea plates that are both subducting beneath Sundaland (Socquet et al. 2006). Between North-West Borneo and the South China Sea, the North-West Borneo Trough (also known as the Palawan Trough) was suggested by Hamilton (1979) to be a trench although he considered the trench to be extinct. Subduction of the Proto-South China plate occurred in this region continued up to the Early Miocene (∼ 20 Ma) (Hesse et al. 2009;Hall 2013). The actual trench involved here is assumed to be located south of the North-West Borneo Trough. Therefore, the North-West Borneo Trough could also be interpreted as the remainder of this original trench (nowadays located below the northwest Sabah basin). This trench may be gradually filled by eroding/gravity sliding sediment material from the Crocker Range that was deposited on the ocean floor (Bol and Hoorn 1980;Pubellier and Morley 2013). The Baram Delta is an active fold-and-thrust belt and situated south of the North-West Borneo Trough (Ingram et al. 2004;King et al. 2010;Fig. 9 Baseline solutions with fixed stations (triangles) MIRI, BEAU, UMSS and BELU. Error ellipses show 99% confidence level Sapin et al. 2013). The Crocker Range dates back from 35-25 Ma, but subsequently underwent uplift and widening from 15-5 Ma (Hall 2013). During this period (8-7 Ma) also the Mount Kinabalu granite was intruded into the already existing Crocker Range (Cottam et al. 2010(Cottam et al. , 2013. Deposition of the Rajang Group and Trusmadi Formation ended 40-35 Ma and predates the Crocker Range. There are (still) active faults that can be found along the Crocker and Trusmadi Ranges. The highest peak of Borneo is sited on the Crocker Range near a few intersecting active faults. Several small earthquakes have occurred in this region with the biggest earthquake of magnitude 6.0 recorded on June 4th, 2015 (Fig. 8). The denser 2009-2011 GPS network data were used to study the present deformation in North Borneo in more detail. Three stations (LAWS, LAB1 and MRDU) were excluded from the deformation analysis due to local subsidence and unstable monument history (vertical velocity >10 mm/year and a significant offset in the horizontal motion (magnitude/direction) compared to nearby stations). Velocities IGS08 (ITRF2008) plots and time series of these stations are available in the auxiliary materials (Fig.  B to E).

Baseline analysis of North-West Borneo
Four stations (time span 2009-2011) along the coastal area of North Borneo were selected as fixed stations (Fig. 9) to further investigate possible deformations in northwest direc-tion toward the North-West Borneo Trough (approximately perpendicular to the costal line). The MIRI-MRDI baseline shows an increasing baseline length of 4 ± 0.8 mm/year in southeast direction. This motion appears compatible with gravity sliding of material from the collapsing Crocker Formation complex sourcing the active fold-and-thrust belt of Baram Delta that is located close to MIRI station. GPS observations and offshore measurements clearly show present activity of the belt (Simons et al. 2007;Hesse et al. 2010;Hutchison 2010;Morley et al. 2011;Sapin et al. 2013). However, stations KENI and TENM show no relative motions with respect to station BEAU. Lengths of these baselines are approximately ∼30 km which might be too short to detect any relative displacements. Meanwhile, the UMSS-TMBN baseline illustrates these stations are closing in on each other at a rate of 2±0.9 mm/year. Further to the north, the BELU-RANA baseline, however, again indicates spreading in north/south direction. RANA, TMBN and KENI are also located in valleys on the opposite site of the Crocker Range. The baseline changes thus present inconsistent magnitudes and directions if attributed to common geophysical processes along and behind the coastal region of North Borneo. If gravity sliding would be the main driver behind the observed (local) motions in North-West Borneo, spreading must appear along the coastal areas. However, some of the selected stations might be affected by local deformations that disguise their spreading motion. Contribution of other deformation processes might also be occurring here resulting from a common WNW (block) motion or from additional (compression) motion in stations near the complex collapse zone. Figure 10 demonstrates the interaction between stations in North Borneo from the 2009-2011 GPS solution. As can be seen in panel (a), there is a clear trend of extension in the coastal areas around MIRI and NIAH at rates of more than 70 nanostrain/year. This result agrees with the baseline analysis in the previous section. No significant local deformation is detected in southern and eastern parts of North Borneo. In North-West Borneo, from north to south, there is a mix of deformations occurring. Inconsistency of strain patterns makes it difficult to conclude what is going on in this area. TMBN might be affected by local deformation, as the station is located close to a fault line. Meanwhile, UMSS also behaves strangely because it is moving faster eastwards [in IGS08 (ITRF2008)] than the other stations and whether unstable monumentation (e.g., slowly tilting of a relatively high pillar) is a plausible cause needs to be verified by DSMM. Further analysis will also be needed to identify any additional corrupted velocities that may have been influenced by other local effects. In panel (b) of Fig. 10, most  (2009-2011 solution). a Red arrows show extension deformation, blue arrows represent shortening deformation, and green arrows indicate shear deformation. 1 nanostrain/year represents 1 mm deformation over a 1000 km baseline per year. b Red wedges represent clockwise rotation rates and blue wedges denote anti-clockwise rotation rates of the rotation rates in North-West Borneo have clockwise directions and vary between 0.2 and 2.5 • /Myr. Meanwhile, along the south coast anti-clockwise rotation can be observed while in the center of Borneo, no significant rotations seem to occur. For comparison purposes strain and rotation rate tensors were also computed for the 1999-2004 GPS solution but these prove to be less useful due to the availability of only six stations on North Borneo in that period (Fig. F in the auxiliary materials). But there the KINA (UMSS), SAND and TAWX triangle also indicates clockwise rotation at a rate of ∼ 1 • /Myr.

North Borneo motion
North Borneo is clearly moving slower than it would if it were part of a rigid Sundaland tectonic block. That is signaling it may be undergoing intra-plate deformation near Sundaland's boundary with the Philippine Sea plate. Figure 8 shows that the (relative) westward motions of North Borneo from the 1999-2004 solutions are slower (1-3 mm/year) than those from the 2009-2011 solution. Three stations are available in both solutions; TAWX, SAND and MIRI with differences in their velocities in the 1-3 mm/year range. KINA (stopped in 2007) is close to UMSS and both appear to have a similar (slower) relative velocity of 5 mm/year. Meanwhile, BINT (stopped in 2007) is nearby BIN1 and displays an offset in both magnitude and direction of, respectively, 2 mm/year and ∼ 30 • . LABU (stopped in 2006) exhibits a velocity pattern similar to KINA/UMSS. Although some small (regional) increase in westward motion (∼1-2 mm/year) during the 2009-2011 epoch might come from remnant post-seismic motion of the 2004 Mw 9.2 Sumatra-Andaman earthquake, comparison of the two solutions however does not indicate any significant intra-Borneo deformation due to differences in post-seismic relaxation (i.e., North-West Borneo spans ∼ 300 km distance in east-west direction and is situated 2000-2300 km from the Mw 9.2 epicenter). The North Borneo relative motion (2009-2011 solution) estimated in this study is slower (1-2 mm/year) than using the Altamimi et al. (2012) predicted Sundaland motion. The motion is also directed a bit more to the south (Fig. 8). Therefore, Altamimi et al. (2012) may have over-estimated the absolute Sundaland motion in ITRF2008 by being able to use only 2 stations on Peninsular Malaysia as discussed before.
Three coastal stations (BIN1, NIAH and MIRI) along the Sarawak basin to the southwest of Brunei indicate that the coastal area here moves faster (than BELA/MRDI land inwards) and is currently spreading at 3±1 mm/year. Thereby NIAH and BIN1 appear to have similar motions to MIRI but BIN1 does move about 2 mm/year slower. The direction of motion of the MIRI station is also oriented more northward toward the Baram Delta. BIN1, NIAH and BELA all move westwards but with different rates. The directions of their motions appear parallel with fault lines to the south of these stations. A strike-slip earthquake also happened near the Bukit Mersing line in the past two decades.
A clockwise rotation pattern can be noticed both in the station velocities (Fig. 8) and the rotation rate tensors (Fig. 10b) in the northern part (above TENM/TAWX) of North Borneo with magnitudes increasing from north to south. The direction of this apparent (uniform) rotation and motion of North Borneo is also directed toward the North-West Borneo Through. Thus, it is plausible that convergence between Sundaland and the Philippine Sea plate at the very high rate of 10 cm/year (illustrated in Fig. 8) is still contributing to the tectonic motion here as a remaining convergence force that needs to be accommodated beyond the Sulu and Celebes seas by North Borneo (Simons et al. 2007;Franke et al. 2008;Hesse et al. 2009). The direction of relative velocity of the Philippine Sea plate toward Sundaland is also perpendicular to the North-West Borneo Trough. The rotation pattern also might be explained by initiation of a collision between the eastern Borneo shelf and the northwestward moving North Sulawesi micro-block (Socquet et al. 2006), although that motion is more likely to be taken up in the offshore West Sulawesi fold-and-thrust belt. A number of earthquake moment tensors in NE Borneo (Fig. 8 near TAWX, SEMP and DATU) also have strike-slip signatures indicating differences in motion along the Sulu arc and the Celebes Sea. However, the findings in this study do not support the suggestion made by Simons et al. (2007) that a trans-Borneo continuation from the Tinjar fault (near the West Baram line) toward the Sangkulirang fault in East Kalimantan, Indonesia might exist to accommodate North Borneo separation from the remainder of Borneo. The slower velocity magnitudes of MRDI and BELA indicate that this western area cannot be part of a northern North Borneo micro-block rotation. Sapin et al. (2013) and Hall (2013) proposed a collapse of the North-West Borneo Wedge (Crocker Range) after the end of convergence as a source for the gravity sliding deformation of North-West Borneo. Therefore, a vertical trend analysis was carried out with the denser 2009-2011 network solution. To date no vertical GPS displacements of this region have been published. Although the available data set just spans 3 years, a first indication of any ongoing vertical motions can already be obtained, whereby also seasonal variations (Blewitt and Lavallée 2002) are taken into account. The results (at 99% confidence level) in Fig. 11 indicate negative vertical motion trends for all but two (RANA and TENM) stations in North-West Borneo ranging from −0.3 ± 1.68 to −2.8 ± 1.4 mm/year at 99% confidence level. The trend rate of 0.1 ± 2.1 mm/year for TENM, although close to zero, can be regarded as an uplift signal. The result for RANA suggests uplift at a trend of 0.9 ± 1.8 mm/year and is located close to a few active faults and Mount Kinabalu (Fig. 8) whereby the latter is assumed to be still growing at a comparable rate of 0.5 mm/year ). These overall negative vertical rates along the western coasts of North-West Borneo are supportive of a present collapse of the Crocker Range. Meanwhile, on the eastern side of North Borneo, vertical trends are small except at TAWX, which demonstrates a clear uplift rate of 2.4±1.5 mm/year. Any separation between North Borneo and Sundaland might be located south of Tawau where the rotation of the northern part of North Borneo is accommodated. Besides, SEMP moves faster at rate 1 mm/year toward TAWX (Fig. 8). To accommodate the impact of convergence, the land in the Tawau area might be lifted up. Nevertheless, most errors (99% confidence level) on the vertical trends are still rather large in comparison with the small magnitudes of the observed trends. A longer time span of data is going to be needed from DSMM to evaluate the vertical component with higher confidence. In general, based on error propagations, a continuous position time series of at least 5 years will be required to have a good estimate of any slow rate vertical motions (1-2 mm/year) at 99% confidence level.

Conclusions/recommendations
Data from a significantly densified GNSS network on North Borneo have been analyzed using both absolute and relative GPS processing strategies. Thereby additional attention was given to optimally map absolute station positions in the IGS08 (ITRF2008) global reference frame. For this regional network in SE Asia, the inclusion of the same consistent (both in size and geometry) IGS08 sub-network of 32 steering stations in both the analyzed data periods (1999-2004 and 2009-2011) was given preference over using downloadable daily global X-file transformation parameters made available by JPL. However, the IGS08 horizontal and vertical velocity estimates obtained with either method do not significantly differ and are mainly biased at sub-millimeter per year level and remain well within their error estimates at 99% confidence level. Hence, the use of either approach will not result in different geophysical interpretations. Independent representation of GPS results in any global reference frame, though more time consuming, does provide an additional mapping quality check especially if reliable vertical velocity estimates are required.
The present day (relative) velocities of GPS stations in North-West Borneo do not illustrate any significant extension that exclusively supports a gravity-driven collapse of the entire Crocker Range into an offshore fold-and-thrust belt. At a number of GPS stations, the observed horizontal velocity magnitudes and directions may be influenced by local motions and although vertical displacements signify downward motion these estimates still come with relatively high uncertainties. However, the likelihood of all coastal stations in North Borneo moving toward the North-West Borneo Trough due to gravity sliding appears not possible due to small strain rates in the middle of North Borneo. Therefore, the findings of this study cannot support the previous research by Sapin et al. (2013) that suggests that at present mainly gravity sliding is occurring in North-West Borneo. Further investigation is required to determine whether the main driving force behind the apparent motion of the northern part of North Borneo could be the result of the high Sundaland-Philippine Sea plates convergence rate either directly as a residual from different tectonic loading and deformation patterns (Simons et al. 2007) or indirectly to accommodate the triple plate junction convergence through a series of interacting tectonic micro-blocks as also present in Sulawesi (Socquet et al. 2006). Therefore, the latest results justify support for denser GPS networks (land inwards) and access to longer time spans of data. In the future new stations from the Indonesian INA-CORS GNSS network (with also stations located more land inwards in the Indonesian part of Borneo) will allow us to directly verify if any remainder of the Philippine Sea/Sundaland plate convergence presently is also accommodated by deformation in North Borneo and model possible elastic loading effects similar to the convergence at the western boundaries of Sundaland (Sumatra). A new GNSS station (MKNB) installed recently (2013) by TU Delft and UiTM on stable bedrock near the summit of Mount Kinabalu might also provide additional insights in the observed horizontal and vertical deformations. Thereby it might be important to also explore the present day interaction between the North-West Borneo Trough, the Sulu fold-and-trust belt, the Sulu arc and the Celebes Sea as the associated deformation patterns are likely to transfer to the surface. Hall (2013) provides a very comprehensible overview on the (possible) tectonic evolution of North Borneo in geological time scale, but a gap on deciphering the neotectonics in current time still needs to be bridged. Thereby quantifying North Borneo as a slowly clockwise rotating and westward moving microblock (associated with low seismicity) cannot be excluded. This study, however, confirms that station velocities in North Borneo differ significantly from those of an undeformed Sundaland and that gravity sliding alone along its western coastal zones cannot explain the presently observed deformation patterns. It is more plausible that ongoing gravity sliding nowadays is actually co-acting with the slow westward motion of North-West Borneo toward the North-West Borneo Trough.
In order to obtain more reliable results on the vertical displacements, more information is needed to accurately determine the small trends of deformation. In a follow-on study, Interferometric Synthetic Aperture Radar (InSAR) can be beneficial to support vertical solutions from GPS. InSAR is capable of measuring vertical displacement with submillimeter accuracy from analyzing SAR image time series (Hooper 2006;Ferretti et al. 2007). More detailed geological information on Borneo and its surroundings is required to evaluate the possible slow rotation of northern North Borneo as a micro-block. Additional (micro) seismic and geological information can help toward tracing active faults and identify possible separation faults and/or widespread deformation zones between Sundaland and the northern tip of North Borneo. Thereby the earthquake mechanism behind the 2015 Mw 6.0 tremor near Mount Kinabalu may provide more insights into the tectonic processes that are currently active in this region.