Outlier Detection and Optimal Anchor Placement for 3D Underwater Optical Wireless Sensor Networks Localization
Abstract
Location is one of the basic information required for underwater optical wireless sensor networks (UOWSNs) for different purposes such as relating the sensing measurements with precise sensor positions, enabling efficient geographic routing techniques, and sustaining link connectivity between the nodes. Even though various twodimensional UOWSNs localization methods have been proposed in the past, the directive nature of optical wireless communications and threedimensional (3D) deployment of sensors require to develop 3D underwater localization methods. Additionally, the localization accuracy of the network strongly depends on the placement of the anchors. Therefore, we propose a robust 3D localization method for partially connected UOWSNs which can accommodate the outliers and optimize the placement of the anchors to improve the localization accuracy. The proposed method formulates the problem of missing pairwise distances and outliers as an optimization problem which is solved through half quadratic minimization. Furthermore, analysis is provided to optimally place the anchors in the network which improves the localization accuracy. The problem of optimal anchor placement is formulated as a combination of Fisher information matrices for the sensor nodes where the condition of Doptimality is satisfied. The numerical results indicate that the proposed method outperforms the literature substantially in the presence of outliers.
I Introduction
The deployment of underwater wireless sensor networks (UWSNs) opens up a wide range of possible applications as the large percentage of Earth’s surface is covered by water in the form of oceans, seas, and rivers. UWSN applications include but not limited to assisted navigation, disaster management, offshore exploration, ocean sampling, and underwater monitoring [2, 3]. However, connecting underwater sensors is a challenging task due to the hostile underwater wireless channel conditions [4]. Radio frequency (RF), acoustic, and optical systems are three main forms of underwater wireless communications and each has virtues and drawbacks. Underwater RF systems can provide average speed data rate, however, they can only operate at surface levels and cannot operate at higher depths due to significantly high absorption. Therefore, acoustic channels are mostly preferred for UWSNs due to its low absorption and long transmission distance. Nonetheless, acoustic systems suffer from bandwidth scarcity and low propagation speed (i.e., 1500 m/s) which yields a considerable latency at large distances. To overcome the problem of low data rate and high latency acoustic signals, underwater optical wireless communication (UOWC) has recently attracted attention by its high data rates and low latency levels [5]. Nevertheless, the low transmission range of UOWC systems necessitates practical networking and control solutions to realize underwater optical wireless sensor networks (UOWSNs) in real life.
The unique properties of UOWC also necessitate for an innovative reexamination of the localization problem as it is also an important task for UOWSNs to geographically tag the data, mitigate the limited transmission range via multihop routing based on sensor locations, and to sustain the link connectivity and quality via pointing, acquisition, and tracking (PAT) mechanisms. Arguably, UOWSN localization is a far more challenging task than terrestrial localization due to the nonavailability of global positioning system (GPS) signals and distinct propagation characteristics of light beams in aquatic medium [6]. In addition, twodimensional localization methods can work well for terrestrial networks but the directivity of optical beams and 3D deployment of sensors in the underwater environment require to develop threedimensional (3D) localization methods which are more challenging due to the sparsity of UOWSN deployment and hostile underwater environment. Besides that, the number of anchor nodes (beacons) is limited for UOWC in comparison to terrestrial communication. All of these considerations make the 3D localization of UOWSNs a significant task where very few choices are available.
Besides the aforementioned challenges, the limited transmission range and other limiting factors of underwater optical channel such as absorption, scattering, turbulence, salinity, and air bubbles leads to develop a multihop network with outliers. Therefore for a centralized localization scheme, it is not only required to estimate all the pairwise distances but also to remove the outliers. Although a number of conventional matrix completion methods have been proposed in the past to estimate the missing pairwise distances, for example, singular value thresholding (SVT) [7], atomic decomposition for minimum rank approximation (Admira) [8], Optspace [9], augmented Lagrange multiplier method [10, 11], and alternating minimization [12], all of these matrix completion methods only estimate the missing elements of a partially filled matrix and do not consider the outliers, and are therefore not robust to the outliers.
Additionally, placement of anchors for localization is a nontrivial task and is of great interest to the localization research community. The literature on the problem of optimal anchor placement mainly considers a twodimensional scenario where a single node is to be localized with a minimum localization error [13, 14, 15, 16, 17, 18, 19]. Threedimensional optimal placement of anchors for a single vehicle localization is investigated in [20]. The analytical characterization of optimal sensor placement for threedimensional scenarios with multiple sensor nodes is still an open research problem.
In this paper, we propose a robust 3D localization method for UOWSNs with limited connectivity. The problem of network localization seeks to find the position of each node given that few anchor nodes and some of the noisy internode distances are available. Because of the limited UOWC transmission distance, we consider a multihop UOWC system to extend the communication range. Furthermore, it is crucial to estimate the missing distances accurately since the internode noisy distances are not always available for all node pairs. Also, the severe UOWC channel conditions can cause outliers to some of the internode distances. Hence, it is important to remove the outliers otherwise the outliers can propagate the error throughout the network and yield a large localization error. Additionally, we provide an analytical expression for the optimal placement of anchors which reduces the network localization error.
The main contributions of the paper can be summarized as follows:

The directive nature and limited transmission distance of UOWC lead to a partially connected network where there are many missing internode distances which are required to be accurately estimated. Hence, we develop a lowrank matrix approximation method which can accurately estimate the missing internode distances.

Some of the internode distance can have a large error and introduces outliers to which the conventional 3D network localization methods are quite susceptible. Consequently, a closedform convergent iterative solution is proposed which can accommodate the outliers.

