Some Mechanical Constraints to the Biomimicry with Peripheral Nerves

Novel high technology devices built to restore impaired peripheral nerves should be biomimetic in both their structure and in the biomolecular environment created around regenerating axons. Nevertheless, the structural biomimicry with peripheral nerves should follow some basic constraints due to their complex mechanical behaviour. However, it is not currently clear how these constraints could be defined. As a consequence, in this work, an explicit, deterministic, and physical-based framework was proposed to describe some mechanical constraints needed to mimic the peripheral nerve behaviour in extension. More specifically, a novel framework was proposed to investigate whether the similarity of the stress/strain curve was enough to replicate the natural nerve behaviour. An original series of computational optimizing procedures was then introduced to further investigate the role of the tangent modulus and of the rate of change of the tangent modulus with strain in better defining the structural biomimicry with peripheral nerves.


Introduction
Biomimetic materials are able to behave similar to materials originating from living organisms and can effectively replicate selected characteristics of natural materials.They could replace biological materials, restoring natural functions when the biological counterparts are unpaired, absent, or unable to correctly perform [1].Biomimetic materials were often used to build scaffolds promoting the regeneration of the natural tissue: biomimetic 3D-printed scaffolds for spinal cord injury repair were built [2], as well as electrically conductive scaffolds mimicking the hierarchical structure of cardiac myofibers [3].Again, in the literature, biomimetic scaffolds for tendon tissue regeneration were described [4][5][6][7].More generally, some research groups focused on biomimetic approaches for bone tissue engineering [8-10], while others obtained a supramolecular biomimetic skin [11].Some scientists succeeded in building a biomimetic nerve platform using cell line as neuronal population [12,13], while others provided some general design principles for a biomimetic artificial nerve [14].More specifically, as "desired biomimetic devices and scaffolds for neurotrophic biomolecules to be implemented in a future PNI-repair scaffold", they described a three-dimensional scaffold architecture designed to improve the regeneration of sensory and motor fibers.The internal channels were described as able to guide a very small group of axons through scaffolds, which were biofunctionalized to bind the receptors expressed at the surface of the growth cone.Interestingly, all the previous structures were enveloped by a material, which formed the outer sheet of the biomimetic artificial nerve, and which was identified as "semipermeable", "semirigid", or "flexible".Nevertheless, the mechanical characteristics of this "material" were not described and are not currently clear, since it should mimic the mechanics of the natural outer sheet of connective tissue enveloping peripheral nerves.In other words, the knowledge of the peripheral nerve mechanical response to axial strain was needed to allow this material to be properly designed.Indeed, peripheral nerves are physical objects [15], which react to elongation by increasing their stiffness to keep the axons integrity and to protect endoneural structures from longitudinal overstretch [16].The study of the mechanical behaviour of peripheral nerves was performed at different scale levels.Thus, the relationship between macroscopic tensile loads and micro-scale deformations [17] was investigated.Again, the nerve was studied as a complex structure showing the characteristics of an isotropic material, providing a deterministic elastic [18,19] or hyperelastic response [20][21][22][23][24] as a function of the strain magnitude.Differently, a stochastic iterative fibril-scale mechanical model was implemented to reproduce the straightening of wavy fibrils and to account for the effects of interfibrillar crosslinks on the overall properties of the tissue [25,26].Recently, an explicit and deterministic framework was proposed to model the behaviour of peripheral nerves under longitudinal stretch [27].This framework was applied here to investigate the main mechanical constraints to the behaviour of a material, which could mimic the axial response of peripheral nerves.The governing equations were proposed in polynomial form together with an original procedure to decrease their overall complexity in physiologically relevant cases.Finally, the minimum number of constraints to achieve an acceptable level of biomimicry with peripheral nerves was presented.

