# Real-Time Nonlinear Characterization of Soft Tissue Mechanical Properties

Authors: Jaehyun Shin; Yongmin Zhong; Chengfan Gu
Journal: Journal of Sensors
 (2020)

---

## Abstract

Online soft tissue characterization is important for robotic-assisted minimally invasive surgery to achieve precise and stable robotic control with haptic feedback. This paper presents a new nonlinear recursive adaptive filtering methodology for online nonlinear soft tissue characterization. An adaptive unscented Kalman filter is developed based on the Hunt-Crossley model by windowing approximation to online estimate system and measurement noise covariances. To improve the accuracy of noise covariance estimations, a recursive formulation is subsequently developed for estimation of system and measurement noise covariances by introducing a weighting factor. This weighting factor is further modified to accommodate noise statistics of large variation which could be caused by rupture events and geometric discontinuities in robotic-assisted surgery. Simulations, experiments, and comparison analyses demonstrate that the proposed nonlinear recursive adaptive filtering methodology can characterize soft tissue parameters in the presence of system or measurement noise statistics in both small and large variations for robotic-assisted surgery. The proposed methodology can effectively estimate soft tissue parameters under system and measurement noises in both small and large variations, leading to improved filtering accuracy and robustness in comparison with UKF.


---

## Body

## 1. Introduction

Soft tissue properties are of great importance in robotic minimally invasive surgery to characterize the interaction between surgical tools and soft tissues for robotic control with force feedback. However, soft tissue properties are time-dependent and patient-specific. They are variable as per various tissue layers and regions, various organs, and various physiological conditions. In order to improve the robustness and stability of haptic control, it is necessary to characterize soft tissue properties dynamically [1].A contact model describing the contact interaction between surgical tools and soft tissues is essential to characterization of soft tissue properties. The existing contact models can be generally classified into two categories. One is the continuum mechanics approach focusing on accurate characterization of the tool-tissue contact interaction. The typical examples in this category include the viscoelastic model [2] and finite element model (FEM) [3], where the mechanical behaviours of a soft tissue are accurately characterized based on continuum mechanics of elasticity. However, this approach is very complex and involves a large amount of computational load. Thus, it is only suitable for offline analysis, while unable to satisfy the real-time requirement of soft tissue characterization. The other is the analytical approach focusing on analytical establishment of the explicit relationship between force and displacement using a spring and damper. In spite of the less accuracy compared to the continuum mechanics approach, the analytical approach has the advantages of real-time performance, suitable for a robotic control purpose. Further, since the spring and damper represent the elasticity and viscosity of soft tissues, the approach still maintains certain accuracy for characterization of soft tissue behaviours. Thus, the existing studies are mainly dominated by the analytical approach, leading to various spring-damper models for online soft tissue characterization. The simple one is made up of a spring, unable to model the viscosity effect [4]. The Maxwell model consists of a serialized spring and damper [5]. The Kelvin-Voigt model consists of a parallel spring and damper [6]. The Kelvin-Boltzmann (KB) model improves the K-V model adding a serialized spring [7]. In general, the above models are a linear model, involving physical inconsistencies at the initial and end stages of contact with soft tissues [4]. Comparing to the above models, the H-C model is a nonlinear contact model and shows consistencies with physical contact. However, the use of such a nonlinear model requires a nonlinear estimation algorithm, which is a challenging research task. Accordingly, the use of the H-C model for accurate and real-time characterization of soft tissue properties has been very limited. A survey of the existing contact models can be found in [8, 9].In order to achieve the real-time performance for soft tissue characterization, it also needs to develop an online estimation algorithm. So far, the existing studies are mainly dominated by a linear estimation algorithm such as the recursive least square (RLS) [10, 11] and Kalman filter [8]. However, a linear estimation algorithm is only suitable for linear contact models. In order to use a linear estimation algorithm with the nonlinear H-C model, it requires a linearization process. Since the linearization process of the nonlinear H-C model causes error for estimation, the resultant solution may be biased or even divergent [12]. Further, the Kalman filter also requires prior knowledge of system noise statistics. A survey on estimation and tuning of system noise statistics for Kalman filtering can be found in [13, 14].The unscented Kalman filter (UKF) is an online nonlinear estimation algorithm, where the state mean and covariance are approximated in third-order accuracy using unscented transformation [15]. However, similar to the Kalman filter, the performance of UKF relies on accurate contact and measurement models. Currently, there have been few studies focusing on development of a UKF algorithm for characterization of nonlinear soft tissue properties. Xi et al. studied a reduced-order UKF based on a cubic Hermite FEM for characterization of myocardial stiffness [3]. However, this FEM is subject to the quasistatic assumption, which does not account for the contact dynamics. The reduced dimension of the state vector also deteriorates the estimation performance. The authors also studied a UKF based on the H-C model for soft tissue characterization [16], but without considering the UKF requirement of accurate system models. Just recently, the authors reported an improved UKF to handle the error of the H-C contact model using the scaling factor [17]. However, this method is based on the assumption that the measurement model is accurate, which does not correspond to the real-world situation with errors involved in both contact and measurement models. Further, it is also based on windowing approximation [18, 19] to use residuals within a small window of historical time points to calculate system noise characteristics at the current time point. In spite of the fast response to dynamic model noise, the estimation accuracy is limited due to the restricted data in the small time window. Recursion [20, 21] is a strategy to assess system characteristics at the present time point using all historical data, leading to higher accuracy than windowing approximation. Since all data are assembled via a recursive process, rather than stored in memory, this strategy also has a small computational load. However, compared to windowing approximation, the recursion strategy is slow in response to model noise due to the use of all available data. As windowing approximation and recursion complement with each other, it is necessary to study how to effectively combine both strategies together to fully leverage their individual advantages for improvement of the UKF performance.This paper presents a new nonlinear recursive adaptive UKF to estimate the parameters of the H-C contact model. It combines windowing approximation with recursion to address the UKF limitation in requirement of accurate contact and measurement models. An adaptive UKF is developed based on windowing approximation to offer the classical UKF with the capability to online estimate the noise covariances of the H-C contact and measurement models. Subsequently, a recursive adaptive UKF is developed by incorporating a weighting factor in the adaptive filtering process to improve the estimation accuracy. This weighting factor is further modified to accommodate the noise statistics of large variation, which are occurring in robotic-assisted surgery due to rupture events and geometric discontinuities. Simulation and experimental results together with comparison analysis demonstrate the efficiency of the proposed method for online soft tissue characterization.


## 2. Hunt-Crossley Model

The Hunt-Crossley (H-C) model describes the dynamic mechanism of mechanical contact between soft tissues and surgical tools with a nonlinear spring and damper [22]:
(1)F=Kdn+Bdnḋp,where F, K, B, d, ḋ, n, and p denote the interaction force, stiffness, damping coefficient, displacement, velocity, and power indices of displacement and velocity, respectively.From (1), the system state equation is defined as
(2)xk=fxk−1+qk−1=dk−1+ḋk−1×Δtḋk−1Kk−1dk−1nk−1+Bk−1dk−1nk−1ḋk−1pk−1Kk−1Bk−1nk−1pk−1+qk−1,where xk=dk,ḋk,Fk,Kk,Bk,nk,pk is the system state at time tk, qk~0,Qk is the system noise which is assumed as a white Gaussian noise with zero mean and covariance Qk, f· is the system function, and Δt is the time step.The measurement function is defined as(3)yk=Hkxk+rk=dkFk+rk,where yk denotes the measurement vector, Hk is the measurement matrix, and rk~0,Rk is the measurement noise which is also assumed as a white Gaussian noise with zero mean and covariance Rk.


## 3. Analysis of Conventional UKF

Let us briefly analyse the problems of the classical UKF. The conventional UKF procedure is as follows:(i)
Initialization.Calculate the initial estimated state and its associated covariance:(4)x0=Ext0,P0=Ex0−xt0x0−xt0T,where xt0 is the predefined initial state.
(ii)
Time update.Select the sigma points of the estimated state at timetk−1 (k=1,2,⋯):
(5)xk−1i=xk−1,i=0,xk−1i=xk−1+aNPk−1i,i=1,⋯,N,xk−1i=xk−1−aNPk−1iT,i=N+1,⋯,2N,where N is the dimension of state vector x, a is the tuning parameter to control the spread of the sigma points around xk−1, and NPk−1 is the ith column of the mean square root of NPk−1.Transform the sigma points with the state function, yielding a set of transformed samples:(6)xk∣k−1i=fxk−1i.Calculate predicted statexk∣k−1 and its associated covariance Pk∣k−1:
(7)xk∣k−1=∑i=12Nwixk∣k−1i,(8)Pk∣k−1=∑i=12Nwixk∣k−1i−xk∣k−1xki−xk∣k−1T+Qk,where
(9)wi=1−1a2,i=0,wi=12Na2,i=1,2,⋯,2N.
(iii)
Measurement update.As the measurement model is linear, the measurement update is performed with the same equations as the Kalman filter.Calculate predicted measurementyk∣k−1 and its associated covariance Pyk∣k−1:
(10)yk∣k−1=Hkxk∣k−1,(11)Pyk∣k−1=HkPk∣k−1HkT+Rk.Calculate the crosscovariance between predicted statexk∣k−1 and predicted measurement yk∣k−1:
(12)Pxk∣k−1yk∣k−1=Pk∣k−1HkT.Determine the Kalman gain:(13)Kk=Pxk∣k−1yk∣k−1Pyk∣k−1−1.Update estimated statexk and its associated covariance Pk:
(14)xk=xk∣k−1+Kkyk−yk∣k−1,(15)Pk=Pk∣k−1−KkPyk∣k−1KkT.
(iv)
Go to step (i) for the next sample until all samples are processedIt can be seen from (8) that if the system noise covariance involves uncertainty, predicted state covariance Pk∣k−1 will become inaccurate, leading crosscovariance Pxk∣k−1yk∣k−1 given by (12) to become inaccurate. The inaccurate Pxk∣k−1yk∣k−1 will further lead Kalman gain Kk to become biased, thus degenerating the state estimate given by (15). In a similar way, if the measurement noise covariance involves uncertainty, predicted measurement covariance Pyk∣k−1 will become inaccurate, leading Kalman gain Kk to be biased. Accordingly, the state estimate given by (14) will be corroded. From above, it is evident that without the accurate system noise covariances, the filtering solution of the conventional UKF will be biased or even divergent.


