Modelación de andamios porosos basados en las estructuras triplemente periódicas P y G

Introducción

Los andamios son estructuras de soporte sólidas utilizadas en injerto de tejidos y que tienen varios requisitos relacionados con sus propiedades mecánicas, biológicas y físico-químicas que les dan la capacidad para facilitar la creación de tejidos. Para diseñar andamios para regeneración de tejido óseo hay que tener en cuenta varios factores como (a) que tengan una red totalmente interconectada de poros que favorezca el transporte de sustancias y nutrientes (Nikolova & Chavali, 2019), (b) que sean biocompatibles y biodegradables, (c) que tengan las propiedades mecánicas apropiadas como la rigidez y la resistencia mecánica teniendo en cuenta el uso que se desee dar al producto (Kanczler et al., 2020) y (d) que la morfología interna de los poros permita el crecimiento celular y la formación de nuevo tejido (Knychala et al., 2013; Zadpoor, 2014).

Existen diversos métodos para diseñar estructuras de andamios porosos. Entre ellos tenemos los métodos CAD para generar la celda unitaria a partir de sólidos elementales (Ahmadi et al., 2015; Naing, Chua, Leong, & Wang, 2005), el método de intersección de imágenes binarias (Hollister, 2005; Saito et al., 2010; Wu, Spanou, Diez-Escudero, & Persson, 2019), la modelación de la trayectoria de entramado en el método de fabricación por filamento fundido o deposición fundida (FDM) (Amirkhani, Bagheri, & Zehtab Yazdi, 2012; Enrique Cuan-Urquizo et al., 2019; E. Cuan-Urquizo, Yang, & Bhaskar, 2015; Giannitelli, Accoto, Trombetta, & Rainer, 2014; Lee et al., 2012; Malinauskas et al., 2015; Moroni, de Wijn, & van Blitterswijk, 2006; Norato & Johnson, 2011) y el método de superficies implícitas que permite diseñar estructuras celulares con alto grado de complejidad (Gandy, Bardhan, Mackay, & Klinowski, 2001) como las superficies minimales triplemente periódicas (SMTP).

Se denominan superficies minimales o extremales mínimas a las que satisfacen el problema de cálculo de variaciones de hallar la superficie orientable cerrada de menor área con frontera limitada por una curva cerrada espacial (Meeks III, 1990). Aunque existen muchas de estas superficies, las más estudiadas por lo interesantes de sus propiedades y de mayores aplicaciones técnicas, son las SMTP que son periódicas en las tres direcciones espaciales independientes, no se autointersecan y dividen al espacio tridimensional donde ellas están en dos sub-volúmenes conexos que no se interpenetran entre sí y que se denominan laberintos (D. Anderson, Davis, Nitsche, & Scriven, 1990). El método implícito de diseño computacional de estas superficies se basa en el truncamiento del desarrollo de Fourier de la superficie que tiene la periodicidad espacial de alguna red cristalina (Gandy et al., 2001), lo que resulta en una ecuación implícita para cada estructura, con un grupo de parámetros que permiten tener mayor control del tamaño de poros y la porosidad. Entre las superficies más conocidas se encuentran la Giroide (G) y la Primitiva de Schwarz (P), ambas con simetría BCC.

El interés por este tipo de superficies matemáticas viene de dos direcciones diferentes pues por un lado se ha descubierto la existencia de las SMTP en la naturaleza (Han & Che, 2018) y por otro lado se ha publicado que la estructura del hueso trabecular presenta una curvatura promedio casi nula, que es una propiedad distintiva de las SMTP (Bidan, Wang, & Dunlop, 2013; Bobbert et al., 2017). Por lo tanto, el diseño de andamios porosos con morfología interna de SMTP es una opción viable si además se cuenta con la posibilidad de construir dichas estructuras mediante fabricación aditiva o impresión 3D utilizando ácido poliláctico (PLA) que es un poliéster alifático sintético biodegradable y biocompatible, muy utilizado en aplicaciones biomédicas (Farah, Anderson, & Langer, 2016). De forma general, las SMTP exhiben las mejores propiedades requeridas en aplicaciones biomédicas. Por ejemplo las estructuras P y G muestran los mayores valores de permeabilidad a fluido y distribución de los esfuerzos bajo carga axial entre las estructuras con simetría cúbica para los mismos valores de tamaño de poros y porosidad (Blanquer et al., 2017; Jung & Torquato, 2005; Melchels et al., 2010). Por otro lado hay estudios experimentales y computacionales que concluyen que la estructura P es la que exhibe las mejores propiedades mecánicas de

