Porous Structure Design in Tissue Engineering Using Anisotropic Radial Basis Function

12/06/2016 ∙ by Ke Liu, et al. ∙ University of Wisconsin-Milwaukee 0

Development of additive manufacturing in last decade greatly improves tissue engineering. During the manufacturing of porous scaffold, simplified but functionally equivalent models are getting focused for practically reasons. Scaffolds can be classified into regular porous scaffolds and irregular porous scaffolds. Several methodologies are developed to design these scaffolds. A novel method is proposed in this paper using anisotropic radial basis function (ARBF) interpolation. This is method uses geometric models such as volumetric meshes as input and proves to be flexible because geometric models are able to capture the characteristics of complex tissues easily. Moreover, this method is straightforward and easy to implement.



There are no comments yet.


page 10

page 11

page 12

page 13

page 14

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1. Introduction

In order to improve biological tissues, tissue engineering (TE) which uses a scaffold to form new tissues for a medical purpose has a wide range of applications. Part of the applications in practice are repairing or replacing portion of or whole tissues and performing specific biochemical functions. Because of the inherent ability to produce customized porous scaffolds with different required architectures, the development of additive manufacturing (AM) techniques during last decade greatly improves tissue engineering. The latest ASTM standards defines AM as ”a process of joining materials to make objects from three-dimensional (3D) model data, usually layer upon layer, as opposed to subtractive manufacturing methodologies” [1]. More specifically, additive manufacturing starts from a 3D computer model and builds the final product by the addition of material, usually from a layer-by-layer fashion. This is a new manufacturing techniques comparing to conventional subtractive processes which removes material from a 3D block. Commercial AM techniques to fabricate scaffolds for tissue engineering applications include selective laser sintering (SLS), stereolithography (SLA), fused deposition modeling (FDM), precision extrusion deposition (PED) and 3D printing (3DP). Interested readers can refer to [2, 3, 4] for details of these techniques.

Because native tissues are inherently heterogeneous and often have complex physiological architectures, literature, in practice, is primarily focused on the manufacturing of models which are simplified but functionally equivalent to the tissue to be repaired in terms of porosity and mechanical properties. Two types of porous scaffolds, namely regular porous scaffolds and irregular porous scaffolds, can be designed to achieve this goal. Numerous methodologies are proposed and categorized to fabricate these two types of scaffolds [5]. Table 1 categorizes these methods.

Periodic porous structures have a limitation that slight local modifications can affect the entire structure globally. Our proposed method provides a new implicit approach to generate porous tissue scaffold through volumetric mesh. Thus scaffold architectures can be adjusted by only modifying the mesh. Our proposed method has three major advantages over other implicit methods like TPMS-based ones. Firstly, local modifications of pore shape, size or distribution is achieved by changing local mesh accordingly, which is easy because there are only geometric changes in the mesh. Secondly, our method is flexible to simulate heterogeneities and discontinuities in natural tissue structures by using purposely-designed mesh as input. Depending on the features of tissue structure, meshes with different type (such as tetrahedron mesh and hexahedron mesh), size and density can be used to represent these characteristics. With different meshes as input, our method is able to build different tissue scaffolds with slight modifications in algorithms. Thirdly, unlike many implicit methods need post-actions like Boolean operations to get the final pieces built, the only post-action in our method is taking iso-surfaces, which is easier, faster and more flexible that user can get different scaffold architectures by taking different iso-values. In general, our method is superior in flexibility and easy to implement.

The rest of paper is organized as the following. Section 2 summarizes progresses so far of methods categorized in Table 1. Section 3.1 introduces porous scaffold reconstruction using conventional RBF interpolation. The proposed scaffold reconstruction is presented in Section 3.2. Section 3.3 shows the overall algorithms of our proposed method. Finally, Section 4 shows some experimental results and discussions. Section 5 concludes this paper.

Type of scaffolds Method
Regular porous scaffolds CAD-based methods
Image-based methods
Implicit surface modeling (ISM)
Space-filling curves
Irregular porous scaffolds An optimization method proposed by [19]
Stochastic methods using Voronoi models [20, 21]
A hybrid Voronoi-spline method [22]
Methods using volumetric meshes [25]
Table 1. Category of methods to design porous scaffolds in tissue engineering.