An analytical optimal condition is provided for the anchor nodes to satisfy which maximizes the localization accuracy for multiple sensor nodes.
The remainder of the paper is organized as follows. In Section II the literature on UOWSNs localization, matrix completion, and optimal anchor placement is presented. In Section III, we introduce the network model and define the single hop pairwise distance measurements for UOWSNs. Section IV presents the closed form solution for the proposed 3D localization method where the missing pairwise distances are estimated and outliers are removed. Furthermore, Section IV also presents the optimal placement of anchors in the 3D underwater environment to improve the localization accuracy. Section V provide the numerical results to validate the performance of the proposed method. Finally, Section VI concludes the paper with a few remarks.
Ii Related Work
Localization of nodes in UOWSNs is of great importance which can enable numerous applications. Conventionally, localization methods are either range based or range free where the range based methods rely on different ranging methods to estimate the distances, and then, estimate position of the node based on the estimated distances. The rangefree localization methods provide coarse position estimation where usually the area containing the node is estimated. Since rangefree localization methods are not yet developed for UOWSNs, we focus on range based methods.
The twodimensional rangebased localization methods for UOWNs can be classified into two categories as distributed and centralized methods. The authors in [6] have proposed for the first time a recieved signal strength (RSS) and time of arrival (ToA) based distributed localization method. The authors have considered an optical base station (OBS) placed in a hexagonal cell which serves as an anchor node for the users. The users are able to estimate the distance to multiple OBSs and then estimate its position by using linear least square estimation. A centralized hybrid acoustic and optical RSS ranging based localization method have been proposed in [21] where the authors have considered a weighting strategy to give more importance to accurate ranging measurements. A centralized RSS based localization method has also been proposed in [22], where the nodes estimate and forward the single hop RSS based distances to the centralized node. The centralized node is then able to estimate the position of each node. However, all of the above UOWSNs localization methods are twodimensional and do not consider the 3D nature of the underwater environment. We refer the interested reader to [23] which presents a comprehensive survey on UOWSNs localization.
As the centralized localization methods are based on dimensionality reduction techniques which require all of the pairwise distance estimations between the nodes, robust matrix completion methods are needed to estimate the missing pairwise distances and remove the outliers. A number of researchers have provided theoretical constraints about the matrix rank, sampling scheme, and missing rate to effectively complete a matrix in low rank. In general, matrix completion schemes can be classified into three classes. The first class consists of matrix factorization based schemes where the missing entries of a partially filled matrix of size and rank are estimated from the product of two matrices and [24, 25]. Matrix completion schemes based on matrix factorization are nonconvex and they require prior information about the rank of the matrix. The authors in [26] have dynamically adjusted the rank of the matrix to estimate the missing elements via matrix factorization. The matrix completion problem in[26] is formulated as a nonlinear successive over relaxation problem which is solved by using linear least square method.
The second class of matrix completion schemes consists of nuclear norm minimization methods which are convex. The authors in [7] proposed a singular value thresholding (SVT) method to estimate the missing entries of a partially filled matrix. In [10, 11], the authors have applied an augmented Lagrange multiplier (ALM) method for nuclear norm minimization which pads the missing entries of the matrix with zeroes and considers the completed matrix as a sum of true complete matrix and error matrix. The alternating direction (AD) method is proposed in [12] to minimize the nuclear norm which is similar to the method proposed in [10, 11] but the error matrix is not explicit. In [27], the authors have proposed truncated nuclear norm (TNN) minimization which is a special case of SVT where the largest singular values are equal to zero. The third class of matrix completion schemes is based on manifold optimization methods. The authors in [9] have proposed a lowrank matrix completion method (also called Optspace) by minimization over the Grassman manifolds [9]. A lowrank matrix completion method based on manifold optimization is also proposed in [28] which uses the Riemannian manifold. However, all of these matrix completion methods only estimate the missing elements of a partially filled matrix and do not consider the outliers, and therefore, are not robust to the outliers.
Outliers can significantly contaminate the ranging measurements, and thus lead to large localization error. A number of outlier detection techniques have been presented in [29] for data mining, predictive modeling, cluster analysis, and wireless sensor networks. Recently, outlier rejection methods have been developed for multilateration based iterative localization [30]. However, in sparse networks, each node may not be directly connected to the anchors and thus cannot perform multilateration for localization. Hence, networks such as UOWSNs with sparsity and hostile underwater environment (absorption, scattering, geometrical losses, turbulence, and air bubbles) require to develop robust localization scheme to encounter the outliers produced by all these phenomena. Therefore, in this paper, we have developed a novel 3D localization algorithm which is not only robust to the ranging errors but also accounts for such outliers.
In addition to the presence of outliers, UOWSNs require to develop 3D localization methods which are more applicable for UOWSNs. Until now, 3D localization methods for UOWSNs do not exist, however, there exists a number of 3D localization methods for underwater acoustic sensor networks (UASNs). The authors in [31] have proposed for the first time a 3D localization scheme for UASNs where a projection based method is used to project the 3D localization problem into its 2D counterpart. The authors in [32] have proposed an anchorfree 3D localization scheme for UASNs where the anchors are able to float along the network. The sensor node computes the ranges to a mobile anchor and computes its 3D location. A 3D multihop localization method has been proposed in [33] where multilateration is used by the sensor nodes to compute their position by getting the measurements from at least three anchor nodes. A topdown approach has been used in [34] for UASNs localization where the sensor nodes near to the surface buoys (anchors) estimate their positions and then act as anchors. Similarly, a 3D localization method has been proposed in [35] with a small number of anchors where once a sensor node estimates its location it becomes a new anchor. Recently, a 3D localization method for UASNs have been proposed in [36] where the uncertainty in the anchor positions is considered.
Furthermore, the placement of anchors for network localization methods is also an important and challenging problem. The literature on optimal placement of anchors consists of two types of mathematical formulations. One is parameter optimization and the other is optimal control. The optimal control based formulation is usually adopted for path planning problems [37, 38] while parameter optimization is mainly used for optimal placement of anchors [13, 14, 15, 16, 17, 18, 19]. Parameter optimization method estimates the position of sensors such that the uncertainty in position estimation is minimized. The uncertainty in position is characterized by Fisher information matrix (FIM) which is the inverse of Cramer Rao lower bound (CRLB). Any unbiased estimator which can achieve the CRLB is efficient, and therefore, the position of anchors which maximizes the determinant of FIM for the position estimation of all nodes is said to be the optimal position of anchors.
The proposed localization method for UOWSNs is fundamentally different and practical than the aforementioned schemes where we have tackled multiple problems such as estimating the missing distances, removal of outliers, and optimizing the anchor positions to improve the localization accuracy. The derived expressions for the proposed 3D localization are general enough to be applicable for any multihop wireless communication network (including terrestrial and underwater wireless communication systems). However, we especially focused on multihop UOWSNs since the UOWC channel is more susceptible to the outliers caused by the intrinsic properties of optical light and underwater aquatic medium. Before going into the details of the proposed localization scheme, we formally introduce the network model.
Iii Network Model
Consider a 3D UOWSN architecture as shown in Fig. 1 where the sensor nodes are deployed at the seabed and their information is delivered to the surface buoys via relay sensors in a multihop fashion. It is assumed that the location of surface buoys are known and they act as anchor nodes in the network. To introduce the notations for seabed sensors, relay sensors, and surface buoys we denote their 3D position as , , and , respectively. Combination of the 3D positions of seabed sensors, relay sensors, and surface buoys yields where , is the three dimension position of node , and the dimensions of matrix is . If surface buoys, who act as anchors, are located at the same plane (i.e., located at the same depth), the 3D localization problem cannot be resolved regardless of how many surface buoys are deployed. Therefore, it is of great interest to find the optimal depth of the anchors which provide the better localization accuracy.
Iiia Optical Ranging for UOWSNs
Every localization technique requires to estimate the distance between the sensor nodes or between the sensor nodes and anchors. The distance can be estimated by using different ranging techniques such as ToA, time difference of arrival (TDoA), angle of arrival (AoA), and RSS. Each of the ranging techniques has its own advantages and disadvantages, e.g., time and angle based ranging methods are more accurate as compared to the RSSbased methods but are more complex and require extra hardware. Whereas, the RSS based methods are simple to implement but provide coarse localization. In this paper, we consider the RSSbased method where the distance between any two neighboring nodes is estimated by using the underwater optical wireless channel model. The underwater aquatic medium consists of different suspended and dissolved elements with different concentrations [39]. These elements cause the propagation of light to suffer from absorption and scattering. Haltran’s model is one of the wellknown models to account for the absorption and scattering by introducing the extinction coefficient given as
(1) 
where is the absorption coefficient, is the scattering coefficient, and is the operating wavelength. The major source for the absorption of optical light in underwater aquatic medium is chlorophyll, therefore is expressed as
(2) 
where and are constants, represents the absorption in pure water, represents the absorption coefficient for chlorophyll, is the absorption coefficient of fulvic acid, is the absorption coefficient of humic acid, is the concentrations of fulvic acid, and represents the concentrations of humic acid. and are given in [39] as
(3) 
and
(4) 
where and . Similarly the scattering coefficient is modeled as
(5) 
where is the scattering coefficient for pure water, represents the scattering coefficient for small particles, is the scattering coefficient for large particles, is the concentration of small particles, and represents the concentration of large particles [39].
Based on the extinction coefficient and the hardware specifications, the received power at node from node is modeled in [40] as follows
(6) 
where is the transmission power of node , and are the optical efficiencies of transmitter and receiver respectively, is the distance between nodes and , is the divergence angle of the transmitter, is the receiver aperture area, and is the trajectory angle between node and . Solving (6) for distance estimation yields
(7) 
where is the real part of Lambart function [41] and is the ranging estimation noise modeled as zero mean Gaussian random variable with variance . Fig. 2 shows the received power as a function of distance for different types of water where the simulation parameters are taken from [42]. It is clear from Fig. 2 that increase in the turbidity of water reduces the received power significantly which can cause outliers in distance estimation.
In Fig. 2 only absorption, scattering and geometrical losses are considered. However, the UOWC is also effected by other phenomena such as turbulence [43], air bubbles [44], and salinity [45]. These degrading phenomena can also introduce outliers into the distance estimation between any two nodes where the estimated distance strongly deviates from its true value. Therefore, , where is the Euclidean distance and is the outlier between any two arbitrary nodes and . Including the outliers into the ranging error deviates its distribution from normal to heavytailed where the outliers are far away from the mean value [46].
IiiB Construction of Pairwise Distance Matrix
The corresponding matrix of measured pairwise distances is denoted as . The noisy single hop pairwise estimated distances are collected at the surface station to form the observation distance matrix which is represented as
(8) 
However, it is very rare to directly obtain all the pairwise distances because of harsh propagation characteristics of UOWC channels. Instead, a subset of pairwise distances are available which are noisy and there is a great possibility to get outliers. Accordingly, matrix can be rewritten as
(9) 
where are the missing pairwise distances and represents the outliers. Hence, goal of the proposed localization method is to determine the 3D positions of seabed sensors and relay sensors, i.e., from matrix . As mentioned before that matrix has missing distances as well as outliers, therefore, it is first required to complete the missing distances and remove the outliers in as illustrated in Fig. 3.
Iv Proposed 3D Localization Method
In this section, first, we propose a novel lowrank matrix completion and outliers removal strategy for 3D localization of UOWSNs with limited connectivity. Matrix completion assigns values to the missing entries of a partially filled matrix and removes the outliers. Secondly, we investigate the impact of anchors placement on the localization accuracy and optimize the placement (especially the depth) of anchors to minimize the average localization error.
Iva Low Rank Matrix Completion and Outliers Removal
A number of conventional matrix completion methods have been proposed in the past, for example, singular value threshold [7], Optspace [9], augmented Lagrange multiplier method [10, 11], Poisson matrix completion [47], and alternating minimization [12]. All of aforementioned matrix completion methods are not robust to the outliers present in matrix and consider a least square loss function which measures the raw stress between and , i.e.,
(10) 
where represents the userdefined nonnegative weights to give importance to the measured distances. In majority of the cases, weights are considered to be equal to one. However, there are several other cases where unequal weights are employed such as Sammon mapping [48] with and elastic scaling [49] with . Error propagation for each element of matrix during the double centering method is given by where is known as the centering operation, is an identity matrix, is the vector of ones, and is the transpose operator [50]. Following from the error propagation of the double centering method, it is obvious that even a single outlier in can severely distort the solution for (10). Here the outliers are referred to the measured distances which are significantly different than the actual Euclidean distances. Therefore, it is necessary to investigate methods which can remove these outliers. Even if the problem defined in (10) is a nonconvex optimization problem without a unique solution, it can be solved by iterative majorization approach which minimizes
(11) 
where is the Frobenius norm. The effect of double centering method can easily be compensated by modifying the error function to norm, i.e., . However, norm is not smooth due to the fact that it has a singularity at its origin. To mitigate this problem, the use of Huber’s loss function can be of great benefit which interpolates between and norms minimizations, i.e.,
(12) 
where is the Huber’s loss function given as,
(13) 
is the argument of the Huber’s loss function which represents the residual error and is the threshold which can be chosen adaptively or arbitrarily from matrix . It is argued by Huber that if then Huber loss function is 95% as efficient as the least square solution [51]. Similarly, the loss function proposed by Tukey provide the same result with where
(14) 
As the estimated distance between any two nodes and is modeled as , where represents the outliers and is a zero mean random variable to model the nominal errors. In the presence of outliers and nominal noise error, the problem defined in (12) follows contaminated distribution and therefore adopts Huber’s function for the residuals [52]. Nonetheless, the outliers are less in number and sparse, therefore inserting the norm minimization problem into (10) is justified and yields
(15)  
where matrix represents the estimated outliers. In order to get to the solution in (15), the following function needs to be minimized
(16)  
The first term in (16) corresponds to the level of fitness between and after removing the outlier and the second term corresponds to the penalty linked to the sparsity of matrix where represents the regularization parameter. The function in (16) is nonconvex and nondifferentiable therefore expanding it yields
(17)  
where . Consider that is a set of indicator vectors with elements if and zero if . Then the term in (17) can be written as
(18)  
where is the trace operator and is a matrix with nondiagonal elements equal to and diagonal elements equal to . The term in (17) is nonconvex and nondifferentiable especially when . Majorization of this term is required to minimize (17), therefore, defining a majorization function . can be split into two terms depending on the value of , i.e., , where
(19) 
and
(20) 
Majorization of is obtained by using CauchySchwartz inequality, given as
(21) 
where is a matrix of auxiliary variables. Majorization of is obtained by using the fact that is upper bounded by for any positive value of . Consider that and , then using the bound leads to
(22) 
where for . Finally adding (21) and (22) yields
where is Laplacian matrix. When equality holds between and , then elements of the Laplacian matrix are defined as [53]
(24) 
where , , and . Substituting (IVA) and (18) in (17) yields
The majorization function in (IVA) is convex and differentiable in terms of both and . Therefore, solving (IVA) by minimization, i.e.,
(26) 
Indeed, each entry of (26) corresponds to a standard least absolute shrinkage and selection operator (LASSO) solution which is expressed as
(27) 
where is the soft thresholding operator, and [54]. Similarly, the position estimates are obtained by minimizing the following function
(28) 
The closedform solution of (28) is obtained by using the first order optimality condition [55], i.e.,
(29) 
where is a MoorePenrose pseudoinverse of matrix . For this iterative process, the initial configuration of is randomly chosen whereas the elements of the outlier matrix are set to zero. Given , the estimation of via (27) is a regularization problem while for a given the estimation of is the least square optimization problem, i.e., . Even though the estimation of accounts for the outliers, it is still a least square solution and is greatly influenced by the outliers. Therefore, it is proposed to apply function on to the residual in (29) where is nonnegative and differentiable with respect to to impose smoothness such that
(30)  
The optimization problem defined in (30) can be solved by using half quadratic minimization (HQM). HQM was pioneered to reconstruct the signal and images from corrupted signal and images [56]. Recently, HQM functions have been widely used in computer vision, machine learning, and feature extraction [57]. Based on the theory of the conjugate function in HQM, the problem in (30) can be written as [56]
(31) 
where is known as potential loss function, , is matrix of auxiliary variables, is the quadratic function, and is the conjugate function of . The quadratic function in (31) is given in [56] as
(32) 
where is a constant parameter. Substituting (32) into (31) and adding the regularization parameter for yields
(33)  
The minimization function in (33) is called resultant cost function of , and is solved by using alternating minimization given in [58] as
(34) 
where is the minimizer function for auxiliary variables and
(35) 
The minimization problem in (35) can be rewritten as
(36)  
where . The closedform solution of (36) with respect to is obtained by using the first order optimality condition, i.e.,
(37) 
where is an matrix which constitutes of the 3D relative position estimation of each node in the network. After getting the relative position estimations , the obtained solution can be transformed into the global position estimation by using the positions of anchors in the network. The common procedure for transformation to the global position estimation is orthogonal Procrustes analysis given by [59]
(38) 
where , , and are the scaling, rotation, and translation elements respectively.
IvB Optimal Placement of Anchors
Optimal anchor placement is crucial for better position estimation. Therefore, in this section, we develop an optimal anchor placement strategy to improve the localization accuracy of the proposed technique. As in the previous sections, the pairwise distances are estimated between every node in the network, now it is possible to optimize the position of anchors to reduce the localization error. Consider that there are number of surface buoys in the network which also work as anchors, then the estimated distance between any arbitrary sensor node and th surface buoy can be written as
(39) 
where is the Euclidean distance to th anchor and is the corresponding range measurement error. The Cramer Rao lower bound (CRLB) of the th parameter of is defined as
(40) 
where represents the Fisher information matrix (FIM). The elements of FIM are given as
(41) 
where E is the expectation operator and is a vector of size containing the estimated distances of a single node to all the anchors. Since the ranging measurements are effected by Gaussian noise, then the probability density function for the range measurements is defined as
(42) 
The submatrices of are derived from (41) and given as
(43) 
(44) 
(45) 
(46) 
(47) 
(48) 
Now the CRLB for the threedimensional position estimation is obtained as
(49)  
(50)  
(51)  
Note that from the above analysis the CRLB is minimized by maximizing the determinant of FIM. Therefore, it is important to manipulate the geometries (position) of anchors to maximize the determinant of FIM. The above analysis is easy to follow for a single node but for multiple nodes in the network, it is possible that changing the position of anchors may reduce the localization error for some nodes but may increase the error for others. Therefore, the optimal position estimation of anchors for multiple nodes in the network can be formulated as a determinant (D)optimality criteria given as
(52) 
where are the optimal anchor positions, are the number of sensor and relay nodes, and stands for the determinant. Geometrically speaking, the analogue of the above maximization problem for twodimension networks result in a regular optimal placement of anchors. The regular placement of anchor means that the anchors are placed in such a way that there are no breaking triangles between them for twodimensional localization. Therefore, here we focus to optimize the depth of anchors to minimize the localization error. To find the optimal depth of anchors, the CRLB for the depth component, needs to be minimized which is formulated as
(53) 
The optimization problem defined in (53) is nonconvex and hard to solve. Alternatively, an approximate solution to this problem is achieved by using iterative gradient method as follows
(54) 
where is the iteration number, is the step size, and represents the gradient. The selection of step size is crucial for the iterative gradient method and there are several step size selection methods among which the most common is constant step size, i.e., . However, constant step size leads to poor performance and require a large number of iterations to acquire certain convergence level. Alternatively, the optimal step size can be achieved as
(55) 
In practice, the minimization problem in (55) is as hard to solve as the original optimization problem in (53). Therefore, a compromise is made between the simplicity of selecting a constant step size and the complexity of selecting the optimal solution by using Armijo rule [60] given by
(56) 
where and . The anchor positions in (54) are updated by using the step size calculated in (56). If then the anchor positions are updated by (54) while if then the CRLB is not reduced, thus the iterative process stops, and are considered to be the optimal position of anchors for the given network setup.
V Numerical Results
In this section, we provide numerical results to validate our proposed solution. In this regard, we consider two scenarios of ten sensor nodes and four relay nodes randomly distributed in a 100 volume as shown in Fig. 4. The transmission range is kept to 80 m to achieve a connected network, ranging error is 0.6 m, and four anchor nodes are considered with their projections at different depths for the global transformation. The performance metric for both of the scenarios is considered to be the root mean square error which is given as
(57) 
In Fig. 3(a) to Fig 3(d) the impact of outliers and depth of anchors are investigated. Following four different cases have been evaluated:

