Automated and Data-driven Computational Design of Patient-Specific Biomechanical Interfaces

Biomechanical interfaces are mechanical structures that form the connection between a device and a tissue region, and, through appropriate load transfer, aim to minimize tissue discomfort and injury. A patient-specific and data-driven computational framework for the automated design of biomechanical interfaces is presented here. Optimization of the design of biomechanical interfaces is complex since it is affected by the interplay of the geometry and mechanical properties of both the tissue and the interface. The proposed framework is presented for the application of transtibial amputee prostheses where the interface is formed by a prosthetic liner and socket. Conventional socket design and manufacturing is largely artisanal, non-standard, and insufficiently data-driven, leading to discrepancies between the quality of sockets produced by different prosthetists. Furthermore, current prosthetic liners are often not patient-specific. The proposed framework involves: A) non-invasive imaging to record patient geometry, B) indentation to assess tissue mechanical properties, C) data-driven and automated creation of patientspecific designs, D) patient-specific finite element analysis (FEA) and design evaluation, and finally E) computer aided manufacturing. Uniquely, the FEA procedure controls both the design and mechanical properties of the devices, and simulates, not only the loading during use, but also the pre-load induced by the donning of both the liner and the socket independently. Through FEA evaluation, detailed information on internal and external tissue loading, which are directly responsible for discomfort and injury, are made available. Further, these provide quantitative evidence on the implications of design choices, e.g. : 1) alterations in the design can be used to locally enhance or reduce tissue loading, 2) compliant features can aid in relieving local surface pressure. The proposed methods form a patient-specific, data-driven and repeatable design framework for biomechanical interfaces, and by enabling FEA-based optimization reduces the requirement for repeated patient involvement in the currently manual and iterative design process.


