11institutetext:  National-Regional Key Technology Engineering Laboratory for
Medical Ultrasound, School of Biomedical Engineering, Health Science Center,
Shenzhen University, China
11email<EMAIL_ADDRESS>22institutetext: Medical Ultrasound Image Computing
(MUSIC) Lab, Shenzhen University, China 33institutetext: Marshall Laboratory
of Biomedical Engineering, Shenzhen University, China 44institutetext:
Shenzhen RayShape Medical Technology Inc., Shenzhen, China

# Multi-IMU with Online Self-Consistency for Freehand 3D Ultrasound
Reconstruction

Mingyuan Luo Mingyuan Luo and Xin Yang contribute equally to this work.112233
Xin Yang0 112233 Zhongnuo Yan 112233 Junyu Li 112233 Yuanji Zhang 112233
Jiongquan Chen 112233 Xindi Hu 44 Jikuan Qian 44 Jun Cheng 112233 Dong Ni(✉)
112233

###### Abstract

Ultrasound (US) imaging is a popular tool in clinical diagnosis, offering
safety, repeatability, and real-time capabilities. Freehand 3D US is a
technique that provides a deeper understanding of scanned regions without
increasing complexity. However, estimating elevation displacement and
accumulation error remains challenging, making it difficult to infer the
relative position using images alone. The addition of external lightweight
sensors has been proposed to enhance reconstruction performance without adding
complexity, which has been shown to be beneficial. We propose a novel online
self-consistency network (OSCNet) using multiple inertial measurement units
(IMUs) to improve reconstruction performance. OSCNet utilizes a modal-level
self-supervised strategy to fuse multiple IMU information and reduce
differences between reconstruction results obtained from each IMU data.
Additionally, a sequence-level self-consistency strategy is proposed to
improve the hierarchical consistency of prediction results among the scanning
sequence and its sub-sequences. Experiments on large-scale arm and carotid
datasets with multiple scanning tactics demonstrate that our OSCNet
outperforms previous methods, achieving state-of-the-art reconstruction
performance.

###### Keywords:

Multiple IMU Online Learning Freehand 3D Ultrasound

Figure 1: Pipeline of freehand 3D US reconstruction with multiple lightweight
inertial measurement unit (IMU) sensors.

##  1 Introduction

Ultrasound (US) imaging has been widely used in clinical diagnosis due to its
advantages of safety, repeatability, and real-time imaging. Compared with 2D
US, 3D US can provide more comprehensive spatial information. Freehand 3D US
can enhance the understanding of physicians about the scanned region of
interest without increasing the complexity of scanning [13, 14, 12]. However,
the difficulty in estimating elevation displacement and accumulation error
makes it very challenging to infer the relative position only from images. In
this regard, it is expected to improve the reconstruction performance with the
help of external lightweight sensors, which will not significantly increase
the scanning complexity.

Sensorless freehand 3D US reconstructs the volume by calculating the relative
transformation of a series of US images. Previous studies were mainly based on
speckle decorrelation [1, 17], which estimates out-of-plane motion through the
correlation of speckle patterns in two successive frames. With the development
of deep learning technology, recent studies were mainly based on convolutional
neural network (CNN). Prevost et al. [15] proposed an end-to-end method based
on CNN to estimate the relative motion of US images. Guo et al. [4] proposed a
deep contextual learning network (DCL-Net) based on 3D CNN to estimate the
trajectory of US probe, and in a more recent study [3], they proposed a deep
contextual-contrastive network (DC2-Net), which introduced a contrastive
learning strategy to improve the reconstruction performance by leveraging the
label efficiently. Luo et al. [10, 11] proposed an online learning framework
(OLF) that improves reconstruction performance by online learning and shape
priors.

Due to the low cost, small size, and low power consumption of micro-electro-
mechanical-systems (MEMS), the sensor called inertial measurement unit (IMU)
has been widely used in navigation systems. Prevost et al. [14] incorporated
IMU orientation into neural network to improve reconstruction performance. Luo
et al. [12] proposed a deep motion network (MoNet) to mine the valuable
information of low signal-to-noise acceleration, and an online self-supervised
strategy was designed to further improve reconstruction performance. However,
the main disadvantage of IMU is that its measurement noise can not be
completely eliminated by calibration. Existing studies have shown that
combining multiple IMUs may help reduce noise and improve accuracy [2, 16].

In this study, we propose a multi-IMU-based online self-consistency network
(OSCNet) for freehand 3D US reconstruction. Our contribution is two-fold.
First, we equip multiple IMUs (see Fig. 1) to reduce the influence of noise in
individual IMU data. We propose a modal-level self-supervised strategy (MSS)
to fuse the information from multiple IMUs. MSS improves reconstruction
performance by reducing the differences between reconstruction results
obtained from each IMU data. Second, to reduce the estimation instability
caused by scanning differences such as frame rates, we propose a sequence-
level self-consistency strategy (SCS), which improves the hierarchical
consistency of prediction results among the scanning sequence and its sub-
sequences based on a consistent context. Experimental results show that the
proposed OSCNet can effectively fuse the information of multiple IMUs and
achieve state-of-the-art reconstruction performance.

##  2 Methodology

Figure 2: Overview of our proposed multi-IMU online self-consistency network
(OSCNet). IMU data diagrams ($U/\hat{U}$) show angle curves
($\Phi/\hat{\Phi}$, top) and acceleration curves ($A/\hat{A}$, bottom).
Relative transformation parameter diagrams ($\theta/\hat{\theta}$) show angle
curves ($\phi/\hat{\phi}$, top) and translation curves ($t/\hat{t}$, bottom).

Fig. 2 illustrates the proposed OSCNet, which consists of two essential
components: backbone and online learning. We construct a backbone using the
temporal and multi-branch structure from [12]. The main branch in the backbone
consists of ResNet [5] for feature extraction and LSTM [6] for processing
temporal information. It aids future estimation by leveraging temporal
contextual information. Additionally, there is an independent motion branch in
the backbone that fuses IMU information from a motion perspective with US
images. For more details, please refer to [12].

In the training phase, we input an $N$-length scanning sequence
$I=\\{I_{i}|i=1,2,\cdots,N\\}$ and corresponding multiple IMU data
$U=\\{U_{i}|i=1,2,\cdots,N-1\\}$ into the backbone to estimate the 3D relative
transformation parameters $\theta=\\{\theta_{i}|i=1,2,\cdots,N-1\\}$, where
$\theta_{i}$ includes 3-axis translations $t_{i}=(t_{x},t_{y},t_{z})_{i}$ and
rotation angles $\phi_{i}=(\phi_{x},\phi_{y},\phi_{z})_{i}$ between image
$I_{i}$ and $I_{i+1}$. The multiple IMU data consists of $M$ independent IMU
data $U_{i}=\\{U^{j}_{i}|j=1,2,\cdots,M\\}$, where $U^{j}_{i}$ consists of
3-axis angles $\Phi^{j}_{i}=(\Phi_{x},\Phi_{y},\Phi_{z})^{j}_{i}$ and
accelerations $A^{j}_{i}=(A_{x},A_{y},A_{z})^{j}_{i}$. The pre-processing
process for $\Phi_{i}$ and $A_{i}$ is consistent with [12]. Compared to
traditional offline inference strategies, online learning can leverage
valuable information from unlabeled data to improve the model’s performance
[10, 12]. In the testing phase, we propose two online self-supervised
strategies based on both the multiple IMU data (modal-level) and the scanning
sequence itself (sequence-level) to improve the performance of the backbone’s
estimations.

###  2.1 Modal-level Self-supervised Strategy

Multiple IMUs mounted in different directions provide diverse measurement
constraints for the model’s estimation, as shown in Fig. 1. This makes it
possible to reduce the influence of noise in individual IMU data while
adaptively optimizing for estimation. We construct an online modal-level self-
supervised strategy (MSS) that leverages the consistency between the
backbone’s estimation and multiple IMU data to improve the reconstruction
performance.

As shown in the top of Fig. 2, during the training phase, we repeatedly input
the US images and $M$ different IMU data into the backbone to obtain $M$
estimated parameters. We use the average of the $M$ estimated parameters
$\hat{\theta}=\frac{1}{M}\sum^{M}_{j=1}\hat{\theta}^{j}$ as the final output
of the backbone. We then calculate training loss between $\hat{\theta}$ and
ground truth $\theta$ using mean absolute error (MAE) and Pearson correlation
loss [4]:

$L=\|\hat{\theta}-\theta\|_{1}+(1-\frac{\textbf{Cov}(\hat{\theta},\theta)}{\sigma(\hat{\theta})\sigma(\theta)}),$
(1)

where Cov, $\sigma$ and $\left\|\cdot\right\|_{1}$ denote the covariance, the
standard deviation, and L1 normalization, respectively.