## 4. Recursive Adaptive Unscented Kalman Filter

In this paper, a new recursive adaptive UKF is developed by combining windowing approximation with recursion to the online estimate system and measurement noise covariances via covariance matching, which is to balance predicted state covariancePk∣k−1 and predicted measurement covariance Pyk∣k−1 with their theoretical values.


### 4.1. Estimation of System Noise Covariance

By introducing an adaptive scaling factorΦ, we readily have the following relationship:
(16)Qk=ΦQk−1.Substituting (16) into (8), the predicted state covariance can be further represented as
(17)Pk∣k−1=∑i=12Nwixk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1T+ΦQk−1.It can be seen from (17) that the adaptive scaling factor Φ enables us to adjust the system noise covariance for online tuning the predicted state covariance.In theory, the predicted state covariance is defined as(18)P~k∣k−1=Ex^k−xk∣k−1x^k−xk∣k−1T.Since the real statex^k is generally unknown, in (18) replacing the real state with the estimated state and applying windowing approximation to calculate the expectation, we readily have a suboptimal estimation for the predicted state covariance:
(19)P~k∣k−1=1m∑j=1mxk−j−xk−j∣k−j−1xk−j−xk−j∣k−j−1T,where m represents the number of the state vectors within the time window (tk−m, tk−1) and is named the window size.Define the innovation vector by(20)Zk∣k−1=yk−yk∣k−1.By (14), we readily have
(21)xk−xk∣k−1=KkZk∣k−1.Substituting (21) into (19) yields
(22)P~k∣k−1=1m∑j=1mKk−jZk−j∣k−j−1Zk−j∣k−j−1TKk−jT.Based on the above, the adaptive scaling factorΦ can be attained by matching the predicted state covariance given by (17) with its theoretical value given by (22):
(23)∑i=12Nwixk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1T+ΦQk−1=P~k∣k−1.Thus, the adaptive scaling factor for the system noise covariance can be determined as(24)Φ=trP~k∣k−1−tr∑i=12Nwixk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1TtrQk−1,where tr· denotes the trace of a matrix.It can be seen from (22) that P~k∣k−1 is approximated using the information within the current window, and thus, it can account for dynamically changing conditions. However, due to the limited data available within the time window, the estimation accuracy is limited. In order to address this problem, the theoretical predicted state covariance is further expressed into a recursive form:
(25)P~k∣k−1∗=CkP~k∣k−1+1−CkP~k−1∣k−2∗,where Ck is called the weighting factor, which is commonly defined by [23, 24]:
(26)Ck=1Jk,where Jk is the number of the entire available time points.It can be seen from (25) that P~k∣k−1∗ is based on the entire available data rather than the limited data available within the current time window, leading to improved accuracy. P~k∣k−1∗ is contributed by two portions. One is at the current time point, which is calculated by windowing approximation to use the data available within the current window to characterize dynamically changing conditions that occurred. The other is at the previous time point, which reflects the historical tendency. The contributions at both current and previous time points are combined via weighting factor Ck. P~k∣k−1∗ takes the full contribution at the current time point when Ck=1, while taking the full contribution at the previous time point when Ck=0. The smaller the value of Ck is, the smaller the contribution at the current time point is and in turn the larger the contribution at the previous time point is.


### 4.2. Modified Weighting Factor

In robotic-assisted minimally invasive surgery, the variation of noise statistics is small when the dynamic interaction environment is relatively stable [25]. However, in case of abrupt changes such as geometrical discontinuities in soft tissue deformation, rupture events occurring in the penetration of different tissue layers, and unexpected contacts with rigid bones, a large variation in noise statistics will occur in the dynamic interaction environment [26].The standard weighting factorCk defined by (26) cannot account for abrupt changes encountered in the interaction with soft tissues. Figure 1 plots the variation of the standard weighting factor Ck in time series where an abrupt change in noise statistics has occurred at tk=50 s. As shown in Figure 1(a), Ck given by (26) starts from the maximum value of one at the beginning. Subsequently, it decreases drastically to a small value and gradually converges to zero in the time series. As Ck remains at such a small value most of the time, the contribution at the present time point to P~k∣k−1∗ is small while in turn the contribution at the previous time point is large. Due to the small contribution at the present time point most of the time, P~k∣k−1∗ is unable to account for an abrupt change in system noise statistics. As shown in Figure 1(a), as Ck still remains at a small value at time point tk=50 s, it leads to the failure to account for the abrupt change of noise statistics occurring at this time point. Accordingly, Ck given by (26) is only suitable for the case that system noise statistics are constant or have a small variation, making the proposed adaptive UKF unable to account for a large variation in system noise statistics.Figure 1
The effects of both standard and modified weighting factors onP~k∗, where an abrupt change occurs at tk=50 s: the black solid line represents the proportion of P~k at the current time point while the red dashed line the proportion of P~k−1∗ at the previous time point.
(a)
Standard weighting factor
(b)
Modified weighting factorIn order to accommodate large variations in system noise statistics, the standard weighting factorCk defined by (26) is modified as
(27)CkR=1Jk−b,(28)b=kr−1gk−kr,where g· is a unit step function and kr is the index of the time point where an abrupt change occurs.The initial value of parameterb is zero. When an abrupt change occurs at the current time point, the parameter b will be updated to kr−1, making CkR become one again. Consequently, P~k∣k−1∗ is calculated by taking the full contribution at the present time point to account for the abrupt change. Figure 1(b) shows the variation of the modified weighting factor CkR where an abrupt change in system noise statistics occurred at tk=50 s. Similar to Ck, CkR starts from the maximum value of one, where P~k∣k−1∗ takes the full contribution at the present time point. Subsequently, the CkR value drops drastically to and remains at a small value until an abrupt change occurs. With such a small value of CkR, P~k∣k−1∗ mainly follows the one at the previous time point to handle small variations in system noise covariance. When an abrupt change occurs at the time point k=50 s (kr=50 s), the parameter b defined by (28) will be updated to kr−1, restoring CkR back to the maximum value of one to account for the abrupt change that occurred. Subsequently, CkR drastically drops to and remains at a small value until a new abrupt change occurs. This variation process of CkR is repeated whenever an abrupt change occurs, allowing the adjustment of P~k∣k−1∗ to account for an abrupt change by taking the full contribution at the present time point.By substituting the weighting factor with the modified weighting factor, (25) becomes
(29)P~k∣k−1∗=CkRP~k∣k−1+1−CkRP~k−1∣k−2∗.Thus, the adaptive scaling factor is calculated by(30)Φ=trP~k∣k−1∗−tr∑i=12Nwicxk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1TtrQk−1.


### 4.3. Estimation of Measurement Noise Covariance

The measurement noise covariance can also be estimated in the similar way as the system noise covariance. By introducing an adaptive scaling factorΓ, the measurement noise covariance can be represented as
(31)Rk=ΓRk−1.Substituting (31) into (11), the predicted measurement covariance can be further expressed as
(32)Pyk∣k−1=HkPk∣k−1HkT+ΓRk−1.It can be seen from (32) that the adaptive scaling factor Γ enables us to adjust the measurement noise covariance for online tuning the predicted measurement covariance.Using windowing approximation [27], the theoretical predicted measurement covariance can be directly calculated as
(33)P~yk∣k−1=1m∑j=1mZk−j∣k−j−1Zk−j∣k−j−1T.Similar to (29), the theoretical predicted measurement covariance P~yk∣k−1 can be further expressed by a recursive form as
(34)P~yk∣k−1∗=CkRP~yk∣k−1+1−CkRP~yk−1∣k−2,where the modified weighting factor CkR is also applied to handle both small and large variations in the measurement noise covariance.By matching the modified predicted measurement covariance given by (32) and its theoretical value given by (34), we readily have
(35)HkPk∣k−1HkT+ΓRk−1=P~yk∣k−1∗.From (35), the adaptive scaling factor Γ for estimation of measurement noise covariance is determined as
(36)Γ=trP~yk∣k−1∗−trHkPk∣k−1HkTtrRk−1.


## 4.1. Estimation of System Noise Covariance

By introducing an adaptive scaling factorΦ, we readily have the following relationship:
(16)Qk=ΦQk−1.Substituting (16) into (8), the predicted state covariance can be further represented as
(17)Pk∣k−1=∑i=12Nwixk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1T+ΦQk−1.It can be seen from (17) that the adaptive scaling factor Φ enables us to adjust the system noise covariance for online tuning the predicted state covariance.In theory, the predicted state covariance is defined as(18)P~k∣k−1=Ex^k−xk∣k−1x^k−xk∣k−1T.Since the real statex^k is generally unknown, in (18) replacing the real state with the estimated state and applying windowing approximation to calculate the expectation, we readily have a suboptimal estimation for the predicted state covariance:
(19)P~k∣k−1=1m∑j=1mxk−j−xk−j∣k−j−1xk−j−xk−j∣k−j−1T,where m represents the number of the state vectors within the time window (tk−m, tk−1) and is named the window size.Define the innovation vector by(20)Zk∣k−1=yk−yk∣k−1.By (14), we readily have
(21)xk−xk∣k−1=KkZk∣k−1.Substituting (21) into (19) yields
(22)P~k∣k−1=1m∑j=1mKk−jZk−j∣k−j−1Zk−j∣k−j−1TKk−jT.Based on the above, the adaptive scaling factorΦ can be attained by matching the predicted state covariance given by (17) with its theoretical value given by (22):
(23)∑i=12Nwixk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1T+ΦQk−1=P~k∣k−1.Thus, the adaptive scaling factor for the system noise covariance can be determined as(24)Φ=trP~k∣k−1−tr∑i=12Nwixk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1TtrQk−1,where tr· denotes the trace of a matrix.It can be seen from (22) that P~k∣k−1 is approximated using the information within the current window, and thus, it can account for dynamically changing conditions. However, due to the limited data available within the time window, the estimation accuracy is limited. In order to address this problem, the theoretical predicted state covariance is further expressed into a recursive form:
(25)P~k∣k−1∗=CkP~k∣k−1+1−CkP~k−1∣k−2∗,where Ck is called the weighting factor, which is commonly defined by [23, 24]:
(26)Ck=1Jk,where Jk is the number of the entire available time points.It can be seen from (25) that P~k∣k−1∗ is based on the entire available data rather than the limited data available within the current time window, leading to improved accuracy. P~k∣k−1∗ is contributed by two portions. One is at the current time point, which is calculated by windowing approximation to use the data available within the current window to characterize dynamically changing conditions that occurred. The other is at the previous time point, which reflects the historical tendency. The contributions at both current and previous time points are combined via weighting factor Ck. P~k∣k−1∗ takes the full contribution at the current time point when Ck=1, while taking the full contribution at the previous time point when Ck=0. The smaller the value of Ck is, the smaller the contribution at the current time point is and in turn the larger the contribution at the previous time point is.


