Needle bevel geometry influences the flexural deflection magnitude in ultrasound-enhanced fine-needle biopsy

It has been recently demonstrated that use of ultrasound increases the tissue yield in ultrasound-enhanced fine-needle aspiration biopsy (USeFNAB) as compared to conventional fine-needle aspiration biopsy (FNAB). To date, the association between bevel geometry and needle tip action has not been widely explored. In this study, we studied the needle resonance characteristics and deflection magnitude of various needle bevel geometries with varying bevel lengths. With a conventional lancet, having a 3.9 mm long bevel, the tip deflection-to-power ratio (DPR) in air and water was 220 and 105 µm/W, respectively. This was higher in comparison to an axi-symmetric tip, having a bevel length of 4 mm, which achieved a DPR of 180 and 80 µm/W in air and water, respectively. This study emphasised the importance of relationship between flexural stiffness of bevel geometry in the context of various insertion media and, thus, could provide understanding on approaches to control post-puncture cutting action by modifying the needle bevel geometry, essential for the USeFNAB application.

As illustrated in Fig. 2b,c, for an infinite (boundless) beam with a cross-sectional area A, and assuming large wavelength with respect to the cross-sectional dimension of the beam, the flexural (or bending) phase velocity c EI was defined 22 : where E was the Young's modulus ( N/m 2 ), ω 0 = 2πf 0 was the excitation angular frequency (rad/s), where f 0 was the linear frequency (1/s or Hz), I was the area moment of inertia (m 4 ) around the axis of interest, and m ′ = ρ 0 A was the mass per unit length (kg/m), where ρ 0 was the density (kg/m 3 ) , and A was the cross-sectional (xy-plane) area of the beam ( m 2 ). Since the force applied in our case, was parallel to the vertical y-axis i.e. F y � j , we were only concerned with the area moment of inertia around the horizontal x-axis i.e. I xx , hence: (1)  36 , where α was the primary bevel angle, θ was the secondary bevel rotation angle, and φ was the secondary bevel angle, when rotated, measured in degrees ( • ). (b) Linear a-symmetric single step bevel (referred to as "standard" in DIN 13097:2019 37 ), and (c) linear axi-symmetric (circumferential) single-step bevel. finite element model (FEM) employed a harmonic point force F y � j to excite the needle tube at the proximal end, a point deflection and velocity ( ũ y � j , ṽ y � j ) was measured at the tip to allow a calculation of transfer mechanical mobility. y was defined as the flexural wavelength associated with the vertical force F y � j . (c) Definitions of the centre of gravity, the cross-sectional area A, and the moments of inertia I xx and I yy , around the x and y axes, respectively. www.nature.com/scientificreports/ where y CG is the y-coordinate of the centre of gravity of the needle tube in the xy-plane.
For the finite element model (FEM), a purely harmonic displacement (m) was assumed, therefore acceleration ( m/s 2 ) was expressed as ∂ 2 � u/∂t 2 = −ω 2 � u , such as that � u(x, y, z, t) := u x � i + u y � j + u z � k was a three-dimensional displacement vector defined in the spatial coordinates. Substituting the latter, the law of balance of momentum in its Lagrangian form for finite deformation 23 , was given according to its implementation in COMSOL Multiphysics software (version 5.4-5.5, COMSOL Inc., Massachusetts, USA), as: k was the tensor divergence operator, and σ was the second Piola-Kirchhoff stress tensor (of second order, N/m 2 ), and � k was the volumetric force vector per deformed volume ( N/m 3 ), and e jφ was the phase of the volumetric force having a phase angle φ (rad). In our case, the volumetric body force was zero, and our model assumed geometrical linearity, and small purely elastic strain i.e. ε el = ε , where ε el and ε were the elastic and total strains, respectively (of second order, dimensionless). The Constitutive Hookean isotropic elasticity tensor C was defined using the Young's modulus E ( N/m 2 ) and the Poisson's ratio v, so that C := C(E, v) (of fourth order). Therefore the calculation of stress becomes σ := C : ε.
Computation was done with 10-node tetrahedral elements with element size of ≤ 8 µm. The needle was simulated in a vacuum, and a magnitude of transfer mechanical mobility (m s −1 N −1 ) was defined as |Ỹ v y F y | = |ṽ y � j|/|F y � j| 24 , where ṽ y � j was the the output complex velocity at the tip, and F y � j was the complex driving force located at the proximal end of the tube, as illustrated in Fig. 2b. Transfer mechanical mobility was expressed in decibels (dB), using the maximum as a reference, i.e. 20 log 10 (|Ỹ |/|Ỹ max |) . All FEM studies were conducted at 29.75 kHz.
Fabrication of needle constructs. The needle constructs (Fig. 3) consisted of a conventional 21 gauge hypodermic needle (catalogue number: 4665643, Sterican , outer diameter 0.8 mm, length 120 mm, stainless chromium nickel steel AISI type 304 grade, B. Braun Melsungen AG, Melsungen, Germany) fitted with a Luer Lock plastic hub made of polypropylene at the proximal end, and modified accordingly at the tip. Needle tubes were soldered to waveguides, as shown in Fig. 3b. The waveguides were 3D printed with stainless steel (EOS Stainless Steel 316L in EOS M 290 3D Printer, 3D Formtech Oy, Jyväskylä, Finland), then fastened via a M4 bolt to a Langevin transducer. The Langevin transducer consisted of 8 piezo ring elements, loaded by two masses at either end.
Three axi-symmetrically bevelled tips were fabricated ( Fig. 3) (TAs Machine Tools Oy) with bevel lengths (BL, as defined in Fig. 2a) of 4.0, 1.2 and 0.5 mm, corresponding to bevel angles (BA) of ≈ 2 • , 7 • , and 18 • , respectively. The masses of waveguides and needles were 3.4 ± 0.017 g (mean ± s.d., n = 4) for bevels L and AX1-3, respectively    Fig. 3b, respectively. For all needle constructs, the length from needle tip to the tip of the waveguide (i.e. soldering region), was 4.3 cm, and the needle tube was orientated so that bevel planes faced upwards (i.e. parallel to the y-axis), as in (Fig. 2).