As shown in the bottom of Fig. 2, during the testing phase, we use each IMU
data $U^{j}$ ($j=1,2,\cdots,M$) as a weak label to constrain the corresponding
estimated parameters $\hat{\theta}^{j}$. We calculate the estimated
acceleration $\hat{A}^{j}$ at the center point of each image using the
estimated $\hat{\theta}^{j}$. To reduce the influence of acceleration noise,
we scale the $\hat{A}^{j}$ to match the mean-zeroed IMU acceleration.

$\hat{A}^{j}_{i}=((\hat{t}^{j}_{i-1})^{-1}+\hat{t}^{j}_{i})-\frac{1}{N-2}\sum_{i}((\hat{t}^{j}_{i-1})^{-1}+\hat{t}^{j}_{i}),\quad
i=2,3,\cdots,N-1,$ (2)

where $(\hat{t}^{j}_{i-1})^{-1}$ represents the translations in the inversion
of $\hat{\theta}^{j}_{i-1}$. Similar to [12], we use Pearson correlation loss
to measure the difference between the estimated and IMU acceleration, while
the angle is measured using MAE. Therefore, the single-IMU consistency
constraint between the estimated parameters and corresponding IMU data can be
expressed as:

$L_{single-
IMU}=\sum^{M}_{j=1}(1-\frac{\textbf{Cov}(\hat{A}^{j},A^{j})}{\sigma(\hat{A}^{j})\sigma(A^{j})})+\|\hat{\phi}^{j}-\Phi^{j}\|_{1}.$
(3)

In addition, the consistency among multiple IMU data itself also provides the
possibility to improve the reconstruction performance. It constrains the
backbone to obtain similar estimated parameters for different IMU data inputs
from the same scan. Specifically, we construct multi-IMU consistency
constraints as:

$L_{multi-
IMU}=\sum^{M}_{j,k=1,2,\cdots,M,j<k}(1-\frac{\textbf{Cov}(\hat{A}^{j},\hat{A}^{k})}{\sigma(\hat{A}^{j})\sigma(\hat{A}^{k})})+\|\hat{\phi}^{j}-\hat{\phi}^{k}\|_{1}.$
(4)

###  2.2 Sequence-level Self-consistency Strategy

Consistent context should lead to consistent parameter estimation, which
constrains the model at the sequence level, reducing the estimation
instability caused by scanning differences such as frame rates. Inspired by
contrastive learning [7], we construct an online sequence-level self-
consistency strategy (SCS). SCS randomly generates sub-sequences with
consistent context for each scan. The hierarchical consistency constraint
among the generated sub-sequences and the original sequence improves the
reconstruction performance of the backbone. Specifically, as shown in Fig. 2,
we randomly interval sample and flip each scanning sequence $I$ and its IMU
data $U$ to generate a sub-sequence $I_{sub}$ ($U_{sub}$) with consistent
context. In the testing phase, we obtain the estimated parameters
$\hat{\theta}_{sub}$ of $I_{sub}$ ($U_{sub}$) using the trained backbone. Then
compare $\hat{\theta}_{sub}$ with the original estimated parameters
$\hat{\theta}$ after the same interval sampling and flipping to construct the
self-consistency constraint:

$\begin{split}L_{self-
consistency}=&\|\hat{\theta}_{sub}-H_{\tau}(\hat{\theta})\|_{1}\\\
=&\|\frac{1}{M}\sum^{M}_{j=1}B(H_{\tau}(I),H_{\tau}(U^{j}))-H_{\tau}(\frac{1}{M}\sum^{M}_{j=1}B(I,U^{j}))\|_{1},\end{split}$
(5)

where $H_{\tau}$ converts the parameters, sequences, or IMU data under
interval sampling and flipping operation $\tau$. $B$ denotes the backbone.

##  3 Experiments

Figure 3: Comparison of multiple IMU data. The abscissa of each subfigure
indicates the image index.

##### Materials and Implementation.

The equipment we used to collect data includes a portable US machine, four IMU
sensors (WT901C-232, WitMotion) and an electromagnetic (EM) positioning
system. The US images were acquired with a linear probe at 10 MHz, and the
depth was set at 4 cm. As shown in Fig. 1, we bound four IMU sensors to the
probe in different orientations (three for 3-axis directions and one for
redundancy) using a 3D-printed bracket, which can compensate for errors and
reduce measurement singularities [9]. The resolutions of the IMU acceleration
and angle are $5\times 10^{-4}$ g/LSB and 0.5∘, respectively. We used the EM
positioning system to trace the scan route accurately. The direction and angle
resolutions of the EM positioning system are 1.4 mm and 0.5∘, respectively. We
calibrated the multiple IMU sensors and the EM positioning system using the
Levenberg–Marquardt algorithm [8] to ensure accurate measurements and minimise
system errors. As shown in Fig. 3, the calibrated IMU data exhibits a
generally consistent overall trend, although differences still exist.

We constructed two datasets, including arm and carotid, from 36 volunteers.
The arm dataset contains 288 scans, with scanning tactics including linear,
curved, loop, and sector scan. The carotid dataset contains 216 scans, with
scanning tactics including linear, loop, and sector scan. The average lengths
of the arm and carotid scans are 323.96 mm and 203.25 mm, respectively. The
size of scanned images is $248\times 260$ pixels, and the image spacing is
$0.15\times 0.15$ mm2. The collection and use of the data are approved by the
local IRB.

The arm and carotid datasets were randomly divided into 200/40/48 and
150/30/36 scans based on volunteer level to construct training/validation/test
set. To prevent overfitting and enhance the model’s robustness, we performed
random augmentations on each scan, including sub-sequence intercepting,
interval sampling, and sequence inversion. We randomly augmented each training
scan to 20 sequences and each test scan to 10 sequences to simulate complex
real-world situations. We used the Adam optimizer to optimize the OSCNet.
During the training phase, the epochs and batch size are set to 200 and 1,
respectively. To avoid overfitting, we set the initial learning rate to
$2\times 10^{-4}$ and used a learning rate decay strategy that halves the
learning rate every 30 epochs. During the online learning phase, the iteration
epoch and learning rate are set to 60 and $2\times 10^{-6}$, respectively. All
code was implemented in PyTorch and executed on an RTX 3090 GPU.

Table 1: The mean (std) results of different models on the arm and carotid
scans. DC2: DC2-Net, Bk: Backbone. The best results are shown in blue.

Models | FDR(%)$\downarrow$ | ADR(%)$\downarrow$ | MD(mm)$\downarrow$ | SD(mm)$\downarrow$ | HD(mm)$\downarrow$ | EA(deg)$\downarrow$
---|---|---|---|---|---|---
| Arm scans
CNN [14] | 23.31(13.0) | 32.54(13.7) | 67.79(30.0) | 2313.08(1852.5) | 62.48(31.5) | 4.35(1.8)
DC2 [3] | 14.02(7.3) | 26.15(10.2) | 45.50(21.3) | 1560.30(1181.4) | 42.25(40.9) | 4.71(2.3)
MoNet [12] | 11.58(6.2) | 20.35(6.8) | 32.21(11.0) | 1205.42(742.6) | 31.03(11.3) | 3.98(1.5)
Bk | 13.32(8.2) | 23.21(9.6) | 36.39(13.6) | 1339.17(822.4) | 34.91(13.7) | 4.32(1.7)
Bk+MSS | 10.78(5.6) | 19.53(6.3) | 30.52(10.5) | 1142.42(636.8) | 29.32(10.8) | 3.18(1.4)
Bk+SCS | 10.56(5.9) | 19.57(6.6) | 29.84(11.1) | 1126.28(614.9) | 28.64(11.6) | 3.65(1.9)
OSCNet | 10.01(5.7) | 18.86(6.5) | 28.61(11.0) | 1064.06(582.5) | 27.38(11.4) | 2.76(1.3)
| Carotid scans
CNN [14] | 25.85(15.0) | 33.95(16.8) | 49.64(25.5) | 1944.72(1485.1) | 39.30(18.7) | 3.73(2.3)
DC2 [3] | 13.54(7.1) | 21.68(9.2) | 26.47(9.6) | 1025.06(622.8) | 24.49(10.3) | 4.30(3.0)
MoNet [12] | 11.80(5.7) | 20.42(8.8) | 23.48(8.7) | 894.39(381.0) | 20.78(9.2) | 3.67(2.1)
Bk | 12.85(6.5) | 21.78(10.5) | 24.65(9.1) | 965.12(466.6) | 21.81(9.5) | 3.83(2.0)
Bk+MSS | 11.31(5.4) | 20.04(8.8) | 22.72(8.1) | 850.68(321.7) | 20.02(8.5) | 3.16(1.8)
Bk+SCS | 11.30(5.4) | 20.16(8.6) | 23.01(8.4) | 863.48(320.9) | 20.56(8.7) | 3.36(1.8)
OSCNet | 10.90(5.3) | 19.61(8.5) | 21.81(7.2) | 804.27(282.8) | 19.30(7.6) | 2.60(1.6)

Figure 4: Metric decline curves (with 95% confidence interval). Red: MoNet,
Blue: OSCNet. The abscissa and ordinate of each subfigure represent the number
of iterations and the value of metrics, respectively.

##### Quantitative and Qualitative Analysis.