Methods
The mechanical response of the peripheral nervous tissue was assumed to be ruled by an SEF (Strain Energy Function, for unit of volume) affected by the strain measure as well as by the initial direction of the collagen fibrils.In particular, the SEF was written as W = W(C, M ⊗ M), where F was the deformation gradient and C = F T F was the right Cauchy-Green strain tensor, while M was a unit vector field accounting for the direction of the collagen fibrils in the reference configuration.The mean Cauchy-Green stress tensor was, then, written as: and J = det(F).In general, the strain energy function could be written as a function of invariants as W = W(I 1 , I 2 , I 3 , I 4 , I 5 ), where, As usual, the peripheral nervous tissue was considered incompressible, and thus, I 3 = det(C) = 1 [16,[22][23][24]28,29].In addition, since its main structural components are the ground matrix [30,31] and the collagen fibrils [32], the strain energy function was rewritten as W = W(I 1 , I 4 ), where the I 1 invariant accounted for the action of the matrix, and the I 4 invariant modeled the reinforcing action of the collagen fibrils.Therefore, the mean Cauchy-Green stress tensor was: where p was the pressure, B = FF T was the left Cauchy-Green tensor, W i = ∂W ∂I i , i = 1, 4, and m = FM.In particular, the energetic contributions, deriving from the ground matrix and from the collagen fibrils were chosen, as in the recent literature [27].
The nerve underwent to axial strain (z-direction), and thus, its lateral surface was stress-free.Then, the determination of the mean axial stress as a function of the axial strain ( zz ) was performed through Equation (2) solved with the boundary conditions σ xx = σ yy = 0 [28,29].A polynomial form of the mean Cauchy stress was, then, written as: the tangent modulus was written as: while the rate of change of the tangent modulus with strain was: In Equations ( 3)-( 5) the indexes p, q ∈ N should be determined in order to minimize the functional: where σ zz (K 1 , A, D, ∞, ∞, zz ) = lim p→∞,q→∞ σ zz (K 1 , A, D, p, q, zz ).The sensitivity of the mean Cauchy stress in Equation (3) to each parameter was evaluated through the following indexes: S A (1 + ∆A, p, q) = zzmax 0 |σ zz (1, 1 + ∆A, 1, p, q, zz ) − σ zz (1, 1, 1, p, q, zz )|d zz ∆A (8) where ∆K 1 = ∆A = ∆D = 99.In addition, the best combination of the parameters p, q was chosen through the maximization of the metric:

Error Evaluation for Different Polynomial Stress Functions
The logarithm (base 10) of the functional defined in Equation ( 6) was used to quantify the error of Equation (3) with respect to the limit values of p → ∞ and q → ∞ for each combination of finite values for p and q.More specifically, in Figure 1

Sensitivity of Different Polynomial Stress Functions to Parameters K 1 , A, D
The sensitivity of Equation (3) to parameters K 1 , A, D was studied through Equations ( 7)-( 9).More specifically, in Figure 4a, the influence of the change of K 1 from 1 to 100 was studied for all combinations of 0 ≤ p ≤ 10 and 0 ≤ q ≤ 10, when both the other parameters were kept constant (i.e., A = 1, D = 1).The value of sensitivity S K 1 (for sake of simplicity all the dependencies from parameters were dropped out) ranged from 1.314 × 10 −2 to 1.954 × 10 −2 .In the first column, there were the minimum values, while in the second column there were the maximum ones.The value of sensitivity was almost constant (i.e., S K 1 oscillated around 1.192 × 10 −2 ) over all the other combinations.Similarly, in Figure 4b, the influence of the change of the parameter A from 1 to 100 was studied, keeping constant both K 1 = 1 and D = 1.In this case, the sensitivity S A ranged between 0 and 3.635 × 10 −4 .The fist row was the minimum one, the value of sensitivity was S A = 3.602 × 10 −4 in the second row, while its value was almost constant (i.e., S A = 3.635 × 10 −4 ) over all the other combinations.Finally, in Figure 4c, the influence of the change of the parameter D from 1 to 100 was studied, while the other parameters were constant K 1 = 1 and A = 1.The sensitivity S D ranged between 0 and 11.503 × 10 −2 .In this case, the sensitivity increased with the value of q.The fist row was the minimum one, while the value of S D progressively increased from 3.638 × 10 −2 to 11.503 × 10 −2 starting from the 2nd row to the 10th one.