2. Prior Works

CAD-based methods, such as constructive solid geometry (CSG) and boundary representation (B-Rep), are used to design regular porous scaffolds. CSG-based tools combine standard solid primitives (cylinders, spheres or cubes) through Boolean operations (e.g. intersection) to design and represent complex models. B-Rep tools describe the solid cell through its boundaries by a set of vertices, edges and loops without explicitly specify relations among them. So a preliminary check is required to verify there are no gaps or overlaps among the boundaries [6]. However, as objects become large or their internal architectures become more complex, their size increases dramatically hence it is hard or impossible to visualize and manipulate them. To overcome this limitation of most CAD-based tools, different solid cells with more bio-inspired features have been introduced [7, 8].

Image-based methods combine imaging, image processing and free-form fabrication techniques to simplify scaffold design. Scaffolds can be described by 3D binary images (i.e. voxel values are Boolean and correspond to ”solid” and ”void”). Image-based methods produce scaffolds by taking the intersection of two 3D binary images, one representing the shape to be reproduced, and the other consisting of stacking of a binary unit cell. Empirically derived geometries are created in the unit cell with basic geometric shapes (cylinders, spheres) to represent regular pores within a scaffold. Randomly arranged pores can be obtained by the use of a random number generator to set voxel state. The topological optimization algorithms has been proved pivotal to obtain scaffolds in image-based methods [9, 10].

Implicit surface modeling is highly flexible and describes scaffold architecture by a single mathematical equation, with freedom to introduce different pore shapes, pore size gradients. Recently, a large class of periodic minimal surfaces methods such as triply periodic minimal surfaces (TPMS) have become attractive for the design of biomorphic scaffold architectures. An early attempt of using TPMS-based method to control tissue fabrication is presented by [11]. TPMSs like Schwartz’s Primitive (type P), Schwartz’s Diamond (type D) and Schoen’s Gyroid (type G) are demonstrated their efficacy in high-precision fabrication of TPMS-based scaffolds [12, 13, 14]. However, all aforementioned works are limited to simple cubic or cylindrical outer shape. An improved method for constructing a pore network within an arbitrary complex anatomical model has been developed and successively optimized by Yoo et al. [15, 16]. In general, the TPMS methodologies, combining the advantages of both traditional CSG and image-based methods, are computationally efficient for modeling and fabrication of scaffolds.

Space-filling curves methodologies coupled to extrusion-based techniques. Such techniques consist of the extrusion of a micro-diameter polymeric filament terminating with a nozzle having an orifice diameter in the hundreds of microns range. The fabrication process involves the deposition of polymeric layers, which adhere to each other by heating temperature while retaining their shape. This process leads to regular repetition of identical pores. Thus these geometries have been named honeycomb-like patterns [17]. More complex patterns can be obtained by changing the deposition angle between adjacent layers. An alternative approach to design scaffold is using fractal space-filling curves, which can be mathematically generated by starting with a simple pattern that grows through the recursive rules.

Periodic porous structures have several advantages. They are easier to model and their structural properties are possible to predicate. Their disadvantages lie in the difficulty of controlling the pore shape, size and distribution since slight modification of the unit cell will pose global changes to the entire structure. Moreover, current CAD tools are not suitable to reproduce the complex natural structures. In scaffold with variational porous architectures, discontinuities of deposition path planning are often found at the interface of two adjacent regions [18]. To design such scaffold architecture, Khoda et al. implemented an optimization method [19]. Stochastic and Voronoi models have been used to generate random pores in scaffold design as well. Heterogeneous pores distributed according to a given porosity level are generated by stochastic methods in scaffold design [20, 21]. To overcome the limitation of only simple spheres can be used to represent pores, a hybrid Voronoi-spline representation combined with a random colloid-aggregation model is proposed [22]. The proposed method has been extended to implement graded pore sizes and pore distributions [23]. Volumetric mesh generators derived from finite element tools are used to create heterogeneous porosity within a solid model as well [25].

3. Methods

3.1. Review of Radial Basis Function (RBF) Interpolation

The conventional radial basis function interpolation is given by