To demonstrate the effectiveness of our OSCNet, we compared it with three
state-of-the-art methods, including CNN [14], DC2-Net [3] and MoNet [12]. All
comparison methods were trained to convergence using the experimental settings
given in the corresponding papers. We adopt six metrics from [12] to evaluate
reconstruction performance: final drift rate (FDR), average drift rate (ADR),
maximum drift (MD), sum of drift (SD), symmetric Hausdorff distance (HD), and
mean error of angle (EA). In addition, ablation experiments are conducted to
validate the effectiveness of MSS and SCS as proposed in our OSCNet.

Table 1 shows that our OSCNet significantly outperforms CNN, DC2-Net, MoNet,
and our Backbone in all metrics for both arm and carotid scans ($t$-test,
$p<0.05$), except for MoNet’s ADR on the carotid scans ($t$-test, $p=0.10$).
Notably, sensor-based methods (MoNet and OSCNet) have exhibited improvements
in all metrics compared to sensorless methods (DC2-Net and CNN). The multi-
IMU-based OSCNet outperforms the single-IMU-based MoNet, verifying the
effectiveness of multiple IMU integration. Moreover, the ablation experiments
further demonstrate that both multiple IMU integration (MSS) and self-
consistency constraint (SCS) greatly improve the reconstruction performance.

Figure 5: Typical reconstruction cases. Top: comparison of different methods,
Bottom: comparison of different scanning tactics. At the bottom, all of the
estimated image positions of OSCNet are marked with red boxes to visualize the
scanning tactics.

In addition, Fig. 4 displays the metric decline curves during the online
learning phase of MoNet and OSCNet on the arm and carotid datasets. All metric
curves exhibit a decreasing trend followed by stabilization. We note that our
OSCNet has achieved further improvements compared to MoNet, with 13.56%
/7.32%/30.65% and 7.62%/4.00%/29.16% improvement in FDR/ADR/EA on the arm and
carotid datasets, respectively. Fig. 5 presents several typical reconstruction
results of all methods. It can be observed that our OSCNet outperforms other
methods in reconstruction results and closely approximates the ground truth
across all scanning tactics.

##  4 Conclusion

In this study, we propose a novel multi-IMU-based online self-consistency
network (OSCNet) to conduct freehand 3D US reconstruction. We propose an
online modal-level self-supervised strategy (MSS) that integrates multiple
IMUs to reduce the influence of single IMU noise and enhance reconstruction
performance. We propose an online sequence-level self-consistency strategy
(SCS) to improve the reconstruction stability using hierarchical consistency
among the generated sub-sequences and the original sequence. The experimental
results on the arm and carotid datasets show that our OSCNet achieves state-
of-the-art reconstruction performance. Future research will focus on exploring
more general reconstruction methods.

#### Acknowledgements

This work was supported by the grant from National Natural Science Foundation
of China (Nos. 62171290, 62101343), Shenzhen-Hong Kong Joint Research Program
(No. SGDX20201103095613036), and Shenzhen Science and Technology Innovations
Committee (No. 20200812143441001).

## References

  * [1] Chen, J.F., Fowlkes, J.B., Carson, P.L., Rubin, J.M.: Determination of scan-plane motion using speckle decorrelation: Theoretical considerations and initial test. International Journal of Imaging Systems and Technology 8(1), 38–44 (1997)
  * [2] Guerrier, S.: Improving accuracy with multiple sensors: Study of redundant mems-imu/gps configurations. In: Proceedings of the 22nd international technical meeting of the Satellite Division of the Institute of Navigation (ION GNSS 2009). pp. 3114–3121 (2009)
  * [3] Guo, H., Chao, H., Xu, S., Wood, B.J., Wang, J., Yan, P.: Ultrasound volume reconstruction from freehand scans without tracking. IEEE Transactions on Biomedical Engineering 70(3), 970–979 (2023)
  * [4] Guo, H., Xu, S., Wood, B., Yan, P.: Sensorless freehand 3d ultrasound reconstruction via deep contextual learning. In: International Conference on Medical Image Computing and Computer-Assisted Intervention. pp. 463–472. Springer (2020)
  * [5] He, K., Zhang, X., Ren, S., Sun, J.: Deep residual learning for image recognition. In: 2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR). pp. 770–778. IEEE (2016)
  * [6] Hochreiter, S., Schmidhuber, J.: Long short-term memory. Neural Computation 9(8), 1735–1780 (1997)
  * [7] Le-Khac, P.H., Healy, G., Smeaton, A.F.: Contrastive representation learning: A framework and review. IEEE Access 8, 193907–193934 (2020)
  * [8] Levenberg, K.: A method for the solution of certain non-linear problems in least squares. Quarterly of applied mathematics 2(2), 164–168 (1944)
  * [9] Liang, S., Dong, X., Guo, T., Zhao, F., Zhang, Y.: Peripheral-free calibration method for redundant imus based on array-based consumer-grade mems information fusion. Micromachines 13(8) (2022)
  * [10] Luo, M., Yang, X., Huang, X., Huang, Y., Zou, Y., Hu, X., Ravikumar, N., Frangi, A.F., Ni, D.: Self context and shape prior for sensorless freehand 3d ultrasound reconstruction. In: Medical Image Computing and Computer Assisted Intervention. pp. 201–210 (2021)
  * [11] Luo, M., Yang, X., Wang, H., Dou, H., Hu, X., Huang, Y., Ravikumar, N., Xu, S., Zhang, Y., Xiong, Y., Xue, W., Frangi, A.F., Ni, D., Sun, L.: Recon: Online learning for sensorless freehand 3d ultrasound reconstruction. Medical Image Analysis 87, 102810 (2023)
  * [12] Luo, M., Yang, X., Wang, H., Du, L., Ni, D.: Deep motion network for freehand 3d ultrasound reconstruction. In: Medical Image Computing and Computer Assisted Intervention. pp. 290–299 (2022)
  * [13] Mohamed, F., Siang, C.V.: A survey on 3d ultrasound reconstruction techniques. In: Aceves-Fernandez, M.A. (ed.) Artificial Intelligence, chap. 4. IntechOpen, Rijeka (2019)
  * [14] Prevost, R., Salehi, M., Jagoda, S., Kumar, N., Sprung, J., Ladikos, A., Bauer, R., Zettinig, O., Wein, W.: 3d freehand ultrasound without external tracking using deep learning. Medical Image Analysis 48, 187–202 (2018)
  * [15] Prevost, R., Salehi, M., Sprung, J., Ladikos, A., Bauer, R., Wein, W.: Deep learning for sensorless 3d freehand ultrasound imaging. In: Medical Image Computing and Computer-Assisted Intervention. pp. 628–636 (2017)
  * [16] Rasoulzadeh, R., Shahri, A.M.: Implementation of a low-cost multi-imu hardware by using a homogenous multi-sensor fusion. In: 2016 4th International Conference on Control, Instrumentation, and Automation (ICCIA). pp. 451–456 (2016)
  * [17] Tuthill, T.A., Krücker, J., Fowlkes, J.B., Carson, P.L.: Automated three-dimensional us frame positioning computed from elevational speckle decorrelation. Radiology 209(2), 575–582 (1998)


# Generalizing the Lehmer’s totient problem

Marius Tărnăuceanu

(October 18, 2021)

###### Abstract

An important unsolved question in number theory is the Lehmer’s totient
problem that asks whether there exists any composite number $n$ such that
$\varphi(n)\mid n-1$, where $\varphi$ is the Euler’s totient function. It is
known that if any such $n$ exists, it must be odd, square-free, greater that
$10^{30}$, and divisible by at least $15$ distinct primes. Such a number must
be also a Carmichael number.

In this short note, we discuss a group-theoretical analogous problem involving
the function that counts the number of automorphisms of a finite group.
Another way to generalize the Lehmer’s totient problem is also proposed.

MSC 2020 : Primary 20D60, 11A25; Secondary 20D99, 11A99.

Key words : Lehmer’s totient problem, Euler’s totient function, finite group,
group automorphism, exponent of a group.

##  1 Introduction

The Euler’s totient function $\varphi$ is one of the most famous functions in
number theory. Recall that the totient $\varphi(n)$ of a positive integer $n$
is defined to be the number of positive integers less than or equal to $n$
that are coprime to $n$. In algebra this function is important mainly because
it gives the order of the group of units in the ring
$(\mathbb{Z}_{n},+,\cdot)$. Also, $\varphi(n)$ can be seen as the number of
generators or as the number of automorphisms of the cyclic group
$(\mathbb{Z}_{n},+)$.

Lehmer’s totient problem [6] asks whether the well-known property
”$\varphi(n)=n-1\Leftrightarrow n$ is a prime” can be generalized to
”$\varphi(n)\mid n-1\Leftrightarrow n$ is a prime”. This problem has been
studied by many mathematicians (see e.g. [2, 3, 4, 6, 7]), but up to now no
counterexample has been found. Such a counterexample is often called a Lehmer
number.