Minimization of the p and q Values
Once studied, the sensitivity of Equation ( 3) to general values of parameters K 1 , A, D, the attention was narrowed to their physiological values.In particular, literature values were used both for peroneal and vagus nerves [27].The metric Ω in Equation ( 10) was used to identify the best combination of parameters (with finite and integer values ranging between 0 and 10) for Equation (3) with physiological parameters for peroneal nerve (K 1 = 3.018 kPa, A = 0.014 kPa, and D = 66.893).The maximization of this metric resulted in Ω = 0.883 for p = 1 and q = 4. Therefore, to further compare the evolution of σ zz (1, 4, ) to the σ zz (∞, ∞, ) (the dependence of parameters was dropped out for sake of simplicity), both their difference and their ratio were studied.More specifically, their difference started from 0 and increased up to 1.57 × 10 −2 kPa for = 0.058, while it decreased up to −1.24 × 10 −1 for = 0.08 (Figure 5b).Similarly, their ratio started form 1 for = 0 and increased up to 1.0000033 for = 0.069, while it decreased up to 0.974 for = 0.08 (Figure 5c).In a similar way, the maximization of the metric in Equation (10) resulted in Ω = 0.612 for p = 7 and q = 7 when the mean physiological value for vagus nerve were used (K 1 = 9.731 kPa, A = 5.277 kPa, D = 10.171).The difference between σ zz (7, 7, ) and σ zz (∞, ∞, ) (the dependence of parameters was dropped out for sake of simplicity) started from 0 for = 0 and increased up to 7.183 × 10 −9 for = 0.074, while it decreased up to −5.955 × 10 −9 for = 0.08 (Figure 5e).Their ratio started from 1 for = 0 and oscillated from 1.0000000006 for = 0.071 and 0.9999999996 for = 0.08 (Figure 5f).(a) The metric Ω identified as best candidate the stress function with p = 1 and q = 4 (in yellow) for the peroneal nerve.(b) Evolution of the difference between stress functions with p = 1, q = 4, and p = ∞, q = ∞.
(c) Evolution of the ratio between stress functions with p = 1, q = 4, and p = ∞, q = ∞.(d) The metric Ω identified as best candidate the stress function with p = 7 and q = 7 (in yellow) for the vagus nerve.(e) Evolution of the difference between stress functions with p = 7, q = 7, and p = ∞, q = ∞.(f) Evolution of the ratio between stress functions with p = 7, q = 7, and p = ∞, q = ∞.In all relevant plots zz = for sake of simplicity

Evolution of Stress, Tangent Modulus, and Rate of Change of the Tangent Modulus with Strain
The analysis performed through Equation ( 10) was able to identify the best combination of finite values of parameters p and q with physiological parameters.In particular, the evolution of Equation (3) with K 1 = 3.018 kPa, A = 0.014 kPa, D = 66.893, p = 1, and q = 4 is shown in Figure 6a for the porcine peroneal nerve.The mean Cauchy stress increased monotonically and non linearly starting from 0 kPa for = 0 up to 4.676 kPa for = 0.08.Similarly, the evolution of tangent modulus was non linearly increasing from 18.018 kPa for = 0 to 229.125 kPa for = 0.08 (Figure 6b), while the rate of change of the tangent modulus with strain started from 12.072 for = 0 and non linearly increased up to 13, 357.327 for = 0.08 (Figure 6c).Furthermore, the evolution of Equation (3) with K 1 = 9.731 kPa, A = 5.277 kPa, D = 10.171,p = 7, and q = 7 is shown in Figure 6d for the canine vagus nerve.The Cauchy stress increased monotonically in a non linear way starting from 0 kPa for = 0 up to 16.300 kPa for = 0.08, while the tangent modulus non linearly increased from 58.386 kPa for = 0 to 555.722 kPa for = 0.08 (Figure 6e).Finally, the rate of change of the tangent modulus with strain started from 0 for = 0 and non linearly increased up to 16, 479.238 for = 0.08 (Figure 6f).

