Paksusuolen malli
Mallien kehittämisessä ja niihin liittyvien simulaatioiden suorittamisessa käytettiin Abaqus-finite elementtiohjelmistoa (versio 6.13, SIMULIA, Providence, RI). Kaikki simuloinnit suoritettiin Dell Workstation T7810 -tietokoneella, jossa oli kaksi Intel Xeon -prosessoria ja 32 gigatavua muistia. Aiemmassa sian kudoskäyttäytymistä koskevassa tutkimuksessamme laadimme anisotrooppisen materiaalin konstitutiivisen mallin sian paksusuolen laskeville ja spiraalialueille, joka perustui samanaikaisiin puhallus- ja venytystesteihin15. Tämän mallin on useissa aiemmissa tutkimuksissa osoitettu kuvaavan asianmukaisesti paksusuolen anisotrooppista käyttäytymistä19,34,35,36. Aiemmassa tutkimuksessamme määritetyt materiaaliparametrit säilytettiin, koska se on tietojemme mukaan ainoa tutkimus, jossa mallinnettiin suurten eläinten paksusuolta inflaatiossa ja ekstensiossa, jotka ovat paksusuolen kaksi tärkeintä muodonmuutosmuotoa. Lisäksi on todettu, että sian paksusuoli soveltuu divertikuloositutkimukseen33,37. Valittu materiaalimalli ja parametrit ovat siten sopivia valintoja tähän tutkimukseen. Koska spiraalipaksusuolen anatomia ei ole merkityksellinen ihmiselle, paksusuolen kudosmateriaalin karakterisoinnissa käytettiin tässä laskevan paksusuolen näytteistä saatuja tuloksia. Olemme raportoineet aiemmassa työssämme, että sikojen laskevan paksusuolen segmentin halkaisijan paikallinen vaihtelu on vähäistä, joten normaali paksusuoli mallinnettiin sylinterimäisenä segmenttinä (kuva 1a). Sisä- ja ulkosäde asetettiin \({R}_{i}=11.5\) mm:ksi ja \({R}_{o}=12.7\) mm:ksi aiemmassa työssämme mitattujen keskiarvojen perusteella. Aksiaaliseksi pituudeksi asetettiin \(L=10\) cm siten, että se on useita kertoja suurempi kuin säde. Paksusuolen kudoksen oletettiin olevan homogeeninen, anisotrooppinen ja kokoonpuristumaton hyperelastinen materiaali, jonka mekaanista käyttäytymistä kuvaa seuraava anisotrooppinen rasitusenergian funktio \(\bar{W}:\)15,35
Ensimmäinen termi kuvaa neo-Hookeanin vastetta, joka luonnehtii paksusuolen kudoksen ei-kollageenisten ainesosien käyttäytymistä. Suure \({I}_{1}={\lambda }_{z}^{2}+{\lambda }_{\theta }^{2}+\frac{1}{{{\lambda }_{z}^{2}{\lambda }_{\theta }^{2}}}\) edustaa ensimmäistä muodonmuutostensorin invariantti, jossa \({\lambda }_{z}\) ja \({\lambda }_{\theta }\) viittaavat aksiaaliseen ja kehän suuntaiseen venytykseen, vastaavasti. Seuraavat termit kuvaavat kollageenikuitujen osuutta, ja \({k}_{1}\) ja \({k}_{2}\) kuvaavat kuitujen jäykkyyttä. Erityisesti toinen termi (jonka yläindeksi on \(l\)) kuvaa pitkittäissuuntaan (aksiaaliseen suuntaan) suuntautuneiden kollageenikuitujen käyttäytymistä, kun taas kolmas termi (jonka yläindeksi on \(s\)) kuvaa sellaisten kollageenikuitujen vastetta, jotka on hajautettu kahteen symmetriseen etusäteissuuntaan kehäsuuntaan nähden. Määrää \({I}_{4}\) kutsutaan pseudoinvariantiksi, ja se kuvaa kuitujen mekaanista vastetta preferenssisuunnissa. Toisen termin \({I}_{4}^{l}={\lambda }_{z}^{2}\) osalta, koska se ottaa huomioon kuidut pituussuunnassa. Kolmas termi \({I}_{4}^{s}\) ilmaistaan seuraavasti:
Tässä, \(\pm {\gamma }^{s}\) tarkoittaa kahden symmetrisen kuidun suuntauskulman kulmaa kehäsuuntaan nähden. Materiaaliparametrien \({C}_{10},\,{k}_{1}^{l},\,{k}_{2}^{l},\,{k}_{1}^{s},\,{k}_{2}^{s},\) ja \({\gamma }^{s}\) materiaaliparametreina käytettiin tässä työssä viittä materiaaliparametrijoukkoa, jotka oli aiemmin johdettu viidelle sian laskevasta paksusuolesta otetulle näytteelle. Tämä rasitusenergiafunktio ei ole suoraan käytettävissä Abaqusissa. Siksi se integroitiin ensimmäistä kertaa toteuttamalla Fortran-pohjainen käyttäjän aliohjelma nimeltä UANISOHYPER_INV, joka mahdollistaa anisotrooppisen hyperelastisen materiaalikäyttäytymisen määrittelyn invarianttiformulaation avulla. Tämä aliohjelma edellytti muodonmuutos-energiafunktion (yhtälö 1) ensimmäisen ja toisen derivaatan määrittämistä skalaaristen invarianttien \({I}_{1}\) ja \({I}_{4}\) suhteen. Varmistaaksemme, että aliohjelma oli toteutettu oikein, ja validoidaksemme paksusuolen kudoksen laskennallisen mallin, vertasimme mallin estimoimia isokoorisen kehäjännityksen \({\bar{\sigma }}_{\theta }^{m}\) ja isokoorisen aksiaalijännityksen \({\bar{\sigma }}_{z}^{m}\) isokoorisia jännitysarvoja paksusuolen ulommalla pinnalla arvoihin, jotka saatiin analyyttisesti rasitusenergiafunktioon perustuvina ja jotka laskettiin seuraavasti:
Määrityskerroin \({R}^{2}\) laskettiin seuraavasti laskennallisten ja analyyttisten arvojen välisen poikkeaman kvantifioimiseksi:
jossa \(A={\bar{\sigma }}_{z},\,{\bar{\sigma }}_{\theta }\), ja alaviite ”\(avg\)” tarkoittaa analyyttisten arvojen keskiarvoa kaikista \(n\) datapisteistä. \({R}^{2}\) arvo, joka on lähellä 1, osoittaa, että globaalisti saadaan hyvä korrelaatio.
Divertikkelimalli
Sairaan paksusuolen malleissa divertikkelin kaltainen rakenne mallinnettiin halkaisijaltaan \(D\) pallon ja normaalin paksusuolen sylinterinmuotoisen rakenteen leikkauspisteenä korkeudella pinnan yläpuolella olevasta pallosta \(H\). Tässä tutkimuksessa otettiin huomioon vain yksi pussi, jotta voitiin keskittyä yhden pussin vaikutukseen ilman lisävaihteluita, joita useat pussit aiheuttaisivat (pussien lukumäärä, sijainti paksusuolen varrella, suhteellinen sijainti toisiinsa nähden, yksilön koon vaihtelu jne.). Tietojemme mukaan divertikkelien geometriaa ei ole systemaattisesti raportoitu yhdessäkään tutkimuksessa. Yleisesti mainitaan, että niiden läpimitta on tyypillisesti 5-10 millimetriä, mutta ei mainita, mistä läpimitta mitataan38,39. Erilaisten pussikokojen kattamiseksi simulaatioissa otettiin huomioon kolme arvoa pussin alkukorkeudelle \(H\) ja pussin alkuhalkaisijalle \(D\): \(H=2\) mm, 4 mm ja 6 mm sekä \(D\) = 8 mm, 10 mm ja 12 mm. Näiden D- ja H-arvojen yhdistelmillä saatiin 9 sairasta mallia, joilla oli erilainen pussin koko. Normaalia paksusuolta edustavan sylinterin ja pussia edustavan pallon väliseen leikkauspisteeseen sisällytettiin pyöreä silaus, jonka säde oli \({r}_{f}=2\) mm, jotta voidaan poistaa terävät reunat, joita ei välttämättä esiinny todellisessa kudoksessa ja jotka voivat muuttaa jännitysarvoja. Divertikkelikudoksen paksuuden oletettiin olevan sama kuin paksusuolen (eli \({R}_{o}-{R}_{i}\)). Paksusuolen mekaniikkaa on yritetty luonnehtia divertikuloosissa23,30,40,41. Materiaaliominaisuuksia ei kuitenkaan ole määritetty erityisesti divertikkelikudokselle (eli eristetyn divertikkelipussin karakterisointia). Näin ollen pussin poikkileikkauksen ominaisuuksina pidettiin normaalia kudosta vastaavia ominaisuuksia.
Rajaolosuhteet
Rajaolosuhteet asetettiin aikaisemman sian paksusuolta koskevan tutkimuksemme perusteella15. Havaitsimme, että yli 1,5 kPa:n luminaalipaineen (kudoksen sisäseinämään kohdistuva paine, kun ulkoinen paine on nolla) arvot johtivat kudoksen pysyvään vaurioitumiseen passiivisissa olosuhteissa. Laskennallisissa simulaatioissa tätä arvoa pidennettiin edelleen, jotta jännitysarvojen kehityssuuntaus saataisiin näkyviin, ja paksusuolen sisäpinnalle, pussi mukaan luettuna, asetettiin jopa 5 kPa:n luminaalipaine. Koska kyseessä on ohutseinäinen kudos (sisähalkaisijan ja paksuuden suhde lähes 20), luminaalipaine voidaan rinnastaa sisäisen paineen (esim. ulosteen tai suolistokaasun aiheuttama paine kudoksen kanavan sisällä) ja ulkoisen paineen (esim. vatsaontelon paine) väliseksi deltaksi. Tyypillisesti tarvitaan pieni paine avaamaan paksusuolen kudos, joka muuten romahtaisi. Tämä paine (~0,2 kPa) on kuitenkin hyvin pieni verrattuna tässä tarkasteltavaan alueeseen, joten sitä ei oteta huomioon. Painetta käytettiin kvasistaattisissa portaissa. Aksiaalinen noin 10 prosentin venytys havaittiin tyypillisesti paksusuolessa in vivo aiemmassa tutkimuksessamme. Näin ollen 10 prosentin aksiaalinen venytys asetettiin \(Z\)-suunnassa ennen paineistamista. Olemme havainneet hyvin pieniä avautumiskulma-arvoja sian paksusuolen kudoksessa, ja jäännösjännitykset kehän suuntaisesti jätettiin näin ollen mallissa huomiotta. Sekä normaaleissa että sairaissa malleissa käytimme symmetriaa \(YZ\)-tason ja \(XY\)-tason suhteen ja ratkaistiin vain neljäsosa mallista suuremman tehokkuuden saavuttamiseksi, ennen kuin koko malli rekonstruoitiin havainnollistamista varten.
Verkottaminen
Kussakin tapauksessa tehtiin sopiva geometrian leikkaus, jotta varmistuttiin siitä, että verkottaminen on mahdollista lineaaristen kuusiulotteisten elementtien avulla. Materiaalin kokoonpuristumattomuuden pakottamiseksi valittiin hybridielementtimuodostus (eli elementti C3D8H Abaqus-kirjastosta). Verkon konvergenssitutkimus tehtiin normaalimallille pienentämällä referenssiverkon kokoa Abaqusin automaattisesti ehdottamasta oletuskoosta (1,3 mm) alkaen, kunnes suurimpien pääjännitysten enimmäisarvon vaihtelu kahden peräkkäisen verkon välillä oli alle 10 %. Lähestyvää verkkokokoa käytettiin tämän jälkeen myös sairaiden malleissa, ja pussin ympärillä tehtiin joitakin hienosäätöjä. Normaalimallissa käytettiin yhteensä 69608 elementtiä ja 88350 solmua hienoimmalla analyysiä varten säilytetyllä verkkokoolla (neljäsosa geometriasta). Yhdeksässä sairaan mallissa käytettiin 69764 elementtiä ja 80656 elementtiä (88605 ja 102440 solmua) (neljännes geometriasta). Paikallinen materiaalin orientaatio määritettiin käyttämällä Abaqusissa käytettävissä olevaa diskreettia muotoiluvaihtoehtoa, jotta kuitujen orientaatio kussakin elementissä voitiin määritellä oikein, kuten venymäenergiafunktiossa edellytetään.
Tulosten analyysi
Kokonaisuudessaan tarkasteltiin kymmentä mallia (yksi normaali ja yhdeksän sairastunutta, joilla oli eri kokoiset pussit). Kutakin mallia simuloitiin viidellä eri materiaaliparametrijoukolla, jolloin saatiin yhteensä 50 simuloitua tapausta. Eri mallien vertailemiseksi tarvitaan jännitystason mittari. Koska emme ole määrittämässä vikaantumiskriteeriä, vaan tarkastelemme ainoastaan vaihtelua eri olosuhteissa, mikä tahansa mittari, kuten suurin pääjännitys tai von Misesin jännitys, olisi riittävä. Tulokset esitetään tässä vain suurimman pääjännityksen osalta turhan tiedon välttämiseksi, koska von Mises -jännityksellä tehtyjen havaintojen todettiin olevan samankaltaisia. Von Mises -jännityksen analyysi ja tulokset sisältyvät tähän tutkimukseen liittyviin tietoihin, jos lukija haluaa tutustua niihin. Kunkin simuloidun tapauksen osalta suurinta pääjännitystä (verkkoelementin keskipisteessä), jota kutsutaan nimellä \({\sigma }_{MPS}\), tarkkailtiin painealueella 0-5 kPa ja 0,25 kPa:n lisäyksellä. Python-skripti kehitettiin etsimään automaattisesti \({\sigma }_{MPS}\):n maksimiarvo, jota kutsutaan nimellä \({\sigma }_{MPS,max}\) eri paineissa kullekin simuloidulle tapaukselle. Kullekin mallille laskettiin kullekin simuloidulle tapaukselle (eli kaikille viidelle materiaaliparametrijoukolle) saatu keskiarvo, joka säilytettiin myöhempää analyysia varten. Näitä arvoja kutsutaan \({\sigma }_{MPS,\,max}^{avg}\) ja \({\sigma }_{MPS,\,max}^{avg,\,\,n}\) sairaiden ja normaalien mallien osalta.
Arvot \({\sigma }_{MPS,\,max}^{avg}\), korreloivat 1 ja 5 kPa:n välisellä paineella (0,25 kPa:n askelin) pussin geometrian eri parametreihin (Kuva. 1c) pussin geometrian ja suuren rasituksen välisen suhteen tunnistamiseksi: pussin kaulan leveys \(x\)-akselilla \({D}_{x}\), pussin kaulan leveys \(z\)-akselilla \({D}_{z}\), pinta-ala pussin kaulan tyvellä \({A}_{b}\) ja pussin pinta lumenin puolella \({S}_{p}\). Kullekin parametrille laskettiin Pearsonin korrelaatiokerroin \({r}^{2}\). Korrelaation merkitsevyyden arvioimiseksi laskettiin myös kaksihaarainen p-arvo.
Laskimme lisäksi paineen \({P}_{d}\) sairaissa malleissa siten, että tietylle paineelle \({P}_{n}\) normaalissa paksusuolessa on:
Selvitettiin sen jälkeen seuraava painekerroin:
Tämä kerroin ilmaisee sen %:n, jolla painetta \({P}_{n}\)on alennettava sairastuneessa paksusuolessa, jotta \({\sigma }_{MPS,max}\) on samanlainen kuin \({\sigma }_{MPS,\,max}^{n}\), i.e pääjännityksen maksimiarvo sairaassa mallissa ei ylitä vastaavaa arvoa normaalissa mallissa. Tällöin alennettua painetta nimitetään \({P}_{d}\). Kunkin simuloidun tapauksen osalta \({P}_{n}\) pienennettiin 0,25 kPa:lla, kunnes yhtälön 6 suhde täyttyi, jolloin saatiin vastaava arvo \({P}_{d}\) ja laskettiin painekerroin \(f\). Tämä tehtiin \({P}_{n}\) arvoille 2-5 kPa (1 kPa:n askelin). Laskettiin kaikkien simuloitujen tapausten malleille saatu keskiarvo, joka säilytettiin myöhempää analyysia varten.
Viimeiseksi määritettiin pussin vaikutusalue eli etäisyys pussin ympärillä, jossa jännitys on koholla. Tätä varten laadimme kullekin sairastuneelle tapaukselle kartan suurimmasta pääjännityksestä paksusuolen luumenin puoleisella pinnalla (koska jännitys on suurimmillaan tällä pinnalla) pitkin verkostosolmuja, jotka sijaitsivat pussin keskipisteen läpi kulkevilla pituussuuntaisilla ja ympäryssuuntaisilla poluilla. Suurimman pääjännityksen kehittymistä näitä polkuja pitkin seurattiin pussin keskipisteeseen nähden mitattujen etäisyyksien \({d}_{l}\) ja \({d}_{c}\) funktiona. Tämän jälkeen laskettiin pisimmät etäisyydet \({d}_{l}^{inf}\) pituussuuntaista polkua pitkin ja \({d}_{c}^{inf}\) ympäryssuuntaista polkua pitkin, joiden jälkeen suurin pääjännitys putoaa 10 %:n päähän normaalimallin vastaavasta arvosta samoja polkuja pitkin. Nämä etäisyydet muodostivat vaikutusalueen mittauksemme. Koska havaittiin, että ne saavuttavat tasaisen arvon tietyn paineen jälkeen, niiden keskiarvot kaikista tapauksista kunkin mallin osalta 5 kPa:n paineessa korreloitiin pussin geometrian eri parametrien kanssa samalla tavalla kuin \(\sigma }_{MPS,\,max}^{avg}\).
Kaikki laskelmat tehtiin Python-ohjelmointikielellä kehittämällä Jupyter-tietokirjoja42. Tiedot tallennettiin ja niitä käsiteltiin Pythonin Pandas-kirjastolla43. Matemaattiset perusoperaatiot tehtiin NumPy-kirjastolla44. Korrelaatiokerroin \({r}^{2}\) ja kaksihaarainen p-arvo laskettiin Scipy-kirjaston ”stats”-funktiolla45. Kaikki kuvaajat toteutettiin Matplotlib-kirjastolla46.