rigidez y resistencia (Maskery et al., 2017; Montazerian, Davoodi, Asadi-Eydivand, Kadkhodapour, & Solati-Hashjin, 2017; Restrepo, Ocampo, Ramírez, Paucar, & García, 2017).

En la fase de diseño del andamio se pueden combinar dos o más de estas estructuras logrando una transición más o menos continua de la morfología (porosidad, tamaño de poros y distribución de estos) en la nueva estructura híbrida. La hibridación puede lograrse mediante la combinación de las ecuaciones implicitas de cada superficie a partir de un código en Mathematica (Yang, Quan, Zhang, & Tian, 2014; Yang & Zhou, 2014) utilizando la función sigmoidea, lo que permite variar el tamaño de poros y la porosidad mediante la variación de las constantes de diseño en cada ecuación para cada estructura SMTP individualmente. Por lo tanto es de interés investigar cómo influyen los parámetros de diseño de cada estructura por separado en las propiedades morfológicas de la estructura híbrida.

El objetivo de este trabajo es estudiar la relación funcional entre las constantes de diseño de las estructuras P y G y las propiedades morfológicas de ellas (porosidad y tamaño de poros) y de estructuras híbridas prismáticas y cilíndricas obtenidas por diseño computacional mediante la función sigmoidea. El artículo se divide en tres secciones: Introducción, Materiales y Métodos, donde se describen el equipamiento y los procedimientos empleados en esta investigación y Análisis de Resultados donde se presenta el análisis de regresión multilineal con la discusión de la pertinencia de los modelos bilineales para caracterizar las propiedades morfológicas de las estructuras híbridas. Finalmente se muestran las conclusiones del trabajo y se trazan las líneas de las futuras investigaciones en la temática.

Marco Teórico

Ecuaciones

En la Tabla 1 se presentan los parámetros de las estructuras P y G donde 𝑋 = 𝑛𝑥𝜋𝑥, 𝑌 = 𝑛𝑦𝜋𝑦, 𝑍 = 𝑛𝑧𝜋𝑧. Los valores de (𝑛𝑥, 𝑛𝑦, 𝑛𝑧) permiten manipular el tamaño de la celda unitaria y por lo tanto el tamaño de poros en cada dirección de los ejes coordenados, mientras que el parámetro C (CP y CG) controla la porosidad.

Tabla 1.

Datos de las SMTP Giroide (G) y P de Schwarz.

De manera general, valores grandes de (𝑛𝑥, 𝑛𝑦, 𝑛𝑧) resultan en menores tamaños de poros sin afectar directamente la porosidad, mientras que la porosidad se incrementa linealmente con C, siendo este el único parámetro que la afecta. En estas ecuaciones 𝜙ℎ𝑖𝑏 ≥ 0 representa el sólido mientras que 𝜙ℎ𝑖𝑏 ≤ 0 representa el espacio ocupado por los poros. La estructura híbrida se define como:

Función de enlace

En la Figura 1 se presenta la función de enlace sigmoidea 𝛼(𝑘, 𝑧) ∈ [0,1] para varios valores del parámetro k (parámetro de enlace), donde notamos en línea de rayas y puntos el valor k=0.5 que es recomendado en la literatura (Yang et al., 2014). Esta función actúa como puente para conectar de manera suave y continua a dos morfologías diferentes como la P y la G. El método de la función sigmoidea constituye una manera eficiente de conectar dos estructuras en una híbrida e incluso puede utilizarse de forma recursiva para conectar varias estructuras y de esta forma cambiar la morfología final de poros (Yang et al., 2014). La función 𝐵(𝑥, 𝑦, 𝑧) = 0 define una frontera entre dos estructuras que, de acuerdo con la literatura, puede ser un cilindro, esfera o plano (Yang et al., 2014). Para el diseño de las estructuras prismáticas se utiliza 𝐵(𝑥, 𝑦, 𝑧) = 𝑧 − 𝑧0 que define un plano perpendicular al eje Z (caso de la figura), mientras que para las probetas cilíndricas 𝐵(𝑥, 𝑦, 𝑧) = 𝑟2 − (𝑥2 + 𝑦2) que define un cilindro con eje longitudinal en Z. Esto permite la unión de dos estructuras prismáticas a lo largo del eje Z (mayor longitud) y la unión de dos estructuras distribuidas radialmente en una estructura cilíndrica.