Analysis of the Correlation among Longitudinal Stress, Tangent Modulus, and Rate of Change of the Tangent Modulus with Strain
In order to investigate possible relationships among the previously described quantities, a further analysis is shown in Figure 7.More specifically, in Figure 7a, the evolution of the mean Cauchy stress with strain was compared to the evolution of the tangent modulus with strain for the porcine peroneal nerve.These quantities were positively correlated (R 2 = 0.9997, p-value = 1.88 × 10 −15 , α = 0.05, two-tailed t test) through the non linear expression: Similarly (Figure 7b), the rate of change of the tangent modulus with strain was positively correlated (R 2 = 0.9999, p-value = 4.34 × 10 −10 , α = 0.05, two-tailed t test) to the mean Cauchy stress through the non linear formula: Again (Figure 7c), the evolution of the rate of change of the tangent modulus with strain was positively correlated (R 2 = 0.9999, p-value = 5.27 × 10 −10 , α = 0.05, two-tailed t test) to the tangent modulus through the linear expression: In Equations ( 11)-( 13), the explicit dependence of K 1 = 3.018 kPa, A = 0.014 kPa, and D = 66.893 was neglected for sake of simplicity.In a similar way (Figure 7d), for the canine vagus nerve, the evolution of the tangent modulus was positively correlated with the longitudinal strain (R 2 = 0.9994, p-value = 4.75 × 10 −21 , α = 0.05, two-tailed t test) through the non linear expression: The rate of change of the tangent modulus (Figure 7e) was positively correlated (R 2 = 0.9997, p-value = 1.11 × 10 −19 , α = 0.05, two-tailed t test) to the longitudinal stress through the non linear expression: Finally, the evolution of the rate of change of the tangent modulus with strain was positively correlated (R 2 = 0.9972, p-value = 1.24 × 10 −19 , α = 0.05, two-tailed t test) with the evolution of the tangent modulus (Figure 7f) with the non linear relationship: In Equations ( 14)-( 16), the explicit dependence of K 1 = 9.731 kPa, A = 5.277 kPa, and D = 10.171 was neglected for the sake of simplicity.