Modal analysis.
A custom-script in MATLAB (R2019a, The MathWorks Inc., Massachusetts, USA), running on a computer (Latitude 7490, Dell Inc., Texas, USA), was used to generate a linear sine-sweep from 25 to 35 kHz for duration of 7 s, which was converted to an analogue signal via a digital-to-analogue (DA) converter (Analog Discovery 2, Digilent Inc., Washington, USA). The analogue signal V 0 (0.5 V pk-pk ) was then amplified using a custom-made radio frequency (RF) amplifier (Mariachi Oy, Turku, Finland). The incident amplified voltage V I was output from the RF amplifier at an output impedance of 50 , to the transformer built into the needle construct, which had an input impedance of 50 . The Langevin transducer (back and front mass-loaded sandwich piezoelectric transducer) was used to generate the mechanical wave. The custom-made RF amplifier was equipped with a dual-channel standing-wave power ratio (SWR) meter, which allowed both the incident V I and reflected amplified voltages V R to be recorded via the analogue-digital (AD) converters (Analog Discovery 2) at sampling frequency of 300 kHz. The excitation signal was amplitude modulated at the beginning and end to prevent signal transients overloading the amplifier's input.
Using a custom script implemented in MATLAB, frequency response functions (FRFs) i.e. H (f ) , were estimated offline using a swept-sine dual channel measurement technique 25 (Fig. 4), which assumed a linear timeinvariant system. In addition, a bandpass filter having passband between 20 and 40 kHz was applied to remove any unwanted frequencies from signal. In reference to transmission line theory, H (f ) in this case was equivalent to the voltage reflection coefficient i.e.
Since the amplifier output impedance Z 0 was matched to input impedance of the transformer built-in with the transducer, the electrical power reflection coefficient In the case when absolute values of electrical power were needed, the incident P I and reflected P R powers (W) were calculated by taking the root-mean-square (r.m.s.) of the corresponding voltages, such as that for a sinusoidally-excited transmission line, P = V 2 /(2Z 0 ) 26 , where Z 0 was 50 . The electrical power transmitted to the load P T (i.e. to the insertion medium) could be calculated as |P I − P R | (W, r.m.s.), and the power transfer efficiency (PTE) could be defined and given as a percentage (%), so that 27 : The FRFs were then used to estimate the modal frequencies f 1−3 (kHz) of the needle construct, and their corresponding power transfer efficiencies, PTE 1−3 . The full-width at half-maxima ( FWHM 1−3 , Hz) were estimated directly from PTE 1−3 , obtained from the one-sided linear frequency spectra at modal frequencies f 1−3 described in Table 1.  For each needle bevel type, we recorded 300 high speed camera frames, measuring 128 × 128 pixels with a spatial resolution of 1/180 mm ( ≈ 5 µm) per pixel, and a time resolution of 310 000 frames per second. As outlined in Fig. 6, each frame (1) was cropped (2) so that needle tip was located in the last row (bottom) of frame, then the histogram of the image was computed (3), so that Canny thresholds 1 and 2 could be determined. Then Canny edge detection 28 with a 3 × 3 Sobel operator was applied (4), and the location was computed for a cavitationfree bevel-edge pixel (marked × ) for all 300 time steps. To determine the peak-to-peak deflection at the tip, the derivative (using a central difference algorithm) was calculated (6), and the frames containing the local extrema (i.e. peaks) of deflection were identified (7). Following a visual inspection for cavitation-free edges, a frame pair (or two frames that are half of the time period apart) was chosen (7), and the deflection at the tip was measured (marked × ). The above was implemented in Python (v3.8, Python Software Foundation, python.org), utilising OpenCV's Canny edge detection algorithm (v4.5.1, Open Source Computer Vision Library, opencv.org). Finally, the deflection-to-power ratio (DPR, µm/W), was calculated as the ratio of the peak-to-peak deflection over the transmitted electrical power P T (W, r.m.s.). Statistical analysis. Since the sample size was small (n = 5), and normality could not be assumed, a twosample, two-sided, Wilcoxon rank sum test was used (R, v4.0.3, R Foundation for Statistical Computing, r-project.org), to compare the tip-deflection magnitudes of the different needle bevels. 3 comparisons were done for each bevel, so a Bonferroni-correction was applied, and the adjusted significance level was 0.017, at 5% error rate.  kHz, the flexural half-wavelength ( y /2 ) for the 21 gauge needle tubing was ≈ 8 mm. The flexural wavelength decreased along bevel when approaching the tip. At the tip, y /2 was ≈ 3, 1, and 7 mm for the conventional lancet (a), a-symmetric (b), and axi-symmetric (c) single step bevels, respectively. Consequently this meant the range of the variation was ≈ 5 mm for the lancet (owing to the two lancet planes generating a single sharp point 29,30 ), 7 mm for the a-symmetric bevel, and 1 mm for the axi-symmetric bevel (where the centre of gravity stayed constant, so effectively only tube wall thickness varied along bevel).
Peaks of the mobility |Ỹ v y F y | indicated optimal tube length (TL) and bevel length (BL) combinations (Figs. 8,9). For the conventional lancet, since its dimensions were fixed, the optimal TL was ≈ 29.1 mm (Fig. 8). For the a-symmetric and axi-symmetric bevels (Fig. 9a,b, respectively), FEM studies included BLs 1 to 7 mm, so the optimal TLs varied from 26.9 to 28.7 mm (range 1.8 mm) and 27.9 to 29.2 mm (range 1.3 mm), respectively. For the a-symmetric bevel (Fig. 9a), optimal TLs increased linearly reaching a plateau at a BL of 4 mm, then steeply declined from BLs 5 to 7 mm. For the axi-symmetric bevel (Fig. 9b), optimal TLs increased linearly with longer BLs, and eventually plateaued at BL of ≈ 6 to 7 mm. An extended study of axi-symmetric bevel (Fig. 9c), showed another set of optimal TLs at ≈ 35.1-37.1 mm. The two sets of optimal TLs were separated by a distance of ≈ 8 mm (equivalent to y /2 ), for all BLs.
Modal behaviour. The needle construct exhibited three natural frequencies f 1−3 , which were categorised into low, middle and high modal regions, as summarised in Table 1. The magnitudes of PTE were recorded as in Fig. 10, and then analysed in Fig. 11. The following gives overview of the findings for each modal region: 1st modal region: f 1 did not vary greatly with type of insertion medium, but varied with changing bevel geometry. f 1 decreased with decreasing bevel length (27.1, 26.2, and 25.9 kHz for AX1-3, in air, respectively). The region-averages of PTE 1 and FWHM 1 were ≈ 81% and 230 Hz, respectively. FWHM 1 was highest in gelatin for the Lancet (L, 473 Hz). Note it was not possible to estimate FWHM 1 for AX2 in gelatin, due to low recorded magnitudes of FRF. 2nd modal region: f 2 varied with type of insertion medium and bevel. In air, water, and gelatin, averages of f 2 were 29.1, 27.9, and 28.5 kHz, respectively. This modal region also exhibited PTE as high as 99%, which was the highest among all measurement groups, with a region-average of 84%. The region-average of FWHM 2 was ≈ 910 Hz. 3rd modal region: f 3 frequencies varied with type of insertion medium and bevel. In air, water, and gelatin, the average values of f 3 were 32.0, 31.0, and 31.3 kHz, respectively. The region-average of PTE 3 was ≈ 74%, which was the lowest among all regions. The region-average of FWHM 3 was ≈ 1085 Hz, which was higher than 1st and 2nd regions.
(6) (7) x y z Figure 6. Needle tip deflection was measured using a sequence of frames captured from high-speed camera at 310 kHz, using a 7-step algorithm (1-7), involving cropping (1-2), Canny edge detection (3-4), computation of edge pixel location (5) and its time derivative (6), and finally measuring the peak-to-peak deflection at tip from a visually inspected frame pair (7). www.nature.com/scientificreports/ Measured deflection. The following refers to Fig. 12 and Table 2. The lancet (L) deflected the most (with high significance to all tips, p < 0.017) in both air and water (Fig. 12a), achieving the highest DPR (up to 220 µm/W in air). In air, AX1 which had higher BL, deflected higher than AX2-3 (with significance, p < 0.017), while AX3 (which had lowest BL) deflected more than AX2 with a DPR of 190 µm/W. In water at 20 mm, no significant differences ( p > 0.017) were found in deflection and PTE for AX1-3. PTE levels were overall higher  Fig. 1a,b,c). Mean y /2 was 5.65, 5.17, and 7.52 mm for the lancet, a-symmetric, and axi-symmetric bevels, respectively. Note the tip thickness of a-symmetric and axi-symmetric bevels was limited to ≈ 50 µm.  www.nature.com/scientificreports/ in water (90.2-98.4%) than air (56-77.5%) (Fig. 12c), noting cavitation events were clearly present in water during experimentation (Fig. 13, also see Supplementary Information).