Figura 1. Función de enlace sigmoidea

Metodología

Esta es una investigación teórica que busca determinar la relación funcional entre los factores de diseño CP y CG con propiedades morfológicas de las estructuras como la porosidad y el tamaño de poros. En el siguiente epígrafe se describe el diseño de experimentos empleado en esta investigación y los procedimientos y métodos que se emplearon para la determinación de las propiedades morfológicas mencionadas.

Diseño de experimentos

En el diseño de las estructuras G y P se toma al parámetro C como factor en un esquema de 5 niveles para determinar las ecuaciones de porosidad y tamaño de poros. En el caso de las estructuras híbridas se propone un diseño factorial 32 tomando como factores la constante C de cada SMTP que se denominaron como CG y CP. Las variables de respuesta son la porosidad (Por) medida en porciento y el tamaño de poros (TP) en micrómetros. El diseño elegido responde a la hipótesis de que no existe interacción entre las variables independientes pues cada una de ellas está relacionada con una SMTP diferente. Los niveles de análisis de las variables se presentan en la Tabla 2 para cada estructura por separado y en la Tabla 3 para las estructuras híbridas. El procesamiento estadístico para determinar la relación entre porosidad y tamaño de poros con los factores CP y CG se realizó con Statgraphics 18-X64 y utilizando las funciones de procesamiento estadístico de Mathematica.

Tabla 2.

Niveles de factores para P y G por separado.

Tabla 3.

Niveles de factores para estructuras híbridas GP prismáticas y cilíndricas.

Normas del diseño computacional

Las probetas fueron diseñadas utilizando el software CAS Wolfram Mathematica v11.2 a partir de un código propio y sus dimensiones son de 12.7 × 12.7 × 25.4 mm de acuerdo con la norma ASTM D695_15 (ASTM-International, 2015) para las probetas prismáticas, mientras que las probetas cilíndricas se diseñaron de 15mm de radio y 10mm de altura. El objeto 3D generado en el software es luego exportado a formato STL y abierto en AutoDesk Meshmixer v3.4.474 (Autodesk, Inc 2017) para ser analizado y corregido debido a posibles errores de mallado. El STL final es obtenido una vez que pasa la prueba de corrección en Meshmixer.

Para lograr un tamaño de poros adecuado en el rango del hueso trabecular (100-900 µm) la celda unitaria debe medir entre 2 y 3 mm. Además, el diseño debe tener en cuenta la periodicidad espacial de estas estructuras ya que una cantidad de celdas par o impar puede afectar la calidad de la impresión, por lo que elegimos una cantidad de celdas de 6 × 6 × 12 para las probetas prismáticas lo que nos da un tamaño de celda unitaria de aproximadamente 2.12 mm de acuerdo con:

En el caso de las probetas cilíndricas se utilizó un tamaño de celdas de 2.3 mm. La calidad del STL generado a partir del diseño se controla en el software con la opción 𝑃𝑙𝑜𝑡𝑃𝑜𝑖𝑛𝑡𝑠 → 120 de la función RegionPlot3D que es la que genera la superficie implícita, aunque no es la única función en Mathematica capaz de generar las estructuras en 3D.

Caracterización microestructural

La caracterización microestructural comprende la determinación de la porosidad y el tamaño de poros en la fase de diseño. La porosidad se determina analizando el STL de la probeta con la herramienta AnalysisStability de Meshmixer que permite determinar el volumen de sólido y el área superficial. El volumen total teórico será 𝑙𝑎𝑟𝑔𝑜 × 𝑎𝑛𝑐ℎ𝑜 × 𝑎𝑙𝑡𝑢𝑟𝑎 para las probetas prismáticas y 𝜋𝑟2ℎ para las probetas cilíndricas, por lo que la porosidad teórica será:

Para determinar el tamaño de poros, se realiza una foto mediante captura de pantalla del STL abierto en Meshmixer y la imagen en JPG es pre-procesada en Microsoft Paint para obtener la región de interés. Esta imagen es luego procesada por un código propio creado en Mathematica que utiliza un algoritmo de segmentación por inundación (Watershed) para seleccionar las áreas de poros y separarlas del fondo de la imagen. Mediante la función ComponentMeasurements de Mathematica es posible realizar una serie de cálculos sobre los objetos segmentados que permiten determinar los estimados del tamaño de poros