I. Introduction
I N the United States, over half a million people live with lower limb loss [1] and 130,000 lower extremity amputations (LEAs) are carried out annually [2]. The lifetime healthcare cost after LEA is estimated to be $649,953 [3], Fig. 1. The biomechanical interface for transtibial amputees. Schematic of main tissue structures (right) and the liner and socket system (left). (modified from [5] with permission) A prosthetic liner is a soft sock-like layer which fits tightly around the residual limb. Despite variations in patient geometries and tissue conditions prosthetic liners are generally not patient-specific. Instead, a particular size and design is simply chosen from a range of commercially available liners. Although prosthetic sockets are patient-specific, their design and manufacturing process (Fig. 2) is presently a largely artisanal procedure (see also [6], [7]). The source of the socket geometry is obtained by wrapping a cast around the residual limb of the patient. A derived positive mold is then modified with the aim to remove load from regions that are deemed vulnerable while enhancing load at regions that are deemed safe. These regions are identified using manual palpation. Finally, a test socket is manufactured from the adjusted mold for evaluation with the patient. The adjustment and test socket evaluation process is then repeated until the patient can tolerate the loads on their limb, after which a final socket is manufactured. The success of this traditional socket design process relies heavily on the experience of a prosthetist, and requires manual and iterative design evaluation demanding repeated patient feedback. Fig. 2. Traditional artisan methods for prosthetic socket design. A plaster cast mold is created for the residual limb, and cut lines are manually defined (A), the mold shape is then manually adjusted to define the socket inner shape (B), after the vertical axis is determined an attachment plate is mounted (C), and then carbon fiber layers are wrapped over the mold (D) to produce the final socket design (E) The manual nature of the process means it is non-repeatable and currently largely non-data-driven, and quantitative data is either not obtained or insufficiently employed. As such, there is a reported discrepancy in the quality of sockets produced by prosthetists [8]. Furthermore, it has been reported that 57% of lower extremity prosthetic users experience moderate to severe pain when wearing their device [9]. Discomfort commonly results in skin problems and tissue damage (e.g. [10]- [13]). In severe cases when loading conditions cause tissue deformation thresholds to be exceeded (see also [14] on thresholds), painful pressure ulcers may occur [15]; in some reports pressure ulcers have occurred in as high as 55% of subjects with major amputations [16]. However, even mild discomfort may be concerning as it could result in an altered posture and gait, which in turn may cause long term musculoskeletal conditions such as back pain [17]. Moreover, any limitation in mobility can further contribute to conditions such as obesity, musculoskeletal pathologies including osteoarthritis, osteopenia, and osteoporosis, as well as cardiovascular disease ( [3], [17], [18]). A need in the field of prosthetics is a design and manufacturing framework for biomechanical interfaces based on a clear scientific rationale to maximize comfort and avoid tissue injury. Such a computational design and manufacturing process would provide an accurate, repeatable and fully patient-specific data-driven process, and can also be combined with virtual prototyping techniques for patient-specific design optimization. Virtual prototyping can be realized through finite element analysis (FEA), allowing for the detailed investigation of tissue pressures and internal deformations. FEA based optimization may potentially reduce the need for repeated test socket manufacturing and iterative patient involvement, and is therefore also able to reduce the overall cost and manufacturing time required. Besides optimization of shape in prosthetic design, the use of compliance may also add to comfort. For the design of comfortable shoes and footwear compliant materials have been an obvious choice. However, given that, in contrast to the human foot, the tissues of the residual limb are unevolved for loadbearing, it is surprising that for prosthetic interfaces, rigid materials (with respect to the soft tissue) have predominantly been explored. Some researchers however, have proposed compliant socket designs to offer relief in vulnerable regions, such as near bony protrusions. For instance, by varying the socket wall thickness and by introducing deformable structures [19], by introducing a variable spacing between a flexible inner and rigid outer socket [20], and by spatially varying the elastic material properties of the socket [21]. The preliminary findings of the latter study were reduced contact pressures for a compliant socket compared to a conventional socket. Advancements in the design and manufacturing process of sockets have been proposed. For instance, through the incorporation of computer-aided design (CAD) (e.g. commercial software [22]- [25]) and computer-aided manufacturing (CAM) technologies (see for instance [19], [21], [26], [27]). However, at present, these tools do not inform the design in a data-driven sense [28] since the actual design process remains a manual and experience based procedure. This may explain the reported preferential indifference among patients who used both a socket made using conventional and CAD/CAM techniques [29], and that design errors remain prevalent [30]. Further, non-invasive imaging has been used to study the geometry of the residual limb, e.g. based on magnetic resonance imaging (MRI) [31], [32] and ultrasound [33]. Some authors have proposed frameworks for socket design and evaluation based on non-invasive imaging and computational modeling. For instance, Papaioannou et al. 2011 [34] presents the use of dynamic roentgen stereogrammetry combined with image based modeling and FEA. Colombo et al. 2013 [32] and subsequent studies by the same group [35], [36] present the most detailed patient-specific socket design and evaluation framework to date. Although patient-specific geometries are derived from MRI, the socket designs are created in a computer aided but manual fashion based on experience and a-priori knowledge of manually inspected vulnerable and load-bearing regions. In addition, the above has been combined with FEA based socket design evaluation [37], [38], and socket evaluation using FEA (i.e. solely evaluation without computational design) is also presented in [39]- [41]. However, in all of these cases the soft tissue material behavior was modeled using linear elasticity which is not appropriate for analysis of large deformations. In addition, linear elasticity does not consider deformation induced stiffness enhancement due to the non-linear elastic nature of soft tissue. In order to accurately evaluate the deformations and loading conditions of the soft tissue inside a socket, the FEA should consider the following three loading effects: 1) Liner donning induced pre-loads, 2) socket donning induced pre-loads, and finally 3) loading occurring during functional use (e.g. standing/walking). Modelling of functional use often simply consists of the application of force (e.g. [19]) or displacement (e.g. [42]) boundary conditions (e.g. resulting in loading equivalent to supporting body weight). However, representation of the liner and socket induced pre-loading due to donning is far less trivial. Since the equilibrium shapes of the liner and socket do not match the undeformed soft tissue they create significant pre-strain and pre-stress. The associated large deformations also alter material stiffness and are capable of perturbing the degree of anisotropy due to the non-linear elastic properties of the soft tissue. The same may hold for the liner and socket materials if non-linear elastic materials are employed. Some researchers have attempted to account for socket donning induced pre-loads using prescribed (radial) displacements (e.g. [43]). However, these displacements likely create unnatural deformations since in reality the materials may displace not only normal but also tangential to each other. Faustini et al. [44] did not simulate pre-loading but aimed to account for deformation induced stiffness changes (due to tissue non-linear elasticity) by increasing the linear elastic stiffness in the undeformed geometry at the patellar ligament. However, with this approach the tissue remains undeformed and stress and strain free which is not realistic. Socket pre-loading effects have also been simulated by resolving socket-tissue overlap after placement using contact algorithms (e.g. [45]). Others have simulated a more complete donning process by using contact algorithms and simulation of gradual insertion of the limb inside the socket (e.g. [37], [46]). The approaches involving contact algorithms are more realistic than the use of prescribed displacements since the tissue is free to displace relative (including tangential) to the socket surface. However, simulation of such large motions combined with non-linear FEA based contact evaluation is computationally intensive. In addition, for these studies the socket material stiffness was several orders of magnitudes higher than the soft tissue. Hence the sockets encounter no or minimal deformations during the simulated donning or loading process. Therefore, to the authors knowledge, the accurate simulation (using large strain formulations and nonlinear FEA) of pre-loading has to date not been combined with the evaluation of significantly deformable and compliant socket designs. In addition, pre-loading of the soft tissue due to both a liner and a socket have to date not been investigated. In all the discussed approaches the socket design process is based on human experience, and design evaluation and optimization is manual, involving iterative refinement with repeated patient involvement. Ideally however, the design process should be driven by patient-specific data and quantitative measurements, and design evaluation and optimization should incorporate patient-specific data and computational modeling based virtual prototyping. There are several challenges to overcome to create such a framework: 1) non-invasive imaging is required to assess both external and internal tissue geometries, 2) realistic biomechanical material behavior needs to be considered, 3) design evaluation should employ detailed computational modeling to predict the patient-specific in-vivo tissue loading conditions (reducing the degree of patient involvement), 4) computational modeling should also include tissue pre-loading induced by both the liner and socket, and 5) production should employ CAM techniques to guarantee the fidelity of the design. To address the discrepancies of the traditional prosthetic interface design and manufacturing process, and explore the use of compliant materials, this paper presents a novel quantitative, data-driven and patient-specific socket design framework, which incorporates: 1) MRI for the generation of accurate patient-specific computation model geometries, 2) the use of non-invasive tissue mechanical property assessment based on indentation tests, 3) automated anatomical landmark and biomechanical behavior driven socket and liner design, 4) spatially varying socket design features such as donning induced pre-load, and socket material stiffness, 5) evaluation of liner and socket induced pre-loading, 6) finite element analysis based patient-specific design evaluation, 7) the 3D printing based manufacture. Using the novel framework 4 design strategies are compared in terms of predicted contact pressures and internal strains.

