This article has been accepted for inclusion in a future issue of this journal. Content is final as presented, with the exception of pagination.
ADHIKARI et al.: PHYSIOLOGICAL TREMOR FILTERING WITHOUT PHASE DISTORTION FOR ROBOTIC MICROSURGERY 3
[33]. These algorithms, however, still relied on IIR prefilters
to isolate the tremor signal from the whole motion. As a
consequence, the isolated tremor was completely distorted.
The distorted tremor signal was then rectified by mapping it
with the nondistorted tremor signal obtained from zero-phase
filtering. Since the tremor signal is nonstationary, for an
accurate rectification, these methods necessitated the tremor
signal to be zero-phase filtered continuously to adapt to its
changing characteristics. However, the zero-phase filtering is
a noncausal system and cannot be performed in real time. This
limitation renders the aforementioned methods that are inade-
quate for real-time implementation. Although frequent off-line
mapping can be performed, the filter’s transient characteristic
restricts the prompt update of the mapping, which degrades
the filtering performance significantly.
In this article, we offer a new paradigm to address this
problem. Instead of prefiltering the tremor signal, we directly
extract the voluntary and the tremor signal from the whole
motion by decomposing it. We first extract the voluntary signal
and then deduct it from the whole motion to approximate
the tremor signal. Our approach comprises two units: 1) an
estimator and 2) a predictor. We achieve zero-phase type
estimation of voluntary motion by adopting the recursive
singular spectrum analysis (RSSA) [34] algorithm. However,
this zero-phase estimation comes at the cost of a fixed time
delay. To shorten this delay, we append a predictor unit, that is,
an ELM [35] to the RSSA estimator. This combined structure
is tested with the real physiological tremor data recorded from
five novices and four microsurgeons, and the results show that
the proposed paradigm: 1) yields zero-phase type filtering of
the voluntary and tremor signals without causing any distortion
and 2) attenuates the unwanted tremor below 10-μm limit that
is desired in microsurgeries. The main contribution of this
article is the novel paradigm that can filter the tremor signal
as accurately as a high-order linear phase filter, yet within
a small processing delay, such as a nonlinear phase filter.
Hence, it offers the best of both conventional filters, without
the need for any prefiltering or zero-phase filtering. This is the
major advantage of the proposed framework compared with
the existing state-of-the-art techniques.
This article is organized as follows. Section II provides
the theoretical overviews of the SSA and ELM algorithms.
In Section III, we present the proposed RSSA-ELM and mov-
ing window RSSA-ELM methods. We evaluate the off-line and
real-time performance of these algorithms using real tremor
data in Section IV. Section V contains a discussion of results.
Finally, we conclude this article in Section VI. Preliminary
results were reported in [36].
II. R
EVIEW OF EXISTING METHODS
A. Singular Spectrum Analysis
The SSA algorithm decomposes a time series into a number
of interpretable components by using singular value decompo-
sition (SVD) of a lag-covariance matrix formed from a time
series. The SSA algorithm is computed in two complementary
stages: decomposition and reconstruction [34].
1) Decomposition: The decomposition stage comprises two
substages, embedding followed by SVD. In the embedding
substage, the original 1-D time series x
N
=[x
1
, x
2
,...,x
N
]
is reconstructed into a multidimensional series. The multidi-
mensional series is a sequence of L-lagged vectors x
i
,which
forms the trajectory matrix, X
L×K
=[x
1
x
2
x
3
··· x
K
],
where x
1
=[x
1
, x
2
,...,x
L
]
T
, x
2
=[x
2
, x
3
,...,x
L+1
]
T
, x
K
=
[x
K
, x
K +1
,...,x
N
]
T
,and(.)
T
denotes the transpose operator.
The length of each vector is called the window length, L,
which can range between [2, N −1]. For a time series of length
N, the total number of lagged vectors is K = N − L + 1.
Following embedding, SVD is performed on the covariance
matrix S = XX
T
. The SVD operation enables the trajectory
matrix X to be represented in terms of its L elementary
matrices: X = X
1
+ X
2
+···+X
L
,whereX
i
= (λ
i
)
1/2
u
i
v
T
i
,
λ
i
is the i th eigenvalue, and u
i
and v
i
are the i th left and right
eigenvectors, respectively.
2) Reconstruction: This step comprises eigentriple grouping
and diagonal averaging to reconstruct elementary matrices X
i
.
First, the trajectory matrix for the desired signal is constructed
from one or more elementary matrices. The SVD splits the
time series into several groups I = (i
1
, i
2
,...,i
g
) with g being
the number of eigenvalues. The resultant trajectory matrix of
the desired signal is
ˆ
X
I
=
i
g
j=i
1
X
j
=
⎡
⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
ˆx
11
ˆx
12
ˆx
13
... ˆx
1K
ˆx
21
ˆx
22
ˆx
23
... ˆx
2K
ˆx
31
ˆx
32
ˆx
33
... ˆx
3K
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
ˆx
(L−1)1
ˆx
(L−1)2
ˆx
(L−1)3
... ˆx
(L−1)K
ˆx
L1
ˆx
L2
ˆx
L3
... ˆx
LK
⎤
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎦
L×K
.
The second step is diagonal averaging. Unlike X, the matrix
ˆ
X
I
is not a Hankel matrix. Hence, the antidiagonal elements
of
ˆ
X
I
are averaged to obtain an accurate estimation of each
sample. The reconstructed desired series is expressed as
ˆx
N
=[ˆx
1
, ˆx
2
, ˆx
3
,...,ˆx
N
],where
⎧
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎨
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎩
ˆx
1
= ˆx
11
ˆx
2
= ( ˆx
12
+ ˆx
21
)/2
.
.
.
ˆx
L−1
= ( ˆx
(L−1)1
+ ˆx
(L−2)2
+···+ ˆx
1(L−1)
)/(L − 1)
ˆx
L
= ( ˆx
L1
+ ˆx
(L−1)2
+···+ ˆx
1L
)/L
.
.
.
ˆx
K
= ( ˆx
L(N −2L+2)
+ ˆx
(L−1)(N −2L+3)
+···+ ˆx
1K
)/L
ˆx
K +1
= ( ˆx
L(N −2L+3)
+ ˆx
(L−1)(N −2L+4)
+···+ ˆx
2K
)/(L − 1)
.
.
.
ˆx
N−1
= ( ˆx
L(K −1)
+ ˆx
(L−1)K
)/2
ˆx
N
= ˆx
LK
.
Any sample can be estimated accurately only when an L
number of antidiagonal elements are averaged. For example,
ˆx
L−1
and ˆx
K +1
are calculated by averaging L − 1elements
Authorized licensed use limited to: Raja Ramanna Centre for Advanced Technology RRCAT. Downloaded on May 18,2021 at 19:08:37 UTC from IEEE Xplore. Restrictions apply.