Colon model
A modellek kidolgozásához és a kapcsolódó szimulációk elvégzéséhez az Abaqus végeselemes szoftver (6.13-as verzió, SIMULIA, Providence, RI) került felhasználásra. Minden szimulációt egy két Intel Xeon processzorral és 32 GByte memóriával rendelkező Dell T7810 munkaállomáson végeztünk. A sertés szöveti viselkedését vizsgáló korábbi tanulmányunkban a sertés vastagbél ereszkedő és spirális régióira egy anizotróp anyagkonstituáló modellt állítottunk fel, amely egyidejű felfúvódási és kitágulási vizsgálatokon alapult15. Ez a modell több korábbi tanulmányban is bizonyította, hogy megfelelően megragadja a vastagbél anizotróp viselkedését19,34,35,36. A korábbi tanulmányunkban meghatározott anyagparamétereket megtartottuk, mivel tudomásunk szerint ez az egyetlen tanulmány, amely a vastagbél két fő deformációs módját, a felfúvódást és a kitágulást modellezi. Ezenkívül megállapították, hogy a sertés vastagbele releváns a divertikulózis tanulmányozásához33,37. A kiválasztott anyagmodell és paraméterek tehát megfelelő választást jelentenek ehhez a tanulmányhoz. Mivel a spirális vastagbél anatómiája nem releváns az emberre nézve, itt a vastagbél szövetanyagának jellemzésére a leszálló vastagbél mintáiból származó eredményeket használtuk. Korábbi munkánkban arról számoltunk be, hogy a sertés ereszkedő vastagbélszakasz hosszában az átmérő lokálisan elhanyagolható mértékben változik, ezért a normál vastagbelet hengeres szakaszként modelleztük (1a. ábra). A belső és külső sugarat \({R}_{i}=11,5\) mm-re, illetve \({R}_{o}=12,7\) mm-re állítottuk be, a korábbi munkánk során mért átlagos értékek alapján. Az axiális hosszúságot \(L=10\) cm-re állítottuk be úgy, hogy az többszörösen nagyobb legyen, mint a sugár. A vastagbélszövetet homogén, anizotróp és inkompresszibilis hiperelasztikus anyagnak feltételeztük, amelynek mechanikai viselkedését a következő anizotróp alakváltozási energiafüggvény \(\bar{W}:\)15,35
Az első tag egy neo-Hookean-választ, amely a vastagbélszövet nem kollagén tartalmú összetevőinek viselkedését jellemzi. A \({I}_{1}={\lambda }_{z}^{2}+{\lambda }_{\theta }^{2}+\frac{1}{{{\lambda }_{z}^{2}{\lambda }_{\theta }^{2}}\) mennyiség az első deformációs tenzor invariánsát, ahol \({\lambda }_{z}\) és \({\lambda }_{\theta }\) a tengelyirányú és a kerületi nyúlásra utal, illetve. A következő kifejezések a kollagénszálak hozzájárulását jellemzik, az \({k}_{1}\) és \({k}_{2}\) a szálak merevségét jelöli. Konkrétan, a második kifejezés (feliratos \(l\)) a hosszanti (axiális) irányba igazított kollagénszálak viselkedését jellemzi, míg a harmadik kifejezés (feliratos \(s\)) a kerületi irányhoz képest két preferenciális szimmetrikus irányban eloszló kollagénszálak válaszát jellemzi. Az \({I}_{4}\) mennyiséget pszeudoinvarianciának nevezzük, és a szálak mechanikai válaszát jellemzi a preferenciális irányok mentén. A második kifejezés \({I}_{4}^{l}={\lambda }_{z}^{2}\), mivel ez a hosszirányú szálakat veszi figyelembe. A harmadik termet \({I}_{4}^{s}\) a következőképpen fejezzük ki:
Itt, \(\pm {\gamma }^{s}\) a két szimmetrikus szálorientáció szögét jelöli a kerületi irányhoz képest. Ebben a munkában öt, korábban öt sertés leszálló vastagbélmintára levezetett anyagparaméter-készletet használtunk az \({C}_{10},\,{k}_{1}^{l},\,{k}_{2}^{l},\,{k}_{1}^{s},\,{k}_{2}^{s},\) és \({\gamma }^{s}\) anyagparaméterekhez. Ez az alakváltozási energiafüggvény közvetlenül nem áll rendelkezésre az Abaqusban. Ezért először egy Fortran alapú felhasználói alprogram, az UANISOHYPER_INV implementálásával integráltuk, amely lehetővé teszi az anizotróp hiperelasztikus anyagviselkedés meghatározását az invariáns formulával. Ehhez az alprogramhoz meg kellett határozni az alakváltozás-energia függvény (1. egyenlet) első és második deriváltját az \({I}_{1}\) és \({I}_{4}\) skaláris invariánsok tekintetében. Annak ellenőrzésére, hogy az alprogramot megfelelően implementálták-e, és a vastagbélszövet számítási modelljének validálására összehasonlítottuk a modell által becsült izokorikus kerületi feszültség \({\bar{\sigma }}_{\theta }^{m}\) és izokorikus axiális feszültség \({\bar{\sigma }}_{z}^{m}\) modell által becsült feszültségértékeit a vastagbél külső felületén a nyúlási energiafüggvény alapján analitikusan kapott és az alábbiak szerint számított értékekkel:
A számítási és az analitikus értékek közötti eltérés számszerűsítésére az \({R}^{2}\) meghatározási együtthatót a következőképpen számoltuk ki:
hol \(A={\bar{\sigma }}_{z},\,{\bar{\sigma }}_{\theta }\), és a “\(avg\)” index az analitikus értékek átlagát jelöli az összes \(n\) adatponton. Az \({R}^{2}\) 1 közeli értéke azt jelzi, hogy globálisan jó korrelációt kapunk.
Divertikulum modell
A beteg vastagbél modellek esetében a divertikulumszerű struktúrát egy \(D\) átmérőjű gömb és a normál vastagbél hengeres alakjának metszéspontjaként modelleztük a felszín feletti \(H\) magasságban. Ebben a tanulmányban csak egy tasakot vettünk figyelembe, hogy egyetlen tasak hatására összpontosítsunk, a több tasakból adódó további variabilitások nélkül (a tasakok száma, elhelyezkedése a vastagbél mentén, egymáshoz viszonyított relatív helyzetük, az egyedi méret változékonysága stb.) Tudomásunk szerint nincs olyan tanulmány, amely szisztematikusan beszámolna a divertikulák geometriájáról. Gyakran említik, hogy átmérőjük jellemzően 5-10 mm között van, anélkül, hogy megfelelően jeleznék, hogy hol mérik az átmérőt38,39. A különböző tasakméretek lefedése érdekében a szimulációkban a kezdeti tasakmagasság \(H\) és a kezdeti tasakátmérő \(D\) három értékét vettük figyelembe: \(H=2\) mm, 4 mm és 6 mm, és \(D\) = 8 mm, 10 mm és 12 mm. A D és H ezen értékeinek kombinációi 9 beteg modellt eredményeztek különböző tasakméretekkel. Egy \({r}_{f}=2\) mm sugarú kör alakú filet került a normál vastagbelet ábrázoló henger és a tasakot ábrázoló gömb metszéspontjába, hogy eltávolítsuk az éles éleket, amelyek nem feltétlenül fordulnak elő a tényleges szövetben, és megváltoztathatják a feszültségértékeket. A divertikulumszövet vastagságát a vastagbélével azonosnak feltételeztük (azaz \({R}_{o}-{R}_{i}\)). Történtek kísérletek a vastagbél mechanikájának jellemzésére divertikulózisban23,30,40,41. Kifejezetten a divertikulumszövetre (azaz az izolált divertikulumzacskó jellemzésére) azonban nem állapítottak meg anyagjellemzőket. Ezért a normál szövethez hasonló tulajdonságokat vettünk figyelembe a tasakszelvényre vonatkozóan.
Rétegkörülmények
A peremfeltételeket a sertés vastagbelén végzett korábbi vizsgálatunk alapján határoztuk meg15. Megfigyeltük, hogy az 1,5 kPa feletti luminális nyomás (a szövet belső falára ható nyomás a nulla külső nyomás mellett) értékek passzív körülmények között a szövet maradandó károsodásához vezettek. A számítógépes szimulációkban ezt az értéket tovább toltuk, hogy láthatóvá tegyük a feszültségértékek alakulásának tendenciáját, és 5 kPa-ig terjedő luminális nyomást alkalmaztunk a vastagbél belső felületére, beleértve a tasakot is. A külső falnyomást 0-val tartottuk egyenlőnek. Mivel vékony falú szövettel dolgozunk (a belső átmérő és a vastagság aránya közel 20), a luminális nyomás így a belső nyomás (pl. a széklet vagy a bélgáz által a szövet vezetékén belül kifejtett nyomás) és a külső nyomás (pl. a hasi nyomás) közötti deltaként volt értelmezhető. Jellemzően kis nyomás szükséges a vastagbélszövet megnyitásához, amely egyébként összeesne. Ez a nyomás (~0,2 kPa) azonban nagyon kicsi az itt vizsgált tartományhoz képest, ezért figyelmen kívül hagyjuk. A nyomást kvázi statikus lépésekben alkalmaztuk. Korábbi vizsgálataink során jellemzően 10% körüli axiális nyúlást tapasztaltunk in vivo a vastagbélben. Ezért a nyomás alá helyezés előtt a \(Z\) irányban 10%-os axiális nyújtást alkalmaztunk. A sertés vastagbélszövetében nagyon kis nyitási szögértékeket találtunk, ezért a maradó kerületi feszültségeket elhanyagoltuk a modellben. Mind a normál, mind a beteg modellek esetében a \(YZ\) sík és az \(XY\) sík tekintetében szimmetriát alkalmaztunk, és a nagyobb hatékonyság elérése érdekében a modellnek csak egynegyedét oldottuk meg, mielőtt a teljes modellt rekonstruáltuk volna vizualizációs céllal.
Hálózás
A geometria megfelelő metszését minden esetben elvégeztük, hogy a hálózás lineáris hatszögletes elemekkel lehetséges legyen. Hibrid elemformulációt választottunk az anyag összenyomhatatlanságának kikényszerítésére (azaz C3D8H elemet az Abaqus könyvtárból). A normál modell esetében a háló konvergencia vizsgálatát úgy végeztük el, hogy a referenciaháló méretét az Abaqus által automatikusan javasolt alapértelmezett méretről (1,3 mm) kiindulva addig csökkentettük, amíg a legnagyobb főfeszültségek maximális értékének változása két egymást követő háló között 10%-nál kisebb nem lett. Ezt követően a konvergáló hálóméretet használtuk a beteg modellekben is, néhány további finomítással a tasak körül. A normál modellhez összesen 69608 elemet és 88350 csomópontot használtak az elemzéshez megtartott legfinomabb hálóméretnél (a geometria negyede). A kilenc beteg modellhez 69764 elemet és 80656 elemet (88605, illetve 102440 csomópontot) használtak (a geometria negyedét). A helyi anyagorientációt az Abaqusban elérhető diszkrét formulázási opció segítségével adtuk meg, hogy megfelelően meghatározzuk a szálorientációt minden egyes elemben, ahogyan az a nyúlási energiafüggvényben szükséges.
Eredmények elemzése
Ez összesen 10 modellt (egy normális és kilenc beteg, különböző méretű tasakkal) vettünk figyelembe. Minden modellt 5 különböző anyagparaméter-készlettel szimuláltunk, ami összesen 50 szimulált esetet eredményezett. A modellek közötti összehasonlításhoz szükség van a feszültségszint mérésére. Mivel nem egy tönkremeneteli kritériumot határozunk meg, hanem egyszerűen csak a különböző körülmények közötti változásokat vizsgáljuk, bármilyen mérőszám, például a maximális főfeszültség vagy a von Mises-feszültség megfelelő lehet. Az eredményeket itt csak a maximális főfeszültség tekintetében mutatjuk be a redundáns információk elkerülése érdekében, mivel a von Mises-feszültséggel végzett általános megfigyelések hasonlónak bizonyultak. A von Mises-feszültséggel végzett elemzések és eredmények az e tanulmányhoz kapcsolódó adatok között találhatók, ha az olvasó szeretné azokat megtekinteni. Minden egyes szimulált esetben a legnagyobb főfeszültséget (a hálóelem középpontjában), amelyre \({\sigma }_{MPS}\) néven hivatkozunk, 0 és 5 kPa közötti nyomástartományban figyeltük meg, 0,25 kPa lépéssel. Egy Python szkriptet fejlesztettünk ki, amely automatikusan megkeresi az \({\sigma }_{MPS}\) maximális értékét, amelyet \({\sigma }_{MPS,max}\) néven jelölünk a különböző nyomásoknál minden egyes szimulált esetre. Minden modell esetében az egyes szimulált esetekre (azaz mind az öt anyagparaméter-készletre) kapott átlagértéket kiszámítottuk és megtartottuk a későbbi elemzéshez. Ezeket az értékeket \({\sigma }_{MPS,\,max}^{avg}\) és \({\sigma }_{MPS,\,max}^{avg,\,n}\) néven nevezzük a beteg és normál modellek esetében.
Az \({\sigma }_{MPS,\,max}^{avg}\) értékek 1 és 5 kPa közötti nyomás esetén (0,25 kPa lépésekben) korreláltak a tasak geometriájának különböző paramétereivel (ábra. 1c) a tasakgeometria és a magas feszültség közötti kapcsolat azonosítása érdekében: a tasaknyak szélessége a \(x\) tengely mentén \({D}_{x}\), a tasaknyak szélessége a \(z\) tengely mentén \({D}_{z}\), a tasaknyak alján lévő terület \({A}_{b}\) és a tasak lumen felőli felülete \({S}_{p}\). A Pearson-féle korrelációs együtthatót \({r}^{2}\) minden paraméterre kiszámítottuk. A korreláció szignifikanciájának értékelésére a kétcsapos p-értéket is kiszámítottuk.
Mellett kiszámítottuk az \({P}_{d}\) nyomást a beteg modellekben úgy, hogy egy adott \({P}_{n}\) nyomáshoz a normális vastagbélben a következőkkel rendelkezünk:
A következő nyomástényezőt ezt követően a következőképpen számítottuk:
Ez a tényező azt a %-ot jelzi, amellyel a nyomást \({P}_{n}\)csökkenteni kell a beteg vastagbélben, hogy \({\sigma }_{MPS,max}\) hasonló legyen \({\sigma }_{MPS,\,max}^{n}\), i.azaz a maximális főfeszültség maximális értéke a beteg modellben nem haladja meg a normál modell megfelelő értékét. A csökkentett nyomást ezután \({P}_{d}\) jelöli. Konkrétan, minden egyes szimulált esetre az \({P}_{n}\) értékét 0,25 kPa értékkel csökkentettük, amíg a 6. egyenletből származó összefüggés nem teljesült, megadva az \({P}_{d}\) megfelelő értékét, és kiszámítottuk az \(f\) nyomástényezőt. Ezt \({P}_{n}\) esetén a 2 és 5 kPa közötti tartományban (1 kPa lépésekben) végeztük el. A modellekre kapott átlagértéket az összes szimulált esetre vonatkozóan kiszámítottuk és megtartottuk a későbbi elemzéshez.
Végül meghatároztuk a tasak befolyási zónáját, azaz a tasak körüli távolságot, ahol a feszültség megemelkedik. Ennek érdekében minden egyes beteg esetre létrehoztuk a vastagbél lumen oldali felületén a maximális főfeszültség térképét (mivel a feszültség ezen a felületen a legnagyobb) a tasak közepén áthaladó hosszanti és kerületi útvonalakon elhelyezkedő hálócsomópontok mentén. A maximális főfeszültség alakulását ezen útvonalak mentén az \({d}_{l}\) és \({d}_{c}\) távolságok függvényében figyeltük meg, az erszény középpontjához képest mérve. Ezután kiszámítottuk azokat a legtávolabbi \({d}_{l}^{inf}\) távolságokat a hosszirányú útvonal mentén és \({d}_{c}^{inf}\) távolságokat a kerületi útvonal mentén, amelyek után a maximális főfeszültség a normál modell megfelelő értékének 10%-án belülre csökken ugyanezen útvonalak mentén. Ezek a távolságok képezték a befolyási zóna mértékét. Megfigyelve, hogy egy bizonyos nyomás után elérik a plató értékét, az összes esetre vonatkozó átlagértéküket minden modell esetében 5 kPa mellett az \({\sigma }_{MPS,\,max}^{avg}\) \(\sigma }_{MPS,\,max}^{avg}\) \(\sigma }_{MPS,\,max}^{avg}\) \(\sigma })
A számításokat Python programozási nyelven, Jupyter notebooks42 fejlesztésével végeztük. Az adatokat a Python Pandas könyvtárával43 tároltuk és kezeltük. Az alapvető matematikai műveleteket a NumPy könyvtárral44 végeztük. A korrelációs együtthatót \({r}^{2}\) és a kétfarkú p-értéket a Scipy könyvtár “stats” funkciójával számoltuk45. Az összes ábrát a Matplotlib könyvtárral46 valósítottuk meg.