## 4.2. Modified Weighting Factor

In robotic-assisted minimally invasive surgery, the variation of noise statistics is small when the dynamic interaction environment is relatively stable [25]. However, in case of abrupt changes such as geometrical discontinuities in soft tissue deformation, rupture events occurring in the penetration of different tissue layers, and unexpected contacts with rigid bones, a large variation in noise statistics will occur in the dynamic interaction environment [26].The standard weighting factorCk defined by (26) cannot account for abrupt changes encountered in the interaction with soft tissues. Figure 1 plots the variation of the standard weighting factor Ck in time series where an abrupt change in noise statistics has occurred at tk=50 s. As shown in Figure 1(a), Ck given by (26) starts from the maximum value of one at the beginning. Subsequently, it decreases drastically to a small value and gradually converges to zero in the time series. As Ck remains at such a small value most of the time, the contribution at the present time point to P~k∣k−1∗ is small while in turn the contribution at the previous time point is large. Due to the small contribution at the present time point most of the time, P~k∣k−1∗ is unable to account for an abrupt change in system noise statistics. As shown in Figure 1(a), as Ck still remains at a small value at time point tk=50 s, it leads to the failure to account for the abrupt change of noise statistics occurring at this time point. Accordingly, Ck given by (26) is only suitable for the case that system noise statistics are constant or have a small variation, making the proposed adaptive UKF unable to account for a large variation in system noise statistics.Figure 1
The effects of both standard and modified weighting factors onP~k∗, where an abrupt change occurs at tk=50 s: the black solid line represents the proportion of P~k at the current time point while the red dashed line the proportion of P~k−1∗ at the previous time point.
(a)
Standard weighting factor
(b)
Modified weighting factorIn order to accommodate large variations in system noise statistics, the standard weighting factorCk defined by (26) is modified as
(27)CkR=1Jk−b,(28)b=kr−1gk−kr,where g· is a unit step function and kr is the index of the time point where an abrupt change occurs.The initial value of parameterb is zero. When an abrupt change occurs at the current time point, the parameter b will be updated to kr−1, making CkR become one again. Consequently, P~k∣k−1∗ is calculated by taking the full contribution at the present time point to account for the abrupt change. Figure 1(b) shows the variation of the modified weighting factor CkR where an abrupt change in system noise statistics occurred at tk=50 s. Similar to Ck, CkR starts from the maximum value of one, where P~k∣k−1∗ takes the full contribution at the present time point. Subsequently, the CkR value drops drastically to and remains at a small value until an abrupt change occurs. With such a small value of CkR, P~k∣k−1∗ mainly follows the one at the previous time point to handle small variations in system noise covariance. When an abrupt change occurs at the time point k=50 s (kr=50 s), the parameter b defined by (28) will be updated to kr−1, restoring CkR back to the maximum value of one to account for the abrupt change that occurred. Subsequently, CkR drastically drops to and remains at a small value until a new abrupt change occurs. This variation process of CkR is repeated whenever an abrupt change occurs, allowing the adjustment of P~k∣k−1∗ to account for an abrupt change by taking the full contribution at the present time point.By substituting the weighting factor with the modified weighting factor, (25) becomes
(29)P~k∣k−1∗=CkRP~k∣k−1+1−CkRP~k−1∣k−2∗.Thus, the adaptive scaling factor is calculated by(30)Φ=trP~k∣k−1∗−tr∑i=12Nwicxk∣k−1i−xk∣k−1xk∣k−1i−xk∣k−1TtrQk−1.


## 4.3. Estimation of Measurement Noise Covariance

The measurement noise covariance can also be estimated in the similar way as the system noise covariance. By introducing an adaptive scaling factorΓ, the measurement noise covariance can be represented as
(31)Rk=ΓRk−1.Substituting (31) into (11), the predicted measurement covariance can be further expressed as
(32)Pyk∣k−1=HkPk∣k−1HkT+ΓRk−1.It can be seen from (32) that the adaptive scaling factor Γ enables us to adjust the measurement noise covariance for online tuning the predicted measurement covariance.Using windowing approximation [27], the theoretical predicted measurement covariance can be directly calculated as
(33)P~yk∣k−1=1m∑j=1mZk−j∣k−j−1Zk−j∣k−j−1T.Similar to (29), the theoretical predicted measurement covariance P~yk∣k−1 can be further expressed by a recursive form as
(34)P~yk∣k−1∗=CkRP~yk∣k−1+1−CkRP~yk−1∣k−2,where the modified weighting factor CkR is also applied to handle both small and large variations in the measurement noise covariance.By matching the modified predicted measurement covariance given by (32) and its theoretical value given by (34), we readily have
(35)HkPk∣k−1HkT+ΓRk−1=P~yk∣k−1∗.From (35), the adaptive scaling factor Γ for estimation of measurement noise covariance is determined as
(36)Γ=trP~yk∣k−1∗−trHkPk∣k−1HkTtrRk−1.


## 5. Performance Evaluation and Discussions

A prototype system has been implemented using the proposed methodology for online soft tissue characterization. This system takes measurements of contact force and displacement as input signals to online estimate the nonlinear parameters of the H-C model. Based on the estimated parameters of the H-C model, the contact forces with soft tissues are reconstructed from displacements and further compared with the input forces to calculate the estimation error. Simulations and experiments together with comparison analysis have been conducted to evaluate the performance of the proposed methodology in terms of (i) system noise covariance estimation, (ii) measurement noise covariance estimation, (iii) the effect of the weighting factor, and (iv) force reconstruction based on the estimated parameters of the H-C model.


### 5.1. Simulations and Analysis

Simulations were conducted to evaluate the performance of the proposed methodology in estimating system and measurement noise covariances together with associated force reconstructions. The time steps for the simulations were set to 10 ms. The initial state was defined asx0=0.1,0.1,0.01,0.0001,0.0001,1,1.


#### 5.1.1. System Noise Covariance Estimation

For the performance evaluation in terms of system noise covariance estimation, it is supposed that the measurement noise statistics are exactly known. Without loss of generality, setRk=1 and Erk=0. Two cases are studied. One is that system noise covariance is constant or has a small variation, and the other is that system noise covariance involves abrupt changes.(1) Constant or Small-Variation System Noise Covariance. Simulation trials were conducted under the same conditions by the adaptive UKF (AUKF) and recursive adaptive UKF (RAUKF) for the case that system noise covariance is constant or has a small variation. The measurement data were obtained from the H-C model by adding a white Gaussian noise with zero mean and the covariance of 0.036 N. The initial estimation value was Q^0=0.004 N, the window size was m=4, the real covariance of system noise was Qk=0.036 N, and Eqk=0. The simulation time was 90 s. As shown in Figure 2, the system noise covariance estimated by AUKF via windowing approximation involves large oscillations during the entire simulation time. This is because the limited data available within the small window restricts the estimation accuracy. In contrast, the system noise covariance estimated by RAUKF via recursion is close to the real value after 20 s. This is because the estimation via recursion is obtained based on the entire available data. The RMSE (Root Mean Square Error) of AUKF is 0.0125 N, while after 20 s, the RMSE of RAUKF is 0.0032 N. Thus, it is evident that the proposed RAUKF can significantly improve the AUKF accuracy for estimation of system noise covariance.Figure 2
Estimations of system noise covariance by both AUKF and RAUKF for the case that system noise covariance is constant or in a small variation: the estimated system noise covariances are indicated by the blue lines and the reference values by the red lines.
(a)
AUKF via windowing approximation
(b)
RAUKF via recursionFigure3 illustrates the estimation errors in terms of force reconstruction by both AUKF and RAUKF, based on the estimated system noise covariances as shown in Figure 2. It is clear that the estimation error of RAUKF is much smaller than that of AUKF. During the entire simulation time, since the system noise covariance estimated by AUKF involves large oscillations (see Figure 2(a)), the resultant estimation error also involves large oscillations. However, for RAUKF, the resultant estimation error becomes small and is further converged after the initial time period of 20 s, following the trend of the estimated system noise covariance. The maximum estimation error is 0.0394 N for RAUKF, but 0.1853 N for AUKF. The mean estimation error of AUKF is 0.0374 N, while that of RAUKF is 0.0053 N. The RMSE of RAUKF is 0.0070 N, while that of AUKF is almost nine times larger, which is 0.0475 N. Table 1 summarizes the estimation errors of both AUKF and RAUKF.Figure 3
Estimation errors in force reconstruction by AUKF and RAUKF corresponding to the estimated system noise covariances as shown in Figure2: the estimation error of AUKF is indicated by the black lines and that of RAUKF by the red lines.Table 1
Estimation errors of both AUKF and RAUKF for the case with constant or small-variation system noise covariance.
Error (N)
AUKF
RAUKF
Mean error
0.0374
0.0053
Max error
0.1853
0.0394
RMSE
0.0475
0.0070(2) System Noise Covariance Involving Abrupt Changes. In order to analyse the effect of the weighting factor on system noise covariance estimation, simulation trials were conducted under the same conditions by RAUKF via both standard and modified weighting factors for the case of system noise covariance involving abrupt changes. The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean. The noise covariance was set to 0.036 N within [0, 36 s] and 0.576 N within (36 s, 70 s], leading to an abrupt change from 0.036 N to 0.576 N at 36 s. The initial estimation value was Q^0=0.004 N, the window size was set to m=14, and the simulation time was 70 s.Figure4 shows the system noise covariances estimated by RAUKF via both standard and modified weighting factors. Although the estimation results with both weighting factors follow the reference value of 0.036 N within [0, 36 s], the estimation under the standard weighting factor does not follow the abrupt change from 0.036 N to 0.576 N at 36 s, leading to the mean error of 0.372 N after 36 s. In contrast, the estimation of system noise covariance under the modified weighting factor rapidly follows the abrupt change from 0.036 N to 0.576 N at 36 s, leading to the mean error of 0.0088 N after 36 s. This demonstrates that with the modified weighting factor, RAUKF can accommodate abrupt changes in system noise covariance.Figure 4
Estimations of system noise covariance by RAUKF with both standard and modified weighting factors for the case that system noise covariance involves abrupt changes: the estimated system noise covariances are indicated by the red lines and the reference values by the blue lines.
(a)
RAUKF with the standard weighting factor
(b)
RAUKF with the modified weighting factorFigure5 shows the estimation errors in force reconstruction based on the estimated system noise covariances as shown in Figure 4. It can be seen that before 36 s, although the RAUKF estimation errors via both standard and modified weighting factors are small and converged to zero, the one via the modified weighting factor is even smaller and has a faster convergence speed. However, after 36 s, due to the biased estimate of the system noise covariance, the RAUKF estimation error via the standard weighting factor involves large oscillations, leading to the RMSE of 0.0749 N. In contrast, the RAUKF estimation error via the modified weighting factor is much smaller and converged to zero, leading to the RMSE of 0.0046 N. This demonstrates that the modified weighting factor enables RAUKF to account for abrupt changes in system noise covariance, leading to improved estimation accuracy. Table 2 summarizes the RAUKF estimation errors via both standard and modified weighting factors.Figure 5
Estimation errors in terms of force reconstruction by RAUKF based on the estimated system noise covariances via both standard and modified weighting factors as shown in Figure4: the estimation error via the standard weighting factor is indicated by the black lines and the one via the modified weighting factor by the red lines.Table 2
RAUKF estimation errors via both standard and modified weighting factors for the case that system noise covariance involves abrupt changes.
Error (N)
Standard weighting factor
Modified weighting factor
Mean error
0.0030
0.0009
Max error
0.4517
0.0389
RMSE
0.0749
0.0046