We observe that an integer $n\geq 2$ is a prime or a Lehmer number if and only
if

$|G|-1\equiv 0\,\,({\rm mod}\,\,|{\rm Aut}(G)|),$ (1)

where $G$ a cyclic group of order $n$. This fact suggests us to consider
arbitrary finite groups $G$ which satisfy the relation (1). Their description
is given by the following theorem.

######  Theorem 1.1.

A finite group $G$ satisfies the relation $(1)$ if and only if it is cyclic
and its order is a prime or a Lehmer number.

Finally, we indicate another way to extend the Lehmer’s totient problem via
group theory.

Open problem. Determine the finite groups $G$ satisfying

$|G|-1\equiv 0\,\,({\rm mod}\,\,\varphi(G)),$ (2)

where

$\varphi(G)=|\\{a\in G\mid o(a)=\exp(G)\\}|$

is the generalization of the Euler’s totient function studied in [8].

Note that since $\varphi(\mathbb{Z}_{n})=\varphi(n)$ for all
$n\in\mathbb{N}^{*}$, cyclic groups of prime or Lehmer order will be solutions
of (2).

##  2 Proof of Theorem 1.1

First of all, we recall the well-known formula for the number of automorphisms
of a finite abelian $p$-group (see e.g. [1, 5]).

######  Theorem 2.1.

Let $G=\prod_{i=1}^{k}\mathbb{Z}_{p^{n_{i}}}$ be a finite abelian$p$-group,
where $1\leq n_{1}\leq n_{2}\leq...\leq n_{k}$. Then

$|{\rm
Aut}(G)|=\prod_{i=1}^{k}(p^{a_{i}}-p^{i-1})\prod_{u=1}^{k}p^{n_{u}(k-a_{u})}\prod_{v=1}^{k}p^{(n_{v}-1)(k-b_{v}+1)}\,,$
(3)

where

$a_{r}=max\\{s\mid n_{s}=n_{r}\\}\mbox{ and }b_{r}=min\\{s\mid
n_{s}=n_{r}\\}\,,\,r=1,2,...,k\,.$

By using Theorem 2.1, we easily get the following corollary.

######  Corollary 2.2.

Let $G$ be a finite abelian $p$-group. If $p\nmid|{\rm Aut}(G)|$, then
$G\cong\mathbb{Z}_{p}$.

###### Proof.

Under the notation in Theorem 2.1, we infer that $n_{1}=n_{2}=\dots=n_{k}=1$
and $a_{1}=a_{2}=\dots=a_{k}=k$. Then the first product in the right side of
(3) is $\prod_{i=1}^{k}(p^{k}-p^{i-1})$. Obviously, this is not divisible by
$p$ if and only if $k=1$. Thus $G\cong\mathbb{Z}_{p}$, as desired. ∎

We are now able to prove our main result.

Proof of Theorem 1.1. Let $G$ be a finite group satisfying (1).

We first prove that $G$ is abelian. If not, then $Z(G)\neq G$ and so there
exists a prime $p$ dividing $|G/Z(G)|$. Since $G/Z(G)$ can be embedded in
${\rm Aut}(G)$, it follows that $p$ divides $|{\rm Aut}(G)|$. Consequently,
$p\mid|G|-1$, contradicting the fact that $p\mid|G|$. Thus $G$ is abelian.

Let $G=\prod_{i=1}^{m}G_{i}$, where $G_{i}$ is a finite abelian $p_{i}$-group,
$i=1,2,\dots,m$. Since

$|{\rm Aut}(G)|=\prod_{i=1}^{m}|{\rm Aut}(G_{i})|\mbox{ and
}|G|=\prod_{i=1}^{m}|G_{i}|,$

by (1) we infer that $p_{i}\nmid|{\rm Aut}(G_{i})|$ for each $i$. Then
Corollary 2.2 implies $G_{i}\cong\mathbb{Z}_{p_{i}}$ and therefore

$G\cong\prod_{i=1}^{m}\mathbb{Z}_{p_{i}}\cong\mathbb{Z}_{p_{1}p_{2}\cdots
p_{m}}$

is cyclic. Moreover, (1) becomes $|G|-1\equiv 0\,\,({\rm
mod}\,\,\varphi(|G|))$, i.e. $|G|$ is a prime or a Lehmer number. This
completes the proof.∎

## References

  * [1] J.N.S. Bidwell, M.J. Curran and D.J. McCaughan, Automorphisms of direct products of finite groups, Arch. Math. 86 (2006), 481-489.
  * [2] P. Burcsi, S. Czirbusz and G. Farkas, Computational investigation of Lehmer’s totient problem, Annales Univ. Sci. Budapest., Sect. Comp. 35 (2011), 43–49.
  * [3] A. Grytczuk and M. Wójtowicz, On a Lehmer problem concerning Euler’s totient function, Proc. Japan Acad. Ser. A Math. Sci. 79 (2003), 136-138.
  * [4] P. Hagis, On the equation $M\varphi(n)=n-1$, Nieuw Arch. Wiskd. IV Series. 6 (1988), 255-261.
  * [5] C.J. Hillar and D.L. Rhea, Automorphisms of finite abelian groups, Amer. Math. Monthly 114 (2007), 917-923.
  * [6] D.H. Lehmer, On Euler’s totient function, Bull. Amer. Math. Soc. 38 (1932), 745-751.
  * [7] F. Luca and C. Pomerance, On composite integers $n$ for which $\varphi(n)\mid n-1$, Bol. Soc. Mat. Mexicana. 17 (2011), 13-21.
  * [8] M. Tărnăuceanu, A generalization of the Euler’s totient function, Asian-Eur. J. Math. 8 (2015), article ID 1550087.

Marius Tărnăuceanu

Faculty of Mathematics

“Al.I. Cuza” University

Iaşi, Romania

e-mail<EMAIL_ADDRESS>


# A Novel Approach to Model the Kinematics of Human Fingers Based on an
Elliptic Multi-Joint Configuration