II. Methods
All data processing and visualization was conducted using   (Fig. 3E), 3D and patient-specific FEA models can be constructed, which can be used to generate custom liner and socket designs. The FEA process starts with the unloaded patient geometry and socket and liner source geometries. These source geometries are not yet referred as designs, similar to a prosthetist's patient cast, these source geometries are simply copies of the unloaded patient geometry and are created by offsetting from the patient skin surface. The source geometries therefore need to be adjusted to create the desired socket and liner designs. An automatic and datadriven adjustment process is proposed here using FEA based morphing and fitting on the virtual patient. During FEA, geometries are simultaneously morphed into designs and donned onto the patient. Morphing of the source geometries takes place using so called fitting pressures. These are pressure fields applied to the skin surface which deform the patient geometry, while the source geometries still lack any mechanical strength and are freely carried along with the skin motion (without developing stresses). Once a desired equilibrium design shape is obtained the geometries are "solidified" during FEA by assigning them with appropriate and stress free mechanical properties. Next, the fitting pressures are ramped down in FEA allowing the soft tissues to push back and relax into the devices. For compliant designs the tissue is now able to settle into and deform the devices until both are at (a donning induced pre-loaded) equilibrium. Using this process the designs are therefore initiated, fitted and donned onto the patient. Next, the designs are evaluated using FEA by applying body weight loading. The liner proposed here is manufactured by 3D printing a mold which can be used to create a patient specific liner, for instance using silicone rubber liner (Dragon Skin R 10 FAST, Smooth-On, Inc., Macungie, PA, USA). The sockets proposed here can be manufactured using 3D printing. For rigid sockets the designs can be printed as a single part, and from a single material. For compliant sockets a spatially varying socket wall material is used. Such designs can be realized using a multi-material printer (Connex 500, Stratasys Ltd., Eden Prairie, Minnesota, USA). For the compliant designs the socket system consists of an inner and an outer socket. The inner socket is compliant while the outer socket is rigid and provides additional strength. The outer socket can be bonded to the inner socket at sites where the inner socket is stiff while it can be offset from the inner socket at compliant sites to allow for socket deformation.

B. Non-invasive imaging
In order to capture the anatomical structure and geometry of the residual limb MRI is used. A male volunteer (age 48, weight 78 kg, activity level K3 [see also [2] on activity levels]) was recruited and placed prone and feet-first inside an MRI scanner (Siemens Magnetom 3 Tesla, Siemens Medical Systems, Erlangen, Germany). Informed consent was obtained and the research protocol was approved by the Committee on the Use of Humans as Experimental Subjects at Massachusetts Institute of Technology. All imaging was performed with a RF body coil and an Ultra-short TE MRI (UTE-MRI) sequence (e.g. [49]) was used, (T R /T E =5.71/0.07, acquisition matrix 192x192, 192 slices, field of view 220x220x220 mm, voxel size 1.145x1.145x1.145 mm) enabling visualization of bone tissue contours despite its short T 2 time. Several slices of the MRI data are visualized in Fig. 4A-B.

C. Obtaining patient specific geometries
In order to construct the detailed computational model, skin and bone contours were segmented (based on GIBBON [47] uiContourSegment function). Segmentation is applied to the raw image data and is based on adjustment, selection and combination of iso-contour lines. Once contours for a specific feature are recorded ( Fig. 4A-B) for each slice, these can be converted to iso-surface descriptions which are resampled geodesically at a desired density and smoothened (shrinkage avoiding smoothening was used [50]) ( Fig. 4C-D). Following segmentation, the geometries are reoriented in two steps. First a rotation is performed such that the femur and tibia are aligned in the feet-head direction (z-axis). This direction is identified using principal component analysis of the bone surfaces coordinates. Next, the model is rotated around the zaxis such that the anterior direction corresponds to the positive y-axis direction. This is achieved by rotating the geometry such that the XY-projection of the vector spanned between the center of the patella and the center of the femur, is most aligned with the y-axis. Fig. 4 shows how the image data and therefore the derived surface geometry does not extend far beyond the top of the patella. Since it was of interest to model beyond this region the surface geometry of the femur and skin was extended by 60 mm in the z-direction providing the extended geometry seen in Fig. 4D and elsewhere.