#### 5.1.2. Measurement Noise Covariance Estimation

In order to evaluate the estimation of measurement noise covariance, it is assumed that the system noise statistics are exactly known. Without loss of generality, chooseQk=1 and Eqk=0. Similar to Section 5.1.1., two cases are studied. One is that the measurement noise covariance is constant or in a small variation, and the other is that the measurement noise covariance involves abrupt changes.(1) Constant or Small-Variation Measurement Noise Covariance. Simulation trials were conducted under the condition that the measurement noise covariance is constant or has a small variation by both AUKF and RAUKF.The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean and the covariance of 0.036 N. The initial estimation value wasR^0=0.004 N. The window size was set to m=4. The simulation time was 70 s. As shown in Figure 6, the estimate of measurement noise covariance by AUKF involves large oscillations during the entire simulation time, while the one by RAUKF is very close to the reference value after 15 s. The RMSE of RAUKF is 0.0051 N and is almost twice smaller than that of AUKF which is 0.0101 N. This demonstrates that the accuracy of RAUKF is significantly higher than that of AUKF for estimation of measurement noise covariance.Figure 6
Estimations of measurement noise covariance by AUKF and RAUKF for the case that measurement noise covariance is constant or small-variation: the estimated measurement noise covariances are indicated by the blue lines and the reference values by the red lines.
(a)
AUKF via windowing approximation
(b)
RAUKF via recursionFigure7 illustrates the estimation errors in force reconstruction by both AUKF and RAUKF based on the estimated measurement noise covariances as shown in Figure 6. It can be seen clearly that the estimation error of AUKF is much larger than that of RAUKF during the entire simulation time. In addition, after 15 s, the estimation error of RAUKF is close to zero, while that of AUKF is slightly smaller than the initial estimation error. This is because the estimated measurement noise covariance by windowing approximation is fluctuated during the entire simulation time (see Figure 6(a)), while that by RAUKF is converged to the reference value (see Figure 6(b)). The RMSE of RAUKF is 0.0286 N and is two-three times smaller than that of AUKF which is 0.0708 N. Further, the maximum and mean estimation errors of AUKF are about two times larger than those of RAUKF. The maximum estimation error is 0.2905 N for AUKF, while 0.145 N for RAUKF. The mean estimation error is 0.0551 N for AUKF, while 0.0102 N for RAUKF. Therefore, the estimation accuracy of RAUKF is significantly higher than that of AUKF. Table 3 compares the estimation errors of both AUKF and RAUKF.Figure 7
Estimation errors in force reconstruction by AUKF and RAUKF based on the estimated measurement noise covariances as shown in Figure6: the estimation errors of AUKF are indicated by the black and red lines, respectively.Table 3
Estimation errors of AUKF and RAUKF for the case with constant or small-variation measurement noise covariance.
Error (N)
AUKF
RAUKF
Mean error
0.0551
0.0102
Max error
0.2950
0.1450
RMSE
0.0708
0.0286(2) Measurement Noise Covariance Involving Abrupt Changes. In order to analyse the effect of the weighting factor on measurement noise covariance estimation, trials were conducted under the same conditions by RAUKF via both standard and modified weighting factors for the case of measurement noise covariance involving abrupt changes. The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean. The covariance of this measurement noise was set to 0.036 N within [0, 32 s] and 0.576 N within (32 s, 70s], leading to an abrupt change at 32 s. The initial estimation value was R^0=0.004 N. The window size was m=14. The simulation time was 70 s. Figure 8 shows the system noise covariances estimated by RAUKF via both standard weighting factor and modified weighting factor. As shown in Figure 8, the estimations of measurement noise covariance under both weighting factors follow the reference value of 0.036 N within [0, 32 s]. However, the estimation under the standard weighting factor does not follow the abrupt change from 0.036 N to 0.576 N at 32 s, leading to the large mean error of 0.2662 N after 32 s. In contrast, the estimation under the modified weighting factor follows closely the abrupt change, leading to the mean error of 0.043 N after 32 s.Figure 8
Measurement noise covariance estimations by RAUKF via both standard and modified weighting factors for the case that measurement noise covariance involves abrupt changes: the estimated measurement noise covariances are indicated by the blues lines and the reference values by the red lines.
(a)
RAUKF with the standard weighting factor
(b)
RAUKF with the modified weighting factorFigure9 shows the estimation errors in force reconstruction based on the estimated measurement noise covariances as shown in Figure 8. Before the abrupt change at 32 s, both standard and modified weighting factors result in the similar estimation error. The mean error is 0.0025 N by the standard weighting factor and 0.0023 N by the modified weighting factor. However, after the abrupt change at 32 s, due to the biased estimate of the measurement noise covariance, the use of the standard weighting factor results in a large magnitude of oscillations in the estimation error curve, leading to the RMSE of 0.0907 N. In contrast, with the modified weighting factor, the magnitude of oscillations in the estimation error curve after 32 s is significantly decreased and gradually converged to zero, leading to the RMSE of 0.0444 N. Therefore, it is proved that with the modified weighting factor, RAUKF is able to accommodate abrupt changes in measurement noise covariance. Table 4 lists the estimation errors of RAUKF via both standard and modified weighting factors, where the mean error and maximum error via the standard weighting factor are 0.0603 N and 0.5493 N, while 0.0284 N and 0.2961 N via the modified weighting factor.Figure 9
Estimation errors in force reconstruction by RAUKF corresponding to the estimated measurement noise covariances via both standard and modified weighting factors as shown in Figure8: the estimation error via the standard weighting factor is indicated by the black lines, and the one via the modified weighting factor is indicated by the red lines.Table 4
Estimation errors of RAUKF via both standard and modified weighting factors.
Error (N)
Standard weighting factor
Modified weighting factor
Mean error
0.0603
0.0284
Max error
0.5493
0.2961
RMSE
0.0907
0.0444


### 5.2. Experiments and Analysis

The performance evaluation is also conducted on two experimental cases in the presence of abrupt changes. One is robotic-assisted needle insertion in the presence of rapture events, where the measurement data was acquired from the literature. The other is the mechanical indentation in the presence of geometrical discontinuities to acquire the data of mechanical load and displacement on a phantom tissue.


#### 5.2.1. Robotic-Assisted Needle Insertion

Trials were conducted to evaluate the performance of the proposed methodology for the case of robotic-assisted needle insertion. The experimental data of force and displacement on biological soft tissues during the process of needle insertion were obtained from the literature [28]. The input force profile involves four peaks at the displacements of around 20 mm, 29 mm, 32 mm, and 36 mm (see Figure 10). These peaks represent the rupture events encountered by the robotic needle during the insertion process. For comparison analysis, trials were conducted under the same conditions by UKF, AUKF, and RAUKF with the modified weighting factor. The initial state and noise covariance were set as x0=0.1,0.1,0.01,0.0001,0.0001,1,1, Q0=1, and R0=1. The window size for both AUKF and RAUKF was set to m=4.Figure 10
Forces reconstructed by UKF, AUKF, and RAUKF: the input forces are indicated by the black dashed lines and the reconstructed forces by the red solid lines.
(a)
UKF
(b)
AUKF
(c)
RAUKFFigure10 shows the reconstructed forces by UKF, AUKF, and RAUKF, respectively. It can be seen that the UKF estimation results in large errors at the four peaks, where the estimation error at the first peak (at displacement of around 20 mm) is most significant, which is 0.738 N. An obvious deviation also remained in the displacement range from around 40 mm to 50 mm after the four peaks, leading to the mean error of 0.068 N and the RMSE of 0.15 N.AUKF improves the estimation accuracy of UKF at the four peaks, leading to the mean error of 0.044 N and RMSE of 0.069 N. However, there are still pronounced errors at and after the four peaks. This is because the estimated noise covariance by AUKF involves large error as shown in the previous simulations. In contrast, the reconstructed forces by RAUKF are very close to the reference values at and after the four peaks, leading to the mean error of 0.020 N and RMSE of 0.039 N. Table5 lists the mean and maximum errors as well as RMSEs of the three methods. It can be seen that although both AUKF and RAUKF outperform UKF, RAUKF has the smallest estimation errors (the mean error of 0.020 N, maximum error of 0.205 N, and RMSE of 0.039 N). The mean and maximum errors of RAUKF are more than three times smaller than those of UKF. The RMSE of RAUKF is almost four times smaller than that of UKF. The mean and maximum errors as well as RMSE of RAUKF are also much smaller than those of AUKF.Table 5
Estimation errors by UKF, AUKF, and RAUKF for the case of robotic-assisted needle insertion.
Error (N)
UKF
AUKF
RAUKF
Mean error
0.683
0.044
0.020
Maximum error
0.738
0.346
0.205
RMSE
0.150
0.069
0.039