where the interpolated function is represented as a weighted sum of radial basis functions , each centered differently at and weighted by . Let . By given conditions , the weights can be solved by


where . Once the unknown weights are solved, the value at an arbitrary voxel in porous architecture can be calculated by


In conventional RBF, the distance between point and center is measured by Euclidean distance. Let , commonly used radial basis functions include:
Multiquadrics (MQ):
Inverse MQ (IMQ):
Thin Plate Spline (TPS):
where is a shape parameter. Shape parameter plays a major role in improving accuracy of numerical solutions. In general, the optimal shape parameter depends on densities, distributions and function values at mesh nodes. Choosing shape parameter has been an active topic in approximation theory [26]. Interested readers can refer to [27, 28, 29, 30, 31] for more details.

After interpolation, iso-surfaces are taken on the interpolated piece to create porous architecture. In detail, for both 2D and 3D meshes, mesh vertices are given values 1. For 2D meshes, edge centers and tile (triangle or quadrangle) centers are given values -1. For 3D meshes, face centers and sub-volume (tetrahedron or hexahedron) centers are given values -1. Fig. 1 (a-c) shows a sample 2D triangle mesh, a 3D tetrahedron mesh and a 3D hexahedron mesh with assigned values, respectively. For simplicity of drawing, 2D RBF interpolation scheme is illustrated by Fig. 2 (a).

Conventional RBF seems viable to construct porous scaffolds. However, since RBF is isotropic in the sense that only geometrical distance is considered, the support domains of underlying basis function always tend to be circular (in 2D) or spherical (in 3D), which sets an limitation to the customization of the internal architecture, especially for complex tissue scaffolds. Therefore, conventional RBF is only viable to generate simple architectures. Openings in desired scaffold architecture should start from sub-volume center and grows towards pores on structure surface. Given the condition that the face centers and sub-volume centers are assigned of value -1, the interpolated voxels should have values close to -1. Therefore, the internal opening directions have to be considered during interpolation and the shape of support domain should thus be anisotropic.

Figure 1. Values assigned to mesh nodes. Red dots represent value of 1. Blue dots represent value of -1. In 3D meshes, interior dots are represented by lighter colors. (a) Sample 2D triangle mesh. Note both edge and triangle centers are given values -1. (b) Sample 3D tetrahedron mesh. Small triangle represent tetrahedron centers. (c) Sample 3D hexahedron mesh. Small triangle represent hexahedron centers.
Figure 2. 2D interpolation schemes. x is the pixel to be interpolated. Dashed circle (for RBF) or ellipses (for anisotropic RBF) are support domains of underlying basis functions. (a) RBF interpolation. (b) Anisotropic RBF interpolation.

3.2. Anisotropic Radial Basis Function (ARBF) Interpolation

The main difference between the isotropic and aniso-tropic RBFs is the definition of distance used. Given distinct line segments the anisotropic radial basis function is defined by


where is defined as the distance between any arbitrary point and a line segment or the distance between two line segments.

To calculate the distance between point x and line segment , there are three cases to consider. For case 1, if point x is on the line segment , distance is 0. For case 2, if point x, and form an acute triangle, the distance is defined as the length of x’s projection to . For case 3, if x, and form an obtuse triangle, the distance is defined as . Fig. 3 (a-c) illustrates the three cases. The distance between any two arbitrary line segments and is defined as
. Fig. 3 (d) illustrates this case. Table LABEL:tab:types_distances summaries these types of distances described above. Geometrically, the new distance changes support domain of basis function, from conventional circular or spherical shape to elliptical or hyper-elliptical shape. In other words, the new distance definition considers both magnitude and direction of line segment, which allows it to be a great tool to control the direction of internal structure. Customization of porous architecture can be achieved by measuring the distance between any arbitrary voxel and different line segments along with the distance between this voxel to other mesh nodes. Table LABEL:tab:distance_definition shows the three types of distance definition used during interpolation. In our experiments, the distances between sub-volume centers and tile centers are used. Fig. 2 (b) illustrates this idea in simplified 2D triangular mesh. Comparing 2D and 3D meshes, sub-volume centers in 3D map to the tile centers in 2D, and tile centers in 3D map to the edge centers in 2D. As the figure shows, x is the pixel whose intensity is unknown. o is the center of a tile. a, b, c are the centers of three edges of a triangle. The distance between x and the three line segments (a, o), (b, o), (c, o) are calculated. Additionally, distances between the three line segments, namely, (a, o), (b, o) and (c, o) themselves are calculated as well. The final interpolation matrix is built on these distances along with the Euclidean distance between x to other mesh nodes. Similar to isotropic RBF, with a modified distance, the ARBF interpolation problem becomes