Discussion
In this work, some basic mechanical constraints were investigated to design a material biomimetic to the connective tissue of peripheral nerves.This material was anticipated as required to build a biomimetic artificial nerve [14].However, no quantitative information was up to now provided about its mechanical characteristics.As a consequence, since the peripheral nerves are mainly exposed to axial elongations, here, the characteristic for a material to behave like peripheral nerves in extension are quantitatively explored.An original polynomial framework was also provided, in order to simplify as much as possible the writing of the governing equations.To this aim, some different numerical analyses of the behaviour of Equations ( 3)-( 5) were performed.
A first analysis of the magnitude of errors, expressed through the logarithm of the functional Φ(K 1 , A, D, p, q) in Equation ( 6), was able to show that different patterns arose for different values of parameters.In other words, the distribution of errors for all the combinations of integers values of p, q ∈ N and 0 ≤ p ≤ 10, 0 ≤ q ≤ 10 changed as a function of K 1 , A, D parameters.More specifically, for K 1 = 1, D = 1, and q ≥ 5 (Figure 1), the magnitude of errors depended only of p and decreased from about −4 to −13 when p increased from 0 to 10.For K 1 = 100 and D = 100, instead, the magnitude of errors practically only depended of q, ranging from −1 to −4, when q ranged between 0 and 10, while A = 1.It is worthy to notice that the error value ranged between 2 and −2 for A = 10 and between 3 and −1 for A = 100.For D = 10 the pattern of errors was hybrid between the two previous cases.Similarly, when K 1 = 10 (Figure 2) for D = 1 the values of the logarithm of the functional in Equation ( 6) were independent of q and ranged between about −2 and −12 when p increased from 0 to 10.For D = 100, the values of error were practically not dependent of p and ranged between 1 and −3 for A = 1, between 2 and −2 for A = 10, while between 3 and −1 for A = 100, when p ranged from 0 to 10. Furthermore, in this case, for D = 10 the error pattern was hybrid between the two previous cases.Finally, for K 1 = 100 the distribution of errors was similar (see Figure 3), while the magnitude of errors (absolute value) was slightly decreased for D = 1 and D = 10.
This kind of analysis was useful to learn about the distribution of errors among different combinations.However, because of the previously described dependency of parameters K 1 , A, D and the possible lack of dependency of p or q, several different combinations resulted in the same amount of error.As a consequence, a further sensitivity analysis was performed (see Figure 4) in order to clarify the influence of each parameter on the amount of error.In particular, the influence of the parameter K 1 was studied through the sensitivity index in Equation ( 7) and underlined a lack of dependency of q for p ≥ 2, where the sensitivity was the same for all combinations of p and q.In a similar way, the sensitivity index defined in Equation (8) showed a lack of dependency of q for p ≥ 1, where the sensitivity pattern was flat.Moreover, the sensitivity index in Equation ( 9) clearly showed the total lack of dependence of p together with the increase of q.In addition, in general, the order of magnitude of S D was 5.56 times S K1 , while it was 333.334 times S A .Thus, the D parameter had the maximum influence on the evolution of Equation (3), followed by K 1 and then by the parameter A. However, also in this case, several combinations had the same amount of sensitivity with respect to these parameters, which means this analysis was not enough to choose the best combination of p and q for Equation (3).
As a consequence, the metric in Equation ( 10) was used to account for at the same time both the value of the functional in Equation ( 6) and the values of p and q.Indeed, from a side, this metric was directly proportional to the absolute values of the logarithm of errors, which means it was able to found the combinations resulting in small errors.From the other side, this metric was inversely proportional to the values of p and q, then it was able to identify good candidates with small values of p and q.The use of this metric resulted also dependent of the values of parameters K 1 , A, D. Therefore, to perform an accurate selection of the best possible combinations related to a given biological tissue, previously identified values [27] were used, both for the porcine peroneal nerve and for the canine vagus nerve.This procedure resulted in the clear selection of the values p = 1 and q = 4 for the porcine peroneal nerve.The maximum difference between the evolution of the mean Cauchy stress in Equation ( 3) with p = 1 and q = 4 with respect to the same equation with p → ∞ and q → ∞ was within the 2.66%, as well as the ratio, which was 0.974.Therefore, Equations (3)- (5) were greatly simplified in their writing, keeping a good level of accuracy.On the contrary, for the canine vagus nerve the use of the Ω metric in Equation (6), resulted in the selection of the values p = 7 and q = 7.In this last case, the errors for the difference and the ratio between the reduced polynomial expression and the complete one were both in the order of 1 × 10 −9 , while together a single best combination, some other suboptimal combinations were found, as shown in Figure 5d (diagonal values).
Once simplified forms for Equations (3)-( 5) were provided, their evolution with strain was explicitly explored (Figure 6).Both for porcine peroneal nerve and canine vagus nerves, the increase of the mean Cauchy stress with strain was highly non linear (see Figure 6a,d), in accordance to the previous literature [20][21][22][23][24]27].This phenomenon was related to the strain stiffening ability of peripheral nerves, which was likely developed to preserve the integrity of their internal structures [33][34][35][36][37][38].However, to provide a more complete picture of the behaviour of peripheral nerves under axial elongations, the evolution of the tangent modulus was investigated.Furthermore, in this case, for both nerves, the tangent modulus evolved in a highly non linear way (see Figure 6b,e).In addition, the rate of change of the tangent modulus with strain was explored.This further analysis was able to underline the highly non linear instantaneous evolution of the tangent modulus with the strain (see Figure 6c,f).It could be worthy to notice that even if the non linear evolution of the stress with strain was quite well known in the literature [20][21][22], few works were provided to describe the evolution of the tangent modulus, as well as the evolution of the rate of change of the tangent modulus with strain [24].
Although, intuitively, a material should behave "like natural peripheral nerves" to provide a good level of structural biomimicry, it is not clear whether only the similarity of the stress/strain response is enough.To better investigate this point, a correlation analysis was performed for both the previously described responses of both the porcine peroneal and canine vagus nerves (see Figure 7).This analysis resulted in a highly non linear positive correlations between the tangent modulus and the mean Cauchy stress both for the porcine peroneal and the canine vagus nerves (see Figure 7a,d), as in Equations ( 11) and ( 14).Similarly, a positive non linear correlation was found between the rate of change of the tangent modulus and the mean Cauchy stress for both nerves (see Figure 7b,e), as in Equations ( 12) and (15).Finally, for the porcine peroneal nerve a linear correlation was found between the rate of change of the tangent modulus and the evolution of the tangent modulus (see Figure 7c), as in Equation ( 13), while for the canine vagus nerve the same correlation was non linear (see Figure 7f), as in Equation ( 16).The previous results supported the general need of a comparison with the complete dynamic state of the material, namely evolution of stress, tangent modulus, and rate of change of the tangent modulus with strain, to provide an acceptable level of biomimicry with peripheral nerves.The linear correlation found between the rate of change of tangent modulus and the evolution of tangent modulus for porcine peroneal nerve was likely due to the particular numeric values of parameters (i.e., K 1 , A, D), which were used in this work, in particular the low value of A. Indeed, different values of these parameters for the canine nerve were enough to show that a comparison with only two quantities (stress and tangent modulus with strain, due to the linear correlation between the rate of change of the tangent modulus and the evolution of the tangent modulus with strain) was not enough in general.