#### 5.2.2. Mechanical Indentation with Geometrical Discontinuities

Mechanical indentation tests were also conducted to evaluate the performance of the proposed methodology for soft tissue characterization in the presence of geometrical discontinuities. The indentation tests were conducted on a phantom tissue sample using the DMA (Dynamic Mechanical Analyser, Seiko Instruments). As shown in Figure11, the phantom tissue is made up of silicone rubbers (Ecoflex 0030), which have the similar characteristics as human soft tissues [29]. The phantom tissue is of cubic shape (1cm×1cm×0.55cm) to fit into the DMA machine. An indenter of 1 cm diameter is used to compress the phantom tissue vertically. In order to present geometrical discontinuities, the indenter is controlled to compress the phantom tissue with a given displacement and then held at the displacement for around 0.9 s. Subsequently, the indenter is controlled with the maximum speed (1,000,000 μm/min) to increase the displacement by 1 mm and then held at the new displacement for around 0.9 s, with the conduction of three times. The data of displacement and load were acquired from the DMA during the compression tests. As shown in Figure 12, the force profile has three abrupt changes (three steps) at around 0.9 s. 1.8 s, and 2.7 s, which represent geometrical discontinuities in robotic-assisted surgery. Trials were conducted based on measured forces and displacements under the same conditions by UKF, AUKF, and RAUKF with the modified weighting factor. The parameter settings for these methods are the same as that in the previous case of robotic-assisted needle insertion.Figure 11
Experimental setup for mechanical indentation on the phantom tissue.
(a)
The indenter and phantom tissue
(b)
The appearance of DMAFigure 12
Reconstructed forces of UKF, AUKF, and RAUKF for the case of mechanical indentation: the input forces are indicated by the black dashed line and the reconstructed forces by the red solid line.
(a)
UKF
(b)
AUKF
(c)
RAUKFFigure12 shows the reconstructed forces by UKF, AUKF, and RAUKF from measured forces and displacements. It can be seen that the UKF estimation has large estimation errors at each abrupt change, leading to the maximum error of 0.6341 N, mean error of 0.0099 N, and RMSE of 0.0157 N, respectively. Comparing to the UKF estimation, the AUKF estimation follows the reference force curve, without involving a large estimation error. The mean error of AUKF is 0.0064 N, which is smaller than that of UKF. However, due to large error involved in noise covariance estimation, the estimation error of AUKF becomes large after around 2 s with a large peak error at 2.5 s where the input force is changed abruptly. However, due to the use of recursion, the estimation of RAUKF closely follows the reference force curve, especially at the three abrupt changes. Its RMSE is 0.0098 N, which is much smaller than those of UKF and AUKF. Table 6 compares the estimation errors of the three methods. It can be seen that in addition to the mean error, the maximum error and RMSE of RAUKF are also much smaller than those of UKF and AUKF.Table 6
Estimation errors of UKF, AUKF, and RAUKF for the case of mechanical indentation with geometrical discontinuities.
Error (N)
UKF
AUKF
RAUKF
Mean error
0.0099
0.0064
0.0062
Maximum error
0.6341
0.0727
0.0638
RMSE
0.0157
0.0109
0.0098


## 5.1. Simulations and Analysis

Simulations were conducted to evaluate the performance of the proposed methodology in estimating system and measurement noise covariances together with associated force reconstructions. The time steps for the simulations were set to 10 ms. The initial state was defined asx0=0.1,0.1,0.01,0.0001,0.0001,1,1.


### 5.1.1. System Noise Covariance Estimation

For the performance evaluation in terms of system noise covariance estimation, it is supposed that the measurement noise statistics are exactly known. Without loss of generality, setRk=1 and Erk=0. Two cases are studied. One is that system noise covariance is constant or has a small variation, and the other is that system noise covariance involves abrupt changes.(1) Constant or Small-Variation System Noise Covariance. Simulation trials were conducted under the same conditions by the adaptive UKF (AUKF) and recursive adaptive UKF (RAUKF) for the case that system noise covariance is constant or has a small variation. The measurement data were obtained from the H-C model by adding a white Gaussian noise with zero mean and the covariance of 0.036 N. The initial estimation value was Q^0=0.004 N, the window size was m=4, the real covariance of system noise was Qk=0.036 N, and Eqk=0. The simulation time was 90 s. As shown in Figure 2, the system noise covariance estimated by AUKF via windowing approximation involves large oscillations during the entire simulation time. This is because the limited data available within the small window restricts the estimation accuracy. In contrast, the system noise covariance estimated by RAUKF via recursion is close to the real value after 20 s. This is because the estimation via recursion is obtained based on the entire available data. The RMSE (Root Mean Square Error) of AUKF is 0.0125 N, while after 20 s, the RMSE of RAUKF is 0.0032 N. Thus, it is evident that the proposed RAUKF can significantly improve the AUKF accuracy for estimation of system noise covariance.Figure 2
Estimations of system noise covariance by both AUKF and RAUKF for the case that system noise covariance is constant or in a small variation: the estimated system noise covariances are indicated by the blue lines and the reference values by the red lines.
(a)
AUKF via windowing approximation
(b)
RAUKF via recursionFigure3 illustrates the estimation errors in terms of force reconstruction by both AUKF and RAUKF, based on the estimated system noise covariances as shown in Figure 2. It is clear that the estimation error of RAUKF is much smaller than that of AUKF. During the entire simulation time, since the system noise covariance estimated by AUKF involves large oscillations (see Figure 2(a)), the resultant estimation error also involves large oscillations. However, for RAUKF, the resultant estimation error becomes small and is further converged after the initial time period of 20 s, following the trend of the estimated system noise covariance. The maximum estimation error is 0.0394 N for RAUKF, but 0.1853 N for AUKF. The mean estimation error of AUKF is 0.0374 N, while that of RAUKF is 0.0053 N. The RMSE of RAUKF is 0.0070 N, while that of AUKF is almost nine times larger, which is 0.0475 N. Table 1 summarizes the estimation errors of both AUKF and RAUKF.Figure 3
Estimation errors in force reconstruction by AUKF and RAUKF corresponding to the estimated system noise covariances as shown in Figure2: the estimation error of AUKF is indicated by the black lines and that of RAUKF by the red lines.Table 1
Estimation errors of both AUKF and RAUKF for the case with constant or small-variation system noise covariance.
Error (N)
AUKF
RAUKF
Mean error
0.0374
0.0053
Max error
0.1853
0.0394
RMSE
0.0475
0.0070(2) System Noise Covariance Involving Abrupt Changes. In order to analyse the effect of the weighting factor on system noise covariance estimation, simulation trials were conducted under the same conditions by RAUKF via both standard and modified weighting factors for the case of system noise covariance involving abrupt changes. The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean. The noise covariance was set to 0.036 N within [0, 36 s] and 0.576 N within (36 s, 70 s], leading to an abrupt change from 0.036 N to 0.576 N at 36 s. The initial estimation value was Q^0=0.004 N, the window size was set to m=14, and the simulation time was 70 s.Figure4 shows the system noise covariances estimated by RAUKF via both standard and modified weighting factors. Although the estimation results with both weighting factors follow the reference value of 0.036 N within [0, 36 s], the estimation under the standard weighting factor does not follow the abrupt change from 0.036 N to 0.576 N at 36 s, leading to the mean error of 0.372 N after 36 s. In contrast, the estimation of system noise covariance under the modified weighting factor rapidly follows the abrupt change from 0.036 N to 0.576 N at 36 s, leading to the mean error of 0.0088 N after 36 s. This demonstrates that with the modified weighting factor, RAUKF can accommodate abrupt changes in system noise covariance.Figure 4
Estimations of system noise covariance by RAUKF with both standard and modified weighting factors for the case that system noise covariance involves abrupt changes: the estimated system noise covariances are indicated by the red lines and the reference values by the blue lines.
(a)
RAUKF with the standard weighting factor
(b)
RAUKF with the modified weighting factorFigure5 shows the estimation errors in force reconstruction based on the estimated system noise covariances as shown in Figure 4. It can be seen that before 36 s, although the RAUKF estimation errors via both standard and modified weighting factors are small and converged to zero, the one via the modified weighting factor is even smaller and has a faster convergence speed. However, after 36 s, due to the biased estimate of the system noise covariance, the RAUKF estimation error via the standard weighting factor involves large oscillations, leading to the RMSE of 0.0749 N. In contrast, the RAUKF estimation error via the modified weighting factor is much smaller and converged to zero, leading to the RMSE of 0.0046 N. This demonstrates that the modified weighting factor enables RAUKF to account for abrupt changes in system noise covariance, leading to improved estimation accuracy. Table 2 summarizes the RAUKF estimation errors via both standard and modified weighting factors.Figure 5
Estimation errors in terms of force reconstruction by RAUKF based on the estimated system noise covariances via both standard and modified weighting factors as shown in Figure4: the estimation error via the standard weighting factor is indicated by the black lines and the one via the modified weighting factor by the red lines.Table 2
RAUKF estimation errors via both standard and modified weighting factors for the case that system noise covariance involves abrupt changes.
Error (N)
Standard weighting factor
Modified weighting factor
Mean error
0.0030
0.0009
Max error
0.4517
0.0389
RMSE
0.0749
0.0046