Resultados

En la Tabla 4 aparecen los resultados de porosidad y tamaño de poros para cada una de las estructuras P y G para los 5 niveles. Con estos datos, los resultados de regresión para la porosidad y el tamaño de poros son los de la Tabla 5.

Tabla 4.

Resultados del experimento para las estructuras P y G por separado.

Tabla 5.

Resultados de regresión para porosidad y tamaño de poros para P y G.

En la Tabla 6 se exponen los resultados de los experimentos de las estructuras híbridas GP prismáticas y cilíndricas. El diseño factorial utilizado permite disponer de mayor cantidad de datos para estimar los efectos de interacción entre los factores y de esta forma establecer el mejor modelo que correlaciona las variables de respuesta de porosidad y tamaño de poros con los factores CP y CG, dado que en esta fase de diseño no se realizan réplicas experimentales.

Tabla 6.

Resultados del experimento para las estructuras GP prismáticas y cilíndricas.

Por otra parte, en la Tabla 7 se observan los resultados más importantes de la regresión multilineal para la porosidad y el tamaño de poros en las probetas prismáticas y cilíndricas sin tener en cuenta los efectos de interacción.

Tabla 7.

Resultados de regresión para las estructuras híbridas prismáticas y cilíndricas

Discusión de resultados

Resultados de la segmentación con Mathematica

El software Mathematica nos permite diseñar las estructuras mediante la función RegionPlot3D y realizar el análisis de segmentación de las fotos mediante un algoritmo de inundación (Watershed) para determinar el tamaño de poros. En la Figura 2 (a y b) se muestran una probeta prismática diseñada en Mathematica junto a la foto frontal procesada por el algoritmo de segmentación que señala en color rojo los poros de la estructura. En la Figura 2 (c y d) se muestra el mismo procesamiento para las probetas cilíndricas.

Figura 2. (a) y (b): Probeta prismática del experimento 3 (CP=0.6 y CG=-0.7) con su correspondiente foto procesada para tamaño de poros. (c) y (d) Idem para la probeta cilíndrica del experimento 3.

Como puede notarse del detalle en la Figura 2, el algoritmo desecha aquellas partes abiertas que están en contacto con el borde. Aunque la foto aparece en colores grises, en realidad la segmentación se lleva a cabo sobre una versión binaria (negro y blanco) de cada foto para lograr el contraste adecuado para el análisis, tomando el fondo más oscuro de la figura como poros abiertos. En realidad lo que se obtiene mediante este método es un estimado de la porosidad o podríamos decir que es una porosidad superficial. No obstante como estas estructuras son totalmente conectadas, no existe porosidad oculta ni oclusiones. Además, aunque los poros en la parte frontal no tienen las mismas dimensiones que en la profundidad de la probeta, esto puede tenerse en cuenta al evaluar esta porosidad superficial en otras caras de la probeta. Un análisis más detallado necesitaría de cortes en la estructura 3D, lo cual puede ser tema de trabajos futuros.

Análisis de regresión de las estructuras P y G

En la Tabla 5 puede observarse que para cada una de las estructuras, en cada variable de respuesta, los p- valores de los coeficientes de regresión son mucho menores que el nivel de significación del 5 %, lo que indica que en todos los casos estos valores son significativamente diferentes de cero para el nivel de confianza de 95 % y no pueden ser eliminados del modelo. Los valores de R2 indican que el modelo lineal escogido explica más del 99 % de la variabilidad en la porosidad y el tamaño de poros lo que se confirma además por los valores de ��2 que son muy próximos al 99 % e indican una buena confiabilidad en el modelo lineal dando una concordancia casi perfecta entre los valores observados o medidos con los predichos por las ecuaciones de regresión lineal.

Otros resultados de test estadísticos realizados sobre estos datos permiten concluir que el modelo de regresión lineal escogido es apropiado para la descripción de la dependencia funcional entre porosidad y tamaño de poros con el parámetro de diseño C. Para la estructura P, el test de Durbin-Watson (González, Liste, & Felpeto, 2013a) da valores de 2.31737 y 1.69889 para la porosidad y tamaño de poros,

