$$\rightleftharpoonup{xx}$$
$$\longleftharp{xx}$$,
$$\longrightharp{xx}$$,
Kwestie etyczne i populacja badawcza
To badanie retrospektywne zostało zatwierdzone przez komisję etyki instytucjonalnej. Świadoma zgoda została zniesiona, ponieważ wszystkie chronione informacje zdrowotne zostały zanonimizowane. Populacja badawcza składała się z pacjentów pobranych z bazy obrazowej Pierwszego Szpitala Ludowego w Nantong, którzy przeszli rezonans magnetyczny odcinka lędźwiowego kręgosłupa z powodu bólu dolnej części pleców między styczniem 2022 a grudniem 2023 roku. Zebrano również cechy kliniczne pacjentów (Tabela 1).
Kryteria włączenia i wykluczenia
Obecny standard diagnostyczny dla DLBP opiera się na metodzie dyskografii Międzynarodowego Stowarzyszenia Badań nad Bólem z 1995 roku, która wywołuje ból poprzez zwiększony nacisk, ale nie jest powszechnie akceptowana ze względu na wysoką inwazjność8. Dlatego badanie przyjęło restrykcyjny intraoperacyjny test prowokacji bólu dysku jako metodę diagnostyczną DLBP, zgodny z zasadami dyskografii, jednocześnie minimalizując dalsze uszkodzenia dysku i powikłania31. Procedura polegała na tym, że pacjenci leżeli w pozycji brzucha pod znieczuleniem, z docelowym krążkiem lędźwiowym lokalizowanym za pomocą rentgenowskiego naprowadzania ramienia C. Igła o masie 18 G (<22 G) została wprowadzona przez podejście posterolateralne do centralnego jądra pulsposus, omijając korzenie nerwowe i worek oponowiony. Po potwierdzeniu wszczepienia igły podano sól fizjologiczną lub kontrast niejonowy (np. joheksol) w tempie nieprzekraczającym 0,5 mL/min, ciśnieniem poniżej 50 psi i całkowitą objętością wstrzyku poniżej 3 mL, aby symulować wzrost ciśnienia wewnątrzdyskowego i wywołać znane objawy bólu dolnej części pleców, z wynikiem skali analogowej wzroku (VAS) ≥ 7.
Kryteria włączenia do grupy DLBP: Pacjenci uwzględnieni tutaj przeszli badanie MRI, mają nawracający ból dolnej części pleców przez ponad 3 miesiące, z niepowodzeniem leczenia zachowawczym, z drętwieniem kończyn dolnych lub bez niego lub promieniującym bólem, a także pozytywnie wychodzą na śródoperacyjny test wywołujący ból dysku.
Kryteria włączenia do grupy niebędącej DLBP: Pacjenci uwzględnieni tutaj przeszli rezonans magnetyczny; nie mają historii bólu dolnej części pleców w ciągu 3 miesięcy lub są zdrowymi osobami poddawanymi badaniu fizykalnym; nie ma nieprawidłowości w rezonanse magnetycznym; oraz posiadają standaryzowaną ocenę z Oswestry Disability Index (ODI) <10 i VAS ≤2.
Kryteria wykluczenia: Wykluczyć pacjentów z innymi przyczynami bólu dolnej części pleców, takimi jak znaczna przepuklina dysku uciskająca nerwy, złamania, infekcje kręgosłupa, spondylolistyka, guzy, osteoporoza lub choroby metaboliczne; historię operacji przed badaniem; niewyraźne lub słabej jakości obrazy; niemożność identyfikacji regionu zainteresowania (ROI).
Ostatecznie uwzględniono 243 pacjentów, w tym 81 pacjentów z DLBP i 162 grupy kontrolnej.
Parametry MRI
Wszyscy pacjenci uwzględnieni w tym badaniu przeszli badania MRI 3,0 T, wykorzystując sekwencje obejmujące sagittalne i osiowe obrazowanie ważone T1 (T1WI) oraz T2-ważone (T2WI). W badaniu zastosowano trzy różne urządzenia MRI: Siemens Verio, Siemens Prisma oraz Philips Ingenia CX. Parametry skanowania zostały ustalone następująco: dla sagittalnego T2WI TR mieściło się od 2000 do 4597 ms, TE od 90 do 120 ms, grubość przekroju od 4,0 do 4,8 mm, z 15 fragmentami, szerokość pasma od 250 Hz do 340 Hz, rozmiar matrycy 384 × 384 lub 512 × 512, procent pola widzenia fazowego 100%, oraz pole widzenia odczytu 300 mm.
Funkcje obrazowania HIZ – pomiary i statystyka
W tym badaniu wykorzystano oprogramowanie 3D Slicer (wersja 5.6.1, https://download.slicer.org/?version=5.6.1) do ręcznej analizy obrazów ważonych T2 rezonansu lędźwiowego u 243 pacjentów i grup kontrolnych. Proces przetwarzania danych wyglądał następująco: Obrazy były importowane za pomocą funkcji Dodaj danych DICOM w 3D Slicerze. Pomiary przeprowadzono przy użyciu funkcji Markups przez dwóch doktorantów badań nad kręgosłupem pod kierunkiem radiologów i chirurgów kręgosłupa z ponad 10-letnim doświadczeniem klinicznym. Rozbieżności w pomiarach zostały rozwiązane dzięki konsultacji z oboma lekarzami. Wszystkie wartości pomiarów były uśredniane przez dwóch statystyków, aby zapewnić dokładność. Na zdjęciach strzałkowych mierzono dwie wzajemnie prostopadłe maksymalne średnice strefy wysokiej intensywności (HIZ), definiowane jako Hizh (kierunek bliski pionowi) i Hizw (kierunek bliski poziomu). Funkcja segmentacji została użyta do wyznaczenia regionu HIZ, a obszar najbardziej wyrazistego HIZ na płaszczyźnie strzałkowej obliczono, oznaczany jako Hizarea. Na obrazach osiowych mierzono maksymalną długość HIZ, nazwaną Hizl. Dodatkowo zdefiniowano trzy zmienne binarne: Hiz (obecność HIZ), Other (obecność wielosegmentalnego HIZ) oraz Position (czy HIZ przecina tylną linię środkową) (Tabela 2).
Ekstrakcja cech radiomiki i standaryzacja
Aby zintegrować specyficzną operację Resampling Scalar Volume w kontekście przetwarzania obrazów ważonych T2 oddzielenia lędźwiowego od 243 pacjentów i kontrolnych za pomocą 3D Slicer, zastosowano następujące szczegółowe kroki: Proces rozpoczął się od importu obrazów ważonych T2 odcinka lędźwiowego do 3D Slicer. Aby zapewnić spójność i zmniejszyć błędy heterogeniczności, wszystkie obrazy zostały ponownie próbkowane do rozmiaru wokselowego 0,6 × 0,6 × 0,6 mm za pomocą modułu Resample Scalar Volume . Konkretne kroki to: Przejście do sekcji Moduły i wybór Resample Scalar Volume. W Ustawieniu parametrów upewnij się, że wybrano objętość skalarną ponownego próbkowania . W sekcji Parametry Resampling ustaw odstępy na 0,6, 0,6, 0,6, aby zdefiniować docelowe wymiary wokseli. Wybierz odpowiednią metodę interpolacji spośród opcji takich jak liniowa, najbliższy sąsiad, bspline, hamming, cosinus, welch, lanczos lub blackman, z liniową jako domyślną. Dla Output Volume wybierz lub stwórz nowy wolumin, aby przechować dane z resamplowaniem. Po weryfikacji ustawień kliknij Zastosowaj , aby uruchomić proces ponownego samplingu.
Po ponownym pobieraniu próbek dwóch doktorantów zajmujących się badaniami kręgosłupa, pod kierunkiem radiologów i chirurgów kręgosłupa z ponad 10-letnim doświadczeniem, przeprowadziło półautomatyczne wyznaczanie obszaru zainteresowania (ROI). Rozbieżności zostały rozwiązane poprzez konsultacje. Konkretne kroki to: Wybranie warstwy (ponowne próbkowanie nowego woluminu), utworzenie nowej warstwy segmentacji oraz użycie funkcji Draw do wyznaczania warstw po warstwach. Użyj funkcji Wypełnienia między fragmentami do wypełniania warstw. Wygładz powstały ROI metodą Mediana o rozmiarze jądra 3,0 mm, 5 x 5 pikseli. Korzystając z biblioteki PyRadiomics, z ponownego próbkowania wyodrębniono 107 cech radiomicznych, w tym cechy kształtu, cechy statystyczne pierwszego rzędu, macierz współwystępowania na poziomie szarości (GLCM), macierz długości długości serii na poziomie szarości (GLRLM), macierz strefy rozmiaru na poziomie szarości (GLSZM), sąsiednią macierz różnic tonów szarych (NGTDM) oraz cechy macierzy zależności na poziomie szarości (GLDM) (Tabela Uzupełniająca 1). Na koniec wyodrębnione dane cech zostały standaryzowane metodą Z-score, aby przekształcić ich naturalny zakres w standaryzowany zakres.
Grupowanie na podstawie cech
Grupy zostały przydzielone następująco:
d0 (grupa bazowa): Uwzględniono cechy kliniczne, n(d0) = 5.
d1 (grupa finostroju bazowego): Uwzględniono cechy kliniczne i cechy obrazowania HIZ, n(d1) = 12.
d2 (grupa modelowa): Uwzględniono cechy kliniczne i radiomiki, n(d2) = 112.
D (Grupa finetune modelu): Uwzględniono cechy kliniczne, cechy obrazowania HIZ oraz cechy radiomiki, n(D) = 119.
Odczyt i wstępne przetwarzanie danych
Dane dla poszczególnych grup były odczytywane za pomocą funkcji read_excel w oprogramowaniu R (wersja 4.3.1, https://www.r-project.org/, platforma: x86_64-w64-mingw32/x64 [64-bit]) z kodowaniem UTF-8 (domyślne w systemie), które obsługuje większość języków pisanych na świecie. Funkcja select służyła do oddzielenia zmiennej docelowej od zmiennych cech.
Wybór funkcji
Aby zapewnić powtarzalność podziału danych, cztery zbiory danych (D, d2, d1, d0) przebiegały tym samym procesem przetwarzania, z ustalonym losowym seedem 80. Do ramki danych dodawano kolumnę etykiet, przekształcano ją na typ czynnika, a losowe etykiety (1 i 0) generowano za pomocą rozkładu dwumianowego, dzieląc dane na zestawy treningowe i testowe w stosunku 8:2 (80% prawdopodobieństwa dla zbioru treningowego). Ze względu na ograniczoną liczbę cech, grupy d0 i d1 nie wymagały wyboru cech, podczas gdy grupy d2 i D przechodziły wybór cech:
Najpierw zastosowano test Mann-Whitney U do zbioru treningowego, aby wybrać cechy o wartościach p <0,05. Następnie wykonano regresję Lasso z 10-krotną walidacją krzyżową i regularizacją L1, aby określić optymalny współczynnik kary i wybrać cechy z jednej iteracji. Ustawiano losowe ziarno od 1 do 100, a powyższe kroki (podział danych, test Mann-Whitney U oraz regresja Lasso) powtarzano w pętli. Cechy pojawiające się ponad 50 razy w 100 iteracjach zostały wybrane do kolejnych modelowania, aby zapewnić wiarygodność konstrukcji modelu i diagnostyki. Cechy modelowania dla każdej grupy wymienione są w Tabeli 3.
Wyszukiwanie siatki i strojenie modeli
Aby zoptymalizować modele, zastosowano wyszukiwanie siatki do systematycznego badania i oceny różnych kombinacji hiperparametrów. Dzięki iteracyjnemu strojeniu i ocenie wydajności badanie to zidentyfikowało zestaw hiperparametrów, które zapewniły najlepsze wartości ROC AUC dla zestawów pociągowych i testowych, optymalizując tym samym wydajność modelu.
Rozwój i ocena modelu
Cały rozwój modeli i analiza danych były prowadzone w R (wersja 4.3.1). Badanie to wykorzystało różne algorytmy uczenia maszynowego, w tym Random Forest (RF), Support Vector Machine (SVM), Decision Tree (TREE), K-Nearest Neighbors (KNN) oraz Logistic Regression (LOG), aby opracować 20 modeli w czterech grupach (d0, d1, d2 i D) do przewidywania diagnozy DLBP. Wydajność modelu oceniano na podstawie następujących metryk: ROC AUC, PR AUC, dokładność, czułość, specyficzność, dodatnia wartość predykcyjna (PPV), ujemna wartość predykcyjna (NPV) oraz wynik F1.
Wartości AUC ROC i PR AUC zostały bezpośrednio wygenerowane z krzywych ROC i PR, natomiast pozostałe metryki (dokładność, czułość, specyficzność, PPV, NPV i wynik F1) zostały obliczone i statystycznie przeanalizowane przy użyciu odpowiadających funkcji w R. Spośród 20 opracowanych modeli wybrano 8 reprezentatywnych modeli do oceny na podstawie powyższych wskaźników wydajności.
Wzory dla każdej metryki są następujące:




Q

TP: Prawdziwy pozytyw; TN: Prawdziwy negatyw; FP: Fałszywy alarm; FN: Fałszywie negatywne.
Analiza interpretowalności SHAP
Model Random Forest w grupie D wykazał najlepsze wyniki w ocenie modelu. Aby lepiej zinterpretować prognozy modelu, badanie to przebudowało model w Pythonie (wersja 3.7.9) i przeprowadziło analizę interpretowalności SHAP. Aby zachować spójność, hiperparametry trenowania modelu Random Forest w Pythonie były identyczne jak te używane w R, a proces przetwarzania danych również był spójny. Konkretnie użyliśmy interpretera modelu drzewa do analizy wkładu każdej cechy w diagnostyczne prognozy modelu w zbiorach treningowych i testowych. Dla lepszej wizualizacji wygenerowaliśmy wykresy rozkładu wartości SHAP, wykresy ważności cech oraz wykresy siły SHAP dla poszczególnych prognoz.
Analiza statystyczna
Wszystkie analizy danych przeprowadzono w R (wersja 4.3.1, https://www.r-project.org/, platforma: x86_64-w64-mingw32/x64 (64-bit)) oraz Pythonie (wersja 3.7.9, https://www.python.org/downloads/release/python-379/). Zmienne ciągłe opisywano jako średnią ± odchylenie standardowe, natomiast zmienne kategoryczne jako częstotliwość i procent. Porównanie między DLBP a Non-DLBP przeanalizowano za pomocą testu U, przy czym wartość P poniżej 0,05 została uznana za istotną statystycznie.