### 5.1.2. Measurement Noise Covariance Estimation

In order to evaluate the estimation of measurement noise covariance, it is assumed that the system noise statistics are exactly known. Without loss of generality, chooseQk=1 and Eqk=0. Similar to Section 5.1.1., two cases are studied. One is that the measurement noise covariance is constant or in a small variation, and the other is that the measurement noise covariance involves abrupt changes.(1) Constant or Small-Variation Measurement Noise Covariance. Simulation trials were conducted under the condition that the measurement noise covariance is constant or has a small variation by both AUKF and RAUKF.The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean and the covariance of 0.036 N. The initial estimation value wasR^0=0.004 N. The window size was set to m=4. The simulation time was 70 s. As shown in Figure 6, the estimate of measurement noise covariance by AUKF involves large oscillations during the entire simulation time, while the one by RAUKF is very close to the reference value after 15 s. The RMSE of RAUKF is 0.0051 N and is almost twice smaller than that of AUKF which is 0.0101 N. This demonstrates that the accuracy of RAUKF is significantly higher than that of AUKF for estimation of measurement noise covariance.Figure 6
Estimations of measurement noise covariance by AUKF and RAUKF for the case that measurement noise covariance is constant or small-variation: the estimated measurement noise covariances are indicated by the blue lines and the reference values by the red lines.
(a)
AUKF via windowing approximation
(b)
RAUKF via recursionFigure7 illustrates the estimation errors in force reconstruction by both AUKF and RAUKF based on the estimated measurement noise covariances as shown in Figure 6. It can be seen clearly that the estimation error of AUKF is much larger than that of RAUKF during the entire simulation time. In addition, after 15 s, the estimation error of RAUKF is close to zero, while that of AUKF is slightly smaller than the initial estimation error. This is because the estimated measurement noise covariance by windowing approximation is fluctuated during the entire simulation time (see Figure 6(a)), while that by RAUKF is converged to the reference value (see Figure 6(b)). The RMSE of RAUKF is 0.0286 N and is two-three times smaller than that of AUKF which is 0.0708 N. Further, the maximum and mean estimation errors of AUKF are about two times larger than those of RAUKF. The maximum estimation error is 0.2905 N for AUKF, while 0.145 N for RAUKF. The mean estimation error is 0.0551 N for AUKF, while 0.0102 N for RAUKF. Therefore, the estimation accuracy of RAUKF is significantly higher than that of AUKF. Table 3 compares the estimation errors of both AUKF and RAUKF.Figure 7
Estimation errors in force reconstruction by AUKF and RAUKF based on the estimated measurement noise covariances as shown in Figure6: the estimation errors of AUKF are indicated by the black and red lines, respectively.Table 3
Estimation errors of AUKF and RAUKF for the case with constant or small-variation measurement noise covariance.
Error (N)
AUKF
RAUKF
Mean error
0.0551
0.0102
Max error
0.2950
0.1450
RMSE
0.0708
0.0286(2) Measurement Noise Covariance Involving Abrupt Changes. In order to analyse the effect of the weighting factor on measurement noise covariance estimation, trials were conducted under the same conditions by RAUKF via both standard and modified weighting factors for the case of measurement noise covariance involving abrupt changes. The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean. The covariance of this measurement noise was set to 0.036 N within [0, 32 s] and 0.576 N within (32 s, 70s], leading to an abrupt change at 32 s. The initial estimation value was R^0=0.004 N. The window size was m=14. The simulation time was 70 s. Figure 8 shows the system noise covariances estimated by RAUKF via both standard weighting factor and modified weighting factor. As shown in Figure 8, the estimations of measurement noise covariance under both weighting factors follow the reference value of 0.036 N within [0, 32 s]. However, the estimation under the standard weighting factor does not follow the abrupt change from 0.036 N to 0.576 N at 32 s, leading to the large mean error of 0.2662 N after 32 s. In contrast, the estimation under the modified weighting factor follows closely the abrupt change, leading to the mean error of 0.043 N after 32 s.Figure 8
Measurement noise covariance estimations by RAUKF via both standard and modified weighting factors for the case that measurement noise covariance involves abrupt changes: the estimated measurement noise covariances are indicated by the blues lines and the reference values by the red lines.
(a)
RAUKF with the standard weighting factor
(b)
RAUKF with the modified weighting factorFigure9 shows the estimation errors in force reconstruction based on the estimated measurement noise covariances as shown in Figure 8. Before the abrupt change at 32 s, both standard and modified weighting factors result in the similar estimation error. The mean error is 0.0025 N by the standard weighting factor and 0.0023 N by the modified weighting factor. However, after the abrupt change at 32 s, due to the biased estimate of the measurement noise covariance, the use of the standard weighting factor results in a large magnitude of oscillations in the estimation error curve, leading to the RMSE of 0.0907 N. In contrast, with the modified weighting factor, the magnitude of oscillations in the estimation error curve after 32 s is significantly decreased and gradually converged to zero, leading to the RMSE of 0.0444 N. Therefore, it is proved that with the modified weighting factor, RAUKF is able to accommodate abrupt changes in measurement noise covariance. Table 4 lists the estimation errors of RAUKF via both standard and modified weighting factors, where the mean error and maximum error via the standard weighting factor are 0.0603 N and 0.5493 N, while 0.0284 N and 0.2961 N via the modified weighting factor.Figure 9
Estimation errors in force reconstruction by RAUKF corresponding to the estimated measurement noise covariances via both standard and modified weighting factors as shown in Figure8: the estimation error via the standard weighting factor is indicated by the black lines, and the one via the modified weighting factor is indicated by the red lines.Table 4
Estimation errors of RAUKF via both standard and modified weighting factors.
Error (N)
Standard weighting factor
Modified weighting factor
Mean error
0.0603
0.0284
Max error
0.5493
0.2961
RMSE
0.0907
0.0444


## 5.1.1. System Noise Covariance Estimation

For the performance evaluation in terms of system noise covariance estimation, it is supposed that the measurement noise statistics are exactly known. Without loss of generality, setRk=1 and Erk=0. Two cases are studied. One is that system noise covariance is constant or has a small variation, and the other is that system noise covariance involves abrupt changes.(1) Constant or Small-Variation System Noise Covariance. Simulation trials were conducted under the same conditions by the adaptive UKF (AUKF) and recursive adaptive UKF (RAUKF) for the case that system noise covariance is constant or has a small variation. The measurement data were obtained from the H-C model by adding a white Gaussian noise with zero mean and the covariance of 0.036 N. The initial estimation value was Q^0=0.004 N, the window size was m=4, the real covariance of system noise was Qk=0.036 N, and Eqk=0. The simulation time was 90 s. As shown in Figure 2, the system noise covariance estimated by AUKF via windowing approximation involves large oscillations during the entire simulation time. This is because the limited data available within the small window restricts the estimation accuracy. In contrast, the system noise covariance estimated by RAUKF via recursion is close to the real value after 20 s. This is because the estimation via recursion is obtained based on the entire available data. The RMSE (Root Mean Square Error) of AUKF is 0.0125 N, while after 20 s, the RMSE of RAUKF is 0.0032 N. Thus, it is evident that the proposed RAUKF can significantly improve the AUKF accuracy for estimation of system noise covariance.Figure 2
Estimations of system noise covariance by both AUKF and RAUKF for the case that system noise covariance is constant or in a small variation: the estimated system noise covariances are indicated by the blue lines and the reference values by the red lines.
(a)
AUKF via windowing approximation
(b)
RAUKF via recursionFigure3 illustrates the estimation errors in terms of force reconstruction by both AUKF and RAUKF, based on the estimated system noise covariances as shown in Figure 2. It is clear that the estimation error of RAUKF is much smaller than that of AUKF. During the entire simulation time, since the system noise covariance estimated by AUKF involves large oscillations (see Figure 2(a)), the resultant estimation error also involves large oscillations. However, for RAUKF, the resultant estimation error becomes small and is further converged after the initial time period of 20 s, following the trend of the estimated system noise covariance. The maximum estimation error is 0.0394 N for RAUKF, but 0.1853 N for AUKF. The mean estimation error of AUKF is 0.0374 N, while that of RAUKF is 0.0053 N. The RMSE of RAUKF is 0.0070 N, while that of AUKF is almost nine times larger, which is 0.0475 N. Table 1 summarizes the estimation errors of both AUKF and RAUKF.Figure 3
Estimation errors in force reconstruction by AUKF and RAUKF corresponding to the estimated system noise covariances as shown in Figure2: the estimation error of AUKF is indicated by the black lines and that of RAUKF by the red lines.Table 1
Estimation errors of both AUKF and RAUKF for the case with constant or small-variation system noise covariance.
Error (N)
AUKF
RAUKF
Mean error
0.0374
0.0053
Max error
0.1853
0.0394
RMSE
0.0475
0.0070(2) System Noise Covariance Involving Abrupt Changes. In order to analyse the effect of the weighting factor on system noise covariance estimation, simulation trials were conducted under the same conditions by RAUKF via both standard and modified weighting factors for the case of system noise covariance involving abrupt changes. The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean. The noise covariance was set to 0.036 N within [0, 36 s] and 0.576 N within (36 s, 70 s], leading to an abrupt change from 0.036 N to 0.576 N at 36 s. The initial estimation value was Q^0=0.004 N, the window size was set to m=14, and the simulation time was 70 s.Figure4 shows the system noise covariances estimated by RAUKF via both standard and modified weighting factors. Although the estimation results with both weighting factors follow the reference value of 0.036 N within [0, 36 s], the estimation under the standard weighting factor does not follow the abrupt change from 0.036 N to 0.576 N at 36 s, leading to the mean error of 0.372 N after 36 s. In contrast, the estimation of system noise covariance under the modified weighting factor rapidly follows the abrupt change from 0.036 N to 0.576 N at 36 s, leading to the mean error of 0.0088 N after 36 s. This demonstrates that with the modified weighting factor, RAUKF can accommodate abrupt changes in system noise covariance.Figure 4
Estimations of system noise covariance by RAUKF with both standard and modified weighting factors for the case that system noise covariance involves abrupt changes: the estimated system noise covariances are indicated by the red lines and the reference values by the blue lines.
(a)
RAUKF with the standard weighting factor
(b)
RAUKF with the modified weighting factorFigure5 shows the estimation errors in force reconstruction based on the estimated system noise covariances as shown in Figure 4. It can be seen that before 36 s, although the RAUKF estimation errors via both standard and modified weighting factors are small and converged to zero, the one via the modified weighting factor is even smaller and has a faster convergence speed. However, after 36 s, due to the biased estimate of the system noise covariance, the RAUKF estimation error via the standard weighting factor involves large oscillations, leading to the RMSE of 0.0749 N. In contrast, the RAUKF estimation error via the modified weighting factor is much smaller and converged to zero, leading to the RMSE of 0.0046 N. This demonstrates that the modified weighting factor enables RAUKF to account for abrupt changes in system noise covariance, leading to improved estimation accuracy. Table 2 summarizes the RAUKF estimation errors via both standard and modified weighting factors.Figure 5
Estimation errors in terms of force reconstruction by RAUKF based on the estimated system noise covariances via both standard and modified weighting factors as shown in Figure4: the estimation error via the standard weighting factor is indicated by the black lines and the one via the modified weighting factor by the red lines.Table 2
RAUKF estimation errors via both standard and modified weighting factors for the case that system noise covariance involves abrupt changes.
Error (N)
Standard weighting factor
Modified weighting factor
Mean error
0.0030
0.0009
Max error
0.4517
0.0389
RMSE
0.0749
0.0046