respectivamente. En la estructura G los resultados son de 1.44619 y 2.07531 para porosidad y tamaño de poros, respectivamente. Estos valores permiten concluir que no existe autocorrelación entre los residuos de regresión para cada variable de respuesta al nivel de significación del 95%. Por otro lado, los test de Breusch-Pagan (González, Liste, & Felpeto, 2013b) realizados en cada variable y para cada estructura revelan que existe homocedasticidad en los residuos de regresión (varianza constante).

Para comprobar la distribución normal de los residuos de regresión se realizaron varias pruebas estadísticas como las de Kolmogorov-Smirnov (Corder & Foreman, 2009), Anderson-Darling (T. W. Anderson & Darling, 1952; D´Agostino & Stephens, 1986; Stephens, 1974) y Shapiro-Wilk (Shapiro, Wilk, & Chen, 1968) que realizan comparaciones entre distribuciones de probabilidad tomando como referencia la distribución normal. La prueba de Kolmogorov-Smirnov tiene mayor sensibilidad para valores próximos a la mediana de la distribución, mientras que la prueba de Anderson-Darling es más sensible a valores extremos cercanos a las colas de las distribuciones y la prueba de Shapiro-Wilk es muy definitiva para muestras pequeñas (n<50). En todos los casos, los p-valores de cada uno de estos estimadores son mayores que el nivel de significación del 5 % indicando que los residuos de regresión tienen distribución normal, lo cual es un indicativo de la pertinencia del modelo lineal en cada caso.

Análisis de regresión de las estructuras híbridas GP

En la Figura 3 podemos observar los Gráficos de Pareto Estandarizados para la porosidad y el tamaño de poros de las estructuras híbridas prismáticas y cilíndricas. En las probetas cilíndricas los p-valores del efecto de interacción (AB) entre los factores son de 0.8422 y 0.2890 para la porosidad y el tamaño de poros respectivamente, por lo que al ser mayores que 0.05 indican que estos efectos pueden despreciarse y el modelo puede simplificarse al eliminar el término de interacción. En el caso de las probetas prismáticas esto también puede realizarse para la porosidad (p=0.137204) pero en el caso del tamaño de poros p=0.02288<0.05 por lo no debería ser eliminado del modelo. Sin embargo este valor es bastante cercano a 0.05 lo que implica que el efecto de interacción es muy débil y podría eliminarse.

Figura 3. Gráficos de Pareto Estandarizados: (a) y (b) para la porosidad y tamaño de poros de las probetas prismáticas, (c) y (d) para las probetas cilíndricas.

Podemos apreciar que aún cuando se eliminan los efectos de interacción en el modelo del tamaño de poros para las probetas prismáticas (caso b de la Figura 3) se tiene que 𝑅2 = 0.984305 lo cual es aún muy buen indicador de linealidad en el modelo simplificado. En todos los casos, el estadístico de Durbin-Watson se obtiene con p-valores muy superiores a 0.05 lo que indica que no existe autocorrelación serial entre los residuos de regresión al 95 % del nivel de confianza, lo que a su vez indica la pertinencia del modelo elegido.

Conclusiones

En el diseño de estructuras híbridas de SMTP mediante software CAS a partir de sus ecuaciones implícitas es necesario tener en cuenta la periodicidad de las SMTP por separado y controlar el tamaño de la celda unitaria de acuerdo a la precisión permitida por la impresora 3D. El software Mathematica es una opción viable para realizar este tipo de diseños dada la precisión y calidad del STL obtenido y la disponibilidad de un arsenal de funciones matemáticas para lograr los resultados.

El análisis de los datos de regresión muestra que los efectos de interacción entre los factores CG y CP, en cada modelo para las variables de respuesta de Porosidad y Tamaño de Poros, pueden ser despreciados y que en ambos casos se obtiene un modelo lineal generalizado que permite establecer una relación funcional entre los factores de diseño de cada estructura SMTP por separado y las variables de respuesta. Esto implica que los factores CG y CP que controlan la porosidad de manera lineal en cada SMTP por separado (G y P) también controlan la porosidad de forma lineal en el modelo híbrido. El tamaño de poros se controla también por los factores C y por el tamaño de la celda unitaria que en este experimento se ha mantenido igual para ambas SMTP y de esta forma se obtiene también un control lineal del tamaño de poros a partir de los factores CG y CP. Los modelos obtenidos son estadísticamente significativos hasta un 95 % de nivel de confianza.