D. Creating the liner and socket source geometries
Similar to how the plaster cast in the traditional approach is used as source geometry to initiate the design process, the MRI derived patient geometry is used here. The source geometry for the liner is created by offsetting a layer from the skin surface. For the current study the liner thickness varied linearly from 4 mm to 6 mm from the top of the model to the distal end. Once the liner source geometry is defined the source geometry of the socket is formed by offsetting from the outer surface of the liner source geometry. However, first the socket upper boundary, known as the cut-line, is defined. Automatic cut-line creation is based on anatomical landmarks on the bones (Fig.  5A). The landmarks define curve control points on the outer liner surface (Fig. 5B), see also appendix A. Next, a smooth cubic spline is fitted to the curve control points and mapped to the liner surface to create the cut-line (Fig. 5B). The source geometry for the socket is then formed by offsetting the surface under the cut-line. A uniform socket wall thickness of 6 mm is used here and the socket cut-line rim was rounded with a rounding diameter matching the socket wall thickness. This creates the socket geometry shown in Fig. 5C.

E. Solid meshing
For FEA the following four material regions are modeled: 1) the bulk soft tissue (which includes skin, adipose and muscle tissue), 2) the patellar ligament, 3) the liner, and 4) the socket. Bones were not modeled as solid materials but were instead represented by rigidly supported voids. The solid material regions were meshed using a total of 146502 trilinear tetrahedral elements (see Fig. 6) using the free and open source meshing code TetGen (version 1.5.0, www.tetgen.org, see [51]) integrated within the GIBBON toolbox. The mesh density was biased based on proximity to the bones.

F. Constitutive modeling
The bulk soft tissue, the patellar ligament, and the liner are modeled as homogeneous materials. The socket is allowed to be spatially varying in mechanical behavior and can therefore be heterogeneous, i.e. each element may have different desired constitutive parameters. The non-linear elastic behavior of Fig. 3. Overview of the data-driven computational design framework. By segmenting MRI data (A), the patient-specific geometry is obtained (B), using anatomical landmarks the socket cut-lines can be automatically created (C), the liner and socket source geometries can be offset from the skin surface and can be meshed with the soft tissue to form a FEA model (D), assigned with indentation derived tissue mechanical properties (E), spatially varying design features, such as socket compliance and fitting pressure, can be defined using FEA based measures of tissue vulnerability (F), fitting pressure fields can be used to morph the liner and socket into their desired shape, while also pre-loading the tissue due to donning (G), the designs can now be evaluated for body weight loading (H), and optimal designs can be exported (I), for 3D printing based manufacturing (J). all materials is modeled using the following isotropic, and coupled hyperelastic strain energy density function [52]: The material parameters c and κ have units of stress and define a shear-modulus like and bulk-modulus like parameter respectively. The unitless parameter m sets the degree of nonlinearity. Finally, λ i are the principal stretches, and J is the volume ratio (determinant of the deformation gradient tensor).
The constitutive parameters used are listed in Table I. For the bulk soft tissue the parameters were identified from dedicated patient-specific indentation tests (see our recent study [53]).
Since viscoelasticity is not considered here only the elastic parameters are used. As is common for constitutive modelling of soft tissue, in Sengeh et al. 2016 [53] near incompressibility was assumed by using the equivalent of κ = 100 · c. However, since in our former study no tissue deformation or shape changes were recorded, the degree of compressibility of the tissue was not sufficiently validated. Further, residual limbs are known to be capable of changing volume due to loading and with use of sockets, across different time scales [54], [55]. Therefore to allow for realistic pressures and deformations κ = 18 · c was used here, which does allow for some volume change of the tissue. For the patient included in this study a cast of the limb while pressurized at 90 kPa was available.
Using FEA and the above parameters the response to such a pressure could also be simulated. The value for κ to use in the current study was estimated by altering it such that a similar degree of pressure induced skin displacement was qualitatively observed. For the patellar ligament the parameters were based on literature data for tensile testing of human patellar ligament tissue [56] (the reported linear elastic Young's modulus E = 660 MPa was used to set c = E 3 and m = 2). The patellar ligament, liner and socket materials were made relatively incompressible by setting κ = 100 · c MPa. For the liner and socket materials the parameters were identified using uniaxial compression experiments (data not shown). For rigid regions the stiffest available material was used with c = 1558 MPa. For compliant regions a spatial variation of socket materials is proposed which depends on the choice of design variation (see also section II-G). As listed in Table I

G. Controlling design features
The liner fitting pressure was set at a homogeneous 90 kPa. This pressure was manually determined by varying it until the mean skin surface pressure was qualitatively observed to be approximately 15 kPa, which was deemed a target donned liner skin surface pressure. For the sockets a more complex procedure is followed. A design map is defined, denoted by  (1) , p (2) , . . . , p (14) }, through which a smooth curve can be fitted (B). The source geometry for the socket is then formed by offsetting the region found under the curve by a desired thickness from the skin surface (C). D, with D ∈ [0, 1], which can be used to set the local socket fitting pressure and local socket element stiffness through a linear mapping. The spatially varying fitting pressures P for each triangular skin surface element can be determined using: Here p max and p min are the desired minimum and maximum fitting pressures. For the spatially varying mechanical properties of the socket the constitutive parameters c for each socket element C can then be derived from: Here c max and c min are the desired minimum and maximum c values. In principle the constitutive parameters can be continuously varied between the minimum and maximum level allowing for the creation of smooth parameter variations. However, as mentioned before, for the current 3D printer only 5 compliant material types are available (see Table I). By using the design map the spatial variation of fitting pressure and socket material parameters can each be controlled with two parameters (a desired minimum and a maximum). For iterative design optimization procedures the design map D can be made to evolve with each iteration and a different design map can be employed for the fitting pressure and the socket materials. For the traditional artisanal approach the socket source geometry (plaster mold) is used to inform the designs, and local shape adjustments are made manually based on knowledge of safe and unsafe regions. These regions are largely identified using palpation of the patient's limb. Fig. 7A-B show typically reported vulnerable and safe regions for loading (see also: [35]). Since palpation assesses a combination of local tissue stiffness and thickness (i.e. distance to bones), this assessment is here termed displaceability, i.e. the ability of the tissue to locally deform when loaded. A prosthetists design map is therefore based on experience and palpation. To create an automated assessment of displaceability FEA is used here. Displaceability is computed as the magnitude of skin surface displacement following the application of a homogenous pressure of 90 kPa (i.e. the response to the liner fitting pressure is used here). Fig. 7C shows FEA derived relative displaceability data (normalized total displacement). Reported vulnerable and safe regions are seen to correlate well with regions with low and high displaceability respectively. Fig. 7D shows that the displaceability data can also be mapped onto the socket elements (based on nearest neighbor interpolation). This mapped data can be used as a displaceability based design map to control socket design features. However, although this data appears informative of most of the safe and vulnerable regions it does not sufficiently highlight the patellar ligament region which is generally deemed safe for loading. The region at the patellar ligament (marked with a white ellipse in Fig.  7D) was therefore enhanced. The design map in Fig. 7D is denoted D 123 . Fig. 7E is a variation of this design map where the design map was reduced for elements close to the fibular head and distal end of the tibia. This reduction was informed by the fact that high pressures are observed here for simulations with the mapping D 123 . This adjusted design map is denoted D 4 , and can be viewed as the result of one manual design iteration with respect to the mapping D 123 . Both design maps have also been nulled at the cut-line rim (10 mm high) to create a comfortable rim (as nulled regions result in compliant materials and low fitting pressures). In order to study the effect of the socket shape (controlled by fitting pressures) and material properties, 4 different design strategies are evaluated: 1) a rigid socket designed using a homogeneous fitting pressure, 2) a rigid socket designed using a spatially varying fitting pressure, 3) a compliant socket designed using a spatially varying fitting pressure and featuring spatially varying socket stiffness, 4) the same as 3 but with added soft features over the distal end and fibular head. In terms of fitting pressures the concept of approach 1 is comparable to a total surface bearing (TSB) design, while approaches 2-4 are comparable to a patellar tendon bearing (PTB) design (see also [57]). The parameters used for these design variations are listed in Table II. A visualization of the resulting socket material and fitting pressure distributions is shown in Fig. 8. The rigid material regions in all designs are fully supported during FEA. For the compliant designs (variation 3 and 4) the rigid material is employed where the design map is >0.25. These regions are highlighted in red in the material visualizations of Fig. 8 (and are to be bonded to a rigid external socket). The compliant materials are then linearly mapped for the remaining unsupported regions. The fitting pressures are nulled at the rims and have one further adjustment for the spatially varying pressure designs, the fitting pressure at the patellar ligament (a load safe region) is further enhanced by multiplying the fitting pressure suggested by the map by a factor f pat (Table II).