Please note that matrix in equation (2) should also be updated accordingly with the new distance metric. Therefore, the new set of weights would be different from the weights in isotropic interpolation. After interpolation, a proper iso-value is applied to the interpolated pieces to get the final result. Fig. 4 and Fig. 5 show the ARBF interpolation results.

Figure 3. Cases to calculate anisotropic distance. (a) Point x is on line segment . (b) Point x and line segment form an acute triangle then distance is defined as the length of . (c) Point x and line segment form an obtuse triangle. (d) Distance between line segment and .

3.3. Algorithms

The following algorithm illustrates major steps to build porous structure from triangular (2D), tetrahedron and hexahedron (3D) meshes. The primary step is ARBF interpolation which consists of three sub-steps. Firstly, values are assigned to mesh nodes as explained in Section 3.1 to build the matrix. Secondly, the weight coefficients are solved by using the new distance definition explained in Section 3.2. Thirdly, the weights are used to interpolate the value at each voxel. After the piece is interpolated, an iso-surface is chosen and applied to get the final porous scaffold.

Algorithm: Porous Scaffold Construction

    loadMesh(); // read mesh
    ARBFInterpolation(); // do ARBF interpolation
    marchCube(); // take iso-surface
    exportResult(); // print result

    assignMeshValues(); // assign mesh nodes and build matrix
    solveCoefficients(); // solve unknown weights
    applyCoefficientsToInterpolation(); // do the actual interpolation

4. Results and Discussion

This section includes results of scaffold shapes using different types of meshes and shows several experimental results. Fig. 4 (a) shows a scaffold obtained from single tetrahedron. The opening grows from tetrahedron center toward the four triangle centers. 4 (b) illustrates an icosahedron mesh comprised of 20 tetrahedrons. 4 (c) shows the scaffold interpolated from 4 (b). Fig. 5 (a) shows a scaffold obtained from single hexahedron. Fig. 5 (b) shows a scaffold obtained from a block-shape hexahedron mesh which is comprised of 8 smaller hexahedrons. Fig. 5 (c) shows a scaffold obtained from a rod-shape hexahedron mesh which is comprised of 4 smaller hexahedrons. Input meshes for these results are regular so the final scaffolds are regular as well. The results are interpolated by ARBF interpolation introduced above using inverse multiquadrics (IMQ) as basis functions and shape parameter of value 0.1. After interpolation, a proper iso-value is applied to the interpolated piece to show the final porous structure. Because the input meshes are regular, the output structures tend to be regular as well. Iso-values can be adjusted to get structure of different size of pores. Results of different iso-values are included in Fig. 6. As the figure shows, pore sizes increase as the iso-values increases. So after interpolation, choosing a proper iso-value is required to obtain the desired porous structure. Besides iso-values, basis function can also affect the porous architecture in terms of mostly the size and shape of pores. Because the shape of basis functions are different, when the final piece is interpolated, different voxels are included in their support domains. Interpolated results using different basis functions are included in Fig. 7. To compare isotropic RBF interpolation and anisotropic RBF interpolation, a 2D isotropic RBF interpolated result (Fig. 8 (a)) is also included. As explained before, because basis functions in isotropic RBF have circular support domains, the circle-shape artifacts can be seen in the result. Additionally, interpolation results based on disturbed meshes are also investigated and shown in Fig. 8 (b-d). These disturbed meshes are generated by moving some mesh vertices in random direction. As the results shown, disturbed porous structures can be obtained from disturbed meshes. Therefore, adding disturbance is a viable way to construct porous structure with randomness, which is very useful to construct tissues with complex physiological architectures. Finally, as a contrast, TPMS results are also included in Fig. 9.