## 5.1.2. Measurement Noise Covariance Estimation

In order to evaluate the estimation of measurement noise covariance, it is assumed that the system noise statistics are exactly known. Without loss of generality, chooseQk=1 and Eqk=0. Similar to Section 5.1.1., two cases are studied. One is that the measurement noise covariance is constant or in a small variation, and the other is that the measurement noise covariance involves abrupt changes.(1) Constant or Small-Variation Measurement Noise Covariance. Simulation trials were conducted under the condition that the measurement noise covariance is constant or has a small variation by both AUKF and RAUKF.The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean and the covariance of 0.036 N. The initial estimation value wasR^0=0.004 N. The window size was set to m=4. The simulation time was 70 s. As shown in Figure 6, the estimate of measurement noise covariance by AUKF involves large oscillations during the entire simulation time, while the one by RAUKF is very close to the reference value after 15 s. The RMSE of RAUKF is 0.0051 N and is almost twice smaller than that of AUKF which is 0.0101 N. This demonstrates that the accuracy of RAUKF is significantly higher than that of AUKF for estimation of measurement noise covariance.Figure 6
Estimations of measurement noise covariance by AUKF and RAUKF for the case that measurement noise covariance is constant or small-variation: the estimated measurement noise covariances are indicated by the blue lines and the reference values by the red lines.
(a)
AUKF via windowing approximation
(b)
RAUKF via recursionFigure7 illustrates the estimation errors in force reconstruction by both AUKF and RAUKF based on the estimated measurement noise covariances as shown in Figure 6. It can be seen clearly that the estimation error of AUKF is much larger than that of RAUKF during the entire simulation time. In addition, after 15 s, the estimation error of RAUKF is close to zero, while that of AUKF is slightly smaller than the initial estimation error. This is because the estimated measurement noise covariance by windowing approximation is fluctuated during the entire simulation time (see Figure 6(a)), while that by RAUKF is converged to the reference value (see Figure 6(b)). The RMSE of RAUKF is 0.0286 N and is two-three times smaller than that of AUKF which is 0.0708 N. Further, the maximum and mean estimation errors of AUKF are about two times larger than those of RAUKF. The maximum estimation error is 0.2905 N for AUKF, while 0.145 N for RAUKF. The mean estimation error is 0.0551 N for AUKF, while 0.0102 N for RAUKF. Therefore, the estimation accuracy of RAUKF is significantly higher than that of AUKF. Table 3 compares the estimation errors of both AUKF and RAUKF.Figure 7
Estimation errors in force reconstruction by AUKF and RAUKF based on the estimated measurement noise covariances as shown in Figure6: the estimation errors of AUKF are indicated by the black and red lines, respectively.Table 3
Estimation errors of AUKF and RAUKF for the case with constant or small-variation measurement noise covariance.
Error (N)
AUKF
RAUKF
Mean error
0.0551
0.0102
Max error
0.2950
0.1450
RMSE
0.0708
0.0286(2) Measurement Noise Covariance Involving Abrupt Changes. In order to analyse the effect of the weighting factor on measurement noise covariance estimation, trials were conducted under the same conditions by RAUKF via both standard and modified weighting factors for the case of measurement noise covariance involving abrupt changes. The input signals were obtained from the H-C model by adding a white Gaussian noise with zero mean. The covariance of this measurement noise was set to 0.036 N within [0, 32 s] and 0.576 N within (32 s, 70s], leading to an abrupt change at 32 s. The initial estimation value was R^0=0.004 N. The window size was m=14. The simulation time was 70 s. Figure 8 shows the system noise covariances estimated by RAUKF via both standard weighting factor and modified weighting factor. As shown in Figure 8, the estimations of measurement noise covariance under both weighting factors follow the reference value of 0.036 N within [0, 32 s]. However, the estimation under the standard weighting factor does not follow the abrupt change from 0.036 N to 0.576 N at 32 s, leading to the large mean error of 0.2662 N after 32 s. In contrast, the estimation under the modified weighting factor follows closely the abrupt change, leading to the mean error of 0.043 N after 32 s.Figure 8
Measurement noise covariance estimations by RAUKF via both standard and modified weighting factors for the case that measurement noise covariance involves abrupt changes: the estimated measurement noise covariances are indicated by the blues lines and the reference values by the red lines.
(a)
RAUKF with the standard weighting factor
(b)
RAUKF with the modified weighting factorFigure9 shows the estimation errors in force reconstruction based on the estimated measurement noise covariances as shown in Figure 8. Before the abrupt change at 32 s, both standard and modified weighting factors result in the similar estimation error. The mean error is 0.0025 N by the standard weighting factor and 0.0023 N by the modified weighting factor. However, after the abrupt change at 32 s, due to the biased estimate of the measurement noise covariance, the use of the standard weighting factor results in a large magnitude of oscillations in the estimation error curve, leading to the RMSE of 0.0907 N. In contrast, with the modified weighting factor, the magnitude of oscillations in the estimation error curve after 32 s is significantly decreased and gradually converged to zero, leading to the RMSE of 0.0444 N. Therefore, it is proved that with the modified weighting factor, RAUKF is able to accommodate abrupt changes in measurement noise covariance. Table 4 lists the estimation errors of RAUKF via both standard and modified weighting factors, where the mean error and maximum error via the standard weighting factor are 0.0603 N and 0.5493 N, while 0.0284 N and 0.2961 N via the modified weighting factor.Figure 9
Estimation errors in force reconstruction by RAUKF corresponding to the estimated measurement noise covariances via both standard and modified weighting factors as shown in Figure8: the estimation error via the standard weighting factor is indicated by the black lines, and the one via the modified weighting factor is indicated by the red lines.Table 4
Estimation errors of RAUKF via both standard and modified weighting factors.
Error (N)
Standard weighting factor
Modified weighting factor
Mean error
0.0603
0.0284
Max error
0.5493
0.2961
RMSE
0.0907
0.0444


## 5.2. Experiments and Analysis

The performance evaluation is also conducted on two experimental cases in the presence of abrupt changes. One is robotic-assisted needle insertion in the presence of rapture events, where the measurement data was acquired from the literature. The other is the mechanical indentation in the presence of geometrical discontinuities to acquire the data of mechanical load and displacement on a phantom tissue.


### 5.2.1. Robotic-Assisted Needle Insertion

Trials were conducted to evaluate the performance of the proposed methodology for the case of robotic-assisted needle insertion. The experimental data of force and displacement on biological soft tissues during the process of needle insertion were obtained from the literature [28]. The input force profile involves four peaks at the displacements of around 20 mm, 29 mm, 32 mm, and 36 mm (see Figure 10). These peaks represent the rupture events encountered by the robotic needle during the insertion process. For comparison analysis, trials were conducted under the same conditions by UKF, AUKF, and RAUKF with the modified weighting factor. The initial state and noise covariance were set as x0=0.1,0.1,0.01,0.0001,0.0001,1,1, Q0=1, and R0=1. The window size for both AUKF and RAUKF was set to m=4.Figure 10
Forces reconstructed by UKF, AUKF, and RAUKF: the input forces are indicated by the black dashed lines and the reconstructed forces by the red solid lines.
(a)
UKF
(b)
AUKF
(c)
RAUKFFigure10 shows the reconstructed forces by UKF, AUKF, and RAUKF, respectively. It can be seen that the UKF estimation results in large errors at the four peaks, where the estimation error at the first peak (at displacement of around 20 mm) is most significant, which is 0.738 N. An obvious deviation also remained in the displacement range from around 40 mm to 50 mm after the four peaks, leading to the mean error of 0.068 N and the RMSE of 0.15 N.AUKF improves the estimation accuracy of UKF at the four peaks, leading to the mean error of 0.044 N and RMSE of 0.069 N. However, there are still pronounced errors at and after the four peaks. This is because the estimated noise covariance by AUKF involves large error as shown in the previous simulations. In contrast, the reconstructed forces by RAUKF are very close to the reference values at and after the four peaks, leading to the mean error of 0.020 N and RMSE of 0.039 N. Table5 lists the mean and maximum errors as well as RMSEs of the three methods. It can be seen that although both AUKF and RAUKF outperform UKF, RAUKF has the smallest estimation errors (the mean error of 0.020 N, maximum error of 0.205 N, and RMSE of 0.039 N). The mean and maximum errors of RAUKF are more than three times smaller than those of UKF. The RMSE of RAUKF is almost four times smaller than that of UKF. The mean and maximum errors as well as RMSE of RAUKF are also much smaller than those of AUKF.Table 5
Estimation errors by UKF, AUKF, and RAUKF for the case of robotic-assisted needle insertion.
Error (N)
UKF
AUKF
RAUKF
Mean error
0.683
0.044
0.020
Maximum error
0.738
0.346
0.205
RMSE
0.150
0.069
0.039