H. FEA based design and evaluation
The FEA based design and evaluation procedure consists of 5 steps which are schematically illustrated in Fig. 9. During all steps the bone nodes were rigidly supported and therefore constrained from moving. In addition, the nodes of the top surface of the model were constrained from moving in the zdirection (but free to move in the x-and y-directions). The socket, liner, and tissue regions share nodes at each interface, simulating high friction tied interfaces. The liner and socket are each designed and donned in separate 2 step procedures. First fitting pressures are used to morph the geometries into desired designs. During the design phase the liner and socket are in a "ghosted" form i.e. they do not have significant stiffness and develop no significant stresses (hence shown as transparent in Fig. 9). Once their desired design is achieved they are assigned with natural mechanical properties in a stress free state. This process of morphing the designs (while the soft tissue is pre-loaded) without developing stresses in the liner or socket regions is achieved by modelling the liner and socket materials as multi-generational materials (see [58]). The liner and socket parameter c is made generation dependent in the following way: Where γ is a generation index (see also Fig. 9), c so f t denotes the c parameter for soft tissue, and c true denotes the true (physically realistic) c parameter. Since during the design phase the material can be made to have negligible stiffness (γ = 1) they remain in an effectively stress free state when the source geometry is morphed into a desired design. However, critically during this deformation, the soft tissue is pre-loaded and does build up material stresses. Hence effectively the multi-generational approach is used here to allow for the switching off of the liner and socket material properties during FEA based morphing, and the subsequent switching on of the liner and socket materials in a stress free state when the desired shape is obtained. Once the liner or socket are designed by the fitting pressures and in their second generation (γ = 2) they are able to develop stresses. The use of multi-generation materials in this way therefore forms a simple means of driving the designs of the liner and socket and simultaneously provides a means to simulate the pre-loading induced by donning.
In step 1 only the soft tissues are developing significant stresses (γ liner = 1, γ socket = 1), liner fitting pressures are applied loading and deforming the limb and freely morphing and carrying the liner source geometry with it to its desired design. In step 2 the liner material is assigned with its natural mechanical properties (γ liner = 2, γ socket = 1). The liner starts out stress free in its equilibrium shape but, as the liner fitting pressures are ramped down, the tissue experiences some relaxation and starts to push against the liner, deforming it until the tissue and liner reach equilibrium. Steps 1 and 2 therefore function to design and don the liner. During steps 1 and 2 the socket source geometry remained in its ghosted form, bonded to the liner and was carried along with it. In step 3 only the liner and soft tissues are able to develop significant stresses (γ liner = 2, γ socket = 1) and the socket fitting pressures are applied to the skin surface. The socket fitting pressures deform the limb and the liner, and morph and carry the socket to its desired design. In step 4 the socket material is assigned with its natural (and potentially spatially varying) mechanical properties (γ liner = 2, γ socket = 2). The socket starts out stress free in its equilibrium shape but, as the socket fitting pressures are ramped down, the tissue experiences some relaxation and starts to push against the liner and socket, deforming both until the tissue, liner and socket reach equilibrium. Steps 3 and 4 therefore function to design and don the socket. Finally, in step 5 the supported socket nodes are moved upward by 3 mm. Displacement controlled simulations are used and the fitting pressures listed in Table II were iteratively adjusted until the reaction force was 765.18 ± 2N (the force due to body weight). Once the 5 step FEA procedure is completed tissue loading measures indicative of tissue comfort or injury risk can be studied. In this study the final stresses and strains are output to derive skin surface pressure and tissue maximum shear strains.