Figure 4. Results based on tetrahedron meshes. (a) Scaffold based on single tetrahedron. (b) Input icosahedron mesh (formed by 20 tetrahedrons). (c) Scaffold using (b) as input.
Figure 5. Results based on hexahedron meshes. (a) Scaffold based on single hexahedron. (b) Scaffold based on 8 hexahedrons arranged to form a large cube. (c) Scaffold based on 4 hexahedrons arranged to form a rod shape.
Figure 6. Results taken by different iso-values. (a) - (d) are results taken by increasing iso-values.
Figure 7. Results obtained by using different basis functions in ARBF interpolation. Shape parameter is 0.1 for all results. (a) Result interpolated by multiquadrics (MQ) basis. b) Result interpolated by inverse multiquadrics (IMQ) basis. (c) Result interpolated by Gaussian basis. (d) Result interpolated by thin plate spline (TPS) basis.
Figure 8. Some experimental results. (a) Result interpolated by isotropic RBF and based on a 2D triangle mesh. (b-d) Results based on hexahedron meshes with disturbance.
Figure 9. Porous structure obtained by TPMS method. (a) Structure obtained by P-type function. (b) Structure obtained by D-type function. (c) Structure obtained by G-type function. (d) Structure obtained by IWP-type function.

5. Conclusions