Conclusions
Peripheral nerves under stretch behave in a complex way, increasing their instantaneous stiffness.As a consequence, in this work, an original and complete study was provided to identify general and specific polynomial equations ruling the evolution of stress, tangent modulus and rate of change of tangent modulus with strain.Different analyses were provided to identify the simplest equations able to model the physiological behaviour of different nerve specimens (porcine peroneal and canine vagus nerves).The complete picture of the "static" and "dynamic" evolution of stress with strain was provided for these two different cases.The comparison with all these three quantities was shown to be the basic reference for an acceptable level of biomimicry with peripheral nerves.Indeed, material mimicking the mechanical behaviour of nerves could be strategic in established medical procedures [39,40], interactions with implantable devices [41,42], as well as in novel technological applications as neural interfaces [43? ], stretchable electronics [45], and direct connection between engineered biomaterials and nerves cells [46][47][48][49][50] or peripheral nerves [51,52].In addition, a good level of biomimicry with peripheral nerves could be crucial to control the immune response of the nervous tissue [53,54] in response to the stiffness mismatch between nervous tissue and biomaterials [45,55].Indeed, this mismatch is likely due to both the "static" actual difference between stress, but also to the "dynamic" evolution of this difference with strain.Therefore, the investigation of the evolution of the tangent modulus, and of the rate of change of the tangent modulus with strain could provide explicit and quantitative mechanical constraints for structural materials biomimetic to peripheral nerves.

Figure 5 .
Figure 5. Suitability assessment of combinations with finite values of p and q.(a)The metric Ω identified as best candidate the stress function with p = 1 and q = 4 (in yellow) for the peroneal nerve.(b) Evolution of the difference between stress functions with p = 1, q = 4, and p = ∞, q = ∞.(c) Evolution of the ratio between stress functions with p = 1, q = 4, and p = ∞, q = ∞.(d) The metric Ω identified as best candidate the stress function with p = 7 and q = 7 (in yellow) for the vagus nerve.(e) Evolution of the difference between stress functions with p = 7, q = 7, and p = ∞, q = ∞.(f) Evolution of the ratio between stress functions with p = 7, q = 7, and p = ∞, q = ∞.In all relevant plots zz = for sake of simplicity

Figure 6 .
Figure 6.Evolution of stress, tangent modulus, and rate of change of tangent modulus with strain.(a) Evolution of the mean Cauchy stress with strain for the porcine peroneal nerve.(b) Evolution of the tangent modulus with strain for the porcine peroneal nerve.(c) Evolution of the rate of change of the tangent modulus with strain for the porcine peroneal nerve.(d) Evolution of the mean Cauchy stress with strain for the canine vagus nerve.(e) Evolution of the tangent modulus with strain for the canine vagus nerve.(f) Evolution of the rate of change of the tangent modulus with strain for the canine vagus nerve. a