III. Results
A. Image-based modeling and data-driven liner and socket design Fig. 3 shows the data-driven and automated design and computational modeling framework for development of patientspecific sockets and liners. For FEA the computational time for design and evaluation is currently 12 minutes (32Gb RAM, Intel Core i7-4910MQ CPU). Given this computational speed it is feasible to do iterative FEA based design optimization. The outcome of the framework is a set of CAD files allowing for computer aided manufacture. Fig. 11 and Fig. 12 illustrate production of liner and socket designs.

B. FEA based evaluation of patient-specific socket design strategies
This paper explores a single liner design and compares 4 socket design variations (see section II-G): Design 1 is a rigid socket with a constant fitting pressure, design 2 is a rigid socket with a spatially varying fitting pressure, design 3 is a compliant socket with a spatially varying fitting pressure, design 4 is the same as design 3 but with reduced fitting pressure and socket material stiffness at the fibular head and the distal end of the tibia where high pressures were observed.

IV. Discussion
This study presents a novel framework for the quantitative design, and computational evaluation, of patient-specific prosthetic liners and sockets for people with lower limb amputations. An overview of the components of the novel framework is presented in Fig. 3 and includes: 1) Non-invasive imaging (i.e. MRI) for the generation of accurate patient-specific computation model geometries, 2) automated anatomical landmark and biomechanical behavior driven socket design, 3) spatially varying socket design features such as fitting pressures and socket material stiffness, 4) FEA based patient-specific design evaluation 5) multi-material 3D printing based manufacture of the compliant and spatially varying stiffness designs. The design process is automated and driven by patient-specific The column on the right shows the model shaded towards total displacement (mm) to visualize the shape changes. Model regions are shown as opaque or transparent, respectively to denotes that they either do, or do not have significant mechanical properties assigned to them. data. The design is generated and evaluated using FEA. Evaluation is based on analysis of tissue loading during simulated standing on one leg, i.e. the application of a force equivalent to body weight.  Fig. 11. Liner manufacturing. The inner surface of the FEA derived liner design (at the end of step 1 of the FEA process) (A) can be exported to a CAD file (B), which can be 3D printed to serve as a liner mold (C) for silicone liner production (D), after donning the liner on (E) its shape qualitatively resembles that of the liner at the end of step 2 in the FEA process (F).
Since FEA-based design and design evaluation takes place in 12 minutes, the framework opens the door to iterative FEA based socket and liner design optimization. For optimization, the design controlling parameters can be updated in an iterative fashion while minimizing measures predictive of comfort, such as skin contact pressure and tissue strain. In order to study the effect of the socket shape (controlled by fitting pressures) and material properties, 4 different design strategies are evaluated: 1) a rigid socket designed using a homogeneous fitting pressure, 2) a rigid socket designed using a spatially varying fitting pressure, 3) a compliant socket designed using a spatially varying fitting pressure and featuring spatially varying socket stiffness, 4) the same as 3 but with reduced fitting pressures and socket material stiffness over the distal end and fibular head. Following body weight loading these approaches have varying outcomes in terms of skin surface pressure and tissue maximum shear strain, as shown in Fig.  13 and Fig. 14 respectively. From Fig. 13 it can be observed that a rigid socket design based on a homogeneous socket fitting pressure (design 1) presents with high pressures (close to 100 kPa) at many known vulnerable regions such as the distal end of the tibia, the front of the tibia and the fibular head. By instead using a spatially varying fitting pressure (design 2) the pressure can be enhanced at safe regions such as the patellar ligament and the calf region. Through these enhancements, relief is obtained for the front of the tibia and the top of the fibular head. However, high pressures remain at the lower portion of the fibular head and the region close to the distal end of the tibia. By utilizing not only spatially varying fitting pressures but also compliant materials such as in design 3, the contact pressure at the front of the tibia can be further reduced. However, high pressures remain at the lower part of the fibular head and at the distal end of the limb. If based on these findings the socket material stiffness and fitting pressures are reduced further for these regions (design 4), these observed pressures can be reduced. These results show the possible benefit, in terms of skin contact pressure, of enhancing loading at safe regions while reducing loading at vulnerable regions. Further, these results show that skin contact pressures can greatly, and selectively, be reduced by incorporating deformable and soft socket materials at vulnerable regions.  is used (design 2) some of the load-safe regions are loaded more increasing deformations in these safe regions. Although the maximum shear strain at the fibular head remains 0.5, a reduction to 0.4 is seen at the distal end of the tibia and fibula. The pattern remains similar for design 3 where compliant materials are employed, although the strains are reduced at the front of the tibia. Design 4 presents with a similar pattern, however the reduced fitting pressures and socket stiffness distally have further reduced maximum shear strains at the fibular head and distal end.
The framework presented provides, for the first time, a complete set of tools for patient-specific, data-driven, and automated design and evaluation of lower limb prosthesis. This forms a significant advancement on current procedures which often are insufficiently data-driven, manual, not repeatable (depends on prosthetist experience), and requires repeated patient involvement for design evaluation. In contrast, the presented framework allows data-driven and automatic procedures and offers the ability to perform virtual iterative design evaluation thereby reducing patient involvement. Since the entire pipeline, from MRI segmentation to FEA and CAD file export for 3D printing, is managed in a single automated framework, repeatability and geometric fidelity are guaranteed.
Current efforts towards non-invasive imaging, data-driven, patient-specific, and FEA based prosthetic socket design include the use of dynamic roentgen stereogrammetry [34] and MRI [32], [35], [36], [59], [60], [61], [38], [62]. However, these studies have all considered soft tissue elastic material as linear, which is not realistic since it does not capture deformation induced stiffness enhancement and anisotropy, and may not be compatible with finite strain formulations. In contrast the current study does feature non-linear elastic and finite strain formulations.
Donning induced pre-loading (whereby the liner and socket cause tissue to be loaded and in a state of stress prior to addition, e.g. standing induced, loading) is often not represented in FEA of prosthetic sockets. Some researchers have used induced displacements (e.g. [40]) or local isotropic stiffness enhancement (e.g. [44]) to account for the effects of donning. However, these approaches are not realistic since the actual displacements and (generally anisotropic) stiffness enhancements are not known. More realistic approaches utilize contact algorithms e.g. to resolve the socket-tissue overlap after socket placement (e.g. [63]) or to compute deformation during gradual insertion of the residual limb inside the socket (e.g. [46], [59], [61]). The latter appears to most closely resemble reality. However, contact simulations combined with non-linear analysis for such large relative motions are computationally intensive, especially if both the socket and tissue are deformable. In addition, even for the gradual insertion approach, the results may vary depending on the contact algorithm, the assumed friction conditions, and, most importantly, on the exact motion path of the limb. Rather than a gradual insertion, in reality the patient might push their limb inside the socket and move the limb in various directions to settle their limb inside the socket. Such settling motions might remove and alter tangential forces that develop during the initial large motion of the insertion. Hence there is no consensus as to what motion history to simulate for these donning simulations. For the donning studies mentioned the socket material was either rigid or the stiffness was several orders of magnitudes higher than the soft tissue. The current study uniquely combines the use of significantly deformable socket designs and considers donning induced pre-loading due to both the liner and the socket. The donning procedure used here is unique since it utilizes multi-generational materials (see also: [58]). Using this approach, the patient specific liner and socket source geometries can be morphed and generated using FEA such that the liner, socket and tissues are each in a preloaded and deformed state following the donning processes. The donning process here follows from the application of fitting pressure fields which are ramped down after the liner or socket layers have been defined fully, allowing the tissue to relax into the fitted liner and socket. This approach provides a computationally efficient means to simulate the automated design and donning procedures.
We will now address several limitations of the presented framework.
Currently the FEA procedures assume isotropic, and hyperelastic constitutive behavior for the soft tissues and engineering materials (socket and liner). In addition, the soft tissues are modelled as two regions only, the patellar ligament, and a grouped soft tissue region. In the future a multi-layered material structure, including a skin and adipose tissue layer, can be modelled. Anisotropy of tissues such as ligament and muscle tissue, can be incorporated through the use of material formulations with fibrous reinforcement (fibre directions for such formulations can be obtained from Diffusion Tensor MRI e.g. [64]).
The liner and socket designs are currently evaluated for static loading, i.e. body weight, along the length direction of the limb. However, during daily use the limb and socket are subjected to dynamic forces in access of body weight and at varying orientations. Investigations of the socket designs for a wider array of load cases is therefore also a topic of future investigations. For dynamic load evaluation the viscoelastic nature of the tissues (and engineering materials) will also need to be considered.
The non-linear elastic, viscoelastic and anisotropic material parameters can be derived from dedicated mechanical experiments. In the current work mechanical properties for soft tissue were based on our previous work [53] where an indentation device was used that recorded indentor force and indentor force boundary conditions. That study assumed incompressible behaviour. However, in the current study we relaxed that constraint and allowed volumetric deformations to occur since this provided the most realistic response. In future work the indentation experiment should be combined with tissue deformation assessment allowing for the evaluation of the volumetric deformations of the tissue. This can be achieved by employing MRI based indentation experiments [65], or fullfield deformation measurement using, for instance, 3D digital image correlation.
Currently bones are represented in FEA as rigid voids that are all kinematically linked and fixed in space. In addition, the patellar ligament is not bonded to patella and tibia. For the static evaluations considered here these assumptions may be acceptable but require further investigation. Further, for dynamic FEA relative bone motions, as well as potential muscle activation, need to be considered.
For the current FEA procedure nodes are shared across all interfaces. This simulates a non-slip and high friction interface. Although high friction is a common assumption it may be beneficial in the future to allow for a sliding contact implementation whereby a friction coefficient can be prescribed. In addition, the elements used at present are 4-node tri-linear tetrahedral elements. More accurate results are obtainable using non-linear element formulations. For instance, 10-node quadratic tetrahedral elements. Mesh refinement and the use of higher order elements is the topic of future research.
In the current framework the compliance of the socket is modulated by spatially varying its mechanical properties, while the liner is defined as isotropic and homogeneous. However, spatially varying materials can also be explored for the liner. Further, the use of anisotropic liner or socket materials may be pursued through, for instance, the use of fiber reinforced materials. The favored directions for such anisotropic materials can for instance potentially be related to local tissue deformation directions thereby offering anisotropic support/relief. In addition to spatial modulation of the elastic behavior of the materials, it is possible to spatially vary the thickness of the liner and or socket.
At present the spatially varying pattern of fitting pressures and socket material properties are based on a so called design map. Using such a mapping is convenient since only two parameters (a desired minimum and maximum value) are required to define each of the patterns. The design map is here based on relative displaceability (the normalized total displacement due to a homogeneous external pressure) which appeared to be reasonable estimate of load-safe and vulnerable regions. Other types of design maps incorporating refined detail of safe and vulnerable regions (e.g. taking into account nerves and other painful regions). may be explored in the future, and for iterative FEA based design optimization the design maps can be made to evolve which each iteration.
At present FEA based design evaluation focused on skin surface pressure and internal maximum shear strain. These have not been compared to experimental measurements for validation. Therefore, future work should implement sensors on realized designs in order to validate the predictability of the FEA model. For instance, pressure sensors can be placed inside the socket (e.g. see [66]) and printable sensor designs have also been proposed (e.g. [67]). Motion capture systems may also be employed to study and capture dynamic loading data which may serve as input for load evaluation and validation. In addition, loading experiments can be conducted within an MRI scanner to allow for in-vivo deformation evaluation [68].
In future work other model outcomes indicative of patient comfort can be studied and combined with design optimization. Here skin surface pressures and tissue maximum shear strains were presented but further research is required into what loading measures should be used for design optimization and are most predictive of patient comfort.
The presented framework will benefit amputees as it offers a fully data-driven and patient-specific design procedure. Due to the computational efficiency achieved, the framework can be combined with iterative design optimization. Through the use of FEA based evaluation, or virtual prototyping. By making the entire design and fabrication process repeatable and datadriven, the prosthetic device manufacturing process can be based on a scientific rationale, whereby the comfort outcomes can be clearly linked to design choices. The applications of the work presented are not limited to prosthetic interfaces for lower limb amputees. The framework can easily be adjusted for other biomechanical interfaces, such as optimization of interfaces for wearable devices, for the design of optimal support structures, pressure ulcer prevention padding, and footwear.