Case 1: Considering the random depth of anchors in the presence of 35 % outliers. In this case, the RMSE performance is 23.70 m as shown in Fig. 3(a) where the performance is greatly influenced by both the outliers and random depths of anchors.

Case 2: In this case, we consider that the outliers are still 35 % but the depth of anchors is optimized. In this case, the RMSE performance is improved to 8.46 m as shown in Fig. 3(b).

Case 3: In this case, the depth of anchors is random and the outliers are removed from the estimated pairwise distances. In this case, the RMSE performance is 11.64 m as shown in Fig. 3(c) where the performance is improved as compared to the first case.

Case 4: In this case, the depth of the anchors is optimized as well as the outliers are removed from the estimated pairwise distances. Therefore, the RMSE performance is greatly improved to 0.659 m as shown in Fig. 3(d). Thus, Fig. 4 concludes that removing the outliers from pairwise estimated distances and optimizing the depth of anchors incredibly reduces the localization error for UOWSNs.
To analyze the robustness of the proposed 3D localization method in the presence of outliers we consider a scenario of 50 seabed sensor nodes, 4 relay nodes, and 4 anchor nodes with their optimal depths. A total of 100 Monte Carlo simulation runs were performed and compared with well known matrix completion methods such as nuclear norm minimization, low rank matrix factorization (LMaFit) [25], atomic decomposition for minimum rank approximation (Admira) [8], low rank matrix completion over Riemannian manifold (LRMCRM) [28], and Optspace [9]. To compare the results, first, we need to find the optimal regularization value which provides the minimum RMSE. Fig. 5 shows the performance of the proposed method in the presence of 35 % outliers where the optimal regularization values for the given network setup is . Therefore, we set the value of for rest of the experiments. Fig. 6 shows the impact of outliers on the RMSE of the proposed 3D localization method. It is clear from Fig. 6 that the proposed method has better RMSE performance as compared to the other approaches because it accounts well for the outliers present in matrix .
Undoubtedly, the ranging errors have a negative effect on the accuracy of every localization method. Here, we examine the performance of the proposed technique in the presence of ranging error only. To examine the impact of ranging error we considered 50 seabed sensor nodes, 4 relay nodes, and 4 anchor nodes with their optimal depths. Assuming that the ranging errors are Gaussian distributed with zero mean and variance , where the values of are set to 01 m. Note that the results are averaged over 100 different network setups. Fig. 7 shows that the proposed localization method have same performance as Optspace in the absence of outliers.
As mentioned earlier that the proposed 3D localization method optimizes the location of anchors to improve the localization accuracy. Therefore, we have evaluated two different scenarios with 14 and 54 sensor nodes respectively as shown in Fig. 8. The number of anchors is set to 4 for both setups and RMSE performance is evaluated in terms of noise variance. As illustrated in Fig. 8 that optimizing the depth information of anchors considerably improves the localization performance in both scenarios. Note that the results are averaged our 100 different network setups for both scenarios. Therefore, the results in Fig. 8 conclude that optimizing the depth of anchors in a network significantly improves the localization accuracy.
Vi Conclusions
Even though twodimensional localization methods for UOWSNs have been investigated in the past, 3D nature of UOWC environment requires to develop 3D localization methods. Therefore, we have proposed a robust 3D localization method for UOWSNs with limited connectivity. As the transmission distance of underwater optical sensors is limited, it leads to a partially connected network and many of internode distances are missing. Hence, we have employed a lowrank matrix approximation method which can accurately estimate the missing internode distances. Additionally, some of the estimated internode distances may have a large error and naturally introduces outliers. The traditional 3D network localization methods are susceptible to these outliers. Moreover, the placement of anchors for network localization methods is also an important and challenging problem. Therefore, a closedform convergent iterative solution is proposed which can accommodate these outliers and optimize the placement of the anchors to improve the localization accuracy. Numerical results validate the performance of the proposed method by showing accurate and robust results to the ranging errors and outliers, respectively.
Vii Acknowledgement
The authors would like to thank the anonymous reviewers for their fruitful comments.
References
 [1] N. Saeed, A. Celik, T. Y. AlNaffouri, and M.S. Alouini, “Robust 3D Localization of Underwater Optical Wireless Sensor Networks via Low Rank Matrix Completion,” in IEEE Int. Works. on Signal Process. Adv. in Wireless Commun. (SPAWC), Jun. 2018, pp. 1–5.
 [2] I. F. Akyildiz, D. Pompili, and T. Melodia, “Underwater acoustic sensor networks: research challenges,” Ad Hoc Networks, vol. 3, no. 3, pp. 257 – 279, May 2005.
 [3] I. F. Akyildiz, “Wireless sensor networks in challenged environments such as underwater and underground,” in Proc. of the 17th ACM Int. Conf. on Modeling, Ana. and Sim. of Wireless and Mobile Systems, 2014, pp. 1–2.
 [4] A. Celik, N. Saeed, T. Y. AlNaffouri, and M.S. Alouini, “Modeling and performance analysis of multihop underwater optical wireless sensor networks,” in IEEE Wireless Commun. and Netw. Conf., (WCNC), Apr. 2018, pp. 1–6.
 [5] Z. Zeng, S. Fu, H. Zhang, Y. Dong, and J. Cheng, “A survey of underwater optical wireless communications,” IEEE Commun. Surveys Tuts., vol. 19, no. 1, pp. 204–238, Oct. 2016.
 [6] F. Akhoundi, A. Minoofar, and J. A. Salehi, “Underwater positioning system based on cellular underwater wireless optical CDMA networks,” in Wireless and Optical Commun. Conf., (WOCC), Apr. 2017, pp. 1–3.
 [7] J.F. Cai, E. J. CandÃ¨s, and Z. Shen, “A singular value thresholding algorithm for matrix completion,” SIAM J. on Optimization, vol. 20, no. 4, pp. 1956–1982, Mar. 2010.
 [8] K. Lee and Y. Bresler, “ADMiRA: Atomic Decomposition for Minimum Rank Approximation,” ArXiv eprints, Apr. 2009.
 [9] R. H. Keshavan, A. Montanari, and S. Oh, “Matrix completion from a few entries,” IEEE Trans. Inf. Theory, vol. 56, no. 6, pp. 2980–2998, Jun. 2010.
 [10] Z. Lin, M. Chen, and Y. Ma, “The augmented lagrange multiplier method for exact recovery of corrupted lowrank matrices,” arXiv preprint arXiv:1009.5055, Oct. 2013.
 [11] C. Wang, C. Li, and J. Wang, “A modified augmented lagrange multiplier algorithm for toeplitz matrix completion,” Adv. in Comput. Math., vol. 42, no. 5, pp. 1209–1224, Oct. 2016.
 [12] D. Gamarnik and S. Misra, “A note on alternating minimization algorithm for the matrix completion problem,” IEEE Signal Process. Lett., vol. 23, no. 10, pp. 1340–1343, Oct. 2016.
 [13] H. Zhang, “Twodimensional optimal sensor placement,” IEEE Trans. Syst., Man, Cybern., vol. 25, no. 5, pp. 781–792, May 1995.
 [14] A. N. Bishop, B. Fidan, B. D. O. Anderson, K. Dogancay, and P. N. Pathirana, “Optimality analysis of sensortarget geometries in passive localization: Part 1  bearingonly localization,” in 3rd Int. Conf. on Intelligent Sensors, Sensor Networks and Information, Dec. 2007, pp. 7–12.
 [15] A. N. Bishop, B. Fidan, B. D. O. Anderson, P. N. Pathirana, and K. Dogancay, “Optimality analysis of sensortarget geometries in passive localization: Part 2 timeofarrival based localization,” in 3rd Int. Conf. on Intelligent Sensors, Sensor Networks and Information, Dec. 2007, pp. 13–18.
 [16] A. N. Bishop and P. Jensfelt, “An optimality analysis of sensortarget geometries for signal strength based localization,” in Int. Conf. on Intelligent Sensors, Sensor Networks and Information Processing, Dec. 2009, pp. 127–132.
 [17] M. Hamdollahzadeh, S. Adelipour, and F. Behnia, “Optimal sensor configuration for two dimensional source localization based on TDOA/FDOA measurements,” in 17th Int. Radar Symposium (IRS), May 2016, pp. 1–6.
 [18] N. Saeed and H. Nam, “Energy efficient localization algorithm with improved accuracy in cognitive radio networks,” IEEE Commun. Lett., vol. 21, no. 9, pp. 2017–2020, Sep. 2017.
 [19] C. Rusu and J. Thompson, “On the use of tight frames for optimal sensor placement in timedifference of arrival localization,” in 25th European Signal Processing. Conf., (EUSIPCO), Aug. 2017, pp. 1415–1419.
 [20] J. PerezRamirez, D. K. Borah, and D. G. Voelz, “Optimal 3D landmark placement for vehicle localization using heterogeneous sensors,” IEEE Trans. Veh. Technol., vol. 62, no. 7, pp. 2987–2999, Sep. 2013.
 [21] N. Saeed, A. Celik, T. Y. AlNaffouri, and M.S. Alouini, “Energy harvesting hybrid acousticoptical underwater wireless sensor networks localization,” Sensors, vol. 18, no. 1, Dec. 2017.
 [22] N. Saeed, A. Celik, T. Y. AlNaffouri, and M.S. Alouini, “Underwater optical sensor networks localization with limited connectivity,” in IEEE Int. Conf. on Acoustics, Speech and Signal Process., (ICASSP), Apr. 2018, pp. 1–5.
 [23] N. Saeed, A. Celik, T. Y. AlNaffouri, and M.S. Alouini, “Underwater Optical Wireless Communications, Networking, and Localization: A Survey,” ArXiv eprints, Feb. 2018.
 [24] N. Srebro, “Learning with matrix factorizations,” Ph.D. dissertation, Massachusetts Institute of Technology (MIT), 2004.
 [25] Y. Shen, Z. Wen, and Y. Zhang, “Augmented lagrangian alternating direction method for matrix separation based on lowrank factorization,” Optimization Methods and Software, vol. 29, no. 2, pp. 239–263, Jul. 2014.
 [26] Z. Wen, W. Yin, and Y. Zhang, “Solving a lowrank factorization model for matrix completion by a nonlinear successive overrelaxation algorithm,” Math. Programm. Comput., vol. 4, no. 4, pp. 333–361, Dec. 2012.
 [27] Y. Hu, D. Zhang, J. Ye, X. Li, and X. He, “Fast and accurate matrix completion via truncated nuclear norm regularization,” IEEE Trans. Pattern Anal. Mach. Intell., vol. 35, no. 9, pp. 2117–2130, Sep. 2013.
 [28] B. Vandereycken, “Lowrank matrix completion by riemannian optimization,” SIAM J. on Optimization, vol. 23, no. 2, pp. 1214–1236, 2013.
 [29] Y. Zhang, N. Meratnia, and P. J. Havinga, “Outlier detection techniques for wireless sensor networks: A survey,” IEEE Commun. Surveys Tuts., vol. 12, no. 2, pp. 159–170, Apr. 2010.
 [30] X. Wang, J. Luo, S. Li, D. Dong, and W. Cheng, “Component based localization in sparse wireless ad hoc and sensor networks,” in IEEE Int. Conf. on Network Protocols, Oct. 2008, pp. 288–297.
 [31] A. Y. Teymorian, W. Cheng, L. Ma, X. Cheng, X. Lu, and Z. Lu, “3D underwater sensor network localization,” IEEE Trans. Mobile Comput., vol. 8, no. 12, pp. 1610–1621, Dec. 2009.
 [32] Y. Guo and Y. Liu, “Localization for anchorfree underwater sensor networks,” Computers & Elect. Eng., vol. 39, no. 6, pp. 1812–1821, 2013.
 [33] Y. Ren, J. Zhong, J. Huang, Y. Song, X. Xin, N. Yu, and R. Feng, “Orthogonal regression based multihop localization algorithm for largescale underwater wireless sensor networks,” Int. J. of Distrib. Sens. Netw., vol. 10, no. 3, p. 596082, 2014.
 [34] S. Zhang, Q. Zhang, M. Liu, and Z. Fan, “A topdown positioning scheme for underwater wireless sensor networks,” Sci. China Inf. Sci., vol. 57, no. 3, pp. 1–10, 2014.
 [35] M. Y. S. Uddin, “Lowoverhead rangebased 3D localization technique for underwater sensor networks,” in IEEE Int. Conf. on Commun., (ICC), 2016, pp. 1–6.
 [36] K. Mridula and P. Ameer, “Localization under anchor node uncertainty for underwater acoustic sensor networks,” Int. J. of Commun. Sys., vol. 31, no. 2, pp. 1–14, 2018.
 [37] X. Xu, “Impulsive control in continuous and discretecontinuous systems,” IEEE Trans. Autom. Control, vol. 48, no. 12, pp. 2288–2289, Dec. 2003.
 [38] J. Ousingsawat and M. E. Campbell, “Optimal cooperative reconnaissance using multiple vehicles,” J. of Guidance, Control, and Dynamics, vol. 30, no. 1, pp. 122–132, Jan. 2007.
 [39] V. I. Haltrin, “Chlorophyll based model of seawater optical properties,” Appl. Opt., vol. 38, no. 33, pp. 6826–6832, Nov. 1999.
 [40] L. K. Gkoura, G. D. Roumelas, H. E. Nistazakis, H. G. Sandalidis, A. Vavoulas, A. D. Tsigopoulos, and G. S. Tombras, “Underwater optical wireless communication systems: A concise review,” in Turbulence Modelling Approaches  Current State, Development Prospects, Applications. InTech, 2017.
 [41] R. M. Corless, G. H. Gonnet, D. E. G. Hare, D. J. Jeffrey, and D. E. Knuth, “On the lambert W function,” Adv. in Comput. Math., vol. 5, no. 1, pp. 329–359, Dec. 1996.
 [42] S. Arnon and D. Kedar, “Nonlineofsight underwater optical wireless communication network,” J. Opt. Soc. Am. A, vol. 26, no. 3, pp. 530–539, Mar. 2009.
 [43] H. M. Oubei, E. Zedini, R. T. ElAfandy, A. Kammoun, M. Abdallah, T. K. Ng, M. Hamdi, M.S. Alouini, and B. S. Ooi, “Simple statistical channel model for weak temperatureinduced turbulence in underwater wireless optical communication systems,” Opt. Lett., vol. 42, no. 13, pp. 2455–2458, Jul. 2017.
 [44] M. V. Jamali, P. Khorramshahi, A. Tashakori, A. Chizari, S. Shahsavari, S. AbdollahRamezani, M. Fazelian, S. Bahrani, and J. A. Salehi, “Statistical distribution of intensity fluctuations for underwater wireless optical channels in the presence of air bubbles,” in Iran Workshop on Commun. and Info. Theory, (IWCIT), May 2016, pp. 1–6.
 [45] M. V. Jamali, A. Mirani, A. Parsay, B. Abolhassani, P. Nabavi, A. Chizari, P. Khorramshahi, S. Abdollahramezani, and J. A. Salehi, “Statistical studies of fading in underwater wireless optical channels in the presence of air bubble, temperature, and salinity random variations,” IEEE Trans. Commun., pp. 1–16, May 2018, Early Access.
 [46] B. M. C., HeavyTailed Distributions. American Cancer Society, 2014.
 [47] Y. Cao and Y. Xie, “Poisson matrix completion,” in IEEE Int. Symp. on Info. Theory, (ISIT), Jun. 2015, pp. 1841–1845.
 [48] J. W. Sammon, “A nonlinear mapping for data structure analysis,” IEEE Trans. on Computers, vol. 18, no. 5, pp. 401–409, May 1969.
 [49] V. E. McGee, “The multidimensional analysis of elastic distances,” British J. of Math. and Stat. Psy., vol. 19, no. 2, pp. 181–196, 1966.
 [50] A. Buja, D. F. Swayne, M. L. Littman, N. Dean, H. Hofmann, and L. Chen, “Data visualization with multidimensional scaling,” J. of Comput. and Graphical Stat., vol. 17, no. 2, pp. 444–472, Dec. 2008.
 [51] M. Udell, C. Horn, R. Zadeh, and S. Boyd, “Generalized low rank models,” Found. Trends Mach. Learn., vol. 9, no. 1, pp. 1–118, Jun. 2016.
 [52] P. J. Huber and E. M. Ronchetti, Robust Statistics. New York: Wiley, 2009.
 [53] Y. Aflalo and R. Kimmel, “Spectral multidimensional scaling,” Proc. of the National Academy of Sci., vol. 110, no. 45, pp. 18 052–18 057, 2013.
 [54] R. Tibshirani, “Regression shrinkage and selection via the lasso,” J. of the Roy. Stat. Soc. Series B (Methodol.), vol. 58, no. 1, pp. 267–288, 1996.
 [55] I. Markovsky, Low Rank Approximation Algorithms, Implementation, Applications. New York: Springer International Publishing, 2012.
 [56] G. H. Golub and C. F. Van Loan, Matrix computations. JHU Press, 2012, vol. 3.
 [57] T. Zhang, “Multistage convex relaxation for learning with sparse regularization,” in Adv. in Neural Info. Processing Sys., 2009, pp. 1929–1936.
 [58] R. He, B. Hu, X. Yuan, and L. Wang, MEstimators and HalfQuadratic Minimization. Cham: Springer International Publishing, 2014, pp. 3–11.
 [59] C. Goodall, “Procrustes methods in the statistical analysis of shape,” J. R. Statist. Soc. B, vol. 53, no. 2, pp. 285–339, Oct. 1991.
 [60] J. N. S. Wright, Numerical Optimization. NY, USA: SpringerVerlag, 2006.