In practical, constructed scaffold may be very complicated to simulate complex physiological tissue architecture in terms of equivalent internal connectivity and mass transportation. Periodic porous structure cannot meet this challenge satisfactorily. Our proposed method uses volumetric mesh as input which can be very complex C thanks to modern mesh modeling techniques. Thus, using complex mesh as input, our proposed method is able to construct complicated porous scaffold structures to meet this challenge. Moreover, modifications to the final structure can be achieved by common mesh operations or changing iso-value, which makes our method very flexible comparing to other period porous manufacturing techniques (like implicit surface methods). Finally, implementation of our proposed method is easy and computational cost is low because the core algorithm of interpolation is calculating distances and solving unknown weights. There are a lot of mature and fast linear algebra libraries available to use. In the future, this method will be tested with adaptive meshes which represent porous architectures with heterogeneous and discontinuous structures. Moreover, a criteria to measure final scaffold will be developed to further test the efficacy and efficiency of this method.


  • [1] ASTM International: ASTM F2792-10: standard terminology for additive manufacturing technologies. West Conshohocken, PA: ASTM International (2010)
  • [2] Melchels FPW., Domingos MAN., Klein TJ., Malda J., Bartolo PJS., Hutmacher DW.: Additive manufacturing of tissues and organs. Progress in Polymer Science 37 (2012) 1079-1104
  • [3] Bartolo PJS., Almeida H., Laoui T.: Rapid prototyping and manufacturing for tissue engineering scaffolds. International Journal of Computer Applications in Technology 36 (2009) 1-9
  • [4] Peltola SM., Melchels FPW., Grijpma DW., Kellomaki M.: A review of rapid prototyping techniques for tissue engineering purposes. Annals of Medicine 40 (2008) 268-280
  • [5] Giannitelli SM., Accoto D., Trombetta M., Rainer A.: Current trends in the design of scaffolds for computer-aided tissue engineering. Acta Biomaterialia 10 (2014) 580-594
  • [6] Chiu WK., Yeung YC., Yu KM.: Toolpath generation for layer manufacturing of fractal objects. Rapid Prototyping Journal 12 (2006) 214-221
  • [7] Sun W., Starly B., Nam J., Darling A.: Bio-CAD modeling and its applications in computer-aided tissue engineering. Computer-Aided Design 37 (2005) 1097-1114
  • [8] Bucklen BS., Wettergreen WA., Yuksel E., Liebschner MAK.: Bone-derived CAD library for assembly of scaffolds in computer-aided tissue engineering. Virtual and Physical Prototyping 3 (2008) 13-23
  • [9] Hollister SJ., Maddonx RD., Taboas JM.: Optimal design and fabrication of scaffolds to mimic tissue properties and satisfy biological constraints. Biomaterials 23 (2002) 4095-4103
  • [10] Hollister SJ.: Porous scaffold design for tissue engineering. Nature Materials 4 (2005) 518-524
  • [11] Rajagopalan S., Robb RA.: Schwarz meets Schwann: design fabrication of biomorphic and durataxic tissue engineering scaffolds. Medical Image Analysis 10 (2006) 693-712
  • [12] Seck TM., Melchels FPW., Feijen J., Grijpma DW.: Designed biodegradable hydrogel structures prepared by stereolithography using poly(ethylene glycol)/poly(D,L-lactide)-based resins. Journal of Controlled Release 148 (2010), 34-41
  • [13] Melchels FPW., Feijen J., Grijpma DW.: A poly(d, l-lactide) resin for the preparation of tissue engineering scaffolds by stereolithography. Biomaterials 30 (2009) 3801-3809
  • [14] Elomaa L., Teixeira S., Hakala R., Korhonen H., Grijpma DW., Seppala JV.: Preparation of poly(e-caprolactone)-based tissue engineering scaffolds by stereolithography. Acta Biomaterialia 7 (2011) 3850-3856
  • [15] Yoo D-J.: Computer-aided porous scaffold design for tissue engineering using triply periodic minimal surfaces. International Journal of Precision Engineering and Manufacturing 12 (2011) 61-71
  • [16] Yoo D-J.: Porous scaffold design using the distance field and triply periodic minimal surface models. Biomaterials 32 (2011) 7741-7754
  • [17] Zein I., Hutmacher DW., Tan KC., Teoh SH.: Fused deposition modeling of novel scaffold architectures for tissue engineering applications. Biomaterials 23 (2002) 1169-1185
  • [18] Kalita SJ.: Development of controlled porosity polymer-ceramic composite scaffolds via fused deposition modeling. Materials Science and Engineering 23 (2003) 611-620
  • [19] Khoda AKMB., Ozbolat IT., Koc B.: Engineered tissue scaffolds with variational porous architecture. Journal of Biomechanical Engineering 133 (2011) 011001
  • [20] Schroeder C., Regli WC., Shokoufandeh A., Sun W.: Computer-aided design of porous artifacts. Computer-Aided Design 37 (2005) 339-353
  • [21] Sogutlu S., Koc B.: Stochatic modeling of tissue engineering scaffolds with varying porosity levels. Compute-Aided Design & Applications 4 (2007) 661-670
  • [22] Schaefer DW., Keefer KD.: Structure of random porous materials: silica aerogel. Physical Review Letters 56 (1986) 2199-2202
  • [23] Kou XY., Tan ST.: Microstructural modeling of functionally graded materials using stochastic Voronoi diagram and B-spline representations. International Journal of Computer Integrated Manufacturing 25 (2012) 177-188
  • [24] Cai S., Xi J.: A control approach for pore size distribution in the bone scaffold based on the hexahedral mesh refinement. Computer-Aided Design 40 (2008) 1040-1050
  • [25] Melchels FPW., Wiggenhauser PS., Warne D., Barry M., Ong FR., Chong WS., et al.: CAD/CAM-assisted breast reconstruction. Biofabrication 3 (2011) 034114
  • [26] Wang, J., Liu, G.: On the optimal shape parameters of radial basis functions used for 2-d meshless methods. Computer Methods in Applied Mechanics and Engineering 191 (2002) 2611–2630
  • [27] Divo, E., Kassab, A.: An efficient localized RBF meshless method for fluid flow and conjugate hear transfer. ASME Journal of Heat Transfer 129 (2007) 124–136
  • [28] Kosec, G., S̃arler, B.: Local RBF collocation method for Darcy flow. Computer Modeling in Engineering and Sciences 25 (2008) 197–208
  • [29] S̃arler, B., Vertnik, R.: Meshfree explicit local radial basis function collocation method for diffusion problems. Computers and Mathematics with Applications 51 (2006) 1269–1282
  • [30] Vertnik, R., S̃arler, B.: Meshless local radial basis function collocation method for convective-diffusive solid-liquid phase change problems. International Journal of Numerical Methods for Heat and Fluid Flow 16 (2006) 617–640
  • [31] Vertnik, R., S̃arler, B.: Solution of incompressible turbulent flow by a mesh-free method. Computer Modeling in Engineering and Sciences 44 (2009) 65–95