### 5.2.2. Mechanical Indentation with Geometrical Discontinuities

Mechanical indentation tests were also conducted to evaluate the performance of the proposed methodology for soft tissue characterization in the presence of geometrical discontinuities. The indentation tests were conducted on a phantom tissue sample using the DMA (Dynamic Mechanical Analyser, Seiko Instruments). As shown in Figure11, the phantom tissue is made up of silicone rubbers (Ecoflex 0030), which have the similar characteristics as human soft tissues [29]. The phantom tissue is of cubic shape (1cm×1cm×0.55cm) to fit into the DMA machine. An indenter of 1 cm diameter is used to compress the phantom tissue vertically. In order to present geometrical discontinuities, the indenter is controlled to compress the phantom tissue with a given displacement and then held at the displacement for around 0.9 s. Subsequently, the indenter is controlled with the maximum speed (1,000,000 μm/min) to increase the displacement by 1 mm and then held at the new displacement for around 0.9 s, with the conduction of three times. The data of displacement and load were acquired from the DMA during the compression tests. As shown in Figure 12, the force profile has three abrupt changes (three steps) at around 0.9 s. 1.8 s, and 2.7 s, which represent geometrical discontinuities in robotic-assisted surgery. Trials were conducted based on measured forces and displacements under the same conditions by UKF, AUKF, and RAUKF with the modified weighting factor. The parameter settings for these methods are the same as that in the previous case of robotic-assisted needle insertion.Figure 11
Experimental setup for mechanical indentation on the phantom tissue.
(a)
The indenter and phantom tissue
(b)
The appearance of DMAFigure 12
Reconstructed forces of UKF, AUKF, and RAUKF for the case of mechanical indentation: the input forces are indicated by the black dashed line and the reconstructed forces by the red solid line.
(a)
UKF
(b)
AUKF
(c)
RAUKFFigure12 shows the reconstructed forces by UKF, AUKF, and RAUKF from measured forces and displacements. It can be seen that the UKF estimation has large estimation errors at each abrupt change, leading to the maximum error of 0.6341 N, mean error of 0.0099 N, and RMSE of 0.0157 N, respectively. Comparing to the UKF estimation, the AUKF estimation follows the reference force curve, without involving a large estimation error. The mean error of AUKF is 0.0064 N, which is smaller than that of UKF. However, due to large error involved in noise covariance estimation, the estimation error of AUKF becomes large after around 2 s with a large peak error at 2.5 s where the input force is changed abruptly. However, due to the use of recursion, the estimation of RAUKF closely follows the reference force curve, especially at the three abrupt changes. Its RMSE is 0.0098 N, which is much smaller than those of UKF and AUKF. Table 6 compares the estimation errors of the three methods. It can be seen that in addition to the mean error, the maximum error and RMSE of RAUKF are also much smaller than those of UKF and AUKF.Table 6
Estimation errors of UKF, AUKF, and RAUKF for the case of mechanical indentation with geometrical discontinuities.
Error (N)
UKF
AUKF
RAUKF
Mean error
0.0099
0.0064
0.0062
Maximum error
0.6341
0.0727
0.0638
RMSE
0.0157
0.0109
0.0098


## 5.2.1. Robotic-Assisted Needle Insertion

Trials were conducted to evaluate the performance of the proposed methodology for the case of robotic-assisted needle insertion. The experimental data of force and displacement on biological soft tissues during the process of needle insertion were obtained from the literature [28]. The input force profile involves four peaks at the displacements of around 20 mm, 29 mm, 32 mm, and 36 mm (see Figure 10). These peaks represent the rupture events encountered by the robotic needle during the insertion process. For comparison analysis, trials were conducted under the same conditions by UKF, AUKF, and RAUKF with the modified weighting factor. The initial state and noise covariance were set as x0=0.1,0.1,0.01,0.0001,0.0001,1,1, Q0=1, and R0=1. The window size for both AUKF and RAUKF was set to m=4.Figure 10
Forces reconstructed by UKF, AUKF, and RAUKF: the input forces are indicated by the black dashed lines and the reconstructed forces by the red solid lines.
(a)
UKF
(b)
AUKF
(c)
RAUKFFigure10 shows the reconstructed forces by UKF, AUKF, and RAUKF, respectively. It can be seen that the UKF estimation results in large errors at the four peaks, where the estimation error at the first peak (at displacement of around 20 mm) is most significant, which is 0.738 N. An obvious deviation also remained in the displacement range from around 40 mm to 50 mm after the four peaks, leading to the mean error of 0.068 N and the RMSE of 0.15 N.AUKF improves the estimation accuracy of UKF at the four peaks, leading to the mean error of 0.044 N and RMSE of 0.069 N. However, there are still pronounced errors at and after the four peaks. This is because the estimated noise covariance by AUKF involves large error as shown in the previous simulations. In contrast, the reconstructed forces by RAUKF are very close to the reference values at and after the four peaks, leading to the mean error of 0.020 N and RMSE of 0.039 N. Table5 lists the mean and maximum errors as well as RMSEs of the three methods. It can be seen that although both AUKF and RAUKF outperform UKF, RAUKF has the smallest estimation errors (the mean error of 0.020 N, maximum error of 0.205 N, and RMSE of 0.039 N). The mean and maximum errors of RAUKF are more than three times smaller than those of UKF. The RMSE of RAUKF is almost four times smaller than that of UKF. The mean and maximum errors as well as RMSE of RAUKF are also much smaller than those of AUKF.Table 5
Estimation errors by UKF, AUKF, and RAUKF for the case of robotic-assisted needle insertion.
Error (N)
UKF
AUKF
RAUKF
Mean error
0.683
0.044
0.020
Maximum error
0.738
0.346
0.205
RMSE
0.150
0.069
0.039


## 5.2.2. Mechanical Indentation with Geometrical Discontinuities

Mechanical indentation tests were also conducted to evaluate the performance of the proposed methodology for soft tissue characterization in the presence of geometrical discontinuities. The indentation tests were conducted on a phantom tissue sample using the DMA (Dynamic Mechanical Analyser, Seiko Instruments). As shown in Figure11, the phantom tissue is made up of silicone rubbers (Ecoflex 0030), which have the similar characteristics as human soft tissues [29]. The phantom tissue is of cubic shape (1cm×1cm×0.55cm) to fit into the DMA machine. An indenter of 1 cm diameter is used to compress the phantom tissue vertically. In order to present geometrical discontinuities, the indenter is controlled to compress the phantom tissue with a given displacement and then held at the displacement for around 0.9 s. Subsequently, the indenter is controlled with the maximum speed (1,000,000 μm/min) to increase the displacement by 1 mm and then held at the new displacement for around 0.9 s, with the conduction of three times. The data of displacement and load were acquired from the DMA during the compression tests. As shown in Figure 12, the force profile has three abrupt changes (three steps) at around 0.9 s. 1.8 s, and 2.7 s, which represent geometrical discontinuities in robotic-assisted surgery. Trials were conducted based on measured forces and displacements under the same conditions by UKF, AUKF, and RAUKF with the modified weighting factor. The parameter settings for these methods are the same as that in the previous case of robotic-assisted needle insertion.Figure 11
Experimental setup for mechanical indentation on the phantom tissue.
(a)
The indenter and phantom tissue
(b)
The appearance of DMAFigure 12
Reconstructed forces of UKF, AUKF, and RAUKF for the case of mechanical indentation: the input forces are indicated by the black dashed line and the reconstructed forces by the red solid line.
(a)
UKF
(b)
AUKF
(c)
RAUKFFigure12 shows the reconstructed forces by UKF, AUKF, and RAUKF from measured forces and displacements. It can be seen that the UKF estimation has large estimation errors at each abrupt change, leading to the maximum error of 0.6341 N, mean error of 0.0099 N, and RMSE of 0.0157 N, respectively. Comparing to the UKF estimation, the AUKF estimation follows the reference force curve, without involving a large estimation error. The mean error of AUKF is 0.0064 N, which is smaller than that of UKF. However, due to large error involved in noise covariance estimation, the estimation error of AUKF becomes large after around 2 s with a large peak error at 2.5 s where the input force is changed abruptly. However, due to the use of recursion, the estimation of RAUKF closely follows the reference force curve, especially at the three abrupt changes. Its RMSE is 0.0098 N, which is much smaller than those of UKF and AUKF. Table 6 compares the estimation errors of the three methods. It can be seen that in addition to the mean error, the maximum error and RMSE of RAUKF are also much smaller than those of UKF and AUKF.Table 6
Estimation errors of UKF, AUKF, and RAUKF for the case of mechanical indentation with geometrical discontinuities.
Error (N)
UKF
AUKF
RAUKF
Mean error
0.0099
0.0064
0.0062
Maximum error
0.6341
0.0727
0.0638
RMSE
0.0157
0.0109
0.0098


## 6. Conclusions

This paper presents a new methodology for online nonlinear soft tissue characterization in robotic-assisted surgery. An adaptive UKF is established based on the nonlinear H-C model for online estimate soft tissue parameters without prior knowledge on noise statistics. Based on this, a recursive adaptive UKF with the modified weighting factor is further developed to improve the estimation accuracy and accommodate large variations in noise statistics. Simulation and experimental results as well as comparison analysis demonstrate that the proposed methodology can effectively estimate soft tissue parameters under system and measurement noises in both small and large variations, leading to improved accuracy and robustness in comparison with the conventional UKF.Future work will focus on two aspects. One is to investigate the stability of the proposed methodology, including the effects of window size and initial state error on the filtering performance. The other is to incorporate the proposed methodology in robotic control and force feedback for robotic-assisted minimally invasive surgery. It is expected that new control strategies will be developed by integrating the proposed methodology into the control loop for precise and stable robotic control with force feedback in the presence of unknown system noises.



---

*Source: 9873410-2020-03-04.xml*