Nasir Saeed(S’14M’2016) received his Bachelors of Telecommunication degree from N.W.F.P University of Engineering and Technology, Peshawar, Pakistan, in 2009 and received Masters degree in satellite navigation from Polito di Torino, Italy, in 2012. He received his Ph.D. degree in electronics and communication engineering from Hanyang University, Seoul, South Korea in 2015. He was an assistant professor at the Department of Electrical Engineering, Gandhara Institute of Science and IT, Peshawar, Pakistan from August 2015 to September 2016. Dr. Saeed worked as an assistant professor at IQRA National University, Peshawar, Pakistan from October 2017 to July 2017. He is currently a postdoctoral research fellow at Communication Theory Lab, King Abdullah University of Science and Technology (KAUST). His current areas of interest include cognitive radio networks, underwater optical sensor networks, localization, and random matrix theory. 
Tareq Y. AlNaffouri (M’10) Tareq AlNaffouri received the B.S. degrees in mathematics and electrical engineering (with first honors) from King Fahd University of Petroleum and Minerals, Dhahran, Saudi Arabia, the M.S. degree in electrical engineering from the Georgia Institute of Technology, Atlanta, in 1998, and the Ph.D. degree in electrical engineering from Stanford University, Stanford, CA, in 2004. He was a visiting scholar at California Institute of Technology, Pasadena, CA, from January to August 2005 and during summer 2006. He was a Fulbright Scholar at the University of Southern California from February to September 2008. He has held internship positions at NEC Research Labs, Tokyo, Japan, in 1998, Adaptive Systems Lab, University of California at Los Angeles in 1999, National Semiconductor, Santa Clara, CA, in 2001 and 2002, and Beceem Communications Santa Clara, CA, in 2004. He is currently an Associate professor at the Electrical Engineering Department, King Abdullah University of Science and Technology (KAUST). His research interests lie in the areas of sparse, adaptive, and statistical signal processing and their applications and in network information theory. He has over 150 publications in journal and conference proceedings, 9 standard contributions, 10 issued patents, and 6 pending. Dr. AlNaffouri is the recipient of the IEEE Education Society Chapter Achievement Award in 2008 and AlMarai Award for innovative research in communication in 2009. Dr. AlNaffouri has also been serving as an Associate Editor of Transactions on Signal Processing since August 2013. 
MohamedSlim Alouini (S’94M’98SM’03F’09) was born in Tunis, Tunisia. He received the Ph.D. degree in Electrical Engineering from the California Institute of Technology (Caltech), Pasadena, CA, USA, in 1998. He served as a faculty member in the University of Minnesota, Minneapolis, MN, USA, then in the Texas A&M University at Qatar, Education City, Doha, Qatar before joining King Abdullah University of Science and Technology (KAUST), Thuwal, Makkah Province, Saudi Arabia as a Professor of Electrical Engineering in 2009. His current research interests include the modeling, design, and performance analysis of wireless communication systems. 