Discussion
To summarise, we modelled the change of the flexural wavelength (Fig. 7), and computed the transfer mechanical mobility for a combination of tube and bevel lengths (Figs. 8,9), for a conventional lancet, a-symmetric and axisymmetric bevel geometries. Based on the latter, we estimated an optimal distance of 43 mm (or ≈ 2.75 y at 29.75 kHz) from tip to soldering region, as illustrated in Fig. 5, and fabricated accordingly three axi-symmetric bevels with varying bevel lengths. We then characterised their frequency behaviour in comparison to the conventional lancet, in air, water, and ballistic gelatin 10% (w/v) (Figs. 10,11), and identified the mode most appropriate for comparing deflection of bevels. Finally, we measured the flexural-wave deflection at the needle tip in air and at 20 mm depth in water, and quantified the electrical power transfer efficiency to insertion medium (PTE, %) and the deflection-to-power ratio (DPR, µm/W) for each bevel type (Fig. 12).
The results show that the needle-bevel geometry affects the deflection amplitude at the needle-tip. The lancet achieved the highest deflection, as well as the highest DPR, in comparison to axi-symmetric bevels, which on average deflected less (Fig. 12). The axi-symmetric 4 mm bevel (AX1) having the longest bevel length, achieved statistically significant highest deflection in air ( p < 0.017 , Table 2), in comparison to other axi-symmetric needles (AX2-3), but no significant differences were observed, when the needle was placed in water. Therefore, in terms of the peak deflection at the tip, there was no clear benefit of having longer bevel lengths. Considering this, the results suggest that the bevel geometries investigated in this study have a greater effect on deflection amplitudes than the bevel length. This may be associated with the flexural stiffness, depending e.g. on the overall thickness of the flexurally-bending material and needle structure.
In the experimental studies, the magnitude of the reflected flexural waves were affected by the boundary conditions at the needle-tip. When the needle-tip was inserted in water and gelatin, the average of PTE 2 was ≈ 95%, compared to an average of 73% and 77% for PTE 1 and PTE 3 , respectively (Fig. 11). This suggested that the greatest transmission of acoustic energy into the embedding medium, i.e. water or gelatin, occurred at f 2 . Similar behaviour was observed in a previous study 31  www.nature.com/scientificreports/ depth 32 and mechanical properties of tissue provide a mechanical load on the needle, and hence are expected to affect resonant behavior of USeFNAB. Therefore, resonance tracking algorithms e.g. 17,18,33 , could be utilised to optimise the delivered acoustic power through the needle. The simulation study of flexural wavelengths (Fig. 7) revealed the axi-symmetric had higher structural rigidity at the tip (i.e. higher flexural stiffness), than both the lancet and a-symmetric bevels. Deduced from (1), and utilising the known velocity-frequency relationship, we estimated the flexural stiffnesses at the tip to be ≈ 200, 20, and 1500 MPa for the lancet, a-symmetric, and axi-symmetric bevels, respectively. This corresponded to y of ≈ 5.3, 1.7, and 14.2 mm at 29.75 kHz, respectively (Fig. 7a-c). Considering clinical safety during USeFNAB procedures, the influence of geometry on structural stiffness 34 of bevel needs to be assessed.
The bevel vs. tube length parametric study (Fig. 9), revealed that the range of optimal TLs was higher for the a-symmetric (1.8 mm), than axi-symmetric bevel (1.3 mm). In addition, the mobilities plateaued at ≈ 4 to 4.5 mm and at 6 to 7 mm, for a-symmetric and axi-symmetric bevels, respectively (Fig. 9a,b). The practical relevance of this finding translates to manufacturing tolerances, e.g. a lower range of optimal TLs may mean a higher precision for the lengths is required. Meanwhile, the plateaus in mobility provide greater tolerance for selecting bevel lengths at a given frequency, without significantly affecting mobility.
The study included the following limitations. Measuring needle deflection directly using edge detection and high-speed imaging (Fig. 12), meant that we were limited to optically-transparent media such as air and water. We would also like to note that we did not use experiments to validate the modelled transfer mobilities, or vice versa, rather, FEM studies were used to determine the optimal lengths for fabrication of the needles. In terms of practical limitations, the length from tip to needle hub was ≈ 0.4 cm longer for the lancet, than the other needles (AX1-3), see Fig. 3b. This could have influenced the modal response of needle construct. In addition, the shape  www.nature.com/scientificreports/ and volume of the soldering at the waveguide-needle termination (see Fig. 3), may have affected the mechanical impedance of the needle construct, introducing uncertainty in mechanical impedance and bending behavior.
To conclude, we have demonstrated experimentally bevel-geometry affects deflection amplitudes in USeFNAB. In the case that higher deflection magnitudes would positively influence the effect of the needle on tissue, e.g. efficacy of post-puncture cutting, the conventional lancet may be recommended for use in USeFNAB, since it achieved the highest deflection magnitude, whilst still maintaining adequate structural rigidity at the tip. In addition, greater tip deflections could enhance bioeffects, e.g. cavitation, as suggested by a recent study 35 , which could be helpful in the development of applications for minimally invasive surgical interventions. Considering that it has already been shown that increasing the total acoustic power could increase biopsy yield in USeFNAB 13 , further quantitative study on sample yield and quality is needed, to assess the detailed clinical benefits of the studied needle geometries.   Figure 11. Analysis of the modal responses shown in Fig. 10 (mean ± s.d., n = 5), for bevels L and AX1-3, in air, water and gelatin 10% (depth 20 mm), featuring (top) three modal regions (low, middle and high), and their corresponding modal frequencies f 1−3 (kHz), (middle) power efficiency PTE 1−3 calculated using Eq. (4), and (bottom) the full-width at half-maximum measurements FWHM 1−3 (Hz), respectively. Note that the measurement of bandwidth was omitted, when low PTE was recorded, i.e. in the case of bevel AX2, FWHM 1 . Mode f 2 was considered to be the most appropriate for comparing deflection of bevels, since it exhibited the highest levels of power transfer efficiency ( PTE 2 ), which were as high as 99%.  4), and (d) the deflection-to-power ratio (DPR, µm/W), which was calculated as the ratio of the peak-to-peak deflection over the transmitted electrical power P T (W, r.m.s.). Table 2. Was there significant difference in tip-deflection between the lancet (L) and the axi-symmetric bevels (AX1-3)? a . a Each pair of bevels (e.g. L vs. AX1) were tested using a two-sample, two-sided, Wilcoxon rank sum test, n = 5, 5% error rate, significance level was 0.017 (Bonferroni-corrected). Bolded entries indicate no statistical significance found between bevel pair (i.e. p ≥ 0.017).

Data availability
The datasets produced during this study are available on reasonable request. (1) (1) (1) (2) Figure 13. Typical high-speed camera shadowgraphs showing the peak-to-peak tip deflection (green and red dashed lines) for the lancet (L) and axi-symmetric tips (AX1-3), in water (depth 20 mm), during a half-cycle, at excitation frequency f 2 (sampling frequency 310 kHz). Captured greyscale images measured 128 × 128 pixels and the pixel size was ≈ 5 µm. A video can be found in Supplementary Information.