Zeyu Wu, Luiza Labazanova, Peng Zhou, and David Navarro-Alarcon This work is
supported in part by: The Research Grants Council of Hong Kong (grants
14203917, F-PolyU503/18), the Jiangsu Industrial Technology Research Institute
Collaborative Research Program Scheme (grant ZG9V), the 2019/20 Belt and Road
Scholarship (Research Postgraduate), and PolyU (grants G-YBYT, 4-ZZHJ).
Corresponding author: Luiza Labazanova.All authors are with the Department of
Mechanical Engineering, The Hong Kong Polytechnic University (PolyU), Kowloon,
Hong Kong. (e-mail<EMAIL_ADDRESS>
###### Abstract

In this paper, we present a novel kinematic model of the human phalanges based
on the elliptical motion of their joints. The presence of the soft elastic
tissues and the general anatomical structure of the hand joints highly affect
the relative movement of the bones. Commonly used assumption of circular
trajectories simplifies the designing process but leads to divergence with the
actual hand behavior. The advantages of the proposed model are demonstrated
through the comparison with the conventional revolute joint model. Conducted
simulations and experiments validate designed forward and inverse kinematic
algorithms. Obtained results show a high performance of the model in mimicking
the human fingertip motion trajectory.

##  I INTRODUCTION

The human hand is an essential tool that is used to interact with the
surroundings. Hand’s dexterity and precision make it highly effective in
numerous grasping, manipulation and communication tasks that occur in people
routines every day. Researchers have been trying to implement human hand
abilities in humanoid robots that demonstrated advanced development in
locomotion but lacked upper limbs dexterity. There are multiple examples of
bionic robotic hands (BRH) aiming to mimic the hand’s shape and general
structure. However, their performance is far from their counterpart skills
capacity. We suppose that deeper investigation and understanding of the arm’s
anatomy and functionality is required in order to fulfil the gap between the
biological and artificial hands.

Existing BRHs can be divided into two types according to the structure of the
joint: Fixed Rotational Joint (FRJ) and Contact Rotational Joint (CRJ). In
FRJ, bionic bones are linked by a rigid shaft; thus, only pure rotation occurs
between them. Such a joint has one degree of freedom (DoF), which determines
its common use in bionic applications [1, 2, 3]. In CRJ, the articulation
surfaces of the bones are connected by soft tissues. Therefore, the bones are
not strictly constrained relative to each other, enabling slight extension or
compression of the joint [4, 5].

Both FRJ and CRJ employ the joints that define the circular motion of the
output bones. However, there are ongoing disputes about the trajectories
generated by the joints within the human hand. In early studies, [6, 7, 8],
the common assumption was that a circular arc is a proper approximation,
whereas more recent anatomical investigations revealed that elliptical paths
of the joints are significant to be neglected [9, 10, 11]. Therefore, in this
work, we suggest developing a general kinematic model that will describe the
rotation and translation of the joints with non-circular trajectories.

There are various studies striving to implement elliptical motion on the other
human joints, such as the shoulder [12] and knee [13]. These are examples of
single joints while human fingers can be considered as kinematic chains and,
thus, their modelling is more complicated. One of the attempts to solve this
issue is presented in [14]. However, this work presumes that ellipsoids are in
contact with each other. On the other hand, the biological finger joints are
composed of one bone with a convex (proximal) head and another one with a
concave (distal) head, which are separated by synovial fluid.

In this paper, we aim to present a novel kinematic model of the human finger
based on the anatomical features of its joints (see Fig. 1). A multi-joint
finger model is designed with the assumption that joints exhibit an elliptic
arc motion. The model was further adjusted by curve-fitting with data
collected during the movement of the human finger. The main contributions of
this work are:

  * •

Analytical design of the finger joint that enables elliptical trajectories.

  * •

Solution of the forward and inverse kinematics equations that describe the
motion of the multi-finger system.

  * •

Simulations and experiments that evaluate the performance of the proposed
model.

Figure 1: The structure of bones of human right hand [15]

##  II MODELING OF A BIOMIMETIC HYBRID ROBOTIC FINGER

###  II-A The Kinematic Model of an Elliptic Joint

This section introduces the designing approach of a kinematic model of the
joints with non-circular motion.

Figure 2: The trajectory of elliptical joint: (a) The workspace of MCP joint.
(b) the trajectory of joint on xz plane and yz plane.

Fig. 2(a) shows the workspace of a metacarpophalangeal (MCP) joint, which is
formed by elliptic trajectories ${{c}_{0}}$ (adduction/abduction) and
${{c}_{1}}$ (flexion/extension). According to the geometric relationships in
Fig. 2(b), the adduction and abduction (Ad/Ab) movements occur on the
${{x}_{i-1}}{{z}_{i-1}}$ plane, and the flexion and extension movements act on
the ${{y}_{i-1}}{{z}_{i-1}}$ plane. We define ${{o}_{i}^{i-1}}$ as the origin
point of the coordinate frame $i$ under coordinate frame $i-1$. Additionally,
${{o}_{i}^{i-1}}$ represents the translation from coordinate frame $i-1$ to
$i$, that is
${{o}_{i}^{i-1}}={{\left[\begin{matrix}{{x}_{i-1}}&{{y}_{i-1}}&{{z}_{i-1}}\\\
\end{matrix}\right]}^{T}}$ $i\in\left[1...4\right]$, which correspond the
three joints of the finger. MCP joint can be considered as a combination of
the two joints: flexion/extension (F/E) joint and adduction/abduction (Ad/Ab)
joint. Thus, the origin point can be determined as follows.

On the ${{x}_{i-1}}{{z}_{i-1}}$ plane:

$o_{i}^{i-1}({{\theta}_{i}})={{\left[\begin{matrix}0&{{z}_{i-1}}\tan({{\theta}_{i}})&{{l}_{i-1}}{{s}_{i-1}}\sqrt{\frac{1}{{{\tan}^{2}}{{\theta}_{i}}l_{i-1}^{2}+s_{i-1}^{2}}}\\\
\end{matrix}\text{ }\right]}^{T}}.$ (1)

On the ${{y}_{i-1}}{{z}_{i-1}}$ plane:

$o_{i}^{i-1}({{\theta}_{i}})\text{=}{{\left[\begin{matrix}{{z}_{i-1}}\tan({{\theta}_{i}})&0&{{l}_{i-1}}{{s}_{i-1}}\sqrt{\frac{1}{{{\tan}^{2}}{{\theta}_{i}}s_{i-1}^{2}+l_{i-1}^{2}}}\\\
\end{matrix}\text{ }\right]}^{T}},$ (2)

where ${{l}_{i-1}}$ and ${{s}_{i-1}}$ are the major and minor axes of the $i$
joint, respectively. We assume that the major axis is on ${{x}_{i-1}}$ and
${{y}_{i-1}}$ according to the anatomy geometry. ${{\theta}_{i}}$ is the
rotation angle of the joint, on ${{x}_{i-1}}{{z}_{i-1}}$ plane,
${{\theta}_{i}}\in(0,\pi/2)$, on ${{y}_{i-1}}{{z}_{i-1}}$ plane,
${{\theta}_{i}}\in(-\pi/2,\pi/2)$.

According to Fig.2(b), ${{\varphi}_{i}}$ is the acute angle formed by the $i$
bone extension line and the horizontal axis, which is the rotation angle of
coordinate frame $i-1$ to $i$. It is worth noting that ${{\varphi}_{i}}$ is a
function of ${{\theta}_{i}}$, they are not equal except the major axis and
minor axis of ellipse are same (at this time, the trajectory is a circle).
Their relationship is:

$\displaystyle\text{on }{{x}_{i-1}}{{z}_{i-1}}~{}\text{plane:
}{{\varphi}_{i}}={{\tan}^{-1}}\left(\frac{s_{i-1}^{2}}{l_{i-1}^{2}}\tan({{\theta}_{i}})\right)$
(3) $\displaystyle\text{on }{{y}_{i-1}}{{z}_{i-1}}~{}\text{plane:
}{{\varphi}_{i}}={{\tan}^{-1}}\left(\frac{l_{i-1}^{2}}{s_{i-1}^{2}}\tan({{\theta}_{i}})\right)$
(4)

The trajectory of the bone has both translation and rotation. We assume a
point on coordinate frame $i-1$ is
${{q}_{i-1}}={{\left[\begin{matrix}{{q}_{{{x}_{i-1}}}}&{{q}_{{{y}_{i-1}}}}&{{q}_{{{z}_{i-1}}}}\\\
\end{matrix}\right]}^{T}}$. It can be expressed as:

$\displaystyle{{\overset{\\_}{\mathop{q}}\,}_{i-1}}=\left[\begin{matrix}{{q}_{i-1}}\\\
1\\\
\end{matrix}\right]=\left[\begin{matrix}R_{i-1}^{i}({{\varphi}_{i}})&p_{i-1}^{i}({{\theta}_{i}})\\\
0&1\\\ \end{matrix}\right]\left[\begin{matrix}{{q}_{i}}\\\ 1\\\
\end{matrix}\right]=\overline{T}_{i-1}^{i}{{\overset{\\_}{\mathop{q}}\,}_{i}}$
(5)

${\overset{\\_}{\mathop{q}}_{i-1}\,}$ is the homogeneous coordinate of
${{q}_{i}}$. $\text{ }R_{i-1}^{i}({{\varphi}_{i}})\in SO(3)\text{ }$ is the
rotation matrix from $i-1$ to $i$. It’s a function of ${{\varphi}_{i}}$.
$p_{i-1}^{i}({{\theta}_{i}})\in{{\mathbb{R}}^{3}}$ is the translation from
coordinate frame $i-1$ to $i$, which is equal to
${{o}_{i}^{i-1}({{\theta}_{i}})}$. ${T}_{i-1}^{i}$ is the transformation
matrix from $i-1$ to $i$.

###  II-B Forward Kinematics of the Finger

Figure 3: The comparison of circular and elliptical joint

In order to develop the forward kinematics equation of the multi-joints
finger, we need to get its transformation matrix. The transformation matrix of
a joint with elliptical trajectory (elliptical joint) is different from that
of a joint with circular trajectory (circular joint). Generally, for getting
the transformation matrix, we need to obtain the rotation angle of the
coordinate frame ($\varphi_{i}$) of each joint relative to the previous joint,
which we name it as the rotation angle of the coordinate frame. In circular
joint, we consider that it equals to the rotation angle of the joint
($\theta_{i}$), which is the angle formed by the intersection of the dashed
line and the vertical axis (Fig.3). Then we input this angle into the joint
motor or use it to calculate the displacement of the tendon in tendon-driven
system. However, the transformation matrix of elliptical joint is different.
It’s rotation angle of joint is not equal to rotation angle of coordinate
frame. This is because the trajectory of its motion is not only pure rotation,
but also having the translation in plane, and the translation is a function of
the rotation angle of joint. This means that the translation will change with
the rotation of the joint. This is the most significant difference between the
elliptical and circular joint. We will introduce how we develop the forward
kinematics model of the multi-joints finger with elliptical trajectory joints
as follows.

Figure 4: The coordinate frame and rotational axes of finger. (a) The
rotational axes (b) The coordinate frame from 0-3.

Fig.4 (b) shows that the coordinate frames of the finger. There are three
condyloid joints with four movements. We choose coordinate frame $0$
represents the Ad/Ab movement, frame $1$ as the F/E of MCP joint, frame $2$
and $3$ as the F/E of PIP and DIP joints. ${l}_{i-1}$ is the length of bone.
As equation (5), we can calculate $\overline{T}_{i-1}^{i}$, the first
transformation matrix is on ${{x}_{i-1}}{{z}_{i-1}}$ plane, the rest are on on
${{y}_{i-1}}{{z}_{i-1}}$ plane. So we have:

$\displaystyle\overline{T}_{0}^{1}={{\left[\begin{matrix}R_{0}^{1}({{\varphi}_{1}})&o_{0}^{1}({{\theta}_{1}})\\\
0&1\\\ \end{matrix}\right]}_{4\times 4}}\text{ }$ (6)
$\displaystyle\overline{T}_{i-1}^{i}={{\left[\begin{matrix}R_{i-1}^{i}({{\varphi}_{i}})&o_{i-1}^{i}({{\theta}_{i}})+\left[\begin{matrix}0\\\
0\\\ {{l}_{i-1}}\\\ \end{matrix}\right]\\\ 0&1\\\
\end{matrix}\right]}_{4\times 4}}\text{(i=2}...\text{4) }$ (7)

After we have all transformation matrix, we can get the forward kinematics
expression:

${{\overline{q}}_{0}}=(\prod\limits_{i=1}^{4}{\overline{T}_{i-1}^{i}){{\overline{q}}_{4}}}\text{
i}\in\text{ }\\!\\![\\!\\!\text{ 1,4 }\\!\\!]\\!\\!\text{ }$ (8)

${{\overline{q}}_{0}}$ and ${{\overline{q}}_{4}}$ are the fingertip
homogeneous coordinate in coordinate frame $0$ and $4$,
${{\overline{q}}_{4}}=[0,0,0,1]^{T}$.

###  II-C Inverse Kinematics of the Finger

Figure 5: The geometry of MCP Ad/Ab trajectory and fingertip. (a) The
projection of finger model on ${{x}_{0}}{{z}_{0}}$ plane. (b) The posture of
finger model. (c) The Ad/Ab movement of human finger. (d) The lateral side of
finger model.

We have applied three methods of inverse kinematics. The first method refers
the work done by Jiewen Lai et al. [16], firstly input random rotation angle
of joints and get the fingertip position by forward kinematics, then save all
the data as Point Cloud Library (PCL) with fingertip position and rotation
angle of joints. When doing inverse kinematics, we need to search the
approximate fingertip position in PCL and return rotation angles.

The second method is to get the analytical solution of inverse kinematics. We
firstly solve the rotation angle of Ad/Ab by using it’s physical properties
(Fig.5). During the movement of finger (except thumb), Ad/Ab movement occurs
on ${{x}_{i-1}}{{z}_{i-1}}$ plane and it’s the only movement on this plane.
Because we have set it on first coordinate frame, we can directly solve it
according to fingertip position ${{\overline{f}}}=[f_{x},f_{y},f_{z}]^{T}$.
${{\theta}_{1}}$ can be expressed as:

$F_{1}({{\theta}_{1}})={{\left[\begin{matrix}{{f}_{x}}\\\ {{f}_{y}}\\\
{{f}_{z}}\\\ 1\\\ \end{matrix}\right]}^{T}}\left[\begin{matrix}1\\\
-\tan({{\varphi}_{1}}({{\theta}_{1}}))\\\ 0\\\
{{z}_{0}}({{\theta}_{1}})(\tan({{\varphi}_{1}}({{\theta}_{1}}))-\tan({{\theta}_{1}}))\\\
\end{matrix}\right]$ (9)

${{z}_{0}}({{\theta}_{1}})$ is the ${{Z}_{0}}$ axis coordinate of $o_{1}^{0}$,
which is a function of ${{\theta}_{1}}$. $F_{1}({{\theta}_{1}})$ equal to zero
and has complex expression which is difficult to solve directly. We use
Newton-Raphson method to get the root ($F_{1}({{\theta}_{1}})=0$) of it. The
iteration process is:

$\Delta{{\theta}_{1}}={{J}^{T}}{{(J{{J}^{T}})}^{-1}}\overset{\to}{\mathop{e}}\,={{J}^{\dagger}}\overset{\to}{\mathop{e}}\
$ (10)

$\Delta{{\theta}_{1}}$ is the change value of ${{\theta}_{1}}$ for each
iteration ($\vartriangle{{\theta}_{1}}=\theta_{1}^{j+1}-\theta_{1}^{j}$),
${j}$ is the increment of iteration. $J$ is the Jacobian matrix,
$J=\frac{\partial F({{\theta}_{1}})}{\partial{{\theta}_{1}}}$. ${J}^{\dagger}$
is the Jacobian Pseudo-inverse and
${\overset{\to}{\mathop{e}}}=F({{\theta}_{1}})-0$.

After calculating ${{\theta}_{1}}$, we need to solve the inverse kinematics of
F/E of three joints. From equation (9), we can get:

${{F}_{2}}(\overline{\theta})=\prod\limits_{i=2}^{4}{T_{i-1}^{i}-{{(\overline{T}_{0}^{1})}^{-1}}{{\overline{q}}_{0}}{{({{\overline{q}}_{4}})}^{-1}}\text{
}i\in[2,4]}$ (11)

${{F}_{2}}$ is function of $\overline{\theta}$ and equals to zero.
$\overline{\theta}$ is a column vector:
${{\left[\begin{matrix}{{\theta}_{2}}&{{\theta}_{3}}&{{\theta}_{4}}\\\
\end{matrix}\right]}^{T}}$. This three movements are all in sagittal plane
which can be regarded as the plane motion, making it have infinite solutions.
In order to determine the finger configuration in space, we add a constraint
according to anatomical structure. In human finger, the movement of DIP and
PIP joint are not independent and have constant ratio when doing flexion
without load [11]. So we assume ${{\theta}_{3}}=\alpha{{\theta}_{4}}$,
$\alpha$ is the constant ratio. We can get the analytical solution of inverse
kinematics from $(12)$, which contain three equations. However, they are
complex implicit equation and difficult to be solved. So we use the numerical
method again to solve this equation. The equation can be expressed as:

$\Delta\overline{\theta}={{\hat{J}}^{T}}{{(\hat{J}{{\hat{J}}^{T}})}^{-1}}\overset{\to}{\mathop{{\hat{e}}}}\,={{\hat{J}}^{\dagger}}\overset{\to}{\mathop{{\hat{e}}}}\,$
(12)

$\hat{J}$ is the Jacobian matrix of ${{F}_{2}}(\overline{\theta})$, which is:

$J{{(\overline{\theta})}_{j}}=\text{
}{{\left(\frac{\partial{{F}_{2}}{{(\overline{\theta})}_{i}}}{\partial({{\overline{\theta}}_{j}})}\right)}_{ij}}\text{
}i,j=[1...3]\text{ }$ (13)

${{F}_{2}}{{(\overline{\theta})}_{i}}$ means the $i_{th}$ equation from
$(12)$, ${{\overline{\theta}}_{j}}$ equals to ${{\theta}_{j+1}}$.
${\hat{J}}^{\dagger}$ represents the Jacobin Pseudo-inverse and
${\hat{{e}}}\,={{F}_{2}}(\overline{\theta})-0$.

The third method is for getting the numerical solution of inverse kinematics
with multiple constraints. In one finger case, it’s similar like the second
method. For Ad/Ab movement of MCP joint, we still use geometry method to solve
it. For the rest three F/E movements, we use optimization-based approach of
numerical inverse kinematics solver as follow:

$\hat{\theta}=\arg\min\left({{\left\|\overline{f}-Fkine\left({\hat{\theta}}\right)\right\|}^{2}}+{{\left\|\Delta\hat{\theta}\right\|}^{2}}+r{{\left\|{{\theta}_{4}}-\alpha{{\theta}_{3}}\right\|}^{2}}\right)$
(14)

${\hat{\theta}}\,\in{{\mathbb{R}}^{4}}$ equals to
${{\left[{{\theta}_{1}}\text{ }{{\theta}_{2}}\text{ }{{\theta}_{3}}\text{
}{{\theta}_{4}}\right]}^{T}}$. $Fkine()\in{{\mathbb{R}}^{3}}$ means the
forward kinematics function, which can output the fingertip position according
to different rotation angle of joints.

##  III RESULTS

In this section, we will present experiment results of the forward and inverse
kinematics of one finger with three ellipse joints. The simulations were done
in Matlab 2019a on a laptop with 32GB RAM and Intel i7-8500H CPU. All code is
available on https://github.com/zhifangsixia/Ellipse_Joint_Finger_Simulation

###  III-A The Forward Kinematics Simulation

In forward kinematics, we have input the rotation angle of joints. For getting
the work space of fingertip, we have used around 8000 sets of
$\hat{{{\theta}}}$, the value of input is random
(${{\theta}_{1}}\in(-{{50}^{\circ}},{{50}^{\circ}}\text{)
}{{\theta}_{2}},{{\theta}_{3}},{{\theta}_{4}}\in({{0}^{\circ}},{{90}^{\circ}})$,
the limits are according to [17]). Using these data we can get the workspace
of fingertip (Fig.6).

Figure 6: The workspace of finger tip: (a) The simulated and the PCL of
finger. (b) (c) (d) The relationship of fingertip position and the rotation
angle of joints.

Workspace Input: $\hat{\theta}\in{{\mathbb{R}}^{4\times 8000}}\text{ }$
$(\text{
}\hat{\theta}=[{{\hat{\theta}}_{1}},{{\hat{\theta}}_{2}}...{{\hat{\theta}}_{i}}]\text{
}i\in[1,8000],{{\hat{\theta}}_{i}}={{[\theta_{1},\theta_{2},\theta_{3},\theta_{4}]}^{T}})$
$\text{ }para\in{{\hat{R}}^{1\times 12}}$

Output: $fingertip=\left[\begin{matrix}\overline{{{f}_{i}}}\\\
{{{\hat{\theta}}}_{i}}\\\ \end{matrix}\right]\text{
}fingertip\in{{\mathbb{R}}^{7\times 8000}}$

Function  _$\varphi\\_T\\_Calculate$(_$\theta$ , $a$, $b$_)_:

According equation (3), (4) and (5) to calculate

the transformation matrix.

return $T$

Function  _$Fkine$(_$\theta$ , $para$_)_:

According to $\varphi\\_T\\_Calculate()$ and equation (9)

to get the collection of fingertip positions.

return $\overline{f}$

for  _$i=1:8000$_ do

$fingertip(1:3,i)={{\hat{\theta}}_{i}}$

$fingertip(4:7,i)=Fkine({{\hat{\theta}}_{i}})$

save $fingertip$

Algorithm 1 Workspace of Fingertip

Algorithm $1$ shows the process of getting the workspace of fingertip. We have
defined two functions, the $\varphi\\_T\\_Calculate()$ is to calculate
rotation angle of coordinate frame by rotation angle of joint $\theta$, major
axis $a$ and minor axis $b$ of elliptical trajectory. This calculation is
according to equation (3), (4) and (5). $T$ is the transformation matrix
between coordinate frames. Another function is $Fkine()$, it’s for calculating
the coordinate of fingertip $\overline{f}$ in coordinate frame $0$. This
calculation is according to equation (9). $fingertip$ is a matrix with seven
rows and eight thousands columns. The first three rows are the coordinate of
fingertip, the rest four rows are the rotation angle of joints.

###  III-B The Inverse Kinematics Simulation

We have applied two methods to solve the inverse kinematics problem of this
finger. The first one is to get the analytical solution of all points on
desired trajectory. The second one is searching in the pre-generated PCL. It
searches the $\hat{\theta}$ which makes fingertip has minimum Euclidean
distance between targets $\overline{f}$ and points in PCL. We have defined two
trajectories: the ’heart’ trajectory and circle trajectory. We have taken 628
points from ’heart’ trajectory and 360 points from circle trajectory, and
using two methods to get the inverse kinematics solutions of trajectories.

In Algorithm 2, we have showed the process of getting analytical solution of
inverse kinematics. The input $T$ means the data set of points that we used to
form the trajectory. The output $\hat{\theta}$ means the collection of the
movement of joint angles. In first function, we have defined a solver for
nonlinear equation. The equation should be the form as $F(\theta)=0$. The
second function is based on equation (10) and the
$NonlinearFunction\\_Solver$, we can calculate the movement angle of Ad/Ab by
using the $para$ (parameters of finger) and $T$. The third function is to
calculate the F/E movement of MCP, PIP and DIP joints based on equation (15).
We input the $T$, $para$ and $\hat{\theta_{1}}$ to get the $\hat{\theta}$.

IK_Analytical Input: $T\in{{\mathbb{R}}^{3\times n}}\text{
}T=[{{\overline{f}}_{1,}}{{\overline{f}}_{2}}...{{\overline{f}}_{n}}]\text{
};{{\overline{f}}_{n}}={{[f_{x}^{n},f_{y}^{n},f_{z}^{n}]}^{T}}$

$\text{ }\text{ }\text{ }\text{ }\text{ }\text{ }\text{ }\text{ }\text{
}\text{ }para\in{{\mathbb{R}}^{1\times 12}}$

Output: $\hat{\theta}\in{{\mathbb{R}}^{4\times n}}\text{
}\hat{\theta}={{\left[{{{\hat{\theta}}}_{1}},{{{\hat{\theta}}}_{2}},{{{\hat{\theta}}}_{3}},{{{\hat{\theta}}}_{4}}\right]}^{T}}$
$\text{
}{{\hat{\theta}}_{i}}=[\theta_{i}^{1},...\theta_{i}^{n}]\in{{\mathbb{R}}^{1\times
n}}\text{ }i\in\left[1...4\right]$

Function  _$NonlinearFunction\\_Solver$(_$F(\theta)$ , $\theta$_)_:

Set the error we can bear: $error$; Using equation

(13) to do iterations until $error$ is enough small.

return $\theta$

Function  _$theta1\\_calculate$(_$T$ , $para$_)_:

Using $NonlinearFunction\\_Solver$ and equation

(10) to calculate the $\hat{\theta}_{1}$

return $\hat{\theta}_{1}$

Function  _$theta234\\_calculate$(_$T$ , $para$,$\hat{\theta}_{1}$_)_:

Substitute $\hat{\theta}_{1}$ into $\hat{\theta}$

Get $\overline{f}_{n}$ from $T$

Using $NonlinearFunction\\_Solver$ to get the $\hat{\theta}$ return
$\hat{\theta}$

$\hat{\theta_{1}}=theta1\\_calculate(T,para)$

$\hat{\theta}=theta234\\_calculate(T,para,\hat{\theta_{1}})$

save $\hat{\theta}$

Algorithm 2 The IK of Analytical Solution

IK_PCL Input: $T\in{{\mathbb{R}}^{3\times n}}\text{
}PCL\in{{\mathbb{R}}^{7\times m}}\text{
}PCL=\left[{{P}_{1}},{{P}_{2}}...{{P}_{m}}\right]$
${{P}_{m}}=\left[\begin{matrix}{{\overline{f}}_{m}}\\\
{{{\hat{\theta}}}_{m}}\\\ \end{matrix}\right]\text{
}{{\overline{f}}_{m}}=({{x}_{m}},{{y}_{m}},{{z}_{m}})\text{
}{{\hat{\theta}}_{m}}=(\theta_{1}^{m}...\theta_{4}^{m})$

Output: $\hat{\theta}\in{{\mathbb{R}}^{4\times n}}$

for  _$i=1:n$_ do

$T\\_error=PCL(1:3,:)-T(:,i)\text{ (}T\\_error\in{{\mathbb{R}}^{3\times m}})$

$error\\_norm=sum(T\\_error.*T\\_error,1)$

$(error\\_norm\in{{\mathbb{R}}^{1\times m}})$

$\left[\ast{\ },Index(i)\right]=\min(error\\_norm)\text{
}Index\in{{\mathbb{R}}^{1\times n}}$
${{\hat{\theta}}_{i}}=PCL\left[4:7,Index(i)\right]$

save $\hat{\theta}$

Algorithm 3 The IK of PCL Solution Figure 7: Inverse kinematics of analytical
solution and PCL solution for ’Heart’ trajectory

In Algorithm 3, we have expressed that how to find the minimum Euclidean
distance between the points in PCL and the target points in $T$. The input are
$T$ and PCL. The PCL has $7$ rows and $m$ columns. The first three rows are
the coordinate of fingertip, the rest four rows are the corresponding movement
angle of joints. $m$ is the amount of points in PCL. For each loop, we
calculate the Euclidean distance between the point on trajectory with all
points in PCL. $sum(A,1)$ is the function that will sum each column of $A$.
$[\ast{\ },Index]=min(A)$ will return the index (position) of minimum value in
$A$.

Then we can get the simulation result as Fig.7 and Fig.8, the blue curve of
’Heart’ and ’Circle’ are the desired trajectory that we have generated.
According to them, we have done the inverse kinematics in two methods and
input the results into $Fkine()$ function to generate the simulation trace.
The red curve is the results of simulation. In fingertip error charts, the
horizontal axis: Points, means the points to form this trace. We have showed
the fingertip error between desired and simulation trace. We can find that for
analytical solution, it has high accuracy, which the error can be neglected.
However, for PCL solution, it has $\pm 2\text{ mm}$ error for ’heart’ trace
and maximum $6\text{ mm}$ error for circle. At bottom of Fig.7 and Fig.8, we
have showed the relationship between the displacement in one axis and the
rotation angle of the joints. The heart and circle trace are on the plane
which parallels to $XZ\text{ plane}$ and $XY\text{ plane}$ respectively.

Figure 8: Inverse kinematics of analytical solution and PCL solution for
circle trajectory

###  III-C Comparison of Fitting Experiment

Figure 9: The experiment setup

In this experiment, we stick reflective markers on index finger, then use
OptiTrack cameras to track the markers. We use five OptiTrack cameras to
detect 4 markers position (Fig.9). One of the markers is fixed on fingertip,
and the other three markers are fixed on Intermediate Phalange, Proximal
Phalange and Metacarpal of index finger (for locating the finger
configuration). The collection of fingertip position is shown in Fig.10 (a).
The blue points are the instantaneous position of the fingertip. It is the
movement of the fingertip relative to the Intermediate Phalange, we name it as
the movement of DIP joint. We project these points onto a plane to facilitate
fitting. In Fig.10 (b), we have used elliptical and circular trajectories to
fit these points. For a single joint, the end (fingertip) trajectory of an
elliptical joint is an ellipse, the end trajectory of a circular joint is a
circle. Both of them have the same Center of Rotation (CoR), The position of
CoR is calculated based on the position of the skin creases of the human index
DIP joint.

Figure 10: The comparison of result fitting by elliptical and circular
trajectory: (a) the movement of DIP joint. (b) The error of fitting. (c) The
fitting trajectory.

Fig.10 (b) shows the result of fitting. We can see that the fitting
performance of the elliptical trajectory is better. Fig.10 (c) shows the error
of fitting. The red curve is the error of circular trajectory and the black
curve is the error of elliptical trajectory. We firstly use the actual
position of fingertip to get the rotation angle of the joint and calculate the
position of the corresponding point on the trajectory, then calculate the
distance between the corresponding point and the actual position, and divide
by the distance between the corresponding point and the CoR to get the error.
The Mean Squared Error (MSE) of elliptical trajectory is 0.0196. The MSE of
circular trajectory is 0.0244. The MSE of elliptical trajectory is about
$81.7\%$ of circular trajectory, which means the fitting result of the
elliptical trajectory has smaller and more stable error.

###  III-D Leap Motion with Index Fingers

Figure 11: The experiment setup: (a) the rest posture. (b) and (c) the process
of flexion. (d) the process of adduction.

After the first experiment, we have conducted the second experiment. In second
experiment, we have collected the data of the entire index finger. Because the
OptiTrack system can only recognize large markers (compared with finger size)
and it is easy to confuse different markers when they are close, we choose
Leap Motion to collect the data of finger.

The Leap Motion is able to detect the hand moving and it can return the start
and end points coordinate of bones of the fingers. From the data collected we
can get the length of bones. We also need to get the length of the major and
minor axis of the elliptical trajectory. We let the fingers to do random
movement to collect the data, and use nonlinear optimization method to fit the
length of the major and minor axis as we have done in previous experiment.

Figure 12: The simulation of finger based on real fingertip position

In this experiment, we try to make our index finger drawing square and circle.
Fig.12 (1a) and (2a) show the projections of trajectories on plane that we
collected from real human index fingertip. We have collected about 300 points
for square trajectory and 160 points for circle trajectory. We have done the
inverse kinematics according to these two trajectories by using the kinematics
model we have fitted in the previous experiment. Then using the result of
inverse kinematics to do the forward kinematics to get the simulated fingertip
position. The error in Fig.12 (1b) and (2b) is from the norm of the difference
between the real fingertip position and the simulated fingertip position. We
can find that for both trajectories, the simulation results have very high
accuracy.

##  IV CONCLUSIONS AND FUTURE WORK

In this article, we have proposed a new kinematics model of the multi-joint
finger based on elliptical joints, which can better describe the movement of
the human finger compared to the conventional circular joints. The error (MSE)
of the elliptical joint is $18.3\%$ smaller than that of a circular joint when
fitting the human finger joint movement. In our kinematics model, there are
four bones with three joints. The Metacarpal is fixed, and the rest bones can
move, which are based on joints. Among the joints, DIP and PIP each have one
degree of freedom, which is F/E. MCP has two degrees of freedom; one is F/E,
the other is Ad/Ab. All joints of the finger are elliptical joints. In our
work, we have also used elliptical trajectory to fit the motion of the joint
and analyzed the geometric relationship between the rotation angle of the
joint and the rotation angle of the coordinate frame in a single joint. Then
we have analyzed the difference between our novel elliptical joint and the
conventional circular joint. The most obvious difference is the relationship
between the rotation angle of the joint and the rotation angle of the
coordinate frame. Finally, we have proposed a method for developing a single
joint kinematics model.

After that, we have established the forward kinematics model of a multi-joint
finger, and then we have developed two methods to get the inverse kinematics
solution. The first method is to derive the analytical solution formula of
inverse kinematics and then use the numerical method to find the solution of
the formula. The second method is to use forward kinematics to generate PCL.
When solving inverse kinematics, the algorithm will search for points that are
the closest to the target point in PCL and then return the joint angles.

In the experimental part, we have done four experiments to verify our proposed
finger kinematics model. The first one is a forward kinematics simulation
experiment. In this experiment, we use forward kinematics to input random
joint angles to generate PCL. The second experiment is a simulation experiment
of inverse kinematics. We first define two target trajectories: a ’heart’
trajectory and the other is a circular trajectory, and then we use two methods
we proposed to solve the inverse kinematics of the trajectory. We found that
the first method has higher accuracy than the PCL searching method.

The third experiment is to fit the movement of the DIP joint by using an
elliptical and circular joint. The experiment results reveal that the
elliptical joint shows high fitting performance, and it has a more minor and
more stable error. The last experiment uses a multi-joint finger kinematics
model based on elliptical joints to fit the movement of the human index
fingertip. We first captured the configuration of random finger movement with
Leap Motion and used the nonlinear optimization method to change the major and
minor axes of the elliptical joints of our model to fit these data. Then we
use a human finger to draw square and circular trajectories and record them by
Leap Motion. Then use our fitted kinematics model to solve the inverse
kinematics of the trajectories. We found that our simulated fingertip
trajectory is similar to the fingertip trajectory collected from the subjects,
and the error is minimal, indicating that our kinematics model has high
accuracy. In the future, we plan to make a prototype of a robotic finger with
elliptical joints to validate our kinematics model.

## References

  * [1] M. Grebenstein, M. Chalon, W. Friedl, S. Haddadin, T. Wimböck, G. Hirzinger, and R. Siegwart, “The hand of the dlr hand arm system: Designed for interaction,” The International Journal of Robotics Research, vol. 31, pp. 1531–1555, 11 2012.
  * [2] P. Bizony, “Robonaut flies in [humanoid space robot],” Engineering & Technology, vol. 5, pp. 30–33, 11 2010.
  * [3] A. D. Deshpande, Z. Xu, M. J. V. Weghe, B. H. Brown, J. Ko, L. Y. Chang, D. D. Wilkinson, S. M. Bidic, and Y. Matsuoka, “Mechanisms of the anatomically correct testbed hand,” IEEE/ASME Transactions on Mechatronics, vol. 18, pp. 238–250, 02 2013.
  * [4] Zhe Xu and E. Todorov, “Design of a highly biomimetic anthropomorphic robotic hand towards artificial limb regeneration,” in 2016 IEEE International Conference on Robotics and Automation (ICRA), pp. 3485–3492, May 2016.
  * [5] Y.-J. Kim, J. Yoon, and Y.-W. Sim, “Fluid lubricated dexterous finger mechanism for human-like impact absorbing capability,” IEEE Robotics and Automation Letters, vol. 4, pp. 3971–3978, 10 2019.
  * [6] W. Khalil, E. Dombre, and M. Nagurka, “Modeling, identification and control of robots,” Applied Mechanics Reviews, vol. 56, pp. B37–B38, 05 2003.
  * [7] K. KUCZYNSKI, “The thumb and the saddle,” Hand, vol. 7, pp. 120–122, 06 1975.
  * [8] Y. Youm, T. Gillespie, A. Flatt, and B. Sprague, “Kinematic investigation of normal mcp joint,” Journal of Biomechanics, vol. 11, pp. 109–118, 01 1978\.
  * [9] M. Truman, L. Focht, and A. Hollister, “Prosthetic joint replacement design methods to restore kinematics and stability while preventing material overload,” in Proceedings of the 1995 Fourteenth Southern Biomedical Engineering Conference, pp. 201–205, 1995.
  * [10] J. Friedman and T. Flash, “Trajectory of the index finger during grasping,” Experimental Brain Research, vol. 196, pp. 497–509, 06 2009.
  * [11] D. G. Kamper, E. G. Cruz, and M. P. Siegel, “Stereotypical fingertip trajectories during grasp,” Journal of Neurophysiology, vol. 90, pp. 3702–3710, 12 2003.
  * [12] A. Seth, M. Sherman, P. Eastman, and S. Delp, “Minimal formulation of joint motion for biomechanisms,” Nonlinear Dynamics, vol. 62, pp. 291–303, 05 2010.
  * [13] K.-M. Lee and J. Guo, “Kinematic and dynamic analysis of an anatomically based knee joint,” Journal of Biomechanics, vol. 43, pp. 1231–1236, 05 2010\.
  * [14] B. Buchholz and T. J. Armstrong, “A kinematic model of the human hand to evaluate its prehensile capabilities,” Journal of Biomechanics, vol. 25, pp. 149–162, 02 1992.
  * [15] M. Ruiz Villarreal, “Category:bones of the human hand - wikimedia commons,” 11 2013.
  * [16] J. Lai, K. Huang, B. Lu, Q. Zhao, and H. Chu, “Verticalized-tip trajectory tracking of a 3d-printable soft continuum robot: Enabling surgical blood suction automation,” IEEE/ASME Transactions on Mechatronics, pp. 1–1, 2021\.
  * [17] C. C. Norkin and D. J. White, Measurement of joint motion : a guide to goniometry.  F.A. Davis Company, 5th ed., 2016.