V. Conclusions
This paper presents a patient-specific and date-driven computational framework for the automated design of biomechanical interfaces. Based on non-invasive image data and segmentation patient specific tissue geometries are captured. These in turn are used to create 3D volumetric meshes suitable for FEA. FEA then enables the assessment of the local deformability of the tissue, which allows for the creation of design maps that can be used to drive the spatially varying features of the socket. After the source geometry for the liner and socket are created these are morphed, using FEA, into a desired design. Uniquely the FEA procedure controls both the design and mechanical properties of the devices, and simulates, not only the loading during use, but also the pre-load induced by the donning of both the liner and the socket independently. Through FEA evaluation detailed information on internal and external tissue loading, which are directly responsible for discomfort and injury, are available. By studying the outcomes of rigid and compliant design variations it is demonstrated that: socket shape and fitting pressures can be used to locally enhance or reduce loading. Further, it is shown that socket compliance can be used to locally relieve surface pressure and internal strain. The patient-specific, data-driven and design framework for biomechanical interfaces, forms a repeatable design process based on a scientific rationale, which through virtual prototyping reduces iterative patient-involvement.

Appendix A Socket cut-line definition
The socket cut-line curve is generated by ray tracing line directions from anatomical landmarks outwards. The notation u = v −y refers to the intersection point u of a ray and a surface, where the ray starts at the position v in the direction of −y, i.e. the negative y-axis direction. w f emur = min(g y , f y ) − e y (9)